TRD/AliTRDCalibraFillHisto.cxx

   1 //**************************************************************************
   2 // * Copyright(c) 1998-1999, ALICE Experiment at CERN, All rights reserved. *
   3 // *                                                                        *
   4 // * Author: The ALICE Off-line Project.                                    *
   5 // * Contributors are mentioned in the code where appropriate.              *
   6 //   *                                                                        *
   7 //   * Permission to use, copy, modify and distribute this software and its   *
   8 //  * documentation strictly for non-commercial purposes is hereby granted   *
   9 //    * without fee, provided that the above copyright notice appears in all   *
  10 //    * copies and that both the copyright notice and this permission notice   *
  11 //* appear in the supporting documentation. The authors make no claims     *
  12 //* about the suitability of this software for any purpose. It is          *
  13 //* provided "as is" without express or implied warranty.                  *
  14 //**************************************************************************/
  15
  16 /* $Id$ */
  17
  18 /////////////////////////////////////////////////////////////////////////////////
  19 //
  20 // AliTRDCalibraFillHisto
  21 //
  22 // This class is for the TRD calibration of the relative gain factor, the drift velocity,
  23 // the time 0 and the pad response function. It fills histos or vectors.
  24 // It can be used for the calibration per chamber but also per group of pads and eventually per pad.
  25 // The user has to choose with the functions SetNz and SetNrphi the precision of the calibration (see AliTRDCalibraMode).
  26 // 2D Histograms (Histo2d) or vectors (Vector2d), then converted in Trees, will be filled
  27 // from RAW DATA in a run or from reconstructed TRD tracks during the offline tracking
  28 // in the function "FollowBackProlongation" (AliTRDtracker)
  29 // Per default the functions to fill are off.
  30 //
  31 // Authors:
  32 //   R. Bailhache (R.Bailhache@gsi.de, rbailhache@ikf.uni-frankfurt.de)
  33 //   J. Book (jbook@ikf.uni-frankfurt.de)
  34 //
  35 //////////////////////////////////////////////////////////////////////////////////////
  36
  37 #include <TProfile2D.h>
  38 #include <TProfile.h>
  39 #include <TFile.h>
  40 #include <TStyle.h>
  41 #include <TCanvas.h>
  42 #include <TObjArray.h>
  43 #include <TObject.h>
  44 #include <TH1F.h>
  45 #include <TH2I.h>
  46 #include <TH2.h>
  47 #include <TStopwatch.h>
  48 #include <TMath.h>
  49 #include <TDirectory.h>
  50 #include <TTreeStream.h>
  51 #include <TVectorD.h>
  52 #include <TLinearFitter.h>
  53
  54 #include "AliLog.h"
  55
  56 #include "AliTRDCalibraFillHisto.h"
  57 #include "AliTRDCalibraMode.h"
  58 #include "AliTRDCalibraVector.h"
  59 #include "AliTRDCalibraVdriftLinearFit.h"
  60 #include "AliTRDcalibDB.h"
  61 #include "AliTRDCommonParam.h"
  62 #include "AliTRDpadPlane.h"
  63 #include "AliTRDcluster.h"
  64 #include "AliTRDtrack.h"
  65 #include "AliTRDtrackV1.h"
  66 #include "AliTRDrawStreamBase.h"
  67 #include "AliRawReader.h"
  68 #include "AliRawReaderDate.h"
  69 #include "AliTRDgeometry.h"
  70 #include "./Cal/AliTRDCalROC.h"
  71 #include "./Cal/AliTRDCalPad.h"
  72 #include "./Cal/AliTRDCalDet.h"
  73
  74 #include "AliTRDrawFastStream.h"
  75 #include "AliTRDdigitsManager.h"
  76 #include "AliTRDdigitsParam.h"
  77 #include "AliTRDSignalIndex.h"
  78 #include "AliTRDarrayADC.h"
  79
  80 #include "AliTRDrawStream.h"
  81
  82 #include "AliCDBEntry.h"
  83 #include "AliCDBManager.h"
  84
  85 #ifdef ALI_DATE
  86 #include "event.h"
  87 #endif
  88
  89
  90 ClassImp(AliTRDCalibraFillHisto)
  91
  92 AliTRDCalibraFillHisto* AliTRDCalibraFillHisto::fgInstance = 0;
  93 Bool_t AliTRDCalibraFillHisto::fgTerminated = kFALSE;
  94
  95 //_____________singleton implementation_________________________________________________
  96 AliTRDCalibraFillHisto *AliTRDCalibraFillHisto::Instance()
  97 {
  98   //
  99   // Singleton implementation
 100   //
 101
 102   if (fgTerminated != kFALSE) {
 103     return 0;
 104   }
 105
 106   if (fgInstance == 0) {
 107     fgInstance = new AliTRDCalibraFillHisto();
 108   }
 109
 110   return fgInstance;
 111
 112 }
 113
 114 //______________________________________________________________________________________
 115 void AliTRDCalibraFillHisto::Terminate()
 116 {
 117   //
 118   // Singleton implementation
 119   // Deletes the instance of this class
 120   //
 121
 122   fgTerminated = kTRUE;
 123
 124   if (fgInstance != 0) {
 125     delete fgInstance;
 126     fgInstance = 0;
 127   }
 128
 129 }
 130
 131 //______________________________________________________________________________________
 132 AliTRDCalibraFillHisto::AliTRDCalibraFillHisto()
 133   :TObject()
 134   ,fGeo(0)
 135   ,fCalibDB(0)
 136   ,fIsHLT(kFALSE)
 137   ,fCH2dOn(kFALSE)
 138   ,fPH2dOn(kFALSE)
 139   ,fPRF2dOn(kFALSE)
 140   ,fHisto2d(kFALSE)
 141   ,fVector2d(kFALSE)
 142   ,fLinearFitterOn(kFALSE)
 143   ,fLinearFitterDebugOn(kFALSE)
 144   ,fRelativeScale(0)
 145   ,fThresholdClusterPRF2(15.0)
 146   ,fLimitChargeIntegration(kFALSE)
 147   ,fFillWithZero(kFALSE)
 148   ,fNormalizeNbOfCluster(kFALSE)
 149   ,fMaxCluster(0)
 150   ,fNbMaxCluster(0)
 151   ,fVersionGainUsed(0)
 152   ,fSubVersionGainUsed(0)
 153   ,fVersionGainLocalUsed(0)
 154   ,fSubVersionGainLocalUsed(0)
 155   ,fVersionVdriftUsed(0)
 156   ,fSubVersionVdriftUsed(0)
 157   ,fCalibraMode(new AliTRDCalibraMode())
 158   ,fDebugStreamer(0)
 159   ,fDebugLevel(0)
 160   ,fDetectorPreviousTrack(-1)
 161   ,fMCMPrevious(-1)
 162   ,fROBPrevious(-1)
 163   ,fNumberClusters(1)
 164   ,fNumberClustersf(30)
 165   ,fNumberClustersProcent(0.5)
 166   ,fThresholdClustersDAQ(120.0)
 167   ,fNumberRowDAQ(2)
 168   ,fNumberColDAQ(4)
 169   ,fProcent(6.0)
 170   ,fDifference(17)
 171   ,fNumberTrack(0)
 172   ,fTimeMax(0)
 173   ,fSf(10.0)
 174   ,fNumberBinCharge(50)
 175   ,fNumberBinPRF(10)
 176   ,fNgroupprf(3)
 177   ,fAmpTotal(0x0)
 178   ,fPHPlace(0x0)
 179   ,fPHValue(0x0)
 180   ,fGoodTracklet(kTRUE)
 181   ,fLinearFitterTracklet(0x0)
 182   ,fEntriesCH(0x0)
 183   ,fEntriesLinearFitter(0x0)
 184   ,fCalibraVector(0x0)
 185   ,fPH2d(0x0)
 186   ,fPRF2d(0x0)
 187   ,fCH2d(0x0)
 188   ,fLinearFitterArray(540)
 189   ,fLinearVdriftFit(0x0)
 190   ,fCalDetGain(0x0)
 191   ,fCalROCGain(0x0)
 192 {
 193   //
 194   // Default constructor
 195   //
 196
 197   //
 198   // Init some default values
 199   //
 200
 201   fNumberUsedCh[0]       = 0;
 202   fNumberUsedCh[1]       = 0;
 203   fNumberUsedPh[0]       = 0;
 204   fNumberUsedPh[1]       = 0;
 205
 206   fGeo = new AliTRDgeometry();
 207   fCalibDB = AliTRDcalibDB::Instance();
 208 }
 209
 210 //______________________________________________________________________________________
 211 AliTRDCalibraFillHisto::AliTRDCalibraFillHisto(const AliTRDCalibraFillHisto &c)
 212   :TObject(c)
 213   ,fGeo(0)
 214   ,fCalibDB(0)
 215   ,fIsHLT(c.fIsHLT)
 216   ,fCH2dOn(c.fCH2dOn)
 217   ,fPH2dOn(c.fPH2dOn)
 218   ,fPRF2dOn(c.fPRF2dOn)
 219   ,fHisto2d(c.fHisto2d)
 220   ,fVector2d(c.fVector2d)
 221   ,fLinearFitterOn(c.fLinearFitterOn)
 222   ,fLinearFitterDebugOn(c.fLinearFitterDebugOn)
 223   ,fRelativeScale(c.fRelativeScale)
 224   ,fThresholdClusterPRF2(c.fThresholdClusterPRF2)
 225   ,fLimitChargeIntegration(c.fLimitChargeIntegration)
 226   ,fFillWithZero(c.fFillWithZero)
 227   ,fNormalizeNbOfCluster(c.fNormalizeNbOfCluster)
 228   ,fMaxCluster(c.fMaxCluster)
 229   ,fNbMaxCluster(c.fNbMaxCluster)
 230   ,fVersionGainUsed(c.fVersionGainUsed)
 231   ,fSubVersionGainUsed(c.fSubVersionGainUsed)
 232   ,fVersionGainLocalUsed(c.fVersionGainLocalUsed)
 233   ,fSubVersionGainLocalUsed(c.fSubVersionGainLocalUsed)
 234   ,fVersionVdriftUsed(c.fVersionVdriftUsed)
 235   ,fSubVersionVdriftUsed(c.fSubVersionVdriftUsed)
 236   ,fCalibraMode(0x0)
 237   ,fDebugStreamer(0)
 238   ,fDebugLevel(c.fDebugLevel)
 239   ,fDetectorPreviousTrack(c.fDetectorPreviousTrack)
 240   ,fMCMPrevious(c.fMCMPrevious)
 241   ,fROBPrevious(c.fROBPrevious)
 242   ,fNumberClusters(c.fNumberClusters)
 243   ,fNumberClustersf(c.fNumberClustersf)
 244   ,fNumberClustersProcent(c.fNumberClustersProcent)
 245   ,fThresholdClustersDAQ(c.fThresholdClustersDAQ)
 246   ,fNumberRowDAQ(c.fNumberRowDAQ)
 247   ,fNumberColDAQ(c.fNumberColDAQ)
 248   ,fProcent(c.fProcent)
 249   ,fDifference(c.fDifference)
 250   ,fNumberTrack(c.fNumberTrack)
 251   ,fTimeMax(c.fTimeMax)
 252   ,fSf(c.fSf)
 253   ,fNumberBinCharge(c.fNumberBinCharge)
 254   ,fNumberBinPRF(c.fNumberBinPRF)
 255   ,fNgroupprf(c.fNgroupprf)
 256   ,fAmpTotal(0x0)
 257   ,fPHPlace(0x0)
 258   ,fPHValue(0x0)
 259   ,fGoodTracklet(c.fGoodTracklet)
 260   ,fLinearFitterTracklet(0x0)
 261   ,fEntriesCH(0x0)
 262   ,fEntriesLinearFitter(0x0)
 263   ,fCalibraVector(0x0)
 264   ,fPH2d(0x0)
 265   ,fPRF2d(0x0)
 266   ,fCH2d(0x0)
 267   ,fLinearFitterArray(540)
 268   ,fLinearVdriftFit(0x0)
 269   ,fCalDetGain(0x0)
 270   ,fCalROCGain(0x0)
 271 {
 272   //
 273   // Copy constructor
 274   //
 275   if(c.fCalibraMode)   fCalibraMode = new AliTRDCalibraMode(*c.fCalibraMode);
 276   if(c.fCalibraVector) fCalibraVector = new AliTRDCalibraVector(*c.fCalibraVector);
 277   if(c.fPH2d) {
 278     fPH2d = new TProfile2D(*c.fPH2d);
 279     fPH2d->SetDirectory(0);
 280   }
 281   if(c.fPRF2d) {
 282     fPRF2d = new TProfile2D(*c.fPRF2d);
 283     fPRF2d->SetDirectory(0);
 284   }
 285   if(c.fCH2d) {
 286     fCH2d = new TH2I(*c.fCH2d);
 287     fCH2d->SetDirectory(0);
 288   }
 289   if(c.fLinearVdriftFit){
 290     fLinearVdriftFit = new AliTRDCalibraVdriftLinearFit(*c.fLinearVdriftFit);
 291   }
 292
 293   if(c.fCalDetGain)  fCalDetGain   = new AliTRDCalDet(*c.fCalDetGain);
 294   if(c.fCalROCGain)  fCalROCGain   = new AliTRDCalROC(*c.fCalROCGain);
 295
 296   if (fGeo) {
 297     delete fGeo;
 298   }
 299   fGeo = new AliTRDgeometry();
 300   fCalibDB = AliTRDcalibDB::Instance();
 301
 302   fNumberUsedCh[0]       = 0;
 303   fNumberUsedCh[1]       = 0;
 304   fNumberUsedPh[0]       = 0;
 305   fNumberUsedPh[1]       = 0;
 306
 307 }
 308
 309 //____________________________________________________________________________________
 310 AliTRDCalibraFillHisto::~AliTRDCalibraFillHisto()
 311 {
 312   //
 313   // AliTRDCalibraFillHisto destructor
 314   //
 315
 316   ClearHistos();
 317   if ( fDebugStreamer ) delete fDebugStreamer;
 318
 319   if ( fCalDetGain )  delete fCalDetGain;
 320   if ( fCalROCGain )  delete fCalROCGain;
 321
 322   if( fLinearFitterTracklet ) { delete fLinearFitterTracklet; }
 323
 324   delete [] fPHPlace;
 325   delete [] fPHValue;
 326   delete [] fEntriesCH;
 327   delete [] fEntriesLinearFitter;
 328   delete [] fAmpTotal;
 329
 330   for(Int_t idet=0; idet<AliTRDgeometry::kNdet; idet++){
 331     TLinearFitter *f = (TLinearFitter*)fLinearFitterArray.At(idet);
 332     if(f) { delete f;}
 333   }
 334   if(fLinearVdriftFit) delete fLinearVdriftFit;
 335   if (fGeo) {
 336     delete fGeo;
 337   }
 338
 339 }
 340 //_____________________________________________________________________________
 341 void AliTRDCalibraFillHisto::Destroy()
 342 {
 343   //
 344   // Delete instance
 345   //
 346
 347   if (fgInstance) {
 348     delete fgInstance;
 349     fgInstance = 0x0;
 350   }
 351 }
 352 //_____________________________________________________________________________
 353 void AliTRDCalibraFillHisto::DestroyDebugStreamer()
 354 {
 355   //
 356   // Delete DebugStreamer
 357   //
 358
 359   if ( fDebugStreamer ) delete fDebugStreamer;
 360   fDebugStreamer = 0x0;
 361
 362 }
 363 //_____________________________________________________________________________
 364 void AliTRDCalibraFillHisto::ClearHistos()
 365 {
 366   //
 367   // Delete the histos
 368   //
 369
 370   if (fPH2d) {
 371     delete fPH2d;
 372     fPH2d  = 0x0;
 373   }
 374   if (fCH2d) {
 375     delete fCH2d;
 376     fCH2d  = 0x0;
 377   }
 378   if (fPRF2d) {
 379     delete fPRF2d;
 380     fPRF2d = 0x0;
 381   }
 382
 383 }
 384 //////////////////////////////////////////////////////////////////////////////////
 385 // calibration with AliTRDtrackV1: Init, Update
 386 //////////////////////////////////////////////////////////////////////////////////
 387 //____________Functions for initialising the AliTRDCalibraFillHisto in the code_________
 388 Bool_t AliTRDCalibraFillHisto::Init2Dhistos(Int_t nboftimebin)
 389 {
 390   //
 391   // Init the histograms and stuff to be filled
 392   //
 393
 394   // DB Setting
 395   // Get cal
 396   AliTRDcalibDB *cal = AliTRDcalibDB::Instance();
 397   if (!cal) {
 398     AliInfo("Could not get calibDB");
 399     return kFALSE;
 400   }
 401   AliTRDCommonParam *parCom = AliTRDCommonParam::Instance();
 402   if (!parCom) {
 403     AliInfo("Could not get CommonParam");
 404     return kFALSE;
 405   }
 406
 407   // Some parameters
 408   if(nboftimebin > 0) fTimeMax = nboftimebin;
 409   else fTimeMax = cal->GetNumberOfTimeBinsDCS();
 410   if(fTimeMax <= 0) fTimeMax = 30;
 411   printf("////////////////////////////////////////////\n");
 412   printf("Number of time bins in calibration component %d\n",fTimeMax);
 413   printf("////////////////////////////////////////////\n");
 414   fSf                 = parCom->GetSamplingFrequency();
 415   if(!fNormalizeNbOfCluster) fRelativeScale = 20.0;
 416   else fRelativeScale = 1.18;
 417   fNumberClustersf    = fTimeMax;
 418   fNumberClusters     = (Int_t)(fNumberClustersProcent*fTimeMax);
 419
 420   // Init linear fitter
 421   if(!fLinearFitterTracklet) {
 422     fLinearFitterTracklet = new TLinearFitter(2,"pol1");
 423     fLinearFitterTracklet->StoreData(kTRUE);
 424   }
 425
 426   // Calcul Xbins Chambd0, Chamb2
 427   Int_t ntotal0 = CalculateTotalNumberOfBins(0);
 428   Int_t ntotal1 = CalculateTotalNumberOfBins(1);
 429   Int_t ntotal2 = CalculateTotalNumberOfBins(2);
 430
 431   // If vector method On initialised all the stuff
 432   if(fVector2d){
 433     fCalibraVector = new AliTRDCalibraVector();
 434     fCalibraVector->SetNumberBinCharge(fNumberBinCharge);
 435     fCalibraVector->SetTimeMax(fTimeMax);
 436     if(fNgroupprf != 0) {
 437       fCalibraVector->SetNumberBinPRF(2*fNgroupprf*fNumberBinPRF);
 438       fCalibraVector->SetPRFRange((Float_t)(3.0*fNgroupprf));
 439     }
 440     else {
 441       fCalibraVector->SetNumberBinPRF(fNumberBinPRF);
 442       fCalibraVector->SetPRFRange(1.5);
 443     }
 444     for(Int_t k = 0; k < 3; k++){
 445       fCalibraVector->SetDetCha0(k,fCalibraMode->GetDetChamb0(k));
 446       fCalibraVector->SetDetCha2(k,fCalibraMode->GetDetChamb2(k));
 447     }
 448     fCalibraVector->SetNzNrphi(0,fCalibraMode->GetNz(0),fCalibraMode->GetNrphi(0));
 449     fCalibraVector->SetNzNrphi(1,fCalibraMode->GetNz(1),fCalibraMode->GetNrphi(1));
 450     fCalibraVector->SetNzNrphi(2,fCalibraMode->GetNz(2),fCalibraMode->GetNrphi(2));
 451     fCalibraVector->SetNbGroupPRF(fNgroupprf);
 452   }
 453
 454   // Create the 2D histos corresponding to the pad groupCalibration mode
 455   if (fCH2dOn) {
 456
 457     AliInfo(Form("The pad calibration mode for the relative gain calibration: Nz %d, and Nrphi %d"
 458                 ,fCalibraMode->GetNz(0)
 459                 ,fCalibraMode->GetNrphi(0)));
 460
 461     // Create the 2D histo
 462     if (fHisto2d) {
 463       CreateCH2d(ntotal0);
 464     }
 465     // Variable
 466     fAmpTotal = new Float_t[TMath::Max(fCalibraMode->GetDetChamb2(0),fCalibraMode->GetDetChamb0(0))];
 467     for (Int_t k = 0; k < TMath::Max(fCalibraMode->GetDetChamb2(0),fCalibraMode->GetDetChamb0(0)); k++) {
 468       fAmpTotal[k] = 0.0;
 469     }
 470     //Statistics
 471     fEntriesCH = new Int_t[ntotal0];
 472     for(Int_t k = 0; k < ntotal0; k++){
 473       fEntriesCH[k] = 0;
 474     }
 475
 476   }
 477   if (fPH2dOn) {
 478
 479     AliInfo(Form("The pad calibration mode for the drift velocity calibration: Nz %d, and Nrphi %d"
 480                 ,fCalibraMode->GetNz(1)
 481                 ,fCalibraMode->GetNrphi(1)));
 482
 483     // Create the 2D histo
 484     if (fHisto2d) {
 485       CreatePH2d(ntotal1);
 486     }
 487     // Variable
 488     fPHPlace = new Short_t[fTimeMax];
 489     for (Int_t k = 0; k < fTimeMax; k++) {
 490       fPHPlace[k] = -1;
 491     }
 492     fPHValue = new Float_t[fTimeMax];
 493     for (Int_t k = 0; k < fTimeMax; k++) {
 494       fPHValue[k] = 0.0;
 495     }
 496   }
 497   if (fLinearFitterOn) {
 498     if(fLinearFitterDebugOn) {
 499       fLinearFitterArray.SetName("ArrayLinearFitters");
 500       fEntriesLinearFitter = new Int_t[540];
 501       for(Int_t k = 0; k < 540; k++){
 502         fEntriesLinearFitter[k] = 0;
 503       }
 504     }
 505     fLinearVdriftFit = new AliTRDCalibraVdriftLinearFit();
 506   }
 507
 508   if (fPRF2dOn) {
 509
 510     AliInfo(Form("The pad calibration mode for the PRF calibration: Nz %d, and Nrphi %d"
 511                 ,fCalibraMode->GetNz(2)
 512                 ,fCalibraMode->GetNrphi(2)));
 513     // Create the 2D histo
 514     if (fHisto2d) {
 515       CreatePRF2d(ntotal2);
 516     }
 517   }
 518
 519   return kTRUE;
 520
 521 }
 522 ///////////////////////////////////////////////////////////////////////////////////////////////////////////////
 523 Bool_t AliTRDCalibraFillHisto::InitCalDet()
 524 {
 525   //
 526   // Init the Gain Cal Det
 527   //
 528
 529   // DB Setting
 530   // Get cal
 531   AliCDBEntry *entry = AliCDBManager::Instance()->Get("TRD/Calib/ChamberGainFactor",AliCDBManager::Instance()->GetRun(),fVersionGainUsed,fSubVersionGainUsed);
 532
 533   if (entry) {
 534
 535
 536     if( fCalDetGain ){
 537       fCalDetGain->~AliTRDCalDet();
 538       new(fCalDetGain) AliTRDCalDet(*(dynamic_cast<AliTRDCalDet *>(entry->GetObject())));
 539     }else fCalDetGain = new AliTRDCalDet(*(dynamic_cast<AliTRDCalDet *>(entry->GetObject())));
 540
 541
 542   }
 543   else {
 544
 545     AliError("No new gain calibration entry found");
 546     return kFALSE;
 547
 548   }
 549
 550   // title CH2d
 551   TString name("Ver");
 552   name += fVersionGainUsed;
 553   name += "Subver";
 554   name += fSubVersionGainUsed;
 555   name += "Nz";
 556   name += fCalibraMode->GetNz(0);
 557   name += "Nrphi";
 558   name += fCalibraMode->GetNrphi(0);
 559
 560   fCH2d->SetTitle(name);
 561
 562   // title PH2d
 563   TString namee("Ver");
 564   namee += fVersionVdriftUsed;
 565   namee += "Subver";
 566   namee += fSubVersionVdriftUsed;
 567   namee += "Nz";
 568   namee += fCalibraMode->GetNz(1);
 569   namee += "Nrphi";
 570   namee += fCalibraMode->GetNrphi(1);
 571
 572   fPH2d->SetTitle(namee);
 573
 574
 575   return kTRUE;
 576
 577 }
 578 ///////////////////////////////////////////////////////////////////////////////////////////////////////////////
 579 Bool_t AliTRDCalibraFillHisto::InitCalPad(Int_t detector)
 580 {
 581   //
 582   // Init the Gain Cal Pad
 583   //
 584
 585   // DB Setting
 586   // Get cal
 587   AliCDBEntry *entry = AliCDBManager::Instance()->Get("TRD/Calib/LocalGainFactor",AliCDBManager::Instance()->GetRun(),fVersionGainLocalUsed,fSubVersionGainLocalUsed);
 588
 589   if (entry) {
 590
 591
 592     if( fCalROCGain ){
 593       fCalROCGain->~AliTRDCalROC();
 594       new(fCalROCGain) AliTRDCalROC(*((dynamic_cast<AliTRDCalPad *>(entry->GetObject()))->GetCalROC(detector)));
 595     }else fCalROCGain = new AliTRDCalROC(*((dynamic_cast<AliTRDCalPad *>(entry->GetObject()))->GetCalROC(detector)));
 596
 597
 598   }
 599   else {
 600
 601     AliError("No new gain calibration entry found");
 602     return kFALSE;
 603
 604   }
 605
 606   return kTRUE;
 607
 608 }
 609 //____________Offline tracking in the AliTRDtracker____________________________
 610 Bool_t AliTRDCalibraFillHisto::UpdateHistograms(const AliTRDtrack *t)
 611 {
 612   //
 613   // Use AliTRDtrack for the calibration
 614   //
 615
 616
 617   AliTRDcluster *cl = 0x0;       // pointeur to cluster
 618   Int_t index0 = 0;              // index of the first cluster in the current chamber
 619   Int_t index1 = 0;              // index of the current cluster in the current chamber
 620
 621   //////////////////////////////
 622   // loop over the clusters
 623   ///////////////////////////////
 624   while((cl = t->GetCluster(index1))){
 625
 626     /////////////////////////////////////////////////////////////////////////
 627     // Fill the infos for the previous clusters if not the same detector
 628     ////////////////////////////////////////////////////////////////////////
 629     Int_t detector = cl->GetDetector();
 630     if ((detector != fDetectorPreviousTrack) &&
 631         (index0 != index1)) {
 632
 633       fNumberTrack++;
 634
 635       //If the same track, then look if the previous detector is in
 636       //the same plane, if yes: not a good track
 637       if ((GetLayer(detector) == GetLayer(fDetectorPreviousTrack))) {
 638         return kFALSE;
 639       }
 640
 641       // Fill only if the track doesn't touch a masked pad or doesn't
 642       if (fGoodTracklet) {
 643
 644         // drift velocity unables to cut bad tracklets
 645         Bool_t  pass = FindP1TrackPHtracklet(t,index0,index1);
 646
 647         // Gain calibration
 648         if (fCH2dOn) {
 649           FillTheInfoOfTheTrackCH(index1-index0);
 650         }
 651
 652         // PH calibration
 653         if (fPH2dOn) {
 654           FillTheInfoOfTheTrackPH();
 655         }
 656
 657         if(pass && fPRF2dOn) HandlePRFtracklet(t,index0,index1);
 658
 659
 660       } // if a good tracklet
 661
 662       // reset stuff
 663       ResetfVariablestracklet();
 664       index0 = index1;
 665
 666     } // Fill at the end the charge
 667
 668
 669     /////////////////////////////////////////////////////////////
 670     // Localise and take the calib gain object for the detector
 671     ////////////////////////////////////////////////////////////
 672     if (detector != fDetectorPreviousTrack) {
 673
 674       //Localise the detector
 675       LocalisationDetectorXbins(detector);
 676
 677       // Get cal
 678       AliTRDcalibDB *cal = AliTRDcalibDB::Instance();
 679       if (!cal) {
 680         AliInfo("Could not get calibDB");
 681         return kFALSE;
 682       }
 683
 684       // Get calib objects
 685       if(!fIsHLT) InitCalPad(detector);
 686
 687     }
 688
 689     // Reset the detectbjobsor
 690     fDetectorPreviousTrack = detector;
 691
 692     ///////////////////////////////////////
 693     // Store the info of the cluster
 694     ///////////////////////////////////////
 695     Int_t row = cl->GetPadRow();
 696     Int_t col = cl->GetPadCol();
 697     CheckGoodTrackletV1(cl);
 698     Int_t     group[2] = {0,0};
 699     if(fCH2dOn)  group[0]  = CalculateCalibrationGroup(0,row,col);
 700     if(fPH2dOn)  group[1]  = CalculateCalibrationGroup(1,row,col);
 701     StoreInfoCHPHtrack(cl,t->GetClusterdQdl(index1),group,row,col);
 702
 703     index1++;
 704
 705   } // while on clusters
 706
 707   ///////////////////////////
 708   // Fill the last plane
 709   //////////////////////////
 710   if( index0 != index1 ){
 711
 712     fNumberTrack++;
 713
 714     if (fGoodTracklet) {
 715
 716       // drift velocity unables to cut bad tracklets
 717       Bool_t  pass = FindP1TrackPHtracklet(t,index0,index1);
 718
 719       // Gain calibration
 720       if (fCH2dOn) {
 721         FillTheInfoOfTheTrackCH(index1-index0);
 722       }
 723
 724       // PH calibration
 725       if (fPH2dOn) {
 726         FillTheInfoOfTheTrackPH();
 727       }
 728
 729       if(pass && fPRF2dOn) HandlePRFtracklet(t,index0,index1);
 730
 731     } // if a good tracklet
 732
 733   }
 734
 735   // reset stuff
 736   ResetfVariablestracklet();
 737
 738   return kTRUE;
 739
 740 }
 741 //____________Offline tracking in the AliTRDtracker____________________________
 742 Bool_t AliTRDCalibraFillHisto::UpdateHistogramsV1(const AliTRDtrackV1 *t)
 743 {
 744   //
 745   // Use AliTRDtrackV1 for the calibration
 746   //
 747
 748
 749   const AliTRDseedV1 *tracklet = 0x0;          // tracklet per plane
 750   AliTRDcluster *cl      = 0x0;                // cluster attached now to the tracklet
 751   AliTRDcluster *cls     = 0x0;                // shared cluster attached now to the tracklet
 752   Bool_t         newtr   = kTRUE;              // new track
 753
 754   // Get cal
 755   //  AliTRDcalibDB *cal = AliTRDcalibDB::Instance();
 756   /*
 757   if (!cal) {
 758     AliInfo("Could not get calibDB");
 759     return kFALSE;
 760   }
 761 */
 762   if (!fCalibDB) {
 763     AliInfo("Could not get calibDB");
 764     return kFALSE;
 765   }
 766
 767
 768   ///////////////////////////
 769   // loop over the tracklet
 770   ///////////////////////////
 771   for(Int_t itr = 0; itr < 6; itr++){
 772
 773     if(!(tracklet = t->GetTracklet(itr))) continue;
 774     if(!tracklet->IsOK()) continue;
 775     fNumberTrack++;
 776     ResetfVariablestracklet();
 777
 778     //////////////////////////////////////////
 779     // localisation of the tracklet and dqdl
 780     //////////////////////////////////////////
 781     Int_t layer    = tracklet->GetPlane();
 782     Int_t ic = 0;
 783     while(!(cl = tracklet->GetClusters(ic++))) continue;
 784     Int_t detector = cl->GetDetector();
 785     if (detector != fDetectorPreviousTrack) {
 786       // if not a new track
 787       if(!newtr){
 788         // don't use the rest of this track if in the same plane
 789         if (layer == GetLayer(fDetectorPreviousTrack)) {
 790           //printf("bad tracklet, same layer for detector %d\n",detector);
 791           break;
 792         }
 793       }
 794       //Localise the detector bin
 795       LocalisationDetectorXbins(detector);
 796       // Get calib objects
 797       if(!fIsHLT) InitCalPad(detector);
 798
 799       // reset
 800       fDetectorPreviousTrack = detector;
 801     }
 802     newtr = kFALSE;
 803
 804     ////////////////////////////
 805     // loop over the clusters
 806     ////////////////////////////
 807     Int_t nbclusters = 0;
 808     for(int jc=0; jc<AliTRDseedV1::kNtb; jc++){
 809       if(!(cl = tracklet->GetClusters(jc))) continue;
 810       nbclusters++;
 811
 812       // Store the info bis of the tracklet
 813       Int_t row = cl->GetPadRow();
 814       Int_t col = cl->GetPadCol();
 815       CheckGoodTrackletV1(cl);
 816       Int_t     group[2] = {0,0};
 817       if(fCH2dOn)  group[0]  = CalculateCalibrationGroup(0,row,col);
 818       if(fPH2dOn)  group[1]  = CalculateCalibrationGroup(1,row,col);
 819       // Add the charge if shared cluster
 820       cls = tracklet->GetClusters(jc+AliTRDseedV1::kNtb);
 821       //
 822       StoreInfoCHPHtrack(cl, tracklet->GetdQdl(jc),group,row,col,cls);
 823     }
 824
 825     ////////////////////////////////////////
 826     // Fill the stuffs if a good tracklet
 827     ////////////////////////////////////////
 828     if (fGoodTracklet) {
 829
 830       // drift velocity unables to cut bad tracklets
 831       Bool_t  pass = FindP1TrackPHtrackletV1(tracklet, nbclusters);
 832
 833       //printf("pass %d and nbclusters %d\n",pass,nbclusters);
 834
 835       // Gain calibration
 836       if (fCH2dOn) {
 837         FillTheInfoOfTheTrackCH(nbclusters);
 838       }
 839
 840       // PH calibration
 841       if (fPH2dOn) {
 842         FillTheInfoOfTheTrackPH();
 843       }
 844
 845       if(pass && fPRF2dOn) HandlePRFtrackletV1(tracklet,nbclusters);
 846
 847     } // if a good tracklet
 848   }
 849
 850   return kTRUE;
 851
 852 }
 853 ///////////////////////////////////////////////////////////////////////////////////
 854 // Routine inside the update with AliTRDtrack
 855 ///////////////////////////////////////////////////////////////////////////////////
 856 //____________Offine tracking in the AliTRDtracker_____________________________
 857 Bool_t AliTRDCalibraFillHisto::FindP1TrackPHtracklet(const AliTRDtrack *t, Int_t index0, Int_t index1)
 858 {
 859   //
 860   // Drift velocity calibration:
 861   // Fit the clusters with a straight line
 862   // From the slope find the drift velocity
 863   //
 864
 865   //Number of points: if less than 3 return kFALSE
 866   Int_t npoints = index1-index0;
 867   if(npoints <= 2) return kFALSE;
 868
 869   /////////////////
 870   //Variables
 871   ////////////////
 872   Int_t    detector                   = ((AliTRDcluster *) t->GetCluster(index0))->GetDetector(); //detector
 873   // parameters of the track
 874   Double_t snp                        = t->GetSnpPlane(GetLayer(detector)); // sin angle in the plan yx track
 875   Double_t tgl                        = t->GetTglPlane(GetLayer(detector)); // dz/dl and not dz/dx!
 876   Double_t tnp                        = 0.0;                                // tan angle in the plan xy track
 877   if( TMath::Abs(snp) <  1.){
 878     tnp = snp / TMath::Sqrt((1.-snp)*(1.+snp));
 879   }
 880   Float_t dzdx                        = tgl*TMath::Sqrt(1+tnp*tnp);         // dz/dx now from dz/dl
 881   // tilting pad and cross row
 882   Int_t    crossrow                   = 0;                                  // if it crosses a pad row
 883   Int_t    rowp                       = -1;                                 // if it crosses a pad row
 884   AliTRDpadPlane *padplane            = fGeo->GetPadPlane(GetLayer(detector),GetStack(detector));
 885   Double_t tiltingangle               = padplane->GetTiltingAngle();        // tiltingangle of the pad
 886   Float_t  tnt                        = TMath::Tan(tiltingangle/180.*TMath::Pi()); // tan tiltingangle
 887   // linear fit
 888   fLinearFitterTracklet->ClearPoints();
 889   Double_t dydt                       = 0.0;                                // dydt tracklet after straight line fit
 890   Double_t errorpar                   = 0.0;                                // error after straight line fit on dy/dt
 891   Double_t pointError                 = 0.0;                                // error after straight line fit
 892   Int_t     nbli                      = 0;                                  // number linear fitter points
 893
 894   //////////////////////////////
 895   // loop over clusters
 896   ////////////////////////////
 897   for(Int_t k = 0; k < npoints; k++){
 898
 899     AliTRDcluster *cl                 = (AliTRDcluster *) t->GetCluster(k+index0);
 900     if((fLimitChargeIntegration) && (!cl->IsInChamber())) continue;
 901     Double_t ycluster                 = cl->GetY();
 902     Int_t time                        = cl->GetPadTime();
 903     Double_t timeis                   = time/fSf;
 904     //Double_t q                        = cl->GetQ();
 905     //See if cross two pad rows
 906     Int_t    row                      = cl->GetPadRow();
 907     if(k==0) rowp                     = row;
 908     if(row != rowp) crossrow          = 1;
 909
 910     fLinearFitterTracklet->AddPoint(&timeis,ycluster,1);
 911     nbli++;
 912
 913   }
 914
 915   //////////////////////////////
 916   // linear fit
 917   //////////////////////////////
 918   if(nbli <= 2){
 919     fLinearFitterTracklet->ClearPoints();
 920     return kFALSE;
 921   }
 922   TVectorD pars;
 923   fLinearFitterTracklet->Eval();
 924   fLinearFitterTracklet->GetParameters(pars);
 925   pointError  =  TMath::Sqrt(fLinearFitterTracklet->GetChisquare()/(nbli-2));
 926   errorpar    =  fLinearFitterTracklet->GetParError(1)*pointError;
 927   dydt        = pars[1];
 928   fLinearFitterTracklet->ClearPoints();
 929
 930   /////////////////////////////
 931   // debug
 932   ////////////////////////////
 933   if(fDebugLevel > 0){
 934     if ( !fDebugStreamer ) {
 935       //debug stream
 936       TDirectory *backup = gDirectory;
 937       fDebugStreamer = new TTreeSRedirector("TRDdebugCalibraFill.root");
 938       if ( backup ) backup->cd();  //we don't want to be cd'd to the debug streamer
 939     }
 940
 941
 942     (* fDebugStreamer) << "FindP1TrackPHtracklet0"<<
 943       //"snpright="<<snpright<<
 944       "npoints="<<npoints<<
 945       "nbli="<<nbli<<
 946       "detector="<<detector<<
 947       "snp="<<snp<<
 948       "tnp="<<tnp<<
 949       "tgl="<<tgl<<
 950       "tnt="<<tnt<<
 951       "dydt="<<dydt<<
 952       "dzdx="<<dzdx<<
 953       "crossrow="<<crossrow<<
 954       "errorpar="<<errorpar<<
 955       "pointError="<<pointError<<
 956       "\n";
 957
 958
 959     Int_t nbclusters = index1-index0;
 960     Int_t layer      = GetLayer(fDetectorPreviousTrack);
 961
 962     (* fDebugStreamer) << "FindP1TrackPHtracklet1"<<
 963       //"snpright="<<snpright<<
 964       "nbclusters="<<nbclusters<<
 965       "detector="<<fDetectorPreviousTrack<<
 966       "layer="<<layer<<
 967       "\n";
 968
 969   }
 970
 971   ///////////////////////////
 972   // quality cuts
 973   ///////////////////////////
 974   if(npoints < fNumberClusters) return kFALSE;
 975   if(npoints > fNumberClustersf) return kFALSE;
 976   if(pointError >= 0.1) return kFALSE;
 977   if(crossrow == 1) return kFALSE;
 978
 979   ////////////////////////////
 980   // fill
 981   ////////////////////////////
 982   if(fLinearFitterOn){
 983     //Add to the linear fitter of the detector
 984     if( TMath::Abs(snp) <  1.){
 985       Double_t x = tnp-dzdx*tnt;
 986       if(fLinearFitterDebugOn) {
 987         (GetLinearFitter(detector,kTRUE))->AddPoint(&x,dydt);
 988         fEntriesLinearFitter[detector]++;
 989       }
 990       fLinearVdriftFit->Update(detector,x,pars[1]);
 991     }
 992   }
 993
 994   return kTRUE;
 995 }
 996 //____________Offine tracking in the AliTRDtracker_____________________________
 997 Bool_t AliTRDCalibraFillHisto::FindP1TrackPHtrackletV1(const AliTRDseedV1 *tracklet, Int_t nbclusters)
 998 {
 999   //
1000   // Drift velocity calibration:
1001   // Fit the clusters with a straight line
1002   // From the slope find the drift velocity
1003   //
1004
1005   ////////////////////////////////////////////////
1006   //Number of points: if less than 3 return kFALSE
1007   /////////////////////////////////////////////////
1008   if(nbclusters <= 2) return kFALSE;
1009
1010   ////////////
1011   //Variables
1012   ////////////
1013   // results of the linear fit
1014   Double_t dydt                       = 0.0;                                // dydt tracklet after straight line fit
1015   Double_t errorpar                   = 0.0;                                // error after straight line fit on dy/dt
1016   Double_t pointError                 = 0.0;                                // error after straight line fit
1017   // pad row problemes: avoid tracklet that cross pad rows, tilting angle in the constant
1018   Int_t    crossrow                   = 0;                                  // if it crosses a pad row
1019   Int_t    rowp                       = -1;                                 // if it crosses a pad row
1020   Float_t  tnt                        = tracklet->GetTilt();                // tan tiltingangle
1021   fLinearFitterTracklet->ClearPoints();
1022
1023
1024   ///////////////////////////////////////////
1025   // Take the parameters of the track
1026   //////////////////////////////////////////
1027   // take now the snp, tnp and tgl from the track
1028   Double_t snp = tracklet->GetSnp();             // sin dy/dx at the end of the chamber
1029   Double_t tnp = 0.0;                            // dy/dx at the end of the chamber
1030   if( TMath::Abs(snp) <  1.){
1031     tnp = snp / TMath::Sqrt((1.-snp)*(1.+snp));
1032   }
1033   Double_t tgl  = tracklet->GetTgl();           // dz/dl
1034   Double_t dzdx = tgl*TMath::Sqrt(1+tnp*tnp);   // dz/dx calculated from dz/dl
1035   // at the entrance
1036   //Double_t tnp = tracklet->GetYref(1);      // dy/dx at the entrance of the chamber
1037   //Double_t tgl = tracklet->GetZref(1);      // dz/dl at the entrance of the chamber
1038   //Double_t dzdx = tgl;                      //*TMath::Sqrt(1+tnp*tnp); // dz/dx from dz/dl
1039   // at the end with correction due to linear fit
1040   //Double_t tnp = tracklet->GetYfit(1);      // dy/dx at the end of the chamber after fit correction
1041   //Double_t tgl = tracklet->GetZfit(1);      // dz/dl at the end of the chamber after fit correction
1042
1043
1044   ////////////////////////////
1045   // loop over the clusters
1046   ////////////////////////////
1047   Int_t  nbli = 0;
1048   AliTRDcluster *cl                   = 0x0;
1049   //////////////////////////////
1050   // Check no shared clusters
1051   //////////////////////////////
1052   for(int icc=AliTRDseedV1::kNtb; icc<AliTRDseedV1::kNclusters; icc++){
1053     if((cl = tracklet->GetClusters(icc)))  crossrow = 1;
1054   }
1055   //////////////////////////////////
1056   // Loop clusters
1057   //////////////////////////////////
1058   for(int ic=0; ic<AliTRDseedV1::kNtb; ic++){
1059     if(!(cl = tracklet->GetClusters(ic))) continue;
1060     if((fLimitChargeIntegration) && (!cl->IsInChamber())) continue;
1061
1062     Double_t ycluster                 = cl->GetY();
1063     Int_t time                        = cl->GetPadTime();
1064     Double_t timeis                   = time/fSf;
1065     //See if cross two pad rows
1066     Int_t    row                      = cl->GetPadRow();
1067     if(rowp==-1) rowp                 = row;
1068     if(row != rowp) crossrow          = 1;
1069
1070     fLinearFitterTracklet->AddPoint(&timeis,ycluster,1);
1071     nbli++;
1072
1073
1074   }
1075
1076   ////////////////////////////////////
1077   // Do the straight line fit now
1078   ///////////////////////////////////
1079   if(nbli <= 2){
1080     fLinearFitterTracklet->ClearPoints();
1081     return kFALSE;
1082   }
1083   TVectorD pars;
1084   fLinearFitterTracklet->Eval();
1085   fLinearFitterTracklet->GetParameters(pars);
1086   pointError  =  TMath::Sqrt(fLinearFitterTracklet->GetChisquare()/(nbli-2));
1087   errorpar    =  fLinearFitterTracklet->GetParError(1)*pointError;
1088   dydt        = pars[1];
1089   //printf("chis %f, nbli %d, pointError %f, parError %f, errorpar %f\n",fLinearFitterTracklet->GetChisquare(),nbli,pointError,fLinearFitterTracklet->GetParError(1),errorpar);
1090   fLinearFitterTracklet->ClearPoints();
1091
1092   ////////////////////////////////
1093   // Debug stuff
1094   ///////////////////////////////
1095
1096
1097   if(fDebugLevel > 0){
1098     if ( !fDebugStreamer ) {
1099       //debug stream
1100       TDirectory *backup = gDirectory;
1101       fDebugStreamer = new TTreeSRedirector("TRDdebugCalibraFill.root");
1102       if ( backup ) backup->cd();  //we don't want to be cd'd to the debug streamer
1103     }
1104
1105
1106     Int_t layer = GetLayer(fDetectorPreviousTrack);
1107
1108     (* fDebugStreamer) << "FindP1TrackPHtrackletV1"<<
1109       //"snpright="<<snpright<<
1110       "nbli="<<nbli<<
1111       "nbclusters="<<nbclusters<<
1112       "detector="<<fDetectorPreviousTrack<<
1113       "layer="<<layer<<
1114       "snp="<<snp<<
1115       "tnp="<<tnp<<
1116       "tgl="<<tgl<<
1117       "tnt="<<tnt<<
1118       "dydt="<<dydt<<
1119       "dzdx="<<dzdx<<
1120       "crossrow="<<crossrow<<
1121       "errorpar="<<errorpar<<
1122       "pointError="<<pointError<<
1123       "\n";
1124
1125   }
1126
1127   /////////////////////////
1128   // Cuts quality
1129   ////////////////////////
1130
1131   if(nbclusters < fNumberClusters) return kFALSE;
1132   if(nbclusters > fNumberClustersf) return kFALSE;
1133   if(pointError >= 0.3) return kFALSE;
1134   if(crossrow == 1) return kFALSE;
1135
1136   ///////////////////////
1137   // Fill
1138   //////////////////////
1139
1140   if(fLinearFitterOn){
1141     //Add to the linear fitter of the detector
1142     if( TMath::Abs(snp) <  1.){
1143       Double_t x = tnp-dzdx*tnt;
1144       if(fLinearFitterDebugOn) {
1145         (GetLinearFitter(fDetectorPreviousTrack,kTRUE))->AddPoint(&x,dydt);
1146         fEntriesLinearFitter[fDetectorPreviousTrack]++;
1147       }
1148       fLinearVdriftFit->Update(fDetectorPreviousTrack,x,pars[1]);
1149     }
1150   }
1151
1152   return kTRUE;
1153 }
1154 //____________Offine tracking in the AliTRDtracker_____________________________
1155 Bool_t AliTRDCalibraFillHisto::HandlePRFtracklet(const AliTRDtrack *t, Int_t index0, Int_t index1)
1156 {
1157   //
1158   // PRF width calibration
1159   // Assume a Gaussian shape: determinate the position of the three pad clusters
1160   // Fit with a straight line
1161   // Take the fitted values for all the clusters (3 or 2 pad clusters)
1162   // Fill the PRF as function of angle of the track
1163   //
1164   //
1165
1166
1167   //////////////////////////
1168   // variables
1169   /////////////////////////
1170   Int_t npoints  = index1-index0;                                           // number of total points
1171   Int_t nb3pc    = 0;                                                       // number of three pads clusters used for fit
1172   Int_t detector = ((AliTRDcluster *) t->GetCluster(index0))->GetDetector(); // detector
1173   // To see the difference due to the fit
1174   Double_t *padPositions;
1175   padPositions = new Double_t[npoints];
1176   for(Int_t k = 0; k < npoints; k++){
1177     padPositions[k] = 0.0;
1178   }
1179   // Take the tgl and snp with the track t now
1180   Double_t  tgl = t->GetTglPlane(GetLayer(detector)); //dz/dl and not dz/dx
1181   Double_t  snp = t->GetSnpPlane(GetLayer(detector)); // sin angle in xy plan
1182   Float_t  dzdx = 0.0;                                // dzdx
1183   Float_t  tnp  = 0.0;
1184   if(TMath::Abs(snp) < 1.0){
1185     tnp = snp / TMath::Sqrt((1.-snp)*(1.+snp));
1186     dzdx = tgl*TMath::Sqrt(1+tnp*tnp);
1187   }
1188   // linear fitter
1189   fLinearFitterTracklet->ClearPoints();
1190
1191   ///////////////////////////
1192   // calcul the tnp group
1193   ///////////////////////////
1194   Bool_t echec   = kFALSE;
1195   Double_t shift = 0.0;
1196   //Calculate the shift in x coresponding to this tnp
1197   if(fNgroupprf != 0.0){
1198     shift      = -3.0*(fNgroupprf-1)-1.5;
1199     Double_t limithigh  = -0.2*(fNgroupprf-1);
1200     if((tnp < (-0.2*fNgroupprf)) || (tnp > (0.2*fNgroupprf))) echec = kTRUE;
1201     else{
1202       while(tnp > limithigh){
1203         limithigh += 0.2;
1204         shift += 3.0;
1205       }
1206     }
1207   }
1208   if(echec) {
1209     delete [] padPositions;
1210     return kFALSE;
1211   }
1212
1213   //////////////////////
1214   // loop clusters
1215   /////////////////////
1216   for(Int_t k = 0;  k < npoints; k++){
1217     //Take the cluster
1218     AliTRDcluster *cl  = (AliTRDcluster *) t->GetCluster(k+index0);
1219     if(!cl) continue;
1220     Short_t  *signals  = cl->GetSignals();
1221     Double_t     time  = cl->GetPadTime();
1222     //Calculate x if possible
1223     Float_t xcenter    = 0.0;
1224     Bool_t  echec1      = kTRUE;
1225     if((time<=7) || (time>=21)) continue;
1226     // Center 3 balanced: position with the center of the pad
1227     if ((((Float_t) signals[3]) > 0.0) &&
1228         (((Float_t) signals[2]) > 0.0) &&
1229         (((Float_t) signals[4]) > 0.0)) {
1230       echec1 = kFALSE;
1231       // Security if the denomiateur is 0
1232       if ((((Float_t) (((Float_t) signals[3]) * ((Float_t) signals[3]))) /
1233            ((Float_t) (((Float_t) signals[2]) * ((Float_t) signals[4])))) != 1.0) {
1234         xcenter = 0.5 * (TMath::Log((Float_t) (((Float_t) signals[4]) / ((Float_t) signals[2]))))
1235           / (TMath::Log(((Float_t) (((Float_t) signals[3]) * ((Float_t) signals[3])))
1236                         / ((Float_t) (((Float_t) signals[2]) * ((Float_t) signals[4])))));
1237       }
1238       else {
1239         echec1 = kTRUE;
1240       }
1241     }
1242     if(TMath::Abs(xcenter) > 0.5) echec = kTRUE;
1243     if(echec) continue;
1244     //if no echec: calculate with the position of the pad
1245     // Position of the cluster
1246     Double_t       padPosition = xcenter +  cl->GetPadCol();
1247     padPositions[k]            = padPosition;
1248     nb3pc++;
1249     fLinearFitterTracklet->AddPoint(&time, padPosition,1);
1250   }//clusters loop
1251
1252
1253   /////////////////////////////
1254   // fit
1255   ////////////////////////////
1256   if(nb3pc < 3) {
1257     delete [] padPositions;
1258     fLinearFitterTracklet->ClearPoints();
1259     return kFALSE;
1260   }
1261   fLinearFitterTracklet->Eval();
1262   TVectorD line(2);
1263   fLinearFitterTracklet->GetParameters(line);
1264   Float_t  pointError  = -1.0;
1265   if( fLinearFitterTracklet->GetChisquare()>=0.0) {
1266     pointError  =  TMath::Sqrt( fLinearFitterTracklet->GetChisquare()/(nb3pc-2));
1267   }
1268   fLinearFitterTracklet->ClearPoints();
1269
1270   /////////////////////////////////////////////////////
1271   // Now fill the PRF: second loop over clusters
1272   /////////////////////////////////////////////////////
1273   for(Int_t k = 0;  k < npoints; k++){
1274     //Take the cluster
1275     AliTRDcluster *cl      = (AliTRDcluster *) t->GetCluster(k+index0);
1276     Short_t *signals       = cl->GetSignals();              // signal
1277     Double_t time          = cl->GetPadTime();              // time bin
1278     Float_t padPosTracklet = line[0]+line[1]*time;          // reconstruct position from fit
1279     Float_t padPos         = cl->GetPadCol();               // middle pad
1280     Double_t dpad          = padPosTracklet - padPos;       // reconstruct position relative to middle pad from fit
1281     Float_t ycenter        = 0.0;                           // relative center charge
1282     Float_t ymin           = 0.0;                           // relative left charge
1283     Float_t ymax           = 0.0;                           // relative right charge
1284
1285     //Requiere simply two pads clusters at least
1286     if(((((Float_t) signals[3]) > 0.0) && (((Float_t) signals[2]) > 0.0)) ||
1287        ((((Float_t) signals[3]) > 0.0) && (((Float_t) signals[4]) > 0.0))){
1288       Float_t sum     = ((Float_t) signals[2]) + ((Float_t) signals[3]) + ((Float_t) signals[4]);
1289       if(sum > 0.0) ycenter = ((Float_t) signals[3])/ sum;
1290       if(sum > 0.0) ymin    = ((Float_t) signals[2])/ sum;
1291       if(sum > 0.0) ymax    = ((Float_t) signals[4])/ sum;
1292     }
1293
1294     //calibration group
1295     Int_t     rowcl        = cl->GetPadRow();                           // row of cluster
1296     Int_t     colcl        = cl->GetPadCol();                           // col of cluster
1297     Int_t     grouplocal   = CalculateCalibrationGroup(2,rowcl,colcl);  // calcul the corresponding group
1298     Int_t     caligroup    = fCalibraMode->GetXbins(2)+ grouplocal;     // calcul the corresponding group
1299     Float_t   xcl          = cl->GetY();                                // y cluster
1300     Float_t   qcl          = cl->GetQ();                                // charge cluster
1301     Int_t     layer        = GetLayer(detector);                        // layer
1302     Int_t     stack        = GetStack(detector);                        // stack
1303     Double_t  xdiff        = dpad;                                      // reconstructed position constant
1304     Double_t  x            = dpad;                                      // reconstructed position moved
1305     Float_t   ep           = pointError;                                // error of fit
1306     Float_t   signal1      = (Float_t)signals[1];                       // signal at the border
1307     Float_t   signal3      = (Float_t)signals[3];                       // signal
1308     Float_t   signal2      = (Float_t)signals[2];                       // signal
1309     Float_t   signal4      = (Float_t)signals[4];                       // signal
1310     Float_t   signal5      = (Float_t)signals[5];                       // signal at the border
1311
1312     //////////////////////////////
1313     // debug
1314     /////////////////////////////
1315     if(fDebugLevel > 0){
1316       if ( !fDebugStreamer ) {
1317         //debug stream
1318         TDirectory *backup = gDirectory;
1319         fDebugStreamer = new TTreeSRedirector("TRDdebugCalibraFill.root");
1320         if ( backup ) backup->cd();  //we don't want to be cd'd to the debug streamer
1321       }
1322
1323
1324       x = xdiff;
1325       Int_t type=0;
1326       Float_t y = ycenter;
1327       (* fDebugStreamer) << "HandlePRFtracklet"<<
1328         "caligroup="<<caligroup<<
1329         "detector="<<detector<<
1330         "layer="<<layer<<
1331         "stack="<< stack <<
1332         "npoints="<<npoints<<
1333         "Np="<<nb3pc<<
1334         "ep="<<ep<<
1335         "type="<<type<<
1336         "snp="<<snp<<
1337         "tnp="<<tnp<<
1338         "tgl="<<tgl<<
1339         "dzdx="<<dzdx<<
1340         "padPos="<<padPos<<
1341         "padPosition="<<padPositions[k]<<
1342         "padPosTracklet="<<padPosTracklet<<
1343         "x="<<x<<
1344         "y="<<y<<
1345         "xcl="<<xcl<<
1346         "qcl="<<qcl<<
1347         "signal1="<<signal1<<
1348         "signal2="<<signal2<<
1349         "signal3="<<signal3<<
1350         "signal4="<<signal4<<
1351         "signal5="<<signal5<<
1352         "time="<<time<<
1353         "\n";
1354       x=-(xdiff+1);
1355       y = ymin;
1356       type=-1;
1357       (* fDebugStreamer) << "HandlePRFtracklet"<<
1358         "caligroup="<<caligroup<<
1359         "detector="<<detector<<
1360         "layer="<<layer<<
1361         "stack="<<stack<<
1362         "npoints="<<npoints<<
1363         "Np="<<nb3pc<<
1364         "ep="<<ep<<
1365         "type="<<type<<
1366         "snp="<<snp<<
1367         "tnp="<<tnp<<
1368         "tgl="<<tgl<<
1369         "dzdx="<<dzdx<<
1370         "padPos="<<padPos<<
1371         "padPosition="<<padPositions[k]<<
1372         "padPosTracklet="<<padPosTracklet<<
1373         "x="<<x<<
1374         "y="<<y<<
1375         "xcl="<<xcl<<
1376         "qcl="<<qcl<<
1377         "signal1="<<signal1<<
1378         "signal2="<<signal2<<
1379         "signal3="<<signal3<<
1380         "signal4="<<signal4<<
1381         "signal5="<<signal5<<
1382         "time="<<time<<
1383         "\n";
1384       x=1-xdiff;
1385       y = ymax;
1386       type=1;
1387       (* fDebugStreamer) << "HandlePRFtracklet"<<
1388         "caligroup="<<caligroup<<
1389         "detector="<<detector<<
1390         "layer="<<layer<<
1391         "stack="<<stack<<
1392         "npoints="<<npoints<<
1393         "Np="<<nb3pc<<
1394         "ep="<<ep<<
1395         "type="<<type<<
1396         "snp="<<snp<<
1397         "tnp="<<tnp<<
1398         "tgl="<<tgl<<
1399         "dzdx="<<dzdx<<
1400         "padPos="<<padPos<<
1401         "padPosition="<<padPositions[k]<<
1402         "padPosTracklet="<<padPosTracklet<<
1403         "x="<<x<<
1404         "y="<<y<<
1405         "xcl="<<xcl<<
1406         "qcl="<<qcl<<
1407         "signal1="<<signal1<<
1408         "signal2="<<signal2<<
1409         "signal3="<<signal3<<
1410         "signal4="<<signal4<<
1411         "signal5="<<signal5<<
1412         "time="<<time<<
1413         "\n";
1414
1415     }
1416
1417     ////////////////////////////
1418     // quality cuts
1419     ///////////////////////////
1420     if(npoints < fNumberClusters) continue;
1421     if(npoints > fNumberClustersf) continue;
1422     if(nb3pc <= 5) continue;
1423     if((time >= 21) || (time < 7)) continue;
1424     if(TMath::Abs(snp) >= 1.0) continue;
1425     if(TMath::Abs(qcl) < 80) continue;
1426
1427     ////////////////////////////
1428     // Fill
1429     ///////////////////////////
1430     if (fHisto2d) {
1431       if(TMath::Abs(dpad) < 1.5) {
1432         fPRF2d->Fill(shift+dpad,(caligroup+0.5),ycenter);
1433         fPRF2d->Fill(shift-dpad,(caligroup+0.5),ycenter);
1434       }
1435       if((ymin > 0.0) && (TMath::Abs(dpad+1.0) < 1.5)) {
1436         fPRF2d->Fill(shift-(dpad+1.0),(caligroup+0.5),ymin);
1437         fPRF2d->Fill(shift+(dpad+1.0),(caligroup+0.5),ymin);
1438       }
1439       if((ymax > 0.0) && (TMath::Abs(dpad-1.0) < 1.5)) {
1440         fPRF2d->Fill(shift+1.0-dpad,(caligroup+0.5),ymax);
1441         fPRF2d->Fill(shift-1.0+dpad,(caligroup+0.5),ymax);
1442       }
1443     }
1444     if (fVector2d) {
1445       if(TMath::Abs(dpad) < 1.5) {
1446         fCalibraVector->UpdateVectorPRF(fDetectorPreviousTrack,grouplocal,shift+dpad,ycenter);
1447         fCalibraVector->UpdateVectorPRF(fDetectorPreviousTrack,grouplocal,shift-dpad,ycenter);
1448       }
1449       if((ymin > 0.0) && (TMath::Abs(dpad+1.0) < 1.5)) {
1450         fCalibraVector->UpdateVectorPRF(fDetectorPreviousTrack,grouplocal,shift-(dpad+1.0),ymin);
1451         fCalibraVector->UpdateVectorPRF(fDetectorPreviousTrack,grouplocal,shift+(dpad+1.0),ymin);
1452       }
1453       if((ymax > 0.0)  && (TMath::Abs(dpad-1.0) < 1.5)) {
1454         fCalibraVector->UpdateVectorPRF(fDetectorPreviousTrack,grouplocal,shift+1.0-dpad,ymax);
1455         fCalibraVector->UpdateVectorPRF(fDetectorPreviousTrack,grouplocal,shift-1.0+dpad,ymax);
1456       }
1457     }
1458   }
1459   delete [] padPositions;
1460   return kTRUE;
1461
1462 }
1463 //____________Offine tracking in the AliTRDtracker_____________________________
1464 Bool_t AliTRDCalibraFillHisto::HandlePRFtrackletV1(const AliTRDseedV1 *tracklet, Int_t nbclusters)
1465 {
1466   //
1467   // PRF width calibration
1468   // Assume a Gaussian shape: determinate the position of the three pad clusters
1469   // Fit with a straight line
1470   // Take the fitted values for all the clusters (3 or 2 pad clusters)
1471   // Fill the PRF as function of angle of the track
1472   //
1473   //
1474
1475   //printf("begin\n");
1476   ///////////////////////////////////////////
1477   // Take the parameters of the track
1478   //////////////////////////////////////////
1479   // take now the snp, tnp and tgl from the track
1480   Double_t snp = tracklet->GetSnp();             // sin dy/dx at the end of the chamber
1481   Double_t tnp = 0.0;                            // dy/dx at the end of the chamber
1482   if( TMath::Abs(snp) <  1.){
1483     tnp = snp / TMath::Sqrt((1.-snp)*(1.+snp));
1484   }
1485   Double_t tgl  = tracklet->GetTgl();           // dz/dl
1486   Double_t dzdx = tgl*TMath::Sqrt(1+tnp*tnp);   // dz/dx calculated from dz/dl
1487   // at the entrance
1488   //Double_t tnp = tracklet->GetYref(1);      // dy/dx at the entrance of the chamber
1489   //Double_t tgl = tracklet->GetZref(1);      // dz/dl at the entrance of the chamber
1490   //Double_t dzdx = tgl;                      //*TMath::Sqrt(1+tnp*tnp); // dz/dx from dz/dl
1491   // at the end with correction due to linear fit
1492   //Double_t tnp = tracklet->GetYfit(1);      // dy/dx at the end of the chamber after fit correction
1493   //Double_t tgl = tracklet->GetZfit(1);      // dz/dl at the end of the chamber after fit correction
1494
1495   ///////////////////////////////
1496   // Calculate tnp group shift
1497   ///////////////////////////////
1498   Bool_t echec   = kFALSE;
1499   Double_t shift = 0.0;
1500   //Calculate the shift in x coresponding to this tnp
1501   if(fNgroupprf != 0.0){
1502     shift      = -3.0*(fNgroupprf-1)-1.5;
1503     Double_t limithigh  = -0.2*(fNgroupprf-1);
1504     if((tnp < (-0.2*fNgroupprf)) || (tnp > (0.2*fNgroupprf))) echec = kTRUE;
1505     else{
1506       while(tnp > limithigh){
1507         limithigh += 0.2;
1508         shift += 3.0;
1509       }
1510     }
1511   }
1512   // do nothing if out of tnp range
1513   //printf("echec %d\n",(Int_t)echec);
1514   if(echec) return kFALSE;
1515
1516   ///////////////////////
1517   // Variables
1518   //////////////////////
1519
1520   Int_t nb3pc    = 0;              // number of three pads clusters used for fit
1521   // to see the difference between the fit and the 3 pad clusters position
1522   Double_t padPositions[AliTRDseedV1::kNtb];
1523   memset(padPositions, 0, AliTRDseedV1::kNtb*sizeof(Double_t));
1524   fLinearFitterTracklet->ClearPoints();
1525
1526   //printf("loop clusters \n");
1527   ////////////////////////////
1528   // loop over the clusters
1529   ////////////////////////////
1530   AliTRDcluster *cl                   = 0x0;
1531   for(int ic=0; ic<AliTRDseedV1::kNtb; ic++){
1532     // reject shared clusters on pad row
1533     if((ic+AliTRDseedV1::kNtb) < AliTRDseedV1::kNclusters) {
1534       if((cl = tracklet->GetClusters(ic+AliTRDseedV1::kNtb))) continue;
1535     }
1536     if(!(cl = tracklet->GetClusters(ic))) continue;
1537     Double_t     time  = cl->GetPadTime();
1538     if((time<=7) || (time>=21)) continue;
1539     Short_t  *signals  = cl->GetSignals();
1540     Float_t xcenter    = 0.0;
1541     Bool_t  echec1      = kTRUE;
1542
1543     /////////////////////////////////////////////////////////////
1544     // Center 3 balanced: position with the center of the pad
1545     /////////////////////////////////////////////////////////////
1546     if ((((Float_t) signals[3]) > 0.0) &&
1547         (((Float_t) signals[2]) > 0.0) &&
1548         (((Float_t) signals[4]) > 0.0)) {
1549       echec1 = kFALSE;
1550       // Security if the denomiateur is 0
1551       if ((((Float_t) (((Float_t) signals[3]) * ((Float_t) signals[3]))) /
1552            ((Float_t) (((Float_t) signals[2]) * ((Float_t) signals[4])))) != 1.0) {
1553         xcenter = 0.5 * (TMath::Log((Float_t) (((Float_t) signals[4]) / ((Float_t) signals[2]))))
1554           / (TMath::Log(((Float_t) (((Float_t) signals[3]) * ((Float_t) signals[3])))
1555                         / ((Float_t) (((Float_t) signals[2]) * ((Float_t) signals[4])))));
1556       }
1557       else {
1558         echec1 = kTRUE;
1559       }
1560     }
1561     if(TMath::Abs(xcenter) > 0.5) echec1 = kTRUE;
1562     if(echec1) continue;
1563
1564     ////////////////////////////////////////////////////////
1565     //if no echec1: calculate with the position of the pad
1566     // Position of the cluster
1567     // fill the linear fitter
1568     ///////////////////////////////////////////////////////
1569     Double_t       padPosition = xcenter +  cl->GetPadCol();
1570     padPositions[ic]            = padPosition;
1571     nb3pc++;
1572     fLinearFitterTracklet->AddPoint(&time, padPosition,1);
1573
1574
1575   }//clusters loop
1576
1577   //printf("Fin loop clusters \n");
1578   //////////////////////////////
1579   // fit with a straight line
1580   /////////////////////////////
1581   if(nb3pc < 3){
1582     fLinearFitterTracklet->ClearPoints();
1583     return kFALSE;
1584   }
1585   fLinearFitterTracklet->Eval();
1586   TVectorD line(2);
1587   fLinearFitterTracklet->GetParameters(line);
1588   Float_t  pointError  = -1.0;
1589   if( fLinearFitterTracklet->GetChisquare()>=0.0) {
1590   pointError  =  TMath::Sqrt( fLinearFitterTracklet->GetChisquare()/(nb3pc-2));
1591   }
1592   fLinearFitterTracklet->ClearPoints();
1593
1594   //printf("PRF second loop \n");
1595   ////////////////////////////////////////////////
1596   // Fill the PRF: Second loop over clusters
1597   //////////////////////////////////////////////
1598   for(int ic=0; ic<AliTRDseedV1::kNtb; ic++){
1599     // reject shared clusters on pad row
1600     if(((ic+AliTRDseedV1::kNtb) < AliTRDseedV1::kNclusters) && (cl = tracklet->GetClusters(ic+AliTRDseedV1::kNtb))) continue;
1601     //
1602     if(!(cl = tracklet->GetClusters(ic))) continue;
1603
1604     Short_t  *signals      = cl->GetSignals();              // signal
1605     Double_t     time      = cl->GetPadTime();         // time bin
1606     Float_t padPosTracklet = line[0]+line[1]*time;          // reconstruct position from fit
1607     Float_t padPos         = cl->GetPadCol();               // middle pad
1608     Double_t dpad          = padPosTracklet - padPos;       // reconstruct position relative to middle pad from fit
1609     Float_t ycenter        = 0.0;                           // relative center charge
1610     Float_t ymin           = 0.0;                           // relative left charge
1611     Float_t ymax           = 0.0;                           // relative right charge
1612
1613     ////////////////////////////////////////////////////////////////
1614     // Calculate ycenter, ymin and ymax even for two pad clusters
1615     ////////////////////////////////////////////////////////////////
1616     if(((((Float_t) signals[3]) > 0.0) && (((Float_t) signals[2]) > 0.0)) ||
1617        ((((Float_t) signals[3]) > 0.0) && (((Float_t) signals[4]) > 0.0))){
1618       Float_t sum     = ((Float_t) signals[2]) + ((Float_t) signals[3]) + ((Float_t) signals[4]);
1619       if(sum > 0.0) ycenter = ((Float_t) signals[3])/ sum;
1620       if(sum > 0.0) ymin    = ((Float_t) signals[2])/ sum;
1621       if(sum > 0.0) ymax    = ((Float_t) signals[4])/ sum;
1622     }
1623
1624     /////////////////////////
1625     // Calibration group
1626     ////////////////////////
1627     Int_t     rowcl        = cl->GetPadRow();                           // row of cluster
1628     Int_t     colcl        = cl->GetPadCol();                           // col of cluster
1629     Int_t     grouplocal   = CalculateCalibrationGroup(2,rowcl,colcl);  // calcul the corresponding group
1630     Int_t     caligroup    = fCalibraMode->GetXbins(2)+ grouplocal;     // calcul the corresponding group
1631     Float_t   xcl          = cl->GetY();                                // y cluster
1632     Float_t   qcl          = cl->GetQ();                                // charge cluster
1633     Int_t     layer        = GetLayer(fDetectorPreviousTrack);          // layer
1634     Int_t     stack        = GetStack(fDetectorPreviousTrack);          // stack
1635     Double_t  xdiff        = dpad;                                      // reconstructed position constant
1636     Double_t  x            = dpad;                                      // reconstructed position moved
1637     Float_t   ep           = pointError;                                // error of fit
1638     Float_t   signal1      = (Float_t)signals[1];                       // signal at the border
1639     Float_t   signal3      = (Float_t)signals[3];                       // signal
1640     Float_t   signal2      = (Float_t)signals[2];                       // signal
1641     Float_t   signal4      = (Float_t)signals[4];                       // signal
1642     Float_t   signal5      = (Float_t)signals[5];                       // signal at the border
1643
1644
1645
1646     /////////////////////
1647     // Debug stuff
1648     ////////////////////
1649
1650     if(fDebugLevel > 0){
1651       if ( !fDebugStreamer ) {
1652         //debug stream
1653         TDirectory *backup = gDirectory;
1654         fDebugStreamer = new TTreeSRedirector("TRDdebugCalibraFill.root");
1655         if ( backup ) backup->cd();  //we don't want to be cd'd to the debug streamer
1656       }
1657
1658       x = xdiff;
1659       Int_t type=0;
1660       Float_t y = ycenter;
1661       (* fDebugStreamer) << "HandlePRFtrackletV1"<<
1662         "caligroup="<<caligroup<<
1663         "detector="<<fDetectorPreviousTrack<<
1664         "layer="<<layer<<
1665         "stack="<<stack<<
1666         "npoints="<<nbclusters<<
1667         "Np="<<nb3pc<<
1668         "ep="<<ep<<
1669         "type="<<type<<
1670         "snp="<<snp<<
1671         "tnp="<<tnp<<
1672         "tgl="<<tgl<<
1673         "dzdx="<<dzdx<<
1674         "padPos="<<padPos<<
1675         "padPosition="<<padPositions[ic]<<
1676         "padPosTracklet="<<padPosTracklet<<
1677         "x="<<x<<
1678         "y="<<y<<
1679         "xcl="<<xcl<<
1680         "qcl="<<qcl<<
1681         "signal1="<<signal1<<
1682         "signal2="<<signal2<<
1683         "signal3="<<signal3<<
1684         "signal4="<<signal4<<
1685         "signal5="<<signal5<<
1686         "time="<<time<<
1687         "\n";
1688       x=-(xdiff+1);
1689       y = ymin;
1690       type=-1;
1691       (* fDebugStreamer) << "HandlePRFtrackletV1"<<
1692         "caligroup="<<caligroup<<
1693         "detector="<<fDetectorPreviousTrack<<
1694         "layer="<<layer<<
1695         "stack="<<stack<<
1696         "npoints="<<nbclusters<<
1697         "Np="<<nb3pc<<
1698         "ep="<<ep<<
1699         "type="<<type<<
1700         "snp="<<snp<<
1701         "tnp="<<tnp<<
1702         "tgl="<<tgl<<
1703         "dzdx="<<dzdx<<
1704         "padPos="<<padPos<<
1705         "padPosition="<<padPositions[ic]<<
1706         "padPosTracklet="<<padPosTracklet<<
1707         "x="<<x<<
1708         "y="<<y<<
1709         "xcl="<<xcl<<
1710         "qcl="<<qcl<<
1711         "signal1="<<signal1<<
1712         "signal2="<<signal2<<
1713         "signal3="<<signal3<<
1714         "signal4="<<signal4<<
1715         "signal5="<<signal5<<
1716         "time="<<time<<
1717         "\n";
1718       x=1-xdiff;
1719       y = ymax;
1720       type=1;
1721       (* fDebugStreamer) << "HandlePRFtrackletV1"<<
1722         "caligroup="<<caligroup<<
1723         "detector="<<fDetectorPreviousTrack<<
1724         "layer="<<layer<<
1725         "stack="<<stack<<
1726         "npoints="<<nbclusters<<
1727         "Np="<<nb3pc<<
1728         "ep="<<ep<<
1729         "type="<<type<<
1730         "snp="<<snp<<
1731         "tnp="<<tnp<<
1732         "tgl="<<tgl<<
1733         "dzdx="<<dzdx<<
1734         "padPos="<<padPos<<
1735         "padPosition="<<padPositions[ic]<<
1736         "padPosTracklet="<<padPosTracklet<<
1737         "x="<<x<<
1738         "y="<<y<<
1739         "xcl="<<xcl<<
1740         "qcl="<<qcl<<
1741         "signal1="<<signal1<<
1742         "signal2="<<signal2<<
1743         "signal3="<<signal3<<
1744         "signal4="<<signal4<<
1745         "signal5="<<signal5<<
1746         "time="<<time<<
1747         "\n";
1748
1749     }
1750
1751     /////////////////////
1752     // Cuts quality
1753     /////////////////////
1754     if(nbclusters < fNumberClusters) continue;
1755     if(nbclusters > fNumberClustersf) continue;
1756     if(nb3pc <= 5) continue;
1757     if((time >= 21) || (time < 7)) continue;
1758     if(TMath::Abs(qcl) < 80) continue;
1759     if( TMath::Abs(snp) >  1.) continue;
1760
1761
1762     ////////////////////////
1763     // Fill the histos
1764     ///////////////////////
1765     if (fHisto2d) {
1766       if(TMath::Abs(dpad) < 1.5) {
1767         fPRF2d->Fill(shift+dpad,(caligroup+0.5),ycenter);
1768         fPRF2d->Fill(shift-dpad,(caligroup+0.5),ycenter);
1769         //printf("place %f, ycenter %f\n",(shift+dpad),ycenter);
1770       }
1771       if((ymin > 0.0) && (TMath::Abs(dpad+1.0) < 1.5)) {
1772         fPRF2d->Fill(shift-(dpad+1.0),(caligroup+0.5),ymin);
1773         fPRF2d->Fill(shift+(dpad+1.0),(caligroup+0.5),ymin);
1774       }
1775       if((ymax > 0.0) && (TMath::Abs(dpad-1.0) < 1.5)) {
1776         fPRF2d->Fill(shift+1.0-dpad,(caligroup+0.5),ymax);
1777         fPRF2d->Fill(shift-1.0+dpad,(caligroup+0.5),ymax);
1778       }
1779     }
1780     // vector method
1781     if (fVector2d) {
1782       if(TMath::Abs(dpad) < 1.5) {
1783         fCalibraVector->UpdateVectorPRF(fDetectorPreviousTrack,grouplocal,shift+dpad,ycenter);
1784         fCalibraVector->UpdateVectorPRF(fDetectorPreviousTrack,grouplocal,shift-dpad,ycenter);
1785       }
1786       if((ymin > 0.0) && (TMath::Abs(dpad+1.0) < 1.5)) {
1787         fCalibraVector->UpdateVectorPRF(fDetectorPreviousTrack,grouplocal,shift-(dpad+1.0),ymin);
1788         fCalibraVector->UpdateVectorPRF(fDetectorPreviousTrack,grouplocal,shift+(dpad+1.0),ymin);
1789       }
1790       if((ymax > 0.0)  && (TMath::Abs(dpad-1.0) < 1.5)) {
1791         fCalibraVector->UpdateVectorPRF(fDetectorPreviousTrack,grouplocal,shift+1.0-dpad,ymax);
1792         fCalibraVector->UpdateVectorPRF(fDetectorPreviousTrack,grouplocal,shift-1.0+dpad,ymax);
1793       }
1794     }
1795   } // second loop over clusters
1796
1797
1798   return kTRUE;
1799 }
1800 ///////////////////////////////////////////////////////////////////////////////////////
1801 // Pad row col stuff: see if masked or not
1802 ///////////////////////////////////////////////////////////////////////////////////////
1803 //_____________________________________________________________________________
1804 void AliTRDCalibraFillHisto::CheckGoodTrackletV1(const AliTRDcluster *cl)
1805 {
1806   //
1807   // See if we are not near a masked pad
1808   //
1809
1810   if(cl->IsMasked()) fGoodTracklet = kFALSE;
1811
1812
1813 }
1814 //_____________________________________________________________________________
1815 void AliTRDCalibraFillHisto::CheckGoodTrackletV0(const Int_t detector,const Int_t row,const Int_t col)
1816 {
1817   //
1818   // See if we are not near a masked pad
1819   //
1820
1821   if (!IsPadOn(detector, col, row)) {
1822     fGoodTracklet = kFALSE;
1823   }
1824
1825   if (col > 0) {
1826     if (!IsPadOn(detector, col-1, row)) {
1827       fGoodTracklet = kFALSE;
1828     }
1829   }
1830
1831   if (col < 143) {
1832     if (!IsPadOn(detector, col+1, row)) {
1833       fGoodTracklet = kFALSE;
1834     }
1835   }
1836
1837 }
1838 //_____________________________________________________________________________
1839 Bool_t AliTRDCalibraFillHisto::IsPadOn(Int_t detector, Int_t row, Int_t col) const
1840 {
1841   //
1842   // Look in the choosen database if the pad is On.
1843   // If no the track will be "not good"
1844   //
1845
1846   // Get the parameter object
1847   AliTRDcalibDB *cal = AliTRDcalibDB::Instance();
1848   if (!cal) {
1849     AliInfo("Could not get calibDB");
1850     return kFALSE;
1851   }
1852
1853   if (!cal->IsChamberInstalled(detector)     ||
1854        cal->IsChamberMasked(detector)        ||
1855        cal->IsPadMasked(detector,col,row)) {
1856     return kFALSE;
1857   }
1858   else {
1859     return kTRUE;
1860   }
1861
1862 }
1863 ///////////////////////////////////////////////////////////////////////////////////////
1864 // Calibration groups: calculate the number of groups, localise...
1865 ////////////////////////////////////////////////////////////////////////////////////////
1866 //_____________________________________________________________________________
1867 Int_t AliTRDCalibraFillHisto::CalculateCalibrationGroup(Int_t i, Int_t row, Int_t col) const
1868 {
1869   //
1870   // Calculate the calibration group number for i
1871   //
1872
1873   // Row of the cluster and position in the pad groups
1874   Int_t posr = 0;
1875   if (fCalibraMode->GetNnZ(i) != 0) {
1876     posr = (Int_t) row / fCalibraMode->GetNnZ(i);
1877   }
1878
1879
1880   // Col of the cluster and position in the pad groups
1881   Int_t posc = 0;
1882   if (fCalibraMode->GetNnRphi(i) != 0) {
1883     posc = (Int_t) col / fCalibraMode->GetNnRphi(i);
1884   }
1885
1886   return posc*fCalibraMode->GetNfragZ(i)+posr;
1887
1888 }
1889 //____________________________________________________________________________________
1890 Int_t AliTRDCalibraFillHisto::CalculateTotalNumberOfBins(Int_t i)
1891 {
1892   //
1893   // Calculate the total number of calibration groups
1894   //
1895
1896   Int_t ntotal = 0;
1897
1898   // All together
1899   if((fCalibraMode->GetNz(i)==100) && (fCalibraMode->GetNrphi(i)==100)){
1900     ntotal = 1;
1901     AliInfo(Form("Total number of Xbins: %d for i %d",ntotal,i));
1902     return ntotal;
1903   }
1904
1905   // Per Supermodule
1906   if((fCalibraMode->GetNz(i)==10) && (fCalibraMode->GetNrphi(i)==10)){
1907     ntotal = 18;
1908     AliInfo(Form("Total number of Xbins: %d for i %d",ntotal,i));
1909     return ntotal;
1910   }
1911
1912   // More
1913   fCalibraMode->ModePadCalibration(2,i);
1914   fCalibraMode->ModePadFragmentation(0,2,0,i);
1915   fCalibraMode->SetDetChamb2(i);
1916   ntotal += 6 * 18 * fCalibraMode->GetDetChamb2(i);
1917   fCalibraMode->ModePadCalibration(0,i);
1918   fCalibraMode->ModePadFragmentation(0,0,0,i);
1919   fCalibraMode->SetDetChamb0(i);
1920   ntotal += 6 * 4 * 18 * fCalibraMode->GetDetChamb0(i);
1921   AliInfo(Form("Total number of Xbins: %d for i %d",ntotal,i));
1922   return ntotal;
1923
1924 }
1925 //_____________________________________________________________________________
1926 void AliTRDCalibraFillHisto::SetNz(Int_t i, Short_t Nz)
1927 {
1928   //
1929   // Set the mode of calibration group in the z direction for the parameter i
1930   //
1931
1932   if ((Nz >= 0) &&
1933       (Nz <  5)) {
1934     fCalibraMode->SetNz(i, Nz);
1935   }
1936   else {
1937     AliInfo("You have to choose between 0 and 4");
1938   }
1939
1940 }
1941 //_____________________________________________________________________________
1942 void AliTRDCalibraFillHisto::SetNrphi(Int_t i, Short_t Nrphi)
1943 {
1944   //
1945   // Set the mode of calibration group in the rphi direction for the parameter i
1946   //
1947
1948   if ((Nrphi >= 0) &&
1949       (Nrphi <  7)) {
1950     fCalibraMode->SetNrphi(i ,Nrphi);
1951   }
1952   else {
1953     AliInfo("You have to choose between 0 and 6");
1954   }
1955
1956 }
1957
1958 //_____________________________________________________________________________
1959 void AliTRDCalibraFillHisto::SetAllTogether(Int_t i)
1960 {
1961   //
1962   // Set the mode of calibration group all together
1963   //
1964   if(fVector2d == kTRUE) {
1965     AliInfo("Can not work with the vector method");
1966     return;
1967   }
1968   fCalibraMode->SetAllTogether(i);
1969
1970 }
1971
1972 //_____________________________________________________________________________
1973 void AliTRDCalibraFillHisto::SetPerSuperModule(Int_t i)
1974 {
1975   //
1976   // Set the mode of calibration group per supermodule
1977   //
1978   if(fVector2d == kTRUE) {
1979     AliInfo("Can not work with the vector method");
1980     return;
1981   }
1982   fCalibraMode->SetPerSuperModule(i);
1983
1984 }
1985
1986 //____________Set the pad calibration variables for the detector_______________
1987 Bool_t AliTRDCalibraFillHisto::LocalisationDetectorXbins(Int_t detector)
1988 {
1989   //
1990   // For the detector calcul the first Xbins and set the number of row
1991   // and col pads per calibration groups, the number of calibration
1992   // groups in the detector.
1993   //
1994
1995   // first Xbins of the detector
1996   if (fCH2dOn) {
1997     fCalibraMode->CalculXBins(detector,0);
1998   }
1999   if (fPH2dOn) {
2000     fCalibraMode->CalculXBins(detector,1);
2001   }
2002   if (fPRF2dOn) {
2003     fCalibraMode->CalculXBins(detector,2);
2004   }
2005
2006   // fragmentation of idect
2007   for (Int_t i = 0; i < 3; i++) {
2008     fCalibraMode->ModePadCalibration((Int_t) GetStack(detector),i);
2009     fCalibraMode->ModePadFragmentation((Int_t) GetLayer(detector)
2010                        , (Int_t) GetStack(detector)
2011                        , (Int_t) GetSector(detector),i);
2012   }
2013
2014   return kTRUE;
2015
2016 }
2017 //_____________________________________________________________________________
2018 void AliTRDCalibraFillHisto::SetNumberGroupsPRF(Short_t numberGroupsPRF)
2019 {
2020   //
2021   // Should be between 0 and 6
2022   //
2023
2024   if ((numberGroupsPRF < 0) || (numberGroupsPRF > 6)) {
2025     AliInfo("The number of groups must be between 0 and 6!");
2026   }
2027   else {
2028     fNgroupprf = numberGroupsPRF;
2029   }
2030
2031 }
2032 ///////////////////////////////////////////////////////////////////////////////////////////
2033 // Per tracklet: store or reset the info, fill the histos with the info
2034 //////////////////////////////////////////////////////////////////////////////////////////
2035 //_____________________________________________________________________________
2036 void AliTRDCalibraFillHisto::StoreInfoCHPHtrack(const AliTRDcluster *cl,const Double_t dqdl,const Int_t *group,const Int_t row,const Int_t col,const AliTRDcluster *cls)
2037 {
2038   //
2039   // Store the infos in fAmpTotal, fPHPlace and fPHValue
2040   // Correct from the gain correction before
2041   // cls is shared cluster if any
2042   //
2043
2044   //printf("StoreInfoCHPHtrack\n");
2045
2046   // time bin of the cluster not corrected
2047   Int_t    time     = cl->GetPadTime();
2048   Float_t  charge   = TMath::Abs(cl->GetQ());
2049   if(cls) {
2050     charge += TMath::Abs(cls->GetQ());
2051     //printf("AliTRDCalibraFillHisto::Add the cluster charge");
2052   }
2053
2054   //printf("Store::time %d, amplitude %f\n",time,dqdl);
2055
2056   //Correct for the gain coefficient used in the database for reconstruction
2057   Float_t correctthegain = 1.0;
2058   if(fIsHLT) correctthegain = fCalDetGain->GetValue(fDetectorPreviousTrack);
2059   else correctthegain = fCalDetGain->GetValue(fDetectorPreviousTrack)*fCalROCGain->GetValue(col,row);
2060   Float_t correction    = 1.0;
2061   Float_t normalisation = 6.67;
2062   // we divide with gain in AliTRDclusterizer::Transform...
2063   if( correctthegain > 0 ) normalisation /= correctthegain;
2064
2065
2066   // take dd/dl corrected from the angle of the track
2067   correction = dqdl / (normalisation);
2068
2069
2070   // Fill the fAmpTotal with the charge
2071   if (fCH2dOn) {
2072     if((!fLimitChargeIntegration) || (cl->IsInChamber())) {
2073       //printf("Store::group %d, amplitude %f\n",group[0],correction);
2074       fAmpTotal[(Int_t) group[0]] += correction;
2075     }
2076   }
2077
2078   // Fill the fPHPlace and value
2079   if (fPH2dOn) {
2080     fPHPlace[time] = group[1];
2081     fPHValue[time] = charge;
2082   }
2083
2084 }
2085 //____________Offine tracking in the AliTRDtracker_____________________________
2086 void AliTRDCalibraFillHisto::ResetfVariablestracklet()
2087 {
2088   //
2089   // Reset values per tracklet
2090   //
2091
2092   //Reset good tracklet
2093   fGoodTracklet = kTRUE;
2094
2095   // Reset the fPHValue
2096   if (fPH2dOn) {
2097     //Reset the fPHValue and fPHPlace
2098     for (Int_t k = 0; k < fTimeMax; k++) {
2099       fPHValue[k] = 0.0;
2100       fPHPlace[k] = -1;
2101     }
2102   }
2103
2104   // Reset the fAmpTotal where we put value
2105   if (fCH2dOn) {
2106     for (Int_t k = 0; k < fCalibraMode->GetNfragZ(0)*fCalibraMode->GetNfragRphi(0); k++) {
2107       fAmpTotal[k] = 0.0;
2108     }
2109   }
2110 }
2111 //____________Offine tracking in the AliTRDtracker_____________________________
2112 void AliTRDCalibraFillHisto::FillTheInfoOfTheTrackCH(Int_t nbclusters)
2113 {
2114   //
2115   // For the offline tracking or mcm tracklets
2116   // This function will be called in the functions UpdateHistogram...
2117   // to fill the info of a track for the relativ gain calibration
2118   //
2119
2120   Int_t nb            =  0;   // Nombre de zones traversees
2121   Int_t fd            = -1;   // Premiere zone non nulle
2122   Float_t totalcharge = 0.0;  // Total charge for the supermodule histo
2123
2124   //printf("CH2d nbclusters %d, fNumberClusters %d, fNumberClustersf %d\n",nbclusters,fNumberClusters,fNumberClustersf);
2125
2126   if(nbclusters < fNumberClusters) return;
2127   if(nbclusters > fNumberClustersf) return;
2128
2129
2130   // Normalize with the number of clusters
2131   Double_t normalizeCst = fRelativeScale;
2132   if(fNormalizeNbOfCluster) normalizeCst = normalizeCst*nbclusters;
2133
2134   //printf("Number of groups in one detector %d\n",fCalibraMode->GetNfragZ(0)*fCalibraMode->GetNfragRphi(0));
2135
2136   // See if the track goes through different zones
2137   for (Int_t k = 0; k < fCalibraMode->GetNfragZ(0)*fCalibraMode->GetNfragRphi(0); k++) {
2138     //printf("fAmpTotal %f for %d\n",fAmpTotal[k],k);
2139     if (fAmpTotal[k] > 0.0) {
2140       totalcharge += fAmpTotal[k];
2141       nb++;
2142       if (nb == 1) {
2143         fd = k;
2144       }
2145     }
2146   }
2147
2148   //printf("CH2d: nb %d, fd %d, calibration group %d, amplitude %f, detector %d\n",nb,fd,fCalibraMode->GetXbins(0),fAmpTotal[fd]/normalizeCst,fDetectorPreviousTrack);
2149
2150   switch (nb)
2151     {
2152     case 1:
2153       fNumberUsedCh[0]++;
2154       fEntriesCH[fCalibraMode->GetXbins(0)+fd]++;
2155       if (fHisto2d) {
2156         FillCH2d(fCalibraMode->GetXbins(0)+fd,fAmpTotal[fd]/normalizeCst);
2157         //fCH2d->Fill(fAmpTotal[fd]/normalizeCst,fCalibraMode->GetXbins(0)+fd+0.5);
2158       }
2159       if (fVector2d) {
2160         fCalibraVector->UpdateVectorCH(fDetectorPreviousTrack,fd,fAmpTotal[fd]/normalizeCst);
2161       }
2162       break;
2163     case 2:
2164       if ((fAmpTotal[fd]   > 0.0) &&
2165           (fAmpTotal[fd+1] > 0.0)) {
2166         // One of the two very big
2167         if (fAmpTotal[fd] > fProcent*fAmpTotal[fd+1]) {
2168           if (fHisto2d) {
2169             FillCH2d(fCalibraMode->GetXbins(0)+fd,fAmpTotal[fd]/normalizeCst);
2170             //fCH2d->Fill(fAmpTotal[fd]/normalizeCst,fCalibraMode->GetXbins(0)+fd+0.5);
2171           }
2172           if (fVector2d) {
2173             fCalibraVector->UpdateVectorCH(fDetectorPreviousTrack,fd,fAmpTotal[fd]/normalizeCst);
2174           }
2175           fNumberUsedCh[1]++;
2176           fEntriesCH[fCalibraMode->GetXbins(0)+fd]++;
2177         }
2178         if (fAmpTotal[fd+1] > fProcent*fAmpTotal[fd]) {
2179           if (fHisto2d) {
2180             FillCH2d(fCalibraMode->GetXbins(0)+fd+1,fAmpTotal[fd+1]/normalizeCst);
2181             //fCH2d->Fill(fAmpTotal[fd+1]/normalizeCst,fCalibraMode->GetXbins(0)+fd+1.5);
2182           }
2183           if (fVector2d) {
2184             fCalibraVector->UpdateVectorCH(fDetectorPreviousTrack,fd+1,fAmpTotal[fd+1]/normalizeCst);
2185           }
2186           fNumberUsedCh[1]++;
2187           fEntriesCH[fCalibraMode->GetXbins(0)+fd+1]++;
2188         }
2189       }
2190       if (fCalibraMode->GetNfragZ(0) > 1) {
2191         if (fAmpTotal[fd] > 0.0) {
2192           if ((fd+fCalibraMode->GetNfragZ(0)) < (fCalibraMode->GetNfragZ(0)*fCalibraMode->GetNfragRphi(0))) {
2193             if (fAmpTotal[fd+fCalibraMode->GetNfragZ(0)] > 0.0) {
2194               // One of the two very big
2195               if (fAmpTotal[fd] > fProcent*fAmpTotal[fd+fCalibraMode->GetNfragZ(0)]) {
2196                 if (fHisto2d) {
2197                   FillCH2d(fCalibraMode->GetXbins(0)+fd,fAmpTotal[fd]/normalizeCst);
2198                   //fCH2d->Fill(fAmpTotal[fd]/normalizeCst,fCalibraMode->GetXbins(0)+fd+0.5);
2199                 }
2200                 if (fVector2d) {
2201                   fCalibraVector->UpdateVectorCH(fDetectorPreviousTrack,fd,fAmpTotal[fd]/normalizeCst);
2202                 }
2203                 fNumberUsedCh[1]++;
2204                 fEntriesCH[fCalibraMode->GetXbins(0)+fd]++;
2205               }
2206               if (fAmpTotal[fd+fCalibraMode->GetNfragZ(0)] > fProcent*fAmpTotal[fd]) {
2207                 if (fHisto2d) {
2208                   FillCH2d(fCalibraMode->GetXbins(0)+fd+fCalibraMode->GetNfragZ(0),fAmpTotal[fd+fCalibraMode->GetNfragZ(0)]/normalizeCst);
2209                   //fCH2d->Fill(fAmpTotal[fd+fCalibraMode->GetNfragZ(0)]/normalizeCst,fCalibraMode->GetXbins(0)+fd+fCalibraMode->GetNfragZ(0)+0.5);
2210                 }
2211                 fNumberUsedCh[1]++;
2212                 fEntriesCH[fCalibraMode->GetXbins(0)+fd+fCalibraMode->GetNfragZ(0)]++;
2213                 if (fVector2d) {
2214                   fCalibraVector->UpdateVectorCH(fDetectorPreviousTrack,fd+fCalibraMode->GetNfragZ(0),fAmpTotal[fd+fCalibraMode->GetNfragZ(0)]/normalizeCst);
2215                 }
2216               }
2217             }
2218           }
2219         }
2220       }
2221       break;
2222     default: break;
2223     }
2224 }
2225 //____________Offine tracking in the AliTRDtracker_____________________________
2226 void AliTRDCalibraFillHisto::FillTheInfoOfTheTrackPH()
2227 {
2228   //
2229   // For the offline tracking or mcm tracklets
2230   // This function will be called in the functions UpdateHistogram...
2231   // to fill the info of a track for the drift velocity  calibration
2232   //
2233
2234   Int_t nb  =  1; // Nombre de zones traversees 1, 2 ou plus de 3
2235   Int_t fd1 = -1; // Premiere zone non nulle
2236   Int_t fd2 = -1; // Deuxieme zone non nulle
2237   Int_t k1  = -1; // Debut de la premiere zone
2238   Int_t k2  = -1; // Debut de la seconde zone
2239   Int_t nbclusters = 0; // number of clusters
2240
2241
2242
2243   // See if the track goes through different zones
2244   for (Int_t k = 0; k < fTimeMax; k++) {
2245     if (fPHValue[k] > 0.0) {
2246       nbclusters++;
2247       if (fd1 == -1) {
2248         fd1 = fPHPlace[k];
2249         k1  = k;
2250       }
2251       if (fPHPlace[k] != fd1) {
2252         if (fd2 == -1) {
2253           k2  = k;
2254           fd2 = fPHPlace[k];
2255           nb  = 2;
2256         }
2257         if (fPHPlace[k] != fd2) {
2258           nb = 3;
2259         }
2260       }
2261     }
2262   }
2263
2264   // See if noise before and after
2265   if(fMaxCluster > 0) {
2266     if(fPHValue[0] > fMaxCluster) return;
2267     if(fTimeMax > fNbMaxCluster) {
2268       for(Int_t k = (fTimeMax-fNbMaxCluster); k < fTimeMax; k++){
2269         if(fPHValue[k] > fMaxCluster) return;
2270       }
2271     }
2272   }
2273
2274   //printf("nbclusters %d, low limit %d, high limit %d\n",nbclusters,fNumberClusters,fNumberClustersf);
2275
2276   if(nbclusters < fNumberClusters) return;
2277   if(nbclusters > fNumberClustersf) return;
2278
2279   switch(nb)
2280     {
2281     case 1:
2282       fNumberUsedPh[0]++;
2283       for (Int_t i = 0; i < fTimeMax; i++) {
2284         if (fHisto2d) {
2285           if(fFillWithZero) fPH2d->Fill((Float_t) i/fSf,(fCalibraMode->GetXbins(1)+fd1)+0.5,(Float_t) fPHValue[i]);
2286           else {
2287             if(((Float_t) fPHValue[i] > 0.0)) fPH2d->Fill((Float_t) i/fSf,(fCalibraMode->GetXbins(1)+fd1)+0.5,(Float_t) fPHValue[i]);
2288               }
2289           //printf("Fill the time bin %d with %f\n",i,fPHValue[i]);
2290         }
2291         if (fVector2d) {
2292           if(fFillWithZero) fCalibraVector->UpdateVectorPH(fDetectorPreviousTrack,fd1,i,fPHValue[i]);
2293           else {
2294             if(((Float_t) fPHValue[i] > 0.0)) fCalibraVector->UpdateVectorPH(fDetectorPreviousTrack,fd1,i,fPHValue[i]);
2295           }
2296         }
2297       }
2298       break;
2299     case 2:
2300       if ((fd1 == fd2+1) ||
2301           (fd2 == fd1+1)) {
2302         // One of the two fast all the think
2303         if (k2 > (k1+fDifference)) {
2304           //we choose to fill the fd1 with all the values
2305           fNumberUsedPh[1]++;
2306           for (Int_t i = 0; i < fTimeMax; i++) {
2307             if (fHisto2d) {
2308               if(fFillWithZero) fPH2d->Fill((Float_t) i/fSf,(fCalibraMode->GetXbins(1)+fd1)+0.5,(Float_t) fPHValue[i]);
2309               else {
2310                 if(((Float_t) fPHValue[i] > 0.0)) fPH2d->Fill((Float_t) i/fSf,(fCalibraMode->GetXbins(1)+fd1)+0.5,(Float_t) fPHValue[i]);
2311                   }
2312             }
2313             if (fVector2d) {
2314               if(fFillWithZero) fCalibraVector->UpdateVectorPH(fDetectorPreviousTrack,fd1,i,fPHValue[i]);
2315               else {
2316                 if(((Float_t) fPHValue[i] > 0.0)) fCalibraVector->UpdateVectorPH(fDetectorPreviousTrack,fd1,i,fPHValue[i]);
2317                   }
2318             }
2319           }
2320         }
2321         if ((k2+fDifference) < fTimeMax) {
2322           //we choose to fill the fd2 with all the values
2323           fNumberUsedPh[1]++;
2324           for (Int_t i = 0; i < fTimeMax; i++) {
2325             if (fHisto2d) {
2326               if(fFillWithZero) fPH2d->Fill((Float_t) i/fSf,(fCalibraMode->GetXbins(1)+fd2)+0.5,(Float_t) fPHValue[i]);
2327               else {
2328                 if(((Float_t) fPHValue[i] > 0.0)) fPH2d->Fill((Float_t) i/fSf,(fCalibraMode->GetXbins(1)+fd2)+0.5,(Float_t) fPHValue[i]);
2329               }
2330             }
2331           if (fVector2d) {
2332             if(fFillWithZero) fCalibraVector->UpdateVectorPH(fDetectorPreviousTrack,fd2,i,fPHValue[i]);
2333             else {
2334               if(((Float_t) fPHValue[i] > 0.0)) fCalibraVector->UpdateVectorPH(fDetectorPreviousTrack,fd2,i,fPHValue[i]);
2335             }
2336           }
2337           }
2338         }
2339       }
2340       // Two zones voisines sinon rien!
2341       if (fCalibraMode->GetNfragZ(1) > 1) {
2342         // Case 2
2343         if ((fd1+fCalibraMode->GetNfragZ(1)) < (fCalibraMode->GetNfragZ(1)*fCalibraMode->GetNfragRphi(1))) {
2344           if (fd2 == (fd1+fCalibraMode->GetNfragZ(1))) {
2345             // One of the two fast all the think
2346             if (k2 > (k1+fDifference)) {
2347               //we choose to fill the fd1 with all the values
2348               fNumberUsedPh[1]++;
2349               for (Int_t i = 0; i < fTimeMax; i++) {
2350                 if (fHisto2d) {
2351                   if(fFillWithZero) fPH2d->Fill((Float_t) i/fSf,(fCalibraMode->GetXbins(1)+fd1)+0.5,(Float_t) fPHValue[i]);
2352                   else {
2353                     if(((Float_t) fPHValue[i] > 0.0)) fPH2d->Fill((Float_t) i/fSf,(fCalibraMode->GetXbins(1)+fd1)+0.5,(Float_t) fPHValue[i]);
2354                   }
2355                 }
2356                 if (fVector2d) {
2357                   if(fFillWithZero) fCalibraVector->UpdateVectorPH(fDetectorPreviousTrack,fd1,i,fPHValue[i]);
2358                   else {
2359                     if(((Float_t) fPHValue[i] > 0.0)) fCalibraVector->UpdateVectorPH(fDetectorPreviousTrack,fd1,i,fPHValue[i]);
2360                   }
2361                 }
2362               }
2363             }
2364             if ((k2+fDifference) < fTimeMax) {
2365               //we choose to fill the fd2 with all the values
2366               fNumberUsedPh[1]++;
2367               for (Int_t i = 0; i < fTimeMax; i++) {
2368                 if (fHisto2d) {
2369                   if(fFillWithZero) fPH2d->Fill((Float_t) i/fSf,(fCalibraMode->GetXbins(1)+fd2)+0.5,(Float_t) fPHValue[i]);
2370                   else {
2371                     if(((Float_t) fPHValue[i] > 0.0)) fPH2d->Fill((Float_t) i/fSf,(fCalibraMode->GetXbins(1)+fd2)+0.5,(Float_t) fPHValue[i]);
2372                   }
2373                 }
2374                 if (fVector2d) {
2375                   if(fFillWithZero) fCalibraVector->UpdateVectorPH(fDetectorPreviousTrack,fd2,i,fPHValue[i]);
2376                   else {
2377                     if(((Float_t) fPHValue[i] > 0.0)) fCalibraVector->UpdateVectorPH(fDetectorPreviousTrack,fd2,i,fPHValue[i]);
2378                   }
2379                 }
2380               }
2381             }
2382           }
2383         }
2384         // Two zones voisines sinon rien!
2385         // Case 3
2386         if ((fd1 - fCalibraMode->GetNfragZ(1)) >= 0) {
2387           if (fd2 == (fd1 - fCalibraMode->GetNfragZ(1))) {
2388             // One of the two fast all the think
2389             if (k2 > (k1 + fDifference)) {
2390               //we choose to fill the fd1 with all the values
2391               fNumberUsedPh[1]++;
2392               for (Int_t i = 0; i < fTimeMax; i++) {
2393                 if (fHisto2d) {
2394                   if(fFillWithZero) fPH2d->Fill((Float_t) i/fSf,(fCalibraMode->GetXbins(1)+fd1)+0.5,(Float_t) fPHValue[i]);
2395                   else {
2396                     if(((Float_t) fPHValue[i] > 0.0)) fPH2d->Fill((Float_t) i/fSf,(fCalibraMode->GetXbins(1)+fd1)+0.5,(Float_t) fPHValue[i]);
2397                   }
2398                 }
2399                 if (fVector2d) {
2400                   if(fFillWithZero) fCalibraVector->UpdateVectorPH(fDetectorPreviousTrack,fd1,i,fPHValue[i]);
2401                   else {
2402                     if(((Float_t) fPHValue[i] > 0.0)) fCalibraVector->UpdateVectorPH(fDetectorPreviousTrack,fd1,i,fPHValue[i]);
2403                   }
2404                 }
2405               }
2406             }
2407             if ((k2+fDifference) < fTimeMax) {
2408               //we choose to fill the fd2 with all the values
2409               fNumberUsedPh[1]++;
2410               for (Int_t i = 0; i < fTimeMax; i++) {
2411                 if (fHisto2d) {
2412                   if(fFillWithZero) fPH2d->Fill((Float_t) i/fSf,(fCalibraMode->GetXbins(1)+fd2)+0.5,(Float_t) fPHValue[i]);
2413                   else {
2414                     if(((Float_t) fPHValue[i] > 0.0)) fPH2d->Fill((Float_t) i/fSf,(fCalibraMode->GetXbins(1)+fd2)+0.5,(Float_t) fPHValue[i]);
2415                   }
2416                 }
2417                 if (fVector2d) {
2418                   if(fFillWithZero) fCalibraVector->UpdateVectorPH(fDetectorPreviousTrack,fd2,i,fPHValue[i]);
2419                   else {
2420                     if(((Float_t) fPHValue[i] > 0.0)) fCalibraVector->UpdateVectorPH(fDetectorPreviousTrack,fd2,i,fPHValue[i]);
2421                   }
2422                 }
2423               }
2424             }
2425           }
2426         }
2427       }
2428       break;
2429     default: break;
2430     }
2431 }
2432 //////////////////////////////////////////////////////////////////////////////////////////
2433 // DAQ process functions
2434 /////////////////////////////////////////////////////////////////////////////////////////
2435 //_____________________________________________________________________
2436 Int_t AliTRDCalibraFillHisto::ProcessEventDAQ(AliTRDrawStreamBase *rawStream, Bool_t nocheck)
2437 {
2438   //
2439   // Event Processing loop - AliTRDrawStreamBase
2440   // TestBeam 2007 version
2441   // 0 timebin problem
2442   // 1 no input
2443   // 2 input
2444   // Same algorithm as TestBeam but different reader
2445   //
2446
2447   rawStream->SetSharedPadReadout(kFALSE);
2448
2449   Int_t withInput = 1;
2450
2451   Double_t phvalue[16][144][36];
2452   for(Int_t k = 0; k < 36; k++){
2453     for(Int_t j = 0; j < 16; j++){
2454       for(Int_t c = 0; c < 144; c++){
2455         phvalue[j][c][k] = 0.0;
2456       }
2457     }
2458   }
2459
2460   fDetectorPreviousTrack = -1;
2461   fMCMPrevious           = -1;
2462   fROBPrevious           = -1;
2463   Int_t nbtimebin = 0;
2464   Int_t baseline  = 10;
2465   //printf("------------Detector\n");
2466
2467   if(!nocheck){
2468
2469     fTimeMax = 0;
2470
2471     while (rawStream->Next()) {
2472
2473       Int_t idetector = rawStream->GetDet();                            //  current detector
2474       Int_t imcm      = rawStream->GetMCM();                            //  current MCM
2475       Int_t irob      = rawStream->GetROB();                            //  current ROB
2476
2477       //printf("Detector %d\n",idetector);
2478
2479       if((fDetectorPreviousTrack != idetector) && (fDetectorPreviousTrack != -1)){
2480
2481         // Fill
2482         withInput = TMath::Max(FillDAQ(phvalue),withInput);
2483
2484
2485         // reset
2486         for(Int_t k = 0; k < 36; k++){
2487           for(Int_t j = 0; j < 16; j++){
2488             for(Int_t c = 0; c < 144; c++){
2489               phvalue[j][c][k] = 0.0;
2490             }
2491           }
2492         }
2493       }
2494
2495       fDetectorPreviousTrack = idetector;
2496       fMCMPrevious           = imcm;
2497       fROBPrevious           = irob;
2498
2499       nbtimebin              = rawStream->GetNumberOfTimeBins();              //  number of time bins read from data
2500       if(nbtimebin == 0) return 0;
2501       if((fTimeMax != 0) && (nbtimebin != fTimeMax)) return 0;
2502       fTimeMax          = nbtimebin;
2503
2504       //baseline          = rawStream->GetCommonAdditive();                // common additive baseline
2505       fNumberClustersf    = fTimeMax;
2506       fNumberClusters     = (Int_t)(fNumberClustersProcent*fTimeMax);
2507
2508
2509       Int_t *signal     = rawStream->GetSignals();                       //  current ADC signal
2510       Int_t  col        = rawStream->GetCol();
2511       Int_t row         = rawStream->GetRow();
2512
2513
2514       // printf("detector %d, signal[0] %d, signal[1] %d, signal[2] %d, baseline %d\n",idetector,signal[0],signal[1],signal[2], baseline);
2515
2516
2517       for(Int_t itime = 0; itime < nbtimebin; itime++){
2518         phvalue[row][col][itime] = signal[itime]-baseline;
2519       }
2520     }
2521
2522     // fill the last one
2523     if(fDetectorPreviousTrack != -1){
2524
2525       // Fill
2526       withInput = TMath::Max(FillDAQ(phvalue),withInput);
2527
2528       // reset
2529       for(Int_t k = 0; k < 36; k++){
2530         for(Int_t j = 0; j < 16; j++){
2531           for(Int_t c = 0; c < 144; c++){
2532             phvalue[j][c][k] = 0.0;
2533           }
2534         }
2535       }
2536     }
2537
2538   }
2539   else{
2540
2541     while (rawStream->Next()) { //iddetecte
2542
2543       Int_t idetector = rawStream->GetDet();                            //  current detector
2544       Int_t imcm      = rawStream->GetMCM();                            //  current MCM
2545       Int_t irob      = rawStream->GetROB();                            //  current ROB
2546
2547       //printf("Detector %d\n",idetector);
2548
2549       if((fDetectorPreviousTrack != idetector) && (fDetectorPreviousTrack != -1)){
2550
2551         // Fill
2552         withInput = TMath::Max(FillDAQ(phvalue),withInput);
2553
2554         // reset
2555         for(Int_t k = 0; k < 36; k++){
2556           for(Int_t j = 0; j < 16; j++){
2557             for(Int_t c = 0; c < 144; c++){
2558               phvalue[j][c][k] = 0.0;
2559             }
2560           }
2561         }
2562       }
2563
2564       fDetectorPreviousTrack = idetector;
2565       fMCMPrevious           = imcm;
2566       fROBPrevious           = irob;
2567
2568       //baseline          = rawStream->GetCommonAdditive();                //  common baseline
2569
2570       fTimeMax          = rawStream->GetNumberOfTimeBins();              //  number of time bins read from data
2571       fNumberClustersf    = fTimeMax;
2572       fNumberClusters     = (Int_t)(fNumberClustersProcent*fTimeMax);
2573       Int_t *signal     = rawStream->GetSignals();                       //  current ADC signal
2574       Int_t col         = rawStream->GetCol();
2575       Int_t row         = rawStream->GetRow();
2576
2577
2578       //printf("detector %d, signal[0] %d, signal[1] %d, signal[2] %d, baseline %d\n",idetector,signal[0],signal[1],signal[2], baseline);
2579
2580       for(Int_t itime = 0; itime < fTimeMax; itime++){
2581         phvalue[row][col][itime] = signal[itime]-baseline;
2582         /*if(phvalue[row][col][itime] >= 20) {
2583                 printf("----------> phvalue[%d][%d][%d] %d  baseline %d \n",
2584                        row,
2585                        col,
2586                        itime,
2587                        signal[itime],
2588                        baseline);
2589                        }*/
2590       }
2591     }
2592
2593     // fill the last one
2594     if(fDetectorPreviousTrack != -1){
2595
2596       // Fill
2597       withInput = TMath::Max(FillDAQ(phvalue),withInput);
2598
2599       // reset
2600       for(Int_t k = 0; k < 36; k++){
2601         for(Int_t j = 0; j < 16; j++){
2602           for(Int_t c = 0; c < 144; c++){
2603             phvalue[j][c][k] = 0.0;
2604           }
2605         }
2606       }
2607     }
2608   }
2609
2610   return withInput;
2611
2612 }
2613 //_____________________________________________________________________
2614 Int_t AliTRDCalibraFillHisto::ProcessEventDAQ(AliRawReader *rawReader, Bool_t nocheck)
2615 {
2616   //
2617   //  Event processing loop - AliRawReader
2618   //  Testbeam 2007 version
2619   //
2620
2621   AliTRDrawStreamBase rawStream(rawReader);
2622
2623   rawReader->Select("TRD");
2624
2625   return ProcessEventDAQ(&rawStream, nocheck);
2626 }
2627
2628 //_________________________________________________________________________
2629 Int_t AliTRDCalibraFillHisto::ProcessEventDAQ(
2630 #ifdef ALI_DATE
2631                                                const eventHeaderStruct *event,
2632                                                Bool_t nocheck
2633 #else
2634                                                const eventHeaderStruct* /*event*/,
2635                                                Bool_t /*nocheck*/
2636
2637 #endif
2638                                    )
2639 {
2640   //
2641   //  process date event
2642   //  Testbeam 2007 version
2643   //
2644 #ifdef ALI_DATE
2645     AliRawReader *rawReader = new AliRawReaderDate((void*)event);
2646     Int_t result=ProcessEventDAQ(rawReader, nocheck);
2647     delete rawReader;
2648     return result;
2649 #else
2650     Fatal("AliTRDCalibraFillHisto", "this class was compiled without DATE");
2651     return 0;
2652 #endif
2653
2654 }
2655 //_____________________________________________________________________
2656 Int_t AliTRDCalibraFillHisto::ProcessEventDAQ2(AliRawReader *rawReader)
2657  { //main
2658   //
2659   // Event Processing loop - AliTRDrawFastStream
2660   //
2661   // 0 timebin problem
2662   // 1 no input
2663   // 2 input
2664   // Same algorithm as TestBeam but different reader
2665   //
2666
2667    //  AliTRDrawFastStream *rawStream = AliTRDrawFastStream::GetRawStream(rawReader);
2668   AliTRDrawFastStream *rawStream = new AliTRDrawFastStream(rawReader);
2669   rawStream->SetNoErrorWarning();
2670   rawStream->SetSharedPadReadout(kFALSE);
2671
2672   AliTRDdigitsManager *digitsManager = new AliTRDdigitsManager(kTRUE);
2673   digitsManager->CreateArrays();
2674
2675   Int_t withInput = 1;
2676
2677   Double_t phvalue[16][144][36];
2678   for(Int_t k = 0; k < 36; k++){
2679     for(Int_t j = 0; j < 16; j++){
2680       for(Int_t c = 0; c < 144; c++){
2681         phvalue[j][c][k] = 0.0;
2682       }
2683     }
2684   }
2685
2686   fDetectorPreviousTrack = -1;
2687   fMCMPrevious           = -1;
2688   fROBPrevious           = -1;
2689
2690   Int_t nbtimebin = 0;
2691   Int_t baseline  = 10;
2692
2693
2694     fTimeMax = 0;
2695
2696     Int_t det    = 0;
2697     while ((det = rawStream->NextChamber(digitsManager, NULL, NULL)) >= 0) { //idetector
2698
2699       if (digitsManager->GetIndexes(det)->HasEntry()) {//QA
2700         //      printf("there is ADC data on this chamber!\n");
2701
2702         AliTRDarrayADC *digits = (AliTRDarrayADC *) digitsManager->GetDigits(det); //mod
2703         if (digits->HasData()) { //array
2704
2705           AliTRDSignalIndex   *indexes = digitsManager->GetIndexes(det);
2706           if (indexes->IsAllocated() == kFALSE) {
2707             AliError("Indexes do not exist!");
2708             break;
2709           }
2710           Int_t iRow  = 0;
2711           Int_t iCol  = 0;
2712           indexes->ResetCounters();
2713
2714           while (indexes->NextRCIndex(iRow, iCol)) { //column,row
2715             //printf(" det %d \t row %d \t col %d \t digit\n",det,iRow,iCol);
2716             //while (rawStream->Next()) {
2717
2718             Int_t idetector = det;                                             //  current detector
2719             //Int_t imcm      = rawStream->GetMCM();                            //  current MCM
2720             //Int_t irob      = rawStream->GetROB();                            //  current ROB
2721
2722
2723             if((fDetectorPreviousTrack != idetector) && (fDetectorPreviousTrack != -1)) {
2724               // Fill
2725               withInput = TMath::Max(FillDAQ(phvalue),withInput);
2726
2727               // reset
2728               for(Int_t k = 0; k < 36; k++){
2729                 for(Int_t j = 0; j < 16; j++){
2730                   for(Int_t c = 0; c < 144; c++){
2731                     phvalue[j][c][k] = 0.0;
2732                   }
2733                 }
2734               }
2735             }
2736
2737             fDetectorPreviousTrack = idetector;
2738             //fMCMPrevious           = imcm;
2739             //fROBPrevious           = irob;
2740
2741             //    nbtimebin              = rawStream->GetNumberOfTimeBins();              //  number of time bins read from data
2742             AliTRDdigitsParam *digitParam = (AliTRDdigitsParam *)digitsManager->GetDigitsParam();
2743             nbtimebin              = digitParam->GetNTimeBins(det);              //  number of time bins read from data
2744             baseline               = digitParam->GetADCbaseline(det);            //  baseline
2745
2746             if(nbtimebin == 0) return 0;
2747             if((fTimeMax != 0) && (nbtimebin != fTimeMax)) return 0;
2748             fTimeMax          = nbtimebin;
2749
2750             fNumberClustersf    = fTimeMax;
2751             fNumberClusters     = (Int_t)(fNumberClustersProcent*fTimeMax);
2752
2753
2754             for(Int_t itime = 0; itime < nbtimebin; itime++) {
2755               //            phvalue[row][col][itime] = signal[itime]-baseline;
2756               phvalue[iRow][iCol][itime] = (Short_t)(digits->GetData(iRow,iCol,itime) - baseline);
2757               /*if(phvalue[iRow][iCol][itime] >= 20) {
2758                  printf("----------> phvalue[%d][%d][%d] %d  baseline %d \n",
2759                        iRow,
2760                        iCol,
2761                        itime,
2762                        (Short_t)(digits->GetData(iRow,iCol,itime)),
2763                        baseline);
2764                        }*/
2765             }
2766
2767           }//column,row
2768
2769           // fill the last one
2770           if(fDetectorPreviousTrack != -1){
2771
2772             // Fill
2773             withInput = TMath::Max(FillDAQ(phvalue),withInput);
2774             //      printf("\n ---> withinput %d\n\n",withInput);
2775             // reset
2776             for(Int_t k = 0; k < 36; k++){
2777               for(Int_t j = 0; j < 16; j++){
2778                 for(Int_t c = 0; c < 144; c++){
2779                   phvalue[j][c][k] = 0.0;
2780                 }
2781               }
2782             }
2783           }
2784
2785         }//array
2786       }//QA
2787       digitsManager->ClearArrays(det);
2788     }//idetector
2789     delete digitsManager;
2790
2791     delete rawStream;
2792     return withInput;
2793  }//main
2794 //_____________________________________________________________________
2795 Int_t AliTRDCalibraFillHisto::ProcessEventDAQ3(AliRawReader *rawReader)
2796  { //main
2797   //
2798   // Event Processing loop - AliTRDrawStream
2799   //
2800   // 0 timebin problem
2801   // 1 no input
2802   // 2 input
2803   // Same algorithm as TestBeam but different reader
2804   //
2805
2806    //  AliTRDrawFastStream *rawStream = AliTRDrawFastStream::GetRawStream(rawReader);
2807   AliTRDrawStream *rawStream = new AliTRDrawStream(rawReader);
2808   rawStream->SetNoErrorWarning();
2809   rawStream->SetSharedPadReadout(kFALSE);
2810
2811   AliTRDdigitsManager *digitsManager = new AliTRDdigitsManager(kTRUE);
2812   digitsManager->CreateArrays();
2813
2814   Int_t withInput = 1;
2815
2816   Double_t phvalue[16][144][36];
2817   for(Int_t k = 0; k < 36; k++){
2818     for(Int_t j = 0; j < 16; j++){
2819       for(Int_t c = 0; c < 144; c++){
2820         phvalue[j][c][k] = 0.0;
2821       }
2822     }
2823   }
2824
2825   fDetectorPreviousTrack = -1;
2826   fMCMPrevious           = -1;
2827   fROBPrevious           = -1;
2828
2829   Int_t nbtimebin = 0;
2830   Int_t baseline  = 10;
2831
2832
2833     fTimeMax = 0;
2834
2835     Int_t det    = 0;
2836     while ((det = rawStream->NextChamber(digitsManager, NULL, NULL)) >= 0) { //idetector
2837
2838       if (digitsManager->GetIndexes(det)->HasEntry()) {//QA
2839         //      printf("there is ADC data on this chamber!\n");
2840
2841         AliTRDarrayADC *digits = (AliTRDarrayADC *) digitsManager->GetDigits(det); //mod
2842         if (digits->HasData()) { //array
2843
2844           AliTRDSignalIndex   *indexes = digitsManager->GetIndexes(det);
2845           if (indexes->IsAllocated() == kFALSE) {
2846             AliError("Indexes do not exist!");
2847             break;
2848           }
2849           Int_t iRow  = 0;
2850           Int_t iCol  = 0;
2851           indexes->ResetCounters();
2852
2853           while (indexes->NextRCIndex(iRow, iCol)) { //column,row
2854             //printf(" det %d \t row %d \t col %d \t digit\n",det,iRow,iCol);
2855             //while (rawStream->Next()) {
2856
2857             Int_t idetector = det;                                             //  current detector
2858             //Int_t imcm      = rawStream->GetMCM();                            //  current MCM
2859             //Int_t irob      = rawStream->GetROB();                            //  current ROB
2860
2861
2862             if((fDetectorPreviousTrack != idetector) && (fDetectorPreviousTrack != -1)) {
2863               // Fill
2864               withInput = TMath::Max(FillDAQ(phvalue),withInput);
2865
2866               // reset
2867               for(Int_t k = 0; k < 36; k++){
2868                 for(Int_t j = 0; j < 16; j++){
2869                   for(Int_t c = 0; c < 144; c++){
2870                     phvalue[j][c][k] = 0.0;
2871                   }
2872                 }
2873               }
2874             }
2875
2876             fDetectorPreviousTrack = idetector;
2877             //fMCMPrevious           = imcm;
2878             //fROBPrevious           = irob;
2879
2880             //    nbtimebin              = rawStream->GetNumberOfTimeBins();              //  number of time bins read from data
2881             AliTRDdigitsParam *digitParam = (AliTRDdigitsParam *)digitsManager->GetDigitsParam();
2882             nbtimebin              = digitParam->GetNTimeBins(det);              //  number of time bins read from data
2883             baseline               = digitParam->GetADCbaseline(det);            //  baseline
2884
2885             if(nbtimebin == 0) return 0;
2886             if((fTimeMax != 0) && (nbtimebin != fTimeMax)) return 0;
2887             fTimeMax          = nbtimebin;
2888
2889             fNumberClustersf    = fTimeMax;
2890             fNumberClusters     = (Int_t)(fNumberClustersProcent*fTimeMax);
2891
2892
2893             for(Int_t itime = 0; itime < nbtimebin; itime++) {
2894               //            phvalue[row][col][itime] = signal[itime]-baseline;
2895               phvalue[iRow][iCol][itime] = (Short_t)(digits->GetData(iRow,iCol,itime) - baseline);
2896               /*if(phvalue[iRow][iCol][itime] >= 20) {
2897                  printf("----------> phvalue[%d][%d][%d] %d  baseline %d \n",
2898                        iRow,
2899                        iCol,
2900                        itime,
2901                        (Short_t)(digits->GetData(iRow,iCol,itime)),
2902                        baseline);
2903                        }*/
2904             }
2905
2906           }//column,row
2907
2908           // fill the last one
2909           if(fDetectorPreviousTrack != -1){
2910
2911             // Fill
2912             withInput = TMath::Max(FillDAQ(phvalue),withInput);
2913             //      printf("\n ---> withinput %d\n\n",withInput);
2914             // reset
2915             for(Int_t k = 0; k < 36; k++){
2916               for(Int_t j = 0; j < 16; j++){
2917                 for(Int_t c = 0; c < 144; c++){
2918                   phvalue[j][c][k] = 0.0;
2919                 }
2920               }
2921             }
2922           }
2923
2924         }//array
2925       }//QA
2926       digitsManager->ClearArrays(det);
2927     }//idetector
2928     delete digitsManager;
2929
2930     delete rawStream;
2931     return withInput;
2932  }//main
2933 //_____________________________________________________________________
2934 //////////////////////////////////////////////////////////////////////////////
2935 // Routine inside the DAQ process
2936 /////////////////////////////////////////////////////////////////////////////
2937 //_______________________________________________________________________
2938 Int_t AliTRDCalibraFillHisto::FillDAQ(Double_t phvalue[16][144][36]){
2939
2940   //
2941   // Look for the maximum by collapsing over the time
2942   // Sum over four pad col and two pad row
2943   //
2944
2945   Int_t used = 0;
2946
2947
2948   Int_t idect = fDetectorPreviousTrack;
2949   //printf("Enter Detector %d\n",fDetectorPreviousTrack);
2950   Double_t sum[36];
2951   for(Int_t tb = 0; tb < 36; tb++){
2952     sum[tb] = 0.0;
2953   }
2954
2955   //fGoodTracklet = kTRUE;
2956   //fDetectorPreviousTrack = 0;
2957
2958
2959   ///////////////////////////
2960   // look for maximum
2961   /////////////////////////
2962
2963   Int_t imaxRow = 0;
2964   Int_t imaxCol = 0;
2965   Double_t integralMax = -1;
2966
2967   for (Int_t ir = 1; ir <= 15; ir++)
2968     {
2969       for (Int_t ic = 2; ic <= 142; ic++)
2970         {
2971           Double_t integral = 0;
2972           for (Int_t ishiftR = 0; ishiftR < fNumberRowDAQ; ishiftR++)
2973             {
2974               for (Int_t ishiftC = -fNumberColDAQ; ishiftC < fNumberColDAQ; ishiftC++)
2975                 {
2976                   if (ir + ishiftR >= 1 && ir + ishiftR <= 16 &&
2977                       ic + ishiftC >= 1 && ic + ishiftC <= 144)
2978                     {
2979
2980                       for(Int_t tb = 0; tb< fTimeMax; tb++){
2981                         integral += phvalue[ir + ishiftR-1][ic + ishiftC-1][tb];
2982                       }// addtb
2983                     } //addsignal
2984                 } //shiftC
2985             } // shiftR
2986           if (integralMax < integral)
2987             {
2988               imaxRow = ir;
2989               imaxCol = ic;
2990               integralMax = integral;
2991
2992             } // check max integral
2993         } //ic
2994     } // ir
2995
2996   //  printf("imaxRow %d, imaxCol %d, fTimeMax %d, integralMax %f\n",imaxRow,imaxCol,fTimeMax, integralMax);
2997   //if((imaxRow == 0) || (imaxRow >= 15) || (imaxCol <= 3) || (imaxCol >= 140)) {
2998   //  used=1;
2999   //  return used;
3000   // }
3001
3002   if(((imaxRow + fNumberRowDAQ) > 16) || (imaxRow == 0) || ((imaxCol - fNumberColDAQ) <= 1) || ((imaxCol + fNumberColDAQ) >= 144)) {
3003     used=1;
3004     return used;
3005   }
3006   //CheckGoodTrackletV0(fDetectorPreviousTrack,imaxRow,imaxCol);
3007   //if(!fGoodTracklet) used = 1;;
3008
3009   //  /////////////////////////////////////////////////////
3010   // sum ober 2 row and 4 pad cols for each time bins
3011   //  ////////////////////////////////////////////////////
3012
3013
3014
3015   for (Int_t ishiftR = 0; ishiftR < fNumberRowDAQ; ishiftR++)
3016     {
3017       for (Int_t ishiftC = -fNumberColDAQ; ishiftC < fNumberColDAQ; ishiftC++)
3018         {
3019           if (imaxRow + ishiftR >= 1 && imaxRow + ishiftR <= 16 &&
3020               imaxCol + ishiftC >= 1 && imaxCol + ishiftC <= 144)
3021             {
3022               for(Int_t it = 0; it < fTimeMax; it++){
3023                 sum[it] += phvalue[imaxRow + ishiftR-1][imaxCol + ishiftC-1][it];
3024               }
3025             }
3026         } // col shift
3027     }// row shift
3028
3029   Int_t nbcl = 0;
3030   Double_t sumcharge = 0.0;
3031   for(Int_t it = 0; it < fTimeMax; it++){
3032     sumcharge += sum[it];
3033     if(sum[it] > fThresholdClustersDAQ)  nbcl++;
3034   }
3035
3036
3037   /////////////////////////////////////////////////////////
3038   // Debug
3039   ////////////////////////////////////////////////////////
3040   if(fDebugLevel > 0){
3041     if ( !fDebugStreamer ) {
3042       //debug stream
3043       TDirectory *backup = gDirectory;
3044       fDebugStreamer = new TTreeSRedirector("TRDdebugCalibraFill.root");
3045       if ( backup ) backup->cd();  //we don't want to be cd'd to the debug streamer
3046     }
3047
3048     Double_t amph0 = sum[0];
3049     Double_t amphlast = sum[fTimeMax-1];
3050     Double_t rms      = TMath::RMS(fTimeMax,sum);
3051     Int_t    goodtracklet = (Int_t) fGoodTracklet;
3052     for(Int_t it = 0; it < fTimeMax; it++){
3053       Double_t clustera = sum[it];
3054
3055     (* fDebugStreamer) << "FillDAQa"<<
3056       "ampTotal="<<sumcharge<<
3057       "row="<<imaxRow<<
3058       "col="<<imaxCol<<
3059       "detector="<<idect<<
3060       "amph0="<<amph0<<
3061       "amphlast="<<amphlast<<
3062       "goodtracklet="<<goodtracklet<<
3063       "clustera="<<clustera<<
3064       "it="<<it<<
3065       "rms="<<rms<<
3066       "nbcl="<<nbcl<<
3067       "\n";
3068     }
3069   }
3070
3071   ////////////////////////////////////////////////////////
3072   // fill
3073   ///////////////////////////////////////////////////////
3074   //printf("fNumberClusters %d, fNumberClustersf %d\n",fNumberClusters,fNumberClustersf);
3075   if(sum[0] > 100.0) used = 1;
3076   if(nbcl < fNumberClusters) used = 1;
3077   if(nbcl > fNumberClustersf) used = 1;
3078
3079   //if(fDetectorPreviousTrack == 15){
3080   //  printf("rms %f and first time bin %f\n",TMath::RMS(fTimeMax,sum),sum[0]);
3081   //}
3082   //if((TMath::RMS(fTimeMax,sum) <= 10.0) && (sum[0] > 200.0)) return 1;
3083   if(used == 0){
3084     for(Int_t it = 0; it < fTimeMax; it++){
3085       if(fFillWithZero) UpdateDAQ(fDetectorPreviousTrack,0,0,it,sum[it],fTimeMax);
3086       else{
3087         if(sum[it] > 0.0) UpdateDAQ(fDetectorPreviousTrack,0,0,it,sum[it],fTimeMax);
3088       }
3089       //if(fFillWithZero) UpdateDAQ(0,0,0,it,sum[it],fTimeMax);
3090       //else{
3091       // if(sum[it] > 0.0) UpdateDAQ(0,0,0,it,sum[it],fTimeMax);
3092       //}
3093     }
3094
3095
3096     //((TH2I *)GetCH2d()->Fill(sumcharge/30.0,fDetectorPreviousTrack));
3097     used = 2;
3098     //printf("Pass Detector %d\n",fDetectorPreviousTrack);
3099
3100   }
3101
3102   return used;
3103
3104 }
3105 //____________Online trackling in AliTRDtrigger________________________________
3106 Bool_t AliTRDCalibraFillHisto::UpdateDAQ(Int_t det, Int_t /*row*/, Int_t /*col*/, Int_t timebin, Float_t signal, Int_t nbtimebins)
3107 {
3108   //
3109   // For the DAQ
3110   // Fill a simple average pulse height
3111   //
3112
3113
3114   ((TProfile2D *)GetPH2d(nbtimebins,fSf))->Fill((Float_t) timebin/fSf,det+0.5,(Float_t) signal);
3115
3116
3117   return kTRUE;
3118
3119 }
3120 //____________Write_____________________________________________________
3121 //_____________________________________________________________________
3122 void AliTRDCalibraFillHisto::Write2d(const Char_t *filename, Bool_t append)
3123 {
3124   //
3125   //  Write infos to file
3126   //
3127
3128   //For debugging
3129   if ( fDebugStreamer ) {
3130     delete fDebugStreamer;
3131     fDebugStreamer = 0x0;
3132   }
3133
3134   AliInfo(Form("Numbertrack: %d Numberusedch[0]: %d, Numberusedch[1]: %d Numberusedph[0]: %d, Numberusedph[1]: %d"
3135                ,fNumberTrack
3136                ,fNumberUsedCh[0]
3137                ,fNumberUsedCh[1]
3138                ,fNumberUsedPh[0]
3139                ,fNumberUsedPh[1]));
3140
3141   TDirectory *backup = gDirectory;
3142   TString option;
3143
3144   if ( append )
3145     option = "update";
3146   else
3147     option = "recreate";
3148
3149   TFile f(filename,option.Data());
3150
3151   TStopwatch stopwatch;
3152   stopwatch.Start();
3153   if(fVector2d) {
3154     f.WriteTObject(fCalibraVector);
3155   }
3156
3157   if (fCH2dOn ) {
3158     if (fHisto2d) {
3159       f.WriteTObject(fCH2d);
3160     }
3161   }
3162   if (fPH2dOn ) {
3163     if (fHisto2d) {
3164       f.WriteTObject(fPH2d);
3165     }
3166   }
3167   if (fPRF2dOn) {
3168     if (fHisto2d) {
3169         f.WriteTObject(fPRF2d);
3170     }
3171   }
3172   if(fLinearFitterOn){
3173     if(fLinearFitterDebugOn) AnalyseLinearFitter();
3174     f.WriteTObject(fLinearVdriftFit);
3175   }
3176
3177   f.Close();
3178
3179   if ( backup ) backup->cd();
3180
3181   AliInfo(Form("Execution time Write2d: R:%.2fs C:%.2fs"
3182                ,stopwatch.RealTime(),stopwatch.CpuTime()));
3183 }
3184 //////////////////////////////////////////////////////////////////////////////////////////////////////////////
3185 // Stats stuff
3186 //////////////////////////////////////////////////////////////////////////////////////////////////////////////
3187 //___________________________________________probe the histos__________________________________________________
3188 Double_t *AliTRDCalibraFillHisto::StatH(TH2 *h, Int_t i)
3189 {
3190   //
3191   // Check the number of stats in h, 0 is TH2I 1 is TProfile2D
3192   // debug mode with 2 for TH2I and 3 for TProfile2D
3193   // It gives a pointer to a Double_t[7] with the info following...
3194   // [0] : number of calibration groups with entries
3195   // [1] : minimal number of entries found
3196   // [2] : calibration group number of the min
3197   // [3] : maximal number of entries found
3198   // [4] : calibration group number of the max
3199   // [5] : mean number of entries found
3200   // [6] : mean relative error
3201   //
3202
3203   Double_t *info = new Double_t[7];
3204
3205   // Number of Xbins (detectors or groups of pads)
3206   Int_t    nbins   = h->GetNbinsY(); //number of calibration groups
3207   Int_t    nxbins  = h->GetNbinsX(); //number of bins per histo
3208
3209   // Initialise
3210   Double_t nbwe = 0; //number of calibration groups with entries
3211   Double_t minentries = 0; //minimal number of entries found
3212   Double_t maxentries = 0; //maximal number of entries found
3213   Double_t placemin = 0; //calibration group number of the min
3214   Double_t placemax = -1; //calibration group number of the max
3215   Double_t meanstats = 0.0; //mean number of entries over the calibration group with at least ome entry
3216   Double_t meanrelativerror = 0.0; //mean relativ error in the TProfile2D
3217
3218   Double_t counter = 0;
3219
3220   //Debug
3221   TH1F *nbEntries = 0x0;//distribution of the number of entries
3222   TH1F *nbEntriesPerGroup = 0x0;//Number of entries per group
3223   TProfile *nbEntriesPerSp = 0x0;//Number of entries for one supermodule
3224
3225   // Beginning of the loop over the calibration groups
3226   for (Int_t idect = 0; idect < nbins; idect++) {
3227
3228     TH1I *projch = (TH1I *) h->ProjectionX("projch",idect+1,idect+1,(Option_t *)"e");
3229     projch->SetDirectory(0);
3230
3231     // Number of entries for this calibration group
3232     Double_t nentries = 0.0;
3233     if((i%2) == 0){
3234       for (Int_t k = 0; k < nxbins; k++) {
3235         nentries += h->GetBinContent(h->GetBin(k+1,idect+1));
3236       }
3237     }
3238     else{
3239       for (Int_t k = 0; k < nxbins; k++) {
3240         nentries += ((TProfile2D *)h)->GetBinEntries(h->GetBin(k+1,idect+1));
3241         if(h->GetBinContent(h->GetBin(k+1,idect+1)) != 0) {
3242           meanrelativerror += (h->GetBinError(h->GetBin(k+1,idect+1))/(TMath::Abs(h->GetBinContent(h->GetBin(k+1,idect+1)))));
3243           counter++;
3244         }
3245       }
3246     }
3247
3248     //Debug
3249     if(i > 1){
3250       if((!((Bool_t)nbEntries)) && (nentries > 0)){
3251         nbEntries = new TH1F("Number of entries","Number of entries"
3252                                ,100,(Int_t)nentries/2,nentries*2);
3253         nbEntries->SetDirectory(0);
3254         nbEntriesPerGroup = new TH1F("Number of entries per group","Number of entries per group"
3255                                ,nbins,0,nbins);
3256         nbEntriesPerGroup->SetDirectory(0);
3257         nbEntriesPerSp = new TProfile("Number of entries per supermodule","Number of entries per supermodule"
3258                                ,(Int_t)(nbins/18),0,(Int_t)(nbins/18));
3259         nbEntriesPerSp->SetDirectory(0);
3260       }
3261       if(nbEntries){
3262         if(nentries > 0) nbEntries->Fill(nentries);
3263         nbEntriesPerGroup->Fill(idect+0.5,nentries);
3264         nbEntriesPerSp->Fill((idect%((Int_t)(nbins/18)))+0.5,nentries);
3265       }
3266     }
3267
3268     //min amd max
3269     if(nentries > maxentries){
3270       maxentries = nentries;
3271       placemax = idect;
3272     }
3273     if(idect == 0) {
3274       minentries = nentries;
3275     }
3276     if(nentries < minentries){
3277       minentries = nentries;
3278       placemin = idect;
3279     }
3280     //nbwe
3281     if(nentries > 0) {
3282       nbwe++;
3283       meanstats += nentries;
3284     }
3285   }//calibration groups loop
3286
3287   if(nbwe > 0) meanstats /= nbwe;
3288   if(counter > 0) meanrelativerror /= counter;
3289
3290   AliInfo(Form("There are %f calibration groups with entries",nbwe));
3291   AliInfo(Form("The minimum number of entries is %f for the group %f",minentries,placemin));
3292   AliInfo(Form("The maximum number of entries is %f for the group %f",maxentries,placemax));
3293   AliInfo(Form("The mean number of entries is %f",meanstats));
3294   if((i%2) == 1) AliInfo(Form("The mean relative error is %f",meanrelativerror));
3295
3296   info[0] = nbwe;
3297   info[1] = minentries;
3298   info[2] = placemin;
3299   info[3] = maxentries;
3300   info[4] = placemax;
3301   info[5] = meanstats;
3302   info[6] = meanrelativerror;
3303
3304   if(i > 1){
3305     gStyle->SetPalette(1);
3306     gStyle->SetOptStat(1111);
3307     gStyle->SetPadBorderMode(0);
3308     gStyle->SetCanvasColor(10);
3309     gStyle->SetPadLeftMargin(0.13);
3310     gStyle->SetPadRightMargin(0.01);
3311     TCanvas *stat = new TCanvas("stat","",50,50,600,800);
3312     stat->Divide(2,1);
3313     stat->cd(1);
3314     nbEntries->Draw("");
3315     stat->cd(2);
3316     nbEntriesPerSp->SetStats(0);
3317     nbEntriesPerSp->Draw("");
3318     TCanvas *stat1 = new TCanvas("stat1","",50,50,600,800);
3319     stat1->cd();
3320     nbEntriesPerGroup->SetStats(0);
3321     nbEntriesPerGroup->Draw("");
3322   }
3323
3324   return info;
3325
3326 }
3327 //____________________________________________________________________________
3328 Double_t *AliTRDCalibraFillHisto::GetMeanMedianRMSNumberCH()
3329 {
3330   //
3331   // Return a Int_t[4] with:
3332   // 0 Mean number of entries
3333   // 1 median of number of entries
3334   // 2 rms of number of entries
3335   // 3 number of group with entries
3336   //
3337
3338   Double_t *stat      = new Double_t[4];
3339   stat[3]             = 0.0;
3340
3341   Int_t    nbofgroups = CalculateTotalNumberOfBins(0);
3342   Double_t *weight    = new Double_t[nbofgroups];
3343   Int_t    *nonul     = new Int_t[nbofgroups];
3344
3345   for(Int_t k = 0; k < nbofgroups; k++){
3346     if(fEntriesCH[k] > 0) {
3347       weight[k] = 1.0;
3348       nonul[(Int_t)stat[3]] = fEntriesCH[k];
3349       stat[3]++;
3350     }
3351     else weight[k] = 0.0;
3352   }
3353   stat[0]          = TMath::Mean(nbofgroups,fEntriesCH,weight);
3354   stat[1]          = TMath::Median(nbofgroups,fEntriesCH,weight);
3355   stat[2]          = TMath::RMS((Int_t)stat[3],nonul);
3356
3357   return stat;
3358
3359 }
3360 //____________________________________________________________________________
3361 Double_t *AliTRDCalibraFillHisto::GetMeanMedianRMSNumberLinearFitter() const
3362 {
3363   //
3364   // Return a Int_t[4] with:
3365   // 0 Mean number of entries
3366   // 1 median of number of entries
3367   // 2 rms of number of entries
3368   // 3 number of group with entries
3369   //
3370
3371   Double_t *stat      = new Double_t[4];
3372   stat[3]             = 0.0;
3373
3374   Int_t    nbofgroups = 540;
3375   Double_t *weight    = new Double_t[nbofgroups];
3376   Int_t    *nonul     = new Int_t[nbofgroups];
3377
3378   for(Int_t k = 0; k < nbofgroups; k++){
3379     if(fEntriesLinearFitter[k] > 0) {
3380       weight[k] = 1.0;
3381       nonul[(Int_t) stat[3]] = fEntriesLinearFitter[k];
3382       stat[3]++;
3383     }
3384     else weight[k] = 0.0;
3385   }
3386   stat[0]          = TMath::Mean(nbofgroups,fEntriesLinearFitter,weight);
3387   stat[1]          = TMath::Median(nbofgroups,fEntriesLinearFitter,weight);
3388   stat[2]          = TMath::RMS((Int_t)stat[3],nonul);
3389
3390   return stat;
3391
3392 }
3393 //////////////////////////////////////////////////////////////////////////////////////
3394 // Create Histos
3395 //////////////////////////////////////////////////////////////////////////////////////
3396 //_____________________________________________________________________________
3397 void AliTRDCalibraFillHisto::CreatePRF2d(Int_t nn)
3398 {
3399   //
3400   // Create the 2D histos: here we have 2*fNgroupprf bins in tnp of 0.2 amplitude each
3401   // If fNgroupprf is zero then no binning in tnp
3402   //
3403
3404   TString name("Nz");
3405   name += fCalibraMode->GetNz(2);
3406   name += "Nrphi";
3407   name += fCalibraMode->GetNrphi(2);
3408   name += "Ngp";
3409   name += fNgroupprf;
3410
3411   if(fNgroupprf != 0){
3412
3413     fPRF2d = new TProfile2D("PRF2d",(const Char_t *) name
3414                             ,2*fNgroupprf*fNumberBinPRF,-3.0*fNgroupprf,3.0*fNgroupprf,nn,0,nn);
3415     fPRF2d->SetYTitle("Det/pad groups");
3416     fPRF2d->SetXTitle("Position x/W [pad width units]");
3417     fPRF2d->SetZTitle("Q_{i}/Q_{total}");
3418     fPRF2d->SetStats(0);
3419   }
3420   else{
3421     fPRF2d = new TProfile2D("PRF2d",(const Char_t *) name
3422                             ,fNumberBinPRF,-1.5,1.5,nn,0,nn);
3423     fPRF2d->SetYTitle("Det/pad groups");
3424     fPRF2d->SetXTitle("Position x/W [pad width units]");
3425     fPRF2d->SetZTitle("Q_{i}/Q_{total}");
3426     fPRF2d->SetStats(0);
3427   }
3428
3429 }
3430
3431 //_____________________________________________________________________________
3432 void AliTRDCalibraFillHisto::CreatePH2d(Int_t nn)
3433 {
3434   //
3435   // Create the 2D histos
3436   //
3437
3438   TString name("Ver");
3439   name += fVersionVdriftUsed;
3440   name += "Subver";
3441   name += fSubVersionVdriftUsed;
3442   name += "Nz";
3443   name += fCalibraMode->GetNz(1);
3444   name += "Nrphi";
3445   name += fCalibraMode->GetNrphi(1);
3446
3447   fPH2d = new TProfile2D("PH2d",(const Char_t *) name
3448                          ,fTimeMax,-0.5/fSf,(Float_t) (fTimeMax-0.5)/fSf
3449                          ,nn,0,nn);
3450   fPH2d->SetYTitle("Det/pad groups");
3451   fPH2d->SetXTitle("time [#mus]");
3452   fPH2d->SetZTitle("<PH> [a.u.]");
3453   fPH2d->SetStats(0);
3454
3455 }
3456 //_____________________________________________________________________________
3457 void AliTRDCalibraFillHisto::CreateCH2d(Int_t nn)
3458 {
3459   //
3460   // Create the 2D histos
3461   //
3462
3463   TString name("Ver");
3464   name += fVersionGainUsed;
3465   name += "Subver";
3466   name += fSubVersionGainUsed;
3467   name += "Nz";
3468   name += fCalibraMode->GetNz(0);
3469   name += "Nrphi";
3470   name += fCalibraMode->GetNrphi(0);
3471
3472   fCH2d = new TH2I("CH2d",(const Char_t *) name
3473                    ,fNumberBinCharge,0,300,nn,0,nn);
3474   fCH2d->SetYTitle("Det/pad groups");
3475   fCH2d->SetXTitle("charge deposit [a.u]");
3476   fCH2d->SetZTitle("counts");
3477   fCH2d->SetStats(0);
3478   fCH2d->Sumw2();
3479
3480 }
3481 //////////////////////////////////////////////////////////////////////////////////
3482 // Set relative scale
3483 /////////////////////////////////////////////////////////////////////////////////
3484 //_____________________________________________________________________________
3485 void AliTRDCalibraFillHisto::SetRelativeScale(Float_t RelativeScale)
3486 {
3487   //
3488   // Set the factor that will divide the deposited charge
3489   // to fit in the histo range [0,300]
3490   //
3491
3492   if (RelativeScale > 0.0) {
3493     fRelativeScale = RelativeScale;
3494   }
3495   else {
3496     AliInfo("RelativeScale must be strict positif!");
3497   }
3498
3499 }
3500 //////////////////////////////////////////////////////////////////////////////////
3501 // Quick way to fill a histo
3502 //////////////////////////////////////////////////////////////////////////////////
3503 //_____________________________________________________________________
3504 void  AliTRDCalibraFillHisto::FillCH2d(Int_t x, Float_t y)
3505 {
3506   //
3507   // FillCH2d: Marian style
3508   //
3509
3510   //skip simply the value out of range
3511   if((y>=300.0) || (y<0.0)) return;
3512
3513   //Calcul the y place
3514   Int_t yplace = (Int_t) (fNumberBinCharge*y/300.0)+1;
3515   Int_t place = (fNumberBinCharge+2)*(x+1)+yplace;
3516
3517   //Fill
3518   fCH2d->GetArray()[place]++;
3519
3520 }
3521
3522 //////////////////////////////////////////////////////////////////////////////////
3523 // Geometrical functions
3524 ///////////////////////////////////////////////////////////////////////////////////
3525 //_____________________________________________________________________________
3526 Int_t AliTRDCalibraFillHisto::GetLayer(Int_t d) const
3527 {
3528   //
3529   // Reconstruct the layer number from the detector number
3530   //
3531
3532   return ((Int_t) (d % 6));
3533
3534 }
3535
3536 //_____________________________________________________________________________
3537 Int_t AliTRDCalibraFillHisto::GetStack(Int_t d) const
3538 {
3539   //
3540   // Reconstruct the stack number from the detector number
3541   //
3542   const Int_t kNlayer = 6;
3543
3544   return ((Int_t) (d % 30) / kNlayer);
3545
3546 }
3547
3548 //_____________________________________________________________________________
3549 Int_t AliTRDCalibraFillHisto::GetSector(Int_t d) const
3550 {
3551   //
3552   // Reconstruct the sector number from the detector number
3553   //
3554   Int_t fg = 30;
3555
3556   return ((Int_t) (d / fg));
3557
3558 }
3559 ///////////////////////////////////////////////////////////////////////////////////
3560 // Getter functions for DAQ of the CH2d and the PH2d
3561 //////////////////////////////////////////////////////////////////////////////////
3562 //_____________________________________________________________________
3563 TProfile2D* AliTRDCalibraFillHisto::GetPH2d(Int_t nbtimebin, Float_t samplefrequency)
3564 {
3565     //
3566     // return pointer to fPH2d TProfile2D
3567     // create a new TProfile2D if it doesn't exist allready
3568     //
3569     if ( fPH2d )
3570         return fPH2d;
3571
3572     // Some parameters
3573     fTimeMax = nbtimebin;
3574     fSf      = samplefrequency;
3575
3576     CreatePH2d(540);
3577
3578     return fPH2d;
3579 }
3580 //_____________________________________________________________________
3581 TH2I* AliTRDCalibraFillHisto::GetCH2d()
3582 {
3583     //
3584     // return pointer to fCH2d TH2I
3585     // create a new TH2I if it doesn't exist allready
3586     //
3587     if ( fCH2d )
3588         return fCH2d;
3589
3590     CreateCH2d(540);
3591
3592     return fCH2d;
3593 }
3594 ////////////////////////////////////////////////////////////////////////////////////////////
3595 // Drift velocity calibration
3596 ///////////////////////////////////////////////////////////////////////////////////////////
3597 //_____________________________________________________________________
3598 TLinearFitter* AliTRDCalibraFillHisto::GetLinearFitter(Int_t detector, Bool_t force)
3599 {
3600     //
3601     // return pointer to TLinearFitter Calibration
3602     // if force is true create a new TLinearFitter if it doesn't exist allready
3603     //
3604
3605   if ((!force) || (fLinearFitterArray.UncheckedAt(detector))){
3606     return (TLinearFitter*)fLinearFitterArray.UncheckedAt(detector);
3607   }
3608
3609   // if we are forced and TLinearFitter doesn't yet exist create it
3610
3611   // new TLinearFitter
3612   TLinearFitter *linearfitter = new TLinearFitter(2,"pol1");
3613   fLinearFitterArray.AddAt(linearfitter,detector);
3614   return linearfitter;
3615 }
3616
3617 //____________________________________________________________________________
3618 void AliTRDCalibraFillHisto::AnalyseLinearFitter()
3619 {
3620   //
3621   // Analyse array of linear fitter because can not be written
3622   // Store two arrays: one with the param the other one with the error param + number of entries
3623   //
3624
3625   for(Int_t k = 0; k < 540; k++){
3626     TLinearFitter *linearfitter = GetLinearFitter(k);
3627     if((linearfitter!=0) && (fEntriesLinearFitter[k]>10)){
3628       TVectorD  *par  = new TVectorD(2);
3629       TVectorD   pare = TVectorD(2);
3630       TVectorD  *parE = new TVectorD(3);
3631       linearfitter->Eval();
3632       linearfitter->GetParameters(*par);
3633       linearfitter->GetErrors(pare);
3634       Float_t  ppointError =  TMath::Sqrt(TMath::Abs(linearfitter->GetChisquare())/fEntriesLinearFitter[k]);
3635       (*parE)[0] = pare[0]*ppointError;
3636       (*parE)[1] = pare[1]*ppointError;
3637       (*parE)[2] = (Double_t) fEntriesLinearFitter[k];
3638       ((TObjArray *)fLinearVdriftFit->GetPArray())->AddAt(par,k);
3639       ((TObjArray *)fLinearVdriftFit->GetEArray())->AddAt(parE,k);
3640     }
3641   }
3642 }
3643
3644