Add global trigger in RecTriggerTrack (Christian)
[u/mrichter/AliRoot.git] / HBTAN / AliHBTAnalysis.cxx
CommitLineData
bfb09ece 1//_________________________________________________________
2////////////////////////////////////////////////////////////////////////////
3//
4// class AliHBTAnalysis
5//
6// Central Object Of HBTAnalyser:
7// This class performs main looping within HBT Analysis
8// User must plug a reader of Type AliHBTReader
9// User plugs in coorelation and monitor functions
10// as well as monitor functions
11//
12// HBT Analysis Tool, which is integral part of AliRoot,
13// ALICE Off-Line framework:
14//
15// Piotr.Skowronski@cern.ch
16// more info: http://alisoft.cern.ch/people/skowron/analyzer/index.html
17//
18////////////////////////////////////////////////////////////////////////////
19//_________________________________________________________
20
1b446896 21#include "AliHBTAnalysis.h"
dc2c3f36 22#include "AliHBTEvent.h"
1b446896 23#include "AliHBTReader.h"
1b446896 24#include "AliHBTPair.h"
1b446896 25#include "AliHBTFunction.h"
5c58441a 26#include "AliHBTMonitorFunction.h"
bed069a4 27#include "AliHBTEventBuffer.h"
9616170a 28#include "AliHBTPairCut.h"
81b7b887 29
9616170a 30#include <TSystem.h>
e4f2b1da 31
1b446896 32ClassImp(AliHBTAnalysis)
33
34const UInt_t AliHBTAnalysis::fgkFctnArraySize = 100;
81b7b887 35const UInt_t AliHBTAnalysis::fgkDefaultMixingInfo = 1000;
36const Int_t AliHBTAnalysis::fgkDefaultBufferSize = 5;
37
38AliHBTAnalysis::AliHBTAnalysis():
2dc7203b 39 fReader(0x0),
40 fNTrackFunctions(0),
41 fNParticleFunctions(0),
42 fNParticleAndTrackFunctions(0),
43 fNTrackMonitorFunctions(0),
44 fNParticleMonitorFunctions(0),
45 fNParticleAndTrackMonitorFunctions(0),
9616170a 46 fTrackFunctions ( new AliHBTOnePairFctn* [fgkFctnArraySize]),
47 fParticleFunctions ( new AliHBTOnePairFctn* [fgkFctnArraySize]),
48 fParticleAndTrackFunctions ( new AliHBTTwoPairFctn* [fgkFctnArraySize]),
49 fParticleMonitorFunctions ( new AliHBTMonOneParticleFctn* [fgkFctnArraySize]),
50 fTrackMonitorFunctions ( new AliHBTMonOneParticleFctn* [fgkFctnArraySize]),
51 fParticleAndTrackMonitorFunctions ( new AliHBTMonTwoParticleFctn* [fgkFctnArraySize]),
52 fPairCut(new AliHBTEmptyPairCut()),//empty cut - accepts all particles
2dc7203b 53 fBufferSize(2),
54 fDisplayMixingInfo(fgkDefaultMixingInfo),
66d1d1a4 55 fIsOwner(kFALSE),
17d74b37 56 fkPass(&AliHBTAnalysis::PassPartAndTrack), //by default perform cut on both track and particle pair
57 fkPass1(&AliHBTAnalysis::PassPartAndTrack1), //used onluy by ProcessTracksAndParticles
58 fkPass2(&AliHBTAnalysis::PassPartAndTrack2),
59 fkPassPairProp(&AliHBTAnalysis::PassPairPropPartAndTrack)
1b446896 60 {
9616170a 61 //default constructor
66d1d1a4 62
1b446896 63 }
491d1b5d 64/*************************************************************************************/
1b446896 65
81b7b887 66AliHBTAnalysis::AliHBTAnalysis(const AliHBTAnalysis& in):
bdec021f 67 TObject(in),
2dc7203b 68 fReader(0x0),
69 fNTrackFunctions(0),
70 fNParticleFunctions(0),
71 fNParticleAndTrackFunctions(0),
72 fNTrackMonitorFunctions(0),
73 fNParticleMonitorFunctions(0),
74 fNParticleAndTrackMonitorFunctions(0),
75 fTrackFunctions(0x0),
76 fParticleFunctions(0x0),
77 fParticleAndTrackFunctions(0x0),
78 fParticleMonitorFunctions(0x0),
79 fTrackMonitorFunctions(0x0),
80 fParticleAndTrackMonitorFunctions(0x0),
81 fPairCut(0x0),
82 fBufferSize(fgkDefaultBufferSize),
83 fDisplayMixingInfo(fgkDefaultMixingInfo),
66d1d1a4 84 fIsOwner(kFALSE),
17d74b37 85 fkPass(&AliHBTAnalysis::PassPartAndTrack), //by default perform cut on both track and particle pair
86 fkPass1(&AliHBTAnalysis::PassPartAndTrack1),
87 fkPass2(&AliHBTAnalysis::PassPartAndTrack2),
88 fkPassPairProp(&AliHBTAnalysis::PassPairPropPartAndTrack)
81b7b887 89 {
2dc7203b 90//copy constructor
81b7b887 91 Fatal("AliHBTAnalysis(const AliHBTAnalysis&)","Sensless");
92 }
93/*************************************************************************************/
34914285 94AliHBTAnalysis& AliHBTAnalysis::operator=(const AliHBTAnalysis& /*right*/)
81b7b887 95 {
2dc7203b 96//operator =
81b7b887 97 Fatal("AliHBTAnalysis(const AliHBTAnalysis&)","Sensless");
98 return *this;
99 }
100/*************************************************************************************/
1b446896 101AliHBTAnalysis::~AliHBTAnalysis()
102 {
103 //destructor
104 //note that we do not delete functions itself
105 // they should be deleted by whom where created
106 //we only store pointers, and we use "new" only for pointers array
e7a04795 107
108 if (fIsOwner)
109 {
110 if (AliHBTParticle::GetDebug()>5)Info("~AliHBTAnalysis","Is Owner: Attempting to delete functions");
111 DeleteFunctions();
112 if (AliHBTParticle::GetDebug()>5)Info("~AliHBTAnalysis","Delete functions done");
113 }
1b446896 114 delete [] fTrackFunctions;
115 delete [] fParticleFunctions;
116 delete [] fParticleAndTrackFunctions;
117
5c58441a 118 delete [] fParticleMonitorFunctions;
119 delete [] fTrackMonitorFunctions;
120 delete [] fParticleAndTrackMonitorFunctions;
121
1b446896 122 delete fPairCut; // always have an copy of an object - we create we dstroy
123 }
124
125/*************************************************************************************/
e4f2b1da 126
81b7b887 127void AliHBTAnalysis::DeleteFunctions()
128{
e4f2b1da 129 //Deletes all functions added to analysis
81b7b887 130 UInt_t ii;
131 for(ii = 0;ii<fNParticleFunctions;ii++)
e7a04795 132 {
133 if (AliHBTParticle::GetDebug()>5)
134 Info("DeleteFunctions","Deleting ParticleFunction %s",fParticleFunctions[ii]->Name());
135 delete fParticleFunctions[ii];
136 }
e4f2b1da 137 fNParticleFunctions = 0;
81b7b887 138
139 for(ii = 0;ii<fNTrackFunctions;ii++)
e7a04795 140 {
141 if (AliHBTParticle::GetDebug()>5)
142 {
143 Info("DeleteFunctions","Deleting TrackFunction %#x",fParticleFunctions[ii]);
144 Info("DeleteFunctions","Deleting TrackFunction %s",fParticleFunctions[ii]->Name());
145 }
146 delete fTrackFunctions[ii];
147 }
e4f2b1da 148 fNTrackFunctions = 0;
149
81b7b887 150 for(ii = 0;ii<fNParticleAndTrackFunctions;ii++)
e7a04795 151 {
152 if (AliHBTParticle::GetDebug()>5)
153 Info("DeleteFunctions","Deleting ParticleAndTrackFunction %s",fParticleFunctions[ii]->Name());
154 delete fParticleAndTrackFunctions[ii];
155 }
e4f2b1da 156 fNParticleAndTrackFunctions = 0;
81b7b887 157
158 for(ii = 0; ii<fNParticleMonitorFunctions; ii++)
e7a04795 159 {
160 if (AliHBTParticle::GetDebug()>5)
161 Info("DeleteFunctions","Deleting ParticleMonitorFunction %s",fParticleFunctions[ii]->Name());
162 delete fParticleMonitorFunctions[ii];
163 }
e4f2b1da 164 fNParticleMonitorFunctions = 0;
81b7b887 165
166 for(ii = 0; ii<fNTrackMonitorFunctions; ii++)
e7a04795 167 {
168 if (AliHBTParticle::GetDebug()>5)
169 Info("DeleteFunctions","Deleting TrackMonitorFunction %s",fParticleFunctions[ii]->Name());
170 delete fTrackMonitorFunctions[ii];
171 }
e4f2b1da 172 fNTrackMonitorFunctions = 0;
81b7b887 173
174 for(ii = 0; ii<fNParticleAndTrackMonitorFunctions; ii++)
e7a04795 175 {
176 if (AliHBTParticle::GetDebug()>5)
177 Info("DeleteFunctions","Deleting ParticleAndTrackMonitorFunction %s",fParticleFunctions[ii]->Name());
178 delete fParticleAndTrackMonitorFunctions[ii];
179 }
e4f2b1da 180 fNParticleAndTrackMonitorFunctions = 0;
bfb09ece 181
81b7b887 182}
e4f2b1da 183/*************************************************************************************/
184
185void AliHBTAnalysis::Init()
186{
2dc7203b 187//Initializeation method
188//calls Init for all functions
e4f2b1da 189 UInt_t ii;
190 for(ii = 0;ii<fNParticleFunctions;ii++)
191 fParticleFunctions[ii]->Init();
192
193 for(ii = 0;ii<fNTrackFunctions;ii++)
194 fTrackFunctions[ii]->Init();
195
196 for(ii = 0;ii<fNParticleAndTrackFunctions;ii++)
197 fParticleAndTrackFunctions[ii]->Init();
198
199 for(ii = 0; ii<fNParticleMonitorFunctions; ii++)
200 fParticleMonitorFunctions[ii]->Init();
201
202 for(ii = 0; ii<fNTrackMonitorFunctions; ii++)
203 fTrackMonitorFunctions[ii]->Init();
204
205 for(ii = 0; ii<fNParticleAndTrackMonitorFunctions; ii++)
206 fParticleAndTrackMonitorFunctions[ii]->Init();
bfb09ece 207
208
e4f2b1da 209}
210/*************************************************************************************/
211
212void AliHBTAnalysis::ResetFunctions()
213{
214//In case fOwner is true, deletes all functions
215//in other case, just set number of analysis to 0
216 if (fIsOwner) DeleteFunctions();
217 else
218 {
219 fNParticleFunctions = 0;
220 fNTrackFunctions = 0;
221 fNParticleAndTrackFunctions = 0;
222 fNParticleMonitorFunctions = 0;
223 fNTrackMonitorFunctions = 0;
224 fNParticleAndTrackMonitorFunctions = 0;
225 }
226}
227/*************************************************************************************/
228
1b446896 229void AliHBTAnalysis::Process(Option_t* option)
230{
231 //default option = "TracksAndParticles"
232 //Main method of the HBT Analysis Package
233 //It triggers reading with the global cut (default is an empty cut)
234 //Than it checks options and data which are read
235 //if everything is OK, then it calls one of the looping methods
236 //depending on tfReaderhe option
237 //These methods differs on what thay are looping on
238 //
239 // METHOD OPTION
240 //--------------------------------------------------------------------
241 //ProcessTracksAndParticles - "TracksAndParticles"
242 // DEFAULT
243 // it loops over both, tracks(reconstructed) and particles(simulated)
244 // all function gethered in all 3 lists are called for each (double)pair
245 //
246 //ProcessTracks - "Tracks"
247 // it loops only on tracks(reconstructed),
248 // functions ONLY from fTrackFunctions list are called
249 //
250 //ProcessParticles - "Particles"
251 // it loops only on particles(simulated),
252 // functions ONLY from fParticleAndTrackFunctions list are called
253 //
254 //
255 if (!fReader)
256 {
01725374 257 Error("Process","The reader is not set");
1b446896 258 return;
259 }
260
1b446896 261 const char *oT = strstr(option,"Tracks");
262 const char *oP = strstr(option,"Particles");
263
dc2c3f36 264 Bool_t nonid = IsNonIdentAnalysis();
265
bed069a4 266 Init();
267
1b446896 268 if(oT && oP)
269 {
dc2c3f36 270 if (nonid) ProcessTracksAndParticlesNonIdentAnal();
271 else ProcessTracksAndParticles();
1b446896 272 return;
273 }
274
275 if(oT)
276 {
dc2c3f36 277 if (nonid) ProcessTracksNonIdentAnal();
278 else ProcessTracks();
1b446896 279 return;
280 }
281
282 if(oP)
283 {
dc2c3f36 284 if (nonid) ProcessParticlesNonIdentAnal();
285 else ProcessParticles();
1b446896 286 return;
287 }
288
289}
1b446896 290/*************************************************************************************/
491d1b5d 291
1b446896 292void AliHBTAnalysis::ProcessTracksAndParticles()
293{
9616170a 294//Makes analysis for both tracks and particles
295//mainly for resolution study and analysies with weighting algirithms
1b446896 296//In order to minimize calling AliRun::GetEvent (we need at one time particles from different events),
297//the loops are splited
298
66d1d1a4 299// cut on particles only -- why?
300// - PID: when we make resolution analysis we want to take only tracks with correct PID
301// We need cut on tracks because there are data characteristic to
1b446896 302
303 AliHBTParticle * part1, * part2;
304 AliHBTParticle * track1, * track2;
305
306 AliHBTEvent * trackEvent, *partEvent;
bed069a4 307 AliHBTEvent * trackEvent1 = 0x0,*partEvent1 = 0x0;
1b446896 308 AliHBTEvent * trackEvent2,*partEvent2;
309
310// Int_t N1, N2, N=0; //number of particles in current event(we prcess two events in one time)
311
bed069a4 312// Int_t nev = fReader->GetNumberOfTrackEvents();
1b446896 313 AliHBTPair * trackpair = new AliHBTPair();
314 AliHBTPair * partpair = new AliHBTPair();
491d1b5d 315
316 AliHBTPair * tmptrackpair;//temprary pointers to pairs
317 AliHBTPair * tmppartpair;
81b7b887 318
bed069a4 319 AliHBTEventBuffer partbuffer(fBufferSize);
320 AliHBTEventBuffer trackbuffer(fBufferSize);
321
81b7b887 322 register UInt_t ii;
7a2c8238 323
bed069a4 324 Bool_t nocorrfctns = (fNParticleFunctions == 0) && (fNTrackFunctions == 0) && (fNParticleAndTrackFunctions == 0);
325
326 fReader->Rewind();
327 Int_t i = -1;
328 while (fReader->Next() == kFALSE)
1b446896 329 {
bed069a4 330 i++;
331 partEvent= fReader->GetParticleEvent();
332 trackEvent = fReader->GetTrackEvent();
1b446896 333
bed069a4 334 if ( !partEvent || !trackEvent )
335 {
336 Error("ProcessTracksAndParticles","Can not get event");
337 return;
338 }
339
340 if ( partEvent->GetNumberOfParticles() != trackEvent->GetNumberOfParticles() )
341 {
342 Fatal("ProcessTracksAndParticles",
343 "Event %d: Number of simulated particles (%d) not equal to number of reconstructed tracks (%d)",
344 i,partEvent->GetNumberOfParticles() , trackEvent->GetNumberOfParticles());
345 return;
346 }
1b446896 347
bed069a4 348 if(partEvent1 == 0x0)
1b446896 349 {
bed069a4 350 partEvent1 = new AliHBTEvent();
351 partEvent1->SetOwner(kTRUE);
81b7b887 352
bed069a4 353 trackEvent1 = new AliHBTEvent();
354 trackEvent1->SetOwner(kTRUE);
355 }
356 else
357 {
358 partEvent1->Reset();
359 trackEvent1->Reset();
360 }
361
362 for (Int_t j = 0; j<partEvent->GetNumberOfParticles() ; j++)
363 {
364 /***************************************/
365 /****** Looping same events ********/
366 /****** filling numerators ********/
367 /***************************************/
368 if ( (j%fDisplayMixingInfo) == 0)
81b7b887 369 Info("ProcessTracksAndParticles",
370 "Mixing particle %d from event %d with particles from event %d",j,i,i);
1b446896 371
372 part1= partEvent->GetParticle(j);
bed069a4 373 track1= trackEvent->GetParticle(j);
4d91c73a 374
375//PID imperfections ???
376// if( part1->GetPdgCode() != track1->GetPdgCode() )
377// {
378// Fatal("ProcessTracksAndParticles",
379// "Event %d: Particle %d: PID of simulated particle (%d) not the same of reconstructed track (%d)",
380// i,j, part1->GetPdgCode(),track1->GetPdgCode() );
381// return;
382// }
bed069a4 383
17d74b37 384 Bool_t firstcut = (this->*fkPass1)(part1,track1);
9277d10a 385 if (fBufferSize != 0)
17d74b37 386 if ( (firstcut == kFALSE) || ( (this->*fkPass2)(part1,track1) == kFALSE ) )
47d9a058 387 {
388 //accepted by any cut
389 // we have to copy because reader keeps only one event
390
391 partEvent1->AddParticle(new AliHBTParticle(*part1));
392 trackEvent1->AddParticle(new AliHBTParticle(*track1));
393 }
81b7b887 394
bed069a4 395 if (firstcut) continue;
396
81b7b887 397 for(ii = 0; ii<fNParticleMonitorFunctions; ii++)
e4f2b1da 398 fParticleMonitorFunctions[ii]->Process(part1);
81b7b887 399 for(ii = 0; ii<fNTrackMonitorFunctions; ii++)
e4f2b1da 400 fTrackMonitorFunctions[ii]->Process(track1);
81b7b887 401 for(ii = 0; ii<fNParticleAndTrackMonitorFunctions; ii++)
e4f2b1da 402 fParticleAndTrackMonitorFunctions[ii]->Process(track1,part1);
5c58441a 403
bed069a4 404 if (nocorrfctns) continue;
5c58441a 405
1b446896 406 for (Int_t k =j+1; k < partEvent->GetNumberOfParticles() ; k++)
407 {
408 part2= partEvent->GetParticle(k);
d74e6a27 409 if (part1->GetUID() == part2->GetUID()) continue;
1b446896 410 partpair->SetParticles(part1,part2);
411
4d91c73a 412 track2= trackEvent->GetParticle(k);
1b446896 413 trackpair->SetParticles(track1,track2);
414
17d74b37 415 if( (this->*fkPass)(partpair,trackpair) ) //check pair cut
bed069a4 416 { //do not meets crietria of the pair cut, try with swapped pairs
17d74b37 417 if( (this->*fkPass)(partpair->GetSwapedPair(),trackpair->GetSwapedPair()) )
1b446896 418 continue; //swaped pairs do not meet criteria of pair cut as well, take next particle
419 else
bed069a4 420 { //swaped pair meets all the criteria
491d1b5d 421 tmppartpair = partpair->GetSwapedPair();
422 tmptrackpair = trackpair->GetSwapedPair();
1b446896 423 }
424 }
491d1b5d 425 else
426 {//meets criteria of the pair cut
427 tmptrackpair = trackpair;
428 tmppartpair = partpair;
429 }
9616170a 430
1b446896 431 for(ii = 0;ii<fNParticleFunctions;ii++)
491d1b5d 432 fParticleFunctions[ii]->ProcessSameEventParticles(tmppartpair);
1b446896 433
434 for(ii = 0;ii<fNTrackFunctions;ii++)
491d1b5d 435 fTrackFunctions[ii]->ProcessSameEventParticles(tmptrackpair);
1b446896 436
437 for(ii = 0;ii<fNParticleAndTrackFunctions;ii++)
491d1b5d 438 fParticleAndTrackFunctions[ii]->ProcessSameEventParticles(tmptrackpair,tmppartpair);
bed069a4 439 //end of 2nd loop over particles from the same event
440 }//for (Int_t k =j+1; k < partEvent->GetNumberOfParticles() ; k++)
441
442 /***************************************/
443 /***** Filling denominators *********/
444 /***************************************/
445 if (fBufferSize == 0) continue;
446
447 partbuffer.ResetIter();
448 trackbuffer.ResetIter();
449 Int_t m = 0;
450 while (( partEvent2 = partbuffer.Next() ))
451 {
452 trackEvent2 = trackbuffer.Next();
5c58441a 453
bed069a4 454 m++;
455 if ( (j%fDisplayMixingInfo) == 0)
456 Info("ProcessTracksAndParticles",
457 "Mixing particle %d from event %d with particles from event %d",j,i,i-m);
458
459 for(Int_t l = 0; l<partEvent2->GetNumberOfParticles();l++) // ... on all particles
1b446896 460 {
1b446896 461 part2= partEvent2->GetParticle(l);
462 partpair->SetParticles(part1,part2);
463
464 track2= trackEvent2->GetParticle(l);
465 trackpair->SetParticles(track1,track2);
466
17d74b37 467 if( (this->*fkPass)(partpair,trackpair) ) //check pair cut
1b446896 468 { //do not meets crietria of the
17d74b37 469 if( (this->*fkPass)(partpair->GetSwapedPair(),trackpair->GetSwapedPair()) )
81b7b887 470 continue;
471 else
472 {
473 tmppartpair = partpair->GetSwapedPair();
474 tmptrackpair = trackpair->GetSwapedPair();
475 }
1b446896 476 }
491d1b5d 477 else
478 {//meets criteria of the pair cut
479 tmptrackpair = trackpair;
480 tmppartpair = partpair;
481 }
9616170a 482
1b446896 483 for(ii = 0;ii<fNParticleFunctions;ii++)
491d1b5d 484 fParticleFunctions[ii]->ProcessDiffEventParticles(tmppartpair);
1b446896 485
486 for(ii = 0;ii<fNTrackFunctions;ii++)
491d1b5d 487 fTrackFunctions[ii]->ProcessDiffEventParticles(tmptrackpair);
1b446896 488
489 for(ii = 0;ii<fNParticleAndTrackFunctions;ii++)
491d1b5d 490 fParticleAndTrackFunctions[ii]->ProcessDiffEventParticles(tmptrackpair,tmppartpair);
1b446896 491 }//for(Int_t l = 0; l<N2;l++) // ... on all particles
bed069a4 492
493 }
494 //end of loop over particles from first event
495 }//for (Int_t j = 0; j<partEvent->GetNumberOfParticles() ; j++)
496 partEvent1 = partbuffer.Push(partEvent1);
497 trackEvent1 = trackbuffer.Push(trackEvent1);
498 //end of loop over events
499 }//while (fReader->Next() == kFALSE)
4cb6e8f7 500
501 delete trackpair;
502 delete partpair;
503
504 partbuffer.SetOwner(kTRUE);
505 trackbuffer.SetOwner(kTRUE);
1b446896 506}
507/*************************************************************************************/
508
509void AliHBTAnalysis::ProcessTracks()
510{
2dc7203b 511//In order to minimize calling AliRun::GetEvent (we need at one time particles from different events),
1b446896 512//the loops are splited
513 AliHBTParticle * track1, * track2;
1b446896 514 AliHBTEvent * trackEvent;
bed069a4 515 AliHBTEvent * trackEvent1 = 0x0;
1b446896 516 AliHBTEvent * trackEvent2;
81b7b887 517
bed069a4 518 register UInt_t ii;
81b7b887 519
520 AliHBTPair * trackpair = new AliHBTPair();
521 AliHBTPair * tmptrackpair; //temporary pointer
1b446896 522
bed069a4 523 AliHBTEventBuffer trackbuffer(fBufferSize);
524
525 fReader->Rewind();
526 Int_t i = -1;
527 while (fReader->Next() == kFALSE)
1b446896 528 {
bed069a4 529 i++;
530 trackEvent = fReader->GetTrackEvent();
1b446896 531 if (!trackEvent) continue;
1b446896 532
bed069a4 533 if(trackEvent1 == 0x0)
534 {
535 trackEvent1 = new AliHBTEvent();
536 trackEvent1->SetOwner(kTRUE);
537 }
538 else
539 {
540 trackEvent1->Reset();
541 }
542
1b446896 543 for (Int_t j = 0; j<trackEvent->GetNumberOfParticles() ; j++)
544 {
bed069a4 545 /***************************************/
546 /****** Looping same events ********/
547 /****** filling numerators ********/
548 /***************************************/
81b7b887 549 if ( (j%fDisplayMixingInfo) == 0)
550 Info("ProcessTracks",
551 "Mixing particle %d from event %d with particles from event %d",j,i,i);
552
bed069a4 553 track1= trackEvent->GetParticle(j);
554 Bool_t firstcut = fPairCut->GetFirstPartCut()->Pass(track1);
555
9277d10a 556 if (fBufferSize != 0)
47d9a058 557 if ( (firstcut == kFALSE) || (fPairCut->GetSecondPartCut()->Pass(track1) == kFALSE) )
558 {
559 //accepted by any cut
560 // we have to copy because reader keeps only one event
561 trackEvent1->AddParticle(new AliHBTParticle(*track1));
562 }
bed069a4 563
564 if (firstcut) continue;
81b7b887 565
566 for(ii = 0; ii<fNTrackMonitorFunctions; ii++)
e4f2b1da 567 fTrackMonitorFunctions[ii]->Process(track1);
5c58441a 568
bed069a4 569 if ( fNTrackFunctions ==0 ) continue;
5c58441a 570
1b446896 571 for (Int_t k =j+1; k < trackEvent->GetNumberOfParticles() ; k++)
572 {
d74e6a27 573 track2= trackEvent->GetParticle(k);
574 if (track1->GetUID() == track2->GetUID()) continue;
575
5c58441a 576 trackpair->SetParticles(track1,track2);
577 if(fPairCut->Pass(trackpair)) //check pair cut
81b7b887 578 { //do not meets crietria of the
5c58441a 579 if( fPairCut->Pass(trackpair->GetSwapedPair()) ) continue;
580 else tmptrackpair = trackpair->GetSwapedPair();
81b7b887 581 }
582 else
583 {
584 tmptrackpair = trackpair;
585 }
586 for(ii = 0;ii<fNTrackFunctions;ii++)
587 fTrackFunctions[ii]->ProcessSameEventParticles(tmptrackpair);
588 }
bed069a4 589 /***************************************/
590 /***** Filling denominators *********/
591 /***************************************/
592
593 if (fBufferSize == 0) continue;
594
595 trackbuffer.ResetIter();
596
40d40d72 597 Int_t m = 0;
bed069a4 598 while (( trackEvent2 = trackbuffer.Next() ))
81b7b887 599 {
40d40d72 600 m++;
601 if ( (j%fDisplayMixingInfo) == 0)
602 Info("ProcessTracks",
603 "Mixing track %d from event %d with tracks from event %d",j,i,i-m);
81b7b887 604 for(Int_t l = 0; l<trackEvent2->GetNumberOfParticles();l++) // ... on all particles
bed069a4 605 {
606
607 track2= trackEvent2->GetParticle(l);
608 trackpair->SetParticles(track1,track2);
609
610 if( fPairCut->Pass(trackpair) ) //check pair cut
611 { //do not meets crietria of the
612 if( fPairCut->Pass(trackpair->GetSwapedPair()) )
613 continue;
614 else
615 {
616 tmptrackpair = trackpair->GetSwapedPair();
617 }
618 }
619 else
620 {//meets criteria of the pair cut
81b7b887 621 tmptrackpair = trackpair;
bed069a4 622 }
623
624 for(ii = 0;ii<fNTrackFunctions;ii++)
625 fTrackFunctions[ii]->ProcessDiffEventParticles(tmptrackpair);
626
81b7b887 627 }//for(Int_t l = 0; l<N2;l++) // ... on all particles
bed069a4 628 }
1b446896 629 }
bed069a4 630 trackEvent1 = trackbuffer.Push(trackEvent1);
631 }//while (fReader->Next() == kFALSE)
4cb6e8f7 632
633 delete trackpair;
634 trackbuffer.SetOwner(kTRUE);
1b446896 635}
636
637/*************************************************************************************/
491d1b5d 638
1b446896 639void AliHBTAnalysis::ProcessParticles()
640{
81b7b887 641//In order to minimize calling AliRun::GetEvent (we need at one time particles from different events),
1b446896 642//the loops are splited
643 AliHBTParticle * part1, * part2;
1b446896 644 AliHBTEvent * partEvent;
bed069a4 645 AliHBTEvent * partEvent1 = 0x0;
1b446896 646 AliHBTEvent * partEvent2;
491d1b5d 647
bed069a4 648 register UInt_t ii;
649
491d1b5d 650 AliHBTPair * partpair = new AliHBTPair();
bed069a4 651 AliHBTPair * tmppartpair; //temporary pointer
1b446896 652
bed069a4 653 AliHBTEventBuffer partbuffer(fBufferSize);
4cb6e8f7 654 partbuffer.SetOwner(kTRUE);
1b446896 655
bed069a4 656 fReader->Rewind();
657 Int_t i = -1;
658 while (fReader->Next() == kFALSE)
1b446896 659 {
bed069a4 660 i++;
661 partEvent = fReader->GetParticleEvent();
1b446896 662 if (!partEvent) continue;
bed069a4 663
664 if(partEvent1 == 0x0)
665 {
666 partEvent1 = new AliHBTEvent();
667 partEvent1->SetOwner(kTRUE);
668 }
669 else
670 {
671 partEvent1->Reset();
672 }
1b446896 673
674 for (Int_t j = 0; j<partEvent->GetNumberOfParticles() ; j++)
675 {
bed069a4 676 /***************************************/
677 /****** Looping same events ********/
678 /****** filling numerators ********/
679 /***************************************/
81b7b887 680 if ( (j%fDisplayMixingInfo) == 0)
40d40d72 681 Info("ProcessParticles",
bed069a4 682 "Mixing particle %d from event %d with particles from event %d",j,i,i);
683
40d40d72 684 part1 = partEvent->GetParticle(j);
bed069a4 685 Bool_t firstcut = fPairCut->GetFirstPartCut()->Pass(part1);
686
9277d10a 687 if (fBufferSize != 0) //useless in case
47d9a058 688 if ( (firstcut == kFALSE) || (fPairCut->GetSecondPartCut()->Pass(part1) == kFALSE) )
689 {
690 //accepted by any cut
691 // we have to copy because reader keeps only one event
692 partEvent1->AddParticle(new AliHBTParticle(*part1));
693 }
1b446896 694
bed069a4 695 if (firstcut) continue;
5c58441a 696
bed069a4 697 for(ii = 0; ii<fNParticleMonitorFunctions; ii++)
698 fParticleMonitorFunctions[ii]->Process(part1);
1b446896 699
bed069a4 700 if ( fNParticleFunctions == 0 ) continue;
701
702 for (Int_t k =j+1; k < partEvent->GetNumberOfParticles() ; k++)
703 {
704 part2= partEvent->GetParticle(k);
705 if (part1->GetUID() == part2->GetUID()) continue;
81b7b887 706
bed069a4 707 partpair->SetParticles(part1,part2);
708 if(fPairCut->Pass(partpair)) //check pair cut
709 { //do not meets crietria of the
710 if( fPairCut->Pass(partpair->GetSwapedPair()) ) continue;
711 else tmppartpair = partpair->GetSwapedPair();
7a2c8238 712 }
bed069a4 713 else
7a2c8238 714 {
bed069a4 715 tmppartpair = partpair;
7a2c8238 716 }
bed069a4 717 for(ii = 0;ii<fNParticleFunctions;ii++)
718 fParticleFunctions[ii]->ProcessSameEventParticles(tmppartpair);
719 }
720 /***************************************/
721 /***** Filling denominators *********/
722 /***************************************/
723
724 if (fBufferSize == 0) continue;
725
726 partbuffer.ResetIter();
727
40d40d72 728 Int_t m = 0;
bed069a4 729 while (( partEvent2 = partbuffer.Next() ))
730 {
40d40d72 731 m++;
732 if ( (j%fDisplayMixingInfo) == 0)
733 Info("ProcessParticles",
734 "Mixing particle %d from event %d with particles from event %d",j,i,i-m);
bed069a4 735 for(Int_t l = 0; l<partEvent2->GetNumberOfParticles();l++) // ... on all particles
1b446896 736 {
bed069a4 737
1b446896 738 part2= partEvent2->GetParticle(l);
739 partpair->SetParticles(part1,part2);
bed069a4 740
741 if( fPairCut->Pass(partpair) ) //check pair cut
1b446896 742 { //do not meets crietria of the
bed069a4 743 if( fPairCut->Pass(partpair->GetSwapedPair()) )
744 continue;
745 else
746 {
747 tmppartpair = partpair->GetSwapedPair();
748 }
491d1b5d 749 }
750 else
bed069a4 751 {//meets criteria of the pair cut
752 tmppartpair = partpair;
753 }
754
1b446896 755 for(ii = 0;ii<fNParticleFunctions;ii++)
491d1b5d 756 fParticleFunctions[ii]->ProcessDiffEventParticles(tmppartpair);
bed069a4 757
758 }//for(Int_t l = 0; l<N2;l++) // ... on all particles
759 }
1b446896 760 }
bed069a4 761 partEvent1 = partbuffer.Push(partEvent1);
762 }//while (fReader->Next() == kFALSE)
4cb6e8f7 763 delete partpair;
764 partbuffer.SetOwner(kTRUE);
1b446896 765}
1b446896 766/*************************************************************************************/
767
491d1b5d 768void AliHBTAnalysis::WriteFunctions()
1b446896 769{
81b7b887 770//Calls Write for all defined functions in analysis
771//== writes all results
1b446896 772 UInt_t ii;
773 for(ii = 0;ii<fNParticleFunctions;ii++)
774 fParticleFunctions[ii]->Write();
775
776 for(ii = 0;ii<fNTrackFunctions;ii++)
777 fTrackFunctions[ii]->Write();
778
779 for(ii = 0;ii<fNParticleAndTrackFunctions;ii++)
780 fParticleAndTrackFunctions[ii]->Write();
5c58441a 781
782 for(ii = 0;ii<fNParticleMonitorFunctions;ii++)
783 fParticleMonitorFunctions[ii]->Write();
784
785 for(ii = 0;ii<fNTrackMonitorFunctions;ii++)
786 fTrackMonitorFunctions[ii]->Write();
787
788 for(ii = 0;ii<fNParticleAndTrackMonitorFunctions;ii++)
789 fParticleAndTrackMonitorFunctions[ii]->Write();
1b446896 790}
791/*************************************************************************************/
491d1b5d 792
1b446896 793void AliHBTAnalysis::SetGlobalPairCut(AliHBTPairCut* cut)
794{
81b7b887 795//Sets the global cut
1b446896 796 if (cut == 0x0)
797 {
798 Error("AliHBTAnalysis::SetGlobalPairCut","Pointer is NULL. Ignoring");
799 }
800 delete fPairCut;
801 fPairCut = (AliHBTPairCut*)cut->Clone();
802}
803
804/*************************************************************************************/
491d1b5d 805
27b3fe5d 806void AliHBTAnalysis::AddTrackFunction(AliHBTOnePairFctn* f)
1b446896 807{
81b7b887 808//Adds track function
1b446896 809 if (f == 0x0) return;
810 if (fNTrackFunctions == fgkFctnArraySize)
811 {
812 Error("AliHBTAnalysis::AddTrackFunction","Can not add this function, not enough place in the array.");
813 }
814 fTrackFunctions[fNTrackFunctions] = f;
815 fNTrackFunctions++;
816}
491d1b5d 817/*************************************************************************************/
818
27b3fe5d 819void AliHBTAnalysis::AddParticleFunction(AliHBTOnePairFctn* f)
1b446896 820{
81b7b887 821//adds particle function
1b446896 822 if (f == 0x0) return;
823
824 if (fNParticleFunctions == fgkFctnArraySize)
825 {
826 Error("AliHBTAnalysis::AddParticleFunction","Can not add this function, not enough place in the array.");
827 }
828 fParticleFunctions[fNParticleFunctions] = f;
829 fNParticleFunctions++;
1b446896 830}
5c58441a 831/*************************************************************************************/
832
27b3fe5d 833void AliHBTAnalysis::AddParticleAndTrackFunction(AliHBTTwoPairFctn* f)
1b446896 834{
81b7b887 835//add resolution function
1b446896 836 if (f == 0x0) return;
837 if (fNParticleAndTrackFunctions == fgkFctnArraySize)
838 {
839 Error("AliHBTAnalysis::AddParticleAndTrackFunction","Can not add this function, not enough place in the array.");
840 }
841 fParticleAndTrackFunctions[fNParticleAndTrackFunctions] = f;
842 fNParticleAndTrackFunctions++;
843}
5c58441a 844/*************************************************************************************/
845
846void AliHBTAnalysis::AddParticleMonitorFunction(AliHBTMonOneParticleFctn* f)
847{
81b7b887 848//add particle monitoring function
5c58441a 849 if (f == 0x0) return;
850
851 if (fNParticleMonitorFunctions == fgkFctnArraySize)
852 {
853 Error("AliHBTAnalysis::AddParticleMonitorFunction","Can not add this function, not enough place in the array.");
854 }
855 fParticleMonitorFunctions[fNParticleMonitorFunctions] = f;
856 fNParticleMonitorFunctions++;
857}
858/*************************************************************************************/
1b446896 859
5c58441a 860void AliHBTAnalysis::AddTrackMonitorFunction(AliHBTMonOneParticleFctn* f)
861{
81b7b887 862//add track monitoring function
5c58441a 863 if (f == 0x0) return;
1b446896 864
5c58441a 865 if (fNTrackMonitorFunctions == fgkFctnArraySize)
866 {
867 Error("AliHBTAnalysis::AddTrackMonitorFunction","Can not add this function, not enough place in the array.");
868 }
869 fTrackMonitorFunctions[fNTrackMonitorFunctions] = f;
870 fNTrackMonitorFunctions++;
871}
1b446896 872/*************************************************************************************/
873
5c58441a 874void AliHBTAnalysis::AddParticleAndTrackMonitorFunction(AliHBTMonTwoParticleFctn* f)
875{
81b7b887 876//add resolution monitoring function
5c58441a 877 if (f == 0x0) return;
878 if (fNParticleAndTrackMonitorFunctions == fgkFctnArraySize)
879 {
880 Error("AliHBTAnalysis::AddParticleAndTrackMonitorFunction","Can not add this function, not enough place in the array.");
881 }
882 fParticleAndTrackMonitorFunctions[fNParticleAndTrackMonitorFunctions] = f;
883 fNParticleAndTrackMonitorFunctions++;
884}
885
1b446896 886
5c58441a 887/*************************************************************************************/
1b446896 888/*************************************************************************************/
491d1b5d 889
1b446896 890Bool_t AliHBTAnalysis::RunCoherencyCheck()
891{
892 //Checks if both HBTRuns are similar
893 //return true if error found
894 //if they seem to be OK return false
895 Int_t i;
81b7b887 896 Info("RunCoherencyCheck","Checking HBT Runs Coherency");
897
898 Info("RunCoherencyCheck","Number of events ...");
1b446896 899 if (fReader->GetNumberOfPartEvents() == fReader->GetNumberOfTrackEvents() ) //check whether there is the same number of events
900 {
81b7b887 901 Info("RunCoherencyCheck","OK. %d found\n",fReader->GetNumberOfTrackEvents());
1b446896 902 }
903 else
904 { //if not the same - ERROR
905 Error("AliHBTAnalysis::RunCoherencyCheck()",
906 "Number of simulated events (%d) is not equal to number of reconstructed events(%d)",
907 fReader->GetNumberOfPartEvents(),fReader->GetNumberOfTrackEvents());
908 return kTRUE;
909 }
910
81b7b887 911 Info("RunCoherencyCheck","Checking number of Particles AND Particles Types in each event ...");
1b446896 912
913 AliHBTEvent *partEvent;
914 AliHBTEvent *trackEvent;
915 for( i = 0; i<fReader->GetNumberOfTrackEvents();i++)
916 {
917 partEvent= fReader->GetParticleEvent(i); //gets the "ith" event
918 trackEvent = fReader->GetTrackEvent(i);
919
920 if ( (partEvent == 0x0) && (partEvent == 0x0) ) continue;
921 if ( (partEvent == 0x0) || (partEvent == 0x0) )
922 {
923 Error("AliHBTAnalysis::RunCoherencyCheck()",
924 "One event is NULL and the other one not. Event Number %d",i);
925 return kTRUE;
926 }
927
928 if ( partEvent->GetNumberOfParticles() != trackEvent->GetNumberOfParticles() )
929 {
930 Error("AliHBTAnalysis::RunCoherencyCheck()",
931 "Event %d: Number of simulated particles (%d) not equal to number of reconstructed tracks (%d)",
932 i,partEvent->GetNumberOfParticles() , trackEvent->GetNumberOfParticles());
933 return kTRUE;
934 }
935 else
936 for (Int_t j = 0; j<partEvent->GetNumberOfParticles(); j++)
937 {
938 if( partEvent->GetParticle(j)->GetPdgCode() != trackEvent->GetParticle(j)->GetPdgCode() )
939 {
940 Error("AliHBTAnalysis::RunCoherencyCheck()",
941 "Event %d: Particle %d: PID of simulated particle (%d) not the same of reconstructed track (%d)",
942 i,j, partEvent->GetParticle(j)->GetPdgCode(),trackEvent->GetParticle(j)->GetPdgCode() );
943 return kTRUE;
944
945 }
946 }
947 }
81b7b887 948 Info("RunCoherencyCheck"," Done");
949 Info("RunCoherencyCheck"," Everything looks OK");
1b446896 950 return kFALSE;
951}
952
dc2c3f36 953/*************************************************************************************/
954
955void AliHBTAnalysis::ProcessTracksAndParticlesNonIdentAnal()
956{
81b7b887 957//Performs analysis for both, tracks and particles
958
dc2c3f36 959 AliHBTParticle * part1, * part2;
960 AliHBTParticle * track1, * track2;
961
aeba88d7 962 AliHBTEvent * trackEvent1=0x0,*partEvent1=0x0;
963 AliHBTEvent * trackEvent2=0x0,*partEvent2=0x0;
964 AliHBTEvent * trackEvent3=0x0,*partEvent3=0x0;
dc2c3f36 965
bed069a4 966 AliHBTEvent * rawtrackEvent, *rawpartEvent;//this we get from Reader
dc2c3f36 967
dc2c3f36 968 AliHBTPair * trackpair = new AliHBTPair();
969 AliHBTPair * partpair = new AliHBTPair();
970
bed069a4 971 AliHBTEventBuffer partbuffer(fBufferSize);
972 AliHBTEventBuffer trackbuffer(fBufferSize);
973
974 register UInt_t ii;
dc2c3f36 975
976 trackEvent1 = new AliHBTEvent();
977 partEvent1 = new AliHBTEvent();
978 trackEvent1->SetOwner(kFALSE);
979 partEvent1->SetOwner(kFALSE);;
980
bed069a4 981 Bool_t nocorrfctns = (fNParticleFunctions == 0) && (fNTrackFunctions == 0) && (fNParticleAndTrackFunctions == 0);
982
983 fReader->Rewind();
984
81b7b887 985 Info("ProcessTracksAndParticlesNonIdentAnal","**************************************");
986 Info("ProcessTracksAndParticlesNonIdentAnal","***** NON IDENT MODE ****************");
987 Info("ProcessTracksAndParticlesNonIdentAnal","**************************************");
988
bed069a4 989 for (Int_t i = 0;;i++)//infinite loop
dc2c3f36 990 {
bed069a4 991 if (fReader->Next()) break; //end when no more events available
992
993 rawpartEvent = fReader->GetParticleEvent();
994 rawtrackEvent = fReader->GetTrackEvent();
dc2c3f36 995 if ( (rawpartEvent == 0x0) || (rawtrackEvent == 0x0) ) continue;//in case of any error
bed069a4 996
997 if ( rawpartEvent->GetNumberOfParticles() != rawtrackEvent->GetNumberOfParticles() )
998 {
999 Fatal("ProcessTracksAndParticlesNonIdentAnal",
1000 "Event %d: Number of simulated particles (%d) not equal to number of reconstructed tracks (%d)",
1001 i,rawpartEvent->GetNumberOfParticles() , rawtrackEvent->GetNumberOfParticles());
1002 return;
1003 }
1004
dc2c3f36 1005 /********************************/
1006 /* Filtering out */
1007 /********************************/
bed069a4 1008 if ( ( (partEvent2==0x0) || (trackEvent2==0x0)) )//in case fBufferSize == 0 and pointers are created do not eneter
1009 {
81b7b887 1010 partEvent2 = new AliHBTEvent();
1011 trackEvent2 = new AliHBTEvent();
bed069a4 1012 partEvent2->SetOwner(kTRUE);
1013 trackEvent2->SetOwner(kTRUE);
1014 }
1015
dc2c3f36 1016 FilterOut(partEvent1, partEvent2, rawpartEvent, trackEvent1, trackEvent2, rawtrackEvent);
1017
1018 for (Int_t j = 0; j<partEvent1->GetNumberOfParticles() ; j++)
1019 {
81b7b887 1020 if ( (j%fDisplayMixingInfo) == 0)
1021 Info("ProcessTracksAndParticlesNonIdentAnal",
1022 "Mixing particle %d from event %d with particles from event %d",j,i,i);
dc2c3f36 1023
1024 part1= partEvent1->GetParticle(j);
1025 track1= trackEvent1->GetParticle(j);
4d91c73a 1026
1027//PID reconstruction imperfections
1028// if( part1->GetPdgCode() != track1->GetPdgCode() )
1029// {
1030// Fatal("ProcessTracksAndParticlesNonIdentAnal",
1031// "Event %d: Particle %d: PID of simulated particle (%d) not the same of reconstructed track (%d)",
1032// i,j, part1->GetPdgCode(),track1->GetPdgCode() );
1033// return;
1034// }
bed069a4 1035
81b7b887 1036 for(ii = 0; ii<fNParticleMonitorFunctions; ii++)
e4f2b1da 1037 fParticleMonitorFunctions[ii]->Process(part1);
81b7b887 1038 for(ii = 0; ii<fNTrackMonitorFunctions; ii++)
e4f2b1da 1039 fTrackMonitorFunctions[ii]->Process(track1);
81b7b887 1040 for(ii = 0; ii<fNParticleAndTrackMonitorFunctions; ii++)
e4f2b1da 1041 fParticleAndTrackMonitorFunctions[ii]->Process(track1,part1);
dc2c3f36 1042
bed069a4 1043 if (nocorrfctns) continue;
1044
dc2c3f36 1045 /***************************************/
1046 /****** filling numerators ********/
1047 /****** (particles from event2) ********/
1048 /***************************************/
bed069a4 1049
1050 for (Int_t k = 0; k < partEvent2->GetNumberOfParticles() ; k++) //partEvent1 and partEvent2 are particles from the same event but separated to two groups
dc2c3f36 1051 {
1052 part2= partEvent2->GetParticle(k);
d74e6a27 1053 if (part1->GetUID() == part2->GetUID()) continue;//this is the same particle but with different PID
dc2c3f36 1054 partpair->SetParticles(part1,part2);
1055
1056 track2= trackEvent2->GetParticle(k);
1057 trackpair->SetParticles(track1,track2);
1058
17d74b37 1059 if( (this->*fkPassPairProp)(partpair,trackpair) ) //check pair cut
81b7b887 1060 { //do not meets crietria of the pair cut
1061 continue;
1062 }
dc2c3f36 1063 else
1064 {//meets criteria of the pair cut
dc2c3f36 1065 for(ii = 0;ii<fNParticleFunctions;ii++)
1066 fParticleFunctions[ii]->ProcessSameEventParticles(partpair);
1067
1068 for(ii = 0;ii<fNTrackFunctions;ii++)
1069 fTrackFunctions[ii]->ProcessSameEventParticles(trackpair);
1070
1071 for(ii = 0;ii<fNParticleAndTrackFunctions;ii++)
1072 fParticleAndTrackFunctions[ii]->ProcessSameEventParticles(trackpair,partpair);
1073 }
1074 }
5c58441a 1075
81b7b887 1076 if ( fBufferSize == 0) continue;//do not mix diff histograms
dc2c3f36 1077 /***************************************/
1078 /***** Filling denominators *********/
1079 /***************************************/
bed069a4 1080 partbuffer.ResetIter();
1081 trackbuffer.ResetIter();
1082
dc2c3f36 1083 Int_t nmonitor = 0;
1084
bed069a4 1085 while ( (partEvent3 = partbuffer.Next() ) != 0x0)
dc2c3f36 1086 {
bed069a4 1087 trackEvent3 = trackbuffer.Next();
1088
81b7b887 1089 if ( (j%fDisplayMixingInfo) == 0)
1090 Info("ProcessTracksAndParticlesNonIdentAnal",
1091 "Mixing particle %d from event %d with particles from event %d",j,i,i-(++nmonitor));
dc2c3f36 1092
1093 for (Int_t k = 0; k < partEvent3->GetNumberOfParticles() ; k++)
1094 {
1095 part2= partEvent3->GetParticle(k);
1096 partpair->SetParticles(part1,part2);
1b446896 1097
dc2c3f36 1098 track2= trackEvent3->GetParticle(k);
1099 trackpair->SetParticles(track1,track2);
1100
17d74b37 1101 if( (this->*fkPassPairProp)(partpair,trackpair) ) //check pair cut
81b7b887 1102 { //do not meets crietria of the pair cut
1103 continue;
1104 }
dc2c3f36 1105 else
1106 {//meets criteria of the pair cut
1107 UInt_t ii;
1108 for(ii = 0;ii<fNParticleFunctions;ii++)
1109 fParticleFunctions[ii]->ProcessDiffEventParticles(partpair);
1110
1111 for(ii = 0;ii<fNTrackFunctions;ii++)
1112 fTrackFunctions[ii]->ProcessDiffEventParticles(trackpair);
1113
1114 for(ii = 0;ii<fNParticleAndTrackFunctions;ii++)
1115 fParticleAndTrackFunctions[ii]->ProcessDiffEventParticles(trackpair,partpair);
1116 }
1117 }// for particles event2
1118 }//while event2
1119 }//for over particles in event1
1120
bed069a4 1121 partEvent2 = partbuffer.Push(partEvent2);
1122 trackEvent2 = trackbuffer.Push(trackEvent2);
dc2c3f36 1123
1124 }//end of loop over events (1)
bed069a4 1125
1126 partbuffer.SetOwner(kTRUE);
1127 trackbuffer.SetOwner(kTRUE);
1128
b70eb506 1129 delete partEvent1;
1130 delete trackEvent1;
bed069a4 1131 delete partEvent2;
1132 delete trackEvent2;
b70eb506 1133 delete partpair;
1134 delete trackpair;
81b7b887 1135
dc2c3f36 1136}
1b446896 1137/*************************************************************************************/
1138
dc2c3f36 1139void AliHBTAnalysis::ProcessTracksNonIdentAnal()
1140{
2dc7203b 1141//Process Tracks only with non identical mode
dc2c3f36 1142 AliHBTParticle * track1, * track2;
1143
aeba88d7 1144 AliHBTEvent * trackEvent1=0x0;
1145 AliHBTEvent * trackEvent2=0x0;
1146 AliHBTEvent * trackEvent3=0x0;
dc2c3f36 1147
1148 AliHBTEvent * rawtrackEvent;
dc2c3f36 1149
dc2c3f36 1150 AliHBTPair * trackpair = new AliHBTPair();
bed069a4 1151 AliHBTEventBuffer trackbuffer(fBufferSize);
dc2c3f36 1152
81b7b887 1153 register UInt_t ii;
dc2c3f36 1154
1155 trackEvent1 = new AliHBTEvent();
1156 trackEvent1->SetOwner(kFALSE);
1157
bed069a4 1158 fReader->Rewind();
1159
81b7b887 1160 Info("ProcessTracksNonIdentAnal","**************************************");
1161 Info("ProcessTracksNonIdentAnal","***** NON IDENT MODE ****************");
1162 Info("ProcessTracksNonIdentAnal","**************************************");
1163
bed069a4 1164
1165 for (Int_t i = 0;;i++)//infinite loop
dc2c3f36 1166 {
bed069a4 1167 if (fReader->Next()) break; //end when no more events available
1168 rawtrackEvent = fReader->GetTrackEvent();
1169
dc2c3f36 1170 if (rawtrackEvent == 0x0) continue;//in case of any error
1171
1172 /********************************/
1173 /* Filtering out */
1174 /********************************/
bed069a4 1175 if ( trackEvent2==0x0 )//in case fBufferSize == 0 and pointers are created do not eneter
1176 {
81b7b887 1177 trackEvent2 = new AliHBTEvent();
bed069a4 1178 trackEvent2->SetOwner(kTRUE);
1179 }
1180
dc2c3f36 1181 FilterOut(trackEvent1, trackEvent2, rawtrackEvent);
1182
1183 for (Int_t j = 0; j<trackEvent1->GetNumberOfParticles() ; j++)
1184 {
81b7b887 1185 if ( (j%fDisplayMixingInfo) == 0)
1186 Info("ProcessTracksNonIdentAnal",
1187 "Mixing particle %d from event %d with particles from event %d",j,i,i);
dc2c3f36 1188
1189 track1= trackEvent1->GetParticle(j);
1190
81b7b887 1191 for(ii = 0; ii<fNTrackMonitorFunctions; ii++)
e4f2b1da 1192 fTrackMonitorFunctions[ii]->Process(track1);
bed069a4 1193
1194 if (fNTrackFunctions == 0x0) continue;
1195
dc2c3f36 1196 /***************************************/
1197 /****** filling numerators ********/
1198 /****** (particles from event2) ********/
1199 /***************************************/
1200 for (Int_t k = 0; k < trackEvent2->GetNumberOfParticles() ; k++)
1201 {
1202 track2= trackEvent2->GetParticle(k);
d74e6a27 1203 if (track1->GetUID() == track2->GetUID()) continue;//this is the same particle but with different PID
dc2c3f36 1204 trackpair->SetParticles(track1,track2);
1205
1206
1207 if( fPairCut->PassPairProp(trackpair)) //check pair cut
1208 { //do not meets crietria of the pair cut
1209 continue;
1210 }
1211 else
1212 {//meets criteria of the pair cut
1213 UInt_t ii;
1214 for(ii = 0;ii<fNTrackFunctions;ii++)
1215 fTrackFunctions[ii]->ProcessSameEventParticles(trackpair);
1216 }
1217 }
81b7b887 1218 if ( fBufferSize == 0) continue;//do not mix diff histograms
dc2c3f36 1219 /***************************************/
1220 /***** Filling denominators *********/
1221 /***************************************/
bed069a4 1222 trackbuffer.ResetIter();
dc2c3f36 1223 Int_t nmonitor = 0;
1224
bed069a4 1225 while ( (trackEvent3 = trackbuffer.Next() ) != 0x0)
dc2c3f36 1226 {
1227
81b7b887 1228 if ( (j%fDisplayMixingInfo) == 0)
1229 Info("ProcessTracksNonIdentAnal",
1230 "Mixing particle %d from event %d with particles from event %d",j,i,i-(++nmonitor));
dc2c3f36 1231
1232 for (Int_t k = 0; k < trackEvent3->GetNumberOfParticles() ; k++)
1233 {
1234
1235 track2= trackEvent3->GetParticle(k);
1236 trackpair->SetParticles(track1,track2);
1237
1238
1239 if( fPairCut->PassPairProp(trackpair)) //check pair cut
81b7b887 1240 { //do not meets crietria of the pair cut
1241 continue;
1242 }
dc2c3f36 1243 else
1244 {//meets criteria of the pair cut
81b7b887 1245 for(ii = 0;ii<fNTrackFunctions;ii++)
1246 fTrackFunctions[ii]->ProcessDiffEventParticles(trackpair);
dc2c3f36 1247 }
1248 }// for particles event2
1249 }//while event2
1250 }//for over particles in event1
1251
bed069a4 1252 trackEvent2 = trackbuffer.Push(trackEvent2);
dc2c3f36 1253
1254 }//end of loop over events (1)
b70eb506 1255
bed069a4 1256 trackbuffer.SetOwner(kTRUE);
b70eb506 1257 delete trackpair;
1258 delete trackEvent1;
bed069a4 1259 delete trackEvent2;
dc2c3f36 1260}
1261/*************************************************************************************/
1262
1263void AliHBTAnalysis::ProcessParticlesNonIdentAnal()
1264{
2dc7203b 1265//process paricles only with non identical mode
dc2c3f36 1266 AliHBTParticle * part1 = 0x0, * part2 = 0x0;
1267
1268 AliHBTEvent * partEvent1 = 0x0;
1269 AliHBTEvent * partEvent2 = 0x0;
1270 AliHBTEvent * partEvent3 = 0x0;
1271
1272 AliHBTEvent * rawpartEvent = 0x0;
1273
dc2c3f36 1274 AliHBTPair * partpair = new AliHBTPair();
bed069a4 1275 AliHBTEventBuffer partbuffer(fBufferSize);
dc2c3f36 1276
bed069a4 1277 register UInt_t ii;
dc2c3f36 1278
1279 partEvent1 = new AliHBTEvent();
bed069a4 1280 partEvent1->SetOwner(kFALSE);
1281
1282 fReader->Rewind();
dc2c3f36 1283
81b7b887 1284 Info("ProcessParticlesNonIdentAnal","**************************************");
1285 Info("ProcessParticlesNonIdentAnal","***** NON IDENT MODE ****************");
1286 Info("ProcessParticlesNonIdentAnal","**************************************");
dc2c3f36 1287
bed069a4 1288 for (Int_t i = 0;;i++)//infinite loop
dc2c3f36 1289 {
bed069a4 1290 if (fReader->Next()) break; //end when no more events available
1291
1292 rawpartEvent = fReader->GetParticleEvent();
dc2c3f36 1293 if ( rawpartEvent == 0x0 ) continue;//in case of any error
1294
1295 /********************************/
1296 /* Filtering out */
1297 /********************************/
bed069a4 1298 if (partEvent2==0x0)//in case fBufferSize == 0 and pointers are created do not eneter
1299 {
b70eb506 1300 partEvent2 = new AliHBTEvent();
bed069a4 1301 partEvent2->SetOwner(kTRUE);
1302 }
1303
dc2c3f36 1304 FilterOut(partEvent1, partEvent2, rawpartEvent);
1305
1306 for (Int_t j = 0; j<partEvent1->GetNumberOfParticles() ; j++)
1307 {
81b7b887 1308 if ( (j%fDisplayMixingInfo) == 0)
1309 Info("ProcessParticlesNonIdentAnal",
1310 "Mixing particle %d from event %d with particles from event %d",j,i,i);
dc2c3f36 1311
1312 part1= partEvent1->GetParticle(j);
1313
bed069a4 1314 for(ii = 0; ii<fNParticleMonitorFunctions; ii++)
1315 fParticleMonitorFunctions[ii]->Process(part1);
1316
1317 if (fNParticleFunctions == 0) continue;
1318
dc2c3f36 1319 /***************************************/
1320 /****** filling numerators ********/
1321 /****** (particles from event2) ********/
1322 /***************************************/
1323 for (Int_t k = 0; k < partEvent2->GetNumberOfParticles() ; k++)
1324 {
1325 part2= partEvent2->GetParticle(k);
d74e6a27 1326 if (part1->GetUID() == part2->GetUID()) continue;//this is the same particle but with different PID
dc2c3f36 1327 partpair->SetParticles(part1,part2);
1328
1329 if(fPairCut->PassPairProp(partpair) ) //check pair cut
1330 { //do not meets crietria of the pair cut
1331 continue;
1332 }
1333 else
1334 {//meets criteria of the pair cut
dc2c3f36 1335 for(ii = 0;ii<fNParticleFunctions;ii++)
5c58441a 1336 fParticleFunctions[ii]->ProcessSameEventParticles(partpair);
dc2c3f36 1337 }
1338 }
b70eb506 1339 if ( fBufferSize == 0) continue;//do not mix diff histograms
dc2c3f36 1340 /***************************************/
1341 /***** Filling denominators *********/
1342 /***************************************/
bed069a4 1343 partbuffer.ResetIter();
dc2c3f36 1344 Int_t nmonitor = 0;
1345
bed069a4 1346 while ( (partEvent3 = partbuffer.Next() ) != 0x0)
dc2c3f36 1347 {
81b7b887 1348 if ( (j%fDisplayMixingInfo) == 0)
1349 Info("ProcessParticlesNonIdentAnal",
1350 "Mixing particle %d from event %d with particles from event %d",j,i,i-(++nmonitor));
1351
dc2c3f36 1352 for (Int_t k = 0; k < partEvent3->GetNumberOfParticles() ; k++)
1353 {
1354 part2= partEvent3->GetParticle(k);
1355 partpair->SetParticles(part1,part2);
1356
dc2c3f36 1357 if(fPairCut->PassPairProp(partpair) ) //check pair cut
1358 { //do not meets crietria of the pair cut
1359 continue;
1360 }
1361 else
1362 {//meets criteria of the pair cut
dc2c3f36 1363 for(ii = 0;ii<fNParticleFunctions;ii++)
5c58441a 1364 {
5c58441a 1365 fParticleFunctions[ii]->ProcessDiffEventParticles(partpair);
1366 }
dc2c3f36 1367 }
1368 }// for particles event2
1369 }//while event2
1370 }//for over particles in event1
bed069a4 1371 partEvent2 = partbuffer.Push(partEvent2);
dc2c3f36 1372 }//end of loop over events (1)
bed069a4 1373
1374 partbuffer.SetOwner(kTRUE);
b70eb506 1375 delete partpair;
1376 delete partEvent1;
bed069a4 1377 delete partEvent2;
dc2c3f36 1378}
1379
1380/*************************************************************************************/
1381void AliHBTAnalysis::FilterOut(AliHBTEvent* outpart1, AliHBTEvent* outpart2, AliHBTEvent* inpart,
17d74b37 1382 AliHBTEvent* outtrack1, AliHBTEvent* outtrack2, AliHBTEvent* intrack) const
dc2c3f36 1383{
1384 //Puts particles accepted as a first particle by global cut in out1
1385 //and as a second particle in out2
1386
1387 AliHBTParticle* part, *track;
1388
1389 outpart1->Reset();
1390 outpart2->Reset();
1391 outtrack1->Reset();
1392 outtrack2->Reset();
1393
dc2c3f36 1394 Bool_t in1, in2;
1395
1396 for (Int_t i = 0; i < inpart->GetNumberOfParticles(); i++)
1397 {
1398 in1 = in2 = kTRUE;
1399 part = inpart->GetParticle(i);
1400 track = intrack->GetParticle(i);
1401
17d74b37 1402 if ( ((this->*fkPass1)(part,track)) ) in1 = kFALSE; //if part is rejected by cut1, in1 is false
1403 if ( ((this->*fkPass2)(part,track)) ) in2 = kFALSE; //if part is rejected by cut2, in2 is false
dc2c3f36 1404
1405 if (gDebug)//to be removed in real analysis
1406 if ( in1 && in2 ) //both cuts accepted, should never happen, just in case
1407 {
1408 //Particle accpted by both cuts
1409 Error("FilterOut","Particle accepted by both cuts");
1410 continue;
1411 }
1412
1413 if (in1)
1414 {
1415 outpart1->AddParticle(part);
1416 outtrack1->AddParticle(track);
1417 continue;
1418 }
1419
1420 if (in2)
1421 {
bed069a4 1422 outpart2->AddParticle(new AliHBTParticle(*part));
1423 outtrack2->AddParticle(new AliHBTParticle(*track));
dc2c3f36 1424 continue;
1425 }
1426 }
dc2c3f36 1427}
1b446896 1428/*************************************************************************************/
81b7b887 1429
17d74b37 1430void AliHBTAnalysis::FilterOut(AliHBTEvent* out1, AliHBTEvent* out2, AliHBTEvent* in) const
dc2c3f36 1431{
1432 //Puts particles accepted as a first particle by global cut in out1
1433 //and as a second particle in out2
1434 AliHBTParticle* part;
1435
1436 out1->Reset();
1437 out2->Reset();
1438
1439 AliHBTParticleCut *cut1 = fPairCut->GetFirstPartCut();
1440 AliHBTParticleCut *cut2 = fPairCut->GetSecondPartCut();
1441
1442 Bool_t in1, in2;
1443
1444 for (Int_t i = 0; i < in->GetNumberOfParticles(); i++)
1445 {
1446 in1 = in2 = kTRUE;
1447 part = in->GetParticle(i);
1448
1449 if ( cut1->Pass(part) ) in1 = kFALSE; //if part is rejected by cut1, in1 is false
1450 if ( cut2->Pass(part) ) in2 = kFALSE; //if part is rejected by cut2, in2 is false
1451
1452 if (gDebug)//to be removed in real analysis
1453 if ( in1 && in2 ) //both cuts accepted, should never happen, just in case
1454 {
1455 //Particle accpted by both cuts
1456 Error("FilterOut","Particle accepted by both cuts");
1457 continue;
1458 }
1b446896 1459
dc2c3f36 1460 if (in1)
1461 {
1462 out1->AddParticle(part);
1463 continue;
1464 }
1465
1466 if (in2)
1467 {
1468 out2->AddParticle(part);
1469 continue;
1470 }
1471 }
1472}
1473/*************************************************************************************/
1b446896 1474
dc2c3f36 1475Bool_t AliHBTAnalysis::IsNonIdentAnalysis()
1476{
1477 //checks if it is possible to use special analysis for non identical particles
1478 //it means - in global pair cut first particle id is different than second one
1479 //and both are different from 0
1480 //in the future is possible to perform more sophisticated check
1481 //if cuts have excluding requirements
1482
5c58441a 1483 if (fPairCut->IsEmpty())
1484 return kFALSE;
1485
1486 if (fPairCut->GetFirstPartCut()->IsEmpty())
1487 return kFALSE;
1488
1489 if (fPairCut->GetSecondPartCut()->IsEmpty())
1490 return kFALSE;
dc2c3f36 1491
1492 Int_t id1 = fPairCut->GetFirstPartCut()->GetPID();
1493 Int_t id2 = fPairCut->GetSecondPartCut()->GetPID();
dc2c3f36 1494
5c58441a 1495 if ( (id1==0) || (id2==0) || (id1==id2) )
1496 return kFALSE;
1497
dc2c3f36 1498 return kTRUE;
1499}
66d1d1a4 1500/*************************************************************************************/
9616170a 1501
66d1d1a4 1502void AliHBTAnalysis::SetCutsOnParticles()
1503{
1504 // -- aplies only to Process("TracksAndParticles")
1505 // (ProcessTracksAndParticles and ProcessTracksAndParticlesNonIdentAnal)
1506 // Only particles properties are checkes against cuts
17d74b37 1507 fkPass = &AliHBTAnalysis::PassPart;
1508 fkPass1 = &AliHBTAnalysis::PassPart1;
1509 fkPass2 = &AliHBTAnalysis::PassPart2;
1510 fkPassPairProp = &AliHBTAnalysis::PassPairPropPart;
66d1d1a4 1511
1512}
1513/*************************************************************************************/
1514
1515void AliHBTAnalysis::SetCutsOnTracks()
1516{
1517 // -- aplies only to Process("TracksAndParticles")
1518 // (ProcessTracksAndParticles and ProcessTracksAndParticlesNonIdentAnal)
1519 // Only tracks properties are checkes against cuts
17d74b37 1520 fkPass = &AliHBTAnalysis::PassTrack;
1521 fkPass1 = &AliHBTAnalysis::PassTrack1;
1522 fkPass2 = &AliHBTAnalysis::PassTrack2;
1523 fkPassPairProp = &AliHBTAnalysis::PassPairPropTrack;
66d1d1a4 1524
1525}
1526/*************************************************************************************/
1527
1528void AliHBTAnalysis::SetCutsOnTracksAndParticles()
1529{
1530 // -- aplies only to Process("TracksAndParticles")
1531 // (ProcessTracksAndParticles and ProcessTracksAndParticlesNonIdentAnal)
1532 // Both, tracks and particles, properties are checked against cuts
17d74b37 1533 fkPass = &AliHBTAnalysis::PassPartAndTrack;
1534 fkPass1 = &AliHBTAnalysis::PassPartAndTrack1;
1535 fkPass2 = &AliHBTAnalysis::PassPartAndTrack2;
1536 fkPassPairProp = &AliHBTAnalysis::PassPairPropPartAndTrack;
66d1d1a4 1537}
9616170a 1538/*************************************************************************************/
1539
1540void AliHBTAnalysis::PressAnyKey()
17d74b37 1541{
1542 //small utility function that helps to make comfortable macros
9616170a 1543 char c;
1544 int nread = -1;
1545 fcntl(0, F_SETFL, O_NONBLOCK);
1546 ::Info("","Press Any Key to continue ...");
1547 while (nread<1)
1548 {
1549 nread = read(0, &c, 1);
1550 gSystem->ProcessEvents();
1551 }
1552}
1553
66d1d1a4 1554/*************************************************************************************/
9616170a 1555