5 #include "AliRunAnalysis.h"
6 #include "AliHBTAnalysis.h"
8 #include "AliReaderESD.h"
9 #include "AliReaderKineTree.h"
10 #include "AliParticleCut.h"
12 #include "AliAODPairCut.h"
13 #include "AliAODParticleCut.h"
14 #include "AliHBTQResolutionFctns.h"
15 #include "AliHBTQDistributionFctns.h"
16 #include "AliHBTMonDistributionFctns.h"
17 #include "AliHBTTwoTrackEffFctn.h"
18 #include "AliHBTCorrelFctn.h"
19 #include "AliHBTWeightFctn.h"
20 #include "AliHBTWeightTheorFctn.h"
21 #include "AliHBTWeights.h"
22 #include "AliHBTLLWeights.h"
24 #include "TObjString.h"
27 #include "AliHBTMonDistributionFctns.h"
28 #include "AliHBTMonResolutionFctns.h"
32 //This macro can not be executed from
33 AliHBTAnalysis* buildanalysis(Int_t pid1, Int_t pid2, Bool_t threeDcuts, Float_t qinvmax, const char* anasuffix = "");
36 void hbtanalysis(Option_t* datatype="Kine", Option_t* processopt="TracksAndParticles",
37 Int_t first = -1,Int_t last = -1,
40 const char* basedir=".";
41 const char* serie=".";
42 const char* field = ".";
49 Bool_t ident = kTRUE;//identical particles
50 Bool_t nonident = kFALSE;//non identical particles analysis
51 // Dynamically link some shared libs
52 // gROOT->LoadMacro("loadlibs.C");
54 gSystem->Load("$(ALICE_ROOT)/lib/tgt_$(ALICE_TARGET)/libANALYSIS");
55 gSystem->Load("$(ALICE_ROOT)/lib/tgt_$(ALICE_TARGET)/libHBTAN");
57 /***********************************************************/
58 //Create analysis and reader
62 Int_t kine = strcmp(datatype,"Kine");
63 Int_t ESD = strcmp(datatype,"ESD");
64 Int_t intern = strcmp(datatype,"Intern");
68 reader = new AliReaderKineTree();
69 processopt="Particles"; //this reader by definition reads only simulated particles
74 reader = new AliReaderESD();
80 // reader = new AliReaderAOD(dataname+".data.root");
84 cerr<<"Option "<<datatype<<" not recognized. Exiting"<<endl;
88 cout<<"hbtanalysis.C: Reader Set, "<<reader->Class()->GetName()<<" dataname: "<<dataname<<endl;
89 cout<<"hbtanalysis.C: Setting dirs\n";
91 if ( (first >= 0) && (last>=0) && ( (last-first)>=0 ) )
92 {//read from many dirs dirs
94 dirs = new TObjArray(last-first+1);
95 for (Int_t i = first; i<=last; i++)
97 // sprintf(buff,"%s/%s/%s/%05.5d",basedir,field,serie,i);
98 sprintf(buff,"%s/%s/%s/%d",basedir,field,serie,i);
99 TObjString *odir= new TObjString(buff);
103 reader->SetDirs(dirs);
104 reader->SetName("READER");
105 reader->SetTitle("Reader Used in HBT Analysis");
106 cout<<"hbtanalysis.C: Dirs set\n";
107 /***********************************************************/
108 /***** R E A D E R ***************************************/
109 /***********************************************************/
111 //we want have only low pt pi+ so set a cut to reader
112 AliAODParticleCut* readerpartcut= new AliAODParticleCut();
113 readerpartcut->SetPtRange(0.0,10000.0);
114 // readerpartcut->SetThetaRange(0.25*TMath::Pi(), 0.75*TMath::Pi());//45->135 deg.
115 readerpartcut->SetYRange(-1.0,1.0);
120 readerpartcut->SetPID(PID[l]);
121 reader->AddParticleCut(readerpartcut);//read this particle type with this cut
125 TString filename = prefname + dataname + ".anal.";
132 /************************************************************/
133 /**** P R O C E S S ***********************/
134 /************************************************************/
135 AliHBTAnalysis* hbtanalysis = 0x0;
141 pdg1 = TDatabasePDG::Instance()->GetParticle(PID[l]);
144 cerr<<"ERROR: hbtanalysis.C: Can not find particle type "<<PID[l]<<" in PDG DataBase\n";
148 p1name= pdg1->GetName();
152 if ( (k==l)&&(ident==kFALSE) ) { k++;continue;}
153 if ( (k!=l)&&(nonident==kFALSE) ) { k++;continue; }
155 pdg2 = TDatabasePDG::Instance()->GetParticle(PID[k]);
158 cerr<<"ERROR: hbtanalysis.C: Can not find particle type "<<PID[k]<<" in PDG DataBase\n";
162 p2name = pdg2->GetName();
163 cout<<"hbtanalysis.C: buildanalysis for "<<PID[l]<<"("<<p1name<<") "
164 <<PID[k]<<"("<<p2name<<")\n";
165 fn = filename + p1name + p2name + ".root";
166 cout<<"hbtanalysis.C: File Name is: "<<fn<<endl;
168 AliRunAnalysis* analysis = new AliRunAnalysis();
169 analysis->SetReader(reader);
171 cout<<"hbtanalysis.C: buildanalysis "<<endl;
172 hbtanalysis = buildanalysis(PID[l],PID[k],kFALSE,1.0,"-a1");
173 cout<<"hbtanalysis.C: buildanalysis Done"<<endl;
174 if (hbtanalysis == 0x0)
176 cout<<"ERROR: hbtanalysis.C: buildanalysis function returned NULL pointer\n";
177 cout<<"ERROR: Deleting Run Analysis ..."<<endl;
179 cout<<"ERROR: Deleting Run Analysis Done"<<endl;
183 TFile* histoOutput = TFile::Open(fn,"recreate");
184 // hbtanalysis->SetOutputFileName(fn);
186 AliHBTLLWeights::Instance()->SetParticlesTypes(PID[l],PID[k]);
187 AliHBTLLWeights::Instance()->SetColoumb(kFALSE);
188 AliHBTLLWeights::Instance()->SetStrongInterSwitch(kFALSE);
189 AliHBTLLWeights::Instance()->Init();
190 AliHBTLLWeights::Instance()->Set();
192 // AliHBTCrab::Instance()->Init(PID[l],PID[k]);
193 // AliHBTCrab::Instance()->Set();
195 hbtanalysis->SetProcessOption(AliHBTAnalysis::kSimulated);
196 analysis->Add(hbtanalysis);
199 cout<<"hbtanalysis.C: Run\n";
202 cout<<"hbtanalysis.C: Write\n";
204 cout<<"hbtanalysis.C: Delete analysis\n";
205 hbtanalysis->SetOwner(kTRUE);
207 histoOutput->Close();
213 delete readerpartcut;
222 AliHBTAnalysis* buildanalysis(Int_t pid1, Int_t pid2, Bool_t threeDcuts, Float_t qinvmax, const char* anasuffix)
224 /************************************************************/
225 /**** Q INV Correlation Function ************************/
226 /************************************************************/
229 AliHBTAnalysis * analysis = new AliHBTAnalysis();
230 analysis->SetBufferSize(1);
233 AliAODPairCut *paircut = new AliAODPairCut();
235 TString anasfx(anasuffix);
237 Float_t qinvmin = 0.0;
238 paircut->SetQInvRange(qinvmin,1.5*qinvmax);
240 TParticlePDG* pdg1 = TDatabasePDG::Instance()->GetParticle(pid1);
241 TString system = pdg1->GetName();
245 cout<<"Setting cuts on particles: "<<pid1<<" "<<pid2<<endl;
246 AliAODParticleCut* partcut = new AliAODParticleCut();
247 partcut->SetPID(pid1);
248 paircut->SetFirstPartCut(partcut);
249 partcut->SetPID(pid2);
250 paircut->SetSecondPartCut(partcut);
252 TParticlePDG* pdg2 = TDatabasePDG::Instance()->GetParticle(pid2);
253 system+=pdg2->GetName();
257 system+=pdg1->GetName();
259 analysis->SetGlobalPairCut(paircut);
260 paircut->SetName("AnalysisGlobalPairCut");
263 cout<<"Building functions"<<endl;
265 /************************************************************/
266 /**** Q INV Correlation Function ************************/
267 /************************************************************/
268 //general cross check
269 AliHBTQInvCorrelFctn * qinvcfP = new AliHBTQInvCorrelFctn(1000,qinvmax);
270 qinvcfP->SetNumberOfBinsToScale(200);
271 qinvcfP->Rename("qinvcfP"+anasfx,"Particle (simulated) Q_{inv} Correlation Function "+ system + anasfx);
272 analysis->AddParticleFunction(qinvcfP);
274 //===========================================================================================
276 AliHBTWeightTheorQInvFctn *wqinvcf = new AliHBTWeightTheorQInvFctn(1000,qinvmax);
277 wqinvcf->SetNumberOfBinsToScale(200);
278 wqinvcf->Rename("wqinvcf","Particle (simulated) Lednicky Q_{inv} Correlation Function "+system);
279 analysis->AddParticleFunction(wqinvcf);
281 AliHBTWeightTheorOSLFctn *wqoslcf = new AliHBTWeightTheorOSLFctn(60,qinvmax,0.0,60,qinvmax,0.0,60,qinvmax,0.0);
282 wqoslcf->Rename("wqoslcf","Particle (simulated) Lednicky OSL Theoret Correlation Function "+system);
283 wqoslcf->SetNumberOfBinsToScale(20,20,20);
284 analysis->AddParticleFunction(wqoslcf);
286 //===========================================================================================
287 AliAODPairCut *cfPairCutKt = new AliAODPairCut();
289 //===========================================================================================
290 cfPairCutKt->SetKtRange(0.0,0.2);
292 AliHBTWeightTheorQInvFctn *wqinvcfK11 = new AliHBTWeightTheorQInvFctn(1000,qinvmax);
293 wqinvcfK11->SetNumberOfBinsToScale(200);
294 wqinvcfK11->SetPairCut(cfPairCutKt);
295 wqinvcfK11->Rename("wqinvcfK11","Particle (simulated) Lednicky Q_{inv} Correlation Function 0.0<p_{t}<0.2 "+system);
296 analysis->AddParticleFunction(wqinvcfK11);
298 AliHBTWeightTheorOSLFctn *wqoslcfK11 = new AliHBTWeightTheorOSLFctn(60,qinvmax,0.0,60,qinvmax,0.0,60,qinvmax,0.0);
299 wqoslcfK11->Rename("wqoslcfK11","Particle (simulated) Lednicky OSL Theoret Correlation Function 0.0<p_{t}<0.2 "+system);
300 wqoslcfK11->SetNumberOfBinsToScale(20,20,20);
301 wqoslcfK11->SetPairCut(cfPairCutKt);
302 analysis->AddParticleFunction(wqoslcfK11);
303 //===========================================================================================
305 cfPairCutKt->SetKtRange(0.2,0.4);
307 AliHBTWeightTheorQInvFctn *wqinvcfK12 = new AliHBTWeightTheorQInvFctn(1000,qinvmax);
308 wqinvcfK12->SetNumberOfBinsToScale(200);
309 wqinvcfK12->SetPairCut(cfPairCutKt);
310 wqinvcfK12->Rename("wqinvcfK12","Particle (simulated) Lednicky Q_{inv} Correlation Function 0.2<p_{t}<0.4 "+system);
311 analysis->AddParticleFunction(wqinvcfK12);
313 AliHBTWeightTheorOSLFctn *wqoslcfK12 = new AliHBTWeightTheorOSLFctn(60,qinvmax,0.0,60,qinvmax,0.0,60,qinvmax,0.0);
314 wqoslcfK12->Rename("wqoslcfK12","Particle (simulated) Lednicky OSL Theoret Correlation Function 0.2<p_{t}<0.4 "+system);
315 wqoslcfK12->SetNumberOfBinsToScale(20,20,20);
316 wqoslcfK12->SetPairCut(cfPairCutKt);
317 analysis->AddParticleFunction(wqoslcfK12);
318 //===========================================================================================
320 cfPairCutKt->SetKtRange(0.4,0.6);
322 AliHBTWeightTheorQInvFctn *wqinvcfK21 = new AliHBTWeightTheorQInvFctn(1000,qinvmax);
323 wqinvcfK21->SetNumberOfBinsToScale(200);
324 wqinvcfK21->SetPairCut(cfPairCutKt);
325 wqinvcfK21->Rename("wqinvcfK21","Particle (simulated) Lednicky Q_{inv} Correlation Function 0.4<p_{t}<0.6 "+system);
326 analysis->AddParticleFunction(wqinvcfK21);
328 AliHBTWeightTheorOSLFctn *wqoslcfK21 = new AliHBTWeightTheorOSLFctn(60,qinvmax,0.0, 60,qinvmax,0.0, 60,qinvmax,0.0);
329 wqoslcfK21->Rename("wqoslcfK21","Particle (simulated) Lednicky OSL Theoret Correlation Function 0.4<p_{t}<0.6 "+system);
330 wqoslcfK21->SetNumberOfBinsToScale(20,20,20);
331 wqoslcfK21->SetPairCut(cfPairCutKt);
332 analysis->AddParticleFunction(wqoslcfK21);
333 //===========================================================================================
335 cfPairCutKt->SetKtRange(0.6,0.8);
337 AliHBTWeightTheorQInvFctn *wqinvcfK22 = new AliHBTWeightTheorQInvFctn(1000,qinvmax);
338 wqinvcfK22->SetNumberOfBinsToScale(200);
339 wqinvcfK22->SetPairCut(cfPairCutKt);
340 wqinvcfK22->Rename("wqinvcfK22","Particle (simulated) Lednicky Q_{inv} Correlation Function 0.6<p_{t}<0.8 "+system);
341 analysis->AddParticleFunction(wqinvcfK22);
343 AliHBTWeightTheorOSLFctn *wqoslcfK22 = new AliHBTWeightTheorOSLFctn(60,qinvmax,0.0, 60,qinvmax,0.0, 60,qinvmax,0.0);
344 wqoslcfK22->Rename("wqoslcfK22","Particle (simulated) Lednicky OSL Theoret Correlation Function 0.6<p_{t}<0.8 "+system);
345 wqoslcfK22->SetNumberOfBinsToScale(20,20,20);
346 wqoslcfK22->SetPairCut(cfPairCutKt);
347 analysis->AddParticleFunction(wqoslcfK22);
348 //===========================================================================================
350 cfPairCutKt->SetKtRange(0.8,1.2);
352 AliHBTWeightTheorQInvFctn *wqinvcfK3 = new AliHBTWeightTheorQInvFctn(1000,qinvmax);
353 wqinvcfK3->SetNumberOfBinsToScale(200);
354 wqinvcfK3->SetPairCut(cfPairCutKt);
355 wqinvcfK3->Rename("wqinvcfK3","Particle (simulated) Lednicky Q_{inv} Correlation Function 0.8<p_{t}<1.2 "+system);
356 analysis->AddParticleFunction(wqinvcfK3);
358 AliHBTWeightTheorOSLFctn *wqoslcfK3 = new AliHBTWeightTheorOSLFctn(60,qinvmax,0.0, 60,qinvmax,0.0, 60,qinvmax,0.0);
359 wqoslcfK3->Rename("wqoslcfK3","Particle (simulated) Lednicky OSL Theoret Correlation Function 0.8<p_{t}<1.2 "+system);
360 wqoslcfK3->SetNumberOfBinsToScale(20,20,20);
361 wqoslcfK3->SetPairCut(cfPairCutKt);
362 analysis->AddParticleFunction(wqoslcfK3);
363 //===========================================================================================
365 cfPairCutKt->SetKtRange(1.2,100000.);
367 AliHBTWeightTheorQInvFctn *wqinvcfK4 = new AliHBTWeightTheorQInvFctn(1000,qinvmax);
368 wqinvcfK4->SetNumberOfBinsToScale(200);
369 wqinvcfK4->SetPairCut(cfPairCutKt);
370 wqinvcfK4->Rename("wqinvcfK4","Particle (simulated) Lednicky Q_{inv} Correlation Function 1.2<p_{t} "+system);
371 analysis->AddParticleFunction(wqinvcfK4);
373 AliHBTWeightTheorOSLFctn *wqoslcfK4 = new AliHBTWeightTheorOSLFctn(60,qinvmax,0.0, 60,qinvmax,0.0, 60,qinvmax,0.0);
374 wqoslcfK4->Rename("wqoslcfK4","Particle (simulated) Lednicky OSL Theoret Correlation Function 1.2<p_{t} "+system);
375 wqoslcfK4->SetNumberOfBinsToScale(20,20,20);
376 wqoslcfK4->SetPairCut(cfPairCutKt);
377 analysis->AddParticleFunction(wqoslcfK4);
378 //===========================================================================================
381 AliHBTMonVyDistributionVsVxFctn* vydistrvsvx = new AliHBTMonVyDistributionVsVxFctn(800,1.0e-11,-1.0e-11,800,1.0e-11,-1.0e-11);
382 AliHBTMonRtDistributionVsVzFctn* rtdistrvsvz = new AliHBTMonRtDistributionVsVzFctn(800,1.0e-11,-1.0e-11,800,1.0e-11,0.0);
383 analysis->AddParticleMonitorFunction(vydistrvsvx);
384 analysis->AddParticleMonitorFunction(rtdistrvsvz);
386 AliHBTMonPhiDistributionVsPtFctn* phidistr = new AliHBTMonPhiDistributionVsPtFctn(800,0,4,200,6.3,-3.2);
387 AliHBTMonThetaDistributionVsPtFctn* thetadistr = new AliHBTMonThetaDistributionVsPtFctn();
388 analysis->AddParticleMonitorFunction(phidistr);
389 analysis->AddParticleMonitorFunction(thetadistr);
391 cout<<"Building functions Done "<<endl;
396 //===========================================================================================
397 // EEEE N N DDD =====================================================================
398 // E NN N D D =====================================================================
399 // EEE N N N D D =====================================================================
400 // E N NN D D =====================================================================
401 // EEEE N N DDD =====================================================================
402 //===========================================================================================
403 /************************************************************/
404 /**** Q OUT Correlation Function ************************/
405 /************************************************************/
406 AliAODPairCut *outPairCut = new AliAODPairCut();
407 outPairCut->SetQSideCMSLRange(-0.01,0.01);
408 outPairCut->SetQLongCMSLRange(-0.01,0.01);
409 outPairCut->SetKtRange(0.0,0.8);
411 AliHBTWeightTheorQOutFctn *wqcfoutP = new AliHBTWeightTheorQOutFctn(1000,qinvmax,-qinvmax);
412 wqcfoutP->SetNumberOfBinsToScale(200);
413 wqcfoutP->SetPairCut(outPairCut);
414 analysis->AddParticleFunction(wqcfoutP);
416 /************************************************************/
417 /**** Q SIDE Correlation Function ************************/
418 /************************************************************/
419 AliAODPairCut *sidePairCut = new AliAODPairCut();
420 sidePairCut->SetQOutCMSLRange(-0.01,0.01);
421 sidePairCut->SetQLongCMSLRange(-0.01,0.01);
422 sidePairCut->SetKtRange(0.0,0.8);
424 AliHBTWeightTheorQSideFctn *wqcfsideP = new AliHBTWeightTheorQSideFctn(1000,qinvmax,-qinvmax);
425 wqcfsideP->SetNumberOfBinsToScale(200);
426 wqcfsideP->SetPairCut(sidePairCut);
427 analysis->AddParticleFunction(wqcfsideP);
429 /************************************************************/
430 /**** Q LONG Correlation Function ************************/
431 /************************************************************/
432 AliAODPairCut *longPairCut = new AliAODPairCut();
433 longPairCut->SetQOutCMSLRange(-0.01,0.01);
434 longPairCut->SetQSideCMSLRange(-0.01,0.01);
435 longPairCut->SetKtRange(0.0,0.8);
437 AliHBTWeightTheorQLongFctn *wqcflongP = new AliHBTWeightTheorQLongFctn(1000,qinvmax,-qinvmax);
438 wqcflongP->SetNumberOfBinsToScale(200);
439 wqcfsideP->SetPairCut(longPairCut);
440 analysis->AddParticleFunction(wqcflongP);
444 AliHBTMonPxDistributionFctn* pxdistr = new AliHBTMonPxDistributionFctn();
445 AliHBTMonPyDistributionFctn* pydistr = new AliHBTMonPyDistributionFctn();
446 AliHBTMonPzDistributionFctn* pzdistr = new AliHBTMonPzDistributionFctn();
447 AliHBTMonPDistributionFctn* pdistr = new AliHBTMonPDistributionFctn();
448 analysis->AddParticleMonitorFunction(pxdistr);
449 analysis->AddParticleMonitorFunction(pydistr);
450 analysis->AddParticleMonitorFunction(pzdistr);
451 analysis->AddParticleMonitorFunction(pdistr);
453 AliHBTMonVxDistributionFctn* vxdistr = new AliHBTMonVxDistributionFctn(800,1.0e-11,-1.0e-11);
454 AliHBTMonVyDistributionFctn* vydistr = new AliHBTMonVyDistributionFctn(800,1.0e-11,-1.0e-11);
455 AliHBTMonVzDistributionFctn* vzdistr = new AliHBTMonVzDistributionFctn(800,1.0e-11,-1.0e-11);
456 AliHBTMonRDistributionFctn* vrdistr = new AliHBTMonRDistributionFctn (800,1.0e-11,-1.0e-11);
457 analysis->AddParticleMonitorFunction(vxdistr);
458 analysis->AddParticleMonitorFunction(vydistr);
459 analysis->AddParticleMonitorFunction(vzdistr);
460 analysis->AddParticleMonitorFunction(vrdistr);
463 AliAODParticleCut* vxycut = new AliAODParticleCut();
464 AliHBTMonVyDistributionVsVxFctn* vydistrvsvx1 = new AliHBTMonVyDistributionVsVxFctn(800,1.0e-11,-1.0e-11,800,1.0e-11,-1.0e-11);
465 vydistrvsvx1->Rename("vydistrvsvx1","Vx:Vy 0.2 < Pt < 1.8");
466 vxycut->SetPtRange(0.2,0.8);
467 vydistrvsvx1->SetParticleCut(vxycut);
468 analysis->AddParticleMonitorFunction(vydistrvsvx1);
470 AliHBTMonVyDistributionVsVxFctn* vydistrvsvx2 = new AliHBTMonVyDistributionVsVxFctn(800,1.0e-11,-1.0e-11,800,1.0e-11,-1.0e-11);
471 vydistrvsvx2->Rename("vydistrvsvx2","Vx:Vy 0.8 < Pt < 1.2");
472 vxycut->SetPtRange(0.8,1.2);
473 vydistrvsvx2->SetParticleCut(vxycut);
474 analysis->AddParticleMonitorFunction(vydistrvsvx2);
476 AliHBTMonVyDistributionVsVxFctn* vydistrvsvx3 = new AliHBTMonVyDistributionVsVxFctn(800,1.0e-11,-1.0e-11,800,1.0e-11,-1.0e-11);
477 vydistrvsvx3->Rename("vydistrvsvx3","Vx:Vy 1.2 < Pt ");
478 vxycut->SetPtRange(1.2,10000.0);
479 vydistrvsvx3->SetParticleCut(vxycut);
480 analysis->AddParticleMonitorFunction(vydistrvsvx3);
482 // AliHBTTwoKStarCorrelFctn* kstarP = new AliHBTTwoKStarCorrelFctn(600,qinvmax);
483 // kstarP->Rename("kstarP","Particle (simulated) 2K^{*} Correlation Function "+system);
484 // analysis->AddParticleFunction(kstarP);
488 // delete outPairCut;
489 // delete sidePairCut;
490 // delete longPairCut;