PWG2/RESONANCES/AliRsnFunction.cxx

   1 //
   2 // Class AliRsnFunction
   3 //
   4 // This class defines a base classe to implement a function
   5 // which uses the internal RSN package event format (AliRsnEvent).
   6 // It contains some default flags which turn out to be useful:
   7 //  - a flag to select only the "true" pairs (tracks from same resonance)
   8 //  - a flag to know if the computation is done over two events (mixing)
   9 //
  10 // Any kind of analysis object should be implemented as inheriting from this
  11 // because the AliRsnAnalyzer which executes the analysis will accept a collection
  12 // of such objects, in order to have a unique format of processing method
  13 //
  14 // The user who implements a kind of computation type should inherit from
  15 // this class and override the virtual functions defined in it, which
  16 // initialize the final output histogram and define how to process data.
  17 //
  18 //
  19 // author: A. Pulvirenti             (email: alberto.pulvirenti@ct.infn.it)
  20 //
  21
  22 #include <Riostream.h>
  23
  24 #include <TH1.h>
  25 #include <TList.h>
  26 #include <TString.h>
  27
  28 #include "AliLog.h"
  29
  30 #include "AliRsnMCInfo.h"
  31 #include "AliRsnDaughter.h"
  32 #include "AliRsnEvent.h"
  33 #include "AliRsnPairDef.h"
  34 #include "AliRsnCut.h"
  35
  36 #include "AliRsnFunction.h"
  37
  38 ClassImp(AliRsnFunction)
  39
  40 //________________________________________________________________________________________
  41 AliRsnFunction::AliRsnFunction() :
  42     fFcnType(kFcnTypes),
  43     fRotAngle(0.0),
  44     fUseBins(kFALSE),
  45     fSkipOutsideInterval(kFALSE),
  46     fBins(0),
  47     fBinningCut(),
  48     fBinningCutType(AliRsnCut::kLastCutType),
  49     fHistoDef(0x0)
  50 {
  51   //
  52   // Constructor.
  53   // The histogram data member cannot be passed externally,
  54   // its initialization MUST be defined inside the Init() method,
  55   // which must be overridden in any derivate implementation.
  56   //
  57
  58   Int_t i;
  59   for (i = 0; i < 100; i++)
  60   {
  61     fHisto[i] = 0x0;
  62   }
  63 }
  64
  65 //________________________________________________________________________________________
  66 AliRsnFunction::AliRsnFunction
  67 (EFcnType type, AliRsnHistoDef *hd, Bool_t skipOut) :
  68     fFcnType(type),
  69     fRotAngle(0.0),
  70     fUseBins(kFALSE),
  71     fSkipOutsideInterval(skipOut),
  72     fBins(0),
  73     fBinningCut(),
  74     fBinningCutType(AliRsnCut::kLastCutType),
  75     fHistoDef(hd)
  76 {
  77   //
  78   // Constructor.
  79   // The histogram data member cannot be passed externally,
  80   // its initialization MUST be defined inside the Init() method,
  81   // which must be overridden in any derivate implementation.
  82   //
  83
  84   Int_t i;
  85   for (i = 0; i < 100; i++)
  86   {
  87     fHisto[i] = 0x0;
  88   }
  89 }
  90
  91 //________________________________________________________________________________________
  92 AliRsnFunction::AliRsnFunction(const AliRsnFunction &copy) :
  93     TObject(copy),
  94     fFcnType(copy.fFcnType),
  95     fRotAngle(copy.fRotAngle),
  96     fUseBins(copy.fUseBins),
  97     fSkipOutsideInterval(copy.fSkipOutsideInterval),
  98     fBins(0),
  99     fBinningCut(),
 100     fBinningCutType(AliRsnCut::kLastCutType),
 101     fHistoDef(copy.fHistoDef)
 102 {
 103   //
 104   // Copy constructor.
 105   // Calls the function to define binning.
 106   //
 107
 108   Int_t i, n = 100;
 109   for (i = 0; i < n; i++)
 110   {
 111     fHisto[i] = 0x0;
 112   }
 113
 114   if (fUseBins)
 115   {
 116     n = copy.fBins.GetSize();
 117     Double_t *array = new Double_t[n];
 118     for (i = 0; i < n; i++) array[i] = copy.fBins[i];
 119     SetBinningCut(copy.fBinningCutType, copy.fBins.GetSize(), array);
 120     delete [] array;
 121   }
 122 }
 123 //________________________________________________________________________________________
 124 const AliRsnFunction& AliRsnFunction::operator=(const AliRsnFunction& /*copy*/)
 125 {
 126   //
 127   // Assignment operator.
 128   // Behaves like copy constructor.
 129   // Also in this case, the histogram is not copied, and,
 130   // if it was present, it is destroyed and will need to be recreated.
 131   //
 132
 133   return (*this);
 134 }
 135 //________________________________________________________________________________________
 136 void AliRsnFunction::Clear(Option_t* /*option*/)
 137 {
 138   //
 139   // Clear arrays and histogram.
 140   // For the sake of security, all pointers are also set explicitly to NULL.
 141   //
 142
 143   Int_t i;
 144   for (i = 0; i < 100; i++)
 145   {
 146     delete fHisto[i];
 147     fHisto[i] = 0x0;
 148   }
 149 }
 150
 151 //________________________________________________________________________________________
 152 TList* AliRsnFunction::Init(const char *histoName, const char *histoTitle)
 153 {
 154   //
 155   // Initialization function.
 156   // By default, it initializes the owned histogram using the method
 157   // from AliRsnHistoDef class, giving the same name and title of this.
 158   // A user can override this behaviour, if necessary.
 159   // Before creating, the HistoDef is checked for proper initialization.
 160   //
 161
 162   Clear();
 163
 164   Int_t i, ibin, nbins = fHistoDef->GetNBins();
 165   Double_t min = fHistoDef->GetMin(), max = fHistoDef->GetMax();
 166
 167   // list is created and named after the general
 168   // settings used for the contained histograms
 169   TList *histos = new TList;
 170   histos->SetName(Form("%s", GetFcnName().Data()));
 171
 172   // a general histogram is always added,
 173   // which overrides the binning and collects everything
 174   fHisto[0] = new TH1D(histoName, histoTitle, nbins, min, max);
 175   histos->AddLast(fHisto[0]);
 176
 177   // if requested a binning w.r. to some cut variable, histograms are added
 178   // for that in this part of the method (one per each bin)
 179   Char_t hName[255];
 180   Char_t hTitle[255];
 181   if (fUseBins)
 182   {
 183     for (ibin = 0, i = 1; ibin < fBins.GetSize() - 1; ibin++, i++)
 184     {
 185       sprintf(hName, "%s[%.2f-%.2f]", histoName, fBins[ibin], fBins[ibin+1]);
 186       sprintf(hTitle, "%s [%.2f-%.2f]", histoTitle, fBins[ibin], fBins[ibin+1]);
 187       fHisto[i] = new TH1D(hName, hTitle, nbins, min, max);
 188       histos->AddLast(fHisto[i]);
 189     }
 190   }
 191
 192   // returns the full list at the end
 193   return histos;
 194 }
 195
 196 //________________________________________________________________________________________
 197 void AliRsnFunction::Init(const char *histoName, const char *histoTitle, TList *histos)
 198 {
 199   //
 200   // Initialization function.
 201   // By default, it initializes the owned histogram using the method
 202   // from AliRsnHistoDef class, giving the same name and title of this.
 203   // A user can override this behaviour, if necessary.
 204   // Before creating, the HistoDef is checked for proper initialization.
 205   //
 206
 207   Clear();
 208
 209   Int_t i, ibin, nbins = fHistoDef->GetNBins();
 210   Double_t min = fHistoDef->GetMin(), max = fHistoDef->GetMax();
 211
 212   // list is created and named after the general
 213   // settings used for the contained histograms
 214   if (!histos)
 215   {
 216     AliError("NULL target list!");
 217     return;
 218   }
 219
 220   // a general histogram is always added,
 221   // which overrides the binning and collects everything
 222   fHisto[0] = new TH1D(histoName, histoTitle, nbins, min, max);
 223   histos->AddLast(fHisto[0]);
 224
 225   // if requested a binning w.r. to some cut variable, histograms are added
 226   // for that in this part of the method (one per each bin)
 227   Char_t hName[255];
 228   Char_t hTitle[255];
 229   if (fUseBins)
 230   {
 231     for (ibin = 0, i = 1; ibin < fBins.GetSize() - 1; ibin++, i++)
 232     {
 233       sprintf(hName, "%s[%.2f-%.2f]", histoName, fBins[ibin], fBins[ibin+1]);
 234       sprintf(hTitle, "%s [%.2f-%.2f]", histoTitle, fBins[ibin], fBins[ibin+1]);
 235       fHisto[i] = new TH1D(hName, hTitle, nbins, min, max);
 236       histos->AddLast(fHisto[i]);
 237     }
 238   }
 239 }
 240
 241 //________________________________________________________________________________________
 242 void AliRsnFunction::SetBinningCut
 243 (AliRsnCut::EType type, Double_t min, Double_t max, Double_t step)
 244 {
 245   //
 246   // Set fixed bins
 247   //
 248
 249   fUseBins = kTRUE;
 250
 251   Int_t i, nBins = (Int_t)((max - min) / step) + 1;
 252   fBinningCutType = type;
 253   fBins.Set(nBins);
 254   for (i = 0; i < nBins; i++)
 255   {
 256     fBins[i] = min + (Double_t)i * step;
 257   }
 258 }
 259
 260 //________________________________________________________________________________________
 261 void AliRsnFunction::SetBinningCut
 262 (AliRsnCut::EType type, Int_t nbins, Double_t *bins)
 263 {
 264   //
 265   // Set variable bins
 266   //
 267
 268   fUseBins = kTRUE;
 269
 270   Int_t i;
 271   fBinningCutType = type;
 272   fBins.Set(nbins);
 273   for (i = 0; i < nbins; i++)
 274   {
 275     fBins[i] = bins[i];
 276   }
 277 }
 278
 279 //________________________________________________________________________________________
 280 TString AliRsnFunction::GetFcnName()
 281 {
 282   //
 283   // Return a string which names the function type
 284   //
 285
 286   TString text("Undef");
 287
 288   switch (fFcnType)
 289   {
 290     case kInvMass:
 291       text = "IM";
 292       break;
 293     case kInvMassMC:
 294       text = "IM_MC";
 295       break;
 296     case kInvMassRotated:
 297       text = Form("IMR%.2f", fRotAngle);
 298       break;
 299     case kResolution:
 300       text = "RES";
 301       break;
 302     case kPtSpectrum:
 303       text = "PT";
 304       break;
 305     case kEtaSpectrum:
 306       text = "ETA";
 307       break;
 308     default:
 309       AliError("Type not defined");
 310   }
 311
 312   return text;
 313 }
 314
 315 //________________________________________________________________________________________
 316 TString AliRsnFunction::GetFcnTitle()
 317 {
 318   //
 319   // Return a string which names the function type
 320   //
 321
 322   TString text("Undef");
 323
 324   switch (fFcnType)
 325   {
 326     case kInvMass:
 327       text = "Invariant mass";
 328       break;
 329     case kInvMassMC:
 330       text = "Invariant mass (MC)";
 331       break;
 332     case kResolution:
 333       text = "Resolution";
 334       break;
 335     case kPtSpectrum:
 336       text = "p_{#perp} distribution";
 337       break;
 338     case kEtaSpectrum:
 339       text = "#eta distribution";
 340       break;
 341     default:
 342       AliError("Type not defined");
 343   }
 344
 345   return text;
 346 }
 347
 348 //________________________________________________________________________________________
 349 Bool_t AliRsnFunction::Fill(AliRsnPairParticle *pair, AliRsnPairDef *ref, Double_t weight)
 350 {
 351   //
 352   // Fillse the histogram with data contained in a defined pair.
 353   // This method must be overidden by an appropriate definition in each inheriting class.
 354   //
 355
 356   Double_t value = FcnValue(pair, ref);
 357   if (fSkipOutsideInterval)
 358   {
 359     if (value < fHistoDef->GetMin()) return kFALSE;
 360     if (value > fHistoDef->GetMax()) return kFALSE;
 361   }
 362
 363   // fill global histogram
 364   if (weight == 0.0) fHisto[0]->Fill(value);
 365   else fHisto[0]->Fill(value, weight);
 366
 367   // if bins are allocated, find right one and fill it
 368   if (fUseBins)
 369   {
 370     Int_t i, ibin;
 371     for (ibin = 0, i = 1; ibin < fBins.GetSize() - 1; ibin++, i++)
 372     {
 373       if (!fHisto[i]) continue;
 374       fBinningCut.SetCutValues(fBinningCutType, fBins[ibin], fBins[ibin+1]);
 375       if (fBinningCut.IsSelected(AliRsnCut::kPair, pair))
 376       {
 377         if (weight == 0.0) fHisto[i]->Fill(value);
 378         else fHisto[i]->Fill(value, weight);
 379         break;
 380       }
 381     }
 382   }
 383
 384   return kTRUE;
 385 }
 386
 387 //________________________________________________________________________________________
 388 Double_t AliRsnFunction::FcnValue(AliRsnPairParticle *pair, AliRsnPairDef *ref)
 389 {
 390   //
 391   // This method must be overridden in all inheritin functions.
 392   // It computes the value which must be used to fill the histogram.
 393   //
 394
 395   switch (fFcnType)
 396   {
 397     case kInvMass:
 398       return pair->GetInvMass(ref->GetMass(0), ref->GetMass(1));
 399     case kInvMassMC:
 400       return pair->GetInvMassMC(ref->GetMass(0), ref->GetMass(1));
 401     case kInvMassRotated:
 402       //AliInfo(Form("*** ROTATION ANGLE = %f ***", fRotAngle));
 403       //AliInfo(Form("UNROTATED INV MASS = %f", pair->GetInvMass(ref->GetMass(0), ref->GetMass(1))));
 404       //pair->GetDaughter(1)->Print("P");
 405       pair->GetDaughter(1)->RotateP(fRotAngle * TMath::DegToRad());
 406       pair->ResetPair();
 407       //AliInfo(Form("  ROTATED INV MASS = %f", pair->GetInvMass(ref->GetMass(0), ref->GetMass(1))));
 408       //pair->GetDaughter(1)->Print("P");
 409       return pair->GetInvMass(ref->GetMass(0), ref->GetMass(1));
 410     case kResolution:
 411       return FcnResolution(pair, ref);
 412     case kPtSpectrum:
 413       return pair->GetPt();
 414     case kEtaSpectrum:
 415       return pair->GetEta();
 416     default:
 417       AliError("Type not defined");
 418   }
 419
 420   return 0.0;
 421 }
 422
 423 //________________________________________________________________________________________
 424 inline Double_t AliRsnFunction::FcnResolution(AliRsnPairParticle *pair, AliRsnPairDef *ref)
 425 {
 426   //
 427   // Invariant mass resolution (compared between reconstructed and montecarlo)
 428   //
 429
 430   Double_t recInvMass = pair->GetInvMass(ref->GetMass(0), ref->GetMass(1));
 431   Double_t simInvMass = pair->GetInvMassMC(ref->GetMass(0), ref->GetMass(1));
 432
 433   return (simInvMass - recInvMass) / simInvMass;
 434 }