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