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