]>
Commit | Line | Data |
---|---|---|
10757ee9 | 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 | //////////////////////////////////////////////////////////////////////////// | |
0d3279d4 | 17 | // |
18 | // | |
19 | // Gain calibration using tracks | |
20 | // | |
21 | // The main goal: | |
22 | // 1.) Inner TPC gain alignement - (parabolic) parameterization (inside of one sector) | |
23 | // 2.) Angular and z-position correction (parabolic) parameterization | |
24 | // 3.) Sector gain alignment | |
25 | // | |
26 | // Following histograms are accumulated | |
27 | // a.) Simple 1D histograms per chamber | |
28 | // b.) Profile histograms per chamber - local x dependence | |
29 | // c.) 2D Profile histograms - local x - fi dependence | |
30 | // | |
31 | // To get the gain map - the simple solution - use the histograms - is not enough | |
32 | // The resulting mean amplitude map depends strongly on the track topology | |
33 | // These dependence can be reduced, taking into account angular effect, and diffusion effect | |
34 | // Using proper fit modeles | |
35 | // | |
36 | // | |
37 | // | |
10757ee9 | 38 | // |
39 | // === Calibration class for gain calibration using tracks === | |
40 | // | |
41 | // 1) Genereal idea | |
42 | // ================ | |
43 | // A 6-parametric parabolic function | |
44 | // | |
45 | // G(x, y) = p0 + p1*x + p2*y + p3*x^2 + p4*y^2 + p5 * x*y | |
46 | // | |
47 | // is fitted to the maximum charge values or total charge values of | |
48 | // all the clusters contained in the tracks that are added to this | |
49 | // object. This fit is performed for each read out chamber, in fact even | |
50 | // for each type of pad sizes (thus for one segment, which consists of | |
51 | // an IROC and an OROC, there are three fitters used, corresponding to | |
52 | // the three pad sizes). The coordinate origin is at the center of the | |
53 | // particular pad size region on each ROC. | |
54 | // | |
55 | // Because of the Landau nature of the charge deposition we use | |
56 | // different "types" of fitters instead of one to minimize the effect | |
57 | // of the long Landau tail. The difference between the fitters is only | |
58 | // the charge value, that is put into them, i.e. the charge is subject | |
59 | // to a transformation. At this point we use three different fit types: | |
60 | // | |
61 | // a) simple: the charge is put in as it is | |
62 | // b) sqrt: the square root of the charge is put into the fitter | |
63 | // c) log: fgkM * Log(1+q/fgkM) is put into the fitter, with | |
64 | // q being the untransformed charge and fgkM=25 | |
65 | // | |
66 | // The results of the fits may be visualized and further used by | |
67 | // creating an AliTPCCalROC or AliTPCCalPad. You may specify to undo | |
68 | // the transformation and/or to normalize to the pad size. | |
69 | // | |
70 | // Not every track you add to this object is actually used for | |
71 | // calibration. There are some cuts and conditions to exclude bad | |
72 | // tracks, e.g. a pt cut to cut out tracks with too much charge | |
73 | // deposition or a cut on edge clusters which are not fully | |
74 | // registered and don't give a usable signal. | |
75 | // | |
76 | // 2) Interface / usage | |
77 | // ==================== | |
78 | // For each track to be added you need to call Process(). | |
79 | // This method expects an AliTPCseed, which contains the necessary | |
80 | // cluster information. At the moment of writing this information | |
81 | // is stored in an AliESDfriend corresponding to an AliESD. | |
82 | // You may also call AddTrack() if you don't want the cuts and | |
83 | // other quality conditions to kick in (thus forcing the object to | |
84 | // accept the track) or AddCluster() for adding single clusters. | |
85 | // Call one of the Evaluate functions to evaluate the fitter(s) and | |
86 | // to retrieve the fit parameters, erros and so on. You can also | |
87 | // do this later on by using the different Getters. | |
88 | // | |
89 | // The visualization methods CreateFitCalPad() and CreateFitCalROC() | |
90 | // are straight forward to use. | |
91 | // | |
92 | // Note: If you plan to write this object to a ROOT file, make sure | |
93 | // you evaluate all the fitters *before* writing, because due | |
94 | // to a bug in the fitter component writing fitters doesn't | |
95 | // work properly (yet). Be aware that you cannot re-evaluate | |
96 | // the fitters after loading this object from file. | |
97 | // (This will be gone for a new ROOT version > v5-17-05) | |
98 | // | |
684602c8 | 99 | // |
100 | // In order to debug some numerical algorithm all data data which are used for | |
101 | // fitters can be stored in the debug streamers. In case of fitting roblems the | |
102 | // errors and tendencies can be checked. | |
103 | // | |
104 | // Debug Streamers: | |
105 | // | |
106 | // | |
107 | // | |
108 | // | |
109 | // | |
10757ee9 | 110 | //////////////////////////////////////////////////////////////////////////// |
111 | ||
684602c8 | 112 | /* |
113 | .x ~/UliStyle.C | |
114 | gSystem->Load("libANALYSIS"); | |
da7d274e | 115 | gSystem->Load("libSTAT"); |
684602c8 | 116 | gSystem->Load("libTPCcalib"); |
da7d274e | 117 | |
684602c8 | 118 | TFile fcalib("CalibObjects.root"); |
119 | TObjArray * array = (TObjArray*)fcalib.Get("TPCCalib"); | |
120 | AliTPCcalibTracksGain * gain = ( AliTPCcalibTracksGain *)array->FindObject("calibTracksGain"); | |
121 | ||
2acad464 | 122 | // |
123 | // Angular and drift correction | |
124 | // | |
125 | AliTPCClusterParam *param = new AliTPCClusterParam;param->SetInstance(param); | |
126 | gain->UpdateClusterParam(param); | |
127 | TF2 fdrty("fdrty","AliTPCClusterParam::SQnorm(0,0,x,y,0)",0,1,0,1) | |
128 | ||
684602c8 | 129 | // |
130 | // Make visual Tree - compare with Kr calibration | |
131 | // | |
132 | AliTPCCalPad * gain010 = gain->CreateFitCalPad(0,kTRUE,0); gain010->SetName("CGain010"); | |
133 | AliTPCCalPad * gain110 = gain->CreateFitCalPad(1,kTRUE,0); gain110->SetName("CGain110"); | |
134 | AliTPCCalPad * gain210 = gain->CreateFitCalPad(2,kTRUE,0); gain210->SetName("CGain210"); | |
2acad464 | 135 | TFile fkr("/u/miranov/GainMap.root"); |
684602c8 | 136 | AliTPCCalPad *gainKr = fkr.Get("GainMap"); fkr->SetName("KrGain"); |
137 | // | |
138 | AliTPCPreprocessorOnline * preprocesor = new AliTPCPreprocessorOnline; | |
139 | preprocesor->AddComponent(gain010); | |
140 | preprocesor->AddComponent(gain110); | |
141 | preprocesor->AddComponent(gain210); | |
142 | preprocesor->AddComponent(gainKr); | |
143 | preprocesor->DumpToFile("cosmicGain.root"); | |
144 | // | |
145 | // | |
146 | // | |
147 | // Simple session using the debug streamers | |
148 | // | |
149 | ||
150 | gSystem->AddIncludePath("-I$ALICE_ROOT/TPC/macros"); | |
151 | gROOT->LoadMacro("$ALICE_ROOT/TPC/macros/AliXRDPROOFtoolkit.cxx+") | |
152 | AliXRDPROOFtoolkit tool; | |
153 | ||
154 | TChain * chain0 = tool.MakeChain("gain.txt","dEdx",0,1000000); | |
155 | TChain * chain1 = tool.MakeChain("gain.txt","Track",0,1000000); | |
156 | TChain * chain2 = tool.MakeChain("gain.txt","TrackG",0,1000000); | |
157 | chain0->Lookup(); | |
158 | chain1->Lookup(); | |
159 | chain2->Lookup(); | |
160 | ||
161 | chain2->SetAlias("k1","1/0.855"); | |
162 | chain2->SetAlias("k0","1/0.9928"); | |
163 | chain2->SetAlias("k2","1/1.152"); | |
164 | ||
684602c8 | 165 | */ |
166 | ||
167 | ||
168 | ||
169 | ||
170 | ||
10757ee9 | 171 | #include <TPDGCode.h> |
172 | #include <TStyle.h> | |
173 | #include "TSystem.h" | |
174 | #include "TMatrixD.h" | |
175 | #include "TTreeStream.h" | |
176 | #include "TF1.h" | |
177 | #include "AliTPCParamSR.h" | |
178 | #include "AliTPCClusterParam.h" | |
179 | #include "AliTrackPointArray.h" | |
180 | #include "TCint.h" | |
181 | #include "AliTPCcalibTracksGain.h" | |
182 | #include <TH1.h> | |
183 | #include <TH3F.h> | |
184 | #include <TLinearFitter.h> | |
185 | #include <TTreeStream.h> | |
186 | #include <TFile.h> | |
187 | #include <TCollection.h> | |
188 | #include <TIterator.h> | |
5b00528f | 189 | #include <TProfile.h> |
190 | #include <TProfile2D.h> | |
c1418a4c | 191 | #include <TProof.h> |
7f37d7e4 | 192 | #include <TStatToolkit.h> |
10757ee9 | 193 | |
194 | // | |
195 | // AliRoot includes | |
196 | // | |
197 | #include "AliMagF.h" | |
198 | #include "AliMathBase.h" | |
199 | // | |
200 | #include "AliTPCROC.h" | |
201 | #include "AliTPCParamSR.h" | |
202 | #include "AliTPCCalROC.h" | |
203 | #include "AliTPCCalPad.h" | |
8076baa0 | 204 | #include "AliTPCClusterParam.h" |
10757ee9 | 205 | // |
206 | #include "AliTracker.h" | |
207 | #include "AliESD.h" | |
208 | #include "AliESDtrack.h" | |
209 | #include "AliESDfriend.h" | |
210 | #include "AliESDfriendTrack.h" | |
211 | #include "AliTPCseed.h" | |
212 | #include "AliTPCclusterMI.h" | |
213 | #include "AliTPCcalibTracksCuts.h" | |
214 | #include "AliTPCFitPad.h" | |
da7d274e | 215 | #include "TStatToolkit.h" |
216 | #include "TString.h" | |
217 | #include "TCut.h" | |
10757ee9 | 218 | |
da7d274e | 219 | // |
10757ee9 | 220 | #include <TTree.h> |
221 | #include "AliESDEvent.h" | |
222 | ||
223 | /* | |
224 | ||
225 | TFile f("TPCCalibTracksGain.root") | |
226 | ||
227 | gSystem->Load("libPWG1.so") | |
228 | AliTreeDraw comp | |
229 | comp.SetTree(dEdx) | |
230 | Double_t chi2; | |
231 | TVectorD vec(3) | |
232 | TMatrixD mat(3,3) | |
233 | TString * str = comp.FitPlane("Cl.fQ/dedxQ.fElements[0]","Cl.fY++Cl.fX","Cl.fDetector<36",chi2,vec,mat) | |
234 | ||
235 | */ | |
236 | ||
237 | ClassImp(AliTPCcalibTracksGain) | |
238 | ||
239 | const Bool_t AliTPCcalibTracksGain::fgkUseTotalCharge = kTRUE; | |
240 | const Double_t AliTPCcalibTracksGain::fgkM = 25.; | |
241 | const char* AliTPCcalibTracksGain::fgkDebugStreamFileName = "TPCCalibTracksGain.root"; | |
242 | AliTPCParamSR* AliTPCcalibTracksGain::fgTPCparam = new AliTPCParamSR(); | |
243 | ||
244 | AliTPCcalibTracksGain::AliTPCcalibTracksGain() : | |
b8601924 | 245 | AliTPCcalibBase(), |
0d3279d4 | 246 | fCuts(0), // cuts that are used for sieving the tracks used for calibration |
2acad464 | 247 | fGainMap(0), // gain map to be applied |
0d3279d4 | 248 | // |
249 | // Simple Array of histograms | |
250 | // | |
251 | fArrayQM(0), // Qmax normalized | |
252 | fArrayQT(0), // Qtot normalized | |
253 | fProfileArrayQM(0), // Qmax normalized versus local X | |
254 | fProfileArrayQT(0), // Qtot normalized versus local X | |
255 | fProfileArrayQM2D(0), // Qmax normalized versus local X and phi | |
256 | fProfileArrayQT2D(0), // Qtot normalized versus local X and phi | |
257 | // | |
258 | // Fitters | |
259 | // | |
260 | fSimpleFitter(0), // simple fitter for short pads | |
261 | fSqrtFitter(0), // sqrt fitter for medium pads | |
262 | fLogFitter(0), // log fitter for long pads | |
684602c8 | 263 | |
0d3279d4 | 264 | fFitter0M(0), // fitting of the atenuation, angular correction, and mean chamber gain |
265 | fFitter1M(0), // fitting of the atenuation, angular correction, and mean chamber gain | |
266 | fFitter2M(0), // fitting of the atenuation, angular correction, and mean chamber gain | |
267 | fFitter0T(0), // fitting of the atenuation, angular correction, and mean chamber gain | |
268 | fFitter1T(0), // fitting of the atenuation, angular correction, and mean chamber gain | |
269 | fFitter2T(0), // fitting of the atenuation, angular correction, and mean chamber gain | |
684602c8 | 270 | // |
271 | fDFitter0M(0), // fitting of the atenuation, angular correction | |
272 | fDFitter1M(0), // fitting of the atenuation, angular correction | |
273 | fDFitter2M(0), // fitting of the atenuation, angular correction | |
274 | fDFitter0T(0), // fitting of the atenuation, angular correction | |
275 | fDFitter1T(0), // fitting of the atenuation, angular correction | |
276 | fDFitter2T(0), // fitting of the atenuation, angular correction | |
277 | // | |
0d3279d4 | 278 | fSingleSectorFitter(0), // just for debugging |
279 | // | |
280 | // Counters | |
281 | // | |
282 | fTotalTracks(0), // just for debugging | |
283 | fAcceptedTracks(0), // just for debugging | |
284 | fDebugCalPadRaw(0), // just for debugging | |
2acad464 | 285 | fDebugCalPadCorr(0) // just for debugging |
0d3279d4 | 286 | |
10757ee9 | 287 | { |
288 | // | |
289 | // Default constructor. | |
290 | // | |
291 | } | |
292 | ||
293 | AliTPCcalibTracksGain::AliTPCcalibTracksGain(const AliTPCcalibTracksGain& obj) : | |
b8601924 | 294 | AliTPCcalibBase(obj), |
0d3279d4 | 295 | fCuts(obj.fCuts), // cuts that are used for sieving the tracks used for calibration |
2acad464 | 296 | fGainMap(new AliTPCCalPad(*(obj.fGainMap))), // gain map to be applied |
f1c2a4a3 | 297 | fArrayQM(0), // Qmax normalized |
298 | fArrayQT(0), // Qtot normalized | |
0d3279d4 | 299 | // |
300 | // Simple Histograms | |
301 | // | |
302 | fProfileArrayQM(obj.fProfileArrayQM), // Qmax normalized versus local X | |
303 | fProfileArrayQT(obj.fProfileArrayQT), // Qtot normalized versus local X | |
304 | fProfileArrayQM2D(obj.fProfileArrayQM2D), // Qmax normalized versus local X and phi | |
305 | fProfileArrayQT2D(obj.fProfileArrayQT2D), // Qtot normalized versus local X and phi | |
306 | // | |
307 | // Fitters | |
308 | // | |
309 | fSimpleFitter(obj.fSimpleFitter), // simple fitter for short pads | |
310 | fSqrtFitter(obj.fSqrtFitter), // sqrt fitter for medium pads | |
311 | fLogFitter(obj.fLogFitter), // log fitter for long pads | |
312 | fFitter0M(obj.fFitter0M), | |
313 | fFitter1M(obj.fFitter1M), | |
314 | fFitter2M(obj.fFitter2M), | |
315 | fFitter0T(obj.fFitter0T), | |
316 | fFitter1T(obj.fFitter1T), | |
317 | fFitter2T(obj.fFitter2T), | |
684602c8 | 318 | // |
319 | fDFitter0M(obj.fDFitter0M), | |
320 | fDFitter1M(obj.fDFitter1M), | |
321 | fDFitter2M(obj.fDFitter2M), | |
322 | fDFitter0T(obj.fDFitter0T), | |
323 | fDFitter1T(obj.fDFitter1T), | |
324 | fDFitter2T(obj.fDFitter2T), | |
0d3279d4 | 325 | fSingleSectorFitter(obj.fSingleSectorFitter), // just for debugging |
326 | // | |
327 | // Counters | |
328 | // | |
329 | fTotalTracks(obj.fTotalTracks), // just for debugging | |
330 | fAcceptedTracks(obj.fAcceptedTracks), // just for debugging | |
331 | fDebugCalPadRaw(obj.fDebugCalPadRaw), // just for debugging | |
2acad464 | 332 | fDebugCalPadCorr(obj.fDebugCalPadCorr) // just for debugging |
0d3279d4 | 333 | |
10757ee9 | 334 | { |
335 | // | |
336 | // Copy constructor. | |
337 | // | |
10757ee9 | 338 | } |
339 | ||
340 | AliTPCcalibTracksGain& AliTPCcalibTracksGain::operator=(const AliTPCcalibTracksGain& rhs) { | |
341 | // | |
342 | // Assignment operator. | |
343 | // | |
344 | ||
345 | if (this != &rhs) { | |
346 | TNamed::operator=(rhs); | |
347 | fDebugCalPadRaw = new AliTPCCalPad(*(rhs.fDebugCalPadRaw)); | |
348 | fDebugCalPadCorr = new AliTPCCalPad(*(rhs.fDebugCalPadCorr)); | |
349 | fSimpleFitter = new AliTPCFitPad(*(rhs.fSimpleFitter)); | |
350 | fSqrtFitter = new AliTPCFitPad(*(rhs.fSqrtFitter)); | |
351 | fLogFitter = new AliTPCFitPad(*(rhs.fLogFitter)); | |
352 | fSingleSectorFitter = new AliTPCFitPad(*(rhs.fSingleSectorFitter)); | |
10757ee9 | 353 | fCuts = new AliTPCcalibTracksCuts(*(rhs.fCuts)); |
2acad464 | 354 | fGainMap = new AliTPCCalPad(*(rhs.fGainMap)); |
10757ee9 | 355 | } |
356 | return *this; | |
357 | } | |
358 | ||
f1c2a4a3 | 359 | AliTPCcalibTracksGain::AliTPCcalibTracksGain(const char* name, const char* title, AliTPCcalibTracksCuts* cuts, TNamed* /*debugStreamPrefix*/, AliTPCcalibTracksGain* prevIter) : |
b8601924 | 360 | AliTPCcalibBase(), |
0d3279d4 | 361 | fCuts(0), // cuts that are used for sieving the tracks used for calibration |
2acad464 | 362 | fGainMap(0), // gain map to be applied |
f1c2a4a3 | 363 | fArrayQM(0), // Qmax normalized |
364 | fArrayQT(0), // Qtot normalized | |
0d3279d4 | 365 | // |
366 | // Simple Histograms | |
367 | // | |
368 | fProfileArrayQM(0), // Qmax normalized versus local X | |
369 | fProfileArrayQT(0), // Qtot normalized versus local X | |
370 | fProfileArrayQM2D(0), // Qmax normalized versus local X and phi | |
371 | fProfileArrayQT2D(0), // Qtot normalized versus local X and phi | |
372 | // | |
373 | // Fitters | |
374 | // | |
375 | fSimpleFitter(0), // simple fitter for short pads | |
376 | fSqrtFitter(0), // sqrt fitter for medium pads | |
377 | fLogFitter(0), // log fitter for long pads | |
378 | fFitter0M(0), // fitting of the atenuation, angular correction, and mean chamber gain | |
379 | fFitter1M(0), // fitting of the atenuation, angular correction, and mean chamber gain | |
380 | fFitter2M(0), // fitting of the atenuation, angular correction, and mean chamber gain | |
381 | fFitter0T(0), // fitting of the atenuation, angular correction, and mean chamber gain | |
382 | fFitter1T(0), // fitting of the atenuation, angular correction, and mean chamber gain | |
383 | fFitter2T(0), // fitting of the atenuation, angular correction, and mean chamber gain | |
684602c8 | 384 | // |
385 | fDFitter0M(0), // fitting of the atenuation, angular correction | |
386 | fDFitter1M(0), // fitting of the atenuation, angular correction | |
387 | fDFitter2M(0), // fitting of the atenuation, angular correction | |
388 | fDFitter0T(0), // fitting of the atenuation, angular correction | |
389 | fDFitter1T(0), // fitting of the atenuation, angular correction | |
390 | fDFitter2T(0), // fitting of the atenuation, angular correction | |
0d3279d4 | 391 | fSingleSectorFitter(0), // just for debugging |
392 | // | |
393 | // Counters | |
394 | // | |
395 | fTotalTracks(0), // just for debugging | |
396 | fAcceptedTracks(0), // just for debugging | |
397 | fDebugCalPadRaw(0), // just for debugging | |
2acad464 | 398 | fDebugCalPadCorr(0) // just for debugging |
0d3279d4 | 399 | |
10757ee9 | 400 | { |
401 | // | |
402 | // Constructor. | |
0d3279d4 | 403 | // |
b8601924 | 404 | this->SetNameTitle(name, title); |
10757ee9 | 405 | fCuts = cuts; |
0d3279d4 | 406 | // |
407 | // Fitter initialization | |
408 | // | |
10757ee9 | 409 | fSimpleFitter = new AliTPCFitPad(8, "hyp7", ""); |
410 | fSqrtFitter = new AliTPCFitPad(8, "hyp7", ""); | |
411 | fLogFitter = new AliTPCFitPad(8, "hyp7", ""); | |
412 | fSingleSectorFitter = new AliTPCFitPad(8, "hyp7", ""); | |
8076baa0 | 413 | // |
414 | fFitter0M = new TLinearFitter(45,"hyp44"); | |
415 | fFitter1M = new TLinearFitter(45,"hyp44"); | |
416 | fFitter2M = new TLinearFitter(45,"hyp44"); | |
417 | fFitter0T = new TLinearFitter(45,"hyp44"); | |
418 | fFitter1T = new TLinearFitter(45,"hyp44"); | |
419 | fFitter2T = new TLinearFitter(45,"hyp44"); | |
0d3279d4 | 420 | // |
684602c8 | 421 | fDFitter0M = new TLinearFitter(10,"hyp9"); |
422 | fDFitter1M = new TLinearFitter(10,"hyp9"); | |
423 | fDFitter2M = new TLinearFitter(10,"hyp9"); | |
424 | fDFitter0T = new TLinearFitter(10,"hyp9"); | |
425 | fDFitter1T = new TLinearFitter(10,"hyp9"); | |
426 | fDFitter2T = new TLinearFitter(10,"hyp9"); | |
427 | // | |
cbc19295 | 428 | // |
429 | fFitter0M->StoreData(kFALSE); | |
430 | fFitter1M->StoreData(kFALSE); | |
431 | fFitter2M->StoreData(kFALSE); | |
432 | fFitter0T->StoreData(kFALSE); | |
433 | fFitter1T->StoreData(kFALSE); | |
434 | fFitter2T->StoreData(kFALSE); | |
435 | // | |
684602c8 | 436 | fDFitter0M->StoreData(kFALSE); |
437 | fDFitter1M->StoreData(kFALSE); | |
438 | fDFitter2M->StoreData(kFALSE); | |
439 | fDFitter0T->StoreData(kFALSE); | |
440 | fDFitter1T->StoreData(kFALSE); | |
441 | fDFitter2T->StoreData(kFALSE); | |
442 | // | |
0d3279d4 | 443 | // |
444 | // Add profile histograms -JUST for visualization - Not used for real calibration | |
445 | // | |
446 | // | |
447 | fArrayQM=new TObjArray(73); // Qmax normalized | |
448 | fArrayQT=new TObjArray(73); // Qtot normalized | |
449 | fProfileArrayQM = new TObjArray(37); // Qmax normalized versus local X | |
450 | fProfileArrayQT = new TObjArray(37); // Qtot normalized versus local X | |
451 | fProfileArrayQM2D = new TObjArray(37); // Qmax normalized versus local X and phi | |
452 | fProfileArrayQT2D = new TObjArray(37); // Qtot normalized versus local X and phi | |
453 | char hname[1000]; | |
454 | for (Int_t i=0; i<73; i++){ | |
455 | sprintf(hname,"QM_%d",i); | |
456 | fArrayQM->AddAt(new TH1F(hname,hname,200,0,1000),i); | |
457 | sprintf(hname,"QT_%d",i); | |
458 | fArrayQT->AddAt(new TH1F(hname,hname,200,0,1000),i); | |
459 | } | |
10757ee9 | 460 | |
0d3279d4 | 461 | for (Int_t i=0; i<37;i++){ |
462 | sprintf(hname,"QMvsx_%d",i); | |
463 | fProfileArrayQM->AddAt(new TProfile(hname,hname,50,89,250),i); | |
464 | sprintf(hname,"QTvsx_%d",i); | |
465 | fProfileArrayQT->AddAt(new TProfile(hname,hname,50,89,250),i); | |
466 | sprintf(hname,"QM2D_%d",i); | |
467 | fProfileArrayQM2D->AddAt(new TProfile2D(hname,hname,50,89,250,10,-0.15,0.15),i); | |
468 | sprintf(hname,"QT2D_%d",i); | |
469 | fProfileArrayQT2D->AddAt(new TProfile2D(hname,hname,50,89,250,10,-0.15,0.15),i); | |
470 | } | |
471 | // | |
472 | // just for debugging -counters | |
473 | // | |
10757ee9 | 474 | fTotalTracks = 0; |
475 | fAcceptedTracks = 0; | |
476 | fDebugCalPadRaw = new AliTPCCalPad("DebugCalPadRaw", "All clusters simply added up before correction"); | |
0d3279d4 | 477 | fDebugCalPadCorr = new AliTPCCalPad("DebugCalPadCorr", "All clusters simply added up after correction"); |
10757ee9 | 478 | // this will be gone for the a new ROOT version > v5-17-05 |
479 | for (UInt_t i = 0; i < 36; i++) { | |
480 | fNShortClusters[i] = 0; | |
481 | fNMediumClusters[i] = 0; | |
482 | fNLongClusters[i] = 0; | |
483 | } | |
484 | } | |
485 | ||
486 | AliTPCcalibTracksGain::~AliTPCcalibTracksGain() { | |
487 | // | |
488 | // Destructor. | |
489 | // | |
490 | ||
491 | Info("Destructor",""); | |
492 | if (fSimpleFitter) delete fSimpleFitter; | |
493 | if (fSqrtFitter) delete fSqrtFitter; | |
494 | if (fLogFitter) delete fLogFitter; | |
495 | if (fSingleSectorFitter) delete fSingleSectorFitter; | |
496 | ||
10757ee9 | 497 | if (fDebugCalPadRaw) delete fDebugCalPadRaw; |
498 | if (fDebugCalPadCorr) delete fDebugCalPadCorr; | |
499 | } | |
500 | ||
501 | void AliTPCcalibTracksGain::Terminate(){ | |
502 | // | |
503 | // Evaluate fitters and close the debug stream. | |
504 | // Also move or copy the debug stream, if a debugStreamPrefix is provided. | |
505 | // | |
506 | ||
507 | Evaluate(); | |
ae28e92e | 508 | AliTPCcalibBase::Terminate(); |
10757ee9 | 509 | } |
510 | ||
10757ee9 | 511 | |
10757ee9 | 512 | |
513 | void AliTPCcalibTracksGain::Process(AliTPCseed* seed) { | |
514 | // | |
515 | // Main method to be called when a new seed is supposed to be processed | |
516 | // and be used for gain calibration. Its quality is checked before it | |
517 | // is added. | |
518 | // | |
519 | ||
c1418a4c | 520 | |
10757ee9 | 521 | fTotalTracks++; |
b8601924 | 522 | if (!fCuts->AcceptTrack(seed)) return; |
c1418a4c | 523 | // |
524 | // reinint on proof | |
525 | // if (gProof){ | |
526 | static Bool_t doinit= kTRUE; | |
527 | if (doinit){ | |
528 | fSimpleFitter = new AliTPCFitPad(8, "hyp7", ""); | |
529 | fSqrtFitter = new AliTPCFitPad(8, "hyp7", ""); | |
530 | fLogFitter = new AliTPCFitPad(8, "hyp7", ""); | |
531 | fSingleSectorFitter = new AliTPCFitPad(8, "hyp7", ""); | |
532 | // | |
533 | fFitter0M = new TLinearFitter(45,"hyp44"); | |
534 | fFitter1M = new TLinearFitter(45,"hyp44"); | |
535 | fFitter2M = new TLinearFitter(45,"hyp44"); | |
536 | fFitter0T = new TLinearFitter(45,"hyp44"); | |
537 | fFitter1T = new TLinearFitter(45,"hyp44"); | |
538 | fFitter2T = new TLinearFitter(45,"hyp44"); | |
684602c8 | 539 | // |
540 | fDFitter0M = new TLinearFitter(10,"hyp9"); | |
541 | fDFitter1M = new TLinearFitter(10,"hyp9"); | |
542 | fDFitter2M = new TLinearFitter(10,"hyp9"); | |
543 | fDFitter0T = new TLinearFitter(10,"hyp9"); | |
544 | fDFitter1T = new TLinearFitter(10,"hyp9"); | |
545 | fDFitter2T = new TLinearFitter(10,"hyp9"); | |
546 | doinit=kFALSE; | |
c1418a4c | 547 | } |
548 | //} | |
549 | ||
10757ee9 | 550 | fAcceptedTracks++; |
551 | AddTrack(seed); | |
552 | } | |
553 | ||
554 | Long64_t AliTPCcalibTracksGain::Merge(TCollection *list) { | |
555 | // | |
556 | // Merge() merges the results of all AliTPCcalibTracksGain objects contained in | |
557 | // list, thus allowing a distributed computation of several files, e.g. on PROOF. | |
558 | // The merged results are merged with the data members of the AliTPCcalibTracksGain | |
559 | // object used for calling the Merge method. | |
560 | // The return value is 0 /*the total number of tracks used for calibration*/ if the merge | |
561 | // is successful, otherwise it is -1. | |
562 | // | |
563 | ||
564 | if (!list || list->IsEmpty()) return -1; | |
565 | ||
566 | if (!fSimpleFitter) fSimpleFitter = new AliTPCFitPad(8, "hyp7", ""); | |
567 | if (!fSqrtFitter) fSqrtFitter = new AliTPCFitPad(8, "hyp7", ""); | |
568 | if (!fLogFitter) fLogFitter = new AliTPCFitPad(8, "hyp7", ""); | |
569 | if (!fSingleSectorFitter) fSingleSectorFitter = new AliTPCFitPad(8, "hyp7", ""); | |
570 | ||
10757ee9 | 571 | |
572 | // just for debugging | |
0d3279d4 | 573 | if (!fDebugCalPadRaw) fDebugCalPadRaw = new AliTPCCalPad("DebugCalPadRaw", "All clusters simply added up before correction"); |
574 | if (!fDebugCalPadCorr) fDebugCalPadCorr = new AliTPCCalPad("DebugCalPadCorr", "All clusters simply added up after correction"); | |
10757ee9 | 575 | |
576 | TIterator* iter = list->MakeIterator(); | |
577 | AliTPCcalibTracksGain* cal = 0; | |
578 | ||
579 | while ((cal = (AliTPCcalibTracksGain*)iter->Next())) { | |
580 | if (!cal->InheritsFrom(AliTPCcalibTracksGain::Class())) { | |
581 | Error("Merge","Attempt to add object of class %s to a %s", cal->ClassName(), this->ClassName()); | |
582 | return -1; | |
583 | } | |
0d3279d4 | 584 | |
10757ee9 | 585 | Add(cal); |
586 | } | |
587 | return 0; | |
588 | } | |
589 | ||
2acad464 | 590 | Float_t AliTPCcalibTracksGain::GetGain(AliTPCclusterMI* cl){ |
591 | // | |
592 | // Return local gain at cluster position | |
593 | // | |
594 | Float_t factor = 1; | |
09bed5bf | 595 | if(!fGainMap) return factor; |
596 | ||
2acad464 | 597 | AliTPCCalROC * roc = fGainMap->GetCalROC(cl->GetDetector()); |
598 | Int_t irow = cl->GetRow(); | |
599 | if (roc){ | |
600 | if (irow < 63) { // IROC | |
601 | factor = roc->GetValue(irow, TMath::Nint(cl->GetPad())); | |
602 | } else { // OROC | |
603 | factor = roc->GetValue(irow - 63, TMath::Nint(cl->GetPad())); | |
604 | } | |
605 | } | |
606 | if (factor<0.1) factor=0.1; | |
607 | return factor; | |
608 | } | |
609 | ||
610 | ||
611 | Float_t AliTPCcalibTracksGain::GetMaxNorm(AliTPCclusterMI * cl){ | |
612 | // | |
613 | // Get normalized amplituded if the gain map provided | |
614 | // | |
615 | return cl->GetMax()/GetGain(cl); | |
616 | } | |
617 | ||
618 | ||
619 | Float_t AliTPCcalibTracksGain::GetQNorm(AliTPCclusterMI * cl){ | |
620 | // | |
621 | // Get normalized amplituded if the gain map provided | |
622 | // | |
623 | return cl->GetQ()/GetGain(cl); | |
624 | } | |
625 | ||
626 | ||
627 | ||
10757ee9 | 628 | void AliTPCcalibTracksGain::Add(AliTPCcalibTracksGain* cal) { |
629 | // | |
630 | // Adds another AliTPCcalibTracksGain object to this object. | |
631 | // | |
632 | ||
2e5bcb67 | 633 | fSimpleFitter->Add(cal->fSimpleFitter); |
634 | fSqrtFitter->Add(cal->fSqrtFitter); | |
635 | fLogFitter->Add(cal->fLogFitter); | |
636 | fSingleSectorFitter->Add(cal->fSingleSectorFitter); | |
637 | // | |
638 | // | |
639 | // | |
640 | if (cal->fFitter0M->GetNpoints()>0) fFitter0M->Add(cal->fFitter0M); | |
641 | if (cal->fFitter1M->GetNpoints()>0) fFitter1M->Add(cal->fFitter1M); | |
642 | if (cal->fFitter2M->GetNpoints()>0) fFitter2M->Add(cal->fFitter2M); | |
643 | if (cal->fFitter0T->GetNpoints()>0) fFitter0T->Add(cal->fFitter0T); | |
644 | if (cal->fFitter1T->GetNpoints()>0) fFitter1T->Add(cal->fFitter1T); | |
645 | if (cal->fFitter2T->GetNpoints()>0) fFitter2T->Add(cal->fFitter2T); | |
646 | // | |
647 | if (cal->fDFitter0M->GetNpoints()>0) fDFitter0M->Add(cal->fDFitter0M); | |
648 | if (cal->fDFitter1M->GetNpoints()>0) fDFitter1M->Add(cal->fDFitter1M); | |
649 | if (cal->fDFitter2M->GetNpoints()>0) fDFitter2M->Add(cal->fDFitter2M); | |
650 | if (cal->fDFitter0T->GetNpoints()>0) fDFitter0T->Add(cal->fDFitter0T); | |
651 | if (cal->fDFitter1T->GetNpoints()>0) fDFitter1T->Add(cal->fDFitter1T); | |
652 | if (cal->fDFitter2T->GetNpoints()>0) fDFitter2T->Add(cal->fDFitter2T); | |
684602c8 | 653 | // |
0d3279d4 | 654 | // |
655 | // histograms | |
656 | // | |
657 | for (Int_t i=0; i<73; i++){ | |
658 | TH1F *his,*hism; | |
659 | his = (TH1F*)fArrayQM->At(i); | |
660 | hism = (TH1F*)cal->fArrayQM->At(i); | |
661 | if (his && hism) his->Add(hism); | |
662 | his = (TH1F*)fArrayQT->At(i); | |
663 | hism = (TH1F*)cal->fArrayQT->At(i); | |
664 | if (his && hism) his->Add(hism); | |
665 | } | |
666 | // | |
667 | // | |
668 | for (Int_t i=0; i<37; i++){ | |
669 | TProfile *his,*hism; | |
670 | his = (TProfile*)fProfileArrayQM->At(i); | |
671 | hism = (TProfile*)cal->fProfileArrayQM->At(i); | |
672 | if (his && hism) his->Add(hism); | |
673 | his = (TProfile*)fProfileArrayQT->At(i); | |
674 | hism = (TProfile*)cal->fProfileArrayQT->At(i); | |
675 | if (his && hism) his->Add(hism); | |
676 | } | |
677 | // | |
678 | // | |
679 | for (Int_t i=0; i<37; i++){ | |
680 | TProfile2D *his,*hism; | |
681 | his = (TProfile2D*)fProfileArrayQM2D->At(i); | |
682 | hism = (TProfile2D*)cal->fProfileArrayQM2D->At(i); | |
683 | if (his && hism) his->Add(hism); | |
684 | his = (TProfile2D*)fProfileArrayQT2D->At(i); | |
685 | hism = (TProfile2D*)cal->fProfileArrayQT2D->At(i); | |
686 | if (his && hism) his->Add(hism); | |
687 | } | |
688 | // | |
10757ee9 | 689 | // this will be gone for the a new ROOT version > v5-17-05 |
690 | for (UInt_t iSegment = 0; iSegment < 36; iSegment++) { | |
691 | fNShortClusters[iSegment] += cal->fNShortClusters[iSegment]; | |
692 | fNMediumClusters[iSegment] += cal->fNMediumClusters[iSegment]; | |
693 | fNLongClusters[iSegment] += cal->fNLongClusters[iSegment]; | |
694 | } | |
695 | ||
696 | // just for debugging, remove me | |
697 | fTotalTracks += cal->fTotalTracks; | |
698 | fAcceptedTracks += cal->fAcceptedTracks; | |
699 | fDebugCalPadRaw->Add(cal->fDebugCalPadRaw); | |
700 | fDebugCalPadCorr->Add(cal->fDebugCalPadCorr); | |
701 | ||
10757ee9 | 702 | } |
703 | ||
704 | void AliTPCcalibTracksGain::AddTrack(AliTPCseed* seed) { | |
705 | // | |
706 | // The clusters making up the track (seed) are added to various fit functions. | |
707 | // See AddCluster(...) for more detail. | |
708 | // | |
709 | ||
0d3279d4 | 710 | DumpTrack(seed); |
711 | // | |
712 | // simple histograming part | |
713 | for (Int_t i=0; i<159; i++){ | |
714 | AliTPCclusterMI* cluster = seed->GetClusterPointer(i); | |
715 | if (cluster) AddCluster(cluster); | |
716 | } | |
10757ee9 | 717 | } |
718 | ||
0d3279d4 | 719 | void AliTPCcalibTracksGain::AddCluster(AliTPCclusterMI* cluster){ |
720 | // | |
721 | // Adding cluster information to the simple histograms | |
2acad464 | 722 | // No correction, fittings are applied |
723 | ||
724 | ||
0d3279d4 | 725 | Float_t kThreshold=5; |
726 | if (cluster->GetX()<=0) return; | |
727 | if (cluster->GetQ()<=kThreshold) return; | |
728 | // | |
0d3279d4 | 729 | |
730 | Int_t sector = cluster->GetDetector(); | |
731 | TH1F * his; | |
732 | his = GetQT(sector); | |
2acad464 | 733 | if (his) his->Fill(GetQNorm(cluster)); |
0d3279d4 | 734 | his = GetQT(-1); |
2acad464 | 735 | if (his) his->Fill(GetQNorm(cluster)); |
0d3279d4 | 736 | his = GetQM(sector); |
2acad464 | 737 | if (his) his->Fill(GetMaxNorm(cluster)); |
0d3279d4 | 738 | his = GetQM(-1); |
2acad464 | 739 | if (his) his->Fill(GetMaxNorm(cluster)); |
0d3279d4 | 740 | // |
741 | sector = sector%36; | |
742 | TProfile *prof; | |
743 | prof = GetProfileQT(sector); | |
2acad464 | 744 | if (prof) prof->Fill(cluster->GetX(),GetQNorm(cluster)); |
0d3279d4 | 745 | prof = GetProfileQT(-1); |
2acad464 | 746 | if (prof) prof->Fill(cluster->GetX(),GetQNorm(cluster)); |
0d3279d4 | 747 | prof = GetProfileQM(sector); |
2acad464 | 748 | if (prof) prof->Fill(cluster->GetX(),GetMaxNorm(cluster)); |
0d3279d4 | 749 | prof = GetProfileQM(-1); |
2acad464 | 750 | if (prof) prof->Fill(cluster->GetX(),GetMaxNorm(cluster)); |
0d3279d4 | 751 | // |
752 | Float_t phi = cluster->GetY()/cluster->GetX(); | |
753 | TProfile2D *prof2; | |
754 | prof2 = GetProfileQT2D(sector); | |
2acad464 | 755 | if (prof2) prof2->Fill(cluster->GetX(),phi,GetQNorm(cluster)); |
0d3279d4 | 756 | prof2 = GetProfileQT2D(-1); |
2acad464 | 757 | if (prof2) prof2->Fill(cluster->GetX(),phi,GetQNorm(cluster)); |
0d3279d4 | 758 | prof2 = GetProfileQM2D(sector); |
2acad464 | 759 | if (prof2) prof2->Fill(cluster->GetX(),phi,GetMaxNorm(cluster)); |
0d3279d4 | 760 | prof2 = GetProfileQM2D(-1); |
2acad464 | 761 | if (prof2) prof2->Fill(cluster->GetX(),phi,GetMaxNorm(cluster)); |
0d3279d4 | 762 | |
763 | // | |
764 | } | |
765 | ||
766 | ||
767 | ||
f1c2a4a3 | 768 | void AliTPCcalibTracksGain::AddCluster(AliTPCclusterMI* cluster, Float_t /*momenta*/, Float_t/* mdedx*/, Int_t padType, |
769 | Float_t xcenter, TVectorD& dedxQ, TVectorD& /*dedxM*/, Float_t /*fraction*/, Float_t fraction2, Float_t dedge, | |
770 | TVectorD& /*parY*/, TVectorD& /*parZ*/, TVectorD& meanPos) { | |
10757ee9 | 771 | // |
772 | // Adds cluster to the appropriate fitter for later analysis. | |
773 | // The charge used for the fit is the maximum charge for this specific cluster or the | |
774 | // accumulated charge per cluster, depending on the value of fgkUseTotalCharge. | |
775 | // Depending on the pad size where the cluster is registered, the value will be put in | |
776 | // the appropriate fitter. Furthermore, for each pad size three different types of fitters | |
777 | // are used. The fit functions are the same for all fitters (parabolic functions), but the value | |
778 | // added to each fitter is different. The simple fitter gets the charge plugged in as is, the sqrt fitter | |
779 | // gets the square root of the charge, and the log fitter gets fgkM*(1+q/fgkM), where q is the original charge | |
780 | // and fgkM==25. | |
781 | // | |
2acad464 | 782 | Float_t kedge = 3; |
783 | Float_t kfraction = 0.7; | |
784 | Int_t kinorm = 2; | |
785 | ||
786 | ||
787 | // Where to put selection on threshold? | |
788 | // Defined by the Q/dEdxT variable - see debug streamer: | |
789 | // | |
790 | // Debug stream variables: (Where tu cut ?) | |
791 | // chain0->Draw("Cl.fQ/dedxQ.fElements[1]>>his(100,0,3)","fraction2<0.6&&dedge>3","",1000000); | |
792 | // mean 1 sigma 0.25 | |
793 | // chain0->Draw("Cl.fMax/dedxM.fElements[1]>>his(100,0,3)","fraction2<0.6&&dedge>3","",1000000) | |
794 | // mean 1 sigma 0.25 | |
795 | // chain0->Draw("Cl.fQ/dedxQ.fElements[2]>>his(100,0,3)","fraction2<0.7&&dedge>3","",1000000) | |
796 | // mean 1 sigma 0.29 | |
797 | // chain0->Draw("Cl.fMax/dedxM.fElements[2]>>his(100,0,3)","fraction2<0.7&&dedge>3","",1000000) | |
798 | // mean 1 sigma 0.27 | |
799 | // chain0->Draw("Cl.fQ/dedxQ.fElements[3]>>his(100,0,3)","fraction2<0.8&&dedge>3","",1000000) | |
800 | // mean 1 sigma 0.32 | |
801 | // | |
802 | // chain0->Draw("Cl.fQ/dedxQ.fElements[4]>>his(100,0,3)","fraction2<0.9&&dedge>3","",1000000) | |
803 | // mean 1 sigma 0.4 | |
804 | ||
805 | // Fraction choosen 0.7 | |
10757ee9 | 806 | |
807 | if (!cluster) { | |
808 | Error("AddCluster", "Cluster not valid."); | |
809 | return; | |
810 | } | |
811 | ||
2acad464 | 812 | if (dedge < kedge) return; |
813 | if (fraction2 > kfraction) return; | |
10757ee9 | 814 | |
815 | //Int_t padType = GetPadType(cluster->GetX()); | |
816 | Double_t xx[7]; | |
817 | //Double_t centerPad[2] = {0}; | |
818 | //AliTPCFitPad::GetPadRegionCenterLocal(padType, centerPad); | |
819 | //xx[0] = cluster->GetX() - centerPad[0]; | |
820 | //xx[1] = cluster->GetY() - centerPad[1]; | |
821 | xx[0] = cluster->GetX() - xcenter; | |
822 | xx[1] = cluster->GetY(); | |
823 | xx[2] = xx[0] * xx[0]; | |
824 | xx[3] = xx[1] * xx[1]; | |
825 | xx[4] = xx[0] * xx[1]; | |
826 | xx[5] = TMath::Abs(cluster->GetZ()) - TMath::Abs(meanPos[4]); | |
827 | xx[6] = xx[5] * xx[5]; | |
0d3279d4 | 828 | // |
829 | // Update profile histograms | |
830 | // | |
10757ee9 | 831 | |
0d3279d4 | 832 | // |
833 | // Update fitters | |
834 | // | |
10757ee9 | 835 | Int_t segment = cluster->GetDetector() % 36; |
2acad464 | 836 | Double_t q = fgkUseTotalCharge ? ((Double_t)(GetQNorm(cluster))) : ((Double_t)(GetMaxNorm(cluster))); // note: no normalization to pad size! |
10757ee9 | 837 | |
838 | // just for debugging | |
839 | Int_t row = 0; | |
840 | Int_t pad = 0; | |
841 | GetRowPad(cluster->GetX(), cluster->GetY(), row, pad); | |
842 | fDebugCalPadRaw->GetCalROC(cluster->GetDetector())->SetValue(row, pad, q + fDebugCalPadRaw->GetCalROC(cluster->GetDetector())->GetValue(row, pad)); | |
843 | ||
844 | // correct charge by normalising to mean charge per track | |
2acad464 | 845 | q /= dedxQ[kinorm]; |
10757ee9 | 846 | |
847 | // just for debugging | |
848 | fDebugCalPadCorr->GetCalROC(cluster->GetDetector())->SetValue(row, pad, q + fDebugCalPadCorr->GetCalROC(cluster->GetDetector())->GetValue(row, pad)); | |
849 | ||
850 | Double_t sqrtQ = TMath::Sqrt(q); | |
851 | Double_t logQ = fgkM * TMath::Log(1 + q / fgkM); | |
684602c8 | 852 | TLinearFitter * fitter =0; |
853 | // | |
854 | fitter = fSimpleFitter->GetFitter(segment, padType); | |
855 | fitter->AddPoint(xx, q); | |
856 | // | |
857 | fitter = fSqrtFitter->GetFitter(segment, padType); | |
858 | fitter->AddPoint(xx, sqrtQ); | |
859 | // | |
860 | fitter = fLogFitter->GetFitter(segment, padType); | |
861 | fitter->AddPoint(xx, logQ); | |
862 | // | |
863 | fitter=fSingleSectorFitter->GetFitter(0, padType); | |
864 | fitter->AddPoint(xx, q); | |
10757ee9 | 865 | |
866 | // this will be gone for the a new ROOT version > v5-17-05 | |
867 | if (padType == kShortPads) | |
868 | fNShortClusters[segment]++; | |
0d3279d4 | 869 | if (padType == kMediumPads) |
10757ee9 | 870 | fNMediumClusters[segment]++; |
0d3279d4 | 871 | if (padType == kLongPads) |
10757ee9 | 872 | fNLongClusters[segment]++; |
873 | } | |
874 | ||
875 | void AliTPCcalibTracksGain::Evaluate(Bool_t robust, Double_t frac) { | |
876 | // | |
877 | // Evaluates all fitters contained in this object. | |
878 | // If the robust option is set to kTRUE a robust fit is performed with frac as | |
879 | // the minimal fraction of good points (see TLinearFitter::EvalRobust for details). | |
880 | // Beware: Robust fitting is much slower! | |
881 | // | |
882 | ||
883 | fSimpleFitter->Evaluate(robust, frac); | |
884 | fSqrtFitter->Evaluate(robust, frac); | |
885 | fLogFitter->Evaluate(robust, frac); | |
886 | fSingleSectorFitter->Evaluate(robust, frac); | |
0d3279d4 | 887 | fFitter0M->Eval(); |
888 | fFitter1M->Eval(); | |
889 | fFitter2M->Eval(); | |
890 | fFitter0T->Eval(); | |
891 | fFitter1T->Eval(); | |
892 | fFitter2T->Eval(); | |
684602c8 | 893 | // |
894 | fDFitter0M->Eval(); | |
895 | fDFitter1M->Eval(); | |
896 | fDFitter2M->Eval(); | |
897 | fDFitter0T->Eval(); | |
898 | fDFitter1T->Eval(); | |
899 | fDFitter2T->Eval(); | |
10757ee9 | 900 | } |
901 | ||
902 | AliTPCCalPad* AliTPCcalibTracksGain::CreateFitCalPad(UInt_t fitType, Bool_t undoTransformation, Bool_t normalizeToPadSize) { | |
903 | // | |
904 | // Creates the calibration object AliTPCcalPad using fitted parameterization | |
905 | // | |
906 | TObjArray tpc(72); | |
907 | for (UInt_t iSector = 0; iSector < 72; iSector++) | |
908 | tpc.Add(CreateFitCalROC(iSector, fitType, undoTransformation, normalizeToPadSize)); | |
909 | return new AliTPCCalPad(&tpc); | |
910 | } | |
911 | ||
912 | AliTPCCalROC* AliTPCcalibTracksGain::CreateFitCalROC(UInt_t sector, UInt_t fitType, Bool_t undoTransformation, Bool_t normalizeToPadSize) { | |
913 | // | |
914 | // Create the AliTPCCalROC with the values per pad | |
915 | // sector - sector of interest | |
916 | // fitType - | |
917 | // | |
918 | ||
919 | TVectorD par(8); | |
920 | if (sector < 36) { | |
921 | GetParameters(sector % 36, 0, fitType, par); | |
922 | return CreateFitCalROC(sector, 0, par, fitType, undoTransformation, normalizeToPadSize); | |
923 | } | |
924 | else { | |
925 | GetParameters(sector % 36, 1, fitType, par); | |
926 | AliTPCCalROC* roc1 = CreateFitCalROC(sector, 1, par, fitType, undoTransformation, normalizeToPadSize); | |
927 | GetParameters(sector % 36, 2, fitType, par); | |
928 | AliTPCCalROC* roc2 = CreateFitCalROC(sector, 2, par, fitType, undoTransformation, normalizeToPadSize); | |
929 | AliTPCCalROC* roc3 = CreateCombinedCalROC(roc1, roc2); | |
930 | delete roc1; | |
931 | delete roc2; | |
932 | return roc3; | |
933 | } | |
934 | } | |
935 | ||
936 | AliTPCCalROC* AliTPCcalibTracksGain::CreateFitCalROC(UInt_t sector, UInt_t padType, TVectorD &fitParam, UInt_t fitType, Bool_t undoTransformation, Bool_t normalizeToPadSize) { | |
937 | // | |
938 | // This function is essentially a copy of AliTPCCalROC::CreateGlobalFitCalROC(...), with the | |
939 | // modifications, that the center of the region of same pad size is used as the origin | |
940 | // of the fit function instead of the center of the ROC. | |
941 | // The possibility of a linear fit is removed as well because it is not needed. | |
942 | // Only values for pads with the given pad size are calculated, the rest is 0. | |
943 | // Set undoTransformation for undoing the transformation that was applied to the | |
944 | // charge values before they were put into the fitter (thus allowing comparison to the original | |
945 | // charge values). For fitType use 0 for the simple fitter, 1 for the sqrt fitter, 2 for the log fitter. | |
946 | // If normalizeToPadSize is true, the values are normalized to the pad size. | |
947 | // Please be aware, that you even need to specify the fitType if you want to normalize to the pad size without | |
948 | // undoing the transformation (because normalizing involves undoing the trafo first, then normalizing, then | |
949 | // applying the trafo again). | |
950 | // Please note: The normalization to the pad size is a simple linear scaling with the pad length, which | |
951 | // actually doesn't describe reality! | |
952 | // | |
953 | ||
954 | Float_t dlx, dly; | |
955 | Double_t centerPad[2] = {0}; | |
956 | Float_t localXY[3] = {0}; | |
957 | AliTPCROC* tpcROC = AliTPCROC::Instance(); | |
958 | if ((padType == 0 && sector >= tpcROC->GetNInnerSector()) || (padType > 0 && sector < tpcROC->GetNInnerSector()) || sector >= tpcROC->GetNSector()) | |
959 | return 0; | |
960 | AliTPCCalROC* lROCfitted = new AliTPCCalROC(sector); | |
961 | //tpcROC->GetPositionLocal(sector, lROCfitted->GetNrows()/2, lROCfitted->GetNPads(lROCfitted->GetNrows()/2)/2, centerPad); // use this instead of the switch statement if you want to calculate the center of the ROC and not the center of the regions with the same pad size | |
962 | UInt_t startRow = 0; | |
963 | UInt_t endRow = 0; | |
964 | switch (padType) { | |
965 | case kShortPads: | |
966 | startRow = 0; | |
967 | endRow = lROCfitted->GetNrows(); | |
968 | break; | |
969 | case kMediumPads: | |
970 | startRow = 0; | |
971 | endRow = 64; | |
972 | break; | |
973 | case kLongPads: | |
974 | startRow = 64; | |
975 | endRow = lROCfitted->GetNrows(); | |
976 | break; | |
977 | } | |
978 | ||
979 | AliTPCFitPad::GetPadRegionCenterLocal(padType, centerPad); | |
980 | Double_t value = 0; | |
981 | for (UInt_t irow = startRow; irow < endRow; irow++) { | |
982 | for (UInt_t ipad = 0; ipad < lROCfitted->GetNPads(irow); ipad++) { | |
983 | tpcROC->GetPositionLocal(sector, irow, ipad, localXY); // calculate position localXY by pad and row number | |
984 | dlx = localXY[0] - centerPad[0]; | |
985 | dly = localXY[1] - centerPad[1]; | |
986 | value = fitParam[0] + fitParam[1]*dlx + fitParam[2]*dly + fitParam[3]*dlx*dlx + fitParam[4]*dly*dly + fitParam[5]*dlx*dly; | |
987 | ||
988 | // Let q' = value be the transformed value without any pad size corrections, | |
989 | // let T be the transformation and let l be the pad size | |
990 | // 1) don't undo transformation, don't normalize: return q' | |
991 | // 2) undo transformation, don't normalize: return T^{-1} q' | |
992 | // 3) undo transformation, normalize: return (T^{-1} q') / l | |
993 | // 4) don't undo transformation, normalize: return T((T^{-1} q') / l) | |
994 | if (!undoTransformation && !normalizeToPadSize) {/* value remains unchanged */} // (1) | |
995 | else { // (2), (3), (4) | |
996 | //calculate T^{-1} | |
997 | switch (fitType) { | |
998 | case 0: /* value remains unchanged */ break; | |
999 | case 1: value = value * value; break; | |
1000 | case 2: value = (TMath::Exp(value / fgkM) - 1) * fgkM; break; | |
1001 | default: Error("CreateFitCalROC", "Wrong fit type."); break; | |
1002 | } | |
1003 | if (normalizeToPadSize) value /= GetPadLength(localXY[0]); // (3) | |
1004 | } | |
1005 | if (!undoTransformation && normalizeToPadSize) { // (4) | |
1006 | // calculate T | |
1007 | switch (fitType) { | |
1008 | case 0: /* value remains unchanged */ break; | |
1009 | case 1: value = TMath::Sqrt(value); break; | |
1010 | case 2: value = fgkM * TMath::Log(1 + value / fgkM); break; | |
1011 | default: Error("CreateFitCalROC", "Wrong fit type."); break; | |
1012 | } | |
1013 | } | |
1014 | lROCfitted->SetValue(irow, ipad, value); | |
1015 | } | |
1016 | } | |
1017 | return lROCfitted; | |
1018 | } | |
1019 | ||
1020 | AliTPCCalROC* AliTPCcalibTracksGain::CreateCombinedCalROC(const AliTPCCalROC* roc1, const AliTPCCalROC* roc2) { | |
1021 | // | |
1022 | // Combines the medium pad size values of roc1 with the long pad size values of roc2 into a new | |
1023 | // AliTPCCalROC. Returns a null pointer if any one of the ROCs is an IROC; issues a warning message | |
1024 | // if the sectors of roc1 and roc2 don't match, but still continue and use the sector of roc1 as the | |
1025 | // sector of the new ROC. | |
1026 | // | |
1027 | ||
1028 | if (!roc1 || !roc2) return 0; | |
1029 | if ((Int_t)(roc1->GetSector()) < fgTPCparam->GetNInnerSector()) return 0; | |
1030 | if ((Int_t)(roc2->GetSector()) < fgTPCparam->GetNInnerSector()) return 0; | |
1031 | if (roc1->GetSector() != roc2->GetSector()) Warning("CreateCombinedCalROC", "Sector number mismatch."); | |
1032 | AliTPCCalROC* roc = new AliTPCCalROC(roc1->GetSector()); | |
1033 | ||
1034 | for (UInt_t iRow = 0; iRow < 64; iRow++) { | |
1035 | for (UInt_t iPad = 0; iPad < roc->GetNPads(iRow); iPad++) | |
1036 | roc->SetValue(iRow, iPad, roc1->GetValue(iRow, iPad)); | |
1037 | } | |
1038 | for (UInt_t iRow = 64; iRow < roc->GetNrows(); iRow++) { | |
1039 | for (UInt_t iPad = 0; iPad < roc->GetNPads(iRow); iPad++) | |
1040 | roc->SetValue(iRow, iPad, roc2->GetValue(iRow, iPad)); | |
1041 | } | |
1042 | return roc; | |
1043 | } | |
1044 | ||
684602c8 | 1045 | Bool_t AliTPCcalibTracksGain::GetParameters(UInt_t segment, UInt_t padType, UInt_t fitType, TVectorD &fitParam) { |
10757ee9 | 1046 | // |
1047 | // Puts the fit parameters for the specified segment (IROC & OROC), padType and fitType | |
1048 | // into the fitParam TVectorD (which should contain 8 elements). | |
1049 | // padType is one of kShortPads, kMediumPads, kLongPads. fitType is one of kSimpleFitter, kSqrtFitter, kLogFitter. | |
1050 | // Note: The fitter has to be evaluated first! | |
1051 | // | |
684602c8 | 1052 | TLinearFitter * fitter = GetFitter(segment, padType, fitType); |
1053 | if (fitter){ | |
1054 | fitter->Eval(); | |
1055 | fitter->GetParameters(fitParam); | |
1056 | return kTRUE; | |
1057 | }else{ | |
1058 | Error("AliTPCcalibTracksGain::GetParameters", | |
1059 | Form("Fitter%d_%d_%d not availble", segment, padType, fitType)); | |
1060 | return kFALSE; | |
1061 | } | |
1062 | return kFALSE; | |
10757ee9 | 1063 | } |
1064 | ||
1065 | void AliTPCcalibTracksGain::GetErrors(UInt_t segment, UInt_t padType, UInt_t fitType, TVectorD &fitError) { | |
1066 | // | |
1067 | // Puts the fit parameter errors for the specified segment (IROC & OROC), padType and fitType | |
1068 | // into the fitError TVectorD (which should contain 8 elements). | |
1069 | // padType is one of kShortPads, kMediumPads, kLongPads. fitType is one of kSimpleFitter, kSqrtFitter, kLogFitter. | |
1070 | // Note: The fitter has to be evaluated first! | |
1071 | // | |
1072 | ||
1073 | GetFitter(segment, padType, fitType)->GetErrors(fitError); | |
1074 | fitError *= TMath::Sqrt(GetRedChi2(segment, padType, fitType)); | |
1075 | } | |
1076 | ||
1077 | Double_t AliTPCcalibTracksGain::GetRedChi2(UInt_t segment, UInt_t padType, UInt_t fitType) { | |
1078 | // | |
1079 | // Returns the reduced chi^2 value for the specified segment, padType and fitType. | |
1080 | // padType is one of kShortPads, kMediumPads, kLongPads. fitType is one of kSimpleFitter, kSqrtFitter, kLogFitter. | |
1081 | // Note: The fitter has to be evaluated first! | |
1082 | // | |
1083 | ||
1084 | // this will be gone for the a new ROOT version > v5-17-05 | |
1085 | Int_t lNClusters = 0; | |
1086 | switch (padType) { | |
1087 | case kShortPads: | |
1088 | lNClusters = fNShortClusters[segment]; | |
1089 | break; | |
1090 | case kMediumPads: | |
1091 | lNClusters = fNMediumClusters[segment]; | |
1092 | break; | |
1093 | case kLongPads: | |
1094 | lNClusters = fNLongClusters[segment]; | |
1095 | break; | |
1096 | } | |
1097 | return GetFitter(segment, padType, fitType)->GetChisquare()/(lNClusters - 8); | |
1098 | } | |
1099 | ||
1100 | void AliTPCcalibTracksGain::GetCovarianceMatrix(UInt_t segment, UInt_t padType, UInt_t fitType, TMatrixD& covMatrix) { | |
1101 | // | |
1102 | // Returns the covariance matrix for the specified segment, padType, fitType. | |
1103 | // padType is one of kShortPads, kMediumPads, kLongPads. fitType is one of kSimpleFitter, kSqrtFitter, kLogFitter. | |
1104 | // | |
1105 | ||
1106 | GetFitter(segment, padType, fitType)->GetCovarianceMatrix(covMatrix); | |
1107 | } | |
1108 | ||
1109 | TLinearFitter* AliTPCcalibTracksGain::GetFitter(UInt_t segment, UInt_t padType, UInt_t fitType) { | |
1110 | // | |
1111 | // Returns the TLinearFitter object for the specified segment, padType, fitType. | |
1112 | // padType is one of kShortPads, kMediumPads, kLongPads. fitType is one of kSimpleFitter, kSqrtFitter, kLogFitter. | |
1113 | // | |
1114 | ||
1115 | switch (fitType) { | |
1116 | case kSimpleFitter: | |
1117 | return fSimpleFitter->GetFitter(segment, padType); | |
1118 | case kSqrtFitter: | |
1119 | return fSqrtFitter->GetFitter(segment, padType); | |
1120 | case kLogFitter: | |
1121 | return fLogFitter->GetFitter(segment, padType); | |
1122 | case 3: | |
1123 | return fSingleSectorFitter->GetFitter(0, padType); | |
1124 | } | |
1125 | return 0; | |
1126 | } | |
1127 | ||
1128 | Double_t AliTPCcalibTracksGain::GetPadLength(Double_t lx) { | |
1129 | // | |
1130 | // The function returns 0.75 for an IROC, 1. for an OROC at medium pad size position, | |
1131 | // 1.5 for an OROC at long pad size position, -1 if out of bounds. | |
1132 | // | |
1133 | ||
1134 | Double_t irocLow = fgTPCparam->GetPadRowRadiiLow(0) - fgTPCparam->GetInnerPadPitchLength()/2; | |
1135 | Double_t irocUp = fgTPCparam->GetPadRowRadiiLow(fgTPCparam->GetNRowLow()-1) + fgTPCparam->GetInnerPadPitchLength()/2; | |
1136 | Double_t orocLow1 = fgTPCparam->GetPadRowRadiiUp(0) - fgTPCparam->GetOuter1PadPitchLength()/2; | |
1137 | Double_t orocUp1 = fgTPCparam->GetPadRowRadiiUp(fgTPCparam->GetNRowUp1()-1) + fgTPCparam->GetOuter1PadPitchLength()/2; | |
1138 | Double_t orocLow2 = fgTPCparam->GetPadRowRadiiUp(fgTPCparam->GetNRowUp1()) - fgTPCparam->GetOuter2PadPitchLength()/2; | |
1139 | Double_t orocUp2 = fgTPCparam->GetPadRowRadiiUp(fgTPCparam->GetNRowUp()-1) + fgTPCparam->GetOuter2PadPitchLength()/2; | |
1140 | ||
1141 | // if IROC | |
1142 | if (lx >= irocLow && lx <= irocUp) return 0.75; | |
1143 | // if OROC medium pads | |
1144 | if (lx >= orocLow1 && lx <= orocUp1) return 1.; | |
1145 | // if OROC long pads | |
1146 | if (lx >= orocLow2 && lx <= orocUp2) return 1.5; | |
1147 | // if out of bounds | |
1148 | return -1; | |
1149 | } | |
1150 | ||
1151 | Int_t AliTPCcalibTracksGain::GetPadType(Double_t lx) { | |
1152 | // | |
1153 | // The function returns 0 for an IROC, 1 for an OROC at medium pad size position, | |
1154 | // 2 for an OROC at long pad size position, -1 if out of bounds. | |
1155 | // | |
1156 | ||
1157 | if (GetPadLength(lx) == 0.75) return 0; | |
1158 | else if (GetPadLength(lx) == 1.) return 1; | |
1159 | else if (GetPadLength(lx) == 1.5) return 2; | |
1160 | return -1; | |
1161 | } | |
1162 | ||
1163 | // ONLY FOR DEBUGGING PURPOSES - REMOVE ME WHEN NOT NEEDED ANYMORE | |
1164 | Bool_t AliTPCcalibTracksGain::GetRowPad(Double_t lx, Double_t ly, Int_t& row, Int_t& pad) { | |
1165 | // | |
1166 | // Calculate the row and pad number when the local coordinates are given. | |
1167 | // Returns kFALSE if the position is out of range, otherwise return kTRUE. | |
1168 | // WARNING: This function is preliminary and probably isn't very accurate!! | |
1169 | // | |
1170 | ||
1171 | Double_t irocLow = fgTPCparam->GetPadRowRadiiLow(0) - fgTPCparam->GetInnerPadPitchLength()/2; | |
1172 | //Double_t irocUp = fgTPCparam->GetPadRowRadiiLow(fgTPCparam->GetNRowLow()-1) + fgTPCparam->GetInnerPadPitchLength()/2; | |
1173 | Double_t orocLow1 = fgTPCparam->GetPadRowRadiiUp(0) - fgTPCparam->GetOuter1PadPitchLength()/2; | |
1174 | //Double_t orocUp1 = fgTPCparam->GetPadRowRadiiUp(fgTPCparam->GetNRowUp1()-1) + fgTPCparam->GetOuter1PadPitchLength()/2; | |
1175 | Double_t orocLow2 = fgTPCparam->GetPadRowRadiiUp(fgTPCparam->GetNRowUp1()) - fgTPCparam->GetOuter2PadPitchLength()/2; | |
1176 | //Double_t orocUp2 = fgTPCparam->GetPadRowRadiiUp(fgTPCparam->GetNRowUp()-1) + fgTPCparam->GetOuter2PadPitchLength()/2; | |
1177 | ||
1178 | if (GetPadType(lx) == 0) { | |
1179 | row = (Int_t)((lx - irocLow) / fgTPCparam->GetInnerPadPitchLength()); | |
1180 | pad = (Int_t)((ly + fgTPCparam->GetYInner(row)) / fgTPCparam->GetInnerPadPitchWidth()); | |
1181 | } else if (GetPadType(lx) == 1) { | |
1182 | row = (Int_t)((lx - orocLow1) / fgTPCparam->GetOuter1PadPitchLength()); | |
1183 | pad = (Int_t)((ly + fgTPCparam->GetYOuter(row)) / fgTPCparam->GetOuterPadPitchWidth()); | |
1184 | } else if (GetPadType(lx) == 2) { | |
1185 | row = fgTPCparam->GetNRowUp1() + (Int_t)((lx - orocLow2) / fgTPCparam->GetOuter2PadPitchLength()); | |
1186 | pad = (Int_t)((ly + fgTPCparam->GetYOuter(row)) / fgTPCparam->GetOuterPadPitchWidth()); | |
1187 | } | |
1188 | else return kFALSE; | |
1189 | return kTRUE; | |
1190 | } | |
1191 | ||
1192 | void AliTPCcalibTracksGain::DumpTrack(AliTPCseed* track) { | |
1193 | // | |
1194 | // Dump track information to the debug stream | |
1195 | // | |
1196 | ||
10757ee9 | 1197 | Int_t rows[200]; |
0d3279d4 | 1198 | Int_t sector[3]; |
1199 | Int_t npoints[3]; | |
1200 | TVectorD dedxM[3]; | |
1201 | TVectorD dedxQ[3]; | |
1202 | TVectorD parY[3]; | |
1203 | TVectorD parZ[3]; | |
1204 | TVectorD meanPos[3]; | |
1205 | // | |
1206 | Int_t count=0; | |
10757ee9 | 1207 | for (Int_t ipad = 0; ipad < 3; ipad++) { |
0d3279d4 | 1208 | dedxM[ipad].ResizeTo(5); |
1209 | dedxQ[ipad].ResizeTo(5); | |
1210 | parY[ipad].ResizeTo(3); | |
1211 | parZ[ipad].ResizeTo(3); | |
1212 | meanPos[ipad].ResizeTo(6); | |
1213 | Bool_t isOK = GetDedx(track, ipad, rows, sector[ipad], npoints[ipad], dedxM[ipad], dedxQ[ipad], parY[ipad], parZ[ipad], meanPos[ipad]); | |
1214 | if (isOK) | |
1215 | AddTracklet(sector[ipad],ipad, dedxQ[ipad], dedxM[ipad], parY[ipad], parZ[ipad], meanPos[ipad] ); | |
1216 | if (isOK) count++; | |
10757ee9 | 1217 | } |
ae28e92e | 1218 | |
1219 | TTreeSRedirector * cstream = GetDebugStreamer(); | |
1220 | if (cstream){ | |
1221 | (*cstream) << "Track" << | |
0d3279d4 | 1222 | "Track.=" << track << // track information |
0d3279d4 | 1223 | "\n"; |
ae28e92e | 1224 | // |
1225 | // | |
1226 | // | |
1227 | if ( GetStreamLevel()>1 && count>1){ | |
1228 | (*cstream) << "TrackG" << | |
1229 | "Track.=" << track << // track information | |
1230 | "ncount="<<count<< | |
1231 | // info for pad type 0 | |
1232 | "sector0="<<sector[0]<< | |
1233 | "npoints0="<<npoints[0]<< | |
1234 | "dedxM0.="<<&dedxM[0]<< | |
1235 | "dedxQ0.="<<&dedxQ[0]<< | |
1236 | "parY0.="<<&parY[0]<< | |
1237 | "parZ0.="<<&parZ[0]<< | |
1238 | "meanPos0.="<<&meanPos[0]<< | |
1239 | // | |
1240 | // info for pad type 1 | |
1241 | "sector1="<<sector[1]<< | |
1242 | "npoints1="<<npoints[1]<< | |
1243 | "dedxM1.="<<&dedxM[1]<< | |
1244 | "dedxQ1.="<<&dedxQ[1]<< | |
1245 | "parY1.="<<&parY[1]<< | |
1246 | "parZ1.="<<&parZ[1]<< | |
1247 | "meanPos1.="<<&meanPos[1]<< | |
1248 | // | |
1249 | // info for pad type 2 | |
1250 | "sector2="<<sector[2]<< | |
1251 | "npoints2="<<npoints[2]<< | |
1252 | "dedxM2.="<<&dedxM[2]<< | |
1253 | "dedxQ2.="<<&dedxQ[2]<< | |
1254 | "parY2.="<<&parY[2]<< | |
1255 | "parZ2.="<<&parZ[2]<< | |
1256 | "meanPos2.="<<&meanPos[2]<< | |
1257 | // | |
1258 | "\n"; | |
1259 | } | |
0d3279d4 | 1260 | } |
0d3279d4 | 1261 | } |
1262 | ||
f1c2a4a3 | 1263 | Bool_t AliTPCcalibTracksGain::GetDedx(AliTPCseed* track, Int_t padType, Int_t* /*rows*/, |
0d3279d4 | 1264 | Int_t §or, Int_t& npoints, |
1265 | TVectorD &dedxM, TVectorD &dedxQ, | |
1266 | TVectorD &parY, TVectorD &parZ, TVectorD&meanPos) | |
1267 | { | |
1268 | // | |
1269 | // GetDedx for given sector for given track | |
1270 | // padType - type of pads | |
1271 | // | |
1272 | ||
1273 | static TLinearFitter fitY(2, "pol1"); | |
1274 | static TLinearFitter fitZ(2, "pol1"); | |
684602c8 | 1275 | fitY.StoreData(kFALSE); |
1276 | fitZ.StoreData(kFALSE); | |
0d3279d4 | 1277 | // |
1278 | // | |
10757ee9 | 1279 | Int_t firstRow = 0, lastRow = 0; |
1280 | Int_t minRow = 100; | |
1281 | Float_t xcenter = 0; | |
1282 | const Float_t ktany = TMath::Tan(TMath::DegToRad() * 10); | |
1283 | const Float_t kedgey = 4.; | |
1284 | if (padType == 0) { | |
1285 | firstRow = 0; | |
1286 | lastRow = fgTPCparam->GetNRowLow(); | |
1287 | xcenter = 108.47; | |
1288 | } | |
1289 | if (padType == 1) { | |
1290 | firstRow = fgTPCparam->GetNRowLow(); | |
1291 | lastRow = fgTPCparam->GetNRowLow() + fgTPCparam->GetNRowUp1(); | |
1292 | xcenter = 166.60; | |
1293 | } | |
1294 | if (padType == 2) { | |
1295 | firstRow = fgTPCparam->GetNRowLow() + fgTPCparam->GetNRowUp1(); | |
1296 | lastRow = fgTPCparam->GetNRowLow() + fgTPCparam->GetNRowUp(); | |
1297 | xcenter = 222.6; | |
1298 | } | |
1299 | minRow = (lastRow - firstRow) / 2; | |
1300 | // | |
1301 | // | |
1302 | Int_t nclusters = 0; | |
1303 | Int_t nclustersNE = 0; // number of not edge clusters | |
1304 | Int_t lastSector = -1; | |
1305 | Float_t amplitudeQ[100]; | |
1306 | Float_t amplitudeM[100]; | |
1307 | Int_t rowIn[100]; | |
1308 | Int_t index[100]; | |
1309 | // | |
0d3279d4 | 1310 | // |
10757ee9 | 1311 | fitY.ClearPoints(); |
1312 | fitZ.ClearPoints(); | |
10757ee9 | 1313 | |
1314 | for (Int_t iCluster = firstRow; iCluster < lastRow; iCluster++) { | |
1315 | AliTPCclusterMI* cluster = track->GetClusterPointer(iCluster); | |
1316 | if (cluster) { | |
1317 | Int_t detector = cluster->GetDetector() ; | |
1318 | if (lastSector == -1) lastSector = detector; | |
1319 | if (lastSector != detector) continue; | |
2acad464 | 1320 | amplitudeQ[nclusters] = GetQNorm(cluster); |
1321 | amplitudeM[nclusters] = GetMaxNorm(cluster); | |
10757ee9 | 1322 | rowIn[nclusters] = iCluster; |
1323 | nclusters++; | |
1324 | Double_t dx = cluster->GetX() - xcenter; | |
1325 | Double_t y = cluster->GetY(); | |
1326 | Double_t z = cluster->GetZ(); | |
1327 | fitY.AddPoint(&dx, y); | |
1328 | fitZ.AddPoint(&dx, z); | |
1329 | meanPos[0] += dx; | |
1330 | meanPos[1] += dx; | |
1331 | meanPos[2] += y; | |
1332 | meanPos[3] += y*y; | |
1333 | meanPos[4] += z; | |
1334 | meanPos[5] += z*z; | |
1335 | if (TMath::Abs(cluster->GetY()) < cluster->GetX()*ktany - kedgey) nclustersNE++; | |
1336 | } | |
1337 | } | |
1338 | ||
1339 | if (nclusters < minRow / 2) return kFALSE; | |
1340 | if (nclustersNE < minRow / 2) return kFALSE; | |
1341 | for (Int_t i = 0; i < 6; i++) meanPos[i] /= Double_t(nclusters); | |
1342 | fitY.Eval(); | |
1343 | fitZ.Eval(); | |
1344 | fitY.GetParameters(parY); | |
1345 | fitZ.GetParameters(parZ); | |
1346 | // | |
1347 | // calculate truncated mean | |
1348 | // | |
1349 | TMath::Sort(nclusters, amplitudeQ, index, kFALSE); | |
0d3279d4 | 1350 | // |
1351 | // | |
1352 | // | |
10757ee9 | 1353 | Float_t ndedx[5]; |
1354 | for (Int_t i = 0; i < 5; i++) { | |
1355 | dedxQ[i] = 0; | |
1356 | dedxM[i] = 0; | |
1357 | ndedx[i] = 0; | |
1358 | } | |
1359 | // | |
1360 | // dedx calculation | |
1361 | // | |
1362 | Int_t inonEdge = 0; | |
1363 | for (Int_t i = 0; i < nclusters; i++) { | |
1364 | Int_t rowSorted = rowIn[index[i]]; | |
1365 | AliTPCclusterMI* cluster = track->GetClusterPointer(rowSorted); | |
1366 | ||
1367 | if (TMath::Abs(cluster->GetY()) > cluster->GetX()*ktany - kedgey) continue; //don't take edge clusters | |
1368 | inonEdge++; | |
1369 | if (inonEdge < nclustersNE * 0.5) { | |
1370 | ndedx[0]++; | |
1371 | dedxQ[0] += amplitudeQ[index[i]]; | |
1372 | dedxM[0] += amplitudeM[index[i]]; | |
1373 | } | |
1374 | if (inonEdge < nclustersNE * 0.6) { | |
1375 | ndedx[1]++; | |
1376 | dedxQ[1] += amplitudeQ[index[i]]; | |
1377 | dedxM[1] += amplitudeM[index[i]]; | |
1378 | } | |
1379 | if (inonEdge < nclustersNE * 0.7) { | |
1380 | ndedx[2]++; | |
1381 | dedxQ[2] += amplitudeQ[index[i]]; | |
1382 | dedxM[2] += amplitudeM[index[i]]; | |
1383 | } | |
1384 | if (inonEdge < nclustersNE * 0.8) { | |
1385 | ndedx[3]++; | |
1386 | dedxQ[3] += amplitudeQ[index[i]]; | |
1387 | dedxM[3] += amplitudeM[index[i]]; | |
1388 | } | |
1389 | if (inonEdge < nclustersNE * 0.9) { | |
1390 | ndedx[4]++; | |
1391 | dedxQ[4] += amplitudeQ[index[i]]; | |
1392 | dedxM[4] += amplitudeM[index[i]]; | |
1393 | } | |
1394 | } | |
1395 | for (Int_t i = 0; i < 5; i++) { | |
1396 | dedxQ[i] /= ndedx[i]; | |
1397 | dedxM[i] /= ndedx[i]; | |
1398 | } | |
ae28e92e | 1399 | TTreeSRedirector * cstream = GetDebugStreamer(); |
10757ee9 | 1400 | inonEdge = 0; |
0d3279d4 | 1401 | Float_t momenta = track->GetP(); |
1402 | Float_t mdedx = track->GetdEdx(); | |
10757ee9 | 1403 | for (Int_t i = 0; i < nclusters; i++) { |
1404 | Int_t rowSorted = rowIn[index[i]]; | |
1405 | AliTPCclusterMI* cluster = track->GetClusterPointer(rowSorted); | |
1406 | if (!cluster) { | |
1407 | printf("Problem\n"); | |
1408 | continue; | |
1409 | } | |
1410 | if (TMath::Abs(cluster->GetY()) < cluster->GetX()*ktany - kedgey) inonEdge++; | |
1411 | Float_t dedge = cluster->GetX()*ktany - TMath::Abs(cluster->GetY()); | |
1412 | Float_t fraction = Float_t(i) / Float_t(nclusters); | |
1413 | Float_t fraction2 = Float_t(inonEdge) / Float_t(nclustersNE); | |
10757ee9 | 1414 | |
1415 | AddCluster(cluster, momenta, mdedx, padType, xcenter, dedxQ, dedxM, fraction, fraction2, dedge, parY, parZ, meanPos); | |
2acad464 | 1416 | Float_t gain = GetGain(cluster); |
ae28e92e | 1417 | if (cstream) (*cstream) << "dEdx" << |
2acad464 | 1418 | "Cl.=" << cluster << // cluster of interest |
1419 | "gain="<<gain<< // gain at cluster position | |
1420 | "P=" << momenta << // track momenta | |
1421 | "dedx=" << mdedx << // mean dedx - corrected for angle | |
1422 | "IPad=" << padType << // pad type 0..2 | |
1423 | "xc=" << xcenter << // x center of chamber | |
1424 | "dedxQ.=" << &dedxQ << // dedxQ - total charge | |
1425 | "dedxM.=" << &dedxM << // dedxM - maximal charge | |
1426 | "fraction=" << fraction << // fraction - order in statistic (0,1) | |
1427 | "fraction2=" << fraction2 << // fraction - order in statistic (0,1) | |
1428 | "dedge=" << dedge << // distance to the edge | |
1429 | "parY.=" << &parY << // line fit | |
1430 | "parZ.=" << &parZ << // line fit | |
1431 | "meanPos.=" << &meanPos << // mean position (dx, dx^2, y,y^2, z, z^2) | |
1432 | "\n"; | |
10757ee9 | 1433 | } |
0d3279d4 | 1434 | |
ae28e92e | 1435 | if (cstream) (*cstream) << "dEdxT" << |
0d3279d4 | 1436 | "P=" << momenta << // track momenta |
1437 | "npoints="<<inonEdge<< // number of points | |
1438 | "sector="<<lastSector<< // sector number | |
1439 | "dedx=" << mdedx << // mean dedx - corrected for angle | |
1440 | "IPad=" << padType << // pad type 0..2 | |
1441 | "xc=" << xcenter << // x center of chamber | |
1442 | "dedxQ.=" << &dedxQ << // dedxQ - total charge | |
1443 | "dedxM.=" << &dedxM << // dedxM - maximal charge | |
1444 | "parY.=" << &parY << // line fit | |
1445 | "parZ.=" << &parZ << // line fit | |
1446 | "meanPos.=" << &meanPos << // mean position (dx, dx^2, y,y^2, z, z^2) | |
1447 | "\n"; | |
1448 | ||
1449 | sector = lastSector; | |
1450 | npoints = inonEdge; | |
10757ee9 | 1451 | return kTRUE; |
1452 | } | |
0d3279d4 | 1453 | |
1454 | void AliTPCcalibTracksGain::AddTracklet(UInt_t sector, UInt_t padType,TVectorD &dedxQ, TVectorD &dedxM,TVectorD& parY, TVectorD& parZ, TVectorD& meanPos){ | |
1455 | // | |
1456 | // Add measured point - dedx to the fitter | |
1457 | // | |
1458 | // | |
8076baa0 | 1459 | //chain->SetAlias("dr","(250-abs(meanPos.fElements[4]))/250"); |
0d3279d4 | 1460 | //chain->SetAlias("tz","(0+abs(parZ.fElements[1]))"); |
1461 | //chain->SetAlias("ty","(0+abs(parY.fElements[1]))"); | |
1462 | //chain->SetAlias("corrg","sqrt((1+ty^2)*(1+tz^2))"); | |
8076baa0 | 1463 | //expession fast - TString *strq0 = toolkit.FitPlane(chain,"dedxQ.fElements[2]","dr++ty++tz++dr*ty++dr*tz++ty*tz++ty^2++tz^2","IPad==0",chi2,npoints,param,covar,0,100000); |
0d3279d4 | 1464 | |
1465 | Double_t xxx[100]; | |
1466 | // | |
1467 | // z and angular part | |
1468 | // | |
8076baa0 | 1469 | |
1470 | xxx[0] = (250.-TMath::Abs(meanPos[4]))/250.; | |
0d3279d4 | 1471 | xxx[1] = TMath::Abs(parY[1]); |
1472 | xxx[2] = TMath::Abs(parZ[1]); | |
1473 | xxx[3] = xxx[0]*xxx[1]; | |
1474 | xxx[4] = xxx[0]*xxx[2]; | |
1475 | xxx[5] = xxx[1]*xxx[2]; | |
1476 | xxx[6] = xxx[0]*xxx[0]; | |
8076baa0 | 1477 | xxx[7] = xxx[1]*xxx[1]; |
1478 | xxx[8] = xxx[2]*xxx[2]; | |
0d3279d4 | 1479 | // |
1480 | // chamber part | |
1481 | // | |
1482 | Int_t tsector = sector%36; | |
1483 | for (Int_t i=0;i<35;i++){ | |
8076baa0 | 1484 | xxx[9+i]=(i==tsector)?1:0; |
0d3279d4 | 1485 | } |
1486 | TLinearFitter *fitterM = fFitter0M; | |
1487 | if (padType==1) fitterM=fFitter1M; | |
1488 | if (padType==2) fitterM=fFitter2M; | |
1489 | fitterM->AddPoint(xxx,dedxM[1]); | |
1490 | // | |
1491 | TLinearFitter *fitterT = fFitter0T; | |
1492 | if (padType==1) fitterT = fFitter1T; | |
1493 | if (padType==2) fitterT = fFitter2T; | |
1494 | fitterT->AddPoint(xxx,dedxQ[1]); | |
684602c8 | 1495 | // |
1496 | TLinearFitter *dfitterM = fDFitter0M; | |
1497 | if (padType==1) dfitterM=fDFitter1M; | |
1498 | if (padType==2) dfitterM=fDFitter2M; | |
1499 | dfitterM->AddPoint(xxx,dedxM[1]); | |
1500 | // | |
1501 | TLinearFitter *dfitterT = fDFitter0T; | |
1502 | if (padType==1) dfitterT = fDFitter1T; | |
1503 | if (padType==2) dfitterT = fDFitter2T; | |
1504 | dfitterT->AddPoint(xxx,dedxQ[1]); | |
0d3279d4 | 1505 | } |
8076baa0 | 1506 | |
1507 | ||
1508 | TGraph *AliTPCcalibTracksGain::CreateAmpGraph(Int_t ipad, Bool_t qmax){ | |
1509 | // | |
1510 | // create the amplitude graph | |
1511 | // The normalized amplitudes are extrapolated to the 0 angle (y,z) and 0 drift length | |
1512 | // | |
1513 | ||
1514 | TVectorD vec; | |
1515 | if (qmax){ | |
1516 | if (ipad==0) fFitter0M->GetParameters(vec); | |
1517 | if (ipad==1) fFitter1M->GetParameters(vec); | |
1518 | if (ipad==2) fFitter2M->GetParameters(vec); | |
1519 | }else{ | |
1520 | if (ipad==0) fFitter0T->GetParameters(vec); | |
1521 | if (ipad==1) fFitter1T->GetParameters(vec); | |
1522 | if (ipad==2) fFitter2T->GetParameters(vec); | |
1523 | } | |
1524 | ||
1525 | Float_t amp[36]; | |
1526 | Float_t sec[36]; | |
1527 | for (Int_t i=0;i<35;i++){ | |
1528 | sec[i]=i; | |
1529 | amp[i]=vec[10+i]+vec[0]; | |
1530 | } | |
1531 | amp[35]=vec[0]; | |
1532 | Float_t mean = TMath::Mean(36,amp); | |
1533 | for (Int_t i=0;i<36;i++){ | |
1534 | sec[i]=i; | |
1535 | amp[i]=(amp[i]-mean)/mean; | |
1536 | } | |
1537 | TGraph *gr = new TGraph(36,sec,amp); | |
1538 | return gr; | |
1539 | } | |
1540 | ||
1541 | ||
1542 | void AliTPCcalibTracksGain::UpdateClusterParam(AliTPCClusterParam* clparam){ | |
1543 | // | |
1544 | // SetQ normalization parameters | |
1545 | // | |
1546 | // void SetQnorm(Int_t ipad, Int_t itype, TVectorD * norm); | |
1547 | ||
1548 | TVectorD vec; | |
684602c8 | 1549 | |
8076baa0 | 1550 | // |
684602c8 | 1551 | fDFitter0T->Eval(); |
1552 | fDFitter1T->Eval(); | |
1553 | fDFitter2T->Eval(); | |
1554 | fDFitter0M->Eval(); | |
1555 | fDFitter1M->Eval(); | |
1556 | fDFitter2M->Eval(); | |
1557 | fDFitter0T->GetParameters(vec); | |
8076baa0 | 1558 | clparam->SetQnorm(0,0,&vec); |
684602c8 | 1559 | fDFitter1T->GetParameters(vec); |
8076baa0 | 1560 | clparam->SetQnorm(1,0,&vec); |
684602c8 | 1561 | fDFitter2T->GetParameters(vec); |
8076baa0 | 1562 | clparam->SetQnorm(2,0,&vec); |
1563 | // | |
684602c8 | 1564 | fDFitter0M->GetParameters(vec); |
8076baa0 | 1565 | clparam->SetQnorm(0,1,&vec); |
684602c8 | 1566 | fDFitter1M->GetParameters(vec); |
8076baa0 | 1567 | clparam->SetQnorm(1,1,&vec); |
684602c8 | 1568 | fDFitter2M->GetParameters(vec); |
8076baa0 | 1569 | clparam->SetQnorm(2,1,&vec); |
1570 | // | |
1571 | ||
1572 | } | |
b8601924 | 1573 | |
1574 | ||
1575 | void AliTPCcalibTracksGain::Analyze(){ | |
1576 | ||
1577 | Evaluate(); | |
1578 | ||
1579 | } | |
1580 | ||
1581 | ||
da7d274e | 1582 | |
1583 | TVectorD * AliTPCcalibTracksGain::MakeQPosNorm(TTree * chain0, Int_t ipad, Bool_t isMax, Int_t maxPoints, Int_t verbose){ | |
1584 | // | |
1585 | // Input parameters | |
1586 | // chain0 - the tree with information -Debug stream | |
1587 | // ipad - 0 IROC | |
1588 | // - 1 OROC medium | |
1589 | // - 2 OROC LONG | |
1590 | // isMax - kFALSE - total charge param | |
1591 | // kTRUE - Max charge param | |
1592 | // | |
1593 | // maxPoints - number of points for fit | |
1594 | // | |
1595 | // verbose - | |
1596 | // | |
1597 | /* e.g | |
1598 | ipad=0 | |
1599 | isMax=kTRUE; | |
1600 | maxPoints=1000000; | |
1601 | */ | |
1602 | // Make Q normalization as function of following parameters | |
1603 | // 1 - dp - relative pad position | |
1604 | // 2 - dt - relative time position | |
1605 | // 3 - di - drift length (norm to 1); | |
1606 | // 4 - dq0 - Tot/Max charge | |
1607 | // 5 - dq1 - Max/Tot charge | |
1608 | // 6 - sy - sigma y - shape | |
1609 | // 7 - sz - sigma z - shape | |
1610 | // | |
1611 | // Coeficient of Taylor expansion fitted | |
1612 | // Fit parameters returned as TVectorD | |
1613 | // Fit parameters to be used in corresponding correction function | |
1614 | // in AliTPCclusterParam | |
1615 | // | |
1616 | // | |
1617 | TStatToolkit toolkit; | |
1618 | Double_t chi2; | |
1619 | TVectorD fitParam; | |
1620 | TMatrixD covMatrix; | |
1621 | Int_t npoints; | |
1622 | TCut cutA("dedge>3&&fraction2<0.7"); | |
1623 | chain0->SetAlias("dp","((Cl.fPad-int(Cl.fPad)-0.5)/0.5)"); | |
1624 | chain0->SetAlias("dt","((Cl.fTimeBin-int(Cl.fTimeBin)-0.5)/0.5)"); | |
1625 | chain0->SetAlias("di","(sqrt(1.-abs(Cl.fZ)/250.))"); | |
1626 | chain0->SetAlias("dq0","(0.2*(Cl.fQ+2)/(Cl.fMax+2))"); | |
1627 | chain0->SetAlias("dq1","(5*(Cl.fMax+2)/(Cl.fQ+2))"); | |
1628 | chain0->SetAlias("sy","(0.32/sqrt(0.01^2+Cl.fSigmaY2))"); | |
1629 | chain0->SetAlias("sz","(0.32/sqrt(0.01^2+Cl.fSigmaZ2))"); | |
1630 | // | |
1631 | TString fstring=""; | |
1632 | fstring+="dp++"; //1 | |
1633 | fstring+="dt++"; //2 | |
1634 | fstring+="dp*dp++"; //3 | |
1635 | fstring+="dt*dt++"; //4 | |
1636 | fstring+="dt*dt*dt++"; //5 | |
1637 | fstring+="dp*dt++"; //6 | |
1638 | fstring+="dp*dt*dt++"; //7 | |
1639 | fstring+="(dq0)++"; //8 | |
1640 | fstring+="(dq1)++"; //9 | |
1641 | // | |
1642 | // | |
1643 | fstring+="dp*dp*(di)++"; //10 | |
1644 | fstring+="dt*dt*(di)++"; //11 | |
1645 | fstring+="dp*dp*sy++"; //12 | |
1646 | fstring+="dt*sz++"; //13 | |
1647 | fstring+="dt*dt*sz++"; //14 | |
1648 | fstring+="dt*dt*dt*sz++"; //15 | |
1649 | // | |
1650 | fstring+="dp*dp*1*sy*sz++"; //16 | |
1651 | fstring+="dt*sy*sz++"; //17 | |
1652 | fstring+="dt*dt*sy*sz++"; //18 | |
1653 | fstring+="dt*dt*dt*sy*sz++"; //19 | |
1654 | // | |
1655 | fstring+="dp*dp*(dq0)++"; //20 | |
1656 | fstring+="dt*1*(dq0)++"; //21 | |
1657 | fstring+="dt*dt*(dq0)++"; //22 | |
1658 | fstring+="dt*dt*dt*(dq0)++"; //23 | |
1659 | // | |
1660 | fstring+="dp*dp*(dq1)++"; //24 | |
1661 | fstring+="dt*(dq1)++"; //25 | |
1662 | fstring+="dt*dt*(dq1)++"; //26 | |
1663 | fstring+="dt*dt*dt*(dq1)++"; //27 | |
1664 | ||
1665 | TString var; | |
1666 | if (isMax) var = "Cl.fMax/gain/dedxM.fElements[2]"; | |
1667 | if (!isMax) var = "Cl.fQ/gain/dedxQ.fElements[2]"; | |
1668 | TString cutP="IPad=="; | |
1669 | cutP+=ipad; | |
1670 | // | |
1671 | TString *strq0 = toolkit.FitPlane(chain0,var.Data(),fstring.Data(), cutP.Data()+cutA, chi2,npoints,fitParam,covMatrix,-1,0,maxPoints); | |
1672 | // | |
1673 | // | |
1674 | if (verbose){ | |
1675 | printf("Chi2/npoints = %f",TMath::Sqrt(chi2/npoints)); | |
1676 | printf("\nFit function\n:%s\n",strq0->Data()); | |
1677 | } | |
1678 | TVectorD *vec = new TVectorD(fitParam); | |
1679 | return vec; | |
1680 | } | |
1681 | ||
1682 | void AliTPCcalibTracksGain::MakeQPosNormAll(TTree * chain, AliTPCClusterParam * param, Int_t maxPoints, Int_t verbose){ | |
1683 | // | |
1684 | // Fill the content of the of the AliTPCclusterParam | |
1685 | // with fitted values of corrections | |
1686 | // | |
2e5bcb67 | 1687 | param->fPosQTnorm[0] = MakeQPosNorm(chain,0,kTRUE,100000,kTRUE); |
1688 | param->fPosQTnorm[1] = MakeQPosNorm(chain,1,kTRUE,100000,kTRUE); | |
1689 | param->fPosQTnorm[2] = MakeQPosNorm(chain,1,kTRUE,100000,kTRUE); | |
1690 | // | |
1691 | param->fPosQMnorm[0] = MakeQPosNorm(chain,0,kFALSE,100000,kTRUE); | |
1692 | param->fPosQMnorm[1] = MakeQPosNorm(chain,1,kFALSE,100000,kTRUE); | |
1693 | param->fPosQMnorm[2] = MakeQPosNorm(chain,2,kFALSE,100000,kTRUE); | |
da7d274e | 1694 | } |
1695 | ||
1696 | ||
1697 | ||
1698 | /* | |
1699 | ||
1700 | Position correction fit: | |
1701 | // | |
1702 | TStatToolkit toolkit; | |
1703 | Double_t chi2; | |
1704 | TVectorD fitParam; | |
1705 | TMatrixD covMatrix; | |
1706 | Int_t npoints; | |
1707 | // | |
1708 | TCut cutA("dedge>3&&fraction2<0.7"); | |
1709 | chain0->SetAlias("dp","((Cl.fPad-int(Cl.fPad)-0.5)/0.5)"); | |
1710 | chain0->SetAlias("dt","((Cl.fTimeBin-int(Cl.fTimeBin)-0.5)/0.5)"); | |
1711 | chain0->SetAlias("di","(sqrt(1.-abs(Cl.fZ)/250.))"); | |
1712 | chain0->SetAlias("dq0","(0.2*(Cl.fQ+2)/(Cl.fMax+2))"); | |
1713 | chain0->SetAlias("dq1","(5*(Cl.fMax+2)/(Cl.fQ+2))"); | |
1714 | chain0->SetAlias("sy","(0.2/sqrt(0.01^2+Cl.fSigmaY2))"); | |
1715 | chain0->SetAlias("sz","(0.2/sqrt(0.01^2+Cl.fSigmaZ2))"); | |
1716 | ||
1717 | TString fstring=""; | |
1718 | ||
1719 | fstring+="dp++"; //1 | |
1720 | fstring+="dt++"; //2 | |
1721 | fstring+="dp*dp++"; //3 | |
1722 | fstring+="dt*dt++"; //4 | |
1723 | fstring+="dt*dt*dt++"; //5 | |
1724 | fstring+="dp*dt++"; //6 | |
1725 | fstring+="dp*dt*dt++"; //7 | |
1726 | fstring+="(dq0)++"; //8 | |
1727 | fstring+="(dq1)++"; //9 | |
1728 | // | |
1729 | // | |
1730 | fstring+="dp*dp*(di)++"; //10 | |
1731 | fstring+="dt*dt*(di)++"; //11 | |
1732 | fstring+="dp*dp*sy++"; //12 | |
1733 | fstring+="dt*sz++"; //13 | |
1734 | fstring+="dt*dt*sz++"; //14 | |
1735 | fstring+="dt*dt*dt*sz++"; //15 | |
1736 | // | |
1737 | fstring+="dp*dp*1*sy*sz++"; //16 | |
1738 | fstring+="dt*sy*sz++"; //17 | |
1739 | fstring+="dt*dt*sy*sz++"; //18 | |
1740 | fstring+="dt*dt*dt*sy*sz++"; //19 | |
1741 | // | |
1742 | fstring+="dp*dp*(dq0)++"; //20 | |
1743 | fstring+="dt*1*(dq0)++"; //21 | |
1744 | fstring+="dt*dt*(dq0)++"; //22 | |
1745 | fstring+="dt*dt*dt*(dq0)++"; //23 | |
1746 | ||
1747 | fstring+="dp*dp*(dq1)++"; //24 | |
1748 | fstring+="dt*(dq1)++"; //25 | |
1749 | fstring+="dt*dt*(dq1)++"; //26 | |
1750 | fstring+="dt*dt*dt*(dq1)++"; //27 | |
1751 | ||
1752 | ||
1753 | TString *strq0 = toolkit.FitPlane(chain0,"Cl.fMax/gain/dedxM.fElements[2]",fstring->Data(), "IPad==0"+cutA, chi2,npoints,fitParam,covMatrix,-1,0,200000); | |
1754 | TString *strqt0 = toolkit.FitPlane(chain0,"Cl.fQ/gain/dedxQ.fElements[2]",fstring->Data(), "IPad==0"+cutA, chi2,npoints,fitParam,covMatrix,-1,0,200000); | |
1755 | ||
1756 | chain0->SetAlias("qcorM0",strq0->Data()); | |
1757 | chain0->SetAlias("qcorT0",strqt0->Data()); | |
1758 | //chain0->SetAlias("mmqcorM0","min(max(qcorM0,0.75),1.15)"); | |
1759 | chain0->Draw("(Cl.fMax/gain/dedxM.fElements[2]):min(max(qcorM0,0.75),1.15)","IPad==0"+cutA,"prof",100000) | |
1760 | ||
1761 | fraction05 - | |
1762 | sigma 0.2419 | |
1763 | sigma fit 0.2302 | |
1764 | sigma fit with shape 0.2257 | |
1765 | fraction 07 | |
1766 | qtot sigma 0.322 | |
1767 | qmax sigma 0.292 | |
1768 | qmax sigma fit 0.2702 | |
1769 | qmax sigma fit+ratio 0.2638 | |
1770 | ||
1771 | */ | |
1772 |