Monitor functions implemented by Z.Chajecki.
[u/mrichter/AliRoot.git] / HBTAN / AliHBTAnalysis.cxx
1
2 #include "AliHBTAnalysis.h"
3
4 #include <iostream.h>
5
6 #include "AliHBTRun.h"
7 #include "AliHBTEvent.h"
8 #include "AliHBTReader.h"
9 #include "AliHBTParticle.h"
10 #include "AliHBTParticleCut.h"
11 #include "AliHBTPair.h"
12 #include "AliHBTPairCut.h"
13 #include "AliHBTFunction.h"
14 #include "AliHBTMonitorFunction.h"
15
16 #include <TBenchmark.h>
17 #include <TList.h>
18
19
20
21 ClassImp(AliHBTAnalysis)
22
23 const UInt_t AliHBTAnalysis::fgkFctnArraySize = 100;
24 const Int_t AliHBTAnalysis::fgkHbtAnalyzeAll = 0;
25
26 AliHBTAnalysis::AliHBTAnalysis()
27  {
28    fReader = 0x0;
29    
30    fTrackFunctions = new AliHBTOnePairFctn* [fgkFctnArraySize];
31    fParticleFunctions = new AliHBTOnePairFctn* [fgkFctnArraySize];
32    fParticleAndTrackFunctions = new AliHBTTwoPairFctn* [fgkFctnArraySize];
33    
34    fParticleMonitorFunctions = new AliHBTMonOneParticleFctn* [fgkFctnArraySize];    
35    fTrackMonitorFunctions = new AliHBTMonOneParticleFctn* [fgkFctnArraySize];    
36    fParticleAndTrackMonitorFunctions = new AliHBTMonTwoParticleFctn* [fgkFctnArraySize];    
37
38    fNTrackFunctions = 0;
39    fNParticleFunctions = 0;
40    fNParticleAndTrackFunctions = 0;
41   
42    fNParticleMonitorFunctions = 0; 
43    fNTrackMonitorFunctions = 0; 
44    fNParticleAndTrackMonitorFunctions = 0; 
45
46    fPairCut = new AliHBTEmptyPairCut();//empty cut - accepts all particles
47    
48    fBufferSize = 2; 
49  }
50 /*************************************************************************************/ 
51
52 AliHBTAnalysis::~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    
62    delete [] fParticleMonitorFunctions; 
63    delete [] fTrackMonitorFunctions; 
64    delete [] fParticleAndTrackMonitorFunctions; 
65
66    delete fPairCut; // always have an copy of an object - we create we dstroy
67  }
68
69 /*************************************************************************************/ 
70
71 void 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   {
99    Error("Process","The reader is not set");
100    return;
101   }
102  
103  
104  const char *oT = strstr(option,"Tracks");
105  const char *oP = strstr(option,"Particles");
106  
107  Bool_t nonid = IsNonIdentAnalysis();
108  // nonid=0;
109  cout<<"nonid = "<<nonid<<endl; 
110  
111  if(oT && oP)
112   { 
113     if (fReader->GetNumberOfPartEvents() <1)
114      {
115        Error("Process","There is no Particles. Maybe change the option?");
116        return;
117      }
118     if (fReader->GetNumberOfTrackEvents() <1)
119      {
120        Error("Process","There is no Tracks. Maybe change the option?");
121        return;
122      }
123     
124     if ( RunCoherencyCheck() )
125       {
126         Error("Process",
127               "Coherency check not passed. Maybe change the option?\n");
128         return;
129       }
130     if (nonid) ProcessTracksAndParticlesNonIdentAnal();
131     else ProcessTracksAndParticles();
132     return;
133   }
134
135  if(oT)
136   {
137     if (fReader->GetNumberOfTrackEvents() <1)
138      {
139        Error("Process","There is no data to analyze.");
140        return;
141      }
142     if (nonid) ProcessTracksNonIdentAnal();
143     else ProcessTracks();
144     return;
145   }
146  
147  if(oP)
148   {
149     if (fReader->GetNumberOfPartEvents() <1)
150      {
151        Error("Process","There is no data to analyze.");
152        return;
153      }
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     
162     return;
163   }
164  
165 }
166
167 /*************************************************************************************/ 
168
169 void 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();
192
193   AliHBTPair * tmptrackpair;//temprary pointers to pairs
194   AliHBTPair * tmppartpair;
195   
196   
197   
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        {
209          if ( (j%100) == 0) cout<<"Mixing particle "<<j<<" from event "<<i<<" with particles from event "<<i<<endl;
210
211          part1= partEvent->GetParticle(j);
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         
225          for (Int_t k =j+1; k < partEvent->GetNumberOfParticles() ; k++)
226           {
227             part2= partEvent->GetParticle(k);
228             partpair->SetParticles(part1,part2);
229            
230             track2= trackEvent->GetParticle(k);       
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 
238                { //swaped pair meets all the criteria
239                    tmppartpair = partpair->GetSwapedPair();
240                    tmptrackpair = trackpair->GetSwapedPair();
241                  }
242               }
243             else
244              {//meets criteria of the pair cut
245                tmptrackpair = trackpair;
246                tmppartpair = partpair;
247              }
248             UInt_t ii;
249             for(ii = 0;ii<fNParticleFunctions;ii++)
250                    fParticleFunctions[ii]->ProcessSameEventParticles(tmppartpair);
251                 
252             for(ii = 0;ii<fNTrackFunctions;ii++)
253                    fTrackFunctions[ii]->ProcessSameEventParticles(tmptrackpair);
254                  
255             for(ii = 0;ii<fNParticleAndTrackFunctions;ii++)
256                    fParticleAndTrackFunctions[ii]->ProcessSameEventParticles(tmptrackpair,tmppartpair);
257            }
258        }
259     }
260
261
262   /***** Filling denominators    *********/
263   /***************************************/
264   for (Int_t i = 0;i<Nev-1;i++)   //In each event (but last) ....
265     {
266   
267       if ((fNParticleFunctions == 0) && (fNTrackFunctions ==0) && (fNParticleAndTrackFunctions == 0))
268         continue;  
269
270       partEvent= fReader->GetParticleEvent(i);
271       if (!partEvent) continue;
272       
273       trackEvent = fReader->GetTrackEvent(i); 
274       
275 //      N=0;
276
277       for (Int_t j = 0; j< partEvent->GetNumberOfParticles(); j++) // ... Loop over all particles ...
278        {
279            
280            part1= partEvent->GetParticle(j);
281
282            track1= trackEvent->GetParticle(j);
283  
284            Int_t NNN;
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             }
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              
305              if ( (j%100) == 0) cout<<"Mixing particle "<<j<<" from event "<<i<<" with particles from event "<<k<<endl;
306             
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()) ) )
321           continue;
322        else 
323         {
324           tmppartpair = partpair->GetSwapedPair();
325           tmptrackpair = trackpair->GetSwapedPair();
326         }
327                   }
328                 else
329                  {//meets criteria of the pair cut
330                   tmptrackpair = trackpair;
331                   tmppartpair = partpair;
332                  }
333                 UInt_t ii;
334                 for(ii = 0;ii<fNParticleFunctions;ii++)
335                   fParticleFunctions[ii]->ProcessDiffEventParticles(tmppartpair);
336                  
337                 for(ii = 0;ii<fNTrackFunctions;ii++)
338                   fTrackFunctions[ii]->ProcessDiffEventParticles(tmptrackpair);
339                  
340                 for(ii = 0;ii<fNParticleAndTrackFunctions;ii++)
341                   fParticleAndTrackFunctions[ii]->ProcessDiffEventParticles(tmptrackpair,tmppartpair);
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
356 void 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();
374   AliHBTPair * tmptrackpair; //temporary pointer 
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        {
384          if ( (j%100) == 0) cout<<"Mixing particle "<<j<<" from event "<<i<<" with particles from event "<<i<<endl;
385
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         
395          for (Int_t k =j+1; k < trackEvent->GetNumberOfParticles() ; k++)
396           {
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();
403               }
404             else
405               {
406                 tmptrackpair = trackpair;
407               }
408             UInt_t ii;
409                 
410             for(ii = 0;ii<fNTrackFunctions;ii++)
411                    fTrackFunctions[ii]->ProcessSameEventParticles(tmptrackpair);
412                  
413            
414            }
415        }
416     }
417
418   /***************************************/
419   /***** Filling diff histogram *********/
420   /***************************************/
421   for (Int_t i = 0;i<Nev-1;i++)   //In each event (but last) ....
422     {
423       if ( fNTrackFunctions ==0 )
424         continue; 
425
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  
436            Int_t NNN;
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             }
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              
455              if ( (j%100) == 0) cout<<"Mixing particle "<<j<<" from event "<<i<<" with particles from event "<<k<<endl;
456             
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;
467        else tmptrackpair = trackpair->GetSwapedPair();
468                   }
469                 else
470                   {
471                     tmptrackpair = trackpair;
472                   }
473                 UInt_t ii;
474                 for(ii = 0;ii<fNTrackFunctions;ii++)
475                   fTrackFunctions[ii]->ProcessDiffEventParticles(tmptrackpair);
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 /*************************************************************************************/
489
490 void 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;
498
499   AliHBTPair * partpair = new AliHBTPair();
500   AliHBTPair * tmppartpair; //temporary pointer to the pair
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   
506  // Nev = 1;
507   /***************************************/
508   /******   Looping same events   ********/
509   /******   filling numerators    ********/
510   /***************************************/
511 //  gBenchmark->Start("id");
512   
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        {
521          if ( (j%100) == 0) cout<<"Mixing particle "<<j<<" from event "<<i<<" with particles from event "<<i<<endl;
522
523          part1= partEvent->GetParticle(j);
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
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
543                    tmppartpair = partpair->GetSwapedPair();
544                  }
545               }
546             else
547               {
548                 tmppartpair = partpair;
549               }
550             
551             UInt_t ii;
552             for(ii = 0;ii<fNParticleFunctions;ii++)
553                    fParticleFunctions[ii]->ProcessSameEventParticles(tmppartpair);
554            }
555        }
556     }
557
558   /***************************************/
559   /***** Filling diff histogram *********/
560   /***************************************/
561   for (Int_t i = 0;i<Nev-1;i++)   //In each event (but last)....
562     {
563       if ( fNParticleFunctions ==0 )
564         continue; 
565
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  
577            Int_t NNN;
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;
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              
597              if ( (j%100) == 0) cout<<"Mixing particle "<<j<<" from event "<<i<<" with particles from event "<<k<<endl;
598             
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;
609        else tmppartpair = partpair->GetSwapedPair();
610                   }
611                 else
612                   {
613                     tmppartpair = partpair;
614                   }
615                 UInt_t ii;
616                 for(ii = 0;ii<fNParticleFunctions;ii++)
617                   fParticleFunctions[ii]->ProcessDiffEventParticles(tmppartpair);
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   
627 //  gBenchmark->Stop("id");
628 //  gBenchmark->Show("id");
629
630 }
631
632 /*************************************************************************************/
633
634
635 void AliHBTAnalysis::WriteFunctions()
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();
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();
655 }
656 /*************************************************************************************/
657
658 void 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 /*************************************************************************************/
669
670 void AliHBTAnalysis::AddTrackFunction(AliHBTOnePairFctn* f)
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 }
680 /*************************************************************************************/ 
681
682 void AliHBTAnalysis::AddParticleFunction(AliHBTOnePairFctn* f)
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++;
692 }
693 /*************************************************************************************/ 
694
695 void AliHBTAnalysis::AddParticleAndTrackFunction(AliHBTTwoPairFctn* f)
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 }
705 /*************************************************************************************/ 
706
707 void 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 /*************************************************************************************/ 
719
720 void AliHBTAnalysis::AddTrackMonitorFunction(AliHBTMonOneParticleFctn* f)
721 {
722   if (f == 0x0) return;
723
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 }
731 /*************************************************************************************/ 
732
733 void 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
744
745 /*************************************************************************************/ 
746 /*************************************************************************************/  
747
748 Bool_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
812 /*************************************************************************************/  
813  
814 void 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       /********************************/
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
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
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);
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            }
913         
914
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);
931
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 }
965 /*************************************************************************************/  
966  
967 void 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       /********************************/
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
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
1018         UInt_t zz;
1019          for(zz = 0; zz<fNTrackMonitorFunctions; zz++)
1020            fTrackMonitorFunctions[zz]->ProcessSameEventParticles(track1);
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
1085 void 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       /********************************/
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
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
1136          UInt_t zz;
1137          for(zz = 0; zz<fNParticleMonitorFunctions; zz++)
1138            fParticleMonitorFunctions[zz]->ProcessSameEventParticles(part1);
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++)
1157                 {
1158                   //  if ((k%100) == 0) cout<<".";
1159                   fParticleFunctions[ii]->ProcessSameEventParticles(partpair);
1160                 }
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++)
1187                {
1188                 //  if ((k%100) == 0) cout<<"*";
1189                  fParticleFunctions[ii]->ProcessDiffEventParticles(partpair);
1190                }
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 /*************************************************************************************/  
1207 void 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
1259
1260 /*************************************************************************************/  
1261 void 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       }
1290
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 /*************************************************************************************/ 
1305
1306 Bool_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  
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;
1322  
1323  Int_t id1 = fPairCut->GetFirstPartCut()->GetPID();
1324  Int_t id2 = fPairCut->GetSecondPartCut()->GetPID();
1325
1326  if ( (id1==0) || (id2==0) || (id1==id2) ) 
1327    return kFALSE;
1328  
1329  return kTRUE;
1330 }