TPC/AliTPCcalibTracksGain.cxx

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