Splitted the PID cuts into three classes for ITS, TPC and TOF, and added the possibil...
[u/mrichter/AliRoot.git] / 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 }