StRoot: StHbtMaker/CorrFctn/BPLCMSFrame3DCorrFctnKt.cxx Source File

00001 /***************************************************************************
00002  *
00003  * $Id: BPLCMSFrame3DCorrFctnKt.cxx,v 1.3 2003/09/02 17:58:20 perev Exp $
00004  *
00005  * Author: Mercedes Lopez Noriega, OSU, mercedes@pacific.mps.ohio-state.edu
00006  ***************************************************************************
00007  *
00008  * Description: part of STAR HBT Framework: StHbtMaker package
00009  *   This one does 3D Bertsch-Pratt decomposition in the LCMS frame,
00010  *   for kt between ktLo and ktHi max for a given number of bins (nCFs)
00011  *   kt is the transverse four-momentum of the pair
00012  *
00013  ***************************************************************************
00014  *
00015  * $Log: BPLCMSFrame3DCorrFctnKt.cxx,v $
00016  * Revision 1.3  2003/09/02 17:58:20  perev
00017  * gcc 3.2 updates + WarnOff
00018  *
00019  * Revision 1.2  2003/01/31 19:21:09  magestro
00020  * Cleared up simple compiler warnings on i386_linux24
00021  *
00022  * Revision 1.1  2002/05/17 14:25:31  mercedes
00023  * N 3D Bertsch-Pratt decomposition (kt bins) between ktmin and ktmax (LCMSFrame)
00024  *
00025  *
00026  **************************************************************************/
00027 
00028 #include "StHbtMaker/CorrFctn/BPLCMSFrame3DCorrFctnKt.h"
00029 #include <cstdio>
00030 
00031 #ifdef __ROOT__ 
00032 ClassImp(BPLCMSFrame3DCorrFctnKt)
00033 #endif
00034 
00035 //____________________________
00036 BPLCMSFrame3DCorrFctnKt::BPLCMSFrame3DCorrFctnKt(char* title, const int& nbins, const float& QLo, const float& QHi,
00037                                                  const int& nCFs, const float& KtLo, const float& KtHi){
00038 
00039   // set some stuff...
00040 
00041   mNumberCFs = nCFs;
00042   mKtMin = KtLo;
00043   mKtMax = KtHi;
00044   mDeltaKt = (mKtMax-mKtMin)/mNumberCFs;
00045 
00046   mQinvNormLo = 0.15;
00047   mQinvNormHi = 0.18;
00048   mNumRealsNorm = 0;
00049   mNumMixedNorm = 0;
00050   mCorrection = 0;  // pointer to Coulomb Correction object
00051 
00052   // set up numerator
00053   mNumerator = new StHbt3DHisto[mNumberCFs];
00054 
00055   // set up denominator
00056   mDenominator = new StHbt3DHisto[mNumberCFs];
00057 
00058   // set up ratio
00059   mRatio = new StHbt3DHisto[mNumberCFs];
00060 
00061   // set up ave qInv
00062   mQinvHisto = new StHbt3DHisto[mNumberCFs];
00063 
00064   char TitNum[100] = "Num";
00065   strcat(TitNum,title);
00066   char TitDen[100] = "Den";
00067   strcat(TitDen,title);
00068   char TitRat[100] = "Rat";
00069   strcat(TitRat,title);
00070   char TitQinv[100] = "Qinv";
00071   strcat(TitQinv,title);
00072   
00073   for (int i=0; i<mNumberCFs; i++){
00074     
00075     sprintf(TitNum,"NumBPLCMSFrameCFKtBin%d",i);
00076     mNumerator[i].SetName(TitNum);
00077     mNumerator[i].SetTitle(title);
00078     mNumerator[i].SetBins(nbins,QLo,QHi,nbins,QLo,QHi,nbins,QLo,QHi);
00079     mNumerator[i].SetDirectory(0);
00080     mNumerator[i].Sumw2();
00081     
00082     sprintf(TitDen,"DenBPLCMSFrameCFKtBin%d",i);
00083     mDenominator[i].SetName(TitDen);
00084     mDenominator[i].SetTitle(title);
00085     mDenominator[i].SetBins(nbins,QLo,QHi,nbins,QLo,QHi,nbins,QLo,QHi);
00086     mDenominator[i].SetDirectory(0);
00087     mDenominator[i].Sumw2();
00088     
00089     sprintf(TitRat,"RatBPLCMSFrameCFKtBin%d",i);
00090     mRatio[i].SetName(TitRat);
00091     mRatio[i].SetTitle(title);
00092     mRatio[i].SetBins(nbins,QLo,QHi,nbins,QLo,QHi,nbins,QLo,QHi);
00093     mRatio[i].SetDirectory(0);
00094     mRatio[i].Sumw2();
00095 
00096     sprintf(TitQinv,"QInvHistoKtBin%d",i);
00097     mQinvHisto[i].SetName(TitRat);
00098     mQinvHisto[i].SetTitle(title);
00099     mQinvHisto[i].SetBins(nbins,QLo,QHi,nbins,QLo,QHi,nbins,QLo,QHi);
00100     mQinvHisto[i].SetDirectory(0);
00101     mQinvHisto[i].Sumw2();
00102   }
00103 }
00104 //____________________________
00105 BPLCMSFrame3DCorrFctnKt::~BPLCMSFrame3DCorrFctnKt(){
00106   delete [] mNumerator;
00107   delete [] mDenominator;
00108   delete [] mRatio;
00109   delete [] mQinvHisto;
00110 }
00111 //_________________________
00112 void BPLCMSFrame3DCorrFctnKt::Finish(){
00113   // here is where we should normalize, fit, etc...
00114 
00115   for (int i=0; i<mNumberCFs; i++){
00116     double NumFact,DenFact;
00117     if ((mNumRealsNorm !=0) && (mNumMixedNorm !=0)){
00118       NumFact = double(mNumRealsNorm);
00119       DenFact = double(mNumMixedNorm);
00120     }
00121     // can happen that the mNumRealsNorm and mNumMixedNorm = 0 if you do non-standard
00122     //   things like making a new CorrFctn and just setting the Numerator and Denominator
00123     //   from OTHER CorrFctns which you read in (like when doing parallel processing) 
00124     else{
00125       cout << "Warning! - no normalization constants defined - I do the best I can..." << endl;
00126       int nbins = mNumerator[i].GetNbinsX();
00127       int half_way = nbins/2;
00128       NumFact = mNumerator[i].Integral(half_way,nbins,half_way,nbins,half_way,nbins);
00129       DenFact = mDenominator[i].Integral(half_way,nbins,half_way,nbins,half_way,nbins);
00130     }
00131     
00132     mRatio[i].Divide(&mNumerator[i],&mDenominator[i],DenFact,NumFact);
00133     mQinvHisto[i].Divide(&mDenominator[i]);
00134   }
00135 }
00136 
00137 //____________________________
00138 StHbtString BPLCMSFrame3DCorrFctnKt::Report(){
00139   int mNumeratorEntries=0;
00140   int mDenominatorEntries=0;
00141   int mRatioEntries=0;
00142   int mQinvHistoEntries=0;
00143 
00144   for (int i=0; i<mNumberCFs; i++){
00145     mNumeratorEntries += (int)mNumerator[i].GetEntries();
00146     mDenominatorEntries += (int)mDenominator[i].GetEntries();
00147     mRatioEntries += (int)mRatio[i].GetEntries();
00148     mQinvHistoEntries += (int)mQinvHisto[i].GetEntries();
00149   }
00150 
00151   string stemp = "LCMS Frame Bertsch-Pratt 3D Correlation Function Report:\n";
00152   char ctemp[100];
00153   sprintf(ctemp,"Number of entries in numerator:\t%i\n",mNumeratorEntries);
00154   stemp += ctemp;
00155   sprintf(ctemp,"Number of entries in denominator:\t%i\n",mDenominatorEntries);
00156   stemp += ctemp;
00157   sprintf(ctemp,"Number of entries in ratio:\t%i\n",mRatioEntries);
00158   stemp += ctemp;
00159   if (mCorrection)
00160     {
00161       float radius = mCorrection->GetRadius();
00162       sprintf(ctemp,"Coulomb correction used radius of\t%E\n",radius);
00163     }
00164   else
00165     {
00166       sprintf(ctemp,"No Coulomb Correction applied to this CorrFctn\n");
00167     }
00168   stemp += ctemp;
00169   
00170   StHbtString returnThis = stemp;
00171   return returnThis;
00172 }
00173 //____________________________
00174 void BPLCMSFrame3DCorrFctnKt::AddRealPair(const StHbtPair* pair){
00175   int mIndex = (int)(((fabs(pair->kT())-mKtMin)/mDeltaKt));
00176   if ((mIndex>=0)&&(mIndex<mNumberCFs)){
00177     double Qinv = fabs(pair->qInv());   // note - qInv() will be negative for identical pairs...
00178     if ((Qinv < mQinvNormHi) && (Qinv > mQinvNormLo)) mNumRealsNorm++;
00179     double qOut = fabs(pair->qOutCMS());
00180     double qSide = fabs(pair->qSideCMS());
00181     double qLong = fabs(pair->qLongCMS());
00182 
00183     mNumerator[mIndex].Fill(qOut,qSide,qLong);
00184   }
00185 }
00186 
00187 //____________________________
00188 void BPLCMSFrame3DCorrFctnKt::AddMixedPair(const StHbtPair* pair){
00189   int  mIndex = (int)(((fabs(pair->kT())-mKtMin)/mDeltaKt));
00190   if ((mIndex>=0)&&(mIndex<mNumberCFs)){
00191     double CoulombWeight = (mCorrection ? mCorrection->CoulombCorrect(pair) : 1.0);
00192 
00193     double Qinv = fabs(pair->qInv());   // note - qInv() will be negative for identical pairs...
00194     if ((Qinv < mQinvNormHi) && (Qinv > mQinvNormLo)) mNumMixedNorm++;
00195     double qOut = fabs(pair->qOutCMS());
00196     double qSide = fabs(pair->qSideCMS());
00197     double qLong = fabs(pair->qLongCMS());
00198 
00199     mDenominator[mIndex].Fill(qOut,qSide,qLong,CoulombWeight);
00200     mQinvHisto[mIndex].Fill(qOut,qSide,qLong,Qinv);
00201   }
00202 }
00203