TRD/AliTRDCalibTask.cxx

   1
   2 /**************************************************************************
   3  * Copyright(c) 1998-1999, ALICE Experiment at CERN, All rights reserved. *
   4  *                                                                        *
   5  * Author: The ALICE Off-line Project.                                    *
   6  * Contributors are mentioned in the code where appropriate.              *
   7  *                                                                        *
   8  * Permission to use, copy, modify and distribute this software and its   *
   9  * documentation strictly for non-commercial purposes is hereby granted   *
  10  * without fee, provided that the above copyright notice appears in all   *
  11  * copies and that both the copyright notice and this permission notice   *
  12  * appear in the supporting documentation. The authors make no claims     *
  13  * about the suitability of this software for any purpose. It is          *
  14  * provided "as is" without express or implied warranty.                  *
  15  **************************************************************************/
  16
  17 /////////////////////////////////////////////////////////////////////////////////
  18 //
  19 // AliTRDCalibTask
  20 //
  21 // Offline TRD calibration task
  22 //
  23 // Author:
  24 //   R. Bailhache (R.Bailhache@gsi.de)
  25 //
  26 //////////////////////////////////////////////////////////////////////////////////////
  27
  28 #include <iostream>
  29 using namespace std;
  30 #include "Riostream.h"
  31 #include "TChain.h"
  32 #include "TTree.h"
  33 #include "TFile.h"
  34 #include "TProfile2D.h"
  35 #include "TH2I.h"
  36 #include "TH1F.h"
  37 #include "TH2F.h"
  38 #include "TH2S.h"
  39 #include "TList.h"
  40 #include "TMath.h"
  41 #include "TObject.h"
  42 #include "TObjArray.h"
  43 #include "TString.h"
  44 #include "TCanvas.h"
  45 #include "TLegend.h"
  46 #include "TStyle.h"
  47 #include "TLine.h"
  48 #include "TIterator.h"
  49 #include "TLinearFitter.h"
  50 #include "TVectorD.h"
  51
  52 #include "AliAnalysisTask.h"
  53 #include "AliAnalysisManager.h"
  54
  55 #include "AliExternalTrackParam.h"
  56 #include "AliESDVertex.h"
  57 #include "AliESDEvent.h"
  58 #include "AliESDfriend.h"
  59 #include "AliCentrality.h"
  60 #include "AliESDInputHandler.h"
  61 #include "AliESDtrack.h"
  62 #include "AliESDfriendTrack.h"
  63 #include "AliTRDtrackV1.h"
  64 #include "AliTRDseedV1.h"
  65 #include "AliTRDcluster.h"
  66 #include "AliTRDgeometry.h"
  67 #include "AliESDtrackCuts.h"
  68 #include "AliESDVertex.h"
  69 #include "AliTRDCalDet.h"
  70
  71 #include "AliTRDCalibraVector.h"
  72 #include "AliTRDCalibraFillHisto.h"
  73 #include "AliTRDCalibraVdriftLinearFit.h"
  74 #include "AliTRDCalibraExbAltFit.h"
  75
  76 #include "AliTRDcalibDB.h"
  77 #include "AliCDBId.h"
  78 #include "AliLog.h"
  79
  80 #include "AliTRDCalibChamberStatus.h"
  81 #include "AliTRDdEdxBaseUtils.h"
  82 #include "AliTRDdEdxCalibHistArray.h"
  83 #include "AliTRDdEdxCalibUtils.h"
  84
  85 #include "AliTRDCalibTask.h"
  86
  87
  88 ClassImp(AliTRDCalibTask)
  89
  90 //________________________________________________________________________
  91   AliTRDCalibTask::AliTRDCalibTask(const char *name)
  92     : AliAnalysisTaskSE(name), fESD(0),
  93       fkEsdTrack(0),
  94       fFriendTrack(0),
  95       fCalibObject(0),
  96       fTrdTrack(0),
  97       fCl(0),
  98       fListHist(0),
  99       fTRDCalibraFillHisto(0),
 100       fTRDChamberStatus(0),
 101       fNEvents(0),
 102       fNEventsInput(0),
 103       fNbTRDTrack(0),
 104       fNbTRDTrackOffline(0),
 105       fNbTRDTrackStandalone(0),
 106       fNbTPCTRDtrack(0),
 107       fNbGoodTracks(0),
 108       fNbTimeBin(0),
 109       fNbTimeBinOffline(0),
 110       fNbTimeBinStandalone(0),
 111       fNbClusters(0),
 112       fNbClustersOffline(0),
 113       fNbClustersStandalone(0),
 114       fNbTracklets(0),
 115       fNbTrackletsOffline(0),
 116       fNbTrackletsStandalone(0),
 117       fAbsoluteGain(0),
 118       fCH2dSum(0),
 119       fPH2dSum(0),
 120       fCH2dSM(0),
 121       fPH2dSM(0),
 122       fCH2dTest(0),
 123       fPH2dTest(0),
 124       fLinearVdriftTest(0),
 125       fOnInstance(kTRUE),
 126       fHisto2d(kTRUE),
 127       fVector2d(kFALSE),
 128       fVdriftLinear(kTRUE),
 129       fExbAlt(kFALSE),
 130       fDebugLevelTRDCalibraFillHisto(0),
 131       fNbTimeBins(0),
 132       fNumberBinCharge(50),
 133       fSelectTrigger(kTRUE),
 134       fSelectedTrigger(new TObjArray()),
 135       fRejected(kTRUE),
 136       fEsdTrackCuts(0),
 137       fRequirePrimaryVertex(kFALSE),
 138       fVtxTPC(kFALSE),
 139       fVtxSPD(kFALSE),
 140       fMinNbContributors(0),
 141       fRangePrimaryVertexZ(9999999.0),
 142       fMinNbTracks(9),
 143       fMaxNbTracks(999999999),
 144       fCutWithVdriftCalib(kFALSE),
 145       fMinNbTRDtracklets(0),
 146       fMinTRDMomentum(0),
 147       fScaleGainWithTPCSignal(kFALSE),
 148       fLow(0),
 149       fHigh(30),
 150       fFillZero(kFALSE),
 151       fNormalizeNbOfCluster(kFALSE),
 152       fRelativeScale(0.0),
 153       fMaxCluster(100.0),
 154       fNbMaxCluster(2),
 155       fOfflineTracks(kFALSE),
 156       fStandaloneTracks(kFALSE),
 157       fFirstRunGain(-1),
 158       fVersionGainUsed(-1),
 159       fSubVersionGainUsed(-1),
 160       fFirstRunGainLocal(-1),
 161       fVersionGainLocalUsed(-1),
 162       fSubVersionGainLocalUsed(-1),
 163       fFirstRunVdrift(-1),
 164       fVersionVdriftUsed(-1),
 165       fSubVersionVdriftUsed(-1),
 166       fFirstRunExB(-1),
 167       fVersionExBUsed(-1),
 168       fSubVersionExBUsed(-1),
 169       fCalDetGain(0x0),
 170       fMaxEvent(0),
 171       fCounter(0),
 172       fPHQon(kTRUE)
 173 {
 174   //
 175   // Default constructor
 176   //
 177
 178   fNz[0] = 0;
 179   fNz[1] = 0;
 180   fNz[2] = 0;
 181
 182   fNrphi[0] = 0;
 183   fNrphi[1] = 0;
 184   fNrphi[2] = 0;
 185
 186
 187   // Define input and output slots here
 188   // Input slot #0 works with a TChain
 189   DefineInput(0, TChain::Class());
 190
 191   // Output slot #0 writes into a TList container
 192   DefineOutput(1, TList::Class());
 193
 194
 195 }
 196 //____________________________________________________________________________________
 197 AliTRDCalibTask::~AliTRDCalibTask()
 198 {
 199   //
 200   // AliTRDCalibTask destructor
 201   //
 202
 203   // Pointeur
 204   if(fNEvents) delete fNEvents;
 205   if(fNEventsInput) delete fNEventsInput;
 206   if(fNbTRDTrack) delete fNbTRDTrack;
 207   if(fNbTRDTrackOffline) delete fNbTRDTrackOffline;
 208   if(fNbTRDTrackStandalone) delete fNbTRDTrackStandalone;
 209   if(fNbTPCTRDtrack) delete fNbTPCTRDtrack;
 210   if(fNbGoodTracks) delete fNbGoodTracks;
 211   if(fNbTimeBin) delete fNbTimeBin;
 212   if(fNbTimeBinOffline) delete fNbTimeBinOffline;
 213   if(fNbTimeBinStandalone) delete fNbTimeBinStandalone;
 214   if(fNbClusters) delete fNbClusters;
 215   if(fNbClustersOffline) delete fNbClustersOffline;
 216   if(fNbClustersStandalone) delete fNbClustersStandalone;
 217   if(fNbTracklets) delete fNbTracklets;
 218   if(fNbTrackletsOffline) delete fNbTrackletsOffline;
 219   if(fNbTrackletsStandalone) delete fNbTrackletsStandalone;
 220   if(fAbsoluteGain) delete fAbsoluteGain;
 221   if(fCH2dSum) delete fCH2dSum;
 222   if(fPH2dSum) delete fPH2dSum;
 223   if(fCH2dSM) delete fCH2dSM;
 224   if(fPH2dSM) delete fPH2dSM;
 225   if(fCH2dTest) delete fCH2dTest;
 226   if(fPH2dTest) delete fPH2dTest;
 227   if(fLinearVdriftTest) delete fLinearVdriftTest;
 228   if(IsPHQon()){
 229     AliTRDdEdxCalibUtils::DeleteHistArray();
 230   }
 231
 232   if(fCalDetGain) delete fCalDetGain;
 233
 234   if(fSelectedTrigger) {
 235     fSelectedTrigger->Delete();
 236     delete fSelectedTrigger;
 237   }
 238   if(fEsdTrackCuts) {
 239     delete fEsdTrackCuts;
 240   }
 241
 242   if(fTRDChamberStatus) delete fTRDChamberStatus;
 243
 244 }
 245
 246 //________________________________________________________________________
 247 void AliTRDCalibTask::UserCreateOutputObjects()
 248 {
 249   //
 250   // Create the histos
 251   //
 252   //cout << "AliTRDCalibTask::CreateOutputObjects() IN" << endl;
 253
 254   // Number of time bins
 255   if(fNbTimeBins==0) {
 256     AliTRDcalibDB *cal = AliTRDcalibDB::Instance();
 257     fNbTimeBins = cal->GetNumberOfTimeBinsDCS();
 258     if(fNbTimeBins <= 0){
 259       AliWarning(Form("No of TimeBins from DB [%d] use default [30]", fNbTimeBins));
 260       fNbTimeBins = 30;
 261     }
 262   }
 263
 264   // output list
 265   fListHist = new TList();
 266   fListHist->SetOwner();
 267
 268   // init chamber status
 269   fTRDChamberStatus = new AliTRDCalibChamberStatus();
 270   fTRDChamberStatus->Init();
 271   fListHist->Add(fTRDChamberStatus->GetSparseI());
 272
 273   // instance calibration
 274   fTRDCalibraFillHisto = AliTRDCalibraFillHisto::Instance();
 275   if(fOnInstance) {
 276     fTRDCalibraFillHisto->SetNumberBinCharge(fNumberBinCharge); // set number of bin of the charge distribution
 277     fTRDCalibraFillHisto->SetCutWithVdriftCalib(fCutWithVdriftCalib); // cut vdrift
 278     fTRDCalibraFillHisto->SetMinNbTRDtracklets(fMinNbTRDtracklets); // min number of TRD tracklets
 279     fTRDCalibraFillHisto->SetMinTRDMomentum(fMinTRDMomentum); // min TRD momentum
 280     fTRDCalibraFillHisto->SetHisto2d(fHisto2d); // choose to use histograms
 281     fTRDCalibraFillHisto->SetVector2d(fVector2d); // choose to use vectors
 282     fTRDCalibraFillHisto->SetCH2dOn();  // choose to calibrate the gain
 283     fTRDCalibraFillHisto->SetPH2dOn();  // choose to calibrate the drift velocity
 284     fTRDCalibraFillHisto->SetPRF2dOn(); // choose to look at the PRF
 285     fTRDCalibraFillHisto->SetLinearFitterOn(fVdriftLinear); // Other possibility vdrift VDRIFT
 286     fTRDCalibraFillHisto->SetLinearFitterDebugOn(fVdriftLinear); // Other possibility vdrift
 287     fTRDCalibraFillHisto->SetExbAltFitOn(fExbAlt); // Alternative method for exb
 288     fTRDCalibraFillHisto->SetScaleWithTPCSignal(fScaleGainWithTPCSignal); // Scale Gain with TPC signal
 289     for(Int_t k = 0; k < 3; k++){
 290       if(((fNz[k] != 10) && (fNrphi[k] != 10)) && ((fNz[k] != 100) && (fNrphi[k] != 100))) {
 291         fTRDCalibraFillHisto->SetNz(k,fNz[k]);                                    // Mode calibration
 292         fTRDCalibraFillHisto->SetNrphi(k,fNrphi[k]);                             // Mode calibration
 293       }
 294       else {
 295         if((fNz[k] == 100) && (fNrphi[k] == 100))  {
 296           if(fVector2d) AliInfo("The mode all together is not supported by the vector method");
 297           fTRDCalibraFillHisto->SetAllTogether(k);
 298         }
 299         if((fNz[k] == 10) && (fNrphi[k] == 10))  {
 300           if(fVector2d) AliInfo("The mode per supermodule is not supported by the vector method");
 301           fTRDCalibraFillHisto->SetPerSuperModule(k);
 302         }
 303       }
 304     }
 305     // Variables for how to fill
 306     fTRDCalibraFillHisto->SetFillWithZero(fFillZero);
 307     fTRDCalibraFillHisto->SetNormalizeNbOfCluster(fNormalizeNbOfCluster);
 308     fTRDCalibraFillHisto->SetMaxCluster(fMaxCluster);
 309     fTRDCalibraFillHisto->SetNbMaxCluster(fNbMaxCluster);
 310
 311     // Init with 30 timebins
 312     fTRDCalibraFillHisto->Init2Dhistos(fNbTimeBins); // initialise the histos
 313     fTRDCalibraFillHisto->SetNumberClusters(fLow); // At least 11 clusters
 314     fTRDCalibraFillHisto->SetNumberClustersf(fHigh); // At least 11 clusters
 315
 316     // For testing only
 317     fTRDCalibraFillHisto->SetDebugLevel(fDebugLevelTRDCalibraFillHisto); //debug stuff
 318
 319     if(fHisto2d) {
 320       fListHist->Add(fTRDCalibraFillHisto->GetCH2d());
 321       fListHist->Add(fTRDCalibraFillHisto->GetPH2d());
 322       fListHist->Add(fTRDCalibraFillHisto->GetPRF2d());
 323     }
 324     if(fVdriftLinear) fListHist->Add((TObject *)fTRDCalibraFillHisto->GetVdriftLinearFit());
 325     if(fVector2d) fListHist->Add((TObject *) fTRDCalibraFillHisto->GetCalibraVector()); //calibra vector
 326     if(fExbAlt) fListHist->Add((TObject *)fTRDCalibraFillHisto->GetExbAltFit());
 327   }
 328   fRelativeScale = fTRDCalibraFillHisto->GetRelativeScale(); // Get the relative scale for the gain
 329
 330   fNEvents = new TH1I(Form("NEvents_%s",(const char*)fName),"NEvents", 2, 0, 2);
 331   fListHist->Add(fNEvents);
 332   fNEventsInput = new TH1I(Form("NEventsInput_%s",(const char*)fName),"NEventsInput", 2, 0, 2);
 333   fListHist->Add(fNEventsInput);
 334
 335   // absolute gain calibration even without AliESDfriend
 336   Int_t nBinsPt = 25;
 337   Double_t minPt = 0.001;
 338   Double_t maxPt = 10.0;
 339
 340   Double_t *binLimLogPt = new Double_t[nBinsPt+1];
 341   Double_t *binLimPt    = new Double_t[nBinsPt+1];
 342   for(Int_t i=0; i<=nBinsPt; i++) binLimLogPt[i]=(Double_t)TMath::Log10(minPt) + (TMath::Log10(maxPt)-TMath::Log10(minPt))/nBinsPt*(Double_t)i ;
 343   for(Int_t i=0; i<=nBinsPt; i++) binLimPt[i]=(Double_t)TMath::Power(10,binLimLogPt[i]);
 344
 345   fAbsoluteGain = new TH2F(Form("AbsoluteGain_%s",(const char*)fName),"AbsoluteGain", 200, 0.0, 700.0, nBinsPt, binLimPt);
 346   fAbsoluteGain->SetYTitle("Momentum at TRD");
 347   fAbsoluteGain->SetXTitle("charge deposit [a.u]");
 348   fAbsoluteGain->SetZTitle("counts");
 349   fAbsoluteGain->SetStats(0);
 350   fAbsoluteGain->Sumw2();
 351   fListHist->Add(fAbsoluteGain);
 352
 353   if(IsPHQon()){
 354     printf("\n        AliTRDCalibTask PHQ is on!!     \n\n");
 355     AliTRDdEdxBaseUtils::PrintControl();
 356     AliTRDdEdxCalibUtils::IniHistArray(fListHist, kTRUE);
 357   }
 358   else{
 359     printf("\n        AliTRDCalibTask PHQ is off!!     \n\n");
 360   }
 361
 362   /////////////////////////////////////////
 363   // First debug level
 364   ///////////////////////////////////////
 365   if(fDebug > 0) {
 366
 367     fLinearVdriftTest = new TH2S(Form("LFDV0testversion_%s",(const char*)fName),"LFDV0testversion",36,-0.9,0.9,48,-1.2,1.2);
 368     fLinearVdriftTest->SetXTitle("tan(phi_{track})");
 369     fLinearVdriftTest->SetYTitle("dy/dt");
 370     fLinearVdriftTest->SetZTitle("Number of tracklets");
 371     fLinearVdriftTest->SetStats(0);
 372     fLinearVdriftTest->SetDirectory(0);
 373
 374     // Standart with AliESDfriend
 375     fPH2dTest = new TProfile2D(Form("PH2dTest_%s",(const char*)fName),"Nz0Nrphi0"
 376                             ,fNbTimeBins,-0.05,(Double_t)((fNbTimeBins-0.5)/10.0)
 377                            ,540,0,540);
 378     fPH2dTest->SetYTitle("Det/pad groups");
 379     fPH2dTest->SetXTitle("time [#mus]");
 380     fPH2dTest->SetZTitle("<PH> [a.u.]");
 381     fPH2dTest->SetStats(0);
 382     //
 383     fCH2dTest = new TH2I(Form("CH2dTest_%s",(const char*)fName),"Nz0Nrphi0",50,0,300,540,0,540);
 384     fCH2dTest->SetYTitle("Det/pad groups");
 385     fCH2dTest->SetXTitle("charge deposit [a.u]");
 386     fCH2dTest->SetZTitle("counts");
 387     fCH2dTest->SetStats(0);
 388     fCH2dTest->Sumw2();
 389
 390     //
 391     fPH2dSM = new TProfile2D(Form("PH2dSM_%s",(const char*)fName),"Nz10Nrphi10"
 392                             ,fNbTimeBins,-0.05,(Double_t)((fNbTimeBins-0.5)/10.0)
 393                            ,18,0,18);
 394     fPH2dSM->SetYTitle("Det/pad groups");
 395     fPH2dSM->SetXTitle("time [#mus]");
 396     fPH2dSM->SetZTitle("<PH> [a.u.]");
 397     fPH2dSM->SetStats(0);
 398     //
 399     fCH2dSM = new TH2I(Form("CH2dSM_%s",(const char*)fName),"Nz10Nrphi10",50,0,300,18,0,18);
 400     fCH2dSM->SetYTitle("Det/pad groups");
 401     fCH2dSM->SetXTitle("charge deposit [a.u]");
 402     fCH2dSM->SetZTitle("counts");
 403     fCH2dSM->SetStats(0);
 404     fCH2dSM->Sumw2();
 405     //
 406     fPH2dSum = new TProfile2D(Form("PH2dSum_%s",(const char*)fName),"Nz100Nrphi100"
 407                             ,fNbTimeBins,-0.05,(Double_t)((fNbTimeBins-0.5)/10.0)
 408                             ,1,0,1);
 409     fPH2dSum->SetYTitle("Det/pad groups");
 410     fPH2dSum->SetXTitle("time [#mus]");
 411     fPH2dSum->SetZTitle("<PH> [a.u.]");
 412     fPH2dSum->SetStats(0);
 413     //
 414     fCH2dSum = new TH2I(Form("CH2dSum_%s",(const char*)fName),"Nz100Nrphi100",50,0,300,1,0,1);
 415     fCH2dSum->SetYTitle("Det/pad groups");
 416     fCH2dSum->SetXTitle("charge deposit [a.u]");
 417     fCH2dSum->SetZTitle("counts");
 418     fCH2dSum->SetStats(0);
 419     fCH2dSum->Sumw2();
 420
 421
 422     // Add them
 423     fListHist->Add(fLinearVdriftTest);
 424     fListHist->Add(fPH2dTest);
 425     fListHist->Add(fCH2dTest);
 426     fListHist->Add(fPH2dSM);
 427     fListHist->Add(fCH2dSM);
 428     fListHist->Add(fPH2dSum);
 429     fListHist->Add(fCH2dSum);
 430
 431   }
 432
 433   /////////////////////////////////////////
 434   // Second debug level
 435   ///////////////////////////////////////
 436   if(fDebug > 1) {
 437
 438     fNbGoodTracks = new TH2F(Form("NbGoodTracks_%s",(const char*)fName),"NbGoodTracks",500,0.0,2500.0,200,0.0,100.0);
 439     fNbGoodTracks->SetXTitle("Nb of good tracks");
 440     fNbGoodTracks->SetYTitle("Centrality");
 441     fNbGoodTracks->SetStats(0);
 442
 443     fNbTRDTrack = new TH1F(Form("TRDTrack_%s",(const char*)fName),"TRDTrack",50,0,50);
 444     fNbTRDTrack->Sumw2();
 445     fNbTRDTrackOffline = new TH1F(Form("TRDTrackOffline_%s",(const char*)fName),"TRDTrackOffline",50,0,50);
 446     fNbTRDTrackOffline->Sumw2();
 447     fNbTRDTrackStandalone = new TH1F(Form("TRDTrackStandalone_%s",(const char*)fName),"TRDTrackStandalone",50,0,50);
 448     fNbTRDTrackStandalone->Sumw2();
 449     fNbTPCTRDtrack = new TH2F(Form("NbTPCTRDtrack_%s",(const char*)fName),"NbTPCTRDtrack",100,0,100,100,0,100);
 450     fNbTPCTRDtrack->Sumw2();
 451     //
 452     fNbTimeBin = new TH1F(Form("NbTimeBin_%s",(const char*)fName),"NbTimeBin",35,0,35);
 453     fNbTimeBin->Sumw2();
 454     fNbTimeBinOffline = new TH1F(Form("NbTimeBinOffline_%s",(const char*)fName),"NbTimeBinOffline",35,0,35);
 455     fNbTimeBinOffline->Sumw2();
 456     fNbTimeBinStandalone = new TH1F(Form("NbTimeBinStandalone_%s",(const char*)fName),"NbTimeBinStandalone",35,0,35);
 457     fNbTimeBinStandalone->Sumw2();
 458     //
 459     fNbClusters = new TH1F(Form("NbClusters_%s",(const char*)fName),"",35,0,35);
 460     fNbClusters->Sumw2();
 461     fNbClustersOffline = new TH1F(Form("NbClustersOffline_%s",(const char*)fName),"",35,0,35);
 462     fNbClustersOffline->Sumw2();
 463     fNbClustersStandalone = new TH1F(Form("NbClustersStandalone_%s",(const char*)fName),"",35,0,35);
 464     fNbClustersStandalone->Sumw2();
 465     //
 466     fNbTracklets = new TH1F(Form("NbTracklets_%s",(const char*)fName),"NbTracklets",540,0.,540.);
 467     fNbTracklets->Sumw2();
 468     fNbTrackletsOffline = new TH1F(Form("NbTrackletsOffline_%s",(const char*)fName),"NbTrackletsOffline",540,0.,540.);
 469     fNbTrackletsOffline->Sumw2();
 470     fNbTrackletsStandalone = new TH1F(Form("NbTrackletsStandalone_%s",(const char*)fName),"NbTrackletsStandalone",540,0.,540.);
 471     fNbTrackletsStandalone->Sumw2();
 472
 473     fListHist->Add(fNbGoodTracks);
 474
 475     fListHist->Add(fNbTRDTrack);
 476     fListHist->Add(fNbTRDTrackOffline);
 477     fListHist->Add(fNbTRDTrackStandalone);
 478     fListHist->Add(fNbTPCTRDtrack);
 479
 480     fListHist->Add(fNbTimeBin);
 481     fListHist->Add(fNbTimeBinOffline);
 482     fListHist->Add(fNbTimeBinStandalone);
 483     fListHist->Add(fNbClusters);
 484     fListHist->Add(fNbClustersOffline);
 485     fListHist->Add(fNbClustersStandalone);
 486     fListHist->Add(fNbTracklets);
 487     fListHist->Add(fNbTrackletsOffline);
 488     fListHist->Add(fNbTrackletsStandalone);
 489
 490   }
 491
 492   delete [] binLimLogPt;
 493   delete [] binLimPt;
 494
 495   PostData(1,fListHist);
 496
 497   //cout << "AliTRDCalibTask::UserCreateOutputObjects() OUT" << endl;
 498
 499 }
 500
 501 //________________________________________________________________________
 502 void AliTRDCalibTask::UserExec(Option_t *)
 503 {
 504   //
 505   // Filling of the histos
 506   //
 507   //cout << "AliTRDCalibTask::Exec() IN" << endl;
 508
 509   // Init Versions and subversions used
 510   if((fFirstRunGain==-1) || (fVersionGainUsed==-1) || (fSubVersionGainUsed==-1) || (fFirstRunGainLocal==-1) || (fVersionGainLocalUsed==-1) || (fSubVersionGainLocalUsed==-1) || (fFirstRunVdrift==-1) || (fVersionVdriftUsed==-1) || (fSubVersionVdriftUsed==-1)) {
 511     if(!SetVersionSubversion()) {
 512       PostData(1, fListHist);
 513       return;
 514     }
 515   }
 516   if(fCounter==0) {
 517     if(fOnInstance) {
 518       fTRDCalibraFillHisto->SetFirstRunGain(fFirstRunGain); // Gain Used
 519       fTRDCalibraFillHisto->SetVersionGainUsed(fVersionGainUsed); // Gain Used
 520       fTRDCalibraFillHisto->SetSubVersionGainUsed(fSubVersionGainUsed); // Gain Used
 521       fTRDCalibraFillHisto->SetFirstRunGainLocal(fFirstRunGainLocal); // Gain Used
 522       fTRDCalibraFillHisto->SetVersionGainLocalUsed(fVersionGainLocalUsed); // Gain Used
 523       fTRDCalibraFillHisto->SetSubVersionGainLocalUsed(fSubVersionGainLocalUsed); // Gain Used
 524       fTRDCalibraFillHisto->SetFirstRunVdrift(fFirstRunVdrift); // Vdrift Used
 525       fTRDCalibraFillHisto->SetVersionVdriftUsed(fVersionVdriftUsed); // Vdrift Used
 526       fTRDCalibraFillHisto->SetSubVersionVdriftUsed(fSubVersionVdriftUsed); // Vdrift Used
 527       if((fFirstRunExB != -1) && (fVersionExBUsed != -1) && (fSubVersionExBUsed != -1)){
 528         fTRDCalibraFillHisto->SetFirstRunExB(fFirstRunExB); // ExB Used
 529         fTRDCalibraFillHisto->SetVersionExBUsed(fVersionExBUsed); // ExB Used
 530         fTRDCalibraFillHisto->SetSubVersionExBUsed(fSubVersionExBUsed); // ExB Used
 531       }
 532       fTRDCalibraFillHisto->InitCalDet();
 533     }
 534     if(fDebug > 1){
 535       // title CH2dTest
 536       TString name("Ver");
 537       name += fVersionGainUsed;
 538       name += "Subver";
 539       name += fSubVersionGainUsed;
 540       name += "FirstRun";
 541       name += fFirstRunGain;
 542       name += "Nz0Nrphi0";
 543       fCH2dTest->SetTitle(name);
 544       // title PH2dTest
 545       TString namee("Ver");
 546       namee += fVersionVdriftUsed;
 547       namee += "Subver";
 548       namee += fSubVersionVdriftUsed;
 549       namee += "FirstRun";
 550       namee += fFirstRunVdrift;
 551       namee += "Nz0Nrphi0";
 552       fPH2dTest->SetTitle(namee);
 553     }
 554   }
 555
 556   //  AliLog::SetGlobalLogLevel(AliLog::kError);
 557   //  cout << "AliTRDCalibTask::Exec() 1" << endl;
 558   fESD = dynamic_cast<AliESDEvent*>(fInputEvent);
 559   if(!fESD){
 560     AliError("ESD Event missing");
 561     PostData(1, fListHist);
 562     return;
 563   }
 564
 565   const char* type = fESD->GetBeamType();
 566
 567
 568   //printf("Counter %d\n",fCounter);
 569
 570   fCounter++;
 571   fNEventsInput->Fill(1);
 572
 573   //cout << "maxEvent = " << fMaxEvent << endl;
 574   //if(fCounter%100==0) cout << "fCounter = " << fCounter << endl;
 575   if((fMaxEvent != 0) && (fMaxEvent < fCounter)) {
 576     PostData(1, fListHist);
 577     return;
 578   }
 579   //if(fCounter%100==0) cout << "fCounter1 = " << fCounter << endl;
 580   //cout << "event = " << fCounter << endl;
 581
 582   //printf("Counter %d\n",fCounter);
 583
 584   ///////////////////
 585   // Check trigger
 586   ///////////////////
 587   Bool_t pass = kTRUE;
 588
 589   if (fSelectTrigger) {
 590
 591     //printf("Will check the triggers\n");
 592
 593     Int_t numberOfTriggerSelected = fSelectedTrigger->GetEntriesFast();
 594     //printf("numberofTriggerSelected %d\n",numberOfTriggerSelected);
 595     if(fRejected) {
 596       pass = kTRUE;
 597       for(Int_t k = 0; k < numberOfTriggerSelected; k++){
 598         const TObjString *const obString=(TObjString*)fSelectedTrigger->At(k);
 599         const TString tString=obString->GetString();
 600         if(fESD->IsTriggerClassFired((const char*)tString)) {
 601           pass = kFALSE;
 602         }
 603       }
 604     }
 605     else {
 606       pass = kFALSE;
 607       for(Int_t k = 0; k < numberOfTriggerSelected; k++){
 608         const TObjString *const obString=(TObjString*)fSelectedTrigger->At(k);
 609         const TString tString=obString->GetString();
 610         if(fESD->IsTriggerClassFired((const char*)tString)) {
 611           pass = kTRUE;
 612         }
 613       }
 614     }
 615     if(!pass) {
 616       PostData(1, fListHist);
 617       return;
 618     }
 619
 620   }
 621
 622   //printf("Class Fired %s\n",(const char*)fESD->GetFiredTriggerClasses());
 623   //printf("Trigger passed\n");
 624
 625   ///////////////////////////////
 626   // Require a primary vertex
 627   //////////////////////////////
 628   if(fRequirePrimaryVertex) {
 629     const AliESDVertex* vtxESD = 0x0;
 630     if      (fVtxTPC) vtxESD = fESD->GetPrimaryVertexTPC() ;
 631     else if (fVtxSPD) vtxESD = fESD->GetPrimaryVertexSPD() ;
 632     else              vtxESD = fESD->GetPrimaryVertexTracks() ;
 633     if(!vtxESD){
 634       PostData(1, fListHist);
 635       return;
 636     }
 637     Int_t nCtrb = vtxESD->GetNContributors();
 638     if(nCtrb < fMinNbContributors) {
 639       PostData(1, fListHist);
 640       return;
 641     }
 642     Double_t zPosition = vtxESD->GetZ();
 643     if(TMath::Abs(zPosition) > fRangePrimaryVertexZ) {
 644       PostData(1, fListHist);
 645       return;
 646     }
 647
 648   }
 649
 650   //printf("Primary vertex passed\n");
 651
 652   //////////////////////////////////////
 653   // Requirement on number of good tracks
 654   //////////////////////////////////////
 655   Int_t nGoodParticles = 0;
 656   Double_t nbTracks = fESD->GetNumberOfTracks();
 657   for(Int_t itrack = 0; itrack < nbTracks; itrack++) {
 658     if(ParticleGood(itrack)) nGoodParticles++;
 659   }
 660   if(fDebug > 1)  {
 661     // Centrality
 662     AliCentrality *esdCentrality = fESD->GetCentrality();
 663     Float_t centrality = esdCentrality->GetCentralityPercentile("V0M");
 664     //Float_t centralityb = esdCentrality->GetCentralityPercentile("CL1");
 665     fNbGoodTracks->Fill(nGoodParticles,centrality);
 666     //printf("centrality %f, centralityb %f\n",centrality,centralityb);
 667   }
 668
 669   //printf("Beam type %s\n",(const char*)type);
 670   if (!strstr(type,"p-p")) {
 671     //if (strstr(type,"A-A")) {
 672     //printf("Will check the number of good tracks %f %f %f\n",nGoodParticles,fMinNbTracks,fMaxNbTracks);
 673     if((nGoodParticles < fMinNbTracks) || (nGoodParticles > fMaxNbTracks)) {
 674       PostData(1, fListHist);
 675       return;
 676     }
 677   }
 678
 679   fNEvents->Fill(1);
 680
 681   // In total
 682   Int_t nbTrdTracks = 0;
 683   // standalone
 684   Int_t nbTrdTracksStandalone = 0;
 685   // offline
 686   Int_t nbTrdTracksOffline = 0;
 687   // TPC
 688   Int_t nbtrackTPC = 0;
 689
 690
 691
 692   if (nbTracks <= 0.0) {
 693
 694     if(fDebug > 1) {
 695       fNbTRDTrack->Fill(nbTrdTracks);
 696       fNbTRDTrackStandalone->Fill(nbTrdTracksStandalone);
 697       fNbTRDTrackOffline->Fill(nbTrdTracksOffline);
 698     }
 699     PostData(1, fListHist);
 700     return;
 701   }
 702
 703
 704   fESDfriend = dynamic_cast<AliESDfriend*> (fESD->FindListObject("AliESDfriend"));
 705   if(!fESDfriend){
 706     AliError("fESDfriend not available");
 707     PostData(1, fListHist);
 708     return;
 709   }
 710
 711   if(fESDfriend->TestSkipBit()) {
 712     PostData(1, fListHist);
 713     return;
 714   }
 715
 716   //printf("has friends\n");
 717
 718   /////////////////////////////////////
 719   // Loop on AliESDtrack
 720   ////////////////////////////////////
 721   //printf("Nb of tracks %f\n",nbTracks);
 722   for(int itrk=0; itrk < nbTracks; ++itrk){
 723
 724     // Get ESD track
 725     fkEsdTrack = fESD->GetTrack(itrk);
 726     if(!fkEsdTrack) continue;
 727     ULong_t status = fkEsdTrack->GetStatus();
 728     if(status&(AliESDtrack::kTPCout)) ++nbtrackTPC;
 729
 730     fFriendTrack = fESDfriend->GetTrack(itrk);
 731     if(!fFriendTrack)  {
 732       //printf("No friend track %d\n",itrk);
 733       continue;
 734     }
 735
 736     // Other cuts
 737     fTrdTrack = 0x0;
 738     Bool_t good = kTRUE;
 739     Bool_t standalonetrack = kFALSE;
 740     Bool_t offlinetrack = kFALSE;
 741     //ULong_t status = fkEsdTrack->GetStatus();
 742
 743     //////////////////////////////////////
 744     // Loop on calibration objects
 745     //////////////////////////////////////
 746     Int_t icalib=0;
 747     Int_t nTRDtrackV1=0;
 748     while((fCalibObject = (TObject *)(fFriendTrack->GetCalibObject(icalib++)))){
 749       //printf("Name %s\n",fCalibObject->IsA()->GetName());
 750       if(strcmp(fCalibObject->IsA()->GetName(), "AliTRDtrackV1") != 0) continue;
 751       //printf("Find the calibration object\n");
 752       ++nTRDtrackV1;
 753
 754       if((status&(AliESDtrack::kTRDout)) && (!(status&(AliESDtrack::kTRDin)))) {
 755         standalonetrack = kTRUE;
 756       }
 757       if((status&(AliESDtrack::kTRDin))) {
 758         offlinetrack = kTRUE;
 759       }
 760       if(fOfflineTracks){
 761         if(!offlinetrack){
 762           good = kFALSE;
 763         }
 764       }
 765       else if(fStandaloneTracks){
 766         if(!standalonetrack){
 767           good = kFALSE;
 768         }
 769       }
 770
 771       fTrdTrack = (AliTRDtrackV1 *)fCalibObject;
 772       // process chamberstatus
 773       fTRDChamberStatus->ProcessTrack(fTrdTrack);
 774     }
 775
 776     // Quality cuts on the AliESDtrack
 777     if((fEsdTrackCuts) && (!fEsdTrackCuts->IsSelected((AliVParticle *)fkEsdTrack))) {
 778       //printf("Not a good track\n");
 779       continue;
 780     }
 781
 782     // First Absolute gain calibration
 783     Int_t trdNTracklets = (Int_t) fkEsdTrack->GetTRDntracklets();
 784     Int_t trdNTrackletsPID = (Int_t) fkEsdTrack->GetTRDntrackletsPID();
 785     //printf("Number of trd tracklets %d and PID trd tracklets %d\n",trdNTracklets,trdNTrackletsPID);
 786     if((trdNTracklets > 0) && (trdNTrackletsPID > 0)) {
 787       for(Int_t iPlane = 0; iPlane < 6; ++iPlane){
 788         //Double_t slide = fkEsdTrack->GetTRDslice(iPlane);
 789         //printf("Number of slide %d\n",fkEsdTrack->GetNumberOfTRDslices());
 790         //Double_t momentum = fkEsdTrack->GetTRDmomentum(iPlane);
 791         //printf("momentum %f, slide %f\n",momentum,slide);
 792         if(fkEsdTrack->GetTRDslice(iPlane) > 0.0)
 793           fAbsoluteGain->Fill(fkEsdTrack->GetTRDslice(iPlane)*8.0/100.0,
 794                               fkEsdTrack->GetTRDmomentum(iPlane));
 795       }
 796     }
 797
 798
 799     if(!fTrdTrack) continue;
 800
 801     if(good && fOnInstance) {
 802       //cout << "good" << endl;
 803       fTRDCalibraFillHisto->UpdateHistogramsV1(fTrdTrack,fkEsdTrack);
 804       //printf("Fill fTRDCalibraFillHisto\n");
 805     }
 806
 807     if(IsPHQon()){
 808       const Double_t mag = AliTRDdEdxBaseUtils::IsExBOn() ? fESD->GetMagneticField() : -1;
 809       const Int_t charge = AliTRDdEdxBaseUtils::IsExBOn() ? fkEsdTrack->Charge() : -1;
 810       const Double_t toTPCscale = AliTRDdEdxCalibUtils::GetCalibTPCscale(fkEsdTrack->GetTPCncls(), fkEsdTrack->GetTPCsignal());
 811       if(toTPCscale>0){
 812         AliTRDdEdxCalibUtils::FillHist(fTrdTrack, 0, mag, charge, toTPCscale);
 813       }
 814     }
 815
 816     //////////////////////////////////
 817     // Debug
 818     ////////////////////////////////
 819
 820     if(fDebug > 0) {
 821
 822       //printf("Enter debug\n");
 823
 824       Int_t nbtracklets = 0;
 825
 826       // Check some stuff
 827       Bool_t standalonetracklet = kFALSE;
 828       const AliTRDseedV1 *tracklet = 0x0;
 829       //////////////////////////////////////
 830       // Loop tracklets
 831       /////////////////////////////////////
 832       Int_t nbclusters=0;
 833       Double_t phtb[AliTRDseedV1::kNtb];
 834       memset(phtb, 0, AliTRDseedV1::kNtb*sizeof(Double_t));
 835       Double_t sum = 0.0;
 836       Float_t normalisation = 1.13;
 837       Int_t detector = 0;
 838       Int_t sector = 0;
 839       for(Int_t itr = 0; itr < 6; ++itr){
 840
 841         if(!(tracklet = fTrdTrack->GetTracklet(itr))) continue;
 842         if(!tracklet->IsOK()) continue;
 843         ++nbtracklets;
 844         standalonetracklet = kFALSE;
 845         if(tracklet->IsStandAlone()) standalonetracklet = kTRUE;
 846
 847         nbclusters = 0;
 848         memset(phtb, 0, AliTRDseedV1::kNtb*sizeof(Double_t));
 849         sum = 0.0;
 850         detector = 0;
 851         sector = 0;
 852         //Int_t crossrow = 0;
 853
 854           // Check no shared clusters
 855           //for(int icc=AliTRDseedV1::kNtb; icc<AliTRDseedV1::kNclusters; icc++){
 856           //  if((fcl = tracklet->GetClusters(icc)))  crossrow = 1;
 857           // }
 858
 859           // Loop on clusters
 860         Int_t time = 0;
 861         Float_t ch = 0;
 862         Float_t qcl = 0;
 863         for(int ic=0; ic<AliTRDseedV1::kNtb; ++ic){
 864
 865           if(!(fCl = tracklet->GetClusters(ic))) continue;
 866           ++nbclusters;
 867           time = fCl->GetPadTime();
 868           //ch =  tracklet->GetdQdl(ic);
 869           ch =  tracklet->GetQperTB(ic);
 870           qcl = TMath::Abs(fCl->GetQ());
 871           detector = fCl->GetDetector();
 872           // Add the charge if shared cluster
 873           if((ic+AliTRDseedV1::kNtb) < AliTRDseedV1::kNclusters) {
 874             if((fCl = tracklet->GetClusters(ic+AliTRDseedV1::kNtb))) {
 875               qcl += TMath::Abs(fCl->GetQ());
 876               //printf("Add the cluster charge\n");
 877             }
 878           }
 879           if((time>-1) && (time<fNbTimeBins)) phtb[time]=qcl;
 880           if((fCalDetGain) && (fCalDetGain->GetValue(detector) > 0.0)) sum += ch*fCalDetGain->GetValue(detector)/normalisation;
 881           else sum += ch/normalisation;
 882
 883           if(fDebug > 1) {
 884             fNbTimeBin->Fill(time);
 885             if(tracklet->IsStandAlone()) fNbTimeBinStandalone->Fill(time);
 886             else fNbTimeBinOffline->Fill(time);
 887           }
 888         }
 889         sector = AliTRDgeometry::GetSector(detector);
 890
 891         if(fDebug > 1) {
 892           fNbTracklets->Fill(detector);
 893           if(tracklet->IsStandAlone()) fNbTrackletsStandalone->Fill(detector);
 894           else fNbTrackletsOffline->Fill(detector);
 895
 896           fNbClusters->Fill(nbclusters);
 897           if(tracklet->IsStandAlone())  fNbClustersStandalone->Fill(nbclusters);
 898           else  fNbClustersOffline->Fill(nbclusters);
 899         }
 900
 901         if((nbclusters > fLow) && (nbclusters < fHigh)){
 902           if(fRelativeScale > 0.0) sum = sum/fRelativeScale;
 903           fCH2dTest->Fill(sum,detector+0.5);
 904           fCH2dSM->Fill(sum,sector+0.5);
 905           fCH2dSum->Fill(sum,0.5);
 906           Bool_t checknoise = kTRUE;
 907           if(fMaxCluster > 0) {
 908             if(phtb[0] > fMaxCluster) checknoise = kFALSE;
 909             if(fNbTimeBins > fNbMaxCluster) {
 910               for(Int_t k = (fNbTimeBins-fNbMaxCluster); k < fNbTimeBins; k++){
 911                 if(phtb[k] > fMaxCluster) checknoise = kFALSE;
 912               }
 913             }
 914           }
 915           if(checknoise) {
 916             for(int ic=0; ic<fNbTimeBins; ic++){
 917               if(fFillZero) {
 918                 fPH2dTest->Fill((Double_t)(ic/10.0),detector+0.5,(Double_t)phtb[ic]);
 919                 fPH2dSum->Fill((Double_t)(ic/10.0),0.5,(Double_t)phtb[ic]);
 920                 fPH2dSM->Fill((Double_t)(ic/10.0),sector+0.5,(Double_t)phtb[ic]);
 921               }
 922               else {
 923                 if(phtb[ic] > 0.0) {
 924                   fPH2dTest->Fill((Double_t)(ic/10.0),detector+0.5,(Double_t)phtb[ic]);
 925                   fPH2dSum->Fill((Double_t)(ic/10.0),0.0,(Double_t)phtb[ic]);
 926                   fPH2dSM->Fill((Double_t)(ic/10.0),sector+0.5,(Double_t)phtb[ic]);
 927                 }
 928               }
 929             }
 930           }
 931         }
 932         if(detector == 0) FindP1TrackPHtrackletV1Test(tracklet,nbclusters);
 933
 934       } // loop on tracklets
 935
 936
 937     } // debug
 938
 939     if(nTRDtrackV1 > 0) {
 940       ++nbTrdTracks;
 941       if((status&(AliESDtrack::kTRDout)) && (!(status&(AliESDtrack::kTRDin)))) {
 942         ++nbTrdTracksStandalone;
 943       }
 944       if((status&(AliESDtrack::kTRDin))) {
 945         ++nbTrdTracksOffline;
 946       }
 947     }
 948     //delete fFriendTrack;
 949   } // loop ESD track
 950
 951   if(fDebug > 1) {
 952     fNbTRDTrack->Fill(nbTrdTracks);
 953     fNbTRDTrackStandalone->Fill(nbTrdTracksStandalone);
 954     fNbTRDTrackOffline->Fill(nbTrdTracksOffline);
 955     fNbTPCTRDtrack->Fill(nbTrdTracks,nbtrackTPC);
 956   }
 957
 958   // Post output data
 959   PostData(1, fListHist);
 960   //cout << "AliTRDCalibTask::Exec() OUT" << endl;
 961 }
 962
 963 //________________________________________________________________________
 964 void AliTRDCalibTask::Terminate(Option_t *)
 965 {
 966   //
 967   // Terminate
 968   //
 969
 970   if(fTRDCalibraFillHisto) fTRDCalibraFillHisto->DestroyDebugStreamer();
 971
 972
 973 }
 974 //_______________________________________________________
 975 Bool_t AliTRDCalibTask::Load(const Char_t *filename)
 976 {
 977   //
 978   // Generic container loader
 979   //
 980
 981   if(!TFile::Open(filename)){
 982     //AliWarning(Form("Couldn't open file %s.", filename));
 983     return kFALSE;
 984   }
 985   TList *o = 0x0;
 986   if(!(o = (TList*)gFile->Get(GetName()))){
 987     //AliWarning("Missing histogram container.");
 988     return kFALSE;
 989   }
 990   fListHist = (TList*)o->Clone(GetName());
 991   gFile->Close();
 992   return kTRUE;
 993 }
 994 //_______________________________________________________
 995 Bool_t AliTRDCalibTask::Load(TList *lister)
 996 {
 997   //
 998   // Generic container loader
 999   //
1000
1001   fListHist = (TList*)lister->Clone(GetName());
1002   return kTRUE;
1003 }
1004 //_______________________________________________________________________________________
1005 void  AliTRDCalibTask::AddTask(const AliTRDCalibTask * calibTask) {
1006
1007   //
1008   // Add stats
1009   //
1010
1011   TList *listcalibTask = calibTask->GetList();
1012   if(!listcalibTask) return;
1013
1014   THnSparseI *histoEntries = (THnSparseI *) listcalibTask->FindObject("NumberOfEntries");
1015
1016   TH1I *nEvents  = (TH1I *) listcalibTask->FindObject(Form("NEvents_%s",(const char*)calibTask->GetName()));
1017   TH1I *nEventsInput  = (TH1I *) listcalibTask->FindObject(Form("NEventsInput_%s",(const char*)calibTask->GetName()));
1018   TH2F *absoluteGain  = (TH2F *) listcalibTask->FindObject(Form("AbsoluteGain_%s",(const char*)calibTask->GetName()));
1019
1020   TH1F *trdTrack = (TH1F *) listcalibTask->FindObject(Form("TRDTrack_%s",(const char*)calibTask->GetName()));
1021   TH1F *trdTrackOffline = (TH1F *) listcalibTask->FindObject(Form("TRDTrackOffline_%s",(const char*)calibTask->GetName()));
1022   TH1F *trdTrackStandalone = (TH1F *) listcalibTask->FindObject(Form("TRDTrackStandalone_%s",(const char*)calibTask->GetName()));
1023
1024   TH2F *tpctrdTrack = (TH2F *) listcalibTask->FindObject(Form("NbTPCTRDtrack_%s",(const char*)calibTask->GetName()));
1025
1026   TH1F *nbTimeBin = (TH1F *) listcalibTask->FindObject(Form("NbTimeBin_%s",(const char*)calibTask->GetName()));
1027   TH1F *nbTimeBinOffline = (TH1F *) listcalibTask->FindObject(Form("NbTimeBinOffline_%s",(const char*)calibTask->GetName()));
1028   TH1F *nbTimeBinStandalone = (TH1F *) listcalibTask->FindObject(Form("NbTimeBinStandalone_%s",(const char*)calibTask->GetName()));
1029
1030   TH1F *nbClusters = (TH1F *) listcalibTask->FindObject(Form("NbClusters_%s",(const char*)calibTask->GetName()));
1031   TH1F *nbClustersOffline = (TH1F *) listcalibTask->FindObject(Form("NbClustersOffline_%s",(const char*)calibTask->GetName()));
1032   TH1F *nbClustersStandalone = (TH1F *) listcalibTask->FindObject(Form("NbClustersStandalone_%s",(const char*)calibTask->GetName()));
1033
1034   TH1F *nbTracklets = (TH1F *) listcalibTask->FindObject(Form("NbTracklets_%s",(const char*)calibTask->GetName()));
1035   TH1F *nbTrackletsOffline = (TH1F *) listcalibTask->FindObject(Form("NbTrackletsOffline_%s",(const char*)calibTask->GetName()));
1036   TH1F *nbTrackletsStandalone = (TH1F *) listcalibTask->FindObject(Form("NbTrackletsStandalone_%s",(const char*)calibTask->GetName()));
1037
1038   TH2I *ch2d = (TH2I *) listcalibTask->FindObject("CH2d");
1039   TProfile2D *ph2d = (TProfile2D *) listcalibTask->FindObject("PH2d");
1040   TProfile2D *prf2d = (TProfile2D *) listcalibTask->FindObject("PRF2d");
1041
1042   TH2I *ch2dSum = (TH2I *) listcalibTask->FindObject(Form("CH2dSum_%s",(const char*)calibTask->GetName()));
1043   TProfile2D *ph2dSum = (TProfile2D *) listcalibTask->FindObject(Form("PH2dSum_%s",(const char*)calibTask->GetName()));
1044
1045   TH2I *ch2dSM = (TH2I *) listcalibTask->FindObject(Form("CH2dSM_%s",(const char*)calibTask->GetName()));
1046   TProfile2D *ph2dSM = (TProfile2D *) listcalibTask->FindObject(Form("PH2dSM_%s",(const char*)calibTask->GetName()));
1047
1048   AliTRDCalibraVdriftLinearFit *linearfit = (AliTRDCalibraVdriftLinearFit *) listcalibTask->FindObject("AliTRDCalibraVdriftLinearFit");
1049   AliTRDCalibraExbAltFit *exbaltfit = (AliTRDCalibraExbAltFit *) listcalibTask->FindObject("AliTRDCalibraExbAltFit");
1050   AliTRDCalibraVector *calibraVector = (AliTRDCalibraVector *) listcalibTask->FindObject("AliTRDCalibraVector");
1051
1052   //
1053
1054   THnSparseI *inhistoEntries = (THnSparseI *) fListHist->FindObject("NumberOfEntries");
1055
1056   TH1I *inEventsInput  = (TH1I *) fListHist->FindObject(Form("NEventsInput_%s",(const char*)fName));
1057   TH1I *inEvents  = (TH1I *) fListHist->FindObject(Form("NEvents_%s",(const char*)fName));
1058   TH2F *iabsoluteGain  = (TH2F *) fListHist->FindObject(Form("AbsoluteGain_%s",(const char*)fName));
1059
1060   TH1F *itrdTrack = (TH1F *) fListHist->FindObject(Form("TRDTrack_%s",(const char*)fName));
1061   TH1F *itrdTrackOffline = (TH1F *) fListHist->FindObject(Form("TRDTrackOffline_%s",(const char*)fName));
1062   TH1F *itrdTrackStandalone = (TH1F *) fListHist->FindObject(Form("TRDTrackStandalone_%s",(const char*)fName));
1063
1064   TH2F *itpctrdTrack = (TH2F *) fListHist->FindObject(Form("NbTPCTRDtrack_%s",(const char*)fName));
1065
1066   TH1F *inbTimeBin = (TH1F *) fListHist->FindObject(Form("NbTimeBin_%s",(const char*)fName));
1067   TH1F *inbTimeBinOffline = (TH1F *) fListHist->FindObject(Form("NbTimeBinOffline_%s",(const char*)fName));
1068   TH1F *inbTimeBinStandalone = (TH1F *) fListHist->FindObject(Form("NbTimeBinStandalone_%s",(const char*)fName));
1069
1070   TH1F *inbClusters = (TH1F *) fListHist->FindObject(Form("NbClusters_%s",(const char*)fName));
1071   TH1F *inbClustersOffline = (TH1F *) fListHist->FindObject(Form("NbClustersOffline_%s",(const char*)fName));
1072   TH1F *inbClustersStandalone = (TH1F *) fListHist->FindObject(Form("NbClustersStandalone_%s",(const char*)fName));
1073
1074   TH1F *inbTracklets = (TH1F *) fListHist->FindObject(Form("NbTracklets_%s",(const char*)fName));
1075   TH1F *inbTrackletsOffline = (TH1F *) fListHist->FindObject(Form("NbTrackletsOffline_%s",(const char*)fName));
1076   TH1F *inbTrackletsStandalone = (TH1F *) fListHist->FindObject(Form("NbTrackletsStandalone_%s",(const char*)fName));
1077
1078   TH2I *ich2d = (TH2I *) fListHist->FindObject("CH2d");
1079   TProfile2D *iph2d = (TProfile2D *) fListHist->FindObject("PH2d");
1080   TProfile2D *iprf2d = (TProfile2D *) fListHist->FindObject("PRF2d");
1081
1082   TH2I *ich2dSum = (TH2I *) fListHist->FindObject(Form("CH2dSum_%s",(const char*)fName));
1083   TProfile2D *iph2dSum = (TProfile2D *) fListHist->FindObject(Form("PH2dSum_%s",(const char*)fName));
1084
1085   TH2I *ich2dSM = (TH2I *) fListHist->FindObject(Form("CH2dSM_%s",(const char*)fName));
1086   TProfile2D *iph2dSM = (TProfile2D *) fListHist->FindObject(Form("PH2dSM_%s",(const char*)fName));
1087
1088   AliTRDCalibraVdriftLinearFit *ilinearfit = (AliTRDCalibraVdriftLinearFit *) fListHist->FindObject("AliTRDCalibraVdriftLinearFit");
1089   AliTRDCalibraExbAltFit *iexbaltfit = (AliTRDCalibraExbAltFit *) fListHist->FindObject("AliTRDCalibraExbAltFit");
1090   AliTRDCalibraVector *icalibraVector = (AliTRDCalibraVector *) fListHist->FindObject("AliTRDCalibraVector");
1091
1092
1093   // Add
1094
1095   if(histoEntries) {
1096     if(inhistoEntries) {
1097       inhistoEntries->Add(histoEntries);
1098       //printf("Add Events\n");
1099     }
1100     else {
1101       //printf("Create new Events\n");
1102       inhistoEntries = (THnSparseI *) histoEntries->Clone();
1103       fListHist->Add(inhistoEntries);
1104     }
1105   }
1106
1107   if(nEventsInput) {
1108     if(inEventsInput) {
1109       inEventsInput->Add(nEventsInput);
1110       //printf("Add Events\n");
1111     }
1112     else {
1113       //printf("Create new Events\n");
1114       inEventsInput = new TH1I(*nEventsInput);
1115       fListHist->Add(inEventsInput);
1116     }
1117   }
1118
1119   if(nEvents) {
1120     if(inEvents) {
1121       inEvents->Add(nEvents);
1122       //printf("Add Events\n");
1123     }
1124     else {
1125       //printf("Create new Events\n");
1126       inEvents = new TH1I(*nEvents);
1127       fListHist->Add(inEvents);
1128     }
1129   }
1130
1131   if(absoluteGain) {
1132     if(iabsoluteGain) iabsoluteGain->Add(absoluteGain);
1133     else {
1134       iabsoluteGain = new TH2F(*absoluteGain);
1135       fListHist->Add(iabsoluteGain);
1136     }
1137   }
1138
1139   if(trdTrack) {
1140     if(itrdTrack) itrdTrack->Add(trdTrack);
1141     else {
1142      itrdTrack = new TH1F(*trdTrack);
1143      fListHist->Add(itrdTrack);
1144     }
1145   }
1146
1147   if(trdTrackOffline) {
1148     if(itrdTrackOffline) itrdTrackOffline->Add(trdTrackOffline);
1149     else {
1150       itrdTrackOffline = new TH1F(*trdTrackOffline);
1151       fListHist->Add(itrdTrackOffline);
1152     }
1153   }
1154
1155   if(trdTrackStandalone) {
1156     if(itrdTrackStandalone) itrdTrackStandalone->Add(trdTrackStandalone);
1157     else {
1158       itrdTrackStandalone = new TH1F(*trdTrackStandalone);
1159       fListHist->Add(itrdTrackStandalone);
1160     }
1161   }
1162
1163   if(tpctrdTrack) {
1164     if(itpctrdTrack) itpctrdTrack->Add(tpctrdTrack);
1165     else {
1166       itpctrdTrack = new TH2F(*tpctrdTrack);
1167       fListHist->Add(itpctrdTrack);
1168     }
1169   }
1170
1171   if(nbTimeBin) {
1172     if(inbTimeBin) inbTimeBin->Add(nbTimeBin);
1173     else {
1174       inbTimeBin = new TH1F(*inbTimeBin);
1175       fListHist->Add(inbTimeBin);
1176     }
1177   }
1178
1179   if(nbTimeBinOffline) {
1180     if(inbTimeBinOffline) inbTimeBinOffline->Add(nbTimeBinOffline);
1181     else {
1182       inbTimeBinOffline = new TH1F(*nbTimeBinOffline);
1183       fListHist->Add(inbTimeBinOffline);
1184     }
1185   }
1186
1187   if(nbTimeBinStandalone) {
1188     if(inbTimeBinStandalone) inbTimeBinStandalone->Add(nbTimeBinStandalone);
1189     else {
1190       inbTimeBinStandalone = new TH1F(*nbTimeBinStandalone);
1191       fListHist->Add(inbTimeBinStandalone);
1192     }
1193   }
1194
1195   if(nbClusters) {
1196     if(inbClusters) inbClusters->Add(nbClusters);
1197     else {
1198       inbClusters = new TH1F(*nbClusters);
1199       fListHist->Add(inbClusters);
1200     }
1201   }
1202
1203   if(nbClustersOffline) {
1204     if(inbClustersOffline) inbClustersOffline->Add(nbClustersOffline);
1205     else {
1206       inbClustersOffline = new TH1F(*nbClustersOffline);
1207       fListHist->Add(inbClustersOffline);
1208     }
1209   }
1210
1211   if(nbClustersStandalone) {
1212     if(inbClustersStandalone) inbClustersStandalone->Add(nbClustersStandalone);
1213     else {
1214       inbClustersStandalone = new TH1F(*nbClustersStandalone);
1215       fListHist->Add(inbClustersStandalone);
1216     }
1217   }
1218
1219   if(nbTracklets) {
1220     if(inbTracklets) inbTracklets->Add(nbTracklets);
1221     else {
1222       inbTracklets = new TH1F(*nbTracklets);
1223       fListHist->Add(inbTracklets);
1224     }
1225   }
1226
1227   if(nbTrackletsOffline) {
1228     if(inbTrackletsOffline) inbTrackletsOffline->Add(nbTrackletsOffline);
1229     else {
1230       inbTrackletsOffline = new TH1F(*nbTrackletsOffline);
1231       fListHist->Add(inbTrackletsOffline);
1232     }
1233   }
1234
1235   if(nbTrackletsStandalone) {
1236     if(inbTrackletsStandalone) inbTrackletsStandalone->Add(nbTrackletsStandalone);
1237     else {
1238       inbTrackletsStandalone = new TH1F(*nbTrackletsStandalone);
1239       fListHist->Add(inbTrackletsStandalone);
1240     }
1241   }
1242
1243   if(ch2d) {
1244     if(ich2d) ich2d->Add(ch2d);
1245     else {
1246       ich2d = new TH2I(*ch2d);
1247       fListHist->Add(ich2d);
1248     }
1249   }
1250
1251   if(ph2d) {
1252     if(iph2d) iph2d->Add(ph2d);
1253     else {
1254       iph2d = new TProfile2D(*ph2d);
1255       fListHist->Add(iph2d);
1256     }
1257   }
1258
1259   if(prf2d) {
1260     if(iprf2d) iprf2d->Add(prf2d);
1261     else {
1262       iprf2d = new TProfile2D(*prf2d);
1263       fListHist->Add(iprf2d);
1264     }
1265   }
1266
1267   if(ch2dSum) {
1268     if(ich2dSum) ich2dSum->Add(ch2dSum);
1269     else {
1270       ich2dSum = new TH2I(*ch2dSum);
1271       fListHist->Add(ich2dSum);
1272     }
1273   }
1274
1275   if(ph2dSum) {
1276     if(iph2dSum) iph2dSum->Add(ph2dSum);
1277     else {
1278       iph2dSum = new TProfile2D(*ph2dSum);
1279       fListHist->Add(iph2dSum);
1280     }
1281   }
1282
1283   if(ch2dSM) {
1284     if(ich2dSM) ich2dSM->Add(ch2dSM);
1285     else {
1286       ich2dSM = new TH2I(*ch2dSM);
1287       fListHist->Add(ich2dSM);
1288     }
1289   }
1290
1291   if(ph2dSM) {
1292     if(iph2dSM) iph2dSM->Add(ph2dSM);
1293     else {
1294       iph2dSM = new TProfile2D(*ph2dSM);
1295       fListHist->Add(iph2dSM);
1296     }
1297   }
1298
1299   if(linearfit) {
1300     if(ilinearfit) ilinearfit->Add(linearfit);
1301     else {
1302       ilinearfit = new AliTRDCalibraVdriftLinearFit(*linearfit);
1303       fListHist->Add(ilinearfit);
1304     }
1305   }
1306
1307   if(exbaltfit) {
1308     if(iexbaltfit) iexbaltfit->Add(exbaltfit);
1309     else {
1310       iexbaltfit = new AliTRDCalibraExbAltFit(*exbaltfit);
1311       fListHist->Add(iexbaltfit);
1312     }
1313   }
1314
1315   if(calibraVector) {
1316     if(icalibraVector) icalibraVector->Add(calibraVector);
1317     else {
1318       icalibraVector = new AliTRDCalibraVector(*calibraVector);
1319       fListHist->Add(icalibraVector);
1320     }
1321   }
1322
1323 }
1324 //________________________________________________________________________________
1325 Long64_t AliTRDCalibTask::Merge(TCollection *li) {
1326
1327   //
1328   // merge component
1329   //
1330
1331   TIterator* iter = li->MakeIterator();
1332   AliTRDCalibTask* cal = 0;
1333
1334   while ((cal = (AliTRDCalibTask*)iter->Next())) {
1335     if (!cal->InheritsFrom(AliTRDCalibTask::Class())) {
1336       Error("Merge","Attempt to add object of class %s to a %s", cal->ClassName(), this->ClassName());
1337       return -1;
1338     }
1339
1340     // add histograms here...
1341     this->AddTask(cal);
1342
1343   }
1344
1345   return 0;
1346
1347 }
1348 //_____________________________________________________
1349 Bool_t AliTRDCalibTask::SetVersionSubversion(){
1350   //
1351   // Load Chamber Gain factors into the Tender supply
1352   //
1353
1354   //printf("SetVersionSubversion\n");
1355
1356   //find previous entry from the UserInfo
1357   TTree *tree=((TChain*)GetInputData(0))->GetTree();
1358   if (!tree) {
1359     AliError("Tree not found in ESDhandler");
1360     return kFALSE;
1361   }
1362
1363   TList *userInfo=(TList*)tree->GetUserInfo();
1364   if (!userInfo) {
1365     AliError("No UserInfo found in tree");
1366     return kFALSE;
1367   }
1368
1369   TList *cdbList=(TList*)userInfo->FindObject("cdbList");
1370   if (!cdbList) {
1371     AliError("No cdbList found in UserInfo");
1372     if (AliLog::GetGlobalLogLevel()>=AliLog::kError) userInfo->Print();
1373     return kFALSE;
1374   }
1375
1376   TIter nextCDB(cdbList);
1377   TObjString *os=0x0;
1378   while ( (os=(TObjString*)nextCDB()) ){
1379     if(os->GetString().Contains("TRD/Calib/ChamberGainFactor")){
1380       // Get Old gain calibration
1381       AliCDBId *id=AliCDBId::MakeFromString(os->GetString());
1382       fFirstRunGain = id->GetFirstRun();
1383       fVersionGainUsed = id->GetVersion();
1384       fSubVersionGainUsed = id->GetSubVersion();
1385     } else if(os->GetString().Contains("TRD/Calib/ChamberVdrift")){
1386       // Get Old drift velocity calibration
1387       AliCDBId *id=AliCDBId::MakeFromString(os->GetString());
1388       fFirstRunVdrift = id->GetFirstRun();
1389       fVersionVdriftUsed = id->GetVersion();
1390       fSubVersionVdriftUsed = id->GetSubVersion();
1391     } else if(os->GetString().Contains("TRD/Calib/LocalGainFactor")){
1392       // Get Old drift velocity calibration
1393       AliCDBId *id=AliCDBId::MakeFromString(os->GetString());
1394       fFirstRunGainLocal = id->GetFirstRun();
1395       fVersionGainLocalUsed = id->GetVersion();
1396       fSubVersionGainLocalUsed = id->GetSubVersion();
1397     } else if((os->GetString().Contains("TRD/Calib/ChamberExB")) && (!os->GetString().Contains("TRD/Calib/ChamberExBAlt"))){
1398       // Get Old drift velocity calibration
1399       AliCDBId *id=AliCDBId::MakeFromString(os->GetString());
1400       fFirstRunExB = id->GetFirstRun();
1401       fVersionExBUsed = id->GetVersion();
1402       fSubVersionExBUsed = id->GetSubVersion();
1403       //printf("Version %d and subversion %d\n",fVersionExBUsed,fSubVersionExBUsed);
1404     }
1405   }
1406
1407   //printf("VersionGain %d, SubversionGain %d, VersionLocalGain %d, Subversionlocalgain %d, Versionvdrift %d, Subversionvdrift %d\n",fVersionGainUsed,fSubVersionGainUsed,fVersionGainLocalUsed,fSubVersionGainLocalUsed,fVersionVdriftUsed,fSubVersionVdriftUsed);
1408
1409   // Check
1410   if((fFirstRunGain < 0)            ||
1411      (fFirstRunGainLocal < 0)       ||
1412      (fFirstRunVdrift < 0)          ||
1413      (fVersionGainUsed < 0)         ||
1414      (fVersionGainLocalUsed < 0)    ||
1415      (fSubVersionGainUsed < 0)      ||
1416      (fSubVersionGainLocalUsed < 0) ||
1417      (fVersionVdriftUsed < 0)       ||
1418      (fSubVersionVdriftUsed < 0)) {
1419     AliError("No recent calibration found");
1420     return kFALSE;
1421   }
1422   else return kTRUE;
1423
1424 }
1425 //_________________________________________________________________________________________________________________________
1426 Bool_t AliTRDCalibTask::ParticleGood(int i) const {
1427
1428   //
1429   // Definition of good tracks
1430   //
1431
1432
1433   AliESDtrack *track = fESD->GetTrack(i);
1434   if (!track->IsOn(AliESDtrack::kTPCrefit)) return 0;        // TPC refit
1435   if (track->GetTPCNcls() < 90) return 0;                    // number of TPC clusters
1436   if (fabs(track->Eta())>0.8) return 0;                         // fiducial pseudorapidity
1437   Float_t r,z;
1438   track->GetImpactParametersTPC(r,z);
1439   if (fabs(z)>2.0) return 0;                          // impact parameter in z
1440   if (fabs(r)>2.0) return 0;                          // impact parameter in xy
1441   if (r==0) return 0;
1442   return 1;
1443
1444
1445 }
1446 //______________________________________________________________________________________________________________________
1447 Bool_t AliTRDCalibTask::FindP1TrackPHtrackletV1Test(const AliTRDseedV1 *tracklet, Int_t nbclusters)
1448 {
1449   //
1450   // Drift velocity calibration:
1451   // Fit the clusters with a straight line
1452   // From the slope find the drift velocity
1453   //
1454
1455   ////////////////////////////////////////////////
1456   //Number of points: if less than 3 return kFALSE
1457   /////////////////////////////////////////////////
1458   if(nbclusters <= 2) return kFALSE;
1459
1460   ////////////
1461   //Variables
1462   ////////////
1463   // results of the linear fit
1464   Double_t dydt                       = 0.0;                                // dydt tracklet after straight line fit
1465   Double_t errorpar                   = 0.0;                                // error after straight line fit on dy/dt
1466   Double_t pointError                 = 0.0;                                // error after straight line fit
1467   // pad row problemes: avoid tracklet that cross pad rows, tilting angle in the constant
1468   Int_t    crossrow                   = 0;                                  // if it crosses a pad row
1469   Int_t    rowp                       = -1;                                 // if it crosses a pad row
1470   Float_t  tnt                        = tracklet->GetTilt();                // tan tiltingangle
1471   TLinearFitter linearFitterTracklet(2,"pol1");
1472   linearFitterTracklet.StoreData(kTRUE);
1473
1474
1475   ///////////////////////////////////////////
1476   // Take the parameters of the track
1477   //////////////////////////////////////////
1478   // take now the snp, tnp and tgl from the track
1479   Double_t snp = tracklet->GetSnp();             // sin dy/dx at the end of the chamber
1480   Double_t tnp = 0.0;                            // dy/dx at the end of the chamber
1481   if( TMath::Abs(snp) <  1.){
1482     tnp = snp / TMath::Sqrt((1.-snp)*(1.+snp));
1483   }
1484   Double_t tgl  = tracklet->GetTgl();           // dz/dl
1485   Double_t dzdx = tgl*TMath::Sqrt(1+tnp*tnp);   // dz/dx calculated from dz/dl
1486   // at the entrance
1487   //Double_t tnp = tracklet->GetYref(1);      // dy/dx at the entrance of the chamber
1488   //Double_t tgl = tracklet->GetZref(1);      // dz/dl at the entrance of the chamber
1489   //Double_t dzdx = tgl;                      //*TMath::Sqrt(1+tnp*tnp); // dz/dx from dz/dl
1490   // at the end with correction due to linear fit
1491   //Double_t tnp = tracklet->GetYfit(1);      // dy/dx at the end of the chamber after fit correction
1492   //Double_t tgl = tracklet->GetZfit(1);      // dz/dl at the end of the chamber after fit correction
1493
1494
1495   ////////////////////////////
1496   // loop over the clusters
1497   ////////////////////////////
1498   Int_t  nbli = 0;
1499   AliTRDcluster *cl                   = 0x0;
1500   //////////////////////////////
1501   // Check no shared clusters
1502   //////////////////////////////
1503   for(int icc=AliTRDseedV1::kNtb; icc<AliTRDseedV1::kNclusters; icc++){
1504     cl = tracklet->GetClusters(icc);
1505     if(cl)  crossrow = 1;
1506   }
1507   //////////////////////////////////
1508   // Loop clusters
1509   //////////////////////////////////
1510   for(int ic=0; ic<AliTRDseedV1::kNtb; ic++){
1511     if(!(cl = tracklet->GetClusters(ic))) continue;
1512     //if((fLimitChargeIntegration) && (!cl->IsInChamber())) continue;
1513
1514     Double_t ycluster                 = cl->GetY();
1515     Int_t time                        = cl->GetPadTime();
1516     Double_t timeis                   = time/10.0;
1517     //See if cross two pad rows
1518     Int_t    row                      = cl->GetPadRow();
1519     if(rowp==-1) rowp                 = row;
1520     if(row != rowp) crossrow          = 1;
1521
1522     linearFitterTracklet.AddPoint(&timeis,ycluster,1);
1523     nbli++;
1524
1525
1526   }
1527
1528   ////////////////////////////////////
1529   // Do the straight line fit now
1530   ///////////////////////////////////
1531   if(nbli <= 2){
1532     linearFitterTracklet.ClearPoints();
1533     return kFALSE;
1534   }
1535   TVectorD pars;
1536   linearFitterTracklet.Eval();
1537   linearFitterTracklet.GetParameters(pars);
1538   pointError  =  TMath::Sqrt(linearFitterTracklet.GetChisquare()/(nbli-2));
1539   errorpar    =  linearFitterTracklet.GetParError(1)*pointError;
1540   dydt        = pars[1];
1541   //printf("chis %f, nbli %d, pointError %f, parError %f, errorpar %f\n",linearFitterTracklet->GetChisquare(),nbli,pointError,linearFitterTracklet->GetParError(1),errorpar);
1542   linearFitterTracklet.ClearPoints();
1543
1544   /////////////////////////
1545   // Cuts quality
1546   ////////////////////////
1547
1548   if(nbclusters < fLow) return kFALSE;
1549   if(nbclusters > fHigh) return kFALSE;
1550   if(pointError >= 0.3) return kFALSE;
1551   if(crossrow == 1) return kTRUE;
1552
1553   ///////////////////////
1554   // Fill
1555   //////////////////////
1556
1557   if(fDebug > 0){
1558     //Add to the linear fitter of the detector
1559     if( TMath::Abs(snp) <  1.){
1560       Double_t x = tnp-dzdx*tnt;
1561       //if(!fLinearVdriftTest) printf("Not there\n");
1562       Double_t nbentries = fLinearVdriftTest->GetEntries();
1563       if(nbentries < (5.0*32767)) fLinearVdriftTest->Fill(x,dydt);
1564     }
1565   }
1566
1567   return kTRUE;
1568 }
1569