]> git.uio.no Git - u/mrichter/AliRoot.git/blobdiff - PWG2/RESONANCES/AliRsnFunction.cxx
Add new version of macros for RSN analysis
[u/mrichter/AliRoot.git] / PWG2 / RESONANCES / AliRsnFunction.cxx
index 73b835882402e5d048131cc1227b7891df23a67a..99064783c6abb0089b33d3f0c79806a0bcb68516 100644 (file)
 
 #include "AliLog.h"
 
-#include "AliRsnDaughter.h"
-#include "AliRsnEvent.h"
-#include "AliRsnPairDef.h"
-#include "AliRsnMother.h"
 #include "AliRsnValue.h"
 
 #include "AliRsnFunction.h"
@@ -36,15 +32,13 @@ ClassImp(AliRsnFunction)
 
 //________________________________________________________________________________________
 AliRsnFunction::AliRsnFunction(Bool_t useTH1) :
-    TObject(),
-    fPairDef(0x0),
-    fAxisList("AliRsnValue", 0),
-    fPair(0x0),
-    fEvent(0x0),
-    fUseTH1(useTH1),
-    fSize(0),
-    fH1(0x0),
-    fHSparse(0x0)
+   TObject(),
+   fAxisList("AliRsnValue", 0),
+   fUseTH1(useTH1),
+   fSize(0),
+   fH1(0x0),
+   fHSparse(0x0),
+   fValues(0)
 {
 //
 // Constructor.
@@ -53,15 +47,13 @@ AliRsnFunction::AliRsnFunction(Bool_t useTH1) :
 
 //________________________________________________________________________________________
 AliRsnFunction::AliRsnFunction(const AliRsnFunction &copy) :
-    TObject(copy),
-    fPairDef(copy.fPairDef),
-    fAxisList(copy.fAxisList),
-    fPair(copy.fPair),
-    fEvent(copy.fEvent),
-    fUseTH1(copy.fUseTH1),
-    fSize(copy.fSize),
-    fH1(0x0),
-    fHSparse(0x0)
+   TObject(copy),
+   fAxisList(copy.fAxisList),
+   fUseTH1(copy.fUseTH1),
+   fSize(copy.fSize),
+   fH1(0x0),
+   fHSparse(0x0),
+   fValues(copy.fValues)
 {
 //
 // Copy constructor.
@@ -75,22 +67,15 @@ const AliRsnFunction& AliRsnFunction::operator=(const AliRsnFunction& copy)
 // Assignment operator.
 //
 
-  //SetName(copy.GetName());
-  //SetTitle(copy.GetTitle());
+   fAxisList = copy.fAxisList;
+   fUseTH1 = copy.fUseTH1;
+   fSize = copy.fSize;
+   fValues = copy.fValues;
 
-  fPairDef = copy.fPairDef;
-  fPair = copy.fPair;
-  fEvent = copy.fEvent;
-  fUseTH1 = copy.fUseTH1;
-  fSize = copy.fSize;
+   if (fH1) delete fH1;
+   if (fHSparse) delete fHSparse;
 
-  if (fH1) delete fH1;
-  fH1 = 0x0;
-  
-  if (fHSparse) delete fHSparse;
-  fHSparse = 0x0;
-
-  return (*this);
+   return (*this);
 }
 
 //________________________________________________________________________________________
@@ -101,32 +86,43 @@ const char* AliRsnFunction::GetName() const
 // the function type and binning
 //
 
-  TString name("");
+   TString name("");
 
-  TObjArrayIter next(&fAxisList);
-  AliRsnValue *axis = 0;
+   TObjArrayIter next(&fAxisList);
+   AliRsnValue *axis = 0;
 
-  while ((axis = (AliRsnValue*)next())) {
-    if (name.Length() > 1) name += '_';
-    name += axis->GetName();
-  }
+   while ((axis = (AliRsnValue*)next())) {
+      if (name.Length() > 1) name += '_';
+      name += axis->GetName();
+   }
 
-  return name.Data();
+   return name.Data();
 }
 
 //________________________________________________________________________________________
-void AliRsnFunction::AddAxis(AliRsnValue *const axis)
+Bool_t AliRsnFunction::AddAxis(AliRsnValue *const axis)
 {
-  AliDebug(AliLog::kDebug+2,"<-");
-  Int_t size = fAxisList.GetEntries();
-  new(fAxisList[size]) AliRsnValue(*axis);
-  AliDebug(AliLog::kDebug+2,"->");
-  
-  if (fAxisList.GetEntries() > 3)
-  {
-    AliWarning("A TH1-type output cannot add more than 3 axes: switching to THnSparse -- THIS COULD CAUSE VERY LARGE FILES!!!");
-    fUseTH1 = kFALSE;
-  }
+//
+// Try to add a new axis to this function.
+// The operation succeeds only if the related value object
+// has a target compatible with the function type:
+// -- 'single'     functions, for tracks/events: AliRsnDaughter or AliRsnEvent
+// -- 'not single' functions, for pairs/events : AliRsnMother or AliRsnEvent
+// otherwise the axis is not added.
+//
+// If more than 3 axes are added, switch to THnSparseF output.
+// NOTE: this can cause large files and high memory occupancy.
+//
+
+   Int_t size = fAxisList.GetEntries();
+   new (fAxisList[size]) AliRsnValue(*axis);
+
+   if (fAxisList.GetEntries() > 3) {
+      AliWarning("Adding more than 3 axes will produce a THnSparseD output.");
+      fUseTH1 = kFALSE;
+   }
+
+   return kTRUE;
 }
 
 //________________________________________________________________________________________
@@ -142,53 +138,49 @@ TH1* AliRsnFunction::CreateHistogram(const char *histoName, const char *histoTit
 // This version produces a THnSparseF.
 //
 
-  fSize = fAxisList.GetEntries();
-  if (!fSize) {
-    AliError("No axes defined");
-    return 0x0;
-  }
-  else if (fSize < 1 || fSize > 3)
-  {
-    AliError("Too few or too many axes defined");
-    return 0x0;
-  }
-
-  // retrieve binnings for main and secondary axes
-  AliRsnValue *fcnAxis;
-  TArrayD      array[3];
-  for (Int_t i = 0; i < fSize; i++) 
-  {
-    fcnAxis = (AliRsnValue*)fAxisList.At(i);
-    if (!fcnAxis) 
-    {
-      AliError("Empty axis");
-      array[i].Set(2);
-      array[i][0] = -1E5;
-      array[i][1] = -1E5;
-      continue;
-    }
-    else
-    {
-      array[i] = fcnAxis->GetArray();
-    }
-  }
-
-  // create histogram depending on the number of axes
-  switch (fSize)
-  {
-    case 1:
-      fH1 = new TH1F(histoName, histoTitle, array[0].GetSize() - 1, array[0].GetArray());
-      break;
-    case 2:
-      fH1 = new TH2F(histoName, histoTitle, array[0].GetSize() - 1, array[0].GetArray(), array[1].GetSize() - 1, array[1].GetArray());
-      break;
-    case 3:
-      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());
-      break;
-  }
-  fH1->Sumw2();
-
-  return fH1;
+   fSize = fAxisList.GetEntries();
+   if (!fSize) {
+      AliError("No axes defined");
+      return 0x0;
+   } else if (fSize > 3) {
+      AliError("Too many axes defined for a TH1 object");
+      return 0x0;
+   }
+   
+   // initialize the size of values container
+   fValues.Set(fSize);
+
+   // retrieve binnings for main and secondary axes
+   AliRsnValue *fcnAxis;
+   TArrayD      array[3];
+   for (Int_t i = 0; i < fSize; i++) {
+      fcnAxis = (AliRsnValue*)fAxisList.At(i);
+      if (!fcnAxis) {
+         AliError("Empty axis");
+         array[i].Set(2);
+         array[i][0] = -1E5;
+         array[i][1] = -1E5;
+         continue;
+      } else {
+         array[i] = fcnAxis->GetArray();
+      }
+   }
+
+   // create histogram depending on the number of axes
+   switch (fSize) {
+      case 1:
+         fH1 = new TH1F(histoName, histoTitle, array[0].GetSize() - 1, array[0].GetArray());
+         break;
+      case 2:
+         fH1 = new TH2F(histoName, histoTitle, array[0].GetSize() - 1, array[0].GetArray(), array[1].GetSize() - 1, array[1].GetArray());
+         break;
+      case 3:
+         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());
+         break;
+   }
+   fH1->Sumw2();
+
+   return fH1;
 }
 
 //________________________________________________________________________________________
@@ -204,109 +196,121 @@ THnSparseF* AliRsnFunction::CreateHistogramSparse(const char *histoName, const c
 // This version produces a THnSparseF.
 //
 
-  fSize = fAxisList.GetEntries();
-  if (!fSize) {
-    AliError("No axes defined");
-    return 0x0;
-  }
-  
-  // create dummy variables to initialize the histogram
-  Int_t    dummyI;
-  Double_t dummyD;
-
-  // create histogram
-  fHSparse = new THnSparseF(histoName, histoTitle, fSize, &dummyI, &dummyD, &dummyD);
-  fHSparse->Sumw2();
-  
-  // update the various axes using the definitions given in the array of axes here
-  AliRsnValue *fcnAxis = 0;
-  TAxis       *axis = 0;
-  for (Int_t i = 0; i < fSize; i++) 
-  {
-    fcnAxis = (AliRsnValue*)fAxisList.At(i);
-    axis    = fHSparse->GetAxis(i);
-    if (!fcnAxis) {
-      axis->Set(1, -1E5, 1E5);
-      AliError("Empty axis: doing unique bin betweeen -100000 and 100000");
-      continue;
-    }
-    axis->Set(fcnAxis->GetArray().GetSize() - 1, fcnAxis->GetArray().GetArray());
-  }
-
-  return fHSparse;
+   fSize = fAxisList.GetEntries();
+   if (!fSize) {
+      AliError("No axes defined");
+      return 0x0;
+   }
+   
+   // initialize the size of values container
+   fValues.Set(fSize);
+
+   // initialize the array of number of bins for each axis
+   // taking it from the stored values, while for the bins
+   // they are set as summied and defined later
+   Int_t       *nbins   = new Int_t[fSize];
+   AliRsnValue *fcnAxis = 0;
+   for (Int_t i = 0; i < fSize; i++) {
+      fcnAxis = (AliRsnValue*)fAxisList.At(i);
+      if (!fcnAxis) {
+         nbins[i] = 1;
+         AliError("Empty axis");
+         continue;
+      }
+      nbins[i] = fcnAxis->GetArray().GetSize() - 1;
+   }
+
+   // create histogram
+   fHSparse = new THnSparseF(histoName, histoTitle, fSize, nbins);
+   fHSparse->Sumw2();
+
+   // update the various axes using the definitions given in the array of axes here
+   for (Int_t i = 0; i < fSize; i++) {
+      fcnAxis = (AliRsnValue*)fAxisList.At(i);
+      if (!fcnAxis) {
+         AliError("Empty axis: doing unique bin betweeen -100000 and 100000");
+         continue;
+      }
+      TAxis* axis = fHSparse->GetAxis(i);
+      axis->Set(nbins[i], fcnAxis->GetArray().GetArray());
+   }
+
+   delete [] nbins;
+
+   return fHSparse;
 }
 
 
 //________________________________________________________________________________________
-Bool_t AliRsnFunction::Fill()
+Bool_t AliRsnFunction::Fill(TObject *object)
 {
 //
-// Fill function histogram with values computed from given input object.
+// Fill function histogram using the passed object
 //
 
-  AliDebug(AliLog::kDebug +2,"->");
-
-  Int_t  i;
-  Double_t *values = new Double_t[fSize];
-
-  AliRsnValue *fcnAxis = 0;
-  for (i = 0; i < fSize; i++) {
-    fcnAxis = (AliRsnValue*)fAxisList.At(i);
-    if (!fcnAxis) {
-      values[i] = 0.0;
-      continue;
-    }
-    if (fcnAxis->Eval(fPair, fPairDef, fEvent)) values[i] = fcnAxis->GetValue();
-  }
-  
-  // fill histogram
-  if (fUseTH1)
-  {
-    // check presence of output histogram
-    if (!fH1) {
-      AliError("Required a TH1 whish is not initialized");
-      return kFALSE;
-    }
-    
-    // fill according to dimensions
-    switch (fSize)
-    {
-      case 1:
-        {
-          TH1F *h1 = (TH1F*)fH1;
-          h1->Fill(values[0]);
-        }
-        break;
-      case 2:
-        {
-          TH2F *h2 = (TH2F*)fH1;
-          h2->Fill(values[0], values[1]);
-        }
-        break;
-      case 3:
-        {
-          TH3F *h3 = (TH3F*)fH1;
-          h3->Fill(values[0], values[1], values[2]);
-        }
-        break;
-      default:
-        AliError(Form("Wrong size : %d", fSize));
-        return kFALSE;
-    }
-  }
-  else
-  {
-    // check presence of output histogram
-    if (!fHSparse) {
-      AliError("Required a THnSparseF which is not initialized");
-      return kFALSE;
-    }
-    
-    fHSparse->Fill(values);
-  }
-  
-  delete [] values;
-
-  AliDebug(AliLog::kDebug +2,"->");
-  return kTRUE;
+   AliDebug(AliLog::kDebug + 2, "->");
+
+   // loop on axes and try to compute values
+   // using this object or, as an alternative
+   // its reference event
+   Int_t  i;
+   Bool_t globalOK = kTRUE, computeOK;
+   AliRsnValue *fcnAxis = 0;
+   for (i = 0; i < fSize; i++) {
+      fValues[i] = 0.0;
+      computeOK = kFALSE;
+      fcnAxis = (AliRsnValue*)fAxisList.At(i);
+      if (fcnAxis) {
+         computeOK = fcnAxis->Eval(object);
+         if (computeOK) fValues[i] = ((Float_t)fcnAxis->GetComputedValue());
+      }
+      if (!computeOK) globalOK = kFALSE;
+   }
+   
+   // if even one of the computations has failes, the histograms are not filled
+   if (!globalOK) return kFALSE;
+
+   // fill histogram
+   if (fUseTH1) {
+      
+      // check presence of output histogram
+      if (!fH1) {
+         AliError("Required a TH1 which is not initialized");
+         return kFALSE;
+      }
+
+      // fill according to dimensions
+      switch (fSize) {
+         case 1: {
+            TH1F *h1 = (TH1F*)fH1;
+            h1->Fill(fValues[0]);
+         }
+         break;
+         case 2: {
+            TH2F *h2 = (TH2F*)fH1;
+            h2->Fill(fValues[0], fValues[1]);
+         }
+         break;
+         case 3: {
+            TH3F *h3 = (TH3F*)fH1;
+            h3->Fill(fValues[0], fValues[1], fValues[2]);
+         }
+         break;
+         default:
+            AliError(Form("Wrong size : %d", fSize));
+            return kFALSE;
+      }
+   } else {
+      
+      // check presence of output histogram
+      if (!fHSparse) {
+         AliError("Required a THnSparseF which is not initialized");
+         return kFALSE;
+      }
+
+      fHSparse->Fill(fValues.GetArray());
+   }
+
+   AliDebug(AliLog::kDebug + 2, "->");
+   return kTRUE;
 }