PWG2/RESONANCES/AliRsnValue.cxx

   1
   2 //
   3 // Class AliRsnValue
   4 //
   5 // Definition of a single value which can be computed
   6 // from any of the defined input objects implemented
   7 // in the resonance package.
   8 //
   9
  10 #include "AliRsnEvent.h"
  11 #include "AliRsnDaughter.h"
  12 #include "AliRsnMother.h"
  13 #include "AliRsnPairDef.h"
  14
  15 #include "AliRsnValue.h"
  16
  17 ClassImp(AliRsnValue)
  18
  19 //_____________________________________________________________________________
  20 AliRsnValue::AliRsnValue() :
  21   TNamed(),
  22   fType(kValueTypes),
  23   fValue(0.0),
  24   fArray(0)
  25 {
  26 //
  27 // Main constructor (version 1)
  28 // This can also be created without any argument.
  29 //
  30 }
  31
  32 //_____________________________________________________________________________
  33 AliRsnValue::AliRsnValue
  34 (const char *name, EValueType type, Int_t nbins, Double_t min, Double_t max) :
  35   TNamed(name, ""),
  36   fType(type),
  37   fValue(0.0),
  38   fArray(0)
  39 {
  40 //
  41 // Main constructor (version 1)
  42 // This can also be created without any argument.
  43 //
  44
  45   SetBins(nbins, min, max);
  46 }
  47
  48 //_____________________________________________________________________________
  49 AliRsnValue::AliRsnValue
  50 (const char *name, EValueType type, Double_t min, Double_t max, Double_t step) :
  51   TNamed(name, ""),
  52   fType(type),
  53   fValue(0.0),
  54   fArray(0)
  55 {
  56 //
  57 // Main constructor (version 2)
  58 //
  59
  60   SetBins(min, max, step);
  61 }
  62
  63 //_____________________________________________________________________________
  64 AliRsnValue::AliRsnValue
  65 (const char *name, EValueType type, Int_t nbins, Double_t *array) :
  66   TNamed(name, ""),
  67   fType(type),
  68   fValue(0.0),
  69   fArray(0)
  70 {
  71 //
  72 // Main constructor (version 2)
  73 //
  74
  75   SetBins(nbins, array);
  76 }
  77
  78 //_____________________________________________________________________________
  79 void AliRsnValue::SetBins(Int_t nbins, Double_t min, Double_t max)
  80 {
  81 //
  82 // Set binning for the axis in equally spaced bins
  83 // where the number of bins, minimum and maximum are given.
  84 //
  85
  86   fArray.Set(nbins + 1);
  87
  88   Double_t mymax = TMath::Max(min, max);
  89   Double_t mymin = TMath::Min(min, max);
  90
  91   Int_t    k = 0;
  92   Double_t binSize = (mymax - mymin) / ((Double_t)nbins);
  93
  94   fArray[0] = mymin;
  95   for (k = 1; k <= nbins; k++) fArray[k] = fArray[k-1] + binSize;
  96   for (k = 0; k < fArray.GetSize() - 1; k++) AliDebug(AliLog::kDebug + 3, Form("Bin #%d: %f - %f", k, fArray[k], fArray[k+1]));
  97 }
  98
  99 //_____________________________________________________________________________
 100 void AliRsnValue::SetBins(Double_t min, Double_t max, Double_t step)
 101 {
 102 //
 103 // Set binning for the axis in equally spaced bins
 104 // where the bin size, minimum and maximum are given.
 105 //
 106
 107   Double_t dblNbins = TMath::Abs(max - min) / step;
 108   Int_t    intNbins = ((Int_t)dblNbins) + 1;
 109
 110   SetBins(intNbins, min, max);
 111 }
 112
 113 //_____________________________________________________________________________
 114 void AliRsnValue::SetBins(Int_t nbins, Double_t *array)
 115 {
 116 //
 117 // Set binning for the axis in unequally spaced bins
 118 // using the same way it is done in TAxis
 119 //
 120
 121   fArray.Adopt(nbins, array);
 122   for (Int_t k = 0; k < fArray.GetSize() - 1; k++) AliDebug(AliLog::kDebug + 3, Form("Bin #%d: %f - %f", k, fArray[k], fArray[k+1]));
 123 }
 124
 125 //_____________________________________________________________________________
 126 Bool_t AliRsnValue::Eval(AliRsnMother * const mother, AliRsnPairDef * const pairDef, AliRsnEvent * const event)
 127 {
 128 //
 129 // Evaluation of the required value.
 130 // Checks that the passed object is of the right type
 131 // and if this check is successful, returns the required value.
 132 // The output of the function tells if it was successful,
 133 // and the values must be taken with GetValue().
 134 //
 135
 136   // avoid segfaults
 137   if (!mother) return kFALSE;
 138   if (!pairDef) return kFALSE;
 139
 140   Double_t mass = pairDef->GetMotherMass();
 141
 142   switch (fType)
 143   {
 144     case kTrack1P:
 145       fValue = mother->GetDaughter(0)->P().Mag();
 146       break;
 147     case kTrack2P:
 148       fValue = mother->GetDaughter(1)->P().Mag();
 149       break;
 150     case kTrack1Pt:
 151       fValue = mother->GetDaughter(0)->P().Perp();
 152       break;
 153     case kTrack2Pt:
 154       fValue = mother->GetDaughter(1)->P().Perp();
 155       break;
 156     case kTrack1Px:
 157       fValue = mother->GetDaughter(0)->P().X();
 158       break;
 159     case kTrack1Py:
 160       fValue = mother->GetDaughter(0)->P().Y();
 161       break;
 162     case kTrack1Pz:
 163       fValue = mother->GetDaughter(0)->P().Z();
 164       break;
 165     case kTrack2Px:
 166       fValue = mother->GetDaughter(1)->P().X();
 167       break;
 168     case kTrack2Py:
 169       fValue = mother->GetDaughter(1)->P().Y();
 170       break;
 171     case kTrack2Pz:
 172       fValue = mother->GetDaughter(1)->P().Z();
 173       break;
 174     case kPairInvMass:
 175       fValue = mother->Sum().M();
 176       break;
 177     case kPairInvMassMC:
 178       fValue = mother->SumMC().M();
 179       break;
 180     case kPairInvMassRes:
 181       fValue = (mother->SumMC().M() - mother->Sum().M()) / mother->SumMC().M();
 182       break;
 183     case kPairPt:
 184       fValue = mother->Sum().Perp();
 185       break;
 186     case kPairEta:
 187       fValue = mother->Sum().Eta();
 188       break;
 189     case kPairMt:
 190       if (TMath::Abs(mass) < 1E-5) AliWarning(Form("Suspicious mass value specified: %f", mass));
 191       fValue = (TMath::Sqrt(mother->Sum().Perp2() + mass*mass) - mass);
 192       break;
 193     case kPairY:
 194       if (TMath::Abs(mass) < 1E-5) AliWarning(Form("Suspicious mass value specified: %f", mass));
 195       mother->SetDefaultMass(mass);
 196       fValue = mother->Ref().Rapidity();
 197       break;
 198     case kPairCosThetaStar:
 199       fValue = mother->CosThetaStar();
 200       break;
 201     case kPairCosThetaStar1:
 202       //fValue = TMath::Cos(mother->ThetaStar(kTRUE, kFALSE));
 203       break;
 204     case kPairCosThetaStar2:
 205       //fValue = TMath::Cos(mother->ThetaStar(kFALSE, kFALSE));
 206       break;
 207     case kPairCosThetaStarMC1:
 208       //fValue = TMath::Cos(mother->ThetaStar(kTRUE, kTRUE));
 209       break;
 210     case kPairCosThetaStarMC2:
 211       //fValue = TMath::Cos(mother->ThetaStar(kFALSE, kTRUE));
 212       break;
 213     case kAngleToLeading:
 214       {
 215           int ID1 = (mother->GetDaughter(0))->GetID();
 216           int ID2 = (mother->GetDaughter(1))->GetID();
 217           int leadingID = event->SelectLeadingParticle(0);
 218           if(leadingID == ID1 || leadingID == ID2) return kFALSE;
 219           AliRsnDaughter  leadingPart = event->GetDaughter(leadingID);
 220           AliVParticle *ref = leadingPart.GetRef();
 221
 222           fValue = ref->Phi() - mother->Sum().Phi();
 223           //return angle w.r.t. leading particle in the range -pi/2, 3/2pi
 224           while(fValue >= TMath::Pi()) fValue -= 2*TMath::Pi();
 225           while(fValue < -0.5*TMath::Pi()) fValue += 2*TMath::Pi();
 226           //Printf("%g", fValue);
 227
 228       }
 229       break;
 230     case kEventMult:
 231       if (!event) fValue = 0.0;
 232       fValue = (Double_t)event->GetMultiplicity();
 233       break;
 234     case kLeadingPt:
 235       {
 236           int leadingID = event->SelectLeadingParticle(0);
 237           if(leadingID >= 0) {
 238                   AliRsnDaughter leadingPart = event->GetDaughter(leadingID);
 239                   AliVParticle *ref = leadingPart.GetRef();
 240                   fValue = ref->Pt();
 241           }
 242           else fValue = 0;
 243       }
 244       break;
 245     case kQInv:
 246       {
 247         TLorentzVector diff = mother->GetDaughter(0)->P() - mother->GetDaughter(1)->P();
 248         fValue = diff.M();
 249       }
 250       break;
 251     default:
 252       AliWarning("Invalid value type");
 253       return kFALSE;
 254   }
 255
 256   return kTRUE;
 257 }