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