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 <Riostream.h>
  11 #include "AliESDtrackCuts.h"
  12 #include "AliRsnEvent.h"
  13 #include "AliRsnDaughter.h"
  14 #include "AliRsnMother.h"
  15 #include "AliRsnPairDef.h"
  16
  17 #include "AliRsnValue.h"
  18
  19 ClassImp(AliRsnValue)
  20
  21 //_____________________________________________________________________________
  22 AliRsnValue::AliRsnValue() :
  23    AliRsnTarget(),
  24    fComputedValue(0),
  25    fValueType(kValueTypes),
  26    fBinArray(0),
  27    fSupportObject(0x0)
  28 {
  29 //
  30 // Default constructor without arguments.
  31 // Initialize data members to meaningless values.
  32 // This method is provided for ROOT streaming,
  33 // but should never be used directly by a user.
  34 //
  35
  36    AssignTarget();
  37 }
  38
  39 //_____________________________________________________________________________
  40 AliRsnValue::AliRsnValue
  41 (const char *name, EValueType type, Int_t nbins, Double_t min, Double_t max) :
  42    AliRsnTarget(name, AliRsnTarget::kTargetTypes),
  43    fComputedValue(0.0),
  44    fValueType(type),
  45    fBinArray(0),
  46    fSupportObject(0x0)
  47 {
  48 //
  49 // Main constructor (version 1).
  50 // This constructor defines in meaningful way all data members,
  51 // and defined a fixed binnings, subdividing the specified interval
  52 // into that many bins as specified in the integer argument.
  53 // ---
  54 // This method is also the entry point for all instances
  55 // of this class which don't need to do binning (e.g.: TNtuple inputs),
  56 // since arguments 3 to 5 have default values which don't create any
  57 // binning array, in order not to allocate memory when this is useless.
  58 //
  59
  60    AssignTarget();
  61    SetBins(nbins, min, max);
  62 }
  63
  64 //_____________________________________________________________________________
  65 AliRsnValue::AliRsnValue
  66 (const char *name, EValueType type, Double_t min, Double_t max, Double_t step) :
  67    AliRsnTarget(name, AliRsnTarget::kTargetTypes),
  68    fComputedValue(0.0),
  69    fValueType(type),
  70    fBinArray(0),
  71    fSupportObject(0x0)
  72 {
  73 //
  74 // Main constructor (version 2).
  75 // This constructor defines in meaningful way all data members
  76 // and creates enough equal bins of the specified size to cover
  77 // the required interval.
  78 //
  79
  80    AssignTarget();
  81    SetBins(min, max, step);
  82 }
  83
  84 //_____________________________________________________________________________
  85 AliRsnValue::AliRsnValue
  86 (const char *name, EValueType type, Int_t nbins, Double_t *array) :
  87    AliRsnTarget(name, AliRsnTarget::kTargetTypes),
  88    fComputedValue(0.0),
  89    fValueType(type),
  90    fBinArray(0),
  91    fSupportObject(0x0)
  92 {
  93 //
  94 // Main constructor (version 3).
  95 // This constructor defines in meaningful way all data members
  96 // and creates a set of variable bins delimited by the passed array.
  97 //
  98
  99    AssignTarget();
 100    SetBins(nbins, array);
 101 }
 102
 103 //_____________________________________________________________________________
 104 AliRsnValue::AliRsnValue(const AliRsnValue& copy) :
 105    AliRsnTarget(copy),
 106    fComputedValue(copy.fComputedValue),
 107    fValueType(copy.fValueType),
 108    fBinArray(copy.fBinArray),
 109    fSupportObject(copy.fSupportObject)
 110 {
 111 //
 112 // Copy constructor.
 113 // Duplicates the binning array and copies all settings.
 114 // Calls also the function that assigns properly the
 115 // expected target, depending on the computation type.
 116 //
 117
 118    AssignTarget();
 119 }
 120
 121 //_____________________________________________________________________________
 122 AliRsnValue& AliRsnValue::operator=(const AliRsnValue& copy)
 123 {
 124 //
 125 // Assignment operator.
 126 // Works like copy constructor.
 127 //
 128
 129    AliRsnTarget::operator=(copy);
 130
 131    fComputedValue = copy.fComputedValue;
 132    fBinArray = copy.fBinArray;
 133    fSupportObject = copy.fSupportObject;
 134
 135    AssignTarget();
 136
 137    return (*this);
 138 }
 139
 140 //_____________________________________________________________________________
 141 void AliRsnValue::SetBins(Int_t nbins, Double_t min, Double_t max)
 142 {
 143 //
 144 // Set binning for the axis in equally spaced bins
 145 // where the number of bins, minimum and maximum are given.
 146 //
 147
 148    if (!nbins) {
 149       fBinArray.Set(0);
 150       return;
 151    }
 152
 153    fBinArray.Set(nbins + 1);
 154
 155    Double_t mymax = TMath::Max(min, max);
 156    Double_t mymin = TMath::Min(min, max);
 157
 158    Int_t    k = 0;
 159    Double_t binSize = (mymax - mymin) / ((Double_t)nbins);
 160
 161    fBinArray[0] = mymin;
 162    for (k = 1; k <= nbins; k++) fBinArray[k] = fBinArray[k - 1] + binSize;
 163 }
 164
 165 //_____________________________________________________________________________
 166 void AliRsnValue::SetBins(Double_t min, Double_t max, Double_t step)
 167 {
 168 //
 169 // Set binning for the axis in equally spaced bins
 170 // where the bin size, minimum and maximum are given.
 171 //
 172
 173    Double_t dblNbins = TMath::Abs(max - min) / step;
 174    Int_t    intNbins = ((Int_t)dblNbins) + 1;
 175
 176    SetBins(intNbins, min, max);
 177 }
 178
 179 //_____________________________________________________________________________
 180 void AliRsnValue::SetBins(Int_t nbins, Double_t *array)
 181 {
 182 //
 183 // Set binning for the axis in unequally spaced bins
 184 // using the same way it is done in TAxis
 185 //
 186
 187    if (!nbins) {
 188       fBinArray.Set(0);
 189       return;
 190    }
 191
 192    fBinArray.Adopt(nbins, array);
 193 }
 194
 195 //_____________________________________________________________________________
 196 const char* AliRsnValue::GetValueTypeName() const
 197 {
 198 //
 199 // This method returns a string to give a name to each possible
 200 // computation value.
 201 //
 202
 203    switch (fValueType) {
 204       case kTrackP:             return "SingleTrackPtot";
 205       case kTrackPt:            return "SingleTrackPt";
 206       case kTrackEta:           return "SingleTrackEta";
 207       case kPairP1:             return "PairPtotDaughter1";
 208       case kPairP2:             return "PairPtotDaughter2";
 209       case kPairP1t:            return "PairPtDaughter1";
 210       case kPairP2t:            return "PairPtDaughter2";
 211       case kPairP1z:            return "PairPzDaughter1";
 212       case kPairP2z:            return "PairPzDaughter2";
 213       case kPairInvMass:        return "PairInvMass";
 214       case kPairInvMassMC:      return "PairInvMassMC";
 215       case kPairInvMassRes:     return "PairInvMassResolution";
 216       case kPairPt:             return "PairPt";
 217       case kPairPz:             return "PairPz";
 218       case kPairEta:            return "PairEta";
 219       case kPairMt:             return "PairMt";
 220       case kPairY:              return "PairY";
 221       case kPairPhi:            return "PairPhi";
 222       case kPairPhiMC:          return "PairPhiMC";
 223       case kPairPtRatio:        return "PairPtRatio";
 224       case kPairDipAngle:       return "PairDipAngle";
 225       case kPairCosThetaStar:   return "PairCosThetaStar";
 226       case kPairQInv:           return "PairQInv";
 227       case kPairAngleToLeading: return "PairAngleToLeading";
 228       case kEventLeadingPt:     return "EventLeadingPt";
 229       case kEventMult:          return "EventMult";
 230       case kEventMultESDCuts:   return "EventMultESDCuts";
 231       case kEventVz:            return "EventVz";
 232       default:                  return "Undefined";
 233    }
 234 }
 235
 236 //_____________________________________________________________________________
 237 void AliRsnValue::AssignTarget()
 238 {
 239 //
 240 // This method assigns the target to be expected by this object
 241 // in the computation, depending on its type chosen in the enum.
 242 //
 243
 244    switch (fValueType) {
 245          // track related values
 246       case kTrackP:
 247       case kTrackPt:
 248       case kTrackEta:
 249          SetTargetType(AliRsnTarget::kDaughter); // end of track-related values
 250          break;
 251          // pair related values
 252       case kPairP1:
 253       case kPairP2:
 254       case kPairP1t:
 255       case kPairP2t:
 256       case kPairP1z:
 257       case kPairP2z:
 258       case kPairInvMass:
 259       case kPairInvMassMC:
 260       case kPairInvMassRes:
 261       case kPairPt:
 262       case kPairPz:
 263       case kPairEta:
 264       case kPairMt:
 265       case kPairY:
 266       case kPairPhi:
 267       case kPairPhiMC:
 268       case kPairPtRatio:
 269       case kPairDipAngle:
 270       case kPairCosThetaStar:
 271       case kPairQInv:
 272       case kPairAngleToLeading:
 273          SetTargetType(AliRsnTarget::kMother); // end of pair-related values
 274          break;
 275          // event related values
 276       case kEventLeadingPt:
 277       case kEventMult:
 278       case kEventMultESDCuts:
 279       case kEventVz:
 280          SetTargetType(AliRsnTarget::kEvent); // end of event-related values
 281          break;
 282          // undefined value
 283       default:
 284          SetTargetType(AliRsnTarget::kTargetTypes); // undefined targets
 285    }
 286 }
 287
 288 //_____________________________________________________________________________
 289 Bool_t AliRsnValue::Eval(TObject *object, Bool_t useMC)
 290 {
 291 //
 292 // Evaluation of the required value.
 293 // Checks that the passed object is of the right type
 294 // and if this check is successful, computes the required value.
 295 // The output of the function tells if computing was successful,
 296 // and the values must be taken with GetValue().
 297 //
 298
 299    // cast the input to the allowed types
 300    AliRsnDaughter *daughter = dynamic_cast<AliRsnDaughter*>(object);
 301    AliRsnMother   *mother   = dynamic_cast<AliRsnMother*>(object);
 302
 303    // common variables
 304    TLorentzVector pRec;   // 4-momentum for single track or pair sum (reco)
 305    TLorentzVector pSim;   // 4-momentum for single track or pair sum (MC)
 306    TLorentzVector pRec0;  // 4-momentum of first daughter (reco)
 307    TLorentzVector pSim0;  // 4-momentum of first daughter (MC)
 308    TLorentzVector pRec1;  // 4-momentum of second daughter (reco)
 309    TLorentzVector pSim1;  // 4-momentum of second daughter (MC)
 310
 311    // check that the input object is the correct class type
 312    switch (fTargetType) {
 313       case AliRsnTarget::kDaughter:
 314          if (daughter) {
 315             pRec = daughter->Psim();
 316             pSim = daughter->Prec();
 317          } else {
 318             AliError(Form("[%s] expected: AliRsnDaughter, passed: [%s]", GetName(), object->ClassName()));
 319             return kFALSE;
 320          }
 321          break;
 322       case AliRsnTarget::kMother:
 323          if (mother) {
 324             pRec  = mother->Sum();
 325             pSim  = mother->SumMC();
 326             pRec0 = mother->GetDaughter(0)->Prec();
 327             pRec1 = mother->GetDaughter(1)->Prec();
 328             pSim0 = mother->GetDaughter(0)->Psim();
 329             pSim1 = mother->GetDaughter(1)->Psim();
 330          } else {
 331             AliError(Form("[%s] expected: AliRsnMother, passed: [%s]", GetName(), object->ClassName()));
 332             return kFALSE;
 333          }
 334          break;
 335       case AliRsnTarget::kEvent:
 336          if (!AliRsnTarget::GetCurrentEvent()) {
 337             AliError(Form("[%s] expected: AliRsnEvent, passed: [%s]", GetName(), object->ClassName()));
 338             return kFALSE;
 339          }
 340          break;
 341       default:
 342          AliError(Form("[%s] Wrong type", GetName()));
 343          return kFALSE;
 344    }
 345
 346    // cast the support object to the types which could be needed
 347    AliESDtrackCuts *esdCuts = dynamic_cast<AliESDtrackCuts*>(fSupportObject);
 348    AliRsnPairDef   *pairDef = dynamic_cast<AliRsnPairDef*>(fSupportObject);
 349
 350    // compute value depending on type
 351    switch (fValueType) {
 352       case kTrackP:
 353          fComputedValue = useMC ? pSim.Mag() : pRec.Mag();
 354          break;
 355       case kTrackPt:
 356          fComputedValue = useMC ? pSim.Perp() : pRec.Perp();
 357          break;
 358       case kTrackEta:
 359          fComputedValue = useMC ? pSim.Eta() : pRec.Eta();
 360          break;
 361       case kPairP1:
 362          fComputedValue = useMC ? pSim0.Mag() : pRec0.Mag();
 363          break;
 364       case kPairP2:
 365          fComputedValue = useMC ? pSim1.Mag() : pRec1.Mag();
 366          break;
 367       case kPairP1t:
 368          fComputedValue = useMC ? pSim0.Perp() : pRec0.Perp();
 369          break;
 370       case kPairP2t:
 371          fComputedValue = useMC ? pSim1.Perp() : pRec1.Perp();
 372          break;
 373       case kPairP1z:
 374          fComputedValue = useMC ? pSim0.Z() : pRec0.Z();
 375          break;
 376       case kPairP2z:
 377          fComputedValue = useMC ? pSim1.Z() : pRec1.Z();
 378          break;
 379       case kPairInvMass:
 380          fComputedValue = useMC ? pSim.M() : pRec.M();
 381          break;
 382       case kPairInvMassRes:
 383          fComputedValue = (pSim.M() - pRec.M()) / pSim.M();
 384          break;
 385       case kPairPt:
 386          fComputedValue = useMC ? pSim.Perp() : pRec.Perp();
 387          break;
 388       case kPairEta:
 389          fComputedValue = useMC ? pSim.Eta() : pRec.Eta();
 390          break;
 391       case kPairMt:
 392          // for this computation, replace the computed mass with the default mass
 393          // for doing this, an initialized pairDef is required to get the mass
 394          if (!pairDef) {
 395             AliError(Form("[%s] Required a correctly initialized PairDef to compute this value", GetName()));
 396             fComputedValue = 1E+10;
 397             return kFALSE;
 398          } else {
 399             pRec.SetXYZM(pRec.X(), pRec.Y(), pRec.Z(), pairDef->GetMotherMass());
 400             pSim.SetXYZM(pSim.X(), pSim.Y(), pSim.Z(), pairDef->GetMotherMass());
 401             fComputedValue = useMC ? pSim.Mt() : pRec.Mt();
 402          }
 403          break;
 404       case kPairY:
 405          // for this computation, replace the computed mass with the default mass
 406          // for doing this, an initialized pairDef is required to get the mass
 407          if (!pairDef) {
 408             AliError(Form("[%s] Required a correctly initialized PairDef to compute this value", GetName()));
 409             fComputedValue = 1E+10;
 410             return kFALSE;
 411          } else {
 412             pRec.SetXYZM(pRec.X(), pRec.Y(), pRec.Z(), pairDef->GetMotherMass());
 413             pSim.SetXYZM(pSim.X(), pSim.Y(), pSim.Z(), pairDef->GetMotherMass());
 414             fComputedValue = useMC ? pSim.Rapidity() : pRec.Rapidity();
 415          }
 416          break;
 417       case kPairPhi:
 418          fComputedValue = useMC ? pSim.Phi() : pRec.Phi();
 419          break;
 420       case kPairPtRatio:
 421          if (useMC) {
 422             fComputedValue  = TMath::Abs(pSim0.Perp() - pSim1.Perp());
 423             fComputedValue /= TMath::Abs(pSim0.Perp() + pSim1.Perp());
 424          } else {
 425             fComputedValue  = TMath::Abs(pRec0.Perp() - pRec1.Perp());
 426             fComputedValue /= TMath::Abs(pRec0.Perp() + pRec1.Perp());
 427          }
 428          break;
 429       case kPairDipAngle:
 430          fComputedValue = useMC ? pSim0.Angle(pSim1.Vect()) : pRec0.Angle(pRec1.Vect());
 431          fComputedValue = TMath::Abs(TMath::Cos(fComputedValue));
 432          break;
 433       case kPairCosThetaStar:
 434          fComputedValue = mother->CosThetaStar(useMC);
 435          break;
 436       case kPairQInv:
 437          pSim0 -= pSim1;
 438          pRec0 -= pRec1;
 439          fComputedValue = useMC ? pSim0.M() : pRec0.M();
 440          break;
 441       case kPairAngleToLeading: {
 442          AliRsnEvent *event = AliRsnTarget::GetCurrentEvent();
 443          int ID1 = (mother->GetDaughter(0))->GetID();
 444          int ID2 = (mother->GetDaughter(1))->GetID();
 445          //int leadingID = event->SelectLeadingParticle(0);
 446          Int_t leadingID = event->GetLeadingParticleID();
 447          if (leadingID == ID1 || leadingID == ID2) return kFALSE;
 448          AliRsnDaughter leadingPart = event->GetDaughter(leadingID);
 449          AliVParticle  *ref = leadingPart.GetRef();
 450          fComputedValue = ref->Phi() - mother->Sum().Phi();
 451          //return angle w.r.t. leading particle in the range -pi/2, 3/2pi
 452          while (fComputedValue >= TMath::Pi()) fComputedValue -= 2 * TMath::Pi();
 453          while (fComputedValue < -0.5 * TMath::Pi()) fComputedValue += 2 * TMath::Pi();
 454          //Printf("%g", fComputedValue);
 455       }
 456       break;
 457       case kEventMult:
 458          fComputedValue = (Double_t)AliRsnTarget::GetCurrentEvent()->GetMultiplicity(0x0);
 459          break;
 460       case kEventMultESDCuts:
 461          // this value requires an initialized ESDtrackCuts
 462          if (!esdCuts) {
 463             AliError(Form("[%s] Required a correctly initialized ESDtrackCuts to compute this value", GetName()));
 464             fComputedValue = 1E+10;
 465             return kFALSE;
 466          }
 467          fComputedValue = (Double_t)AliRsnTarget::GetCurrentEvent()->GetMultiplicity(esdCuts);
 468          break;
 469       case kEventLeadingPt: {
 470          int leadingID = AliRsnTarget::GetCurrentEvent()->GetLeadingParticleID(); //fEvent->SelectLeadingParticle(0);
 471          if (leadingID >= 0) {
 472             AliRsnDaughter leadingPart = AliRsnTarget::GetCurrentEvent()->GetDaughter(leadingID);
 473             AliVParticle *ref = leadingPart.GetRef();
 474             fComputedValue = ref->Pt();
 475          } else fComputedValue = 0;
 476       }
 477       break;
 478       case kEventVz:
 479          fComputedValue = AliRsnTarget::GetCurrentEvent()->GetRef()->GetPrimaryVertex()->GetZ();
 480          break;
 481       default:
 482          AliError(Form("[%s] Invalid value type for this computation", GetName()));
 483          return kFALSE;
 484    }
 485
 486    return kTRUE;
 487 }
 488
 489 //_____________________________________________________________________________
 490 void AliRsnValue::Print(Option_t * /*option */) const
 491 {
 492 //
 493 // Print informations about this object
 494 //
 495
 496    AliInfo("=== VALUE INFO =================================================");
 497    AliInfo(Form(" Name                  : %s", GetName()));
 498    AliInfo(Form(" Type                  : %s", GetValueTypeName()));
 499    AliInfo(Form(" Current computed value: %f", fComputedValue));
 500    Int_t i;
 501    for (i = 0; i < fBinArray.GetSize(); i++) {
 502       AliInfo(Form(" Bin limit #%d         = %f", i, fBinArray[i]));
 503    }
 504    AliInfo(Form(" Support object        : %s", (fSupportObject ? fSupportObject->ClassName() : " NO SUPPORT")));
 505    AliInfo("=== END VALUE INFO =============================================");
 506 }