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   Int_t nbTracksfriends = fESDfriend->GetNumberOfTracks();
 723   for(int itrk=0; itrk < nbTracksfriends; ++itrk){
 724
 725     // Get ESD track
 726     fkEsdTrack = fESD->GetTrack(itrk);
 727     if(!fkEsdTrack) continue;
 728     ULong_t status = fkEsdTrack->GetStatus();
 729     if(status&(AliESDtrack::kTPCout)) ++nbtrackTPC;
 730
 731     // Fix suggested by Alex
 732     fFriendTrack = fESDfriend->GetTrack(itrk);
 733     //printf("itrk %d\n",itrk);
 734     //fFriendTrack = (fESDfriend->GetNumberOfTracks()>itrk)?fESDfriend->GetTrack(itrk):NULL;
 735     if(!fFriendTrack)  {
 736       //printf("No friend track %d\n",itrk);
 737       continue;
 738     }
 739
 740     // Other cuts
 741     fTrdTrack = 0x0;
 742     Bool_t good = kTRUE;
 743     Bool_t standalonetrack = kFALSE;
 744     Bool_t offlinetrack = kFALSE;
 745     //ULong_t status = fkEsdTrack->GetStatus();
 746
 747     //////////////////////////////////////
 748     // Loop on calibration objects
 749     //////////////////////////////////////
 750     Int_t icalib=0;
 751     Int_t nTRDtrackV1=0;
 752     while((fCalibObject = (TObject *)(fFriendTrack->GetCalibObject(icalib++)))){
 753       //printf("Name %s\n",fCalibObject->IsA()->GetName());
 754       if(strcmp(fCalibObject->IsA()->GetName(), "AliTRDtrackV1") != 0) continue;
 755       //printf("Find the calibration object\n");
 756       ++nTRDtrackV1;
 757
 758       if((status&(AliESDtrack::kTRDout)) && (!(status&(AliESDtrack::kTRDin)))) {
 759         standalonetrack = kTRUE;
 760       }
 761       if((status&(AliESDtrack::kTRDin))) {
 762         offlinetrack = kTRUE;
 763       }
 764       if(fOfflineTracks){
 765         if(!offlinetrack){
 766           good = kFALSE;
 767         }
 768       }
 769       else if(fStandaloneTracks){
 770         if(!standalonetrack){
 771           good = kFALSE;
 772         }
 773       }
 774
 775       fTrdTrack = (AliTRDtrackV1 *)fCalibObject;
 776       // process chamberstatus
 777       fTRDChamberStatus->ProcessTrack(fTrdTrack);
 778     }
 779
 780     // Quality cuts on the AliESDtrack
 781     if((fEsdTrackCuts) && (!fEsdTrackCuts->IsSelected((AliVParticle *)fkEsdTrack))) {
 782       //printf("Not a good track\n");
 783       continue;
 784     }
 785
 786     // First Absolute gain calibration
 787     Int_t trdNTracklets = (Int_t) fkEsdTrack->GetTRDntracklets();
 788     Int_t trdNTrackletsPID = (Int_t) fkEsdTrack->GetTRDntrackletsPID();
 789     //printf("Number of trd tracklets %d and PID trd tracklets %d\n",trdNTracklets,trdNTrackletsPID);
 790     if((trdNTracklets > 0) && (trdNTrackletsPID > 0)) {
 791       for(Int_t iPlane = 0; iPlane < 6; ++iPlane){
 792         //Double_t slide = fkEsdTrack->GetTRDslice(iPlane);
 793         //printf("Number of slide %d\n",fkEsdTrack->GetNumberOfTRDslices());
 794         //Double_t momentum = fkEsdTrack->GetTRDmomentum(iPlane);
 795         //printf("momentum %f, slide %f\n",momentum,slide);
 796         if(fkEsdTrack->GetTRDslice(iPlane) > 0.0)
 797           fAbsoluteGain->Fill(fkEsdTrack->GetTRDslice(iPlane)*8.0/100.0,
 798                               fkEsdTrack->GetTRDmomentum(iPlane));
 799       }
 800     }
 801
 802
 803     if(!fTrdTrack) continue;
 804
 805     if(good && fOnInstance) {
 806       //cout << "good" << endl;
 807       fTRDCalibraFillHisto->UpdateHistogramsV1(fTrdTrack,fkEsdTrack);
 808       //printf("Fill fTRDCalibraFillHisto\n");
 809     }
 810
 811     if(IsPHQon()){
 812       const Double_t mag = AliTRDdEdxBaseUtils::IsExBOn() ? fESD->GetMagneticField() : -1;
 813       const Int_t charge = AliTRDdEdxBaseUtils::IsExBOn() ? fkEsdTrack->Charge() : -1;
 814       const Double_t toTPCscale = AliTRDdEdxCalibUtils::GetCalibTPCscale(fkEsdTrack->GetTPCncls(), fkEsdTrack->GetTPCsignal());
 815       if(toTPCscale>0){
 816         AliTRDdEdxCalibUtils::FillHist(fTrdTrack, 0, mag, charge, toTPCscale);
 817       }
 818     }
 819
 820     //////////////////////////////////
 821     // Debug
 822     ////////////////////////////////
 823
 824     if(fDebug > 0) {
 825
 826       //printf("Enter debug\n");
 827
 828       Int_t nbtracklets = 0;
 829
 830       // Check some stuff
 831       Bool_t standalonetracklet = kFALSE;
 832       const AliTRDseedV1 *tracklet = 0x0;
 833       //////////////////////////////////////
 834       // Loop tracklets
 835       /////////////////////////////////////
 836       Int_t nbclusters=0;
 837       Double_t phtb[AliTRDseedV1::kNtb];
 838       memset(phtb, 0, AliTRDseedV1::kNtb*sizeof(Double_t));
 839       Double_t sum = 0.0;
 840       Float_t normalisation = 1.13;
 841       Int_t detector = 0;
 842       Int_t sector = 0;
 843       for(Int_t itr = 0; itr < 6; ++itr){
 844
 845         if(!(tracklet = fTrdTrack->GetTracklet(itr))) continue;
 846         if(!tracklet->IsOK()) continue;
 847         ++nbtracklets;
 848         standalonetracklet = kFALSE;
 849         if(tracklet->IsStandAlone()) standalonetracklet = kTRUE;
 850
 851         nbclusters = 0;
 852         memset(phtb, 0, AliTRDseedV1::kNtb*sizeof(Double_t));
 853         sum = 0.0;
 854         detector = 0;
 855         sector = 0;
 856         //Int_t crossrow = 0;
 857
 858           // Check no shared clusters
 859           //for(int icc=AliTRDseedV1::kNtb; icc<AliTRDseedV1::kNclusters; icc++){
 860           //  if((fcl = tracklet->GetClusters(icc)))  crossrow = 1;
 861           // }
 862
 863           // Loop on clusters
 864         Int_t time = 0;
 865         Float_t ch = 0;
 866         Float_t qcl = 0;
 867         for(int ic=0; ic<AliTRDseedV1::kNtb; ++ic){
 868
 869           if(!(fCl = tracklet->GetClusters(ic))) continue;
 870           ++nbclusters;
 871           time = fCl->GetPadTime();
 872           //ch =  tracklet->GetdQdl(ic);
 873           ch =  tracklet->GetQperTB(ic);
 874           qcl = TMath::Abs(fCl->GetQ());
 875           detector = fCl->GetDetector();
 876           // Add the charge if shared cluster
 877           if((ic+AliTRDseedV1::kNtb) < AliTRDseedV1::kNclusters) {
 878             if((fCl = tracklet->GetClusters(ic+AliTRDseedV1::kNtb))) {
 879               qcl += TMath::Abs(fCl->GetQ());
 880               //printf("Add the cluster charge\n");
 881             }
 882           }
 883           if((time>-1) && (time<fNbTimeBins)) phtb[time]=qcl;
 884           if((fCalDetGain) && (fCalDetGain->GetValue(detector) > 0.0)) sum += ch*fCalDetGain->GetValue(detector)/normalisation;
 885           else sum += ch/normalisation;
 886
 887           if(fDebug > 1) {
 888             fNbTimeBin->Fill(time);
 889             if(tracklet->IsStandAlone()) fNbTimeBinStandalone->Fill(time);
 890             else fNbTimeBinOffline->Fill(time);
 891           }
 892         }
 893         sector = AliTRDgeometry::GetSector(detector);
 894
 895         if(fDebug > 1) {
 896           fNbTracklets->Fill(detector);
 897           if(tracklet->IsStandAlone()) fNbTrackletsStandalone->Fill(detector);
 898           else fNbTrackletsOffline->Fill(detector);
 899
 900           fNbClusters->Fill(nbclusters);
 901           if(tracklet->IsStandAlone())  fNbClustersStandalone->Fill(nbclusters);
 902           else  fNbClustersOffline->Fill(nbclusters);
 903         }
 904
 905         if((nbclusters > fLow) && (nbclusters < fHigh)){
 906           if(fRelativeScale > 0.0) sum = sum/fRelativeScale;
 907           fCH2dTest->Fill(sum,detector+0.5);
 908           fCH2dSM->Fill(sum,sector+0.5);
 909           fCH2dSum->Fill(sum,0.5);
 910           Bool_t checknoise = kTRUE;
 911           if(fMaxCluster > 0) {
 912             if(phtb[0] > fMaxCluster) checknoise = kFALSE;
 913             if(fNbTimeBins > fNbMaxCluster) {
 914               for(Int_t k = (fNbTimeBins-fNbMaxCluster); k < fNbTimeBins; k++){
 915                 if(phtb[k] > fMaxCluster) checknoise = kFALSE;
 916               }
 917             }
 918           }
 919           if(checknoise) {
 920             for(int ic=0; ic<fNbTimeBins; ic++){
 921               if(fFillZero) {
 922                 fPH2dTest->Fill((Double_t)(ic/10.0),detector+0.5,(Double_t)phtb[ic]);
 923                 fPH2dSum->Fill((Double_t)(ic/10.0),0.5,(Double_t)phtb[ic]);
 924                 fPH2dSM->Fill((Double_t)(ic/10.0),sector+0.5,(Double_t)phtb[ic]);
 925               }
 926               else {
 927                 if(phtb[ic] > 0.0) {
 928                   fPH2dTest->Fill((Double_t)(ic/10.0),detector+0.5,(Double_t)phtb[ic]);
 929                   fPH2dSum->Fill((Double_t)(ic/10.0),0.0,(Double_t)phtb[ic]);
 930                   fPH2dSM->Fill((Double_t)(ic/10.0),sector+0.5,(Double_t)phtb[ic]);
 931                 }
 932               }
 933             }
 934           }
 935         }
 936         if(detector == 0) FindP1TrackPHtrackletV1Test(tracklet,nbclusters);
 937
 938       } // loop on tracklets
 939
 940
 941     } // debug
 942
 943     if(nTRDtrackV1 > 0) {
 944       ++nbTrdTracks;
 945       if((status&(AliESDtrack::kTRDout)) && (!(status&(AliESDtrack::kTRDin)))) {
 946         ++nbTrdTracksStandalone;
 947       }
 948       if((status&(AliESDtrack::kTRDin))) {
 949         ++nbTrdTracksOffline;
 950       }
 951     }
 952     //delete fFriendTrack;
 953   } // loop ESD track
 954
 955   if(fDebug > 1) {
 956     fNbTRDTrack->Fill(nbTrdTracks);
 957     fNbTRDTrackStandalone->Fill(nbTrdTracksStandalone);
 958     fNbTRDTrackOffline->Fill(nbTrdTracksOffline);
 959     fNbTPCTRDtrack->Fill(nbTrdTracks,nbtrackTPC);
 960   }
 961
 962   // Post output data
 963   PostData(1, fListHist);
 964   //cout << "AliTRDCalibTask::Exec() OUT" << endl;
 965 }
 966
 967 //________________________________________________________________________
 968 void AliTRDCalibTask::Terminate(Option_t *)
 969 {
 970   //
 971   // Terminate
 972   //
 973
 974   if(fTRDCalibraFillHisto) fTRDCalibraFillHisto->DestroyDebugStreamer();
 975
 976
 977 }
 978 //_______________________________________________________
 979 Bool_t AliTRDCalibTask::Load(const Char_t *filename)
 980 {
 981   //
 982   // Generic container loader
 983   //
 984
 985   if(!TFile::Open(filename)){
 986     //AliWarning(Form("Couldn't open file %s.", filename));
 987     return kFALSE;
 988   }
 989   TList *o = 0x0;
 990   if(!(o = (TList*)gFile->Get(GetName()))){
 991     //AliWarning("Missing histogram container.");
 992     return kFALSE;
 993   }
 994   fListHist = (TList*)o->Clone(GetName());
 995   gFile->Close();
 996   return kTRUE;
 997 }
 998 //_______________________________________________________
 999 Bool_t AliTRDCalibTask::Load(TList *lister)
1000 {
1001   //
1002   // Generic container loader
1003   //
1004
1005   fListHist = (TList*)lister->Clone(GetName());
1006   return kTRUE;
1007 }
1008 //_______________________________________________________________________________________
1009 void  AliTRDCalibTask::AddTask(const AliTRDCalibTask * calibTask) {
1010
1011   //
1012   // Add stats
1013   //
1014
1015   TList *listcalibTask = calibTask->GetList();
1016   if(!listcalibTask) return;
1017
1018   THnSparseI *histoEntries = (THnSparseI *) listcalibTask->FindObject("NumberOfEntries");
1019
1020   TH1I *nEvents  = (TH1I *) listcalibTask->FindObject(Form("NEvents_%s",(const char*)calibTask->GetName()));
1021   TH1I *nEventsInput  = (TH1I *) listcalibTask->FindObject(Form("NEventsInput_%s",(const char*)calibTask->GetName()));
1022   TH2F *absoluteGain  = (TH2F *) listcalibTask->FindObject(Form("AbsoluteGain_%s",(const char*)calibTask->GetName()));
1023
1024   TH1F *trdTrack = (TH1F *) listcalibTask->FindObject(Form("TRDTrack_%s",(const char*)calibTask->GetName()));
1025   TH1F *trdTrackOffline = (TH1F *) listcalibTask->FindObject(Form("TRDTrackOffline_%s",(const char*)calibTask->GetName()));
1026   TH1F *trdTrackStandalone = (TH1F *) listcalibTask->FindObject(Form("TRDTrackStandalone_%s",(const char*)calibTask->GetName()));
1027
1028   TH2F *tpctrdTrack = (TH2F *) listcalibTask->FindObject(Form("NbTPCTRDtrack_%s",(const char*)calibTask->GetName()));
1029
1030   TH1F *nbTimeBin = (TH1F *) listcalibTask->FindObject(Form("NbTimeBin_%s",(const char*)calibTask->GetName()));
1031   TH1F *nbTimeBinOffline = (TH1F *) listcalibTask->FindObject(Form("NbTimeBinOffline_%s",(const char*)calibTask->GetName()));
1032   TH1F *nbTimeBinStandalone = (TH1F *) listcalibTask->FindObject(Form("NbTimeBinStandalone_%s",(const char*)calibTask->GetName()));
1033
1034   TH1F *nbClusters = (TH1F *) listcalibTask->FindObject(Form("NbClusters_%s",(const char*)calibTask->GetName()));
1035   TH1F *nbClustersOffline = (TH1F *) listcalibTask->FindObject(Form("NbClustersOffline_%s",(const char*)calibTask->GetName()));
1036   TH1F *nbClustersStandalone = (TH1F *) listcalibTask->FindObject(Form("NbClustersStandalone_%s",(const char*)calibTask->GetName()));
1037
1038   TH1F *nbTracklets = (TH1F *) listcalibTask->FindObject(Form("NbTracklets_%s",(const char*)calibTask->GetName()));
1039   TH1F *nbTrackletsOffline = (TH1F *) listcalibTask->FindObject(Form("NbTrackletsOffline_%s",(const char*)calibTask->GetName()));
1040   TH1F *nbTrackletsStandalone = (TH1F *) listcalibTask->FindObject(Form("NbTrackletsStandalone_%s",(const char*)calibTask->GetName()));
1041
1042   TH2I *ch2d = (TH2I *) listcalibTask->FindObject("CH2d");
1043   TProfile2D *ph2d = (TProfile2D *) listcalibTask->FindObject("PH2d");
1044   TProfile2D *prf2d = (TProfile2D *) listcalibTask->FindObject("PRF2d");
1045
1046   TH2I *ch2dSum = (TH2I *) listcalibTask->FindObject(Form("CH2dSum_%s",(const char*)calibTask->GetName()));
1047   TProfile2D *ph2dSum = (TProfile2D *) listcalibTask->FindObject(Form("PH2dSum_%s",(const char*)calibTask->GetName()));
1048
1049   TH2I *ch2dSM = (TH2I *) listcalibTask->FindObject(Form("CH2dSM_%s",(const char*)calibTask->GetName()));
1050   TProfile2D *ph2dSM = (TProfile2D *) listcalibTask->FindObject(Form("PH2dSM_%s",(const char*)calibTask->GetName()));
1051
1052   AliTRDCalibraVdriftLinearFit *linearfit = (AliTRDCalibraVdriftLinearFit *) listcalibTask->FindObject("AliTRDCalibraVdriftLinearFit");
1053   AliTRDCalibraExbAltFit *exbaltfit = (AliTRDCalibraExbAltFit *) listcalibTask->FindObject("AliTRDCalibraExbAltFit");
1054   AliTRDCalibraVector *calibraVector = (AliTRDCalibraVector *) listcalibTask->FindObject("AliTRDCalibraVector");
1055
1056   //
1057
1058   THnSparseI *inhistoEntries = (THnSparseI *) fListHist->FindObject("NumberOfEntries");
1059
1060   TH1I *inEventsInput  = (TH1I *) fListHist->FindObject(Form("NEventsInput_%s",(const char*)fName));
1061   TH1I *inEvents  = (TH1I *) fListHist->FindObject(Form("NEvents_%s",(const char*)fName));
1062   TH2F *iabsoluteGain  = (TH2F *) fListHist->FindObject(Form("AbsoluteGain_%s",(const char*)fName));
1063
1064   TH1F *itrdTrack = (TH1F *) fListHist->FindObject(Form("TRDTrack_%s",(const char*)fName));
1065   TH1F *itrdTrackOffline = (TH1F *) fListHist->FindObject(Form("TRDTrackOffline_%s",(const char*)fName));
1066   TH1F *itrdTrackStandalone = (TH1F *) fListHist->FindObject(Form("TRDTrackStandalone_%s",(const char*)fName));
1067
1068   TH2F *itpctrdTrack = (TH2F *) fListHist->FindObject(Form("NbTPCTRDtrack_%s",(const char*)fName));
1069
1070   TH1F *inbTimeBin = (TH1F *) fListHist->FindObject(Form("NbTimeBin_%s",(const char*)fName));
1071   TH1F *inbTimeBinOffline = (TH1F *) fListHist->FindObject(Form("NbTimeBinOffline_%s",(const char*)fName));
1072   TH1F *inbTimeBinStandalone = (TH1F *) fListHist->FindObject(Form("NbTimeBinStandalone_%s",(const char*)fName));
1073
1074   TH1F *inbClusters = (TH1F *) fListHist->FindObject(Form("NbClusters_%s",(const char*)fName));
1075   TH1F *inbClustersOffline = (TH1F *) fListHist->FindObject(Form("NbClustersOffline_%s",(const char*)fName));
1076   TH1F *inbClustersStandalone = (TH1F *) fListHist->FindObject(Form("NbClustersStandalone_%s",(const char*)fName));
1077
1078   TH1F *inbTracklets = (TH1F *) fListHist->FindObject(Form("NbTracklets_%s",(const char*)fName));
1079   TH1F *inbTrackletsOffline = (TH1F *) fListHist->FindObject(Form("NbTrackletsOffline_%s",(const char*)fName));
1080   TH1F *inbTrackletsStandalone = (TH1F *) fListHist->FindObject(Form("NbTrackletsStandalone_%s",(const char*)fName));
1081
1082   TH2I *ich2d = (TH2I *) fListHist->FindObject("CH2d");
1083   TProfile2D *iph2d = (TProfile2D *) fListHist->FindObject("PH2d");
1084   TProfile2D *iprf2d = (TProfile2D *) fListHist->FindObject("PRF2d");
1085
1086   TH2I *ich2dSum = (TH2I *) fListHist->FindObject(Form("CH2dSum_%s",(const char*)fName));
1087   TProfile2D *iph2dSum = (TProfile2D *) fListHist->FindObject(Form("PH2dSum_%s",(const char*)fName));
1088
1089   TH2I *ich2dSM = (TH2I *) fListHist->FindObject(Form("CH2dSM_%s",(const char*)fName));
1090   TProfile2D *iph2dSM = (TProfile2D *) fListHist->FindObject(Form("PH2dSM_%s",(const char*)fName));
1091
1092   AliTRDCalibraVdriftLinearFit *ilinearfit = (AliTRDCalibraVdriftLinearFit *) fListHist->FindObject("AliTRDCalibraVdriftLinearFit");
1093   AliTRDCalibraExbAltFit *iexbaltfit = (AliTRDCalibraExbAltFit *) fListHist->FindObject("AliTRDCalibraExbAltFit");
1094   AliTRDCalibraVector *icalibraVector = (AliTRDCalibraVector *) fListHist->FindObject("AliTRDCalibraVector");
1095
1096
1097   // Add
1098
1099   if(histoEntries) {
1100     if(inhistoEntries) {
1101       inhistoEntries->Add(histoEntries);
1102       //printf("Add Events\n");
1103     }
1104     else {
1105       //printf("Create new Events\n");
1106       inhistoEntries = (THnSparseI *) histoEntries->Clone();
1107       fListHist->Add(inhistoEntries);
1108     }
1109   }
1110
1111   if(nEventsInput) {
1112     if(inEventsInput) {
1113       inEventsInput->Add(nEventsInput);
1114       //printf("Add Events\n");
1115     }
1116     else {
1117       //printf("Create new Events\n");
1118       inEventsInput = new TH1I(*nEventsInput);
1119       fListHist->Add(inEventsInput);
1120     }
1121   }
1122
1123   if(nEvents) {
1124     if(inEvents) {
1125       inEvents->Add(nEvents);
1126       //printf("Add Events\n");
1127     }
1128     else {
1129       //printf("Create new Events\n");
1130       inEvents = new TH1I(*nEvents);
1131       fListHist->Add(inEvents);
1132     }
1133   }
1134
1135   if(absoluteGain) {
1136     if(iabsoluteGain) iabsoluteGain->Add(absoluteGain);
1137     else {
1138       iabsoluteGain = new TH2F(*absoluteGain);
1139       fListHist->Add(iabsoluteGain);
1140     }
1141   }
1142
1143   if(trdTrack) {
1144     if(itrdTrack) itrdTrack->Add(trdTrack);
1145     else {
1146      itrdTrack = new TH1F(*trdTrack);
1147      fListHist->Add(itrdTrack);
1148     }
1149   }
1150
1151   if(trdTrackOffline) {
1152     if(itrdTrackOffline) itrdTrackOffline->Add(trdTrackOffline);
1153     else {
1154       itrdTrackOffline = new TH1F(*trdTrackOffline);
1155       fListHist->Add(itrdTrackOffline);
1156     }
1157   }
1158
1159   if(trdTrackStandalone) {
1160     if(itrdTrackStandalone) itrdTrackStandalone->Add(trdTrackStandalone);
1161     else {
1162       itrdTrackStandalone = new TH1F(*trdTrackStandalone);
1163       fListHist->Add(itrdTrackStandalone);
1164     }
1165   }
1166
1167   if(tpctrdTrack) {
1168     if(itpctrdTrack) itpctrdTrack->Add(tpctrdTrack);
1169     else {
1170       itpctrdTrack = new TH2F(*tpctrdTrack);
1171       fListHist->Add(itpctrdTrack);
1172     }
1173   }
1174
1175   if(nbTimeBin) {
1176     if(inbTimeBin) inbTimeBin->Add(nbTimeBin);
1177     else {
1178       inbTimeBin = new TH1F(*inbTimeBin);
1179       fListHist->Add(inbTimeBin);
1180     }
1181   }
1182
1183   if(nbTimeBinOffline) {
1184     if(inbTimeBinOffline) inbTimeBinOffline->Add(nbTimeBinOffline);
1185     else {
1186       inbTimeBinOffline = new TH1F(*nbTimeBinOffline);
1187       fListHist->Add(inbTimeBinOffline);
1188     }
1189   }
1190
1191   if(nbTimeBinStandalone) {
1192     if(inbTimeBinStandalone) inbTimeBinStandalone->Add(nbTimeBinStandalone);
1193     else {
1194       inbTimeBinStandalone = new TH1F(*nbTimeBinStandalone);
1195       fListHist->Add(inbTimeBinStandalone);
1196     }
1197   }
1198
1199   if(nbClusters) {
1200     if(inbClusters) inbClusters->Add(nbClusters);
1201     else {
1202       inbClusters = new TH1F(*nbClusters);
1203       fListHist->Add(inbClusters);
1204     }
1205   }
1206
1207   if(nbClustersOffline) {
1208     if(inbClustersOffline) inbClustersOffline->Add(nbClustersOffline);
1209     else {
1210       inbClustersOffline = new TH1F(*nbClustersOffline);
1211       fListHist->Add(inbClustersOffline);
1212     }
1213   }
1214
1215   if(nbClustersStandalone) {
1216     if(inbClustersStandalone) inbClustersStandalone->Add(nbClustersStandalone);
1217     else {
1218       inbClustersStandalone = new TH1F(*nbClustersStandalone);
1219       fListHist->Add(inbClustersStandalone);
1220     }
1221   }
1222
1223   if(nbTracklets) {
1224     if(inbTracklets) inbTracklets->Add(nbTracklets);
1225     else {
1226       inbTracklets = new TH1F(*nbTracklets);
1227       fListHist->Add(inbTracklets);
1228     }
1229   }
1230
1231   if(nbTrackletsOffline) {
1232     if(inbTrackletsOffline) inbTrackletsOffline->Add(nbTrackletsOffline);
1233     else {
1234       inbTrackletsOffline = new TH1F(*nbTrackletsOffline);
1235       fListHist->Add(inbTrackletsOffline);
1236     }
1237   }
1238
1239   if(nbTrackletsStandalone) {
1240     if(inbTrackletsStandalone) inbTrackletsStandalone->Add(nbTrackletsStandalone);
1241     else {
1242       inbTrackletsStandalone = new TH1F(*nbTrackletsStandalone);
1243       fListHist->Add(inbTrackletsStandalone);
1244     }
1245   }
1246
1247   if(ch2d) {
1248     if(ich2d) ich2d->Add(ch2d);
1249     else {
1250       ich2d = new TH2I(*ch2d);
1251       fListHist->Add(ich2d);
1252     }
1253   }
1254
1255   if(ph2d) {
1256     if(iph2d) iph2d->Add(ph2d);
1257     else {
1258       iph2d = new TProfile2D(*ph2d);
1259       fListHist->Add(iph2d);
1260     }
1261   }
1262
1263   if(prf2d) {
1264     if(iprf2d) iprf2d->Add(prf2d);
1265     else {
1266       iprf2d = new TProfile2D(*prf2d);
1267       fListHist->Add(iprf2d);
1268     }
1269   }
1270
1271   if(ch2dSum) {
1272     if(ich2dSum) ich2dSum->Add(ch2dSum);
1273     else {
1274       ich2dSum = new TH2I(*ch2dSum);
1275       fListHist->Add(ich2dSum);
1276     }
1277   }
1278
1279   if(ph2dSum) {
1280     if(iph2dSum) iph2dSum->Add(ph2dSum);
1281     else {
1282       iph2dSum = new TProfile2D(*ph2dSum);
1283       fListHist->Add(iph2dSum);
1284     }
1285   }
1286
1287   if(ch2dSM) {
1288     if(ich2dSM) ich2dSM->Add(ch2dSM);
1289     else {
1290       ich2dSM = new TH2I(*ch2dSM);
1291       fListHist->Add(ich2dSM);
1292     }
1293   }
1294
1295   if(ph2dSM) {
1296     if(iph2dSM) iph2dSM->Add(ph2dSM);
1297     else {
1298       iph2dSM = new TProfile2D(*ph2dSM);
1299       fListHist->Add(iph2dSM);
1300     }
1301   }
1302
1303   if(linearfit) {
1304     if(ilinearfit) ilinearfit->Add(linearfit);
1305     else {
1306       ilinearfit = new AliTRDCalibraVdriftLinearFit(*linearfit);
1307       fListHist->Add(ilinearfit);
1308     }
1309   }
1310
1311   if(exbaltfit) {
1312     if(iexbaltfit) iexbaltfit->Add(exbaltfit);
1313     else {
1314       iexbaltfit = new AliTRDCalibraExbAltFit(*exbaltfit);
1315       fListHist->Add(iexbaltfit);
1316     }
1317   }
1318
1319   if(calibraVector) {
1320     if(icalibraVector) icalibraVector->Add(calibraVector);
1321     else {
1322       icalibraVector = new AliTRDCalibraVector(*calibraVector);
1323       fListHist->Add(icalibraVector);
1324     }
1325   }
1326
1327 }
1328 //________________________________________________________________________________
1329 Long64_t AliTRDCalibTask::Merge(TCollection *li) {
1330
1331   //
1332   // merge component
1333   //
1334
1335   TIterator* iter = li->MakeIterator();
1336   AliTRDCalibTask* cal = 0;
1337
1338   while ((cal = (AliTRDCalibTask*)iter->Next())) {
1339     if (!cal->InheritsFrom(AliTRDCalibTask::Class())) {
1340       Error("Merge","Attempt to add object of class %s to a %s", cal->ClassName(), this->ClassName());
1341       return -1;
1342     }
1343
1344     // add histograms here...
1345     this->AddTask(cal);
1346
1347   }
1348
1349   return 0;
1350
1351 }
1352 //_____________________________________________________
1353 Bool_t AliTRDCalibTask::SetVersionSubversion(){
1354   //
1355   // Load Chamber Gain factors into the Tender supply
1356   //
1357
1358   //printf("SetVersionSubversion\n");
1359
1360   //find previous entry from the UserInfo
1361   TTree *tree=((TChain*)GetInputData(0))->GetTree();
1362   if (!tree) {
1363     AliError("Tree not found in ESDhandler");
1364     return kFALSE;
1365   }
1366
1367   TList *userInfo=(TList*)tree->GetUserInfo();
1368   if (!userInfo) {
1369     AliError("No UserInfo found in tree");
1370     return kFALSE;
1371   }
1372
1373   TList *cdbList=(TList*)userInfo->FindObject("cdbList");
1374   if (!cdbList) {
1375     AliError("No cdbList found in UserInfo");
1376     if (AliLog::GetGlobalLogLevel()>=AliLog::kError) userInfo->Print();
1377     return kFALSE;
1378   }
1379
1380   TIter nextCDB(cdbList);
1381   TObjString *os=0x0;
1382   while ( (os=(TObjString*)nextCDB()) ){
1383     if(os->GetString().Contains("TRD/Calib/ChamberGainFactor")){
1384       // Get Old gain calibration
1385       AliCDBId *id=AliCDBId::MakeFromString(os->GetString());
1386       fFirstRunGain = id->GetFirstRun();
1387       fVersionGainUsed = id->GetVersion();
1388       fSubVersionGainUsed = id->GetSubVersion();
1389     } else if(os->GetString().Contains("TRD/Calib/ChamberVdrift")){
1390       // Get Old drift velocity calibration
1391       AliCDBId *id=AliCDBId::MakeFromString(os->GetString());
1392       fFirstRunVdrift = id->GetFirstRun();
1393       fVersionVdriftUsed = id->GetVersion();
1394       fSubVersionVdriftUsed = id->GetSubVersion();
1395     } else if(os->GetString().Contains("TRD/Calib/LocalGainFactor")){
1396       // Get Old drift velocity calibration
1397       AliCDBId *id=AliCDBId::MakeFromString(os->GetString());
1398       fFirstRunGainLocal = id->GetFirstRun();
1399       fVersionGainLocalUsed = id->GetVersion();
1400       fSubVersionGainLocalUsed = id->GetSubVersion();
1401     } else if((os->GetString().Contains("TRD/Calib/ChamberExB")) && (!os->GetString().Contains("TRD/Calib/ChamberExBAlt"))){
1402       // Get Old drift velocity calibration
1403       AliCDBId *id=AliCDBId::MakeFromString(os->GetString());
1404       fFirstRunExB = id->GetFirstRun();
1405       fVersionExBUsed = id->GetVersion();
1406       fSubVersionExBUsed = id->GetSubVersion();
1407       //printf("Version %d and subversion %d\n",fVersionExBUsed,fSubVersionExBUsed);
1408     }
1409   }
1410
1411   //printf("VersionGain %d, SubversionGain %d, VersionLocalGain %d, Subversionlocalgain %d, Versionvdrift %d, Subversionvdrift %d\n",fVersionGainUsed,fSubVersionGainUsed,fVersionGainLocalUsed,fSubVersionGainLocalUsed,fVersionVdriftUsed,fSubVersionVdriftUsed);
1412
1413   // Check
1414   if((fFirstRunGain < 0)            ||
1415      (fFirstRunGainLocal < 0)       ||
1416      (fFirstRunVdrift < 0)          ||
1417      (fVersionGainUsed < 0)         ||
1418      (fVersionGainLocalUsed < 0)    ||
1419      (fSubVersionGainUsed < 0)      ||
1420      (fSubVersionGainLocalUsed < 0) ||
1421      (fVersionVdriftUsed < 0)       ||
1422      (fSubVersionVdriftUsed < 0)) {
1423     AliError("No recent calibration found");
1424     return kFALSE;
1425   }
1426   else return kTRUE;
1427
1428 }
1429 //_________________________________________________________________________________________________________________________
1430 Bool_t AliTRDCalibTask::ParticleGood(int i) const {
1431
1432   //
1433   // Definition of good tracks
1434   //
1435
1436
1437   AliESDtrack *track = fESD->GetTrack(i);
1438   if (!track->IsOn(AliESDtrack::kTPCrefit)) return 0;        // TPC refit
1439   if (track->GetTPCNcls() < 90) return 0;                    // number of TPC clusters
1440   if (fabs(track->Eta())>0.8) return 0;                         // fiducial pseudorapidity
1441   Float_t r,z;
1442   track->GetImpactParametersTPC(r,z);
1443   if (fabs(z)>2.0) return 0;                          // impact parameter in z
1444   if (fabs(r)>2.0) return 0;                          // impact parameter in xy
1445   if (r==0) return 0;
1446   return 1;
1447
1448
1449 }
1450 //______________________________________________________________________________________________________________________
1451 Bool_t AliTRDCalibTask::FindP1TrackPHtrackletV1Test(const AliTRDseedV1 *tracklet, Int_t nbclusters)
1452 {
1453   //
1454   // Drift velocity calibration:
1455   // Fit the clusters with a straight line
1456   // From the slope find the drift velocity
1457   //
1458
1459   ////////////////////////////////////////////////
1460   //Number of points: if less than 3 return kFALSE
1461   /////////////////////////////////////////////////
1462   if(nbclusters <= 2) return kFALSE;
1463
1464   ////////////
1465   //Variables
1466   ////////////
1467   // results of the linear fit
1468   Double_t dydt                       = 0.0;                                // dydt tracklet after straight line fit
1469   Double_t errorpar                   = 0.0;                                // error after straight line fit on dy/dt
1470   Double_t pointError                 = 0.0;                                // error after straight line fit
1471   // pad row problemes: avoid tracklet that cross pad rows, tilting angle in the constant
1472   Int_t    crossrow                   = 0;                                  // if it crosses a pad row
1473   Int_t    rowp                       = -1;                                 // if it crosses a pad row
1474   Float_t  tnt                        = tracklet->GetTilt();                // tan tiltingangle
1475   TLinearFitter linearFitterTracklet(2,"pol1");
1476   linearFitterTracklet.StoreData(kTRUE);
1477
1478
1479   ///////////////////////////////////////////
1480   // Take the parameters of the track
1481   //////////////////////////////////////////
1482   // take now the snp, tnp and tgl from the track
1483   Double_t snp = tracklet->GetSnp();             // sin dy/dx at the end of the chamber
1484   Double_t tnp = 0.0;                            // dy/dx at the end of the chamber
1485   if( TMath::Abs(snp) <  1.){
1486     tnp = snp / TMath::Sqrt((1.-snp)*(1.+snp));
1487   }
1488   Double_t tgl  = tracklet->GetTgl();           // dz/dl
1489   Double_t dzdx = tgl*TMath::Sqrt(1+tnp*tnp);   // dz/dx calculated from dz/dl
1490   // at the entrance
1491   //Double_t tnp = tracklet->GetYref(1);      // dy/dx at the entrance of the chamber
1492   //Double_t tgl = tracklet->GetZref(1);      // dz/dl at the entrance of the chamber
1493   //Double_t dzdx = tgl;                      //*TMath::Sqrt(1+tnp*tnp); // dz/dx from dz/dl
1494   // at the end with correction due to linear fit
1495   //Double_t tnp = tracklet->GetYfit(1);      // dy/dx at the end of the chamber after fit correction
1496   //Double_t tgl = tracklet->GetZfit(1);      // dz/dl at the end of the chamber after fit correction
1497
1498
1499   ////////////////////////////
1500   // loop over the clusters
1501   ////////////////////////////
1502   Int_t  nbli = 0;
1503   AliTRDcluster *cl                   = 0x0;
1504   //////////////////////////////
1505   // Check no shared clusters
1506   //////////////////////////////
1507   for(int icc=AliTRDseedV1::kNtb; icc<AliTRDseedV1::kNclusters; icc++){
1508     cl = tracklet->GetClusters(icc);
1509     if(cl)  crossrow = 1;
1510   }
1511   //////////////////////////////////
1512   // Loop clusters
1513   //////////////////////////////////
1514   for(int ic=0; ic<AliTRDseedV1::kNtb; ic++){
1515     if(!(cl = tracklet->GetClusters(ic))) continue;
1516     //if((fLimitChargeIntegration) && (!cl->IsInChamber())) continue;
1517
1518     Double_t ycluster                 = cl->GetY();
1519     Int_t time                        = cl->GetPadTime();
1520     Double_t timeis                   = time/10.0;
1521     //See if cross two pad rows
1522     Int_t    row                      = cl->GetPadRow();
1523     if(rowp==-1) rowp                 = row;
1524     if(row != rowp) crossrow          = 1;
1525
1526     linearFitterTracklet.AddPoint(&timeis,ycluster,1);
1527     nbli++;
1528
1529
1530   }
1531
1532   ////////////////////////////////////
1533   // Do the straight line fit now
1534   ///////////////////////////////////
1535   if(nbli <= 2){
1536     linearFitterTracklet.ClearPoints();
1537     return kFALSE;
1538   }
1539   TVectorD pars;
1540   linearFitterTracklet.Eval();
1541   linearFitterTracklet.GetParameters(pars);
1542   pointError  =  TMath::Sqrt(linearFitterTracklet.GetChisquare()/(nbli-2));
1543   errorpar    =  linearFitterTracklet.GetParError(1)*pointError;
1544   dydt        = pars[1];
1545   //printf("chis %f, nbli %d, pointError %f, parError %f, errorpar %f\n",linearFitterTracklet->GetChisquare(),nbli,pointError,linearFitterTracklet->GetParError(1),errorpar);
1546   linearFitterTracklet.ClearPoints();
1547
1548   /////////////////////////
1549   // Cuts quality
1550   ////////////////////////
1551
1552   if(nbclusters < fLow) return kFALSE;
1553   if(nbclusters > fHigh) return kFALSE;
1554   if(pointError >= 0.3) return kFALSE;
1555   if(crossrow == 1) return kTRUE;
1556
1557   ///////////////////////
1558   // Fill
1559   //////////////////////
1560
1561   if(fDebug > 0){
1562     //Add to the linear fitter of the detector
1563     if( TMath::Abs(snp) <  1.){
1564       Double_t x = tnp-dzdx*tnt;
1565       //if(!fLinearVdriftTest) printf("Not there\n");
1566       Double_t nbentries = fLinearVdriftTest->GetEntries();
1567       if(nbentries < (5.0*32767)) fLinearVdriftTest->Fill(x,dydt);
1568     }
1569   }
1570
1571   return kTRUE;
1572 }
1573