]> git.uio.no Git - u/mrichter/AliRoot.git/blame - PWGLF/RESONANCES/AliRsnMiniOutput.cxx
git://git.uio.no / u / mrichter / AliRoot.git / blame
? search:
summary | shortlog | log | commit | commitdiff | tree
blob | blame (incremental) | history | HEAD
Fix for coverity 24399, 24400, 24401
[u/mrichter/AliRoot.git] / PWGLF / RESONANCES / AliRsnMiniOutput.cxx
CommitLineData
03d23846 1//
2// Mini-Output
3// All the definitions needed for building a RSN histogram
4// including:
5// -- properties of resonance (mass, PDG code if needed)
6// -- properties of daughters (assigned mass, charges)
7// -- definition of output histogram
61f275d1 8//
03d23846 9
10#include "Riostream.h"
11
12#include "TH1.h"
13#include "TH2.h"
14#include "TH3.h"
15#include "TList.h"
45aa62b9 16#include "TMath.h"
03d23846 17#include "THnSparse.h"
18#include "TString.h"
19#include "TClonesArray.h"
20
21#include "AliRsnMiniParticle.h"
22#include "AliRsnMiniPair.h"
23#include "AliRsnMiniEvent.h"
088ca370 24#include "AliAODEvent.h"
03d23846 25
26#include "AliLog.h"
27#include "AliRsnCutSet.h"
28#include "AliRsnMiniAxis.h"
29#include "AliRsnMiniOutput.h"
30#include "AliRsnMiniValue.h"
31
32ClassImp(AliRsnMiniOutput)
33
34//__________________________________________________________________________________________________
35AliRsnMiniOutput::AliRsnMiniOutput() :
36 TNamed(),
37 fOutputType(kTypes),
38 fComputation(kComputations),
39 fMotherPDG(0),
40 fMotherMass(0.0),
41 fPairCuts(0x0),
42 fOutputID(-1),
43 fAxes("AliRsnMiniAxis", 0),
44 fComputed(0),
45 fPair(),
a2455d2a 46 fList(0x0),
47 fSel1(0),
74d60285 48 fSel2(0),
088ca370 49 fMaxNSisters(-1),
50 fCheckP(kFALSE),
31dbef4e 51 fCheckFeedDown(kFALSE),
52 fOriginDselection(kFALSE),
53 fKeepDfromB(kFALSE),
54 fKeepDfromBOnly(kFALSE),
55 fRejectIfNoQuark(kFALSE),
74d60285 56 fCheckHistRange(kTRUE)
03d23846 57{
58//
59// Constructor
60//
61
62 fCutID[0] = fCutID[1] = -1;
63 fDaughter[0] = fDaughter[1] = AliRsnDaughter::kUnknown;
64 fCharge[0] = fCharge[1] = 0;
65}
66
67//__________________________________________________________________________________________________
68AliRsnMiniOutput::AliRsnMiniOutput(const char *name, EOutputType type, EComputation src) :
69 TNamed(name, ""),
70 fOutputType(type),
71 fComputation(src),
72 fMotherPDG(0),
73 fMotherMass(0.0),
74 fPairCuts(0x0),
75 fOutputID(-1),
76 fAxes("AliRsnMiniAxis", 0),
77 fComputed(0),
78 fPair(),
a2455d2a 79 fList(0x0),
80 fSel1(0),
74d60285 81 fSel2(0),
088ca370 82 fMaxNSisters(-1),
83 fCheckP(kFALSE),
31dbef4e 84 fCheckFeedDown(kFALSE),
85 fOriginDselection(kFALSE),
86 fKeepDfromB(kFALSE),
87 fKeepDfromBOnly(kFALSE),
88 fRejectIfNoQuark(kFALSE),
74d60285 89 fCheckHistRange(kTRUE)
03d23846 90{
91//
92// Constructor
93//
94
95 fCutID[0] = fCutID[1] = -1;
96 fDaughter[0] = fDaughter[1] = AliRsnDaughter::kUnknown;
97 fCharge[0] = fCharge[1] = 0;
98}
99
100//__________________________________________________________________________________________________
101AliRsnMiniOutput::AliRsnMiniOutput(const char *name, const char *outType, const char *compType) :
102 TNamed(name, ""),
103 fOutputType(kTypes),
104 fComputation(kComputations),
105 fMotherPDG(0),
106 fMotherMass(0.0),
107 fPairCuts(0x0),
108 fOutputID(-1),
109 fAxes("AliRsnMiniAxis", 0),
110 fComputed(0),
111 fPair(),
a2455d2a 112 fList(0x0),
113 fSel1(0),
74d60285 114 fSel2(0),
088ca370 115 fMaxNSisters(-1),
116 fCheckP(kFALSE),
31dbef4e 117 fCheckFeedDown(kFALSE),
118 fOriginDselection(kFALSE),
119 fKeepDfromB(kFALSE),
120 fKeepDfromBOnly(kFALSE),
121 fRejectIfNoQuark(kFALSE),
74d60285 122 fCheckHistRange(kTRUE)
03d23846 123{
124//
125// Constructor, with a more user friendly implementation, where
126// the user sets the type of output and computations through conventional strings:
127//
128// Output:
6aaeb33c 129// -- "HIST" --> common histogram (up to 3 dimensions)
130// -- "SPARSE" --> sparse histogram
61f275d1 131//
132// Computation:
6aaeb33c 133// -- "EVENT" --> event-only computations
134// -- "PAIR" --> track pair computations (default)
135// -- "MIX" --> event mixing (like track pair, but different events)
136// -- "ROTATE1" --> rotated background (rotate first track)
137// -- "ROTATE2" --> rotated background (rotate second track)
138// -- "TRUE" --> true pairs (like track pair, but checking that come from same mother)
139// -- "MOTHER" --> mother (loop on MC directly for mothers --> denominator of efficiency)
03d23846 140//
141
142 TString input;
61f275d1 143
03d23846 144 // understand output type
145 input = outType;
146 input.ToUpper();
147 if (!input.CompareTo("HIST"))
148 fOutputType = kHistogram;
149 else if (!input.CompareTo("SPARSE"))
150 fOutputType = kHistogramSparse;
151 else
152 AliWarning(Form("String '%s' does not define a meaningful output type", outType));
61f275d1 153
03d23846 154 // understand computation type
155 input = compType;
156 input.ToUpper();
157 if (!input.CompareTo("EVENT"))
158 fComputation = kEventOnly;
159 else if (!input.CompareTo("PAIR"))
160 fComputation = kTrackPair;
161 else if (!input.CompareTo("MIX"))
162 fComputation = kTrackPairMix;
6aaeb33c 163 else if (!input.CompareTo("ROTATE1"))
164 fComputation = kTrackPairRotated1;
165 else if (!input.CompareTo("ROTATE2"))
166 fComputation = kTrackPairRotated2;
03d23846 167 else if (!input.CompareTo("TRUE"))
168 fComputation = kTruePair;
169 else if (!input.CompareTo("MOTHER"))
170 fComputation = kMother;
171 else
172 AliWarning(Form("String '%s' does not define a meaningful computation type", compType));
61f275d1 173
03d23846 174 fCutID[0] = fCutID[1] = -1;
175 fDaughter[0] = fDaughter[1] = AliRsnDaughter::kUnknown;
176 fCharge[0] = fCharge[1] = 0;
177}
178
179//__________________________________________________________________________________________________
180AliRsnMiniOutput::AliRsnMiniOutput(const AliRsnMiniOutput &copy) :
181 TNamed(copy),
182 fOutputType(copy.fOutputType),
183 fComputation(copy.fComputation),
184 fMotherPDG(copy.fMotherPDG),
185 fMotherMass(copy.fMotherMass),
186 fPairCuts(copy.fPairCuts),
187 fOutputID(copy.fOutputID),
188 fAxes(copy.fAxes),
189 fComputed(copy.fComputed),
190 fPair(),
a2455d2a 191 fList(copy.fList),
192 fSel1(0),
74d60285 193 fSel2(0),
088ca370 194 fMaxNSisters(-1),
195 fCheckP(kFALSE),
31dbef4e 196 fCheckFeedDown(kFALSE),
197 fOriginDselection(kFALSE),
198 fKeepDfromB(kFALSE),
199 fKeepDfromBOnly(kFALSE),
200 fRejectIfNoQuark(kFALSE),
74d60285 201 fCheckHistRange(copy.fCheckHistRange)
03d23846 202{
203//
204// Copy constructor
205//
206
207 Int_t i;
208 for (i = 0; i < 2; i++) {
209 fCutID[i] = copy.fCutID[i];
210 fDaughter[i] = copy.fDaughter[i];
211 fCharge[i] = copy.fCharge[i];
212 }
213}
214
215//__________________________________________________________________________________________________
61f275d1 216AliRsnMiniOutput &AliRsnMiniOutput::operator=(const AliRsnMiniOutput &copy)
03d23846 217{
218//
219// Assignment operator
220//
61f275d1 221 if (this == &copy)
222 return *this;
223 fOutputType = copy.fOutputType;
224 fComputation = copy.fComputation;
225 fMotherPDG = copy.fMotherPDG;
226 fMotherMass = copy.fMotherMass;
227 fPairCuts = copy.fPairCuts;
228 fOutputID = copy.fOutputID;
229 fAxes = copy.fAxes;
230 fComputed = copy.fComputed;
231 fList = copy.fList;
232
233 Int_t i;
234 for (i = 0; i < 2; i++) {
235 fCutID[i] = copy.fCutID[i];
236 fDaughter[i] = copy.fDaughter[i];
237 fCharge[i] = copy.fCharge[i];
238 }
239
240 fSel1.Set(0);
241 fSel2.Set(0);
088ca370 242 fMaxNSisters = copy.fMaxNSisters;
243 fCheckP = copy.fCheckP;
31dbef4e 244 fCheckFeedDown = copy.fCheckFeedDown;
245 fOriginDselection = copy.fOriginDselection;
246 fKeepDfromB = copy.fOriginDselection;
247 fKeepDfromBOnly = copy.fKeepDfromBOnly;
248 fRejectIfNoQuark = copy.fRejectIfNoQuark;
74d60285 249 fCheckHistRange = copy.fCheckHistRange;
61f275d1 250
03d23846 251 return (*this);
252}
253
254
255//__________________________________________________________________________________________________
256void AliRsnMiniOutput::AddAxis(Int_t i, Int_t nbins, Double_t min, Double_t max)
257{
258//
259// Create a new axis reference
260//
261
61f275d1 262 Int_t size = fAxes.GetEntries();
03d23846 263 new (fAxes[size]) AliRsnMiniAxis(i, nbins, min, max);
264}
265
266//__________________________________________________________________________________________________
267void AliRsnMiniOutput::AddAxis(Int_t i, Double_t min, Double_t max, Double_t step)
268{
269//
270// Create a new axis reference
271//
272
61f275d1 273 Int_t size = fAxes.GetEntries();
03d23846 274 new (fAxes[size]) AliRsnMiniAxis(i, min, max, step);
275}
276
277//__________________________________________________________________________________________________
278void AliRsnMiniOutput::AddAxis(Int_t i, Int_t nbins, Double_t *values)
279{
280//
281// Create a new axis reference
282//
283
61f275d1 284 Int_t size = fAxes.GetEntries();
03d23846 285 new (fAxes[size]) AliRsnMiniAxis(i, nbins, values);
286}
287
288//__________________________________________________________________________________________________
289Bool_t AliRsnMiniOutput::Init(const char *prefix, TList *list)
290{
291//
292// Initialize properly the histogram and add it to the argument list
293//
294
295 if (!list) {
296 AliError("Required an output list");
297 return kFALSE;
298 }
61f275d1 299
03d23846 300 fList = list;
45aa62b9 301 Int_t size = fAxes.GetEntries();
302 if (size < 1) {
303 AliWarning(Form("[%s] Cannot initialize histogram with less than 1 axis", GetName()));
304 return kFALSE;
305 }
03d23846 306
307 switch (fOutputType) {
308 case kHistogram:
d573d2fb 309 if (size <= 3) {
310 CreateHistogram(Form("%s_%s", prefix, GetName()));
311 } else {
312 AliInfo(Form("[%s] Added %d > 3 axes. Creating a sparse histogram", GetName(), size));
313 fOutputType = kHistogramSparse;
314 CreateHistogramSparse(Form("%s_%s", prefix, GetName()));
315 }
03d23846 316 return kTRUE;
317 case kHistogramSparse:
6aaeb33c 318 CreateHistogramSparse(Form("%s_%s", prefix, GetName()));
03d23846 319 return kTRUE;
320 default:
321 AliError("Wrong output histogram definition");
322 return kFALSE;
323 }
324}
325
326//__________________________________________________________________________________________________
327void AliRsnMiniOutput::CreateHistogram(const char *name)
328{
329//
330// Initialize the 'default' TH1 output object.
331// In case one of the expected axes is NULL, the initialization fails.
332//
333
334 Int_t size = fAxes.GetEntries();
335 AliInfo(Form("Histogram name = '%s', with %d axes", name, size));
336
337 // we expect to have maximum 3 axes in this case
338 AliRsnMiniAxis *xAxis = 0x0, *yAxis = 0x0, *zAxis = 0x0;
61f275d1 339 if (size >= 1) xAxis = (AliRsnMiniAxis *)fAxes[0];
340 if (size >= 2) yAxis = (AliRsnMiniAxis *)fAxes[1];
341 if (size >= 3) zAxis = (AliRsnMiniAxis *)fAxes[2];
342
03d23846 343 // create histogram depending on the number of axes
344 TH1 *h1 = 0x0;
345 if (xAxis && yAxis && zAxis) {
346 h1 = new TH3F(name, "", xAxis->NBins(), xAxis->BinArray(), yAxis->NBins(), yAxis->BinArray(), zAxis->NBins(), zAxis->BinArray());
347 } else if (xAxis && yAxis) {
348 h1 = new TH2F(name, "", xAxis->NBins(), xAxis->BinArray(), yAxis->NBins(), yAxis->BinArray());
349 } else if (xAxis) {
350 h1 = new TH1F(name, "", xAxis->NBins(), xAxis->BinArray());
351 } else {
352 AliError("No axis was initialized");
353 return;
354 }
61f275d1 355
03d23846 356 // switch the correct computation of errors
357 if (h1 && fList) {
358 h1->Sumw2();
359 fList->Add(h1);
360 fOutputID = fList->IndexOf(h1);
361 }
362}
363
364//________________________________________________________________________________________
365void AliRsnMiniOutput::CreateHistogramSparse(const char *name)
366{
367//
368// Initialize the THnSparse output object.
369// In case one of the expected axes is NULL, the initialization fails.
370//
371
d573d2fb 372 Int_t size = fAxes.GetEntries();
373 AliInfo(Form("Sparse histogram name = '%s', with %d axes", name, size));
61f275d1 374
03d23846 375 // retrieve binnings and sizes of all axes
376 // since the check for null values is done in Init(),
377 // we assume that here they must all be well defined
03d23846 378 Int_t i, *nbins = new Int_t[size];
379 for (i = 0; i < size; i++) {
61f275d1 380 AliRsnMiniAxis *axis = (AliRsnMiniAxis *)fAxes[i];
03d23846 381 nbins[i] = axis->NBins();
382 }
383
384 // create fHSparseogram
385 THnSparseF *h1 = new THnSparseF(name, "", size, nbins);
386
387 // update the various axes using the definitions given in the array of axes here
388 for (i = 0; i < size; i++) {
61f275d1 389 AliRsnMiniAxis *axis = (AliRsnMiniAxis *)fAxes[i];
03d23846 390 h1->GetAxis(i)->Set(nbins[i], axis->BinArray());
391 }
392
393 // clear heap
394 delete [] nbins;
61f275d1 395
03d23846 396 // add to list
397 if (h1 && fList) {
398 h1->Sumw2();
399 fList->Add(h1);
400 fOutputID = fList->IndexOf(h1);
401 }
402}
403
404//________________________________________________________________________________________
45aa62b9 405Bool_t AliRsnMiniOutput::FillEvent(AliRsnMiniEvent *event, TClonesArray *valueList)
03d23846 406{
407//
408// Compute values for event-based computations (does not use the pair)
409//
410
411 // check computation type
412 if (fComputation != kEventOnly) {
413 AliError("This method can be called only for event-based computations");
414 return kFALSE;
415 }
416
45aa62b9 417 // compute & fill
418 ComputeValues(event, valueList);
419 FillHistogram();
420 return kTRUE;
03d23846 421}
422
423//________________________________________________________________________________________
45aa62b9 424Bool_t AliRsnMiniOutput::FillMother(const AliRsnMiniPair *pair, AliRsnMiniEvent *event, TClonesArray *valueList)
03d23846 425{
426//
427// Compute values for mother-based computations
428//
429
430 // check computation type
431 if (fComputation != kMother) {
432 AliError("This method can be called only for mother-based computations");
433 return kFALSE;
434 }
61f275d1 435
03d23846 436 // copy passed pair info
437 fPair = (*pair);
61f275d1 438
03d23846 439 // check pair against cuts
440 if (fPairCuts) if (!fPairCuts->IsSelected(&fPair)) return kFALSE;
441
45aa62b9 442 // compute & fill
443 ComputeValues(event, valueList);
444 FillHistogram();
445 return kTRUE;
03d23846 446}
447
448//________________________________________________________________________________________
d573d2fb 449Int_t AliRsnMiniOutput::FillPair(AliRsnMiniEvent *event1, AliRsnMiniEvent *event2, TClonesArray *valueList, Bool_t refFirst)
03d23846 450{
451//
45aa62b9 452// Loops on the passed mini-event, and for each pair of particles
453// which satisfy the charge and cut requirements defined here, add an entry.
454// Returns the number of successful fillings.
d573d2fb 455// Last argument tells if the reference event for event-based values is the first or the second.
03d23846 456//
457
458 // check computation type
6aaeb33c 459 Bool_t okComp = kFALSE;
45aa62b9 460 if (fComputation == kTrackPair) okComp = kTRUE;
461 if (fComputation == kTrackPairMix) okComp = kTRUE;
462 if (fComputation == kTrackPairRotated1) okComp = kTRUE;
463 if (fComputation == kTrackPairRotated2) okComp = kTRUE;
464 if (fComputation == kTruePair) okComp = kTRUE;
6aaeb33c 465 if (!okComp) {
45aa62b9 466 AliError(Form("[%s] This method can be called only for pair-based computations", GetName()));
03d23846 467 return kFALSE;
468 }
61f275d1 469
45aa62b9 470 // loop variables
471 Int_t i1, i2, start, nadded = 0;
45aa62b9 472 AliRsnMiniParticle *p1, *p2;
61f275d1 473
45aa62b9 474 // it is necessary to know if criteria for the two daughters are the same
475 // and if the two events are the same or not (mixing)
7196ee4f 476 //Bool_t sameCriteria = ((fCharge[0] == fCharge[1]) && (fCutID[0] == fCutID[1]));
477 Bool_t sameCriteria = ((fCharge[0] == fCharge[1]) && (fDaughter[0] == fDaughter[1]));
45aa62b9 478 Bool_t sameEvent = (event1->ID() == event2->ID());
61f275d1 479
9e3a9020 480 TString selList1 = "";
481 TString selList2 = "";
a2455d2a 482 Int_t n1 = event1->CountParticles(fSel1, fCharge[0], fCutID[0]);
483 Int_t n2 = event2->CountParticles(fSel2, fCharge[1], fCutID[1]);
484 for (i1 = 0; i1 < n1; i1++) selList1.Append(Form("%d ", fSel1[i1]));
485 for (i2 = 0; i2 < n2; i2++) selList2.Append(Form("%d ", fSel2[i2]));
486 AliDebugClass(1, Form("[%10s] Part #1: [%s] -- evID %6d -- charge = %c -- cut ID = %d --> %4d tracks (%s)", GetName(), (event1 == event2 ? "def" : "mix"), event1->ID(), fCharge[0], fCutID[0], n1, selList1.Data()));
487 AliDebugClass(1, Form("[%10s] Part #2: [%s] -- evID %6d -- charge = %c -- cut ID = %d --> %4d tracks (%s)", GetName(), (event1 == event2 ? "def" : "mix"), event2->ID(), fCharge[1], fCutID[1], n2, selList2.Data()));
9e3a9020 488 if (!n1 || !n2) {
d573d2fb 489 AliDebugClass(1, "No pairs to mix");
490 return 0;
491 }
61f275d1 492
45aa62b9 493 // external loop
494 for (i1 = 0; i1 < n1; i1++) {
a2455d2a 495 p1 = event1->GetParticle(fSel1[i1]);
9e3a9020 496 //p1 = event1->GetParticle(i1);
497 //if (p1->Charge() != fCharge[0]) continue;
498 //if (!p1->HasCutBit(fCutID[0])) continue;
45aa62b9 499 // define starting point for inner loop
500 // if daughter selection criteria (charge, cuts) are the same
501 // and the two events coincide, internal loop must start from
502 // the first track *after* current one;
503 // otherwise it starts from the beginning
504 start = ((sameEvent && sameCriteria) ? i1 + 1 : 0);
7196ee4f 505 AliDebugClass(2, Form("Start point = %d", start));
45aa62b9 506 // internal loop
507 for (i2 = start; i2 < n2; i2++) {
a2455d2a 508 p2 = event2->GetParticle(fSel2[i2]);
9e3a9020 509 //p2 = event2->GetParticle(i2);
510 //if (p2->Charge() != fCharge[1]) continue;
511 //if (!p2->HasCutBit(fCutID[1])) continue;
45aa62b9 512 // avoid to mix a particle with itself
513 if (sameEvent && (p1->Index() == p2->Index())) {
d573d2fb 514 AliDebugClass(2, "Skipping same index");
45aa62b9 515 continue;
516 }
517 // sum momenta
518 fPair.Fill(p1, p2, GetMass(0), GetMass(1), fMotherMass);
519 // do rotation if needed
520 if (fComputation == kTrackPairRotated1) fPair.InvertP(kTRUE);
521 if (fComputation == kTrackPairRotated2) fPair.InvertP(kFALSE);
522 // if required, check that this is a true pair
523 if (fComputation == kTruePair) {
524 if (fPair.Mother() < 0) {
525 continue;
55de729a 526 } else if (fPair.MotherPDG() != fMotherPDG) {
45aa62b9 527 continue;
528 }
529 Bool_t decayMatch = kFALSE;
530 if (p1->PDGAbs() == AliRsnDaughter::SpeciesPDG(fDaughter[0]) && p2->PDGAbs() == AliRsnDaughter::SpeciesPDG(fDaughter[1]))
531 decayMatch = kTRUE;
532 if (p2->PDGAbs() == AliRsnDaughter::SpeciesPDG(fDaughter[0]) && p1->PDGAbs() == AliRsnDaughter::SpeciesPDG(fDaughter[1]))
533 decayMatch = kTRUE;
534 if (!decayMatch) continue;
35b36c39 535 if ( (fMaxNSisters>0) && (p1->NTotSisters()==p2->NTotSisters()) && (p1->NTotSisters()>fMaxNSisters)) continue;
536 if ( fCheckP &&(TMath::Abs(fPair.PmotherX()-(p1->Px(1)+p2->Px(1)))/(TMath::Abs(fPair.PmotherX())+1.e-13)) > 0.00001 &&
537 (TMath::Abs(fPair.PmotherY()-(p1->Py(1)+p2->Py(1)))/(TMath::Abs(fPair.PmotherY())+1.e-13)) > 0.00001 &&
538 (TMath::Abs(fPair.PmotherZ()-(p1->Pz(1)+p2->Pz(1)))/(TMath::Abs(fPair.PmotherZ())+1.e-13)) > 0.00001 ) continue;
31dbef4e 539 if ( fCheckFeedDown ){
540 Int_t pdgGranma = 0;
541 Bool_t isFromB=kFALSE;
542 Bool_t isQuarkFound=kFALSE;
543
544 if(fPair.IsFromB() == kTRUE) isFromB = kTRUE;
545 if(fPair.IsQuarkFound() == kTRUE) isQuarkFound = kTRUE;
546 if(fRejectIfNoQuark && !isQuarkFound) pdgGranma = -99999;
547 if(isFromB){
548 if (!fKeepDfromB) pdgGranma = -9999; //skip particle if come from a B meson.
549 }
550 else{
551 if (fKeepDfromBOnly) pdgGranma = -999;
552 }
553 if (pdgGranma == -99999){
554 AliDebug(2,"This particle does not have a quark in his genealogy\n");
555 continue;
556 }
557 if (pdgGranma == -9999){
558 AliDebug(2,"This particle come from a B decay channel but according to the settings of the task, we keep only the prompt charm particles\n");
559 continue;
560 }
561
562 if (pdgGranma == -999){
563 AliDebug(2,"This particle come from a prompt charm particles but according to the settings of the task, we want only the ones coming from B\n");
564 continue;
565 }
566 }
45aa62b9 567 }
568 // check pair against cuts
569 if (fPairCuts) {
9e7b94f5 570 if (!fPairCuts->IsSelected(&fPair)) continue;
45aa62b9 571 }
572 // get computed values & fill histogram
9e3a9020 573 nadded++;
61f275d1 574 if (refFirst) ComputeValues(event1, valueList); else ComputeValues(event2, valueList);
45aa62b9 575 FillHistogram();
45aa62b9 576 } // end internal loop
577 } // end external loop
61f275d1 578
d573d2fb 579 AliDebugClass(1, Form("Pairs added in total = %4d", nadded));
45aa62b9 580 return nadded;
03d23846 581}
31dbef4e 582//___________________________________________________________
583void AliRsnMiniOutput::SetDselection(UShort_t originDselection)
584{
585 // setting the way the D0 will be selected
586 // 0 --> only from c quarks
587 // 1 --> only from b quarks
588 // 2 --> from both c quarks and b quarks
589
590 fOriginDselection = originDselection;
591
592 if (fOriginDselection == 0) {
593 fKeepDfromB = kFALSE;
594 fKeepDfromBOnly = kFALSE;
595 }
596
597 if (fOriginDselection == 1) {
598 fKeepDfromB = kTRUE;
599 fKeepDfromBOnly = kTRUE;
600 }
601
602 if (fOriginDselection == 2) {
603 fKeepDfromB = kTRUE;
604 fKeepDfromBOnly = kFALSE;
605 }
606
607 return;
608}
03d23846 609//________________________________________________________________________________________
45aa62b9 610void AliRsnMiniOutput::ComputeValues(AliRsnMiniEvent *event, TClonesArray *valueList)
03d23846 611{
612//
613// Using the arguments and the internal 'fPair' data member,
614// compute all values to be stored in the histogram
615//
616
617 // check size of computed array
618 Int_t size = fAxes.GetEntries();
619 if (fComputed.GetSize() != size) fComputed.Set(size);
620
621 Int_t i, ival, nval = valueList->GetEntries();
61f275d1 622
03d23846 623 for (i = 0; i < size; i++) {
624 fComputed[i] = 1E20;
61f275d1 625 AliRsnMiniAxis *axis = (AliRsnMiniAxis *)fAxes[i];
03d23846 626 if (!axis) {
627 AliError("Null axis");
628 continue;
629 }
630 ival = axis->GetValueID();
631 if (ival < 0 || ival >= nval) {
632 AliError(Form("Required value #%d, while maximum is %d", ival, nval));
633 continue;
634 }
61f275d1 635 AliRsnMiniValue *val = (AliRsnMiniValue *)valueList->At(ival);
03d23846 636 if (!val) {
637 AliError(Form("Value in position #%d is NULL", ival));
638 continue;
61f275d1 639 }
03d23846 640 // if none of the above exit points is taken, compute value
641 fComputed[i] = val->Eval(&fPair, event);
642 }
03d23846 643}
644
645//________________________________________________________________________________________
45aa62b9 646void AliRsnMiniOutput::FillHistogram()
03d23846 647{
648//
649// Fills the internal histogram using the current values stored in the
650// 'fComputed' array, in the order as they are stored, up to the max
651// dimension of the initialized histogram itself.
652//
653
654 // retrieve object from list
655 if (!fList) {
656 AliError("List pointer is NULL");
45aa62b9 657 return;
03d23846 658 }
659 TObject *obj = fList->At(fOutputID);
660
661 if (obj->InheritsFrom(TH1F::Class())) {
61f275d1 662 ((TH1F *)obj)->Fill(fComputed[0]);
03d23846 663 } else if (obj->InheritsFrom(TH2F::Class())) {
61f275d1 664 ((TH2F *)obj)->Fill(fComputed[0], fComputed[1]);
03d23846 665 } else if (obj->InheritsFrom(TH3F::Class())) {
61f275d1 666 ((TH3F *)obj)->Fill(fComputed[0], fComputed[1], fComputed[2]);
03d23846 667 } else if (obj->InheritsFrom(THnSparseF::Class())) {
74d60285 668 THnSparseF *h = (THnSparseF *)obj;
669 if (fCheckHistRange) {
670 for (Int_t iAxis = 0; iAxis<h->GetNdimensions(); iAxis++) {
671 if (fComputed.At(iAxis)>h->GetAxis(iAxis)->GetXmax() || fComputed.At(iAxis)<h->GetAxis(iAxis)->GetXmin()) return;
672 }
673 }
674 h->Fill(fComputed.GetArray());
03d23846 675 } else {
676 AliError("No output initialized");
03d23846 677 }
03d23846 678}