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