]>
Commit | Line | Data |
---|---|---|
6951a056 | 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 | // | |
17 | // TRD dEdx base utils | |
18 | // xx | |
19 | // xx | |
20 | // xx | |
21 | // xx | |
22 | // | |
23 | // Xianguo Lu | |
24 | // lu@physi.uni-heidelberg.de | |
25 | // Xianguo.Lu@cern.ch | |
26 | // | |
27 | // | |
28 | ||
29 | #include "TF1.h" | |
30 | #include "TFile.h" | |
31 | #include "TH1D.h" | |
32 | #include "TH2D.h" | |
33 | #include "THnSparse.h" | |
34 | #include "TMath.h" | |
35 | #include "TMatrixD.h" | |
36 | #include "TMinuit.h" | |
37 | #include "TObjArray.h" | |
38 | #include "TRandom3.h" | |
39 | #include "TStopwatch.h" | |
40 | #include "TVectorD.h" | |
41 | ||
42 | #include "TTreeStream.h" | |
43 | ||
6951a056 | 44 | #include "AliESDEvent.h" |
45 | #include "AliESDfriendTrack.h" | |
46 | #include "AliESDtrack.h" | |
6951a056 | 47 | #include "AliTRDtrackV1.h" |
48 | ||
49 | #include "AliTRDdEdxBaseUtils.h" | |
50 | ||
6fa94550 | 51 | #define EPSILON 1e-8 |
6951a056 | 52 | |
687aa844 | 53 | Double_t AliTRDdEdxBaseUtils::fgQ0Frac = 0.55; |
6951a056 | 54 | Double_t AliTRDdEdxBaseUtils::fgQ1Frac = 0.5; |
55 | Double_t AliTRDdEdxBaseUtils::fgTimeBinCountCut = 0.0; | |
56 | Int_t AliTRDdEdxBaseUtils::fgCalibTPCnclsCut = 70; | |
57 | Bool_t AliTRDdEdxBaseUtils::fgExBOn = kTRUE; | |
58 | Bool_t AliTRDdEdxBaseUtils::fgPadGainOn = kTRUE; | |
687aa844 | 59 | Double_t AliTRDdEdxBaseUtils::fgQScale = 50; |
6951a056 | 60 | |
61 | //=================================================================================== | |
62 | // Math and Histogram | |
63 | //=================================================================================== | |
6856b6a8 | 64 | void AliTRDdEdxBaseUtils::BinLogX(TAxis *axis) |
65 | { | |
66 | // | |
67 | // Method for the correct logarithmic binning of histograms | |
68 | // copied and modified from AliTPCcalibBase | |
69 | ||
70 | const Int_t bins = axis->GetNbins(); | |
71 | ||
72 | const Double_t from = axis->GetXmin(); | |
73 | const Double_t to = axis->GetXmax(); | |
74 | if (from<EPSILON) return; | |
75 | Double_t *new_bins = new Double_t[bins + 1]; | |
76 | ||
77 | new_bins[0] = from; | |
78 | const Double_t factor = pow(to/from, 1./bins); | |
79 | ||
80 | for (int i = 1; i <= bins; i++) { | |
81 | new_bins[i] = factor * new_bins[i-1]; | |
82 | } | |
83 | axis->Set(bins, new_bins); | |
84 | delete [] new_bins; | |
85 | } | |
86 | ||
6951a056 | 87 | void AliTRDdEdxBaseUtils::GetCDFCuts(const TH1D *hh, const Int_t ncut, Double_t cuts[], const Double_t cdfs[], const Double_t thres) |
88 | { | |
89 | // | |
90 | //counts of hh is sorted | |
91 | // | |
92 | ||
93 | for(Int_t ii=0; ii<ncut; ii++){ | |
94 | cuts[ii] = -999; | |
95 | } | |
96 | ||
97 | Int_t nsel = 0; | |
98 | const Int_t nbin = hh->GetNbinsX(); | |
99 | Double_t datas[nbin]; | |
100 | for(Int_t ii=1; ii<=nbin; ii++){ | |
101 | const Double_t res = hh->GetBinContent(ii); | |
102 | if(res<thres){ | |
103 | continue; | |
104 | } | |
105 | ||
106 | datas[nsel] = res; | |
107 | nsel++; | |
108 | } | |
109 | if(!nsel) | |
110 | return; | |
111 | ||
112 | Int_t id[nsel]; | |
113 | TMath::Sort(nsel, datas, id, kFALSE); | |
114 | ||
115 | for(Int_t ii=0; ii<ncut; ii++){ | |
116 | const Double_t icdf = cdfs[ii]; | |
117 | if(icdf<0 || icdf>1){ | |
118 | printf("AliTRDdEdxBaseUtils::GetCDFCuts error cdfs[%d] %15f out of range!\n", ii, icdf); exit(1); | |
119 | } | |
120 | cuts[ii] = datas[id[Int_t(icdf*nsel)]]; | |
121 | } | |
122 | } | |
123 | ||
124 | Double_t AliTRDdEdxBaseUtils::GetMeanRMS(const Double_t nn, const Double_t sum, const Double_t w2s, Double_t * grms, Double_t * gerr) | |
125 | { | |
126 | // | |
127 | //calculate mean (with error) and rms from sum, w2s, nn | |
128 | //if nn=0, mean, error, and rms all = 0 | |
129 | // | |
130 | ||
131 | Double_t tmean = 0, trms = 0, terr = 0; | |
132 | ||
133 | if(nn>EPSILON){ | |
134 | tmean = sum/nn; | |
135 | ||
136 | const Double_t arg = w2s/nn-tmean*tmean; | |
137 | if(TMath::Abs(arg)<EPSILON){ | |
138 | trms = 0; | |
139 | } | |
140 | else{ | |
141 | if( arg <0 ){ | |
142 | printf("AliTRDdEdxBaseUtils::GetMeanRMS error negative sqrt argument!! %e -- %e %e %f\n", arg, w2s, sum, nn); exit(1); | |
143 | } | |
144 | ||
145 | trms = TMath::Sqrt(arg); | |
146 | } | |
147 | ||
148 | terr = trms/TMath::Sqrt(nn); | |
149 | } | |
150 | ||
151 | if(grms){ | |
152 | (*grms) = trms; | |
153 | } | |
154 | ||
155 | if(gerr){ | |
156 | (*gerr) = terr; | |
157 | } | |
158 | ||
159 | return tmean; | |
160 | } | |
161 | ||
162 | Double_t AliTRDdEdxBaseUtils::TruncatedMean(const Int_t nx, const Double_t xdata[], const Double_t lowfrac, const Double_t highfrac, Double_t * grms, Double_t * gerr, Double_t *wws) | |
163 | { | |
164 | // | |
165 | //calculate truncated mean | |
166 | //return <x*w>_{low-high according to x} | |
167 | // | |
168 | ||
169 | /* | |
170 | //test-> | |
171 | for(Int_t ii=0; ii<nx; ii++){ | |
172 | printf("test %d/%d %f\n", ii, nx, xdata[ii]); | |
173 | } | |
174 | //<--test | |
175 | */ | |
176 | ||
177 | Int_t index[nx]; | |
178 | TMath::Sort(nx, xdata, index, kFALSE); | |
179 | ||
180 | Int_t nused = 0; | |
181 | Double_t sum = 0; | |
182 | Double_t w2s = 0; | |
183 | const Int_t istart = Int_t (nx*lowfrac); | |
184 | const Int_t istop = Int_t (nx*highfrac); | |
185 | ||
186 | //=,< correct, because when low=0, high=1 it is correct | |
187 | for(Int_t ii=istart; ii<istop; ii++){ | |
188 | Double_t weight = 1; | |
189 | if(wws){ | |
190 | weight = wws[index[ii]]; | |
191 | } | |
192 | const Double_t sx = xdata[index[ii]]*weight; | |
193 | ||
194 | sum += sx; | |
195 | w2s += sx*sx; | |
196 | ||
197 | nused++; | |
198 | //printf("test in loop %d/%d %f %f %f\n", ii, nused, sx, sum, w2s); | |
199 | ||
200 | } | |
201 | ||
202 | return GetMeanRMS(nused, sum, w2s, grms, gerr); | |
203 | } | |
204 | ||
205 | Double_t AliTRDdEdxBaseUtils::TruncatedMean(const TH1 *hh, const Double_t lowfrac, const Double_t highfrac, Double_t * grms, Double_t * gerr) | |
206 | { | |
207 | // | |
208 | //do truncation on histogram | |
209 | // | |
210 | //if hh is scaled, be sure Sumw2 is called before scaling!! then mean, rms and err will all be correct | |
211 | ||
212 | //with under-/over-flow | |
213 | Double_t npreTrunc = 0; | |
214 | for(Int_t itmp=0; itmp<=hh->GetNbinsX()+1; itmp++){ | |
215 | const Double_t be = hh->GetBinError(itmp); | |
216 | const Double_t bc = hh->GetBinContent(itmp); | |
217 | if(be<EPSILON){ | |
218 | if(bc>EPSILON){ | |
219 | printf("AliTRDdEdxBaseUtils::TruncatedMean (hist) error %e %e %d\n", bc, be, itmp); exit(1); | |
220 | } | |
221 | continue; | |
222 | } | |
223 | npreTrunc += bc*bc/be/be; | |
224 | } | |
225 | ||
226 | const Double_t nstart = npreTrunc*lowfrac; | |
227 | const Double_t nstop = npreTrunc*highfrac; | |
228 | ||
229 | //with Double_t this should also handle normalized hist | |
230 | Double_t ntot = 0; | |
231 | Double_t nused = 0; | |
232 | Double_t sum = 0; | |
233 | Double_t w2s = 0; | |
234 | for(Int_t itmp=0; itmp<=hh->GetNbinsX()+1; itmp++){ | |
235 | const Double_t be = hh->GetBinError(itmp); | |
236 | const Double_t bc = hh->GetBinContent(itmp); | |
237 | if(be<EPSILON){ | |
238 | if(bc>EPSILON){ | |
239 | printf("AliTRDdEdxBaseUtils::TruncatedMean (hist) error %e %e %d\n", bc, be, itmp); exit(1); | |
240 | } | |
241 | continue; | |
242 | } | |
243 | const Double_t weight = bc*bc/be/be; | |
244 | ntot+=weight; | |
245 | //<= correct, because when high=1, nstop has to be included | |
246 | if(ntot>nstart && ntot<=nstop){ | |
247 | ||
248 | const Double_t val = hh->GetBinCenter(itmp); | |
249 | sum += weight*val; | |
250 | w2s += weight*val*val; | |
251 | ||
252 | nused += weight; | |
253 | ||
254 | //printf("test %d %f %f --- %f %f -- %f %f\n", itmp, weight, val, sum, w2s, nused, nsample); | |
255 | } | |
256 | else if(ntot>nstop){ | |
257 | if(itmp>=hh->GetNbinsX()){ | |
258 | printf("AliTRDdEdxBaseUtils::TruncatedMean warning hist range too small %s %f %f %d %d, %15f %15f %15f; nused w2s sum set to 0\n", hh->GetName(), hh->GetBinLowEdge(1), hh->GetBinLowEdge(itmp), itmp, hh->GetNbinsX(), hh->GetBinContent(hh->GetNbinsX())/hh->Integral(0,hh->GetNbinsX()+1), hh->GetBinContent(hh->GetNbinsX()), hh->Integral(0,hh->GetNbinsX()+1)); //exit(1); | |
259 | nused = 0; | |
260 | w2s = sum = 0; | |
261 | } | |
262 | break; | |
263 | } | |
264 | } | |
265 | ||
266 | return GetMeanRMS(nused, sum, w2s, grms, gerr); | |
267 | } | |
268 | ||
269 | void AliTRDdEdxBaseUtils::FitSlicesY(const TH2D *hh, TH1D *&hnor, TH1D *&hmpv, TH1D *&hwid, TH1D *&hres, const Double_t thres, const Double_t lowfrac, const Double_t highfrac) | |
270 | { | |
271 | // | |
272 | //fit slices of hh using truncation | |
273 | // | |
274 | ||
275 | const Int_t x0 = hh->GetXaxis()->GetFirst(); | |
276 | const Int_t x1 = hh->GetXaxis()->GetLast(); | |
277 | const Int_t y0 = hh->GetYaxis()->GetFirst(); | |
278 | const Int_t y1 = hh->GetYaxis()->GetLast(); | |
279 | ||
280 | const Int_t nx = hh->GetNbinsX(); | |
281 | const Int_t ny = hh->GetNbinsY(); | |
282 | const Double_t xmin = hh->GetXaxis()->GetXmin(); | |
283 | const Double_t xmax = hh->GetXaxis()->GetXmax(); | |
284 | const Double_t ymin = hh->GetYaxis()->GetXmin(); | |
285 | const Double_t ymax = hh->GetYaxis()->GetXmax(); | |
286 | ||
287 | hnor = new TH1D(Form("%s_amp",hh->GetName()), "", nx, xmin, xmax); | |
288 | hmpv = new TH1D(Form("%s_mpv",hh->GetName()), "", nx, xmin, xmax); | |
289 | hwid = new TH1D(Form("%s_wid",hh->GetName()), "", nx, xmin, xmax); | |
290 | hres = new TH1D(Form("%s_res",hh->GetName()), "", nx, xmin, xmax); | |
291 | ||
687aa844 | 292 | Double_t newbins[ny+1]; |
293 | for(Int_t ii=1; ii<=ny+1; ii++){ | |
294 | const Double_t lowx= hh->GetYaxis()->GetBinLowEdge(ii); | |
295 | newbins[ii-1]=lowx; | |
296 | } | |
297 | ||
6951a056 | 298 | for(Int_t ix=x0; ix<=x1; ix++){ |
299 | //to speed up | |
300 | const Double_t rawcount = hh->Integral(ix,ix,0, ny+1); | |
301 | if(rawcount<EPSILON){ | |
302 | continue; | |
303 | } | |
304 | ||
305 | TH1D *htmp = new TH1D(Form("FitSlicesY_%s_%d", hh->GetName(), ix),"",ny, ymin, ymax); | |
687aa844 | 306 | htmp->GetXaxis()->Set(ny, newbins); |
6951a056 | 307 | Double_t ntot = 0; |
308 | for(Int_t iy=y0; iy<=y1; iy++){ | |
309 | const Double_t be = hh->GetBinError(ix,iy); | |
310 | const Double_t bc = hh->GetBinContent(ix, iy); | |
311 | ||
312 | if(be<EPSILON){ | |
313 | if(bc>EPSILON){ | |
314 | printf("AliTRDdEdxBaseUtils::FitSlicesY error %d %d %e %e\n", ix, iy, be, bc); exit(1); | |
315 | } | |
316 | continue; | |
317 | } | |
318 | ||
319 | htmp->SetBinContent(iy, bc); | |
320 | htmp->SetBinError(iy, be); | |
321 | ||
322 | ntot += (bc/be)*(bc/be); | |
323 | ||
324 | //if(be) printf("test %d %d : %f %f %f\n", ix, iy, bc, be, pow(bc/be,2)); | |
325 | } | |
326 | ||
327 | hnor->SetBinContent(ix, ntot); | |
328 | hnor->SetBinError( ix, 0); | |
329 | ||
330 | if(ntot<thres || htmp->GetRMS()<EPSILON){ | |
331 | delete htmp; | |
332 | continue; | |
333 | } | |
334 | ||
335 | //test htmp->Draw(); | |
336 | Double_t trms = -999, terr = -999; | |
337 | const Double_t tmean = TruncatedMean(htmp, lowfrac, highfrac, &trms, &terr); | |
338 | ||
339 | hmpv->SetBinContent(ix, tmean); | |
340 | hmpv->SetBinError( ix, terr); | |
341 | ||
342 | hwid->SetBinContent(ix, trms); | |
343 | hwid->SetBinError( ix, 0); | |
344 | ||
345 | hres->SetBinContent(ix, tmean>EPSILON ? trms/tmean:0); | |
346 | hres->SetBinError( ix, 0); | |
347 | ||
348 | delete htmp; | |
349 | } | |
350 | ||
351 | TH1 *hhs[]={hnor, hmpv, hwid, hres}; | |
352 | const TString yt[]={"N", "MPV", "#sigma", "#sigma/MPV"}; | |
353 | const Int_t nh = sizeof(hhs)/sizeof(TH1*); | |
354 | for(Int_t ii=0; ii<nh; ii++){ | |
355 | hhs[ii]->SetYTitle(Form("%s of %s", yt[ii].Data(), hh->GetYaxis()->GetTitle())); | |
356 | hhs[ii]->SetXTitle(hh->GetXaxis()->GetTitle()); | |
357 | hhs[ii]->GetYaxis()->SetTitleOffset(hh->GetYaxis()->GetTitleOffset()); | |
358 | hhs[ii]->SetTitle(hh->GetTitle()); | |
359 | } | |
360 | } | |
361 | ||
362 | //=================================================================================== | |
363 | // TRD Analysis Fast Tool | |
364 | //=================================================================================== | |
365 | ||
366 | Int_t AliTRDdEdxBaseUtils::GetNtracklet(const AliESDEvent *esd) | |
367 | { | |
368 | // | |
369 | //number of trd tracklet in one esd event | |
370 | // | |
371 | const Int_t ntrk0 = esd->GetNumberOfTracks(); | |
372 | Int_t ntrdv1=0; | |
373 | for(Int_t ii=0; ii<ntrk0; ii++){ | |
374 | ntrdv1 += esd->GetTrack(ii)->GetTRDntracklets(); | |
375 | } | |
376 | return ntrdv1; | |
377 | } | |
378 | ||
379 | AliTRDtrackV1 * AliTRDdEdxBaseUtils::GetTRDtrackV1(const AliESDtrack * esdtrack) | |
380 | { | |
381 | // | |
382 | //Get TRD friend track | |
383 | // | |
384 | ||
385 | AliESDfriendTrack * friendtrk = (AliESDfriendTrack *)esdtrack->GetFriendTrack(); | |
386 | if(!friendtrk){ | |
387 | //printf("xlulog AliAnalysisTaskCosmicTRD::GetTRDtrack no friend!!\n"); exit(1); | |
388 | return 0x0; | |
389 | } | |
390 | ||
391 | TObject *calibObject=0x0; | |
392 | AliTRDtrackV1 * trdtrack=0x0; | |
393 | for(Int_t l=0; (calibObject=friendtrk->GetCalibObject(l)); l++) { | |
394 | if( (trdtrack=dynamic_cast<AliTRDtrackV1*>(calibObject)) ) | |
395 | break; | |
396 | } | |
397 | ||
398 | return trdtrack; | |
399 | } | |
400 | ||
401 | Bool_t AliTRDdEdxBaseUtils::IsInSameStack(const AliTRDtrackV1 *trdtrack) | |
402 | { | |
403 | // | |
404 | // to check if all tracklets are in the same stack, useful in cosmic | |
405 | // | |
406 | ||
407 | TVectorD secs(18), stks(5); | |
408 | ||
409 | for(Int_t ilayer = 0; ilayer < 6; ilayer++){ | |
410 | AliTRDseedV1 *tracklet=trdtrack->GetTracklet(ilayer); | |
411 | if(!tracklet) | |
412 | continue; | |
413 | ||
414 | const Int_t det = tracklet->GetDetector(); | |
415 | const Int_t isector = AliTRDgeometry::GetSector(det); | |
416 | const Int_t istack = AliTRDgeometry::GetStack(det); | |
417 | ||
418 | secs[isector] = 1; | |
419 | stks[istack] = 1; | |
420 | } | |
421 | ||
422 | if(secs.Sum()!=1 || stks.Sum()!=1){ | |
423 | return kFALSE; | |
424 | } | |
425 | else | |
426 | return kTRUE; | |
427 | } | |
428 | ||
fe829aa0 | 429 | Double_t AliTRDdEdxBaseUtils::GetDeltaPhi(const AliTRDseedV1 *tracklet) |
430 | { | |
431 | // | |
432 | //get phi difference to normal incidence | |
433 | // | |
434 | if(tracklet) | |
435 | return TMath::ATan(tracklet->GetYref(1)); | |
436 | else | |
437 | return -999; | |
438 | } | |
439 | ||
440 | AliTRDseedV1 * AliTRDdEdxBaseUtils::GetLastTracklet(const AliTRDtrackV1 *trdtrack) | |
441 | { | |
442 | // | |
443 | //get last tracklet | |
444 | // | |
445 | AliTRDseedV1 *tracklet = 0x0; | |
446 | ||
447 | for(Int_t ilayer = 5; ilayer >= 0; ilayer--){ | |
448 | tracklet=trdtrack->GetTracklet(ilayer); | |
449 | if(!tracklet) | |
450 | continue; | |
451 | ||
452 | break; | |
453 | } | |
454 | ||
455 | return tracklet; | |
456 | } | |
457 | ||
458 | AliTRDseedV1 * AliTRDdEdxBaseUtils::GetFirstSectorStackMomentum(const AliTRDtrackV1 *trdtrack, Int_t & isec, Int_t & istk, Double_t & mom) | |
6951a056 | 459 | { |
460 | // | |
461 | //as function name | |
462 | // | |
463 | isec = istk = -999; | |
464 | mom = -999; | |
465 | ||
fe829aa0 | 466 | AliTRDseedV1 *tracklet = 0x0; |
6951a056 | 467 | for(Int_t ilayer = 0; ilayer < 6; ilayer++){ |
fe829aa0 | 468 | tracklet=trdtrack->GetTracklet(ilayer); |
6951a056 | 469 | if(!tracklet) |
470 | continue; | |
471 | ||
472 | const Int_t det = tracklet->GetDetector(); | |
473 | isec = AliTRDgeometry::GetSector(det); | |
474 | istk = AliTRDgeometry::GetStack(det); | |
475 | ||
476 | mom = tracklet->GetMomentum(); | |
477 | ||
478 | break; | |
479 | } | |
480 | ||
fe829aa0 | 481 | return tracklet; |
6951a056 | 482 | } |
483 | ||
484 | //=================================================================================== | |
485 | // Detector and Data Constant | |
486 | //=================================================================================== | |
487 | Int_t AliTRDdEdxBaseUtils::ToDetector(const Int_t gtb) | |
488 | { | |
489 | // | |
490 | //gtb = det*Ntb+itb | |
491 | // | |
492 | return gtb/AliTRDseedV1::kNtb; | |
493 | } | |
494 | ||
495 | Int_t AliTRDdEdxBaseUtils::ToTimeBin(const Int_t gtb) | |
496 | { | |
497 | // | |
498 | //gtb = det*Ntb+itb | |
499 | // | |
500 | return gtb%AliTRDseedV1::kNtb; | |
501 | } | |
502 | ||
503 | Int_t AliTRDdEdxBaseUtils::ToSector(const Int_t gtb) | |
504 | { | |
505 | // | |
506 | //return sector | |
507 | // | |
508 | return AliTRDgeometry::GetSector(ToDetector(gtb)); | |
509 | } | |
510 | ||
511 | Int_t AliTRDdEdxBaseUtils::ToStack(const Int_t gtb) | |
512 | { | |
513 | // | |
514 | //return stack | |
515 | // | |
516 | return AliTRDgeometry::GetStack(ToDetector(gtb)); | |
517 | } | |
518 | ||
519 | Int_t AliTRDdEdxBaseUtils::ToLayer(const Int_t gtb) | |
520 | { | |
521 | // | |
522 | //return layer | |
523 | // | |
524 | return AliTRDgeometry::GetLayer(ToDetector(gtb)); | |
525 | } | |
526 | ||
527 | TString AliTRDdEdxBaseUtils::GetRunType(const Int_t run) | |
528 | { | |
529 | // | |
530 | //return run type | |
531 | // | |
532 | ||
533 | TString type; | |
534 | if(run>=121527 && run<= 126460)//LHC10d | |
535 | type="pp2010LHC10d"; | |
536 | else if(run>=126461 && run<= 130930)//LHC10e | |
537 | type="pp2010LHC10e"; | |
538 | else if(run>=136782 && run <= 139846)//LHC10h | |
539 | type="PbPb2010LHC10h"; | |
540 | else if(run>= 143224 && run<= 143237)//2011Feb | |
541 | type="cosmic2011Feb"; | |
542 | else if(run>= 150587 && run<= 154930){ | |
543 | type="cosmic2011MayJun"; | |
544 | ||
545 | TString runstr(Form("%d",run)); | |
546 | const TString listrun1kg("154601 154602 154629 154634 154636 154639 154643"); | |
547 | if(listrun1kg.Contains(runstr)){ | |
548 | type+="1kG"; | |
549 | } | |
550 | else{ | |
551 | type+="5kG"; | |
552 | } | |
553 | } | |
554 | else{ | |
555 | type="unknown"; | |
556 | } | |
557 | ||
558 | type.ToUpper(); | |
559 | return type; | |
560 | } | |
561 | ||
562 | void AliTRDdEdxBaseUtils::PrintControl() | |
563 | { | |
564 | // | |
565 | //print out control variable | |
566 | // | |
687aa844 | 567 | printf("\nAliTRDdEdxBaseUtils::PrintControl Q0Frac %.2f, Q1Frac %.2f, TimeBinCountCut %.2f, CalibTPCnclsCut %d, IsExBOn %d, IsPadGainOn %d, QScale %.2f\n\n", Q0Frac(), Q1Frac(), TimeBinCountCut(), CalibTPCnclsCut(), IsExBOn(), IsPadGainOn(), QScale()); |
6951a056 | 568 | } |
569 | ||
570 | //=================================================================================== | |
571 | // dEdx Parameterization | |
572 | //=================================================================================== | |
fe829aa0 | 573 | void AliTRDdEdxBaseUtils::FastFitdEdxTR(TH1 * hh) |
574 | { | |
575 | // | |
576 | //fit dedx tr from mean | |
577 | // | |
578 | ||
579 | TF1 *ff=new TF1("ff", MeandEdxTRLogx, -0.5, 4.5, 8); | |
580 | Double_t par[8]; | |
581 | par[0]= 2.397001e-01; | |
582 | par[1]= 1.334697e+00; | |
583 | par[2]= 6.967470e+00; | |
584 | par[3]= 9.055289e-02; | |
585 | par[4]= 9.388760e+00; | |
586 | par[5]= 9.452742e-04; | |
587 | par[6]= -1.866091e+00; | |
588 | par[7]= 1.403545e+00; | |
589 | ||
590 | ff->SetParameters(par); | |
591 | hh->Fit(ff,"N"); | |
592 | ||
593 | for(int ii=0; ii<8; ii++){ | |
594 | printf("par[%d]=%e;\n", ii, ff->GetParameter(ii)); | |
595 | } | |
596 | ||
597 | TFile *fout=new TFile("fastfit.root","recreate"); | |
598 | hh->Write(); | |
599 | ff->Write(); | |
600 | fout->Save(); | |
601 | fout->Close(); | |
602 | ||
603 | delete fout; | |
604 | delete ff; | |
605 | } | |
6951a056 | 606 | |
607 | Double_t AliTRDdEdxBaseUtils::Q0MeanTRDpp(const Double_t bg) | |
608 | { | |
609 | // | |
610 | //truncated Mean Q_{xx} in TRD | |
611 | // | |
612 | ||
613 | Double_t par[8]; | |
fe829aa0 | 614 | Double_t scale = 1; |
615 | ||
616 | //2012-0605 /u/xlu/.task/CommondEdx/myAnaData/newTune_r56595/inuse/plot/fastFit | |
617 | scale = 0.9257; | |
618 | par[0]=8.316476e-01; | |
619 | par[1]=1.600907e+00; | |
620 | par[2]=7.728447e+00; | |
621 | par[3]=6.235622e-02; | |
622 | par[4]=1.136209e+01; | |
623 | par[5]=-1.495764e-06; | |
624 | par[6]=-2.536119e+00; | |
625 | par[7]=3.416031e+00; | |
626 | ||
627 | return scale*MeandEdxTR(&bg, par); | |
6951a056 | 628 | } |
629 | ||
630 | Double_t AliTRDdEdxBaseUtils::Q0MeanTRDPbPb(const Double_t bg) | |
631 | { | |
632 | // | |
633 | //truncated Mean Q_{xx} in TRD | |
634 | // | |
635 | ||
636 | Double_t par[8]; | |
637 | ||
638 | //03132012161259 | |
639 | //opt: PbPbQ0 | |
640 | par[0]= 1.844912e-01; | |
641 | par[1]= 2.509702e+00; | |
642 | par[2]= 6.744031e+00; | |
643 | par[3]= 7.355123e-02; | |
644 | par[4]= 1.166023e+01; | |
645 | par[5]= 1.736186e-04; | |
646 | par[6]= -1.716063e+00; | |
647 | par[7]= 1.611366e+00; | |
648 | ||
649 | ///u/xlu/.task/CommondEdx/myAnaData/Optimum/check11/Merged/LHC10e_plot/Separation/see4.log:hhtype4Q0b2c2 scale 0.460994 at ltbg 0.650000 | |
650 | return 0.460994*MeandEdxTR(&bg, par); | |
651 | } | |
652 | ||
653 | Double_t AliTRDdEdxBaseUtils::Q1MeanTRDpp(const Double_t bg) | |
654 | { | |
655 | // | |
656 | //truncated Mean 1/(1/Q)_{xx} in TRD | |
657 | // | |
658 | ||
659 | Double_t par[8]; | |
660 | ||
661 |