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