]> git.uio.no Git - u/mrichter/AliRoot.git/blobdiff - PWG2/RESONANCES/AliRsnPair.cxx
bugfix
[u/mrichter/AliRoot.git] / PWG2 / RESONANCES / AliRsnPair.cxx
index 063e8d54e43e850bfbce2a866e186dc86eb9424e..c8a0a8eeb0d86a9b2c0283432550d843eccf1bf5 100644 (file)
-/**************************************************************************
- * Copyright(c) 1998-1999, ALICE Experiment at CERN, All rights reserved. *
- *                                                                        *
- * Author: The ALICE Off-line Project.                                    *
- * Contributors are mentioned in the code where appropriate.              *
- *                                                                        *
- * Permission to use, copy, modify and distribute this software and its   *
- * documentation strictly for non-commercial purposes is hereby granted   *
- * without fee, provided that the above copyright notice appears in all   *
- * copies and that both the copyright notice and this permission notice   *
- * appear in the supporting documentation. The authors make no claims     *
- * about the suitability of this software for any purpose. It is          *
- * provided "as is" without express or implied warranty.                  *
- **************************************************************************/
-//-------------------------------------------------------------------------
-//                     Class AliRsnPair
-//-------------------------------------------------------------------------
-// This class computes the invariant mass spectrum of a specified pair of
-// particles, throughout a list of AliRsnEvents, and returns it as a TH1D.
-// This object is not supposed to be used directly: an AliRsnAnalysis
-// should be initialized in a macro and filled with one or more AliRsnPair's
-// which are then processed with a given sample of events.
-//   
-// author: A. Pulvirenti
-// email : alberto.pulvirenti@ct.infn.it
-//-------------------------------------------------------------------------
+//
+// *** Class AliRsnPair ***
+//
+// "Core" method for defining the work on a pari of particles.
+// For one analysis, one must setup one of this for each pair he wants to analyze,
+// adding to it all analysis which he desires to do.
+// Here he defines the cuts, and the particle types and charges, and can add
+// functions which do different operations on the same pair, and some binning
+// with respect to some kinematic variables (eta, momentum)
+//
+// authors: A. Pulvirenti (email: alberto.pulvirenti@ct.infn.it)
+//          M. Vala (email: martin.vala@cern.ch)
+//
 
 #include <Riostream.h>
-
-#include <TH1.h>
-#include <TString.h>
-#include <TRefArray.h>
-#include <TClonesArray.h>
+#include <TList.h>
 
 #include "AliLog.h"
-#include "AliRsnParticle.h"
-#include "AliRsnDaughter.h"
-#include "AliRsnDaughterCut.h"
-#include "AliRsnDaughterCutPair.h"
+
+#include "AliRsnMother.h"
 #include "AliRsnEvent.h"
+#include "AliRsnFunction.h"
+#include "AliRsnCutSet.h"
+#include "AliRsnCutStd.h"
+#include "AliRsnValue.h"
+#include "AliRsnCutManager.h"
 
 #include "AliRsnPair.h"
 
 ClassImp(AliRsnPair)
-//--------------------------------------------------------------------------------------------------------
-AliRsnPair::AliRsnPair() :
-  TNamed(),
-  fForMixing(kFALSE),
-  fStoreOnlyTrue(kFALSE),
-  fTrueMotherPDG(0),
-  fPtMin(0.0),
-  fPtMax(0.0),
-  fVtMax(0.0),
-  fCutsPair(0x0),
-  fHistogram(0x0)
-{
-//
-// Empty constructor.
-// Initializes the data members to default values:
-//  - default (empty string) name and title for the object
-//  - switch off the 'fStoreOnlyTrue' flag;
-//  - assume working with real data (no PDG code of the mother);
-//  - no cuts of any kind;
-//  - no definition of particles in the pair;
-//  - histogram undefined.
-// When using this constructor, all analysis elements (particles, histogram)
-// must be defined before starting event processing.
-//
 
-       Int_t i;
-       for (i = 0; i < 2; i++) {
-               fMass[i] = 0.0;
-               fCharge[i] = '0';
-               fType[i] = AliRsnPID::kUnknown;
-               fCutsSingle[i] = 0x0;
-       }
-}
-//--------------------------------------------------------------------------------------------------------
-AliRsnPair::AliRsnPair
-(const char *name, const char *title, Int_t nbins, Double_t min, Double_t max, 
- Double_t ptmin, Double_t ptmax, Double_t dmax) :
-  TNamed(name, title),
-  fForMixing(kFALSE),
-  fStoreOnlyTrue(kFALSE),
-  fTrueMotherPDG(0),
-  fPtMin(ptmin),
-  fPtMax(ptmax),
-  fVtMax(dmax),
-  fCutsPair(0x0),
-  fHistogram(0x0)
+//_____________________________________________________________________________
+AliRsnPair::AliRsnPair(const char *name, AliRsnPairDef *def) :
+  TNamed(name, ""),
+  fOnlyTrue(kFALSE),
+  fIsMixed(kFALSE),
+  fPairDef(def),
+  fCutManager(),
+  fMother(),
+  fEvent(0x0)
 {
 //
-// Constructor with arguments.
-// This constructor allows to define some of the initialization values:
-//  - name and title of the object
-//  - histogram binning and edges
-// The other parameters are initialized as in the default constructor.
+// Default constructor
 //
 
-       Int_t i;
-       for (i = 0; i < 2; i++) {
-               fMass[i] = 0.0;
-               fCharge[i] = '0';
-               fType[i] = AliRsnPID::kUnknown;
-               fCutsSingle[i] = 0x0;
-       }
-       fHistogram = new TH1D(Form("histo_%s", name), "AliRsnPair::fHistogram", nbins, min, max);
+  AliDebug(AliLog::kDebug+2,"<-");
+  AliDebug(AliLog::kDebug+2,"->");
 }
-//--------------------------------------------------------------------------------------------------------
-AliRsnPair::AliRsnPair(const AliRsnPair &copy) :
+
+//_____________________________________________________________________________
+AliRsnPair::AliRsnPair(const AliRsnPair& copy) :
   TNamed(copy),
-  fForMixing(copy.fForMixing),
-  fStoreOnlyTrue(copy.fStoreOnlyTrue),
-  fTrueMotherPDG(copy.fTrueMotherPDG),
-  fPtMin(copy.fPtMin),
-  fPtMax(copy.fPtMax),
-  fVtMax(copy.fVtMax),
-  fCutsPair(0x0),
-  fHistogram(0x0)
+  fOnlyTrue(copy.fOnlyTrue),
+  fIsMixed(copy.fIsMixed),
+  fPairDef(copy.fPairDef),
+  fCutManager(copy.fCutManager),
+  fMother(copy.fMother),
+  fEvent(0x0)
 {
 //
-// Copy constructor.
-// Default behavior as a copy constructor for what concerns non-array data-members.
-// The arrays are cloned if they are not NULL.
+// Default constructor
 //
 
-       if (copy.fHistogram) fHistogram = (TH1D*)(copy.fHistogram->Clone());
-       if (copy.fCutsPair) fCutsPair = (TObjArray*)(copy.fCutsPair->Clone());
-       
-       Int_t i;
-       for (i = 0; i < 2; i++) {
-               fMass[i] = copy.fMass[i];
-               fCharge[i] = copy.fCharge[i];
-               fType[i] = copy.fType[i];
-               if (copy.fCutsSingle[i]) fCutsSingle[i] = (TObjArray*)copy.fCutsSingle[i]->Clone();
-       }
+  AliDebug(AliLog::kDebug+2,"<-");
+  AliDebug(AliLog::kDebug+2,"->");
 }
-//--------------------------------------------------------------------------------------------------------
-const AliRsnPair& AliRsnPair::operator=(const AliRsnPair &copy)
-{
-//
-// Assignment operator.
-// Default behavior like copy constructor.
-//
-
-       fHistogram = 0x0;
-       fCutsPair = 0x0;
-       fCutsSingle[0] = fCutsSingle[1] = 0x0;
 
-    fForMixing = copy.fForMixing;
-       fStoreOnlyTrue = copy.fStoreOnlyTrue;
-       fTrueMotherPDG = copy.fTrueMotherPDG;
-       if (copy.fHistogram) fHistogram = (TH1D*)(copy.fHistogram->Clone());
-       if (copy.fCutsPair) fCutsPair = (TObjArray*)(copy.fCutsPair->Clone());
-    
-    fPtMin = copy.fPtMin;
-    fPtMax = copy.fPtMax;
-    fVtMax = copy.fVtMax;
-       
-       Int_t i;
-       for (i = 0; i < 2; i++) {
-               fMass[i] = copy.fMass[i];
-               fCharge[i] = copy.fCharge[i];
-               fType[i] = copy.fType[i];
-               if (copy.fCutsSingle[i]) fCutsSingle[i] = (TObjArray*)copy.fCutsSingle[i]->Clone();
-       }
-       
-       return (*this);
-}
-//--------------------------------------------------------------------------------------------------------
-void AliRsnPair::Clear(Option_t* /*option*/)
+//_____________________________________________________________________________
+AliRsnPair& AliRsnPair::operator=(const AliRsnPair& copy)
 {
-//
-// Clear arrays and histogram.
-// For the sake of security, all pointers are also set explicitly to NULL.
-//
-       fCutsSingle[0]->Delete();
-       fCutsSingle[1]->Delete();
-       fCutsPair->Delete();
-       
-       delete fCutsSingle[0];
-       delete fCutsSingle[1];
-       delete fCutsPair;
-       delete fHistogram;
-       
-       fCutsSingle[0] = fCutsSingle[1] = fCutsPair = 0;
-       fHistogram = 0;
+  fOnlyTrue = copy.fOnlyTrue;
+  fIsMixed = copy.fIsMixed;
+  fPairDef = copy.fPairDef;
+  fMother = copy.fMother;
+  fCutManager = copy.fCutManager;
+  fEvent = 0x0;
+
+  return (*this);
 }
-//--------------------------------------------------------------------------------------------------------
-void AliRsnPair::SetPair
-(Char_t charge1, AliRsnPID::EType type1, Char_t charge2, AliRsnPID::EType type2)
-{
-//
-// This method allows to set at once all the parameters of the particles in the pair.
-// The mass must not be specified, and it is retrieved from TDatabasePDG,
-// using a static method defined in AliRsnDaughter class.
-//
-
-       fCharge[0] = charge1;
-       fType[0] = type1;
-       SetMass(0, AliRsnPID::ParticleMass(type1));
 
-       fCharge[1] = charge2;
-       fType[1] = type2;
-       SetMass(1, AliRsnPID::ParticleMass(type2));
-}
-//--------------------------------------------------------------------------------------------------------
-void AliRsnPair::AddCutPair(AliRsnDaughterCutPair *cut)
+//_____________________________________________________________________________
+AliRsnPair::~AliRsnPair()
 {
 //
-// Add a pair cut.
-// If the cut array is NULL, it is initialized here.
+// Destructor
 //
 
-       if (!fCutsPair) fCutsPair = new TObjArray(0);
-       fCutsPair->AddLast(cut);
+  AliDebug(AliLog::kDebug+2,"<-");
+  AliDebug(AliLog::kDebug+2,"->");
 }
-//--------------------------------------------------------------------------------------------------------
-void AliRsnPair::AddCutSingle(Int_t i, AliRsnDaughterCut *cut)
+
+//_____________________________________________________________________________
+void AliRsnPair::Print(Option_t* /*option*/) const
 {
 //
-// Add a single particle cut.
-// If the cut array is NULL, it is initialized here.
+// Prints info about pair
 //
 
-       if (i < 0 || i > 1) return;
-       if (!fCutsSingle[i]) fCutsSingle[i] = new TObjArray(0);
-       fCutsSingle[i]->AddLast(cut);
+  AliDebug(AliLog::kDebug+2,"<-");
+  AliInfo(Form("PDG %d %d", AliPID::ParticleCode(fPairDef->GetPID(0)), AliPID::ParticleCode(fPairDef->GetPID(1))));
+  AliInfo(Form("Masses %f %f", fPairDef->GetMass(0), fPairDef->GetMass(1)));
+
+  AliDebug(AliLog::kDebug+2,"->");
 }
-//--------------------------------------------------------------------------------------------------------
-Stat_t AliRsnPair::Process(AliRsnEvent *event1, AliRsnEvent *event2, Bool_t usePID)
+
+//_____________________________________________________________________________
+Bool_t AliRsnPair::Fill
+(AliRsnDaughter *daughter0, AliRsnDaughter *daughter1, AliRsnEvent *ev0, AliRsnEvent *ev1)
 {
 //
-// Scans the two events specified in argument to fill the histogram.
-// This method essentially calls the AliRsnPair::Fill() method one or many times.
-// When the "noPID" argument is kFALSE, the analysis is done with identified particles
-// and this causes the Fill() method to be called only once, for the two lists of
-// identified particles of the two kinds specified in AliRsnPair datamembers.
-// When the "noPID" argument is kTRUE, the analysis is done with all collections
-// of particles of the same sign as specified in the two arguments of the pair.
-// ---
-// Particles of type #1 are taken in 'event1', and particles of type #2 are taken in 'event2'.
-// When doing single-event analysis (for resonance signal or like-sign background),
-// the second argument can be simply skipped.
-// When doing event mixing, the two arguments must be not null and different.
-// If argument #1 is NULL, an error is raised, while if argument #2 is NULL, no error is raised,
-// and 'event2' argument is set equal to 'event1' (= single event processing).
-// ---
-// Return value is the total number of pairs processed.
-//
-
-       // preliminary checks
-       if (!event1) {
-               // argument #1 cannot be NULL
-               AliError("Argument #1 cannot be NULL.");
-               return 0.0;
-       }
-       if (!event2) {
-               // if argument #2 is NULL, it is put equal to argument #1
-               event2 = event1;
-       }
+// Sets the two passed daughters to the AliRsnMother data member of this object
+// which is used to perform all computations to fill the value list.
+// This operation is done successfully only when the first passed object matches
+// the required object type (track/V0) and the required charge for first element in pair def,
+// and the second passed object does the same w.r. to the second element in pair def.
+// Moreover, all cuts are checked and the operation fails if a cut check is unsuccessful.
+// Finally, if a true pair is required, this is checked at the end.
+//
+
+  AliDebug(AliLog::kDebug+2,"<-");
+  
+  // check for correct type-charge match for first element
+  if (daughter0->RefType() != fPairDef->GetDaughterType(0)) return kFALSE;
+  if (daughter0->ChargeChar() != fPairDef->GetCharge(0)) return kFALSE;
+  
+  // check for correct type-charge match for second element
+  if (daughter1->RefType() != fPairDef->GetDaughterType(1)) return kFALSE;
+  if (daughter1->ChargeChar() != fPairDef->GetCharge(1)) return kFALSE;
     
-    // define if same indexes must be summed or not depending if
-    // the two events are the same or not
-    Bool_t skipSameIndex = (event1 == event2);
+  // cuts on track #1 & common
+  fCutManager.SetEvent(ev0);
+  if (!fCutManager.PassDaughter1Cuts(daughter0)) 
+  {
+    AliDebug(AliLog::kDebug+2, "Specific cuts for track #1 not passed");
+    return kFALSE;
+  }
+  if (!fCutManager.PassCommonDaughterCuts(daughter0))
+  {
+    AliDebug(AliLog::kDebug+2, "Common cuts for track #1 not passed");
+    return kFALSE;
+  }
+  
+  // cuts on track #2 & common
+  fCutManager.SetEvent(ev1);
+  if (!fCutManager.PassDaughter2Cuts(daughter1))
+  {
+    AliDebug(AliLog::kDebug+2, "Specific cuts for track #2 not passed");
+    return kFALSE;
+  }
+  if (!fCutManager.PassCommonDaughterCuts(daughter1))
+  {
+    AliDebug(AliLog::kDebug+2, "Common cuts for track #2 not passed");
+    return kFALSE;
+  }
+  
+  // point to first event as reference
+  fEvent = ev0;
+  
+  // define pair & check
+  fMother.SetDaughters(daughter0, fPairDef->GetMass(0), daughter1, fPairDef->GetMass(1));
+  fCutManager.SetEvent(fEvent);
+  if (!fCutManager.PassMotherCuts(&fMother)) return kFALSE;
+  
+  // if required a true pair, check this here and eventually return a fail message
+  if (fOnlyTrue)
+  {
+    // are we in a MonteCarlo?
+    if (!daughter0->GetParticle() || !daughter1->GetParticle()) return kFALSE;
     
-    TRefArray *list1, *list2;
-    if (usePID) {
-        // analysis with PID: only two collections are processed
-       list1 = event1->GetTracks(fCharge[0], fType[0]);
-       list2 = event2->GetTracks(fCharge[1], fType[1]);
-    }
-    else {
-        // analysis without PID: directly take the two arrays with all particles of a given sign
-        list1 = event1->GetCharged(fCharge[0]);
-        list2 = event2->GetCharged(fCharge[1]);
-    }
+    // are the daughters really secondaries (in MC)?
+    Int_t m0 = daughter0->GetParticle()->GetFirstMother();
+    Int_t m1 = daughter1->GetParticle()->GetFirstMother();
+    if (m0 < 0 || m1 < 0) return kFALSE;
     
-    /*
-    TRefArray *list[2];
-    Short_t i;
-    for (i = 0; i < 2; i++) {
-        if (usePID) {
-            // analysis with PID: only two collections are processed
-          list[i] = event1->GetTracks(fCharge[i], fType[i]);
-        }
-        else {
-            // analysis without PID: directly take the two arrays with all particles of a given sign
-          list[i] = event1->GetCharged(fCharge[i]);
-        }
-    }
-    Stat_t nPairs = Fill(list[0], list[1], skipSameIndex);
-    */
+    // if they are, do they come from the same mother?
+    if (m0 != m1) return kFALSE;
     
-    Stat_t nPairs = Fill(list1, list2, skipSameIndex);
-    return nPairs;
-}
-//--------------------------------------------------------------------------------------------------------
-Stat_t AliRsnPair::Fill(TRefArray *list1, TRefArray *list2, Bool_t skipSameIndex)
-{
-//
-// [PRIVATE METHOD]
-// This is the core of the pair work flow.
-// It loops on all particles contained in the two lists passed as arguments,
-// and for each pair it computes the invariant mass and fills its histogram.
-// Third argument (skipSameIndex) is a security flag: 
-// when the two lists come from the same event, two tracks with the same index 
-// must not be summed, and setting to kTRUE this flag this is ensured.
-// When it is kFALSE, this check is not done (event mixing).
-// ---
-// Before starting the operation, it checks that the two arguments are meaningful.
-// ---
-// Return value is the number of pairs processed.
-//
-
-       // define two 'cursor' objects
-       AliRsnDaughter *track1 = 0, *track2 = 0;
-       
-       // preliminary checks
-       if (!list1) {
-               AliError("List #1 cannot be NULL.");
-               return 0.0;
-       }
-    if (!list2) {
-               AliError("List #2 cannot be NULL.");
-               return 0.0;
-       }
+    // if they do, is this mother the correct type?
+    if (daughter0->GetMotherPDG() != fPairDef->GetMotherPDG()) return kFALSE;
+    if (daughter1->GetMotherPDG() != fPairDef->GetMotherPDG()) return kFALSE;
+    
+    // do they match the expected decay channel (that is, are they the expected types)?
+    Int_t pdg0 = TMath::Abs(daughter0->GetParticle()->GetPdgCode());
+    Int_t pdg1 = TMath::Abs(daughter1->GetParticle()->GetPdgCode());
+    if (AliPID::ParticleCode(fPairDef->GetPID(0)) != pdg0) return kFALSE;
+    if (AliPID::ParticleCode(fPairDef->GetPID(1)) != pdg1) return kFALSE;
     
-    // create histogram if it is not present
-       if (!fHistogram) {
-               AliError("Histogram undefined.");
-               return 0.0;
-       }
-       
-       // define iterators for the two collections
-       TRefArrayIter iter1(list1);
-       TRefArrayIter iter2(list2);
-       
-       // define temporary variables for better code readability
-       Stat_t nPairs = 0;
-       Int_t  pdgRef = TMath::Abs(fTrueMotherPDG);
-       
-       // loop on particle of type 1 (in 'event1')
-       while ( (track1 = (AliRsnDaughter*)iter1.Next()) ) {
-               // check against impact parameter cut
-        if (fVtMax > 0.0 && track1->Vt() > fVtMax) continue;
-        // check against 1-particle cuts (particle type #1)
-        if (!SingleCutCheck(0, track1)) continue;
-               // loop on particles of type #2 (in 'event2')
-               iter2.Reset();
-               while ( (track2 = (AliRsnDaughter*)iter2.Next()) ) {
-                       // skip the case when particle 2 is the same as particle 1
-                       if (skipSameIndex && (track1->Index() == track2->Index())) continue;
-            // check against impact parameter cut
-            if (fVtMax > 0.0 && track2->Vt() > fVtMax) continue;
-                       // check against 1-particle cuts (particle type #2)
-                       if (!SingleCutCheck(1, track2)) continue;
-                       // check against 2-particle cuts
-                       if (!PairCutCheck(track1, track2)) continue;
-                       // compute total 4-momentum
-                       track1->SetM(fMass[0]);
-                       track2->SetM(fMass[1]);
-                       AliRsnDaughter sum = AliRsnDaughter::Sum(*track1, *track2);
-            // reject wrong pairs with mass = 0
-            if (sum.M() <= 1E-4) continue;
-            // check transverse momentum bin
-            if (fPtMin > 0.0 && sum.Pt() < fPtMin) continue;
-            if (fPtMax > 0.0 && sum.Pt() > fPtMax) continue;
-                       // if the choice to get ONLY true pairs is selected, a check is made on the mothers
-                        if (fStoreOnlyTrue) {
-                            AliRsnParticle *part = sum.GetParticle();
-                            if (!part) continue;
-                            if (TMath::Abs(part->PDG()) != pdgRef) continue;
-            }
-                       // fill histogram
-                       fHistogram->Fill(sum.M());
-                       nPairs++;
-               } // end loop 2
-       } // end loop 1
-       
-       return nPairs;
+    // ok... if we arrive here that must really be a true pair! :-)
+  }
+
+  AliDebug(AliLog::kDebug+2,"->");
+  
+  return kTRUE;
 }
-//--------------------------------------------------------------------------------------------------------
-Bool_t AliRsnPair::SingleCutCheck(Int_t ipart, AliRsnDaughter *track) const
+
+//_____________________________________________________________________________
+void AliRsnPair::Compute()
 {
 //
-// [PRIVATE METHOD]
-// Checks a track against single particle cuts (if defined)
+// Virtual method to compute pair quantities of interest
 //
 
-       if (ipart < 0 || ipart > 1) {
-               AliError(Form("Required cuts collection for particle index %d [allowed 0 or 1]", ipart));
-               return kFALSE;
-       }
-       if (!fCutsSingle[ipart]) return kTRUE;
-    if (fCutsSingle[ipart]->IsEmpty()) return kTRUE;
-               
-       TObjArrayIter iter(fCutsSingle[ipart]);
-       AliRsnDaughterCut *cut = 0;
-       while ( (cut = (AliRsnDaughterCut*)iter.Next()) ) {
-               if (!cut->Pass(track)) return kFALSE;
-       }
-       
-       return kTRUE;
+  AliWarning("Implement this method in derived classes");
 }
-//--------------------------------------------------------------------------------------------------------
-Bool_t AliRsnPair::PairCutCheck(AliRsnDaughter *track1, AliRsnDaughter *track2) const
+
+//_____________________________________________________________________________
+void AliRsnPair::Init(const char* /*prefix*/, TList* /*list*/)
 {
 //
-// [PRIVATE METHOD]
-// Checks a pair against pair cuts (if defined)
+// Virtual method to compute pair quantities of interest
 //
-       if (fCutsPair == 0x0) return kTRUE;
-    if (fCutsPair->IsEmpty()) return kTRUE;
-       
-       TObjArrayIter iter(fCutsPair);
-       AliRsnDaughterCutPair *cut = 0;
-       while ( (cut = (AliRsnDaughterCutPair*)iter.Next()) ) {
-               if (!cut->Pass(track1, track2)) return kFALSE;
-       }
-       
-       return kTRUE;
+
+  AliWarning("Implement this method in derived classes");
 }