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