1 /**************************************************************************
2 * Copyright(c) 1998-1999, ALICE Experiment at CERN, All rights reserved. *
4 * Author: The ALICE Off-line Project. *
5 * Contributors are mentioned in the code where appropriate. *
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 **************************************************************************/
16 ////////////////////////////////////////////////////////////////////////////
19 // Gain calibration using tracks
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
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
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
39 // === Calibration class for gain calibration using tracks ===
43 // A 6-parametric parabolic function
45 // G(x, y) = p0 + p1*x + p2*y + p3*x^2 + p4*y^2 + p5 * x*y
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.
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:
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
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.
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.
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.
89 // The visualization methods CreateFitCalPad() and CreateFitCalROC()
90 // are straight forward to use.
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)
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.
110 ////////////////////////////////////////////////////////////////////////////
114 gSystem->Load("libANALYSIS");
115 gSystem->Load("libSTAT");
116 gSystem->Load("libTPCcalib");
118 TFile fcalib("CalibObjects.root");
119 TObjArray * array = (TObjArray*)fcalib.Get("TPCCalib");
120 AliTPCcalibTracksGain * gain = ( AliTPCcalibTracksGain *)array->FindObject("calibTracksGain");
123 // Angular and drift correction
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)
130 // Make visual Tree - compare with Kr calibration
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");
135 TFile fkr("/u/miranov/GainMap.root");
136 AliTPCCalPad *gainKr = fkr.Get("GainMap"); fkr->SetName("KrGain");
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");
147 // Simple session using the debug streamers
150 gSystem->AddIncludePath("-I$ALICE_ROOT/TPC/macros");
151 gROOT->LoadMacro("$ALICE_ROOT/TPC/macros/AliXRDPROOFtoolkit.cxx+")
152 AliXRDPROOFtoolkit tool;
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);
161 chain2->SetAlias("k1","1/0.855");
162 chain2->SetAlias("k0","1/0.9928");
163 chain2->SetAlias("k2","1/1.152");
169 #include "AliTPCcalibTracksGain.h"
172 #include <TPDGCode.h>
175 #include "TMatrixD.h"
176 #include "TTreeStream.h"
178 #include "AliTPCParamSR.h"
179 #include "AliTPCClusterParam.h"
180 #include "AliTrackPointArray.h"
183 #include <TLinearFitter.h>
184 #include <TTreeStream.h>
186 #include <TCollection.h>
187 #include <TIterator.h>
188 #include <TProfile.h>
189 #include <TProfile2D.h>
191 #include <TStatToolkit.h>
197 #include "AliMathBase.h"
199 #include "AliTPCROC.h"
200 #include "AliTPCParamSR.h"
201 #include "AliTPCCalROC.h"
202 #include "AliTPCCalPad.h"
203 #include "AliTPCClusterParam.h"
204 #include "AliTPCcalibDB.h"
206 #include "AliTracker.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"
215 #include "TStatToolkit.h"
221 #include "AliESDEvent.h"
225 TFile f("TPCCalibTracksGain.root")
227 gSystem->Load("libPWGPP.so")
233 TString * str = comp.FitPlane("Cl.fQ/dedxQ.fElements[0]","Cl.fY++Cl.fX","Cl.fDetector<36",chi2,vec,mat)
237 ClassImp(AliTPCcalibTracksGain)
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();
244 AliTPCcalibTracksGain::AliTPCcalibTracksGain() :
246 fCuts(0), // cuts that are used for sieving the tracks used for calibration
250 fSimpleFitter(0), // simple fitter for short pads
251 fSqrtFitter(0), // sqrt fitter for medium pads
252 fLogFitter(0), // log fitter for long pads
254 fFitter0M(0), // fitting of the atenuation, angular correction, and mean chamber gain
255 fFitter1M(0), // fitting of the atenuation, angular correction, and mean chamber gain
256 fFitter2M(0), // fitting of the atenuation, angular correction, and mean chamber gain
257 fFitter0T(0), // fitting of the atenuation, angular correction, and mean chamber gain
258 fFitter1T(0), // fitting of the atenuation, angular correction, and mean chamber gain
259 fFitter2T(0), // fitting of the atenuation, angular correction, and mean chamber gain
261 fDFitter0M(0), // fitting of the atenuation, angular correction
262 fDFitter1M(0), // fitting of the atenuation, angular correction
263 fDFitter2M(0), // fitting of the atenuation, angular correction
264 fDFitter0T(0), // fitting of the atenuation, angular correction
265 fDFitter1T(0), // fitting of the atenuation, angular correction
266 fDFitter2T(0), // fitting of the atenuation, angular correction
268 fSingleSectorFitter(0), // just for debugging
272 fTotalTracks(0), // just for debugging
273 fAcceptedTracks(0) // just for debugging
277 // Default constructor.
281 AliTPCcalibTracksGain::AliTPCcalibTracksGain(const AliTPCcalibTracksGain& obj) :
282 AliTPCcalibBase(obj),
283 fCuts(obj.fCuts), // cuts that are used for sieving the tracks used for calibration
287 fSimpleFitter(obj.fSimpleFitter), // simple fitter for short pads
288 fSqrtFitter(obj.fSqrtFitter), // sqrt fitter for medium pads
289 fLogFitter(obj.fLogFitter), // log fitter for long pads
290 fFitter0M(obj.fFitter0M),
291 fFitter1M(obj.fFitter1M),
292 fFitter2M(obj.fFitter2M),
293 fFitter0T(obj.fFitter0T),
294 fFitter1T(obj.fFitter1T),
295 fFitter2T(obj.fFitter2T),
297 fDFitter0M(obj.fDFitter0M),
298 fDFitter1M(obj.fDFitter1M),
299 fDFitter2M(obj.fDFitter2M),
300 fDFitter0T(obj.fDFitter0T),
301 fDFitter1T(obj.fDFitter1T),
302 fDFitter2T(obj.fDFitter2T),
303 fSingleSectorFitter(obj.fSingleSectorFitter), // just for debugging
307 fTotalTracks(obj.fTotalTracks), // just for debugging
308 fAcceptedTracks(obj.fAcceptedTracks) // just for debugging
316 AliTPCcalibTracksGain& AliTPCcalibTracksGain::operator=(const AliTPCcalibTracksGain& rhs) {
318 // Assignment operator.
322 TNamed::operator=(rhs);
323 fSimpleFitter = new AliTPCFitPad(*(rhs.fSimpleFitter));
324 fSqrtFitter = new AliTPCFitPad(*(rhs.fSqrtFitter));
325 fLogFitter = new AliTPCFitPad(*(rhs.fLogFitter));
326 fSingleSectorFitter = new AliTPCFitPad(*(rhs.fSingleSectorFitter));
327 fCuts = new AliTPCcalibTracksCuts(*(rhs.fCuts));
332 AliTPCcalibTracksGain::AliTPCcalibTracksGain(const char* name, const char* title, AliTPCcalibTracksCuts* cuts) :
334 fCuts(0), // cuts that are used for sieving the tracks used for calibration
338 fSimpleFitter(0), // simple fitter for short pads
339 fSqrtFitter(0), // sqrt fitter for medium pads
340 fLogFitter(0), // log fitter for long pads
341 fFitter0M(0), // fitting of the atenuation, angular correction, and mean chamber gain
342 fFitter1M(0), // fitting of the atenuation, angular correction, and mean chamber gain
343 fFitter2M(0), // fitting of the atenuation, angular correction, and mean chamber gain
344 fFitter0T(0), // fitting of the atenuation, angular correction, and mean chamber gain
345 fFitter1T(0), // fitting of the atenuation, angular correction, and mean chamber gain
346 fFitter2T(0), // fitting of the atenuation, angular correction, and mean chamber gain
348 fDFitter0M(0), // fitting of the atenuation, angular correction
349 fDFitter1M(0), // fitting of the atenuation, angular correction
350 fDFitter2M(0), // fitting of the atenuation, angular correction
351 fDFitter0T(0), // fitting of the atenuation, angular correction
352 fDFitter1T(0), // fitting of the atenuation, angular correction
353 fDFitter2T(0), // fitting of the atenuation, angular correction
354 fSingleSectorFitter(0), // just for debugging
358 fTotalTracks(0), // just for debugging
359 fAcceptedTracks(0) // just for debugging
365 this->SetNameTitle(name, title);
368 // Fitter initialization
370 fSimpleFitter = new AliTPCFitPad(8, "hyp7", "");
371 fSqrtFitter = new AliTPCFitPad(8, "hyp7", "");
372 fLogFitter = new AliTPCFitPad(8, "hyp7", "");
373 fSingleSectorFitter = new AliTPCFitPad(8, "hyp7", "");
375 fFitter0M = new TLinearFitter(45,"hyp44");
376 fFitter1M = new TLinearFitter(45,"hyp44");
377 fFitter2M = new TLinearFitter(45,"hyp44");
378 fFitter0T = new TLinearFitter(45,"hyp44");
379 fFitter1T = new TLinearFitter(45,"hyp44");
380 fFitter2T = new TLinearFitter(45,"hyp44");
382 fDFitter0M = new TLinearFitter(10,"hyp9");
383 fDFitter1M = new TLinearFitter(10,"hyp9");
384 fDFitter2M = new TLinearFitter(10,"hyp9");
385 fDFitter0T = new TLinearFitter(10,"hyp9");
386 fDFitter1T = new TLinearFitter(10,"hyp9");
387 fDFitter2T = new TLinearFitter(10,"hyp9");
390 fFitter0M->StoreData(kFALSE);
391 fFitter1M->StoreData(kFALSE);
392 fFitter2M->StoreData(kFALSE);
393 fFitter0T->StoreData(kFALSE);
394 fFitter1T->StoreData(kFALSE);
395 fFitter2T->StoreData(kFALSE);
397 fDFitter0M->StoreData(kFALSE);
398 fDFitter1M->StoreData(kFALSE);
399 fDFitter2M->StoreData(kFALSE);
400 fDFitter0T->StoreData(kFALSE);
401 fDFitter1T->StoreData(kFALSE);
402 fDFitter2T->StoreData(kFALSE);
405 // just for debugging -counters
411 AliTPCcalibTracksGain::~AliTPCcalibTracksGain() {
416 Info("Destructor",":");
417 if (fSimpleFitter) delete fSimpleFitter;
418 if (fSqrtFitter) delete fSqrtFitter;
419 if (fLogFitter) delete fLogFitter;
420 if (fSingleSectorFitter) delete fSingleSectorFitter;
427 void AliTPCcalibTracksGain::Terminate(){
429 // Evaluate fitters and close the debug stream.
430 // Also move or copy the debug stream, if a debugStreamPrefix is provided.
434 AliTPCcalibBase::Terminate();
439 void AliTPCcalibTracksGain::Process(AliTPCseed* seed) {
441 // Main method to be called when a new seed is supposed to be processed
442 // and be used for gain calibration. Its quality is checked before it
448 if (!fCuts->AcceptTrack(seed)) return;
452 static Bool_t doinit= kTRUE;
454 fSimpleFitter = new AliTPCFitPad(8, "hyp7", "");
455 fSqrtFitter = new AliTPCFitPad(8, "hyp7", "");
456 fLogFitter = new AliTPCFitPad(8, "hyp7", "");
457 fSingleSectorFitter = new AliTPCFitPad(8, "hyp7", "");
459 fFitter0M = new TLinearFitter(45,"hyp44");
460 fFitter1M = new TLinearFitter(45,"hyp44");
461 fFitter2M = new TLinearFitter(45,"hyp44");
462 fFitter0T = new TLinearFitter(45,"hyp44");
463 fFitter1T = new TLinearFitter(45,"hyp44");
464 fFitter2T = new TLinearFitter(45,"hyp44");
466 fDFitter0M = new TLinearFitter(10,"hyp9");
467 fDFitter1M = new TLinearFitter(10,"hyp9");
468 fDFitter2M = new TLinearFitter(10,"hyp9");
469 fDFitter0T = new TLinearFitter(10,"hyp9");
470 fDFitter1T = new TLinearFitter(10,"hyp9");
471 fDFitter2T = new TLinearFitter(10,"hyp9");
480 Long64_t AliTPCcalibTracksGain::Merge(TCollection *list) {
482 // Merge() merges the results of all AliTPCcalibTracksGain objects contained in
483 // list, thus allowing a distributed computation of several files, e.g. on PROOF.
484 // The merged results are merged with the data members of the AliTPCcalibTracksGain
485 // object used for calling the Merge method.
486 // The return value is 0 /*the total number of tracks used for calibration*/ if the merge
487 // is successful, otherwise it is -1.
490 if (!list || list->IsEmpty()) return -1;
492 if (!fSimpleFitter) fSimpleFitter = new AliTPCFitPad(8, "hyp7", "");
493 if (!fSqrtFitter) fSqrtFitter = new AliTPCFitPad(8, "hyp7", "");
494 if (!fLogFitter) fLogFitter = new AliTPCFitPad(8, "hyp7", "");
495 if (!fSingleSectorFitter) fSingleSectorFitter = new AliTPCFitPad(8, "hyp7", "");
499 TIterator* iter = list->MakeIterator();
500 AliTPCcalibTracksGain* cal = 0;
502 while ((cal = (AliTPCcalibTracksGain*)iter->Next())) {
503 if (!cal->InheritsFrom(AliTPCcalibTracksGain::Class())) {
504 //Error("Merge","Attempt to add object of class %s to a %s", cal->ClassName(), this->ClassName());
513 Float_t AliTPCcalibTracksGain::GetGain(AliTPCclusterMI *cluster){
518 AliTPCCalPad * gainMap = AliTPCcalibDB::Instance()->GetDedxGainFactor();
520 AliTPCCalROC * roc = gainMap->GetCalROC(cluster->GetDetector());
521 gainPad = roc->GetValue(cluster->GetRow(), TMath::Nint(cluster->GetPad()));
526 Float_t AliTPCcalibTracksGain::GetMaxNorm(AliTPCclusterMI * cluster, Float_t ky, Float_t kz){
528 // Get normalized amplituded if the gain map provided
530 AliTPCClusterParam * parcl = AliTPCcalibDB::Instance()->GetClusterParam();
532 parcl->QmaxCorrection(cluster->GetDetector(), cluster->GetRow(),cluster->GetPad(),
533 cluster->GetTimeBin(),ky,kz,0.5,0.5,1.6);
535 return GetGain(cluster)*maxNorm;
539 Float_t AliTPCcalibTracksGain::GetQNorm(AliTPCclusterMI * cluster, Float_t ky, Float_t kz){
541 // Get normalized amplituded if the gain map provided
543 AliTPCClusterParam * parcl = AliTPCcalibDB::Instance()->GetClusterParam();
544 Float_t totNorm = parcl->QtotCorrection(cluster->GetDetector(), cluster->GetRow(),cluster->GetPad(), cluster->GetTimeBin(),ky,kz,0.5,0.5,cluster->GetQ(),2.5,1.6);
545 return GetGain(cluster)*totNorm;
550 void AliTPCcalibTracksGain::Add(AliTPCcalibTracksGain* cal) {
552 // Adds another AliTPCcalibTracksGain object to this object.
555 fSimpleFitter->Add(cal->fSimpleFitter);
556 fSqrtFitter->Add(cal->fSqrtFitter);
557 fLogFitter->Add(cal->fLogFitter);
558 fSingleSectorFitter->Add(cal->fSingleSectorFitter);
562 if (cal->fFitter0M->GetNpoints()>0) fFitter0M->Add(cal->fFitter0M);
563 if (cal->fFitter1M->GetNpoints()>0) fFitter1M->Add(cal->fFitter1M);
564 if (cal->fFitter2M->GetNpoints()>0) fFitter2M->Add(cal->fFitter2M);
565 if (cal->fFitter0T->GetNpoints()>0) fFitter0T->Add(cal->fFitter0T);
566 if (cal->fFitter1T->GetNpoints()>0) fFitter1T->Add(cal->fFitter1T);
567 if (cal->fFitter2T->GetNpoints()>0) fFitter2T->Add(cal->fFitter2T);
569 if (cal->fDFitter0M->GetNpoints()>0) fDFitter0M->Add(cal->fDFitter0M);
570 if (cal->fDFitter1M->GetNpoints()>0) fDFitter1M->Add(cal->fDFitter1M);
571 if (cal->fDFitter2M->GetNpoints()>0) fDFitter2M->Add(cal->fDFitter2M);
572 if (cal->fDFitter0T->GetNpoints()>0) fDFitter0T->Add(cal->fDFitter0T);
573 if (cal->fDFitter1T->GetNpoints()>0) fDFitter1T->Add(cal->fDFitter1T);
574 if (cal->fDFitter2T->GetNpoints()>0) fDFitter2T->Add(cal->fDFitter2T);
577 // just for debugging, remove me
578 fTotalTracks += cal->fTotalTracks;
579 fAcceptedTracks += cal->fAcceptedTracks;
583 void AliTPCcalibTracksGain::AddTrack(AliTPCseed* seed) {
585 // The clusters making up the track (seed) are added to various fit functions.
586 // See AddCluster(...) for more detail.
595 void AliTPCcalibTracksGain::AddCluster(AliTPCclusterMI* cluster, Float_t /*momenta*/, Float_t/* mdedx*/, Int_t padType,
596 Float_t xcenter, TVectorD& dedxQ, TVectorD& /*dedxM*/, Float_t /*fraction*/, Float_t fraction2, Float_t dedge,
597 TVectorD& parY, TVectorD& parZ, TVectorD& meanPos) {
599 // Adds cluster to the appropriate fitter for later analysis.
600 // The charge used for the fit is the maximum charge for this specific cluster or the
601 // accumulated charge per cluster, depending on the value of fgkUseTotalCharge.
602 // Depending on the pad size where the cluster is registered, the value will be put in
603 // the appropriate fitter. Furthermore, for each pad size three different types of fitters
604 // are used. The fit functions are the same for all fitters (parabolic functions), but the value
605 // added to each fitter is different. The simple fitter gets the charge plugged in as is, the sqrt fitter
606 // gets the square root of the charge, and the log fitter gets fgkM*(1+q/fgkM), where q is the original charge
610 Float_t kfraction = 0.7;
614 // Where to put selection on threshold?
615 // Defined by the Q/dEdxT variable - see debug streamer:
617 // Debug stream variables: (Where tu cut ?)
618 // chain0->Draw("Cl.fQ/dedxQ.fElements[1]>>his(100,0,3)","fraction2<0.6&&dedge>3","",1000000);
620 // chain0->Draw("Cl.fMax/dedxM.fElements[1]>>his(100,0,3)","fraction2<0.6&&dedge>3","",1000000)
622 // chain0->Draw("Cl.fQ/dedxQ.fElements[2]>>his(100,0,3)","fraction2<0.7&&dedge>3","",1000000)
624 // chain0->Draw("Cl.fMax/dedxM.fElements[2]>>his(100,0,3)","fraction2<0.7&&dedge>3","",1000000)
626 // chain0->Draw("Cl.fQ/dedxQ.fElements[3]>>his(100,0,3)","fraction2<0.8&&dedge>3","",1000000)
629 // chain0->Draw("Cl.fQ/dedxQ.fElements[4]>>his(100,0,3)","fraction2<0.9&&dedge>3","",1000000)
632 // Fraction choosen 0.7
635 Error("AddCluster", "Cluster not valid.");
639 if (dedge < kedge) return;
640 if (fraction2 > kfraction) return;
642 //Int_t padType = GetPadType(cluster->GetX());
644 //Double_t centerPad[2] = {0};
645 //AliTPCFitPad::GetPadRegionCenterLocal(padType, centerPad);
646 //xx[0] = cluster->GetX() - centerPad[0];
647 //xx[1] = cluster->GetY() - centerPad[1];
648 xx[0] = cluster->GetX() - xcenter;
649 xx[1] = cluster->GetY();
650 xx[2] = xx[0] * xx[0];
651 xx[3] = xx[1] * xx[1];
652 xx[4] = xx[0] * xx[1];
653 xx[5] = TMath::Abs(cluster->GetZ()) - TMath::Abs(meanPos[4]);
654 xx[6] = xx[5] * xx[5];
659 Int_t segment = cluster->GetDetector() % 36;
660 Double_t q = fgkUseTotalCharge ?
661 ((Double_t)(cluster->GetQ()/GetQNorm(cluster,parY[1], parZ[1]))) : ((Double_t)(cluster->GetMax()/GetMaxNorm(cluster,parY[1], parZ[1])));
663 // correct charge by normalising to mean charge per track
666 // just for debugging
668 Double_t sqrtQ = TMath::Sqrt(q);
669 Double_t logQ = fgkM * TMath::Log(1 + q / fgkM);
670 TLinearFitter * fitter =0;
672 fitter = fSimpleFitter->GetFitter(segment, padType);
673 fitter->AddPoint(xx, q);
675 fitter = fSqrtFitter->GetFitter(segment, padType);
676 fitter->AddPoint(xx, sqrtQ);
678 fitter = fLogFitter->GetFitter(segment, padType);
679 fitter->AddPoint(xx, logQ);
681 fitter=fSingleSectorFitter->GetFitter(0, padType);
682 fitter->AddPoint(xx, q);
686 void AliTPCcalibTracksGain::Evaluate(Bool_t robust, Double_t frac) {
688 // Evaluates all fitters contained in this object.
689 // If the robust option is set to kTRUE a robust fit is performed with frac as
690 // the minimal fraction of good points (see TLinearFitter::EvalRobust for details).
691 // Beware: Robust fitting is much slower!
694 fSimpleFitter->Evaluate(robust, frac);
695 fSqrtFitter->Evaluate(robust, frac);
696 fLogFitter->Evaluate(robust, frac);
697 fSingleSectorFitter->Evaluate(robust, frac);
713 AliTPCCalPad* AliTPCcalibTracksGain::CreateFitCalPad(UInt_t fitType, Bool_t undoTransformation, Bool_t normalizeToPadSize) {
715 // Creates the calibration object AliTPCcalPad using fitted parameterization
718 for (UInt_t iSector = 0; iSector < 72; iSector++)
719 tpc.Add(CreateFitCalROC(iSector, fitType, undoTransformation, normalizeToPadSize));
720 return new AliTPCCalPad(&tpc);
723 AliTPCCalROC* AliTPCcalibTracksGain::CreateFitCalROC(UInt_t sector, UInt_t fitType, Bool_t undoTransformation, Bool_t normalizeToPadSize) {
725 // Create the AliTPCCalROC with the values per pad
726 // sector - sector of interest
732 GetParameters(sector % 36, 0, fitType, par);
733 return CreateFitCalROC(sector, 0, par, fitType, undoTransformation, normalizeToPadSize);
736 GetParameters(sector % 36, 1, fitType, par);
737 AliTPCCalROC* roc1 = CreateFitCalROC(sector, 1, par, fitType, undoTransformation, normalizeToPadSize);
738 GetParameters(sector % 36, 2, fitType, par);
739 AliTPCCalROC* roc2 = CreateFitCalROC(sector, 2, par, fitType, undoTransformation, normalizeToPadSize);
740 AliTPCCalROC* roc3 = CreateCombinedCalROC(roc1, roc2);
747 AliTPCCalROC* AliTPCcalibTracksGain::CreateFitCalROC(UInt_t sector, UInt_t padType, TVectorD &fitParam, UInt_t fitType, Bool_t undoTransformation, Bool_t normalizeToPadSize) {
749 // This function is essentially a copy of AliTPCCalROC::CreateGlobalFitCalROC(...), with the
750 // modifications, that the center of the region of same pad size is used as the origin
751 // of the fit function instead of the center of the ROC.
752 // The possibility of a linear fit is removed as well because it is not needed.
753 // Only values for pads with the given pad size are calculated, the rest is 0.
754 // Set undoTransformation for undoing the transformation that was applied to the
755 // charge values before they were put into the fitter (thus allowing comparison to the original
756 // charge values). For fitType use 0 for the simple fitter, 1 for the sqrt fitter, 2 for the log fitter.
757 // If normalizeToPadSize is true, the values are normalized to the pad size.
758 // Please be aware, that you even need to specify the fitType if you want to normalize to the pad size without
759 // undoing the transformation (because normalizing involves undoing the trafo first, then normalizing, then
760 // applying the trafo again).
761 // Please note: The normalization to the pad size is a simple linear scaling with the pad length, which
762 // actually doesn't describe reality!
766 Double_t centerPad[2] = {0};
767 Float_t localXY[3] = {0};
768 AliTPCROC* tpcROC = AliTPCROC::Instance();
769 if ((padType == 0 && sector >= tpcROC->GetNInnerSector()) || (padType > 0 && sector < tpcROC->GetNInnerSector()) || sector >= tpcROC->GetNSector())
771 AliTPCCalROC* lROCfitted = new AliTPCCalROC(sector);
772 //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
778 endRow = lROCfitted->GetNrows();
786 endRow = lROCfitted->GetNrows();
790 AliTPCFitPad::GetPadRegionCenterLocal(padType, centerPad);
792 for (UInt_t irow = startRow; irow < endRow; irow++) {
793 for (UInt_t ipad = 0; ipad < lROCfitted->GetNPads(irow); ipad++) {
794 tpcROC->GetPositionLocal(sector, irow, ipad, localXY); // calculate position localXY by pad and row number
795 dlx = localXY[0] - centerPad[0];
796 dly = localXY[1] - centerPad[1];
797 value = fitParam[0] + fitParam[1]*dlx + fitParam[2]*dly + fitParam[3]*dlx*dlx + fitParam[4]*dly*dly + fitParam[5]*dlx*dly;
799 // Let q' = value be the transformed value without any pad size corrections,
800 // let T be the transformation and let l be the pad size
801 // 1) don't undo transformation, don't normalize: return q'
802 // 2) undo transformation, don't normalize: return T^{-1} q'
803 // 3) undo transformation, normalize: return (T^{-1} q') / l
804 // 4) don't undo transformation, normalize: return T((T^{-1} q') / l)
805 if (!undoTransformation && !normalizeToPadSize) {/* value remains unchanged */} // (1)
806 else { // (2), (3), (4)
809 case 0: /* value remains unchanged */ break;
810 case 1: value = value * value; break;
811 case 2: value = (TMath::Exp(value / fgkM) - 1) * fgkM; break;
812 default: Error("CreateFitCalROC", "Wrong fit type."); break;
814 if (normalizeToPadSize) value /= GetPadLength(localXY[0]); // (3)
816 if (!undoTransformation && normalizeToPadSize) { // (4)
819 case 0: /* value remains unchanged */ break;
820 case 1: value = TMath::Sqrt(value); break;
821 case 2: value = fgkM * TMath::Log(1 + value / fgkM); break;
822 default: Error("CreateFitCalROC", "Wrong fit type."); break;
825 lROCfitted->SetValue(irow, ipad, value);
831 AliTPCCalROC* AliTPCcalibTracksGain::CreateCombinedCalROC(const AliTPCCalROC* roc1, const AliTPCCalROC* roc2) {
833 // Combines the medium pad size values of roc1 with the long pad size values of roc2 into a new
834 // AliTPCCalROC. Returns a null pointer if any one of the ROCs is an IROC; issues a warning message
835 // if the sectors of roc1 and roc2 don't match, but still continue and use the sector of roc1 as the
836 // sector of the new ROC.
839 if (!roc1 || !roc2) return 0;
840 if ((Int_t)(roc1->GetSector()) < fgTPCparam->GetNInnerSector()) return 0;
841 if ((Int_t)(roc2->GetSector()) < fgTPCparam->GetNInnerSector()) return 0;
842 if (roc1->GetSector() != roc2->GetSector()) Warning("CreateCombinedCalROC", "Sector number mismatch.");
843 AliTPCCalROC* roc = new AliTPCCalROC(roc1->GetSector());
845 for (UInt_t iRow = 0; iRow < 64; iRow++) {
846 for (UInt_t iPad = 0; iPad < roc->GetNPads(iRow); iPad++)
847 roc->SetValue(iRow, iPad, roc1->GetValue(iRow, iPad));
849 for (UInt_t iRow = 64; iRow < roc->GetNrows(); iRow++) {
850 for (UInt_t iPad = 0; iPad < roc->GetNPads(iRow); iPad++)
851 roc->SetValue(iRow, iPad, roc2->GetValue(iRow, iPad));
856 Bool_t AliTPCcalibTracksGain::GetParameters(UInt_t segment, UInt_t padType, UInt_t fitType, TVectorD &fitParam) {
858 // Puts the fit parameters for the specified segment (IROC & OROC), padType and fitType
859 // into the fitParam TVectorD (which should contain 8 elements).
860 // padType is one of kShortPads, kMediumPads, kLongPads. fitType is one of kSimpleFitter, kSqrtFitter, kLogFitter.
861 // Note: The fitter has to be evaluated first!
863 TLinearFitter * fitter = GetFitter(segment, padType, fitType);
866 fitter->GetParameters(fitParam);
869 Error("AliTPCcalibTracksGain::GetParameters","Fitter%d_%d_%d not available", segment, padType, fitType);
875 void AliTPCcalibTracksGain::GetErrors(UInt_t segment, UInt_t padType, UInt_t fitType, TVectorD &fitError) {
877 // Puts the fit parameter errors for the specified segment (IROC & OROC), padType and fitType
878 // into the fitError TVectorD (which should contain 8 elements).
879 // padType is one of kShortPads, kMediumPads, kLongPads. fitType is one of kSimpleFitter, kSqrtFitter, kLogFitter.
880 // Note: The fitter has to be evaluated first!
883 GetFitter(segment, padType, fitType)->GetErrors(fitError);
884 //fitError *= TMath::Sqrt(GetRedChi2(segment, padType, fitType));
888 void AliTPCcalibTracksGain::GetCovarianceMatrix(UInt_t segment, UInt_t padType, UInt_t fitType, TMatrixD& covMatrix) {
890 // Returns the covariance matrix for the specified segment, padType, fitType.
891 // padType is one of kShortPads, kMediumPads, kLongPads. fitType is one of kSimpleFitter, kSqrtFitter, kLogFitter.
894 GetFitter(segment, padType, fitType)->GetCovarianceMatrix(covMatrix);
897 TLinearFitter* AliTPCcalibTracksGain::GetFitter(UInt_t segment, UInt_t padType, UInt_t fitType) {
899 // Returns the TLinearFitter object for the specified segment, padType, fitType.
900 // padType is one of kShortPads, kMediumPads, kLongPads. fitType is one of kSimpleFitter, kSqrtFitter, kLogFitter.
905 return fSimpleFitter->GetFitter(segment, padType);
907 return fSqrtFitter->GetFitter(segment, padType);
909 return fLogFitter->GetFitter(segment, padType);
911 return fSingleSectorFitter->GetFitter(0, padType);
916 Double_t AliTPCcalibTracksGain::GetPadLength(Double_t lx) {
918 // The function returns 0.75 for an IROC, 1. for an OROC at medium pad size position,
919 // 1.5 for an OROC at long pad size position, -1 if out of bounds.
922 Double_t irocLow = fgTPCparam->GetPadRowRadiiLow(0) - fgTPCparam->GetInnerPadPitchLength()/2;
923 Double_t irocUp = fgTPCparam->GetPadRowRadiiLow(fgTPCparam->GetNRowLow()-1) + fgTPCparam->GetInnerPadPitchLength()/2;
924 Double_t orocLow1 = fgTPCparam->GetPadRowRadiiUp(0) - fgTPCparam->GetOuter1PadPitchLength()/2;
925 Double_t orocUp1 = fgTPCparam->GetPadRowRadiiUp(fgTPCparam->GetNRowUp1()-1) + fgTPCparam->GetOuter1PadPitchLength()/2;
926 Double_t orocLow2 = fgTPCparam->GetPadRowRadiiUp(fgTPCparam->GetNRowUp1()) - fgTPCparam->GetOuter2PadPitchLength()/2;
927 Double_t orocUp2 = fgTPCparam->GetPadRowRadiiUp(fgTPCparam->GetNRowUp()-1) + fgTPCparam->GetOuter2PadPitchLength()/2;
930 if (lx >= irocLow && lx <= irocUp) return 0.75;
931 // if OROC medium pads
932 if (lx >= orocLow1 && lx <= orocUp1) return 1.;
934 if (lx >= orocLow2 && lx <= orocUp2) return 1.5;
939 Int_t AliTPCcalibTracksGain::GetPadType(Double_t lx) {
941 // The function returns 0 for an IROC, 1 for an OROC at medium pad size position,
942 // 2 for an OROC at long pad size position, -1 if out of bounds.
945 if (GetPadLength(lx) == 0.75) return 0;
946 else if (GetPadLength(lx) == 1.) return 1;
947 else if (GetPadLength(lx) == 1.5) return 2;
951 void AliTPCcalibTracksGain::DumpTrack(AliTPCseed* track) {
953 // Dump track information to the debug stream
966 for (Int_t ipad = 0; ipad < 3; ipad++) {
967 dedxM[ipad].ResizeTo(5);
968 dedxQ[ipad].ResizeTo(5);
969 parY[ipad].ResizeTo(3);
970 parZ[ipad].ResizeTo(3);
971 meanPos[ipad].ResizeTo(6);
972 Bool_t isOK = GetDedx(track, ipad, rows, sector[ipad], npoints[ipad], dedxM[ipad], dedxQ[ipad], parY[ipad], parZ[ipad], meanPos[ipad]);
974 AddTracklet(sector[ipad],ipad, dedxQ[ipad], dedxM[ipad], parY[ipad], parZ[ipad], meanPos[ipad] );
978 TTreeSRedirector * cstream = GetDebugStreamer();
980 (*cstream) << "Track" <<
981 "run="<<fRun<< // run number
982 "event="<<fEvent<< // event number
983 "time="<<fTime<< // time stamp of event
984 "trigger="<<fTrigger<< // trigger
985 "mag="<<fMagF<< // magnetic field
986 "Track.=" << track << // track information
991 if ( GetStreamLevel()>1 && count>1){
992 (*cstream) << "TrackG" <<
993 "run="<<fRun<< // run number
994 "event="<<fEvent<< // event number
995 "time="<<fTime<< // time stamp of event
996 "trigger="<<fTrigger<< // trigger
997 "mag="<<fMagF<< // magnetic field
998 "Track.=" << track << // track information
1000 // info for pad type 0
1001 "sector0="<<sector[0]<<
1002 "npoints0="<<npoints[0]<<
1003 "dedxM0.="<<&dedxM[0]<<
1004 "dedxQ0.="<<&dedxQ[0]<<
1005 "parY0.="<<&parY[0]<<
1006 "parZ0.="<<&parZ[0]<<
1007 "meanPos0.="<<&meanPos[0]<<
1009 // info for pad type 1
1010 "sector1="<<sector[1]<<
1011 "npoints1="<<npoints[1]<<
1012 "dedxM1.="<<&dedxM[1]<<
1013 "dedxQ1.="<<&dedxQ[1]<<
1014 "parY1.="<<&parY[1]<<
1015 "parZ1.="<<&parZ[1]<<
1016 "meanPos1.="<<&meanPos[1]<<
1018 // info for pad type 2
1019 "sector2="<<sector[2]<<
1020 "npoints2="<<npoints[2]<<
1021 "dedxM2.="<<&dedxM[2]<<
1022 "dedxQ2.="<<&dedxQ[2]<<
1023 "parY2.="<<&parY[2]<<
1024 "parZ2.="<<&parZ[2]<<
1025 "meanPos2.="<<&meanPos[2]<<
1032 Bool_t AliTPCcalibTracksGain::GetDedx(AliTPCseed* track, Int_t padType, Int_t* /*rows*/,
1033 Int_t §or, Int_t& npoints,
1034 TVectorD &dedxM, TVectorD &dedxQ,
1035 TVectorD &parY, TVectorD &parZ, TVectorD&meanPos)
1038 // GetDedx for given sector for given track
1039 // padType - type of pads
1042 static TLinearFitter fitY(2, "pol1");
1043 static TLinearFitter fitZ(2, "pol1");
1044 fitY.StoreData(kFALSE);
1045 fitZ.StoreData(kFALSE);
1048 Int_t firstRow = 0, lastRow = 0;
1050 Float_t xcenter = 0;
1051 const Float_t ktany = TMath::Tan(TMath::DegToRad() * 10);
1052 const Float_t kedgey = 4.;
1055 lastRow = fgTPCparam->GetNRowLow();
1059 firstRow = fgTPCparam->GetNRowLow();
1060 lastRow = fgTPCparam->GetNRowLow() + fgTPCparam->GetNRowUp1();
1064 firstRow = fgTPCparam->GetNRowLow() + fgTPCparam->GetNRowUp1();
1065 lastRow = fgTPCparam->GetNRowLow() + fgTPCparam->GetNRowUp();
1068 minRow = (lastRow - firstRow) / 2;
1071 Int_t nclusters = 0;
1072 Int_t nclustersNE = 0; // number of not edge clusters
1073 Int_t lastSector = -1;
1074 Float_t amplitudeQ[100];
1075 Float_t amplitudeM[100];
1083 for (Int_t iCluster = firstRow; iCluster < lastRow; iCluster++) {
1084 AliTPCclusterMI* cluster = track->GetClusterPointer(iCluster);
1086 Int_t detector = cluster->GetDetector() ;
1087 if (lastSector == -1) lastSector = detector;
1088 if (lastSector != detector) continue;
1089 amplitudeQ[nclusters] = cluster->GetQ()/GetQNorm(cluster,parY[1], parZ[1]);
1090 amplitudeM[nclusters] = cluster->GetMax()/GetMaxNorm(cluster,parY[1], parZ[1]);
1091 rowIn[nclusters] = iCluster;
1093 Double_t dx = cluster->GetX() - xcenter;
1094 Double_t y = cluster->GetY();
1095 Double_t z = cluster->GetZ();
1096 fitY.AddPoint(&dx, y);
1097 fitZ.AddPoint(&dx, z);
1104 if (TMath::Abs(cluster->GetY()) < cluster->GetX()*ktany - kedgey) nclustersNE++;
1108 if (nclusters < minRow / 2) return kFALSE;
1109 if (nclustersNE < minRow / 2) return kFALSE;
1110 for (Int_t i = 0; i < 6; i++) meanPos[i] /= Double_t(nclusters);
1113 fitY.GetParameters(parY);
1114 fitZ.GetParameters(parZ);
1116 // calculate truncated mean
1118 TMath::Sort(nclusters, amplitudeQ, index, kFALSE);
1123 for (Int_t i = 0; i < 5; i++) {
1132 for (Int_t i = 0; i < nclusters; i++) {
1133 Int_t rowSorted = rowIn[index[i]];
1134 AliTPCclusterMI* cluster = track->GetClusterPointer(rowSorted);
1136 if (TMath::Abs(cluster->GetY()) > cluster->GetX()*ktany - kedgey) continue; //don't take edge clusters
1138 if (inonEdge < nclustersNE * 0.5) {
1140 dedxQ[0] += amplitudeQ[index[i]];
1141 dedxM[0] += amplitudeM[index[i]];
1143 if (inonEdge < nclustersNE * 0.6) {
1145 dedxQ[1] += amplitudeQ[index[i]];
1146 dedxM[1] += amplitudeM[index[i]];
1148 if (inonEdge < nclustersNE * 0.7) {
1150 dedxQ[2] += amplitudeQ[index[i]];
1151 dedxM[2] += amplitudeM[index[i]];
1153 if (inonEdge < nclustersNE * 0.8) {
1155 dedxQ[3] += amplitudeQ[index[i]];
1156 dedxM[3] += amplitudeM[index[i]];
1158 if (inonEdge < nclustersNE * 0.9) {
1160 dedxQ[4] += amplitudeQ[index[i]];
1161 dedxM[4] += amplitudeM[index[i]];
1164 for (Int_t i = 0; i < 5; i++) {
1165 dedxQ[i] /= ndedx[i];
1166 dedxM[i] /= ndedx[i];
1168 TTreeSRedirector * cstream = GetDebugStreamer();
1170 Float_t momenta = track->GetP();
1171 Float_t mdedx = track->GetdEdx();
1172 for (Int_t i = 0; i < nclusters; i++) {
1173 Int_t rowSorted = rowIn[index[i]];
1174 AliTPCclusterMI* cluster = track->GetClusterPointer(rowSorted);
1176 printf("Problem\n");
1179 if (TMath::Abs(cluster->GetY()) < cluster->GetX()*ktany - kedgey) inonEdge++;
1180 Float_t dedge = cluster->GetX()*ktany - TMath::Abs(cluster->GetY());
1181 Float_t fraction = Float_t(i) / Float_t(nclusters);
1182 Float_t fraction2 = Float_t(inonEdge) / Float_t(nclustersNE);
1184 AddCluster(cluster, momenta, mdedx, padType, xcenter, dedxQ, dedxM, fraction, fraction2, dedge, parY, parZ, meanPos);
1186 Double_t qNorm = GetQNorm(cluster,parY[1], parZ[1]);
1187 Double_t mNorm = GetMaxNorm(cluster,parY[1], parZ[1]);
1190 Float_t gain = GetGain(cluster);
1191 if (cstream) (*cstream) << "dEdx" <<
1192 "run="<<fRun<< // run number
1193 "event="<<fEvent<< // event number
1194 "time="<<fTime<< // time stamp of event
1195 "trigger="<<fTrigger<< // trigger
1196 "mag="<<fMagF<< // magnetic field
1198 "Cl.=" << cluster << // cluster of interest
1199 "gain="<<gain<< // gain at cluster position
1200 "mNorm="<<mNorm<< // Q max normalization
1201 "qNorm="<<qNorm<< // Q tot normalization
1202 "P=" << momenta << // track momenta
1203 "dedx=" << mdedx << // mean dedx - corrected for angle
1204 "IPad=" << padType << // pad type 0..2
1205 "xc=" << xcenter << // x center of chamber
1206 "dedxQ.=" << &dedxQ << // dedxQ - total charge
1207 "dedxM.=" << &dedxM << // dedxM - maximal charge
1208 "fraction=" << fraction << // fraction - order in statistic (0,1)
1209 "fraction2=" << fraction2 << // fraction - order in statistic (0,1)
1210 "dedge=" << dedge << // distance to the edge
1211 "parY.=" << &parY << // line fit
1212 "parZ.=" << &parZ << // line fit
1213 "meanPos.=" << &meanPos << // mean position (dx, dx^2, y,y^2, z, z^2)
1217 if (cstream) (*cstream) << "dEdxT" <<
1218 "run="<<fRun<< // run number
1219 "event="<<fEvent<< // event number
1220 "time="<<fTime<< // time stamp of event
1221 "trigger="<<fTrigger<< // trigger
1222 "mag="<<fMagF<< // magnetic field
1223 "P=" << momenta << // track momenta
1224 "npoints="<<inonEdge<< // number of points
1225 "sector="<<lastSector<< // sector number
1226 "dedx=" << mdedx << // mean dedx - corrected for angle
1227 "IPad=" << padType << // pad type 0..2
1228 "xc=" << xcenter << // x center of chamber
1229 "dedxQ.=" << &dedxQ << // dedxQ - total charge
1230 "dedxM.=" << &dedxM << // dedxM - maximal charge
1231 "parY.=" << &parY << // line fit
1232 "parZ.=" << &parZ << // line fit
1233 "meanPos.=" << &meanPos << // mean position (dx, dx^2, y,y^2, z, z^2)
1236 sector = lastSector;
1241 void AliTPCcalibTracksGain::AddTracklet(UInt_t sector, UInt_t padType,TVectorD &dedxQ, TVectorD &dedxM,TVectorD& parY, TVectorD& parZ, TVectorD& meanPos){
1243 // Add measured point - dedx to the fitter
1246 //chain->SetAlias("dr","(250-abs(meanPos.fElements[4]))/250");
1247 //chain->SetAlias("tz","(0+abs(parZ.fElements[1]))");
1248 //chain->SetAlias("ty","(0+abs(parY.fElements[1]))");
1249 //chain->SetAlias("corrg","sqrt((1+ty^2)*(1+tz^2))");
1250 //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);
1254 // z and angular part
1257 xxx[0] = (250.-TMath::Abs(meanPos[4]))/250.;
1258 xxx[1] = TMath::Abs(parY[1]);
1259 xxx[2] = TMath::Abs(parZ[1]);
1260 xxx[3] = xxx[0]*xxx[1];
1261 xxx[4] = xxx[0]*xxx[2];
1262 xxx[5] = xxx[1]*xxx[2];
1263 xxx[6] = xxx[0]*xxx[0];
1264 xxx[7] = xxx[1]*xxx[1];
1265 xxx[8] = xxx[2]*xxx[2];
1269 Int_t tsector = sector%36;
1270 for (Int_t i=0;i<35;i++){
1271 xxx[9+i]=(i==tsector)?1:0;
1273 TLinearFitter *fitterM = fFitter0M;
1274 if (padType==1) fitterM=fFitter1M;
1275 if (padType==2) fitterM=fFitter2M;
1276 fitterM->AddPoint(xxx,dedxM[1]);
1278 TLinearFitter *fitterT = fFitter0T;
1279 if (padType==1) fitterT = fFitter1T;
1280 if (padType==2) fitterT = fFitter2T;
1281 fitterT->AddPoint(xxx,dedxQ[1]);
1283 TLinearFitter *dfitterM = fDFitter0M;
1284 if (padType==1) dfitterM=fDFitter1M;
1285 if (padType==2) dfitterM=fDFitter2M;
1286 dfitterM->AddPoint(xxx,dedxM[1]);
1288 TLinearFitter *dfitterT = fDFitter0T;
1289 if (padType==1) dfitterT = fDFitter1T;
1290 if (padType==2) dfitterT = fDFitter2T;
1291 dfitterT->AddPoint(xxx,dedxQ[1]);
1295 TGraph *AliTPCcalibTracksGain::CreateAmpGraph(Int_t ipad, Bool_t qmax){
1297 // create the amplitude graph
1298 // The normalized amplitudes are extrapolated to the 0 angle (y,z) and 0 drift length
1303 if (ipad==0) fFitter0M->GetParameters(vec);
1304 if (ipad==1) fFitter1M->GetParameters(vec);
1305 if (ipad==2) fFitter2M->GetParameters(vec);
1307 if (ipad==0) fFitter0T->GetParameters(vec);
1308 if (ipad==1) fFitter1T->GetParameters(vec);
1309 if (ipad==2) fFitter2T->GetParameters(vec);
1314 for (Int_t i=0;i<35;i++){
1316 amp[i]=vec[10+i]+vec[0];
1319 Float_t mean = TMath::Mean(36,amp);
1320 for (Int_t i=0;i<36;i++){
1322 amp[i]=(amp[i]-mean)/mean;
1324 TGraph *gr = new TGraph(36,sec,amp);
1329 void AliTPCcalibTracksGain::UpdateClusterParam(AliTPCClusterParam* clparam){
1331 // SetQ normalization parameters
1333 // void SetQnorm(Int_t ipad, Int_t itype, TVectorD * norm);
1344 fDFitter0T->GetParameters(vec);
1345 clparam->SetQnorm(0,0,&vec);
1346 fDFitter1T->GetParameters(vec);
1347 clparam->SetQnorm(1,0,&vec);
1348 fDFitter2T->GetParameters(vec);
1349 clparam->SetQnorm(2,0,&vec);
1351 fDFitter0M->GetParameters(vec);
1352 clparam->SetQnorm(0,1,&vec);
1353 fDFitter1M->GetParameters(vec);
1354 clparam->SetQnorm(1,1,&vec);
1355 fDFitter2M->GetParameters(vec);
1356 clparam->SetQnorm(2,1,&vec);
1362 void AliTPCcalibTracksGain::Analyze(){