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 <TString.h>
  23 #include <TAxis.h>
  24
  25 #include "AliLog.h"
  26
  27 #include "AliRsnDaughter.h"
  28 #include "AliRsnEvent.h"
  29 #include "AliRsnPairDef.h"
  30 #include "AliRsnMother.h"
  31 #include "AliRsnValue.h"
  32
  33 #include "AliRsnFunction.h"
  34
  35 ClassImp(AliRsnFunction)
  36
  37 //________________________________________________________________________________________
  38 AliRsnFunction::AliRsnFunction(Bool_t useTH1) :
  39   TObject(),
  40   fPairDef(0x0),
  41   fAxisList("AliRsnValue", 0),
  42   fPair(0x0),
  43   fEvent(0x0),
  44   fUseTH1(useTH1),
  45   fSize(0),
  46   fH1(0x0),
  47   fHSparse(0x0)
  48 {
  49 //
  50 // Constructor.
  51 //
  52 }
  53
  54 //________________________________________________________________________________________
  55 AliRsnFunction::AliRsnFunction(const AliRsnFunction &copy) :
  56   TObject(copy),
  57   fPairDef(copy.fPairDef),
  58   fAxisList(copy.fAxisList),
  59   fPair(copy.fPair),
  60   fEvent(copy.fEvent),
  61   fUseTH1(copy.fUseTH1),
  62   fSize(copy.fSize),
  63   fH1(0x0),
  64   fHSparse(0x0)
  65 {
  66 //
  67 // Copy constructor.
  68 //
  69 }
  70
  71 //________________________________________________________________________________________
  72 const AliRsnFunction& AliRsnFunction::operator=(const AliRsnFunction& copy)
  73 {
  74 //
  75 // Assignment operator.
  76 //
  77
  78   fPairDef = copy.fPairDef;
  79   fPair = copy.fPair;
  80   fEvent = copy.fEvent;
  81   fUseTH1 = copy.fUseTH1;
  82   fSize = copy.fSize;
  83
  84   if (fH1) delete fH1;
  85   fH1 = 0x0;
  86
  87   if (fHSparse) delete fHSparse;
  88   fHSparse = 0x0;
  89
  90   return (*this);
  91 }
  92
  93 //________________________________________________________________________________________
  94 const char* AliRsnFunction::GetName() const
  95 {
  96 //
  97 // Defines the name of this object according to
  98 // the function type and binning
  99 //
 100
 101   TString name("");
 102
 103   TObjArrayIter next(&fAxisList);
 104   AliRsnValue *axis = 0;
 105
 106   while ((axis = (AliRsnValue*)next()))
 107   {
 108     if (name.Length() > 1) name += '_';
 109     name += axis->GetName();
 110   }
 111
 112   return name.Data();
 113 }
 114
 115 //________________________________________________________________________________________
 116 Bool_t AliRsnFunction::AddAxis(AliRsnValue *const axis)
 117 {
 118 //
 119 // Try to add a new axis to this function.
 120 // The operation succeeds only if the related value object
 121 // has a target amon those allowed here (AliRsnMother, AliRsnEvent),
 122 // otherwise the axis is not added.
 123 //
 124 // If more than 3 axes are added, switch to THnSparseF output.
 125 // NOTE: this can cause large files and high memory occupancy.
 126 //
 127
 128   RSNTARGET target = axis->GetTargetType();
 129   if (target != AliRsnTarget::kMother && target != AliRsnTarget::kEvent)
 130   {
 131     AliError(Form("Allowed targets are mothers and events; cannot use axis '%s' which has target '%s'", axis->GetName(), axis->GetTargetTypeName()));
 132     return kFALSE;
 133   }
 134
 135   Int_t size = fAxisList.GetEntries();
 136   new (fAxisList[size]) AliRsnValue(*axis);
 137
 138   if (fAxisList.GetEntries() > 3)
 139   {
 140     AliWarning("Adding more than 3 axes will produce a THnSparseD output.");
 141     fUseTH1 = kFALSE;
 142   }
 143
 144   return kTRUE;
 145 }
 146
 147 //________________________________________________________________________________________
 148 TH1* AliRsnFunction::CreateHistogram(const char *histoName, const char *histoTitle)
 149 {
 150 //
 151 // Creates and returns the histogram defined using
 152 // arguments fo name and title, and the first histoDef for binning.
 153 // Variable-sized histogram binning is always called, due to use of histoDef,
 154 // even if the bins are equal, since they are defined in this class.
 155 // Eventually present histoDef's in other slots of array (1, 2) are ignored.
 156 //
 157 // This version produces a THnSparseF.
 158 //
 159
 160   fSize = fAxisList.GetEntries();
 161   if (!fSize)
 162   {
 163     AliError("No axes defined");
 164     return 0x0;
 165   }
 166   else if (fSize > 3)
 167   {
 168     AliError("Too many axes defined for a TH1 object");
 169     return 0x0;
 170   }
 171
 172   // retrieve binnings for main and secondary axes
 173   AliRsnValue *fcnAxis;
 174   TArrayD      array[3];
 175   for (Int_t i = 0; i < fSize; i++)
 176   {
 177     fcnAxis = (AliRsnValue*)fAxisList.At(i);
 178     if (!fcnAxis)
 179     {
 180       AliError("Empty axis");
 181       array[i].Set(2);
 182       array[i][0] = -1E5;
 183       array[i][1] = -1E5;
 184       continue;
 185     }
 186     else
 187     {
 188       array[i] = fcnAxis->GetArray();
 189     }
 190   }
 191
 192   // create histogram depending on the number of axes
 193   switch (fSize)
 194   {
 195     case 1:
 196       fH1 = new TH1F(histoName, histoTitle, array[0].GetSize() - 1, array[0].GetArray());
 197       break;
 198     case 2:
 199       fH1 = new TH2F(histoName, histoTitle, array[0].GetSize() - 1, array[0].GetArray(), array[1].GetSize() - 1, array[1].GetArray());
 200       break;
 201     case 3:
 202       fH1 = new TH3F(histoName, histoTitle, array[0].GetSize() - 1, array[0].GetArray(), array[1].GetSize() - 1, array[1].GetArray(), array[2].GetSize() - 1, array[2].GetArray());
 203       break;
 204   }
 205   fH1->Sumw2();
 206
 207   return fH1;
 208 }
 209
 210 //________________________________________________________________________________________
 211 THnSparseF* AliRsnFunction::CreateHistogramSparse(const char *histoName, const char *histoTitle)
 212 {
 213 //
 214 // Creates and returns the histogram defined using
 215 // arguments fo name and title, and the first histoDef for binning.
 216 // Variable-sized histogram binning is always called, due to use of histoDef,
 217 // even if the bins are equal, since they are defined in this class.
 218 // Eventually present histoDef's in other slots of array (1, 2) are ignored.
 219 //
 220 // This version produces a THnSparseF.
 221 //
 222
 223   fSize = fAxisList.GetEntries();
 224   if (!fSize)
 225   {
 226     AliError("No axes defined");
 227     return 0x0;
 228   }
 229
 230   // initialize the array of number of bins for each axis
 231   // taking it from the stored values, while for the bins
 232   // they are set as summied and defined later
 233   Double_t     dummyD;
 234   Int_t       *nbins   = new Int_t[fSize];
 235   AliRsnValue *fcnAxis = 0;
 236   for (Int_t i = 0; i < fSize; i++)
 237   {
 238     fcnAxis = (AliRsnValue*)fAxisList.At(i);
 239     if (!fcnAxis)
 240     {
 241       nbins[i] = 1;
 242       AliError("Empty axis");
 243       continue;
 244     }
 245     nbins[i] = fcnAxis->GetArray().GetSize() - 1;
 246   }
 247
 248   // create histogram
 249   fHSparse = new THnSparseF(histoName, histoTitle, fSize, nbins, &dummyD, &dummyD);
 250   fHSparse->Sumw2();
 251
 252   // update the various axes using the definitions given in the array of axes here
 253   for (Int_t i = 0; i < fSize; i++)
 254   {
 255     fcnAxis = (AliRsnValue*)fAxisList.At(i);
 256     if (!fcnAxis) {
 257       AliError("Empty axis: doing unique bin betweeen -100000 and 100000");
 258       continue;
 259     }
 260     fHSparse->SetBinEdges(i, fcnAxis->GetArray().GetArray());
 261   }
 262
 263   delete [] nbins;
 264
 265   return fHSparse;
 266 }
 267
 268
 269 //________________________________________________________________________________________
 270 Bool_t AliRsnFunction::Fill()
 271 {
 272 //
 273 // Fill function histogram with values computed from given input object.
 274 //
 275
 276   AliDebug(AliLog::kDebug +2,"->");
 277
 278   Int_t  i;
 279   Double_t *values = new Double_t[fSize];
 280   Bool_t    computeOK = kFALSE;
 281
 282   AliRsnValue *fcnAxis = 0;
 283   for (i = 0; i < fSize; i++)
 284   {
 285     values[i] = 0.0;
 286     fcnAxis = (AliRsnValue*)fAxisList.At(i);
 287     if (!fcnAxis) continue;
 288     switch (fcnAxis->GetTargetType())
 289     {
 290       case AliRsnTarget::kMother:
 291         fcnAxis->SetSupportObject(fPairDef);
 292         computeOK = fcnAxis->Eval(fPair);
 293         break;
 294       case AliRsnTarget::kEvent:
 295         computeOK = fcnAxis->Eval(fEvent);
 296         break;
 297       default:
 298         AliError(Form("Allowed targets are mothers and events; cannot use axis '%s' which has target '%s'", fcnAxis->GetName(), fcnAxis->GetTargetTypeName()));
 299         computeOK = kFALSE;
 300     }
 301     if (computeOK) values[i] = ((Float_t)fcnAxis->GetComputedValue());
 302   }
 303
 304   // fill histogram
 305   if (fUseTH1)
 306   {
 307     // check presence of output histogram
 308     if (!fH1)
 309     {
 310       AliError("Required a TH1 which is not initialized");
 311       delete [] values;
 312       return kFALSE;
 313     }
 314
 315     // fill according to dimensions
 316     switch (fSize)
 317     {
 318       case 1:
 319         {
 320           TH1F *h1 = (TH1F*)fH1;
 321           h1->Fill(values[0]);
 322         }
 323         break;
 324       case 2:
 325         {
 326           TH2F *h2 = (TH2F*)fH1;
 327           h2->Fill(values[0], values[1]);
 328         }
 329         break;
 330       case 3:
 331         {
 332           TH3F *h3 = (TH3F*)fH1;
 333           h3->Fill(values[0], values[1], values[2]);
 334         }
 335         break;
 336       default:
 337         AliError(Form("Wrong size : %d", fSize));
 338         delete [] values;
 339         return kFALSE;
 340     }
 341   }
 342   else
 343   {
 344     // check presence of output histogram
 345     if (!fHSparse)
 346     {
 347       AliError("Required a THnSparseF which is not initialized");
 348       delete [] values;
 349       return kFALSE;
 350     }
 351
 352     fHSparse->Fill(values);
 353   }
 354
 355   delete [] values;
 356
 357   AliDebug(AliLog::kDebug +2,"->");
 358   return kTRUE;
 359 }