1 ///////////////////////////////////////////////////////////////////////////
3 // AliFemtoModelCorrFctn3DLCMSSpherical: a class to calculate 3D correlation //
4 // for pairs of identical particles, binned in spherical coordinates. //
5 // In analysis the function should be first created in a macro, then //
6 // added to the analysis, and at the end of the macro the procedure to //
7 // write out histograms should be called. //
9 ///////////////////////////////////////////////////////////////////////////
11 #include "AliFemtoModelCorrFctn3DLCMSSpherical.h"
12 #include "AliFemtoModelManager.h"
17 ClassImp(AliFemtoModelCorrFctn3DLCMSSpherical)
20 //____________________________
21 AliFemtoModelCorrFctn3DLCMSSpherical::AliFemtoModelCorrFctn3DLCMSSpherical(char* title, const int& nqbins, const float& QLo, const float& QHi, const int& nphibins, const int& ncthetabins):
28 char tTitNum[101] = "NumTrue";
29 strncat(tTitNum,title, 100);
30 fTrueNumeratorSph = new TH3D(tTitNum,title,nqbins,QLo,QHi,nphibins,-TMath::Pi(),TMath::Pi(),ncthetabins,-1.0,1.0);
32 char tTitNumF[101] = "NumFake";
33 strncat(tTitNumF,title, 100);
34 fFakeNumeratorSph = new TH3D(tTitNumF,title,nqbins,QLo,QHi,nphibins,-TMath::Pi(),TMath::Pi(),ncthetabins,-1.0,1.0);
36 char tTitDen[101] = "Den";
37 strncat(tTitDen,title, 100);
38 fDenominatorSph = new TH3D(tTitDen,title,nqbins,QLo,QHi,nphibins,-TMath::Pi(),TMath::Pi(),ncthetabins,-1.0,1.0);
40 // to enable error bar calculation...
41 fTrueNumeratorSph->Sumw2();
42 fFakeNumeratorSph->Sumw2();
43 fDenominatorSph->Sumw2();
46 AliFemtoModelCorrFctn3DLCMSSpherical::AliFemtoModelCorrFctn3DLCMSSpherical(const AliFemtoModelCorrFctn3DLCMSSpherical& aCorrFctn) :
47 AliFemtoModelCorrFctn(),
54 fTrueNumeratorSph = new TH3D(*aCorrFctn.fTrueNumeratorSph);
55 fFakeNumeratorSph = new TH3D(*aCorrFctn.fFakeNumeratorSph);
56 fDenominatorSph = new TH3D(*aCorrFctn.fDenominatorSph);
57 fPairCut = aCorrFctn.fPairCut;
59 //____________________________
60 AliFemtoModelCorrFctn3DLCMSSpherical::~AliFemtoModelCorrFctn3DLCMSSpherical(){
62 delete fTrueNumeratorSph;
63 delete fFakeNumeratorSph;
64 delete fDenominatorSph;
66 //_________________________
67 AliFemtoModelCorrFctn3DLCMSSpherical& AliFemtoModelCorrFctn3DLCMSSpherical::operator=(const AliFemtoModelCorrFctn3DLCMSSpherical& aCorrFctn)
69 // assignment operator
70 if (this == &aCorrFctn)
73 if (fTrueNumeratorSph) delete fTrueNumeratorSph;
74 fTrueNumeratorSph = new TH3D(*aCorrFctn.fTrueNumeratorSph);
75 if (fFakeNumeratorSph) delete fFakeNumeratorSph;
76 fFakeNumeratorSph = new TH3D(*aCorrFctn.fFakeNumeratorSph);
77 if (fDenominatorSph) delete fDenominatorSph;
78 fDenominatorSph = new TH3D(*aCorrFctn.fDenominatorSph);
80 fPairCut = aCorrFctn.fPairCut;
85 //_________________________
86 void AliFemtoModelCorrFctn3DLCMSSpherical::WriteOutHistos(){
87 // Write out all histograms to file
88 fTrueNumeratorSph->Write();
89 fFakeNumeratorSph->Write();
90 fDenominatorSph->Write();
92 //______________________________
93 TList* AliFemtoModelCorrFctn3DLCMSSpherical::GetOutputList()
95 // Prepare the list of objects to be written to the output
96 TList *tOutputList = new TList();
98 tOutputList->Add(fTrueNumeratorSph);
99 tOutputList->Add(fFakeNumeratorSph);
100 tOutputList->Add(fDenominatorSph);
105 //_________________________
106 void AliFemtoModelCorrFctn3DLCMSSpherical::Finish(){
107 // here is where we should normalize, fit, etc...
110 //____________________________
111 AliFemtoString AliFemtoModelCorrFctn3DLCMSSpherical::Report(){
112 // Construct the report
113 string stemp = "LCMS Frame Spherical 3D Model Correlation Function Report:\n";
115 snprintf(ctemp , 100, "Number of entries in numerator:\t%E\n",fTrueNumeratorSph->GetEntries());
117 snprintf(ctemp , 100, "Number of entries in denominator:\t%E\n",fDenominatorSph->GetEntries());
121 snprintf(ctemp , 100, "Here is the PairCut specific to this CorrFctn\n");
123 stemp += fPairCut->Report();
126 snprintf(ctemp , 100, "No PairCut specific to this CorrFctn\n");
131 AliFemtoString returnThis = stemp;
134 //____________________________
135 void AliFemtoModelCorrFctn3DLCMSSpherical::AddRealPair( AliFemtoPair* pair){
136 // perform operations on real pairs
138 if (!(fPairCut->Pass(pair))) return;
141 Double_t weight = fManager->GetWeight(pair);
144 double tKO = (pair->QOutCMS());
145 double tKS = (pair->QSideCMS());
146 double tKL = (pair->QLongCMS());
148 double tKR = sqrt(tKO*tKO + tKS*tKS + tKL*tKL);
150 if ( fabs(tKR) < 1e-10 ) tKC = 0.0;
152 double tKP=atan2(tKS,tKO);
154 fTrueNumeratorSph->Fill(tKR,tKP,tKC,weight);
156 //____________________________
157 void AliFemtoModelCorrFctn3DLCMSSpherical::AddMixedPair( AliFemtoPair* pair){
158 // perform operations on mixed pairs
160 if (!(fPairCut->Pass(pair))) return;
163 Double_t weight = fManager->GetWeight(pair);
165 double tKO = (pair->QOutCMS());
166 double tKS = (pair->QSideCMS());
167 double tKL = (pair->QLongCMS());
169 double tKR = sqrt(tKO*tKO + tKS*tKS + tKL*tKL);
171 if ( fabs(tKR) < 1e-10 ) tKC = 0.0;
173 double tKP=atan2(tKS,tKO);
175 fFakeNumeratorSph->Fill(tKR,tKP,tKC,weight);
176 fDenominatorSph->Fill(tKR,tKP,tKC);