]> git.uio.no Git - u/mrichter/AliRoot.git/blame - PWGLF/SPECTRA/Nuclei/B2/AliLnSecondaries.cxx
git://git.uio.no / u / mrichter / AliRoot.git / blame
? search:
summary | shortlog | log | commit | commitdiff | tree
blob | blame (incremental) | history | HEAD
Merge branch 'feature-movesplit'
[u/mrichter/AliRoot.git] / PWGLF / SPECTRA / Nuclei / B2 / AliLnSecondaries.cxx
CommitLineData
2403d402 1/**************************************************************************
2 * Copyright(c) 1998-1999, ALICE Experiment at CERN, All rights reserved. *
3 * *
4 * Author: The ALICE Off-line Project. *
5 * Contributors are mentioned in the code where appropriate. *
6 * *
7 * Permission to use, copy, modify and distribute this software and its *
8 * documentation strictly for non-commercial purposes is hereby granted *
9 * without fee, provided that the above copyright notice appears in all *
10 * copies and that both the copyright notice and this permission notice *
11 * appear in the supporting documentation. The authors make no claims *
12 * about the suitability of this software for any purpose. It is *
13 * provided "as is" without express or implied warranty. *
14 **************************************************************************/
15
16// removal of secondaries using DCA distributions
17// author: Eulogio Serradilla <eulogio.serradilla@cern.ch>
18
19#include <Riostream.h>
20#include <TFile.h>
21#include <TString.h>
22#include <TF1.h>
23#include <TMath.h>
24#include <TFractionFitter.h>
25#include <TObjArray.h>
26#include <TH1D.h>
27#include <TH2D.h>
1cb68411 28#include <TBackCompFitter.h>
2403d402 29
30#include "AliLnSecondaries.h"
31#include "B2.h"
32
33ClassImp(AliLnSecondaries)
34
35AliLnSecondaries::AliLnSecondaries(const TString& particle, const TString& dataFilename, const TString& simuFilename, const TString& outputFilename, const TString& otag)
36: TObject()
37, fParticle(particle)
38, fDataFilename(dataFilename)
39, fSimuFilename(simuFilename)
40, fOutputFilename(outputFilename)
41, fOutputTag(otag)
aa54def0 42, fPtMin(0.5)
43, fPtMax(3.0)
2403d402 44, fNbin(1)
45, fMinDCAxy(-1.5)
46, fMaxDCAxy(1.5)
47, fFracProc(0)
48, fMatDCAxyMod(AliLnSecondaries::kFlatDCAxy)
fe25d981 49, fANucTemplate(0)
2403d402 50, fScMat(1)
51, fScFd(1)
0f539a2b 52, fAddFakeTracks(1)
2403d402 53{
54//
55// constructor
56//
57 TH1::SetDefaultSumw2(); // switch on histogram errors
2403d402 58}
59
60AliLnSecondaries::~AliLnSecondaries()
61{
62//
63// destructor
64//
65}
66
67Int_t AliLnSecondaries::Exec()
68{
69//
70// Fractions of primaries and secondaries using templates from simulations
71// Secondaries are from feed-down and/or from materials
72//
73 using namespace std;
74
75 TFile* fdata = new TFile(fDataFilename.Data(),"read");
76 if(fdata->IsZombie()) exit(1);
77
78 TFile* fsimu = new TFile(fSimuFilename.Data(),"read");
79 if(fsimu->IsZombie()) exit(1);
80
81 TFile* fdebug = new TFile(Form("debug-%s",fOutputFilename.Data()),"recreate");
82
83 if(fOutputTag != "") fdebug->mkdir(fOutputTag.Data());
84 if(fOutputTag != "") fdebug->cd(fOutputTag.Data());
85
86 // --------- ideal values: primaries + no secondaries --------
87
aa54def0 88 TH1D* hPidPt = FindObj<TH1D>(fdata, fParticle + "_PID_Pt");
2403d402 89 hPidPt->Write();
90
91 const char* contrib[] = { "prim", "mat", "fdwn" };
92
93 TH1D* hFracPt[3];
94 for(Int_t i=0; i<3; ++i)
95 {
96 hFracPt[i] = this->ZeroClone(hPidPt, Form("%s_P%s_Pt",fParticle.Data(),contrib[i]));
97 hFracPt[i]->SetYTitle(Form("P_{%s}",contrib[i]));
98 }
99
100 // fill
101
102 if(fFracProc == kMonteCarlo)
103 {
104 // compute the fractions from the montecarlo simulation
aa54def0 105 TH1D* hAll = FindObj<TH1D>(fsimu, fParticle + "_Sim_Pt");
106 TH1D* hPrim = FindObj<TH1D>(fsimu, fParticle + "_Sim_Prim_Pt");
107 TH1D* hMat = FindObj<TH1D>(fsimu, fParticle + "_Sim_Mat_Pt");
108 TH1D* hFdwn = FindObj<TH1D>(fsimu, fParticle + "_Sim_Fdwn_Pt");
2403d402 109
110 delete hFracPt[0];
111 delete hFracPt[1];
112 delete hFracPt[2];
113
114 hFracPt[0] = Divide(hPrim, hAll, fParticle + "_Pprim_Pt");
115 hFracPt[1] = Divide(hMat, hAll, fParticle + "_Pmat_Pt");
116 hFracPt[2] = Divide(hFdwn, hAll, fParticle + "_Pfdwn_Pt");
117
118 hFracPt[1]->Scale(fScMat);
119 hFracPt[2]->Scale(fScFd);
120 }
121 else if(fParticle == "AntiProton")
122 {
123 // assume only secondaries from feed-down
aa54def0 124 TH2D* hDataDCAxyPt = FindObj<TH2D>(fdata, fParticle + "_PID_DCAxy_Pt");
125 TH2D* hDataMCDCAxyPt = FindObj<TH2D>(fsimu, fParticle + "_PID_DCAxy_Pt");
126 TH2D* hPrimDCAxyPt = FindObj<TH2D>(fsimu, fParticle + "_Sim_PID_Prim_DCAxy_Pt");
127 TH2D* hFdwnDCAxyPt = FindObj<TH2D>(fsimu, fParticle + "_Sim_PID_Fdwn_DCAxy_Pt");
128 TH2D* hFakePrimDCAxyPt = FindObj<TH2D>(fsimu, fParticle + "_Sim_PID_Fake_Prim_DCAxy_Pt");
129 TH2D* hFakeFdwnDCAxyPt = FindObj<TH2D>(fsimu, fParticle + "_Sim_PID_Fake_Fdwn_DCAxy_Pt");
0f539a2b 130
131 if(fAddFakeTracks)
132 {
133 hPrimDCAxyPt->Add(hFakePrimDCAxyPt);
134 hFdwnDCAxyPt->Add(hFakeFdwnDCAxyPt);
135 }
2403d402 136
137 hDataDCAxyPt->Rebin2D(1,fNbin);
138 hDataMCDCAxyPt->Rebin2D(1,fNbin);
139 hPrimDCAxyPt->Rebin2D(1,fNbin);
140 hFdwnDCAxyPt->Rebin2D(1,fNbin);
141
142 this->GetFraction(hFracPt[0], hFracPt[2], hDataDCAxyPt, hDataMCDCAxyPt, hPrimDCAxyPt, hFdwnDCAxyPt, "Fdwn");
143 }
144 else if(fParticle.Contains("Anti"))
145 {
146 // assume there are no secondaries for antinuclei
147 this->GetFraction(hFracPt[0]);
148 }
149 else if(fParticle == "Proton")
150 {
151 // secondaries from materials and feed-down
aa54def0 152 TH2D* hDataDCAxyPt = FindObj<TH2D>(fdata, fParticle + "_PID_DCAxy_Pt");
153 TH2D* hDataMCDCAxyPt = FindObj<TH2D>(fsimu, fParticle + "_PID_DCAxy_Pt");
154 TH2D* hPrimDCAxyPt = FindObj<TH2D>(fsimu, fParticle + "_Sim_PID_Prim_DCAxy_Pt");
155 TH2D* hMatDCAxyPt = FindObj<TH2D>(fsimu, fParticle + "_Sim_PID_Mat_DCAxy_Pt");
156 TH2D* hFdwnDCAxyPt = FindObj<TH2D>(fsimu, fParticle + "_Sim_PID_Fdwn_DCAxy_Pt");
157 TH2D* hFakePrimDCAxyPt = FindObj<TH2D>(fsimu, fParticle + "_Sim_PID_Fake_Prim_DCAxy_Pt");
158 TH2D* hFakeMatDCAxyPt = FindObj<TH2D>(fsimu, fParticle + "_Sim_PID_Fake_Mat_DCAxy_Pt");
159 TH2D* hFakeFdwnDCAxyPt = FindObj<TH2D>(fsimu, fParticle + "_Sim_PID_Fake_Fdwn_DCAxy_Pt");
0f539a2b 160
161 if(fAddFakeTracks)
162 {
163 hPrimDCAxyPt->Add(hFakePrimDCAxyPt);
164 hMatDCAxyPt->Add(hFakeMatDCAxyPt);
165 hFdwnDCAxyPt->Add(hFakeFdwnDCAxyPt);
166 }
2403d402 167
168 if(fMatDCAxyMod == kFlatDCAxy)
169 {
170 hMatDCAxyPt = this->GetFlatDCAxyPt(10, hPrimDCAxyPt, fParticle + "_Sim_PID_Mat_DCAxy_Pt");
171 }
172
173 hDataDCAxyPt->Rebin2D(1,fNbin);
174 hDataMCDCAxyPt->Rebin2D(1,fNbin);
175 hPrimDCAxyPt->Rebin2D(1,fNbin);
176 hMatDCAxyPt->Rebin2D(1,fNbin);
177 hFdwnDCAxyPt->Rebin2D(1,fNbin);
178
179 this->GetFraction(hFracPt, hDataDCAxyPt, hDataMCDCAxyPt, hPrimDCAxyPt, hMatDCAxyPt, hFdwnDCAxyPt);
180 }
181 else
182 {
183 // only secondaries from materials for nuclei
aa54def0 184 TH2D* hDataDCAxyPt = FindObj<TH2D>(fdata, fParticle + "_PID_DCAxy_Pt");
185 TH2D* hDataMCDCAxyPt = FindObj<TH2D>(fsimu, fParticle + "_PID_DCAxy_Pt");
186 TH2D* hPrimDCAxyPt = FindObj<TH2D>(fsimu, fParticle + "_Sim_PID_Prim_DCAxy_Pt");
187 TH2D* hFakePrimDCAxyPt = FindObj<TH2D>(fsimu, fParticle + "_Sim_PID_Fake_Prim_DCAxy_Pt");
0f539a2b 188
189 if(fAddFakeTracks) hPrimDCAxyPt->Add(hFakePrimDCAxyPt);
2403d402 190
fe25d981 191 if(fANucTemplate)
192 {
aa54def0 193 hPrimDCAxyPt = FindObj<TH2D>(fdata, Form("Anti%s_PID_DCAxy_Pt",fParticle.Data()));
fe25d981 194 }
195
aa54def0 196 TH2D* hMatDCAxyPt = FindObj<TH2D>(fsimu, fParticle + "_Sim_PID_Mat_DCAxy_Pt");
2403d402 197 if(fMatDCAxyMod == kFlatDCAxy)
198 {
199 hMatDCAxyPt = this->GetFlatDCAxyPt(10, hPrimDCAxyPt, fParticle + "_Sim_PID_Mat_DCAxy_Pt");
200 }
201
202 hDataDCAxyPt->Rebin2D(1,fNbin);
203 hDataMCDCAxyPt->Rebin2D(1,fNbin);
204 hPrimDCAxyPt->Rebin2D(1,fNbin);
205 hMatDCAxyPt->Rebin2D(1,fNbin);
206
207 this->GetFraction(hFracPt[0], hFracPt[1], hDataDCAxyPt, hDataMCDCAxyPt, hPrimDCAxyPt, hMatDCAxyPt, "Mat");
208 }
209
210 // --------- save --------------
211
212 TFile* foutput = new TFile(fOutputFilename.Data(),"recreate");
213 if(fOutputTag != "")
214 {
215 foutput->mkdir(fOutputTag.Data());
216 foutput->cd(fOutputTag.Data());
217 }
218
219 // ----- fractions w.r.t. primaries ------------
220
221 TH1D* hMatFracPt = (TH1D*)hFracPt[1]->Clone(Form("%s_Frac_Mat_Pt", fParticle.Data()));
222 hMatFracPt->SetYTitle("P_{1} / P_{0}");
223 hMatFracPt->Divide(hFracPt[0]);
224 hMatFracPt->Write();
225
226 TH1D* hFdwnFracPt = (TH1D*)hFracPt[2]->Clone( Form("%s_Frac_Fdwn_Pt",fParticle.Data()));
227 hFdwnFracPt->SetYTitle("P_{2} / P_{0}");
228 hFdwnFracPt->Divide(hFracPt[0]);
229 hFdwnFracPt->Write();
230
231 // --------- optionally fit fractions to extrapolate at high pt -----------
232
233 TF1* fncMat = this->GetMatFraction(fParticle + "_Frac_Mat_Fit_Pt");
234 if(fParticle.Contains("Anti"))
235 {
236 fncMat->SetParameters(0,0,0);
237 }
238 else if(fParticle=="Proton")
239 {
240 hMatFracPt->Fit(fncMat, "NQ", "", 0.,2);
241 }
242 else
243 {
aa54def0 244 hMatFracPt->Fit(fncMat, "NQ", "", 0.,1.6);
2403d402 245 }
246
247 TF1* fncFdwn = this->GetFdwnFraction(fParticle + "_Frac_Fdwn_Fit_Pt");
248 if(!fParticle.Contains("Proton"))
249 {
250 fncFdwn->SetParameters(0,0,0);
251 }
252 else
253 {
254 hFdwnFracPt->Fit(fncFdwn, "NQ", "", 0.5,3.);
255 }
256
257 fncFdwn->Write();
258 fncMat->Write();
259
260 // ---------- clean ---------
261
262 delete fncFdwn;
263 delete fncMat;
264
265 delete hMatFracPt;
266 delete hFdwnFracPt;
267
268 for(Int_t i=0; i<3; ++i)
269 {
270 hFracPt[i]->Write();
271 delete hFracPt[i];
272 }
273
274 delete foutput;
275 delete fdebug;
276 delete fsimu;
277 delete fdata;
278
1cb68411 279 delete dynamic_cast<TBackCompFitter*>(TVirtualFitter::GetFitter());
280
2403d402 281 return 0;
282}
283
284void AliLnSecondaries::GetFraction(TH1D* hPrimPt) const
285{
286//
287// No secondaries only primaries
288//
289 TF1* one = new TF1("one", "1", hPrimPt->GetXaxis()->GetXmin(), hPrimPt->GetXaxis()->GetXmax());
290 hPrimPt->Reset();
fe25d981 291 hPrimPt->Sumw2();
2403d402 292 hPrimPt->Add(one);
293 delete one;
294}
295
296void AliLnSecondaries::GetFraction(TH1D* hPrimPt, TH1D* hSecPt, const TH2D* hDCAxyPt, const TH2D* hMCDCAxyPt, const TH2D* hPrimDCAxyPt, const TH2D* hSecDCAxyPt, const TString& sec) const
297{
298//
299// slice the DCA distributions and compute the fractions
300// 2 contributions (primaries and secondaries)
301//
302 TString nosec = (sec == "Fdwn") ? "Mat" : "Fdwn";
303
aa54def0 304 Int_t lowbin = hDCAxyPt->GetXaxis()->FindFixBin(fPtMin);
305 Int_t hibin = hDCAxyPt->GetXaxis()->FindFixBin(fPtMax);
306
307 for(Int_t i=lowbin; i<hibin; ++i)
2403d402 308 {
309 TH1D* hDCAxy = hDCAxyPt->ProjectionY(Form("%s_Data_DCAxy_%02d",fParticle.Data(),i),i,i);
310 TH1D* hMCDCAxy = hMCDCAxyPt->ProjectionY(Form("%s_SimData_DCAxy_%02d",fParticle.Data(),i),i,i);
311 TH1D* hPrimDCAxy = hPrimDCAxyPt->ProjectionY(Form("%s_Prim_DCAxy_%02d",fParticle.Data(),i),i,i);
312 TH1D* hSecDCAxy = hSecDCAxyPt->ProjectionY(Form("%s_%s_DCAxy_%02d",fParticle.Data(),sec.Data(),i),i,i);
313
314 // fractions
aa54def0 315 Double_t frac[2] = {0};
316 Double_t err[2] = {0};
2403d402 317
aa54def0 318 this->GetTFFfractions(frac, err, hDCAxy, hPrimDCAxy, hSecDCAxy, i, sec);
2403d402 319
320 // fill
321 TH1D* hFracPt[] = { hPrimPt, hSecPt};
322 for(Int_t j=0; j<2; ++j)
323 {
324 hFracPt[j]->SetBinContent(i, frac[j]);
325 hFracPt[j]->SetBinError(i, err[j]);
326 }
327
328 // ------------- debug histograms ------------
329
330 hDCAxy->Write();
331 hMCDCAxy->Write();
332 hPrimDCAxy->Write();
333 hSecDCAxy->Write();
334
335 TH1D* hNoSecDCAxy = this->ZeroClone(hMCDCAxy, Form("%s_%s_DCAxy_%02d",fParticle.Data(),nosec.Data(),i));
336 TH1D* hNoSecFit = this->ZeroClone(hMCDCAxy, Form("%s_Fit_%s_DCAxy_%02d",fParticle.Data(),nosec.Data(),i));
337
338 hNoSecDCAxy->Write();
339 hNoSecFit->Write();
340
341 delete hNoSecDCAxy;
342 delete hNoSecFit;
343
344 // ----------------------- end debug ----------------------
345
346 delete hDCAxy;
347 delete hMCDCAxy;
348 delete hPrimDCAxy;
349 delete hSecDCAxy;
350 }
351}
352
353void AliLnSecondaries::GetFraction(TH1D* hFracPt[3], const TH2D* hDCAxyPt, const TH2D* hMCDCAxyPt, const TH2D* hPrimDCAxyPt, const TH2D* hMatDCAxyPt, const TH2D* hFdwnDCAxyPt) const
354{
355//
356// slice the DCA distribution and get the fractions for each pt bin
357// (3 contributions)
358//
aa54def0 359 Int_t lowbin = hDCAxyPt->GetXaxis()->FindFixBin(fPtMin);
360 Int_t hibin = hDCAxyPt->GetXaxis()->FindFixBin(fPtMax);
361
362 for(Int_t i=lowbin; i<hibin; ++i)
2403d402 363 {
364 // slices
365 TH1D* hDCAxy = hDCAxyPt->ProjectionY(Form("%s_Data_DCAxy_%02d",fParticle.Data(),i),i,i);
366 TH1D* hMCDCAxy = hMCDCAxyPt->ProjectionY(Form("%s_SimData_DCAxy_%02d",fParticle.Data(),i),i,i);
367 TH1D* hPrimDCAxy = hPrimDCAxyPt->ProjectionY(Form("%s_Prim_DCAxy_%02d",fParticle.Data(),i),i,i);
368 TH1D* hMatDCAxy = hMatDCAxyPt->ProjectionY(Form("%s_Mat_DCAxy_%02d",fParticle.Data(),i),i,i);
369 TH1D* hFdwnDCAxy = hFdwnDCAxyPt->ProjectionY(Form("%s_Fdwn_DCAxy_%02d",fParticle.Data(),i),i,i);
370
371 // fractions
372 Double_t frac[3] = {0};
373 Double_t err[3] = {0};
374
aa54def0 375 this->GetTFFfractions(frac, err, hDCAxy, hPrimDCAxy, hMatDCAxy, hFdwnDCAxy, i);
2403d402 376
377 // fill
378 for(Int_t k=0; k<3; ++k)
379 {
380 hFracPt[k]->SetBinContent(i, frac[k]);
381 hFracPt[k]->SetBinError(i, err[k]);
382 }
383
384 // -------------- debug ---------------------------
385
386 hDCAxy->Write();
387 hMCDCAxy->Write();
388 hPrimDCAxy->Write();
389 hMatDCAxy->Write();
390 hFdwnDCAxy->Write();
391
392 // --------------------- end debug ----------------
393
394 delete hDCAxy;
395 delete hMCDCAxy;
396 delete hPrimDCAxy;
397 delete hMatDCAxy;
398 delete hFdwnDCAxy;
399 }
400}
401
fe25d981 402Int_t AliLnSecondaries::GetTFFfractions(Double_t* frac, Double_t* err, TH1D* hData, TH1D* hPrim, TH1D* hSec, Int_t ibin, const TString& sec) const
2403d402 403{
404//
405// find the fractions of 2 contributions
406//
407 TObjArray* mc = new TObjArray(2);
408
409 mc->Add(hPrim);
410 mc->Add(hSec);
411
412 TFractionFitter* fit = new TFractionFitter(hData, mc, "Q"); // initialise
413 fit->Constrain(1,0.0001,1.);
414 fit->Constrain(2,0.0001,1.);
415
416 Int_t status = fit->Fit();
aa54def0 417 status = fit->Fit();
418 status = fit->Fit();
2403d402 419
420 if (status == 0) // get fractions
421 {
422 fit->GetResult(0, frac[0], err[0]);
423 fit->GetResult(1, frac[1], err[1]);
424 }
425
fe25d981 426 const char* contrib[] = {"Prim", sec.Data() };
2403d402 427 this->WriteTFFdebug(hData, fit, status, ibin, contrib, frac, 2);
428
429 delete mc;
430 delete fit;
431
432 return status;
433}
434
435Int_t AliLnSecondaries::GetTFFfractions(Double_t* frac, Double_t* err, TH1D* hData, TH1D* hPrim, TH1D* hMat, TH1D* hFdwn, Int_t ibin) const
436{
437//
438// find the fractions of 3 contributions
439//
440 TObjArray* mc = new TObjArray(3);
441
442 mc->Add(hPrim);
443 mc->Add(hMat);
444 mc->Add(hFdwn);
445
446 TFractionFitter* fit = new TFractionFitter(hData, mc, "Q");
447 fit->Constrain(1,0.0001,1.);
448 fit->Constrain(2,0.0001,1.);
449 fit->Constrain(3,0.0001,1.);
450 //fit->SetRangeX(1,50);
451
452 Int_t status = fit->Fit();
aa54def0 453 status = fit->Fit();
454 status = fit->Fit();
2403d402 455
456 if (status == 0) // get fractions
457 {
458 fit->GetResult(0, frac[0], err[0]);
459 fit->GetResult(1, frac[1], err[1]);
460 fit->GetResult(2, frac[2], err[2]);
461 }
462
463 const char* contrib[] = {"Prim", "Mat", "Fdwn"};
464 this->WriteTFFdebug(hData, fit, status, ibin, contrib, frac, 3);
465
466 delete mc;
467 delete fit;
468
469 return status;
470}
471
1cb68411 472void AliLnSecondaries::WriteTFFdebug(const TH1D* hData, TFractionFitter* fit, Int_t status, Int_t ibin, const char* contrib[], Double_t* frac, Int_t kmax) const
2403d402 473{
474//
475// Write TFractionFitter debug histograms
476//
477 if (status == 0)
478 {
479 TH1D* hResult = (TH1D*) fit->GetPlot();
480
481 hResult->SetName(Form("%s_Fit_Data_DCAxy_%02d",fParticle.Data(),ibin));
482 hResult->Write();
483
484 for(Int_t k=0; k<kmax; ++k)
485 {
486 TH1D* hMCpred = (TH1D*) fit->GetMCPrediction(k);
487 hMCpred->SetName(Form("%s_Fit_%s_DCAxy_%02d",fParticle.Data(),contrib[k],ibin));
488
489 hMCpred->Sumw2();
490 hMCpred->Scale(frac[k]*hResult->Integral()/hMCpred->Integral());
491
492 hMCpred->Write();
493 }
494 }
495 else // write void histograms
496 {
497 TH1D* hZero = this->ZeroClone(hData,Form("%s_Fit_Data_DCAxy_%02d",fParticle.Data(),ibin));
498 hZero->Write();
499
500 for(Int_t k=0; k<kmax; ++k)
501 {
502 hZero->SetName(Form("%s_Fit_%s_DCAxy_%02d",fParticle.Data(),contrib[k],ibin));
503 hZero->Write();
504 }
505 delete hZero;
506 }
507}
508
fe25d981 509TH1D* AliLnSecondaries::ZeroClone(const TH1D* h, const TString& name) const
2403d402 510{
511//
512// clone histogram and reset to zero
513//
aa54def0 514 TH1D* clone = dynamic_cast<TH1D*>(h->Clone(name.Data()));
515 if(clone == 0)
516 {
517 this->Error("ZeroClone", "could not clone %s", h->GetName());
518 exit(1);
519 }
2403d402 520 clone->Reset();
fe25d981 521 clone->Sumw2();
2403d402 522
523 return clone;
524}
525
2403d402 526TH2D* AliLnSecondaries::GetFlatDCAxyPt(Double_t max, const TH2D* hDCAxyPt, const TString& name) const
527{
528//
529// generate a flat dcaxy-pt distribution
530//
531 TF1* fnc = new TF1("gen","1", fMinDCAxy, fMaxDCAxy);
532 TH2D* h = (TH2D*)hDCAxyPt->Clone(name.Data());
533 h->Reset();
fe25d981 534 h->Sumw2();
2403d402 535
536 TAxis* xAxis = h->GetXaxis();
537 TAxis* yAxis = h->GetYaxis();
538
539 for(Int_t i=1; i<xAxis->GetNbins()+1; ++i)
540 {
541 Double_t pt = xAxis->GetBinCenter(i);
542
543 for(Int_t j=1; j<yAxis->GetNbins()+1; ++j)
544 {
545 for(Int_t k=0; k<max; ++k) h->Fill(pt,fnc->GetRandom());
546 }
547 }
548
549 delete fnc;
550 return h;
551}
552
553TF1* AliLnSecondaries::GetMatFraction(const TString& name) const
554{
555//
556// model the material fraction as an exponential + a constant
557//
558 TF1* fnc = new TF1(name.Data(), "[0]*exp(-[1]*x)+[2]",0,10);
559 fnc->SetParameters(1.7,5.5,0.01);
560
561 fnc->SetParLimits(2,0,0.02); // asymptotic value
562
563 return fnc;
564}
565
566TF1* AliLnSecondaries::GetFdwnFraction(const TString& name) const
567{
568//
569// model the feed-down fraction as an exponential + constant
570//
571 TF1* fnc = new TF1(name, "[0]*exp(-[1]*x)+[2]",0,10);
572 fnc->SetParameters(0.5,1,0.02);
573
574 fnc->SetParLimits(2,0,0.3); // asymptotic value
575
576 return fnc;
577}