1 #include "AliHBTAnalysis.h"
2 //_________________________________________________________
3 ////////////////////////////////////////////////////////////////////////////
5 // class AliHBTAnalysis
7 // Central Object Of HBTAnalyser:
8 // This class performs main looping within HBT Analysis
9 // User must plug a reader of Type AliReader
10 // User plugs in coorelation and monitor functions
11 // as well as monitor functions
13 // HBT Analysis Tool, which is integral part of AliRoot,
14 // ALICE Off-Line framework:
16 // Piotr.Skowronski@cern.ch
17 // more info: http://alisoft.cern.ch/people/skowron/analyzer/index.html
19 ////////////////////////////////////////////////////////////////////////////
20 //_________________________________________________________
24 #include "AliAODParticle.h"
25 #include "AliAODPairCut.h"
27 #include "AliEventBuffer.h"
29 #include "AliReader.h"
30 #include "AliHBTPair.h"
31 #include "AliHBTFunction.h"
32 #include "AliHBTMonitorFunction.h"
36 ClassImp(AliHBTAnalysis)
38 const UInt_t AliHBTAnalysis::fgkFctnArraySize = 100;
39 const UInt_t AliHBTAnalysis::fgkDefaultMixingInfo = 1000;
40 const Int_t AliHBTAnalysis::fgkDefaultBufferSize = 5;
42 AliHBTAnalysis::AliHBTAnalysis():
46 fNParticleFunctions(0),
47 fNParticleAndTrackFunctions(0),
48 fNTrackMonitorFunctions(0),
49 fNParticleMonitorFunctions(0),
50 fNParticleAndTrackMonitorFunctions(0),
51 fTrackFunctions ( new AliHBTOnePairFctn* [fgkFctnArraySize]),
52 fParticleFunctions ( new AliHBTOnePairFctn* [fgkFctnArraySize]),
53 fParticleAndTrackFunctions ( new AliHBTTwoPairFctn* [fgkFctnArraySize]),
54 fParticleMonitorFunctions ( new AliHBTMonOneParticleFctn* [fgkFctnArraySize]),
55 fTrackMonitorFunctions ( new AliHBTMonOneParticleFctn* [fgkFctnArraySize]),
56 fParticleAndTrackMonitorFunctions ( new AliHBTMonTwoParticleFctn* [fgkFctnArraySize]),
57 fPairCut(new AliAODEmptyPairCut()),//empty cut - accepts all particles
59 fDisplayMixingInfo(fgkDefaultMixingInfo),
63 fProcessOption(kSimulatedAndReconstructed),
65 fkPass(&AliHBTAnalysis::PassPartAndTrack), //by default perform cut on both track and particle pair
66 fkPass1(&AliHBTAnalysis::PassPartAndTrack1), //used onluy by ProcessTracksAndParticles
67 fkPass2(&AliHBTAnalysis::PassPartAndTrack2),
68 fkPassPairProp(&AliHBTAnalysis::PassPairPropPartAndTrack)
73 /*************************************************************************************/
75 AliHBTAnalysis::AliHBTAnalysis(const AliHBTAnalysis& in):
80 fNParticleFunctions(0),
81 fNParticleAndTrackFunctions(0),
82 fNTrackMonitorFunctions(0),
83 fNParticleMonitorFunctions(0),
84 fNParticleAndTrackMonitorFunctions(0),
86 fParticleFunctions(0x0),
87 fParticleAndTrackFunctions(0x0),
88 fParticleMonitorFunctions(0x0),
89 fTrackMonitorFunctions(0x0),
90 fParticleAndTrackMonitorFunctions(0x0),
92 fBufferSize(fgkDefaultBufferSize),
93 fDisplayMixingInfo(fgkDefaultMixingInfo),
97 fProcessOption(kSimulatedAndReconstructed),
99 fkPass(&AliHBTAnalysis::PassPartAndTrack), //by default perform cut on both track and particle pair
100 fkPass1(&AliHBTAnalysis::PassPartAndTrack1),
101 fkPass2(&AliHBTAnalysis::PassPartAndTrack2),
102 fkPassPairProp(&AliHBTAnalysis::PassPairPropPartAndTrack)
105 Fatal("AliHBTAnalysis(const AliHBTAnalysis&)","Sensless");
107 /*************************************************************************************/
108 AliHBTAnalysis& AliHBTAnalysis::operator=(const AliHBTAnalysis& /*right*/)
111 Fatal("AliHBTAnalysis(const AliHBTAnalysis&)","Sensless");
114 /*************************************************************************************/
115 AliHBTAnalysis::~AliHBTAnalysis()
118 //note that we do not delete functions itself
119 // they should be deleted by whom where created
120 //we only store pointers, and we use "new" only for pointers array
124 if (AliVAODParticle::GetDebug()>5)Info("~AliHBTAnalysis","Is Owner: Attempting to delete functions");
126 if (AliVAODParticle::GetDebug()>5)Info("~AliHBTAnalysis","Delete functions done");
128 delete [] fTrackFunctions;
129 delete [] fParticleFunctions;
130 delete [] fParticleAndTrackFunctions;
132 delete [] fParticleMonitorFunctions;
133 delete [] fTrackMonitorFunctions;
134 delete [] fParticleAndTrackMonitorFunctions;
136 delete fPairCut; // always have an copy of an object - we create we dstroy
139 /*************************************************************************************/
140 Int_t AliHBTAnalysis::ProcessEvent(AliAOD* aodrec, AliAOD* aodsim)
142 //Processes one event
143 if (fProcEvent == 0x0)
145 Error("ProcessEvent","Analysis option not specified");
148 if ( Pass(aodrec,aodsim) ) return 0;
150 return (this->*fProcEvent)(aodrec,aodsim);
152 /*************************************************************************************/
154 Int_t AliHBTAnalysis::Finish()
160 /*************************************************************************************/
162 void AliHBTAnalysis::DeleteFunctions()
164 //Deletes all functions added to analysis
166 for(ii = 0;ii<fNParticleFunctions;ii++)
168 if (AliVAODParticle::GetDebug()>5)
170 Info("DeleteFunctions","Deleting ParticleFunction %#x",fParticleFunctions[ii]);
171 Info("DeleteFunctions","Deleting ParticleFunction %s",fParticleFunctions[ii]->Name());
173 delete fParticleFunctions[ii];
175 fNParticleFunctions = 0;
177 for(ii = 0;ii<fNTrackFunctions;ii++)
179 if (AliVAODParticle::GetDebug()>5)
181 Info("DeleteFunctions","Deleting TrackFunction %#x",fTrackFunctions[ii]);
182 Info("DeleteFunctions","Deleting TrackFunction %s",fTrackFunctions[ii]->Name());
184 delete fTrackFunctions[ii];
186 fNTrackFunctions = 0;
188 for(ii = 0;ii<fNParticleAndTrackFunctions;ii++)
190 if (AliVAODParticle::GetDebug()>5)
192 Info("DeleteFunctions","Deleting ParticleAndTrackFunction %#x",fParticleAndTrackFunctions[ii]);
193 Info("DeleteFunctions","Deleting ParticleAndTrackFunction %s",fParticleAndTrackFunctions[ii]->Name());
195 delete fParticleAndTrackFunctions[ii];
197 fNParticleAndTrackFunctions = 0;
199 for(ii = 0; ii<fNParticleMonitorFunctions; ii++)
201 if (AliVAODParticle::GetDebug()>5)
203 Info("DeleteFunctions","Deleting ParticleMonitorFunction %#x",fParticleMonitorFunctions[ii]);
204 Info("DeleteFunctions","Deleting ParticleMonitorFunction %s",fParticleMonitorFunctions[ii]->Name());
206 delete fParticleMonitorFunctions[ii];
208 fNParticleMonitorFunctions = 0;
210 for(ii = 0; ii<fNTrackMonitorFunctions; ii++)
212 if (AliVAODParticle::GetDebug()>5)
214 Info("DeleteFunctions","Deleting TrackMonitorFunction %#x",fTrackMonitorFunctions[ii]);
215 Info("DeleteFunctions","Deleting TrackMonitorFunction %s",fTrackMonitorFunctions[ii]->Name());
217 delete fTrackMonitorFunctions[ii];
219 fNTrackMonitorFunctions = 0;
221 for(ii = 0; ii<fNParticleAndTrackMonitorFunctions; ii++)
223 if (AliVAODParticle::GetDebug()>5)
225 Info("DeleteFunctions","Deleting ParticleAndTrackMonitorFunction %#x",fParticleAndTrackMonitorFunctions[ii]);
226 Info("DeleteFunctions","Deleting ParticleAndTrackMonitorFunction %s",fParticleAndTrackMonitorFunctions[ii]->Name());
228 delete fParticleAndTrackMonitorFunctions[ii];
230 fNParticleAndTrackMonitorFunctions = 0;
233 /*************************************************************************************/
235 Int_t AliHBTAnalysis::Init()
237 //Initializeation method
238 //calls Init for all functions
240 //Depending on option and pair cut assigns proper analysis method
241 Bool_t nonid = IsNonIdentAnalysis();
242 switch(fProcessOption)
245 if (nonid) fProcEvent = &AliHBTAnalysis::ProcessRecNonId;
246 else fProcEvent = &AliHBTAnalysis::ProcessRec;
250 if (nonid) fProcEvent = &AliHBTAnalysis::ProcessSimNonId;
251 else fProcEvent = &AliHBTAnalysis::ProcessSim;
254 case kSimulatedAndReconstructed:
255 if (nonid) fProcEvent = &AliHBTAnalysis::ProcessRecAndSimNonId;
256 else fProcEvent = &AliHBTAnalysis::ProcessRecAndSim;
260 if (fPartBuffer == 0x0) fPartBuffer = new AliEventBuffer (fBufferSize);
261 else fPartBuffer->Reset();
263 if (fTrackBuffer == 0x0) fTrackBuffer = new AliEventBuffer (fBufferSize);
264 else fTrackBuffer->Reset();
267 fNoCorrfctns = (fNParticleFunctions == 0) && (fNTrackFunctions == 0) && (fNParticleAndTrackFunctions == 0);
270 for(ii = 0;ii<fNParticleFunctions;ii++)
271 fParticleFunctions[ii]->Init();
273 for(ii = 0;ii<fNTrackFunctions;ii++)
274 fTrackFunctions[ii]->Init();
276 for(ii = 0;ii<fNParticleAndTrackFunctions;ii++)
277 fParticleAndTrackFunctions[ii]->Init();
279 for(ii = 0; ii<fNParticleMonitorFunctions; ii++)
280 fParticleMonitorFunctions[ii]->Init();
282 for(ii = 0; ii<fNTrackMonitorFunctions; ii++)
283 fTrackMonitorFunctions[ii]->Init();
285 for(ii = 0; ii<fNParticleAndTrackMonitorFunctions; ii++)
286 fParticleAndTrackMonitorFunctions[ii]->Init();
290 /*************************************************************************************/
292 void AliHBTAnalysis::ResetFunctions()
294 //In case fOwner is true, deletes all functions
295 //in other case, just set number of analysis to 0
296 if (fIsOwner) DeleteFunctions();
299 fNParticleFunctions = 0;
300 fNTrackFunctions = 0;
301 fNParticleAndTrackFunctions = 0;
302 fNParticleMonitorFunctions = 0;
303 fNTrackMonitorFunctions = 0;
304 fNParticleAndTrackMonitorFunctions = 0;
307 /*************************************************************************************/
308 Int_t AliHBTAnalysis::ProcessRecAndSim(AliAOD* aodrec, AliAOD* aodsim)
310 //Does analysis for both tracks and particles
311 //mainly for resolution study and analysies with weighting algirithms
313 // cut on particles only -- why?
314 // - PID: when we make resolution analysis we want to take only tracks with correct PID
315 // We need cut on tracks because there are data characteristic
317 AliVAODParticle * part1, * part2;
318 AliVAODParticle * track1, * track2;
320 AliAOD * trackEvent = aodrec, *partEvent = aodsim;
321 AliAOD* trackEvent1 = new AliAOD();
322 AliAOD* partEvent1 = new AliAOD();
323 partEvent1->SetOwner(kTRUE);
324 trackEvent1->SetOwner(kTRUE);
326 AliAOD * trackEvent2,*partEvent2;
328 // Int_t N1, N2, N=0; //number of particles in current event(we prcess two events in one time)
330 // Int_t nev = fReader->GetNumberOfTrackEvents();
331 static AliHBTPair tpair;
332 static AliHBTPair ppair;
334 AliHBTPair* trackpair = &tpair;
335 AliHBTPair* partpair = &ppair;
337 AliHBTPair * tmptrackpair;//temprary pointers to pairs
338 AliHBTPair * tmppartpair;
344 if ( !partEvent || !trackEvent )
346 Error("ProcessRecAndSim","Can not get event");
350 if ( partEvent->GetNumberOfParticles() != trackEvent->GetNumberOfParticles() )
352 Error("ProcessRecAndSim",
353 "Number of simulated particles (%d) not equal to number of reconstructed tracks (%d). Skipping Event.",
354 partEvent->GetNumberOfParticles() , trackEvent->GetNumberOfParticles());
359 for (Int_t j = 0; j<partEvent->GetNumberOfParticles() ; j++)
361 /***************************************/
362 /****** Looping same events ********/
363 /****** filling numerators ********/
364 /***************************************/
365 if ( (j%fDisplayMixingInfo) == 0)
366 Info("ProcessTracksAndParticles",
367 "Mixing particle %d with particles from the same event",j);
369 part1= partEvent->GetParticle(j);
370 track1= trackEvent->GetParticle(j);
372 Bool_t firstcut = (this->*fkPass1)(part1,track1);
373 if (fBufferSize != 0)
374 if ( (firstcut == kFALSE) || ( (this->*fkPass2)(part1,track1) == kFALSE ) )
376 //accepted by any cut
377 // we have to copy because reader keeps only one event
379 partEvent1->AddParticle(new AliAODParticle(*part1));
380 trackEvent1->AddParticle(new AliAODParticle(*track1));
383 if (firstcut) continue;
385 for(ii = 0; ii<fNParticleMonitorFunctions; ii++)
386 fParticleMonitorFunctions[ii]->Process(part1);
387 for(ii = 0; ii<fNTrackMonitorFunctions; ii++)
388 fTrackMonitorFunctions[ii]->Process(track1);
389 for(ii = 0; ii<fNParticleAndTrackMonitorFunctions; ii++)
390 fParticleAndTrackMonitorFunctions[ii]->Process(track1,part1);
392 if (fNoCorrfctns) continue;
394 for (Int_t k =j+1; k < partEvent->GetNumberOfParticles() ; k++)
396 part2= partEvent->GetParticle(k);
397 if (part1->GetUID() == part2->GetUID()) continue;
398 partpair->SetParticles(part1,part2);
400 track2= trackEvent->GetParticle(k);
401 trackpair->SetParticles(track1,track2);
403 if( (this->*fkPass)(partpair,trackpair) ) //check pair cut
404 { //do not meets crietria of the pair cut, try with swapped pairs
405 if( (this->*fkPass)((AliHBTPair*)partpair->GetSwappedPair(),(AliHBTPair*)trackpair->GetSwappedPair()) )
406 continue; //swaped pairs do not meet criteria of pair cut as well, take next particle
408 { //swaped pair meets all the criteria
409 tmppartpair = (AliHBTPair*)partpair->GetSwappedPair();
410 tmptrackpair = (AliHBTPair*)trackpair->GetSwappedPair();
414 {//meets criteria of the pair cut
415 tmptrackpair = trackpair;
416 tmppartpair = partpair;
419 for(ii = 0;ii<fNParticleFunctions;ii++)
420 fParticleFunctions[ii]->ProcessSameEventParticles(tmppartpair);
422 for(ii = 0;ii<fNTrackFunctions;ii++)
423 fTrackFunctions[ii]->ProcessSameEventParticles(tmptrackpair);
425 for(ii = 0;ii<fNParticleAndTrackFunctions;ii++)
426 fParticleAndTrackFunctions[ii]->ProcessSameEventParticles(tmptrackpair,tmppartpair);
427 //end of 2nd loop over particles from the same event
428 }//for (Int_t k =j+1; k < partEvent->GetNumberOfParticles() ; k++)
430 /***************************************/
431 /***** Filling denominators *********/
432 /***************************************/
433 if (fBufferSize == 0) continue;
435 fPartBuffer->ResetIter();
436 fTrackBuffer->ResetIter();
438 while (( partEvent2 = fPartBuffer->Next() ))
440 trackEvent2 = fTrackBuffer->Next();
443 if ( (j%fDisplayMixingInfo) == 0)
444 Info("ProcessTracksAndParticles",
445 "Mixing particle %d from current event with particles from event %d",j,-m);
447 for(Int_t l = 0; l<partEvent2->GetNumberOfParticles();l++) // ... on all particles
449 part2= partEvent2->GetParticle(l);
450 partpair->SetParticles(part1,part2);
452 track2= trackEvent2->GetParticle(l);
453 trackpair->SetParticles(track1,track2);
455 if( (this->*fkPass)(partpair,trackpair) ) //check pair cut
456 { //do not meets crietria of the
457 if( (this->*fkPass)((AliHBTPair*)partpair->GetSwappedPair(),(AliHBTPair*)trackpair->GetSwappedPair()) )
461 tmppartpair = (AliHBTPair*)partpair->GetSwappedPair();
462 tmptrackpair = (AliHBTPair*)trackpair->GetSwappedPair();
466 {//meets criteria of the pair cut
467 tmptrackpair = trackpair;
468 tmppartpair = partpair;
471 for(ii = 0;ii<fNParticleFunctions;ii++)
472 fParticleFunctions[ii]->ProcessDiffEventParticles(tmppartpair);
474 for(ii = 0;ii<fNTrackFunctions;ii++)
475 fTrackFunctions[ii]->ProcessDiffEventParticles(tmptrackpair);
477 for(ii = 0;ii<fNParticleAndTrackFunctions;ii++)
478 fParticleAndTrackFunctions[ii]->ProcessDiffEventParticles(tmptrackpair,tmppartpair);
479 }//for(Int_t l = 0; l<N2;l++) // ... on all particles
482 //end of loop over particles from first event
483 }//for (Int_t j = 0; j<partEvent->GetNumberOfParticles() ; j++)
484 delete fPartBuffer->Push(partEvent1);
485 delete fTrackBuffer->Push(trackEvent1);
486 //end of loop over events
489 /*************************************************************************************/
491 Int_t AliHBTAnalysis::ProcessSim(AliAOD* /*aodrec*/, AliAOD* aodsim)
493 //Does analysis of simulated data
494 AliVAODParticle * part1, * part2;
496 AliAOD* partEvent = aodsim;
497 AliAOD* partEvent1 = new AliAOD();
498 partEvent1->SetOwner(kTRUE);
504 AliHBTPair* partpair = &ppair;
506 AliHBTPair * tmppartpair;
514 Error("ProcessRecAndSim","Can not get event");
519 for (Int_t j = 0; j<partEvent->GetNumberOfParticles() ; j++)
521 /***************************************/
522 /****** Looping same events ********/
523 /****** filling numerators ********/
524 /***************************************/
525 if ( (j%fDisplayMixingInfo) == 0)
526 Info("ProcessTracksAndParticles",
527 "Mixing particle %d with particles from the same event",j);
529 part1= partEvent->GetParticle(j);
531 Bool_t firstcut = fPairCut->GetFirstPartCut()->Pass(part1);
533 if (fBufferSize != 0)
534 if ( (firstcut == kFALSE) || ( fPairCut->GetSecondPartCut()->Pass(part1) == kFALSE ) )
536 //accepted by any cut
537 // we have to copy because reader keeps only one event
539 partEvent1->AddParticle(new AliAODParticle(*part1));
542 if (firstcut) continue;
544 for(ii = 0; ii<fNParticleMonitorFunctions; ii++)
545 fParticleMonitorFunctions[ii]->Process(part1);
547 if ( fNParticleFunctions == 0 ) continue;
549 for (Int_t k =j+1; k < partEvent->GetNumberOfParticles() ; k++)
551 part2= partEvent->GetParticle(k);
552 if (part1->GetUID() == part2->GetUID()) continue;
553 partpair->SetParticles(part1,part2);
555 if(fPairCut->Pass(partpair)) //check pair cut
556 { //do not meets crietria of the
557 if( fPairCut->Pass((AliHBTPair*)partpair->GetSwappedPair()) ) continue;
558 else tmppartpair = (AliHBTPair*)partpair->GetSwappedPair();
562 tmppartpair = partpair;
565 for(ii = 0;ii<fNParticleFunctions;ii++)
566 fParticleFunctions[ii]->ProcessSameEventParticles(tmppartpair);
568 //end of 2nd loop over particles from the same event
569 }//for (Int_t k =j+1; k < partEvent->GetNumberOfParticles() ; k++)
571 /***************************************/
572 /***** Filling denominators *********/
573 /***************************************/
574 if (fBufferSize == 0) continue;
576 fPartBuffer->ResetIter();
578 while (( partEvent2 = fPartBuffer->Next() ))
581 if ( (j%fDisplayMixingInfo) == 0)
582 Info("ProcessParticles",
583 "Mixing particle %d from current event with particles from event %d",j,-m);
584 for(Int_t l = 0; l<partEvent2->GetNumberOfParticles();l++) // ... on all particles
587 part2= partEvent2->GetParticle(l);
588 partpair->SetParticles(part1,part2);
590 if( fPairCut->Pass(partpair) ) //check pair cut
591 { //do not meets crietria of the
592 if( fPairCut->Pass((AliHBTPair*)partpair->GetSwappedPair()) )
596 tmppartpair = (AliHBTPair*)partpair->GetSwappedPair();
600 {//meets criteria of the pair cut
601 tmppartpair = partpair;
604 for(ii = 0;ii<fNParticleFunctions;ii++)
605 fParticleFunctions[ii]->ProcessDiffEventParticles(tmppartpair);
607 }//for(Int_t l = 0; l<N2;l++) // ... on all particles
610 delete fPartBuffer->Push(partEvent1);
611 //end of loop over events
614 /*************************************************************************************/
615 Int_t AliHBTAnalysis::ProcessRec(AliAOD* aodrec, AliAOD* /*aodsim*/)
617 //Does analysis of reconstructed data
618 AliVAODParticle * track1, * track2;
620 AliAOD* trackEvent = aodrec;
621 AliAOD* trackEvent1 = new AliAOD();
622 trackEvent1->SetOwner(kTRUE);
628 AliHBTPair* trackpair = &tpair;
630 AliHBTPair * tmptrackpair;
637 Error("ProcessRecAndSim","Can not get event");
642 for (Int_t j = 0; j<trackEvent->GetNumberOfParticles() ; j++)
644 /***************************************/
645 /****** Looping same events ********/
646 /****** filling numerators ********/
647 /***************************************/
648 if ( (j%fDisplayMixingInfo) == 0)
649 Info("ProcessTracksAndParticles",
650 "Mixing Particle %d with Particles from the same event",j);
652 track1= trackEvent->GetParticle(j);
654 Bool_t firstcut = fPairCut->GetFirstPartCut()->Pass(track1);
656 if (fBufferSize != 0)
657 if ( (firstcut == kFALSE) || ( fPairCut->GetSecondPartCut()->Pass(track1) == kFALSE ) )
659 //accepted by any cut
660 // we have to copy because reader keeps only one event
662 trackEvent1->AddParticle(new AliAODParticle(*track1));
665 if (firstcut) continue;
667 for(ii = 0; ii<fNParticleMonitorFunctions; ii++)
668 fParticleMonitorFunctions[ii]->Process(track1);
670 if ( fNParticleFunctions == 0 ) continue;
672 for (Int_t k =j+1; k < trackEvent->GetNumberOfParticles() ; k++)
674 track2= trackEvent->GetParticle(k);
675 if (track1->GetUID() == track2->GetUID()) continue;
676 trackpair->SetParticles(track1,track2);
678 if(fPairCut->Pass(trackpair)) //check pair cut
679 { //do not meets crietria of the
680 if( fPairCut->Pass((AliHBTPair*)trackpair->GetSwappedPair()) ) continue;
681 else tmptrackpair = (AliHBTPair*)trackpair->GetSwappedPair();
685 tmptrackpair = trackpair;
688 for(ii = 0;ii<fNTrackFunctions;ii++)
689 fParticleFunctions[ii]->ProcessSameEventParticles(tmptrackpair);
691 //end of 2nd loop over Particles from the same event
692 }//for (Int_t k =j+1; k < trackEvent->GetNumberOfParticles() ; k++)
694 /***************************************/
695 /***** Filling denominators *********/
696 /***************************************/
697 if (fBufferSize == 0) continue;
699 fTrackBuffer->ResetIter();
701 while (( trackEvent2 = fTrackBuffer->Next() ))
704 if ( (j%fDisplayMixingInfo) == 0)
705 Info("ProcessParticles",
706 "Mixing Particle %d from current event with Particles from event %d",j,-m);
707 for(Int_t l = 0; l<trackEvent2->GetNumberOfParticles();l++) // ... on all Particles
710 track2= trackEvent2->GetParticle(l);
711 trackpair->SetParticles(track1,track2);
713 if( fPairCut->Pass(trackpair) ) //check pair cut
714 { //do not meets crietria of the
715 if( fPairCut->Pass((AliHBTPair*)trackpair->GetSwappedPair()) )
719 tmptrackpair = (AliHBTPair*)trackpair->GetSwappedPair();
723 {//meets criteria of the pair cut
724 tmptrackpair = trackpair;
727 for(ii = 0;ii<fNTrackFunctions;ii++)
728 fTrackFunctions[ii]->ProcessDiffEventParticles(tmptrackpair);
730 }//for(Int_t l = 0; l<N2;l++) // ... on all Particles
733 delete fTrackBuffer->Push(trackEvent1);
734 //end of loop over events
737 /*************************************************************************************/
739 Int_t AliHBTAnalysis::ProcessRecAndSimNonId(AliAOD* aodrec, AliAOD* aodsim)
741 //Analyzes both reconstructed and simulated data
744 Error("ProcessTracksAndParticlesNonIdentAnal","Reconstructed event is NULL");
750 Error("ProcessTracksAndParticlesNonIdentAnal","Simulated event is NULL");
754 if ( aodrec->GetNumberOfParticles() != aodsim->GetNumberOfParticles() )
756 Error("ProcessTracksAndParticlesNonIdentAnal",
757 "Number of simulated particles (%d) not equal to number of reconstructed tracks (%d)",
758 aodsim->GetNumberOfParticles() , aodrec->GetNumberOfParticles());
763 AliVAODParticle * part1, * part2;
764 AliVAODParticle * track1, * track2;
766 static AliAOD aodrec1;
767 static AliAOD aodsim1;
769 AliAOD * trackEvent1=&aodrec1,*partEvent1=&aodsim1;//Particle that passes first particle cut, this event
770 AliAOD * trackEvent2=0x0,*partEvent2=0x0;//Particle that passes second particle cut, this event
771 AliAOD * trackEvent3=0x0,*partEvent3=0x0;//Particle that passes second particle cut, events from buffer
773 AliAOD* rawtrackEvent = aodrec;//this we get from Reader
774 AliAOD* rawpartEvent = aodsim;//this we get from Reader
776 static AliHBTPair tpair;
777 static AliHBTPair ppair;
779 AliHBTPair* trackpair = &tpair;
780 AliHBTPair* partpair = &ppair;
785 trackEvent1 = new AliAOD();
786 partEvent1 = new AliAOD();
787 trackEvent1->SetOwner(kFALSE);
788 partEvent1->SetOwner(kFALSE);;
790 /********************************/
792 /********************************/
793 if ( ( (partEvent2==0x0) || (trackEvent2==0x0)) )//in case fBufferSize == 0 and pointers are created do not eneter
795 partEvent2 = new AliAOD();
796 trackEvent2 = new AliAOD();
797 partEvent2->SetOwner(kTRUE);
798 trackEvent2->SetOwner(kTRUE);
801 FilterOut(partEvent1, partEvent2, rawpartEvent, trackEvent1, trackEvent2, rawtrackEvent);
803 for (Int_t j = 0; j<partEvent1->GetNumberOfParticles() ; j++)
805 if ( (j%fDisplayMixingInfo) == 0)
806 Info("ProcessTracksAndParticlesNonIdentAnal",
807 "Mixing particle %d from current event with particles from current event",j);
809 part1= partEvent1->GetParticle(j);
810 track1= trackEvent1->GetParticle(j);
813 for(ii = 0; ii<fNParticleMonitorFunctions; ii++)
814 fParticleMonitorFunctions[ii]->Process(part1);
815 for(ii = 0; ii<fNTrackMonitorFunctions; ii++)
816 fTrackMonitorFunctions[ii]->Process(track1);
817 for(ii = 0; ii<fNParticleAndTrackMonitorFunctions; ii++)
818 fParticleAndTrackMonitorFunctions[ii]->Process(track1,part1);
820 if (fNoCorrfctns) continue;
822 /***************************************/
823 /****** filling numerators ********/
824 /****** (particles from event2) ********/
825 /***************************************/
827 for (Int_t k = 0; k < partEvent2->GetNumberOfParticles() ; k++) //partEvent1 and partEvent2 are particles from the same event but separated to two groups
829 part2= partEvent2->GetParticle(k);
830 if (part1->GetUID() == part2->GetUID()) continue;//this is the same particle but with different PID
831 partpair->SetParticles(part1,part2);
833 track2= trackEvent2->GetParticle(k);
834 trackpair->SetParticles(track1,track2);
836 if( (this->*fkPassPairProp)(partpair,trackpair) ) //check pair cut
837 { //do not meets crietria of the pair cut
841 {//meets criteria of the pair cut
842 for(ii = 0;ii<fNParticleFunctions;ii++)
843 fParticleFunctions[ii]->ProcessSameEventParticles(partpair);
845 for(ii = 0;ii<fNTrackFunctions;ii++)
846 fTrackFunctions[ii]->ProcessSameEventParticles(trackpair);
848 for(ii = 0;ii<fNParticleAndTrackFunctions;ii++)
849 fParticleAndTrackFunctions[ii]->ProcessSameEventParticles(trackpair,partpair);
853 if ( fBufferSize == 0) continue;//do not mix diff histograms
854 /***************************************/
855 /***** Filling denominators *********/
856 /***************************************/
857 fPartBuffer->ResetIter();
858 fTrackBuffer->ResetIter();
862 while ( (partEvent3 = fPartBuffer->Next() ) != 0x0)
864 trackEvent3 = fTrackBuffer->Next();
866 if ( (j%fDisplayMixingInfo) == 0)
867 Info("ProcessTracksAndParticlesNonIdentAnal",
868 "Mixing particle %d from current event with particles from event%d",j,-(++nmonitor));
870 for (Int_t k = 0; k < partEvent3->GetNumberOfParticles() ; k++)
872 part2= partEvent3->GetParticle(k);
873 partpair->SetParticles(part1,part2);
875 track2= trackEvent3->GetParticle(k);
876 trackpair->SetParticles(track1,track2);
878 if( (this->*fkPassPairProp)(partpair,trackpair) ) //check pair cut
879 { //do not meets crietria of the pair cut
883 {//meets criteria of the pair cut
885 for(ii = 0;ii<fNParticleFunctions;ii++)
886 fParticleFunctions[ii]->ProcessDiffEventParticles(partpair);
888 for(ii = 0;ii<fNTrackFunctions;ii++)
889 fTrackFunctions[ii]->ProcessDiffEventParticles(trackpair);
891 for(ii = 0;ii<fNParticleAndTrackFunctions;ii++)
892 fParticleAndTrackFunctions[ii]->ProcessDiffEventParticles(trackpair,partpair);
894 }// for particles event2
896 }//for over particles in event1
898 delete fPartBuffer->Push(partEvent2);
899 delete fTrackBuffer->Push(trackEvent2);
903 /*************************************************************************************/
904 Int_t AliHBTAnalysis::ProcessSimNonId(AliAOD* /*aodrec*/, AliAOD* aodsim)
906 //does analysis of simulated (MC) data in non-identical mode
907 //i.e. when particles selected by first part. cut are a disjunctive set than particles
908 //passed by the second part. cut
915 AliVAODParticle * part1, * part2;
917 static AliAOD aodsim1;
919 AliAOD* partEvent1=&aodsim1;//Particle that passes first particle cut, this event
920 AliAOD* partEvent2=0x0;//Particle that passes second particle cut, this event
921 AliAOD* partEvent3=0x0;//Particle that passes second particle cut, events from buffer
923 AliAOD* rawpartEvent = aodsim;//this we get from Reader
925 static AliHBTPair ppair;
927 AliHBTPair* partpair = &ppair;
932 partEvent1 = new AliAOD();
933 partEvent1->SetOwner(kFALSE);;
936 /********************************/
938 /********************************/
939 if (partEvent2==0x0)//in case fBufferSize == 0 and pointers are created do not eneter
941 partEvent2 = new AliAOD();
942 partEvent2->SetOwner(kTRUE);
945 FilterOut(partEvent1, partEvent2, rawpartEvent);
947 for (Int_t j = 0; j<partEvent1->GetNumberOfParticles() ; j++)
949 if ( (j%fDisplayMixingInfo) == 0)
950 Info("ProcessParticlesNonIdentAnal",
951 "Mixing particle %d from current event with particles from current event",j);
953 part1= partEvent1->GetParticle(j);
956 for(ii = 0; ii<fNParticleMonitorFunctions; ii++)
957 fParticleMonitorFunctions[ii]->Process(part1);
959 if (fNParticleFunctions == 0) continue;
961 /***************************************/
962 /****** filling numerators ********/
963 /****** (particles from event2) ********/
964 /***************************************/
966 for (Int_t k = 0; k < partEvent2->GetNumberOfParticles() ; k++) //partEvent1 and partEvent2 are particles from the same event but separated to two groups
968 part2= partEvent2->GetParticle(k);
969 if (part1->GetUID() == part2->GetUID()) continue;//this is the same particle but with different PID
970 partpair->SetParticles(part1,part2);
973 if(fPairCut->PassPairProp(partpair) ) //check pair cut
974 { //do not meets crietria of the pair cut
978 {//meets criteria of the pair cut
979 for(ii = 0;ii<fNParticleFunctions;ii++)
980 fParticleFunctions[ii]->ProcessSameEventParticles(partpair);
984 if ( fBufferSize == 0) continue;//do not mix diff histograms
985 /***************************************/
986 /***** Filling denominators *********/
987 /***************************************/
988 fPartBuffer->ResetIter();
992 while ( (partEvent3 = fPartBuffer->Next() ) != 0x0)
995 if ( (j%fDisplayMixingInfo) == 0)
996 Info("ProcessParticlesNonIdentAnal",
997 "Mixing particle %d from current event with particles from event%d",j,-(++nmonitor));
999 for (Int_t k = 0; k < partEvent3->GetNumberOfParticles() ; k++)
1001 part2= partEvent3->GetParticle(k);
1002 partpair->SetParticles(part1,part2);
1005 if(fPairCut->PassPairProp(partpair) ) //check pair cut
1006 { //do not meets crietria of the pair cut
1010 {//meets criteria of the pair cut
1011 for(ii = 0;ii<fNParticleFunctions;ii++)
1013 fParticleFunctions[ii]->ProcessDiffEventParticles(partpair);
1016 }// for particles event2
1018 }//for over particles in event1
1020 delete fPartBuffer->Push(partEvent2);
1024 /*************************************************************************************/
1025 Int_t AliHBTAnalysis::ProcessRecNonId(AliAOD* aodrec, AliAOD* /*aodsim*/)
1027 //Analyzes both reconstructed and simulated data
1033 AliVAODParticle * track1, * track2;
1035 static AliAOD aodrec1;
1037 AliAOD * trackEvent1=&aodrec1;//Particle that passes first particle cut, this event
1038 AliAOD * trackEvent2=0x0;//Particle that passes second particle cut, this event
1039 AliAOD * trackEvent3=0x0;//Particle that passes second particle cut, events from buffer
1041 AliAOD* rawtrackEvent = aodrec;//this we get from Reader
1043 static AliHBTPair tpair;
1045 AliHBTPair* trackpair = &tpair;
1050 trackEvent1 = new AliAOD();
1051 trackEvent1->SetOwner(kFALSE);
1053 /********************************/
1055 /********************************/
1056 if ( trackEvent2==0x0 )//in case fBufferSize == 0 and pointers are created do not eneter
1058 trackEvent2 = new AliAOD();
1059 trackEvent2->SetOwner(kTRUE);
1062 FilterOut(trackEvent1, trackEvent2, rawtrackEvent);
1064 for (Int_t j = 0; j<trackEvent1->GetNumberOfParticles() ; j++)
1066 if ( (j%fDisplayMixingInfo) == 0)
1067 Info("ProcessTracksNonIdentAnal",
1068 "Mixing particle %d from current event with particles from current event",j);
1070 track1= trackEvent1->GetParticle(j);
1073 for(ii = 0; ii<fNTrackMonitorFunctions; ii++)
1074 fTrackMonitorFunctions[ii]->Process(track1);
1076 if (fNTrackFunctions == 0x0) continue;
1078 /***************************************/
1079 /****** filling numerators ********/
1080 /****** (particles from event2) ********/
1081 /***************************************/
1083 for (Int_t k = 0; k < trackEvent2->GetNumberOfParticles() ; k++) //partEvent1 and partEvent2 are particles from the same event but separated to two groups
1085 track2= trackEvent2->GetParticle(k);
1086 if (track1->GetUID() == track2->GetUID()) continue;//this is the same particle but with different PID
1087 trackpair->SetParticles(track1,track2);
1090 if( fPairCut->PassPairProp(trackpair)) //check pair cut
1091 { //do not meets crietria of the pair cut
1095 {//meets criteria of the pair cut
1097 for(ii = 0;ii<fNTrackFunctions;ii++)
1098 fTrackFunctions[ii]->ProcessSameEventParticles(trackpair);
1102 if ( fBufferSize == 0) continue;//do not mix diff histograms
1103 /***************************************/
1104 /***** Filling denominators *********/
1105 /***************************************/
1106 fTrackBuffer->ResetIter();
1110 while ( (trackEvent3 = fTrackBuffer->Next() ) != 0x0)
1112 if ( (j%fDisplayMixingInfo) == 0)
1113 Info("ProcessTracksNonIdentAnal",
1114 "Mixing particle %d from current event with particles from event%d",j,-(++nmonitor));
1116 for (Int_t k = 0; k < trackEvent3->GetNumberOfParticles() ; k++)
1118 track2= trackEvent3->GetParticle(k);
1119 trackpair->SetParticles(track1,track2);
1121 if( fPairCut->PassPairProp(trackpair)) //check pair cut
1122 { //do not meets crietria of the pair cut
1126 {//meets criteria of the pair cut
1127 for(ii = 0;ii<fNTrackFunctions;ii++)
1128 fTrackFunctions[ii]->ProcessDiffEventParticles(trackpair);
1130 }// for particles event2
1132 }//for over particles in event1
1134 delete fTrackBuffer->Push(trackEvent2);
1138 /*************************************************************************************/
1140 void AliHBTAnalysis::Process(Option_t* option)
1142 //default option = "TracksAndParticles"
1143 //Main method of the HBT Analysis Package
1144 //It triggers reading with the global cut (default is an empty cut)
1145 //Than it checks options and data which are read
1146 //if everything is OK, then it calls one of the looping methods
1147 //depending on tfReaderhe option
1148 //These methods differs on what thay are looping on
1151 //--------------------------------------------------------------------
1152 //ProcessTracksAndParticles - "TracksAndParticles"
1154 // it loops over both, tracks(reconstructed) and particles(simulated)
1155 // all function gethered in all 3 lists are called for each (double)pair
1157 //ProcessTracks - "Tracks"
1158 // it loops only on tracks(reconstructed),
1159 // functions ONLY from fTrackFunctions list are called
1161 //ProcessParticles - "Particles"
1162 // it loops only on particles(simulated),
1163 // functions ONLY from fParticleAndTrackFunctions list are called
1168 Error("Process","The reader is not set");
1172 const char *oT = strstr(option,"Tracks");
1173 const char *oP = strstr(option,"Particles");
1175 Bool_t nonid = IsNonIdentAnalysis();
1181 if (nonid) ProcessTracksAndParticlesNonIdentAnal();
1182 else ProcessTracksAndParticles();
1188 if (nonid) ProcessTracksNonIdentAnal();
1189 else ProcessTracks();
1195 if (nonid) ProcessParticlesNonIdentAnal();
1196 else ProcessParticles();
1201 /*************************************************************************************/
1203 void AliHBTAnalysis::ProcessTracksAndParticles()
1205 //Makes analysis for both tracks and particles
1206 //mainly for resolution study and analysies with weighting algirithms
1207 //In order to minimize calling AliRun::GetEvent (we need at one time particles from different events),
1208 //the loops are splited
1210 // cut on particles only -- why?
1211 // - PID: when we make resolution analysis we want to take only tracks with correct PID
1212 // We need cut on tracks because there are data characteristic to
1214 AliAOD * trackEvent, *partEvent;
1217 while (fReader->Next() == kFALSE)
1219 partEvent = fReader->GetEventSim();
1220 trackEvent = fReader->GetEventRec();
1221 ProcessRecAndSim(trackEvent,partEvent);
1223 }//while (fReader->Next() == kFALSE)
1226 /*************************************************************************************/
1228 void AliHBTAnalysis::ProcessTracks()
1230 //In order to minimize calling AliRun::GetEvent (we need at one time particles from different events),
1231 //the loops are splited
1232 AliAOD * trackEvent;
1234 while (fReader->Next() == kFALSE)
1236 trackEvent = fReader->GetEventRec();
1237 ProcessRec(trackEvent,0x0);
1238 }//while (fReader->Next() == kFALSE)
1241 /*************************************************************************************/
1243 void AliHBTAnalysis::ProcessParticles()
1245 //In order to minimize calling AliRun::GetEvent (we need at one time particles from different events),
1246 //the loops are splited
1249 while (fReader->Next() == kFALSE)
1251 partEvent = fReader->GetEventSim();
1252 ProcessSim(0x0,partEvent);
1253 }//while (fReader->Next() == kFALSE)
1255 /*************************************************************************************/
1257 void AliHBTAnalysis::WriteFunctions()
1259 //Calls Write for all defined functions in analysis
1260 //== writes all results
1262 for(ii = 0;ii<fNParticleFunctions;ii++)
1264 if (AliVAODParticle::GetDebug()>5)
1266 Info("WriteFunctions","Writing ParticleFunction %#x",fParticleFunctions[ii]);
1267 Info("WriteFunctions","Writing ParticleFunction %s",fParticleFunctions[ii]->Name());
1269 fParticleFunctions[ii]->Write();
1272 for(ii = 0;ii<fNTrackFunctions;ii++)
1274 if (AliVAODParticle::GetDebug()>5)
1276 Info("WriteFunctions","Writing TrackFunction %#x",fTrackFunctions[ii]);
1277 Info("WriteFunctions","Writing TrackFunction %s",fTrackFunctions[ii]->Name());
1279 fTrackFunctions[ii]->Write();
1282 for(ii = 0;ii<fNParticleAndTrackFunctions;ii++)
1284 if (AliVAODParticle::GetDebug()>5)
1286 Info("WriteFunctions","Writing ParticleAndTrackFunction %#x",fParticleAndTrackFunctions[ii]);
1287 Info("WriteFunctions","Writing ParticleAndTrackFunction %s",fParticleAndTrackFunctions[ii]->Name());
1289 fParticleAndTrackFunctions[ii]->Write();
1292 for(ii = 0;ii<fNParticleMonitorFunctions;ii++)
1294 if (AliVAODParticle::GetDebug()>5)
1296 Info("WriteFunctions","Writing ParticleMonitorFunction %#x",fParticleMonitorFunctions[ii]);
1297 Info("WriteFunctions","Writing ParticleMonitorFunction %s",fParticleMonitorFunctions[ii]->Name());
1299 fParticleMonitorFunctions[ii]->Write();
1302 for(ii = 0;ii<fNTrackMonitorFunctions;ii++)
1304 if (AliVAODParticle::GetDebug()>5)
1306 Info("WriteFunctions","Writing TrackMonitorFunction %#x",fTrackMonitorFunctions[ii]);
1307 Info("WriteFunctions","Writing TrackMonitorFunction %s",fTrackMonitorFunctions[ii]->Name());
1309 fTrackMonitorFunctions[ii]->Write();
1312 for(ii = 0;ii<fNParticleAndTrackMonitorFunctions;ii++)
1314 if (AliVAODParticle::GetDebug()>5)
1316 Info("WriteFunctions","Writing ParticleAndTrackMonitorFunction %#x",fParticleAndTrackMonitorFunctions[ii]);
1317 Info("WriteFunctions","Writing ParticleAndTrackMonitorFunction %s",fParticleAndTrackMonitorFunctions[ii]->Name());
1319 fParticleAndTrackMonitorFunctions[ii]->Write();
1322 /*************************************************************************************/
1324 void AliHBTAnalysis::SetGlobalPairCut(AliAODPairCut* cut)
1326 //Sets the global cut
1329 Error("AliHBTAnalysis::SetGlobalPairCut","Pointer is NULL. Ignoring");
1332 fPairCut = (AliAODPairCut*)cut->Clone();
1335 /*************************************************************************************/
1337 void AliHBTAnalysis::AddTrackFunction(AliHBTOnePairFctn* f)
1339 //Adds track function
1340 if (f == 0x0) return;
1341 if (fNTrackFunctions == fgkFctnArraySize)
1343 Error("AliHBTAnalysis::AddTrackFunction","Can not add this function, not enough place in the array.");
1345 fTrackFunctions[fNTrackFunctions] = f;
1348 /*************************************************************************************/
1350 void AliHBTAnalysis::AddParticleFunction(AliHBTOnePairFctn* f)
1352 //adds particle function
1353 if (f == 0x0) return;
1355 if (fNParticleFunctions == fgkFctnArraySize)
1357 Error("AliHBTAnalysis::AddParticleFunction","Can not add this function, not enough place in the array.");
1359 fParticleFunctions[fNParticleFunctions] = f;
1360 fNParticleFunctions++;
1362 /*************************************************************************************/
1364 void AliHBTAnalysis::AddParticleAndTrackFunction(AliHBTTwoPairFctn* f)
1366 //add resolution function
1367 if (f == 0x0) return;
1368 if (fNParticleAndTrackFunctions == fgkFctnArraySize)
1370 Error("AliHBTAnalysis::AddParticleAndTrackFunction","Can not add this function, not enough place in the array.");
1372 fParticleAndTrackFunctions[fNParticleAndTrackFunctions] = f;
1373 fNParticleAndTrackFunctions++;
1375 /*************************************************************************************/
1377 void AliHBTAnalysis::AddParticleMonitorFunction(AliHBTMonOneParticleFctn* f)
1379 //add particle monitoring function
1380 if (f == 0x0) return;
1382 if (fNParticleMonitorFunctions == fgkFctnArraySize)
1384 Error("AliHBTAnalysis::AddParticleMonitorFunction","Can not add this function, not enough place in the array.");
1386 fParticleMonitorFunctions[fNParticleMonitorFunctions] = f;
1387 fNParticleMonitorFunctions++;
1389 /*************************************************************************************/
1391 void AliHBTAnalysis::AddTrackMonitorFunction(AliHBTMonOneParticleFctn* f)
1393 //add track monitoring function
1394 if (f == 0x0) return;
1396 if (fNTrackMonitorFunctions == fgkFctnArraySize)
1398 Error("AliHBTAnalysis::AddTrackMonitorFunction","Can not add this function, not enough place in the array.");
1400 fTrackMonitorFunctions[fNTrackMonitorFunctions] = f;
1401 fNTrackMonitorFunctions++;
1403 /*************************************************************************************/
1405 void AliHBTAnalysis::AddParticleAndTrackMonitorFunction(AliHBTMonTwoParticleFctn* f)
1407 //add resolution monitoring function
1408 if (f == 0x0) return;
1409 if (fNParticleAndTrackMonitorFunctions == fgkFctnArraySize)
1411 Error("AliHBTAnalysis::AddParticleAndTrackMonitorFunction","Can not add this function, not enough place in the array.");
1413 fParticleAndTrackMonitorFunctions[fNParticleAndTrackMonitorFunctions] = f;
1414 fNParticleAndTrackMonitorFunctions++;
1418 /*************************************************************************************/
1419 /*************************************************************************************/
1421 Bool_t AliHBTAnalysis::RunCoherencyCheck()
1423 //Checks if both HBTRuns are similar
1424 //return true if error found
1425 //if they seem to be OK return false
1429 Info("RunCoherencyCheck","Checking HBT Runs Coherency");
1431 //When we use non-buffering reader this is a big waste of time -> We need to read all data to check it
1432 //and reader is implemented safe in this case anyway
1433 // Info("RunCoherencyCheck","Number of events ...");
1434 // if (fReader->GetNumberOfPartEvents() == fReader->GetNumberOfTrackEvents() ) //check whether there is the same number of events
1436 // Info("RunCoherencyCheck","OK. %d found\n",fReader->GetNumberOfTrackEvents());
1439 // { //if not the same - ERROR
1440 // Error("RunCoherencyCheck",
1441 // "Number of simulated events (%d) is not equal to number of reconstructed events(%d)",
1442 // fReader->GetNumberOfPartEvents(),fReader->GetNumberOfTrackEvents());
1446 Info("RunCoherencyCheck","Checking number of Particles AND Particles Types in each event ...");
1450 for( i = 0; i<fReader->GetNumberOfTrackEvents();i++)
1452 partEvent= fReader->GetEventSim(i); //gets the "ith" event
1453 trackEvent = fReader->GetEventRec(i);
1455 if ( (partEvent == 0x0) && (partEvent == 0x0) ) continue;
1456 if ( (partEvent == 0x0) || (partEvent == 0x0) )
1458 Error("RunCoherencyCheck",
1459 "One event is NULL and the other one not. Event Number %d",i);
1463 if ( partEvent->GetNumberOfParticles() != trackEvent->GetNumberOfParticles() )
1465 Error("RunCoherencyCheck",
1466 "Event %d: Number of simulated particles (%d) not equal to number of reconstructed tracks (%d)",
1467 i,partEvent->GetNumberOfParticles() , trackEvent->GetNumberOfParticles());
1471 for (Int_t j = 0; j<partEvent->GetNumberOfParticles(); j++)
1473 if( partEvent->GetParticle(j)->GetPdgCode() != trackEvent->GetParticle(j)->GetPdgCode() )
1475 Error("RunCoherencyCheck",
1476 "Event %d: Particle %d: PID of simulated particle (%d) not the same of reconstructed track (%d)",
1477 i,j, partEvent->GetParticle(j)->GetPdgCode(),trackEvent->GetParticle(j)->GetPdgCode() );
1483 Info("RunCoherencyCheck"," Done");
1484 Info("RunCoherencyCheck"," Everything looks OK");
1489 /*************************************************************************************/
1491 void AliHBTAnalysis::ProcessTracksAndParticlesNonIdentAnal()
1493 //Performs analysis for both, tracks and particles
1494 AliAOD* rawtrackEvent, * rawpartEvent;
1497 Info("ProcessTracksAndParticlesNonIdentAnal","**************************************");
1498 Info("ProcessTracksAndParticlesNonIdentAnal","***** NON IDENT MODE ****************");
1499 Info("ProcessTracksAndParticlesNonIdentAnal","**************************************");
1501 for (Int_t i = 0;;i++)//infinite loop
1503 if (fReader->Next()) break; //end when no more events available
1505 rawpartEvent = fReader->GetEventSim();
1506 rawtrackEvent = fReader->GetEventRec();
1508 ProcessRecAndSimNonId(rawtrackEvent,rawpartEvent);
1509 }//end of loop over events (1)
1511 /*************************************************************************************/
1513 void AliHBTAnalysis::ProcessTracksNonIdentAnal()
1515 //Process Tracks only with non identical mode
1516 AliAOD * rawtrackEvent;
1519 Info("ProcessTracksNonIdentAnal","**************************************");
1520 Info("ProcessTracksNonIdentAnal","***** NON IDENT MODE ****************");
1521 Info("ProcessTracksNonIdentAnal","**************************************");
1523 for (Int_t i = 0;;i++)//infinite loop
1525 if (fReader->Next()) break; //end when no more events available
1526 rawtrackEvent = fReader->GetEventRec();
1527 ProcessRecNonId(rawtrackEvent,0x0);
1528 }//end of loop over events (1)
1530 /*************************************************************************************/
1532 void AliHBTAnalysis::ProcessParticlesNonIdentAnal()
1534 //process paricles only with non identical mode
1535 AliAOD * rawpartEvent = 0x0;
1538 Info("ProcessParticlesNonIdentAnal","**************************************");
1539 Info("ProcessParticlesNonIdentAnal","***** NON IDENT MODE ****************");
1540 Info("ProcessParticlesNonIdentAnal","**************************************");
1542 for (Int_t i = 0;;i++)//infinite loop
1544 if (fReader->Next()) break; //end when no more events available
1546 rawpartEvent = fReader->GetEventSim();
1547 ProcessSimNonId(0x0,rawpartEvent);
1548 }//end of loop over events (1)
1551 /*************************************************************************************/
1552 void AliHBTAnalysis::FilterOut(AliAOD* outpart1, AliAOD* outpart2, AliAOD* inpart,
1553 AliAOD* outtrack1, AliAOD* outtrack2, AliAOD* intrack) const
1555 //Puts particles accepted as a first particle by global cut in out1
1556 //and as a second particle in out2
1558 AliVAODParticle* part, *track;
1567 for (Int_t i = 0; i < inpart->GetNumberOfParticles(); i++)
1570 part = inpart->GetParticle(i);
1571 track = intrack->GetParticle(i);
1573 if ( ((this->*fkPass1)(part,track)) ) in1 = kFALSE; //if part is rejected by cut1, in1 is false
1574 if ( ((this->*fkPass2)(part,track)) ) in2 = kFALSE; //if part is rejected by cut2, in2 is false
1576 if (gDebug)//to be removed in real analysis
1577 if ( in1 && in2 ) //both cuts accepted, should never happen, just in case
1579 //Particle accpted by both cuts
1580 Error("FilterOut","Particle accepted by both cuts");
1586 outpart1->AddParticle(part);
1587 outtrack1->AddParticle(track);
1593 outpart2->AddParticle(new AliAODParticle(*part));
1594 outtrack2->AddParticle(new AliAODParticle(*track));
1599 /*************************************************************************************/
1601 void AliHBTAnalysis::FilterOut(AliAOD* out1, AliAOD* out2, AliAOD* in) const
1603 //Puts particles accepted as a first particle by global cut in out1
1604 //and as a second particle in out2
1605 AliVAODParticle* part;
1610 AliAODParticleCut *cut1 = fPairCut->GetFirstPartCut();
1611 AliAODParticleCut *cut2 = fPairCut->GetSecondPartCut();
1615 for (Int_t i = 0; i < in->GetNumberOfParticles(); i++)
1618 part = in->GetParticle(i);
1620 if ( cut1->Pass(part) ) in1 = kFALSE; //if part is rejected by cut1, in1 is false
1621 if ( cut2->Pass(part) ) in2 = kFALSE; //if part is rejected by cut2, in2 is false
1623 if (gDebug)//to be removed in real analysis
1624 if ( in1 && in2 ) //both cuts accepted, should never happen, just in case
1626 //Particle accpted by both cuts
1627 Error("FilterOut","Particle accepted by both cuts");
1633 out1->AddParticle(part);
1639 out2->AddParticle(part);
1644 /*************************************************************************************/
1646 Bool_t AliHBTAnalysis::IsNonIdentAnalysis()
1648 //checks if it is possible to use special analysis for non identical particles
1649 //it means - in global pair cut first particle id is different than second one
1650 //and both are different from 0
1651 //in the future is possible to perform more sophisticated check
1652 //if cuts have excluding requirements
1654 if (fPairCut->IsEmpty())
1657 if (fPairCut->GetFirstPartCut()->IsEmpty())
1660 if (fPairCut->GetSecondPartCut()->IsEmpty())
1663 Int_t id1 = fPairCut->GetFirstPartCut()->GetPID();
1664 Int_t id2 = fPairCut->GetSecondPartCut()->GetPID();
1666 if ( (id1==0) || (id2==0) || (id1==id2) )
1671 /*************************************************************************************/
1673 void AliHBTAnalysis::SetCutsOnParticles()
1675 // -- aplies only to Process("TracksAndParticles")
1676 // (ProcessTracksAndParticles and ProcessTracksAndParticlesNonIdentAnal)
1677 // Only particles properties are checkes against cuts
1678 fkPass = &AliHBTAnalysis::PassPart;
1679 fkPass1 = &AliHBTAnalysis::PassPart1;
1680 fkPass2 = &AliHBTAnalysis::PassPart2;
1681 fkPassPairProp = &AliHBTAnalysis::PassPairPropPart;
1684 /*************************************************************************************/
1686 void AliHBTAnalysis::SetCutsOnTracks()
1688 // -- aplies only to Process("TracksAndParticles")
1689 // (ProcessTracksAndParticles and ProcessTracksAndParticlesNonIdentAnal)
1690 // Only tracks properties are checkes against cuts
1691 Info("SetCutsOnTracks","Only reconstructed particles will be checked");
1692 fkPass = &AliHBTAnalysis::PassTrack;
1693 fkPass1 = &AliHBTAnalysis::PassTrack1;
1694 fkPass2 = &AliHBTAnalysis::PassTrack2;
1695 fkPassPairProp = &AliHBTAnalysis::PassPairPropTrack;
1698 /*************************************************************************************/
1700 void AliHBTAnalysis::SetCutsOnTracksAndParticles()
1702 // -- aplies only to Process("TracksAndParticles")
1703 // (ProcessTracksAndParticles and ProcessTracksAndParticlesNonIdentAnal)
1704 // Both, tracks and particles, properties are checked against cuts
1705 fkPass = &AliHBTAnalysis::PassPartAndTrack;
1706 fkPass1 = &AliHBTAnalysis::PassPartAndTrack1;
1707 fkPass2 = &AliHBTAnalysis::PassPartAndTrack2;
1708 fkPassPairProp = &AliHBTAnalysis::PassPairPropPartAndTrack;
1710 /*************************************************************************************/
1712 void AliHBTAnalysis::PressAnyKey()
1714 //small utility function that helps to make comfortable macros
1717 fcntl(0, F_SETFL, O_NONBLOCK);
1718 ::Info("","Press Any Key to continue ...");
1721 nread = read(0, &c, 1);
1722 gSystem->ProcessEvents();
1726 /*************************************************************************************/