TRD/qaRec/AliTRDcalibration.cxx

   1
   2 /**************************************************************************
   3 * Copyright(c) 1998-1999, ALICE Experiment at CERN, All rights reserved. *
   4 *                                                                        *
   5 * Author: The ALICE Off-line Project.                                    *
   6 * Contributors are mentioned in the code where appropriate.              *
   7 *                                                                        *
   8 * Permission to use, copy, modify and distribute this software and its   *
   9 * documentation strictly for non-commercial purposes is hereby granted   *
  10 * without fee, provided that the above copyright notice appears in all   *
  11 * copies and that both the copyright notice and this permission notice   *
  12 * appear in the supporting documentation. The authors make no claims     *
  13 * about the suitability of this software for any purpose. It is          *
  14 * provided "as is" without express or implied warranty.                  *
  15 **************************************************************************/
  16
  17 /////////////////////////////////////////////////////////////////////////////////
  18 //
  19 // AliTRDcalibration
  20 //
  21 // Task to run the calibration offline.
  22 // Author:
  23 //   R. Bailhache (rbailhache@ikf.uni-frankfurt.de, R.Bailhache@gsi.de)
  24 //
  25 //////////////////////////////////////////////////////////////////////////////////
  26
  27
  28 #include "Riostream.h"
  29 #include "TChain.h"
  30 #include "TTree.h"
  31 #include "TProfile2D.h"
  32 #include "TH2I.h"
  33 #include "TH1F.h"
  34 #include "TList.h"
  35 #include "TMath.h"
  36 #include "TCanvas.h"
  37 #include "TObject.h"
  38 #include "TFile.h"
  39 #include "TObjArray.h"
  40 #include "TGraph.h"
  41 #include "TStyle.h"
  42 #include "TLegend.h"
  43 #include "TGraphErrors.h"
  44
  45 #include "AliTRDrecoTask.h"
  46 #include "AliAnalysisManager.h"
  47
  48 #include "AliESDInputHandler.h"
  49 #include "AliTRDtrackV1.h"
  50 #include "AliTRDseedV1.h"
  51 #include "AliTRDcluster.h"
  52 #include "info/AliTRDtrackInfo.h"
  53 #include "AliTRDcalibDB.h"
  54
  55 #include "AliTRDCalibraFillHisto.h"
  56 #include "AliTRDCalibraFit.h"
  57 #include "AliTRDCalibraVdriftLinearFit.h"
  58 #include "AliTRDCalibraMode.h"
  59 #include "AliTRDCalibraVector.h"
  60 #include "./Cal/AliTRDCalPad.h"
  61 #include "./Cal/AliTRDCalDet.h"
  62 #include "AliCDBMetaData.h"
  63 #include "AliCDBManager.h"
  64 #include "AliCDBStorage.h"
  65
  66 #include "AliLog.h"
  67
  68 #include "AliTRDcalibration.h"
  69
  70
  71 ClassImp(AliTRDcalibration)
  72
  73 //________________________________________________________________________
  74 AliTRDcalibration::AliTRDcalibration()
  75   :AliTRDrecoTask("Calibration", "Calibration on tracks")
  76   ,fTrackInfo(0)
  77   ,ftrdTrack(0)
  78   ,fcl(0)
  79   ,fTRDCalibraFillHisto(0)
  80   ,fNbTRDTrack(0)
  81   ,fNbTRDTrackOffline(0)
  82   ,fNbTRDTrackStandalone(0)
  83   ,fNbTRDTracklet(0)
  84   ,fNbTRDTrackletOffline(0)
  85   ,fNbTRDTrackletStandalone(0)
  86   ,fNbTimeBin(0x0)
  87   ,fNbTimeBinOffline(0x0)
  88   ,fNbTimeBinStandalone(0x0)
  89   ,fNbClusters(0)
  90   ,fNbClustersOffline(0)
  91   ,fNbClustersStandalone(0)
  92   ,fPHSum(0)
  93   ,fCHSum(0)
  94   ,fDetSum(0)
  95   ,fDetSumVector(0)
  96   ,fHisto2d(kTRUE)
  97   ,fVector2d(kFALSE)
  98   ,fVdriftLinear(kTRUE)
  99   ,flow(0)
 100   ,fhigh(30)
 101   ,fNbTimeBins(30)
 102   ,ffillZero(kFALSE)
 103   ,fnormalizeNbOfCluster(kFALSE)
 104   ,fmaxCluster(0)
 105   ,fOfflineTracks(kFALSE)
 106   ,fStandaloneTracks(kFALSE)
 107   ,fCompressPerDetector(kFALSE)
 108   ,fRunNumber(0)
 109   ,fkNameDirectory("local://.")
 110   ,fGraph(0x0)
 111   ,fPostProcess(kFALSE)
 112 {
 113   // Constructor
 114
 115   fNRefFigures = 17;
 116
 117   for(Int_t k = 0; k < 3; k++)
 118     {
 119       fNz[k]=0;
 120       fNrphi[k]=0;
 121     }
 122
 123 }
 124 //________________________________________________________________________
 125 AliTRDcalibration::~AliTRDcalibration()
 126 {
 127   // Default destructor
 128
 129   if(fNbTRDTrack) delete fNbTRDTrack;
 130   if(fNbTRDTrackOffline) delete fNbTRDTrackOffline;
 131   if(fNbTRDTrackStandalone) delete fNbTRDTrackStandalone;
 132   if(fNbTRDTracklet) delete fNbTRDTracklet;
 133   if(fNbTRDTrackletOffline) delete fNbTRDTrackletOffline;
 134   if(fNbTRDTrackletStandalone) delete fNbTRDTrackletStandalone;
 135   if(fNbTimeBin) delete fNbTimeBin;
 136   if(fNbTimeBinOffline) delete fNbTimeBinOffline;
 137   if(fNbTimeBinStandalone) delete fNbTimeBinStandalone;
 138   if(fNbClusters) delete fNbClusters;
 139   if(fNbClustersOffline) delete fNbClustersOffline;
 140   if(fNbClustersStandalone) delete fNbClustersStandalone;
 141   if(fPHSum) delete fPHSum;
 142   if(fCHSum) delete fCHSum;
 143   if(fDetSum) delete fDetSum;
 144   if(fDetSumVector) delete fDetSumVector;
 145   if(fkNameDirectory) delete fkNameDirectory;
 146   if(fGraph){fGraph->Delete(); delete fGraph;}
 147
 148 }
 149 //________________________________________________________________________
 150 void AliTRDcalibration::CreateOutputObjects()
 151 {
 152   // Create output objects
 153
 154   OpenFile(0, "RECREATE");
 155
 156   // Number of time bins
 157   AliTRDcalibDB *cal = AliTRDcalibDB::Instance();
 158   fNbTimeBins = cal->GetNumberOfTimeBins();
 159
 160   // instance calibration: what to calibrate
 161   fTRDCalibraFillHisto = AliTRDCalibraFillHisto::Instance();
 162   fTRDCalibraFillHisto->SetHisto2d(fHisto2d); // choose to use histograms
 163   fTRDCalibraFillHisto->SetVector2d(fVector2d); // choose to use vectors
 164   fTRDCalibraFillHisto->SetCH2dOn();  // choose to calibrate the gain
 165   fTRDCalibraFillHisto->SetPH2dOn();  // choose to calibrate the drift velocity
 166   fTRDCalibraFillHisto->SetPRF2dOn(); // choose to look at the PRF
 167   fTRDCalibraFillHisto->SetLinearFitterOn(fVdriftLinear); // Other possibility vdrift VDRIFT
 168   fTRDCalibraFillHisto->SetLinearFitterDebugOn(fVdriftLinear); // Other possibility vdrift
 169
 170   // segmentation (should be per default the max and add at the end)
 171   for(Int_t k = 0; k < 3; k++){
 172     if(((fNz[k] != 10) && (fNrphi[k] != 10)) && ((fNz[k] != 100) && (fNrphi[k] != 100))) {
 173       fTRDCalibraFillHisto->SetNz(k,fNz[k]);                                    // Mode calibration
 174       fTRDCalibraFillHisto->SetNrphi(k,fNrphi[k]);                             // Mode calibration
 175     }
 176     else {
 177       if((fNz[k] == 100) && (fNrphi[k] == 100))  {
 178         if(fVector2d) AliInfo("The mode all together is not supported by the vector method");
 179         fTRDCalibraFillHisto->SetAllTogether(k);
 180       }
 181       if((fNz[k] == 10) && (fNrphi[k] == 10))  {
 182         if(fVector2d) AliInfo("The mode per supermodule is not supported by the vector method");
 183         fTRDCalibraFillHisto->SetPerSuperModule(k);
 184       }
 185     }
 186   }
 187
 188   // Debug level
 189   fTRDCalibraFillHisto->SetDebugLevel(fDebugLevel); //debug stuff
 190
 191   // Init the stuff
 192   fTRDCalibraFillHisto->Init2Dhistos(); // initialise the histos
 193
 194   // cuts
 195   fTRDCalibraFillHisto->SetNumberClusters(flow); // At least flow clusters
 196   fTRDCalibraFillHisto->SetNumberClustersf(fhigh); // The more fhigh clusters
 197   fTRDCalibraFillHisto->SetFillWithZero(ffillZero); // Fill zeros
 198   fTRDCalibraFillHisto->SetNormalizeNbOfCluster(fnormalizeNbOfCluster); // For iterations
 199
 200   // Add them to the container
 201   fContainer = new TObjArray();
 202   if(fHisto2d) {
 203     fContainer->Add(fTRDCalibraFillHisto->GetCH2d()); //TH2I
 204     fContainer->Add(fTRDCalibraFillHisto->GetPH2d()); //TProfile2D
 205     fContainer->Add(fTRDCalibraFillHisto->GetPRF2d()); //TProfile2D
 206   }
 207   if(fVdriftLinear) fContainer->Add(fTRDCalibraFillHisto->GetVdriftLinearFit()); // Other drift velocity
 208   if(fVector2d) fContainer->Add(fTRDCalibraFillHisto->GetCalibraVector()); //calibra vector
 209
 210   if(fDebugLevel) {
 211
 212     // Init the debug histos
 213     fNbTRDTrack = new TH1F("TRDTrack","TRDTrack",500,0,500);
 214     fNbTRDTrack->Sumw2();
 215     fNbTRDTrackOffline = new TH1F("TRDTrackOffline","TRDTrackOffline",500,0,500);
 216     fNbTRDTrackOffline->Sumw2();
 217     fNbTRDTrackStandalone = new TH1F("TRDTrackStandalone","TRDTrackStandalone",500,0,500);
 218     fNbTRDTrackStandalone->Sumw2();
 219     //
 220     fNbTRDTracklet = new TH1F("TRDTracklet","TRDTracklet",540,0.,540.);
 221     fNbTRDTracklet->Sumw2();
 222     fNbTRDTrackletOffline = new TH1F("TRDTrackletOffline","TRDTrackletOffline",540,0.,540.);
 223     fNbTRDTrackletOffline->Sumw2();
 224     fNbTRDTrackletStandalone = new TH1F("TRDTrackletStandalone","TRDTrackletStandalone",540,0.,540.);
 225     fNbTRDTrackletStandalone->Sumw2();
 226     //
 227     fNbTimeBin = new TH1F("TimeBin","TimeBin",35,0,35);
 228     fNbTimeBin->Sumw2();
 229     fNbTimeBinOffline = new TH1F("TimeBinOffline","TimeBinOffline",35,0,35);
 230     fNbTimeBinOffline->Sumw2();
 231     fNbTimeBinStandalone = new TH1F("TimeBinStandalone","TimeBinStandalone",35,0,35);
 232     fNbTimeBinStandalone->Sumw2();
 233     //
 234     fNbClusters = new TH1F("NbClusters","",35,0,35);
 235     fNbClusters->Sumw2();
 236     fNbClustersOffline = new TH1F("NbClustersOffline","",35,0,35);
 237     fNbClustersOffline->Sumw2();
 238     fNbClustersStandalone = new TH1F("NbClustersStandalone","",35,0,35);
 239     fNbClustersStandalone->Sumw2();
 240     //
 241     fPHSum = new TProfile2D("PH2dSum","Nz0Nrphi0"
 242                             ,fNbTimeBins,-0.05,(Double_t)(fNbTimeBins/10.0)-0.05
 243                             ,540,0,540);
 244     fPHSum->SetYTitle("Det/pad groups");
 245     fPHSum->SetXTitle("time [#mus]");
 246     fPHSum->SetZTitle("<PH> [a.u.]");
 247     fPHSum->SetStats(0);
 248     //
 249     fCHSum = new TH2I("CH2dSum","Nz0Nrphi0",100,0,300,540,0,540);
 250     fCHSum->SetYTitle("Det/pad groups");
 251     fCHSum->SetXTitle("charge deposit [a.u]");
 252     fCHSum->SetZTitle("counts");
 253     fCHSum->SetStats(0);
 254     fCHSum->Sumw2();
 255
 256     // Add them
 257     fContainer->Add(fNbTRDTrack);
 258     fContainer->Add(fNbTRDTrackOffline);
 259     fContainer->Add(fNbTRDTrackStandalone);
 260     fContainer->Add(fNbTRDTracklet);
 261     fContainer->Add(fNbTRDTrackletOffline);
 262     fContainer->Add(fNbTRDTrackletStandalone);
 263     fContainer->Add(fNbTimeBin);
 264     fContainer->Add(fNbTimeBinOffline);
 265     fContainer->Add(fNbTimeBinStandalone);
 266     fContainer->Add(fNbClusters);
 267     fContainer->Add(fNbClustersOffline);
 268     fContainer->Add(fNbClustersStandalone);
 269     fContainer->Add(fPHSum);
 270     fContainer->Add(fCHSum);
 271
 272   }
 273
 274 }
 275
 276 //________________________________________________________________________
 277 void AliTRDcalibration::Exec(Option_t *)
 278 {
 279   //
 280   // Execute function where the reference data are filled
 281   //
 282
 283   if(!fTracks) return;
 284
 285   // In total
 286   Int_t nbTrdTracks = 0;
 287   // standalone
 288   Int_t nbTrdTracksStandalone = 0;
 289   // offline
 290   Int_t nbTrdTracksOffline = 0;
 291
 292
 293   //
 294   // Loop on track in the event
 295   //
 296   //printf("Total of %d\n",fTracks->GetEntriesFast());
 297   for(Int_t itrk=0; itrk < fTracks->GetEntriesFast(); itrk++){
 298
 299     //printf("itrk %d\n",itrk);
 300
 301     fTrackInfo = (AliTRDtrackInfo*)fTracks->UncheckedAt(itrk);
 302     ftrdTrack = fTrackInfo->GetTrack();
 303     if(!ftrdTrack) continue;
 304
 305     nbTrdTracks++;
 306
 307     fTRDCalibraFillHisto->UpdateHistogramsV1(ftrdTrack);
 308
 309     if(fDebugLevel) {
 310
 311       Bool_t standalonetracklet = kFALSE;
 312       const AliTRDseedV1 *tracklet = 0x0;
 313       //
 314       // Loop on tracklet in the event
 315       //
 316       for(Int_t itr = 0; itr < 6; itr++){
 317         //printf("itr %d\n",itr);
 318         if(!(tracklet = ftrdTrack->GetTracklet(itr))) continue;
 319         if(!tracklet->IsOK()) continue;
 320         // standalone
 321         if(tracklet->IsStandAlone()) standalonetracklet = kTRUE;
 322         Int_t nbclusters = 0;
 323         // For PH
 324         Double_t phtb[AliTRDseedV1::kNtb];
 325         memset(phtb, 0, AliTRDseedV1::kNtb*sizeof(Double_t));
 326         // For CH
 327         Double_t sum = 0.0;
 328         // normalisation
 329         Float_t normalisation = 6.67;
 330         Int_t detector = 0;
 331         Int_t crossrow = 0;
 332         for(int ic=0; ic<AliTRDseedV1::kNtb; ic++){
 333           // Check no shared clusters
 334           for(int icc=AliTRDseedV1::kNtb; icc<AliTRDseedV1::kNclusters; icc++){
 335             if((fcl = tracklet->GetClusters(icc)))  crossrow = 1;
 336           }
 337           if(!(fcl = tracklet->GetClusters(ic))) continue;
 338           nbclusters++;
 339           Int_t time = fcl->GetPadTime();
 340           Float_t ch =  tracklet->GetdQdl(ic);
 341           Float_t qcl = TMath::Abs(fcl->GetQ());
 342           detector = fcl->GetDetector();
 343           if((time>-1) && (time<fNbTimeBins)) phtb[time]=qcl;
 344           sum += ch/normalisation;
 345           fNbTimeBin->Fill(time);
 346           if(tracklet->IsStandAlone()) fNbTimeBinStandalone->Fill(time);
 347           else fNbTimeBinOffline->Fill(time);
 348         }
 349
 350         fNbTRDTracklet->Fill(detector);
 351         if(tracklet->IsStandAlone()) fNbTRDTrackletStandalone->Fill(detector);
 352         else fNbTRDTrackletOffline->Fill(detector);
 353
 354         fNbClusters->Fill(nbclusters);
 355         if(tracklet->IsStandAlone())  fNbClustersStandalone->Fill(nbclusters);
 356         else  fNbClustersOffline->Fill(nbclusters);
 357
 358         if((nbclusters > flow) && (nbclusters < fhigh)){
 359           fCHSum->Fill(sum/20.0,0.0);
 360           for(int ic=0; ic<fNbTimeBins; ic++){
 361             //printf("ic %d and time %f and cluster %f \n",ic,(Double_t)(ic/10.0),(Double_t)phtb[ic]);
 362             if(ffillZero) fPHSum->Fill((Double_t)(ic/10.0),0.0,(Double_t)phtb[ic]);
 363             else {
 364               if(phtb[ic] > 0.0) fPHSum->Fill((Double_t)(ic/10.0),0.0,(Double_t)phtb[ic]);
 365             }
 366           }
 367         }
 368       }
 369
 370     if(standalonetracklet) nbTrdTracksStandalone++;
 371     else nbTrdTracksOffline++;
 372
 373     }
 374
 375   }
 376
 377   if(fDebugLevel) {
 378
 379     //Fill Histos
 380     fNbTRDTrack->Fill(nbTrdTracks);
 381     fNbTRDTrackStandalone->Fill(nbTrdTracksStandalone);
 382     fNbTRDTrackOffline->Fill(nbTrdTracksOffline);
 383
 384   }
 385
 386   //printf("Nbof tracks %d\n",nbTrdTracks);
 387
 388   TFile *file = TFile::Open("test.root","RECREATE");
 389   fContainer->Write();
 390   file->Close();
 391
 392   // Post output data
 393   PostData(0, fContainer);
 394
 395   //printf("post container\n");
 396
 397 }
 398 //________________________________________________________________________
 399 void AliTRDcalibration::Terminate(Option_t *)
 400 {
 401   // Draw result to the screen
 402   // Called once at the end of the query
 403
 404   //printf("terminate\n");
 405
 406   if(fTRDCalibraFillHisto) fTRDCalibraFillHisto->DestroyDebugStreamer();
 407
 408 }
 409 //________________________________________________________
 410 Bool_t AliTRDcalibration::GetRefFigure(Int_t ifig)
 411 {
 412   //
 413   // Draw filled histos
 414   //
 415
 416   gStyle->SetPalette(1);
 417   gStyle->SetOptStat(1111);
 418   gStyle->SetPadBorderMode(0);
 419   gStyle->SetCanvasColor(10);
 420   gStyle->SetPadLeftMargin(0.13);
 421   gStyle->SetPadRightMargin(0.13);
 422
 423   if(!fContainer) return kFALSE;
 424
 425   switch(ifig){
 426   case kNbTrack:{
 427     TCanvas *c0 = new TCanvas("c0","c0",10,10,510,510);
 428     TLegend *legNbTrack = new TLegend(.7, .7, .98, .98);
 429     legNbTrack->SetBorderSize(1);
 430     TH1F *h  = 0x0;
 431     TH1F *ho = 0x0;
 432     TH1F *hs = 0x0;
 433     if(!(h = (TH1F *)fContainer->FindObject("TRDTrack"))) break;
 434     if(!(ho = (TH1F *)fContainer->FindObject("TRDTrackOffline"))) break;
 435     if(!(hs = (TH1F *)fContainer->FindObject("TRDTrackStandalone"))) break;
 436     c0->cd();
 437     //gPad->SetLogy();
 438     gPad->SetGridy();
 439     gPad->SetGridx();
 440     h->Draw();
 441     ho->Draw("same");
 442     hs->Draw("same");
 443     legNbTrack->AddEntry(h, "all", "p");
 444     legNbTrack->AddEntry(ho, "offline", "p");
 445     legNbTrack->AddEntry(hs, "standalone", "p");
 446     legNbTrack->Draw("same");
 447     return kTRUE;
 448   }
 449   case kNbTracklet:{
 450     TLegend *legNbTracklet = new TLegend(.7, .7, .98, .98);
 451     legNbTracklet->SetBorderSize(1);
 452     TH1F *h = 0x0;
 453     TH1F *ho = 0x0;
 454     TH1F *hs = 0x0;
 455     if(!(h = (TH1F *)fContainer->FindObject("TRDTracklet"))) break;
 456     if(!(ho = (TH1F *)fContainer->FindObject("TRDTrackletOffline"))) break;
 457     if(!(hs = (TH1F *)fContainer->FindObject("TRDTrackletStandalone"))) break;
 458     h->Draw();
 459     ho->Draw("same");
 460     hs->Draw("same");
 461     legNbTracklet->AddEntry(h, "all", "p");
 462     legNbTracklet->AddEntry(ho, "offline", "p");
 463     legNbTracklet->AddEntry(hs, "standalone", "p");
 464     legNbTracklet->Draw("same");
 465     gPad->SetLogy();
 466     gPad->SetGridy();
 467     gPad->SetGridx();
 468     return kTRUE;
 469   }
 470   case kNbTimeBin:{
 471     TLegend *legNbTimeBin = new TLegend(.7, .7, .98, .98);
 472     legNbTimeBin->SetBorderSize(1);
 473     TH1F *h = 0x0;
 474     TH1F *ho = 0x0;
 475     TH1F *hs = 0x0;
 476     if(!(h = (TH1F *)fContainer->FindObject("TimeBin"))) break;
 477     if(!(ho = (TH1F *)fContainer->FindObject("TimeBinOffline"))) break;
 478     if(!(hs = (TH1F *)fContainer->FindObject("TimeBinStandalone"))) break;
 479     h->Draw();
 480     ho->Draw("same");
 481     hs->Draw("same");
 482     legNbTimeBin->AddEntry(h, "all", "p");
 483     legNbTimeBin->AddEntry(ho, "offline", "p");
 484     legNbTimeBin->AddEntry(hs, "standalone", "p");
 485     legNbTimeBin->Draw("same");
 486     //gPad->SetLogy();
 487     gPad->SetGridy();
 488     gPad->SetGridx();
 489     return kTRUE;
 490   }
 491   case kNbClusters:{
 492     TLegend *legNbClusters = new TLegend(.7, .7, .98, .98);
 493     legNbClusters->SetBorderSize(1);
 494     TH1F *h = 0x0;
 495     TH1F *ho = 0x0;
 496     TH1F *hs = 0x0;
 497     if(!(h = (TH1F *)fContainer->FindObject("NbClusters"))) break;
 498     if(!(ho = (TH1F *)fContainer->FindObject("NbClustersOffline"))) break;
 499     if(!(hs = (TH1F *)fContainer->FindObject("NbClustersStandalone"))) break;
 500     h->Draw();
 501     ho->Draw("same");
 502     hs->Draw("same");
 503     legNbClusters->AddEntry(h, "all", "p");
 504     legNbClusters->AddEntry(ho, "offline", "p");
 505     legNbClusters->AddEntry(hs, "standalone", "p");
 506     legNbClusters->Draw("same");
 507     gPad->SetLogy();
 508     gPad->SetGridy();
 509     gPad->SetGridx();
 510     return kTRUE;
 511   }
 512   case kPHSum:{
 513     TProfile2D *h = 0x0;
 514     if(!(h = (TProfile2D *)fContainer->FindObject("PH2dSum"))) break;
 515     TH1D *projh = h->ProjectionX("projh",1,1,"e");
 516     projh->Draw();
 517     gPad->SetGridy();
 518     gPad->SetGridx();
 519     return kTRUE;
 520   }
 521   case kCHSum:{
 522     TH2I *h = 0x0;
 523     if(!(h = (TH2I *)fContainer->FindObject("CH2dSum"))) break;
 524     TH1D *projh = h->ProjectionX("projhh",1,1,"e");
 525     projh->Draw();
 526     gPad->SetGridy();
 527     gPad->SetGridx();
 528     return kTRUE;
 529   }
 530   case kPH2D:{
 531     if(!fHisto2d) {
 532       AliInfo("Histo was not filled!");
 533       break;
 534     }
 535     TProfile2D *h = 0x0;
 536     if(!(h = (TProfile2D *)fContainer->FindObject("PH2d"))) break;
 537     h->Draw("lego");
 538     return kTRUE;
 539   }
 540   case kCH2D:{
 541     if(!fHisto2d) {
 542       AliInfo("Histo was not filled!");
 543       break;
 544     }
 545     TH2I *h = 0x0;
 546     if(!(h = (TH2I *)fContainer->FindObject("CH2d"))) break;
 547     h->Draw("lego");
 548     return kTRUE;
 549   }
 550   case kPRF2D:{
 551     if(!fHisto2d) {
 552       AliInfo("Histo was not filled!");
 553       break;
 554     }
 555     TProfile2D *h = 0x0;
 556     if(!(h = (TProfile2D *)fContainer->FindObject("PRF2d"))) break;
 557     h->Draw("lego");
 558     return kTRUE;
 559   }
 560   case kPH2DVector:{
 561     if(!fVector2d) {
 562       AliInfo("vector was not filled!");
 563       break;
 564     }
 565     AliTRDCalibraVector *v = 0x0;
 566     TGraphErrors *vdet = 0x0;
 567     if(!(v = (AliTRDCalibraVector *)fContainer->FindObject("AliTRDCalibraVector"))) break;
 568     Int_t detectormax = -1;
 569     Int_t groupmax    = -1;
 570     if(!v->FindTheMaxEntries(1,detectormax,groupmax)) break;
 571     if(!(vdet = v->ConvertVectorPHTGraphErrors((Int_t)detectormax,groupmax,"plotPH2dVector"))) break;
 572     Int_t nbeentries = 0;
 573     TH1F *ko = v->CorrectTheError(vdet,nbeentries);
 574     ko->Draw();
 575     AliInfo(Form("There are %d entries in the detector %d and group %d",nbeentries,detectormax,groupmax));
 576     return kTRUE;
 577   }
 578 case kCH2DVector:{
 579     if(!fVector2d) {
 580       AliInfo("vector was not filled!");
 581       break;
 582     }
 583     AliTRDCalibraVector *v = 0x0;
 584     TH1F *vdet = 0x0;
 585     if(!(v = (AliTRDCalibraVector *)fContainer->FindObject("AliTRDCalibraVector"))) break;
 586     Int_t detectormax = -1;
 587     Int_t groupmax    = -1;
 588     if(!v->FindTheMaxEntries(0,detectormax,groupmax)) break;
 589     if(!(vdet = v->ConvertVectorCHHisto((Int_t)detectormax,groupmax,"plotCH2dVector"))) break;
 590     vdet->Draw();
 591     AliInfo(Form("The detectormax and groupmax are %d and %d",detectormax,groupmax));
 592     return kTRUE;
 593   }
 594   case kPRF2DVector:{
 595     if(!fVector2d) {
 596       AliInfo("vector was not filled!");
 597       break;
 598     }
 599     AliTRDCalibraVector *v = 0x0;
 600     TGraphErrors *vdet = 0x0;
 601     if(!(v = (AliTRDCalibraVector *)fContainer->FindObject("AliTRDCalibraVector"))) break;
 602     Int_t detectormax  = -1;
 603     Int_t groupmax     = -1;
 604     Int_t nbeentries   = 0;
 605     if(!v->FindTheMaxEntries(2,detectormax,groupmax)) break;
 606     if(!(vdet = v->ConvertVectorPRFTGraphErrors((Int_t)detectormax,groupmax,"plotPRF2dVector"))) break;
 607     TH1F *ko = v->CorrectTheError(vdet,nbeentries);
 608     ko->Draw();
 609     AliInfo(Form("The detectormax and groupmax are %d and %d",detectormax,groupmax));
 610     return kTRUE;
 611   }
 612   case kLinearFitter:{
 613     if(!fVdriftLinear) {
 614       AliInfo("vdrift linear was not filled!");
 615       break;
 616     }
 617     AliTRDCalibraVdriftLinearFit *h = 0x0;
 618     TH2F *hdetector = 0x0;
 619     if(!(h = (AliTRDCalibraVdriftLinearFit *)fContainer->FindObject("AliTRDCalibraVdriftLinearFit"))) break;
 620     Double_t entries[540];
 621     for(Int_t k = 0; k < 540; k++){
 622       entries[k] = 0.0;
 623       hdetector = 0x0;
 624       if(!(hdetector = h->GetLinearFitterHisto(k,kFALSE))) continue;
 625       entries[k] = hdetector->GetEntries();
 626     }
 627     Double_t max = -10.0;
 628     Double_t detectormax = -1;
 629     for(Int_t k = 0; k < 540; k++){
 630       if(entries[k] > max) {
 631         max = entries[k];
 632         detectormax = k;
 633       }
 634     }
 635     hdetector = 0x0;
 636     if((max == 0.0) || (detectormax <0.0) || (detectormax >=540.0)) break;
 637     if(!(hdetector = h->GetLinearFitterHisto((Int_t)detectormax,kFALSE))) break;
 638     AliInfo(Form("The detector with the maximum of entries is %d",detectormax));
 639     hdetector->Draw();
 640     return kTRUE;
 641   }
 642   case kGainFactor:{
 643     if(!fPostProcess){
 644       if(!PostProcess()) break;
 645     }
 646     TGraph *fgain = (TGraph *) fGraph->At(0);
 647     if(!fgain) break;
 648     fgain->Draw("ALP");
 649     return kTRUE;
 650   }
 651   case kVdriftT0Factor:{
 652     if(!fPostProcess){
 653       if(!PostProcess()) break;
 654     }
 655     TCanvas *c = new TCanvas("c","c",10,10,510,510);
 656     c->Divide(2,1);
 657     TGraph *fvd = (TGraph *) fGraph->At(1);
 658     if(fvd){
 659       c->cd(1);
 660       fvd->Draw("ALP");
 661     }
 662     TGraph *ft0 = (TGraph *) fGraph->At(2);
 663     if(ft0){
 664       c->cd(2);
 665       ft0->Draw("ALP");
 666     }
 667     return kTRUE;
 668   }
 669   case kVdriftLorentzAngleFactor:{
 670     if(!fVdriftLinear) {
 671       AliInfo("vdrift linear was not filled!");
 672       break;
 673     }
 674     if(!fPostProcess){
 675       if(!PostProcess()) break;
 676     }
 677     TCanvas *c = new TCanvas("c","c",10,10,510,510);
 678     c->Divide(2,1);
 679     TGraph *fvdl = (TGraph *) fGraph->At(3);
 680     if(fvdl){
 681       c->cd(1);
 682       fvdl->Draw("ALP");
 683     }
 684     TGraph *flal = (TGraph *) fGraph->At(4);
 685     if(flal){
 686       c->cd(2);
 687       flal->Draw("ALP");
 688     }
 689     return kTRUE;
 690   }
 691   case kPRFFactor:{
 692     if(!fPostProcess){
 693       if(!PostProcess()) break;
 694     }
 695     TGraph *fprf = (TGraph *) fGraph->At(5);
 696     if(!fprf) break;
 697     fprf->Draw("ALP");
 698     return kTRUE;
 699   }
 700   }
 701
 702   return kFALSE;
 703
 704 }
 705 //________________________________________________________________________
 706 Bool_t AliTRDcalibration::PostProcess()
 707 {
 708   //
 709   // Fit the filled histos
 710   //
 711
 712   if(!fGraph){
 713     fGraph = new TObjArray(6);
 714     fGraph->SetOwner();
 715   }
 716   else {
 717     delete fGraph;
 718     PostProcess();
 719   }
 720
 721   Bool_t storage[3] = {kFALSE,kFALSE,kFALSE};
 722
 723   // storage element
 724   AliCDBManager *man = AliCDBManager::Instance();
 725   man->SetDefaultStorage("local://$ALICE_ROOT");
 726   AliCDBStorage* storLoc = man->GetStorage(fkNameDirectory);
 727   if (!storLoc)
 728     return kFALSE;
 729   man->SetRun(fRunNumber);
 730
 731   // MetaData
 732   AliCDBMetaData mdDet;
 733   mdDet.SetObjectClassName("AliTRDCalDet");
 734   AliCDBMetaData mdPad;
 735   mdPad.SetObjectClassName("AliTRDCalPad");
 736
 737   // Objects for fitting
 738   AliTRDCalibraFit *calibra = AliTRDCalibraFit::Instance();
 739   calibra->SetDebugLevel(2); // 0 rien, 1 fitvoir, 2 debug files, 3 one detector
 740
 741   // Take the stuff
 742   if (!fContainer) {
 743     Printf("ERROR: list not available");
 744     return kFALSE;
 745   }
 746
 747   if(fHisto2d && fVector2d) AliInfo("We will only look at histos. Set fHisto2d off if you don't want");
 748   AliTRDCalibraVector *calibraVector = 0x0;
 749   if(fVector2d) calibraVector = (AliTRDCalibraVector *) fContainer->FindObject("CalibraVector");
 750   //
 751   // GAIN TH2I
 752   //
 753   Bool_t pass = kFALSE;
 754   AliTRDCalibraVector *vvect = 0x0;
 755   if(fHisto2d) {
 756     TH2I *histogain = (TH2I *) fContainer->FindObject("CH2d");
 757     if(histogain) {
 758       histogain->SetDirectory(0);
 759       calibra->SetMinEntries(20);
 760       if(fCompressPerDetector){
 761         if(AddStatsPerDetector(histogain)) pass = calibra->AnalyseCH(fCHSum);
 762       }
 763       else pass = calibra->AnalyseCH(histogain);
 764     }
 765   }
 766   else {
 767     if(fVector2d && calibraVector) {
 768       calibra->SetMinEntries(20);
 769       if(fCompressPerDetector){
 770         if(!(vvect = calibraVector->AddStatsPerDetectorCH())) return kFALSE;
 771         pass = calibra->AnalyseCH(vvect);
 772       }
 773       else pass = calibra->AnalyseCH(calibraVector);
 774     }
 775   }
 776
 777   if(pass)
 778     {
 779       Int_t nbtg = 6*4*18*((Int_t) ((AliTRDCalibraMode *)calibra->GetCalibraMode())->GetDetChamb0(0))
 780         + 6*  18*((Int_t) ((AliTRDCalibraMode *)calibra->GetCalibraMode())->GetDetChamb2(0));
 781       Int_t nbfit = calibra->GetNumberFit();  //Number of fits
 782       Int_t nbE   = calibra->GetNumberEnt();  //Number of detector mit entries
 783       // enough statistics
 784       if ((nbtg >                  0) &&
 785           (nbfit        >= 0.001*nbE))
 786         {
 787           // create the cal objects
 788           calibra->PutMeanValueOtherVectorFit(1,kTRUE);
 789           TObjArray object           = calibra->GetVectorFit();
 790           AliTRDCalDet *objgaindet   = calibra->CreateDetObjectGain(&object);
 791           TObject *objgainpad        = calibra->CreatePadObjectGain();
 792           // store
 793           AliCDBId id1("TRD/Calib/ChamberGainFactor",fRunNumber, AliCDBRunRange::Infinity());
 794           storLoc->Put((TObject *)objgaindet, id1, &mdDet);
 795           AliCDBId id2("TRD/Calib/LocalGainFactor",fRunNumber, AliCDBRunRange::Infinity());
 796           storLoc->Put((TObject *)objgainpad, id2, &mdPad);
 797           storage[0] = kTRUE;
 798           // Make graph
 799           TGraph *graph = 0x0;
 800           if(FillGraphIndex(&object,graph)){
 801             fGraph->AddAt(graph,0);
 802           }
 803         }//if(enough statistics?)
 804       calibra->ResetVectorFit();
 805     }
 806   else return kFALSE;
 807
 808   //
 809   // VDRIFT average pulse height
 810   //
 811   pass = kFALSE;
 812   if(fHisto2d) {
 813     TProfile2D *histodriftvelocity = (TProfile2D *) fContainer->FindObject("PH2d");
 814     if(histodriftvelocity) {
 815       histodriftvelocity->SetDirectory(0);
 816       calibra->SetMinEntries(20*20);
 817       if(fCompressPerDetector){
 818         if(AddStatsPerDetector(histodriftvelocity,1)) {
 819           pass = calibra->AnalysePH(fDetSumVector);
 820         }
 821       }
 822       else pass = calibra->AnalysePH(histodriftvelocity);
 823     }
 824   }
 825   else {
 826     if(fVector2d && calibraVector) {
 827       calibra->SetMinEntries(20*20);
 828       if(fCompressPerDetector){
 829         if(!(vvect = calibraVector->AddStatsPerDetectorPH())) return kFALSE;
 830         pass = calibra->AnalysePH(vvect);
 831       }
 832       else pass = calibra->AnalysePH(calibraVector);
 833     }
 834   }
 835
 836   if(pass) {
 837     Int_t nbtg = 6*4*18*((Int_t) ((AliTRDCalibraMode *)calibra->GetCalibraMode())->GetDetChamb0(1))
 838       + 6*  18*((Int_t) ((AliTRDCalibraMode *)calibra->GetCalibraMode())->GetDetChamb2(1));
 839     Int_t nbfit  = calibra->GetNumberFit();
 840     Int_t nbE    = calibra->GetNumberEnt();
 841     // enough statistics
 842     if ((nbtg >                  0) &&
 843         (nbfit        >= 0.001*nbE))
 844       {
 845         // create the cal objects
 846         calibra->PutMeanValueOtherVectorFit(1,kTRUE);
 847         calibra->PutMeanValueOtherVectorFit2(1,kTRUE);
 848         TObjArray object  = calibra->GetVectorFit();
 849         AliTRDCalDet *objdriftvelocitydet = calibra->CreateDetObjectVdrift(&object,kTRUE);
 850         TObject *objdriftvelocitypad      = calibra->CreatePadObjectVdrift();
 851         TObjArray objectt          = calibra->GetVectorFit2();
 852         AliTRDCalDet *objtime0det  = calibra->CreateDetObjectT0(&object,kTRUE);
 853         TObject *objtime0pad       = calibra->CreatePadObjectT0();
 854         // store
 855         AliCDBId id1("TRD/Calib/ChamberVdrift",fRunNumber, AliCDBRunRange::Infinity());
 856         storLoc->Put((TObject *)objdriftvelocitydet, id1, &mdDet);
 857         AliCDBId id2("TRD/Calib/LocalVdrift",fRunNumber, AliCDBRunRange::Infinity());
 858         storLoc->Put((TObject *)objdriftvelocitypad, id2, &mdPad);
 859         //
 860         AliCDBId idd1("TRD/Calib/ChamberT0",fRunNumber, AliCDBRunRange::Infinity());
 861         storLoc->Put((TObject *)objtime0det, idd1, &mdDet);
 862         AliCDBId idd2("TRD/Calib/LocalT0",fRunNumber, AliCDBRunRange::Infinity());
 863         storLoc->Put((TObject *)objtime0pad, idd2, &mdPad);
 864         // Make graph
 865         TGraph *graph = 0x0;
 866         if(FillGraphIndex(&object,graph)){
 867           fGraph->AddAt(graph,1);
 868         }
 869         TGraph *graphh = 0x0;
 870         if(FillGraphIndex(&objectt,graphh)){
 871           fGraph->AddAt(graphh,2);
 872         }
 873       }//if(enough statistics)
 874     calibra->ResetVectorFit();
 875   }
 876   else return kFALSE;
 877
 878   //
 879   // PRF
 880   //
 881   pass = kFALSE;
 882   if(fHisto2d) {
 883     TProfile2D *histoprf = (TProfile2D *) fContainer->FindObject("PRF2d");
 884     if (histoprf) {
 885       histoprf->SetDirectory(0);
 886       calibra->SetMinEntries(600);
 887       if(fCompressPerDetector){
 888         if(AddStatsPerDetector(histoprf,2)) pass = calibra->AnalysePRFMarianFit(fDetSumVector);
 889       }
 890       else pass = calibra->AnalysePRFMarianFit(histoprf);
 891     }
 892   }
 893   else {
 894     if(fVector2d && calibraVector) {
 895       calibra->SetMinEntries(600);
 896       if(fCompressPerDetector){
 897         if(!(vvect =calibraVector->AddStatsPerDetectorPRF())) return kFALSE;
 898         pass = calibra->AnalysePRFMarianFit(vvect);
 899       }
 900       else pass = calibra->AnalysePRFMarianFit(calibraVector);
 901     }
 902   }
 903
 904   if(pass){
 905     Int_t nbtg = 6*4*18*((Int_t) ((AliTRDCalibraMode *)calibra->GetCalibraMode())->GetDetChamb0(2))
 906       + 6*  18*((Int_t) ((AliTRDCalibraMode *)calibra->GetCalibraMode())->GetDetChamb2(2));
 907     Int_t nbfit        = calibra->GetNumberFit();
 908     Int_t nbE          = calibra->GetNumberEnt();
 909     // enough statistics
 910     if ((nbtg >                  0) &&
 911         (nbfit        >= 0.001*nbE)) {
 912       TObjArray object            = calibra->GetVectorFit();
 913       TObject *objPRFpad          = calibra->CreatePadObjectPRF(&object);
 914       // store
 915       AliCDBId id2("TRD/Calib/PRFWidth",fRunNumber, AliCDBRunRange::Infinity());
 916       storLoc->Put((TObject *)objPRFpad, id2, &mdPad);
 917       // Make graph
 918       TGraph *graph = 0x0;
 919       if(FillGraphIndex(&object,graph)){
 920         fGraph->AddAt(graph,5);
 921       }
 922     }
 923     calibra->ResetVectorFit();
 924   }
 925   else return kFALSE;
 926
 927   //
 928   // VDRIFT linear fit
 929   //
 930   AliTRDCalibraVdriftLinearFit *vlinearfit = (AliTRDCalibraVdriftLinearFit *) fContainer->FindObject("LinearVdriftFit");
 931   if (vlinearfit) {
 932     calibra->SetMinEntries(20*20);
 933     if(calibra->AnalyseLinearFitters(vlinearfit)) {
 934
 935       Int_t nbtg = 6*4*18*((Int_t) ((AliTRDCalibraMode *)calibra->GetCalibraMode())->GetDetChamb0(2))
 936         + 6*  18*((Int_t) ((AliTRDCalibraMode *)calibra->GetCalibraMode())->GetDetChamb2(2));
 937       Int_t nbfit  = calibra->GetNumberFit();
 938       Int_t nbE    = calibra->GetNumberEnt();
 939       // enough statistics
 940       if ((nbtg >                  0) &&
 941           (nbfit        >= 0.001*nbE))
 942         {
 943           // create the cal objects
 944           calibra->PutMeanValueOtherVectorFit(1,kTRUE);
 945           calibra->PutMeanValueOtherVectorFit2(1,kTRUE);
 946           TObjArray object  = calibra->GetVectorFit();
 947           AliTRDCalDet *objdriftvelocitydet = calibra->CreateDetObjectVdrift(&object,kTRUE);
 948           TObject *objdriftvelocitypad      = calibra->CreatePadObjectVdrift();
 949           TObjArray objectt          = calibra->GetVectorFit2();
 950           AliTRDCalDet *objtime0det  = calibra->CreateDetObjectT0(&object,kTRUE);
 951           TObject *objtime0pad       = calibra->CreatePadObjectT0();
 952           // store dummy
 953           AliCDBId id1("TRD/Calib/ChamberVdriftLinear",fRunNumber, AliCDBRunRange::Infinity());
 954           storLoc->Put((TObject *)objdriftvelocitydet, id1, &mdDet);
 955           AliCDBId id2("TRD/Calib/LocalVdriftLinear",fRunNumber, AliCDBRunRange::Infinity());
 956           storLoc->Put((TObject *)objdriftvelocitypad, id2, &mdPad);
 957           //
 958           AliCDBId idd1("TRD/Calib/ChamberLorentzAngle",fRunNumber, AliCDBRunRange::Infinity());
 959           storLoc->Put((TObject *)objtime0det, idd1, &mdDet);
 960           AliCDBId idd2("TRD/Calib/LocalLorentzAngle",fRunNumber, AliCDBRunRange::Infinity());
 961           storLoc->Put((TObject *)objtime0pad, idd2, &mdPad);
 962           // Make graph
 963           TGraph *graph = 0x0;
 964           if(FillGraphIndex(&object,graph)){
 965             fGraph->AddAt(graph,3);
 966           }
 967           TGraph *graphh = 0x0;
 968           if(FillGraphIndex(&objectt,graphh)){
 969             fGraph->AddAt(graphh,4);
 970           }
 971         }//if(enough statistics)
 972     }// if fit
 973     calibra->ResetVectorFit();
 974   }
 975   else return kFALSE;
 976
 977   fPostProcess = kTRUE;
 978
 979   return kTRUE;
 980
 981 }
 982
 983 //________________________________________________________________________
 984 Bool_t AliTRDcalibration::FillGraphIndex(const TObjArray *vectora,TGraph *graph) const
 985 {
 986   //
 987   // Fill one value (the first one) per detector
 988   //
 989
 990   Int_t loop = (Int_t) vectora->GetEntriesFast();
 991   if(loop != 540) {
 992     AliInfo("The Vector Fit is not complete!");
 993     return kFALSE;
 994   }
 995
 996   Double_t x[540];
 997   Double_t y[540];
 998   for (Int_t k = 0; k < loop; k++) {
 999     if(!vectora->At(k)){
1000       x[k] = -1.0;
1001       y[k] = -1.0;
1002       continue;
1003     }
1004     x[k]  = ((AliTRDCalibraFit::AliTRDFitInfo *) vectora->At(k))->GetDetector();
1005     y[k]  = ((Float_t *)((AliTRDCalibraFit::AliTRDFitInfo *) vectora->At(k))->GetCoef())[0];
1006   }
1007
1008   if(!graph) graph = new TGraph(540,&x[0],&y[0]);
1009   else{
1010     graph->~TGraph();
1011     new(graph) TGraph(540,&x[0],&y[0]);
1012   }
1013
1014   return kTRUE;
1015
1016 }
1017 //________________________________________________________________________
1018 Bool_t AliTRDcalibration::AddStatsPerDetector(const TH2I *ch)
1019 {
1020   //
1021   // Add statistic per detector
1022   //
1023
1024   AliTRDCalibraMode calibMode = AliTRDCalibraMode();
1025   const char *name = ch->GetTitle();
1026   if(!SetNzFromTObject(name,0,&calibMode)) return 0x0;
1027   if(!SetNrphiFromTObject(name,0,&calibMode)) return 0x0;
1028   if(((calibMode.GetNz(0) == 100) && (calibMode.GetNrphi(0) == 100)) || ((calibMode.GetNz(0) == 10) && (calibMode.GetNrphi(0) == 10))) return kFALSE;
1029
1030   Int_t    nybins  = ch->GetNbinsY();// groups number
1031   Int_t    nxbins  = ch->GetNbinsX();// number of bins X
1032   TAxis   *xaxis   = ch->GetXaxis();
1033   Double_t lowedge  = xaxis->GetBinLowEdge(1);
1034   Double_t upedge   = xaxis->GetBinUpEdge(nxbins);
1035
1036   // number per chamber 2
1037   calibMode.ModePadCalibration(2,0);
1038   calibMode.ModePadFragmentation(0,2,0,0);
1039   calibMode.SetDetChamb2(0);
1040   Int_t perChamber2 = (Int_t) calibMode.GetDetChamb2(0);
1041
1042   // number per other chamber
1043   calibMode.ModePadCalibration(0,0);
1044   calibMode.ModePadFragmentation(0,0,0,0);
1045   calibMode.SetDetChamb0(0);
1046   Int_t perChamber0 = (Int_t) calibMode.GetDetChamb0(0);
1047
1048   if(nybins != (6*18*perChamber2+6*4*18*perChamber0)) return kFALSE;
1049
1050   // Create Histo
1051   TString nname((const char *)ch->GetName());
1052   nname  += "PerDetector";
1053   TString title("Nz");
1054   title += 0;
1055   title += "Nrphi";
1056   title += 0;
1057   if(!fCHSum) fCHSum = new TH2I((const char *)nname,(const char *)title
1058                                 ,nxbins,lowedge,upedge,540,0,540);
1059   else{
1060     fCHSum->~TH2I();
1061     new(fCHSum) TH2I((const Char_t *) nname,(const char *)title
1062                      ,nxbins,lowedge,upedge,540,0,540);
1063   }
1064   fCHSum->SetYTitle("Detector number");
1065   fCHSum->SetXTitle("charge deposit [a.u]");
1066   fCHSum->SetZTitle("counts");
1067   fCHSum->SetStats(0);
1068   fCHSum->Sumw2();
1069
1070   Int_t counter = 0;
1071
1072   for(Int_t det = 0; det < 540; det++){
1073
1074     Int_t numberofgroup = 0;
1075     if(AliTRDgeometry::GetStack(det) == 2) numberofgroup = perChamber2;
1076     else numberofgroup = perChamber0;
1077     TH1I *projch = (TH1I *) ch->ProjectionX("projch",counter+1,counter+numberofgroup,(Option_t *)"e");
1078     projch->SetDirectory(0);
1079
1080     for(Int_t nx = 0; nx <= nxbins; nx++) {
1081       fCHSum->SetBinContent(nx,det+1,projch->GetBinContent(nx));
1082       fCHSum->SetBinError(nx,det+1,projch->GetBinError(nx));
1083     }
1084
1085     counter += numberofgroup;
1086
1087     delete projch;
1088
1089   }
1090
1091   return kTRUE;
1092
1093 }
1094 //_____________________________________________________________________________________________________________________
1095 Bool_t AliTRDcalibration::AddStatsPerDetector(const TProfile2D *ph,Int_t i)
1096 {
1097   //
1098   // Add statistic per detector
1099   //
1100
1101   AliTRDCalibraMode calibMode = AliTRDCalibraMode();
1102   const char *name = ph->GetTitle();
1103   //printf("name %s\n",name);
1104   if(!SetNzFromTObject(name,0,&calibMode)) return kFALSE;
1105   if(!SetNrphiFromTObject(name,0,&calibMode)) return kFALSE;
1106   if(((calibMode.GetNz(0) == 100) && (calibMode.GetNrphi(0) == 100)) || ((calibMode.GetNz(0) == 10) && (calibMode.GetNrphi(0) == 10))) return kFALSE;
1107   //printf("Found mode Mz %d, Nrphi %d\n",calibMode.GetNz(0),calibMode.GetNrphi(0));
1108
1109
1110   Int_t    nybins  = ph->GetNbinsY();// groups number
1111   Int_t    nxbins  = ph->GetNbinsX();// number of bins X
1112   TAxis   *xaxis = ph->GetXaxis();
1113   Double_t lowedge  = xaxis->GetBinLowEdge(1);
1114   Double_t upedge   = xaxis->GetBinUpEdge(nxbins);
1115
1116   // number per chamber 2
1117   calibMode.ModePadCalibration(2,0);
1118   calibMode.ModePadFragmentation(0,2,0,0);
1119   calibMode.SetDetChamb2(0);
1120   Int_t perChamber2 = (Int_t) calibMode.GetDetChamb2(0);
1121
1122   // number per other chamber
1123   calibMode.ModePadCalibration(0,0);
1124   calibMode.ModePadFragmentation(0,0,0,0);
1125   calibMode.SetDetChamb0(0);
1126   Int_t perChamber0 = (Int_t) calibMode.GetDetChamb0(0);
1127
1128   if(nybins != (6*18*perChamber2+6*4*18*perChamber0)) return kFALSE;
1129
1130   // Create calvector
1131   TString nbname((const char *)ph->GetName());
1132   nbname  += "PerDetectorVector";
1133   if(!fDetSumVector) fDetSumVector = new AliTRDCalibraVector();
1134   else{
1135     fDetSumVector->~AliTRDCalibraVector();
1136     new(fDetSumVector) AliTRDCalibraVector();
1137   }
1138   if(i==1){
1139     fDetSumVector->SetTimeMax(nxbins);
1140   }
1141   if(i==2){
1142     fDetSumVector->SetNumberBinPRF(nxbins);
1143     fDetSumVector->SetPRFRange(TMath::Abs(lowedge));
1144   }
1145   fDetSumVector->SetDetCha0(i,1);
1146   fDetSumVector->SetDetCha2(i,1);
1147   fDetSumVector->SetNzNrphi(i,0,0);
1148   if(i==2) {
1149     Int_t nbg = GetNumberOfGroupsPRF((const char *)name);
1150     fDetSumVector->SetNbGroupPRF(nbg);
1151   }
1152
1153   // Create Histo
1154   TString nname((const char *)ph->GetName());
1155   nname  += "PerDetector";
1156   TString title("Nz");
1157   title += 0;
1158   title += "Nrphi";
1159   title += 0;
1160   if(!fDetSum) fDetSum = new TH2D((const char *)nname,(const Char_t *) title
1161                                 ,nxbins,lowedge,upedge,540,0,540);
1162   else{
1163     fDetSum->~TH2D();
1164     new(fDetSum) TH2D((const Char_t *) nname,(const Char_t *) title
1165                      ,nxbins,lowedge,upedge,540,0,540);
1166   }
1167   fDetSum->SetYTitle("Detector number");
1168   fDetSum->SetXTitle(xaxis->GetTitle());
1169   fDetSum->SetStats(0);
1170
1171   Int_t counter = 0;
1172
1173   for(Int_t det = 0; det < 540; det++){
1174
1175     Int_t numberofgroup = 0;
1176     if(AliTRDgeometry::GetStack(det) == 2) numberofgroup = perChamber2;
1177     else numberofgroup = perChamber0;
1178
1179     for(Int_t nx = 1; nx <= nxbins; nx++) {
1180
1181       Double_t entries = 0.0;
1182       Double_t sumw2 = 0.0;
1183       Double_t sumw = 0.0;
1184
1185       for(Int_t k = counter+1; k <= (counter+numberofgroup); k++){
1186         Int_t  binnumber = ph->GetBin(nx,k);
1187         entries += ph->GetBinEntries(binnumber);
1188         sumw2 += (ph->GetBinError(binnumber)*ph->GetBinError(binnumber)+ph->GetBinContent(binnumber)*ph->GetBinContent(binnumber))*ph->GetBinEntries(binnumber);
1189         sumw += ph->GetBinContent(binnumber)*ph->GetBinEntries(binnumber);
1190       }
1191
1192       Double_t mean = 0.0;
1193       if(entries > 0.0) mean = sumw/entries;
1194       Double_t squaremean = 0.0;
1195       if(entries > 0.0) squaremean = sumw2/entries;
1196       Double_t errorf = squaremean - mean*mean;
1197       Double_t error = 0.0;
1198       if(entries > 0.0) error = TMath::Sqrt(TMath::Abs(errorf)/entries);
1199
1200       fDetSum->SetBinContent(nx,det+1,mean);
1201       fDetSum->SetBinError(nx,det+1,error);
1202
1203       if(i==1) fDetSumVector->FillVectorPH(det,0,nx-1,(Int_t)entries,(Float_t)mean,(Float_t)squaremean);
1204       if(i==2) fDetSumVector->FillVectorPRF(det,0,nx-1,(Int_t)entries,(Float_t)mean,(Float_t)squaremean);
1205
1206     }
1207
1208     counter += numberofgroup;
1209
1210   }
1211
1212   return kTRUE;
1213
1214
1215 }
1216 //_____________________________________________________________________________
1217 Bool_t AliTRDcalibration::SetNrphiFromTObject(const char *name, Int_t i, AliTRDCalibraMode *calibMode) const
1218 {
1219   //
1220   // Set the granularity from object
1221   //
1222
1223   const Char_t *patternrphi0 = "Nrphi0";
1224   const Char_t *patternrphi1 = "Nrphi1";
1225   const Char_t *patternrphi2 = "Nrphi2";
1226   const Char_t *patternrphi3 = "Nrphi3";
1227   const Char_t *patternrphi4 = "Nrphi4";
1228   const Char_t *patternrphi5 = "Nrphi5";
1229   const Char_t *patternrphi6 = "Nrphi6";
1230
1231
1232   const Char_t *patternrphi10 = "Nrphi10";
1233   const Char_t *patternrphi100 = "Nrphi100";
1234   const Char_t *patternz10 = "Nz10";
1235   const Char_t *patternz100 = "Nz100";
1236
1237   // Nrphi mode
1238   if ((strstr(name,patternrphi100)) && (strstr(name,patternz100))) {
1239     calibMode->SetAllTogether(i);
1240     return kTRUE;
1241   }
1242   if ((strstr(name,patternrphi10)) && (strstr(name,patternz10))) {
1243     calibMode->SetPerSuperModule(i);
1244     return kTRUE;
1245   }
1246
1247   if (strstr(name,patternrphi0)) {
1248     calibMode->SetNrphi(i ,0);
1249     return kTRUE;
1250   }
1251   if (strstr(name,patternrphi1)) {
1252     calibMode->SetNrphi(i, 1);
1253     return kTRUE;
1254   }
1255   if (strstr(name,patternrphi2)) {
1256     calibMode->SetNrphi(i, 2);
1257     return kTRUE;
1258   }
1259   if (strstr(name,patternrphi3)) {
1260     calibMode->SetNrphi(i, 3);
1261     return kTRUE;
1262   }
1263   if (strstr(name,patternrphi4)) {
1264     calibMode->SetNrphi(i, 4);
1265     return kTRUE;
1266   }
1267   if (strstr(name,patternrphi5)) {
1268     calibMode->SetNrphi(i, 5);
1269     return kTRUE;
1270   }
1271   if (strstr(name,patternrphi6)) {
1272     calibMode->SetNrphi(i, 6);
1273     return kTRUE;
1274   }
1275
1276   calibMode->SetNrphi(i ,0);
1277   return kFALSE;
1278
1279 }
1280 //_____________________________________________________________________________
1281 Bool_t AliTRDcalibration::SetNzFromTObject(const char *name, Int_t i, AliTRDCalibraMode *calibMode) const
1282 {
1283   //
1284   // Set fNz[i] of the AliTRDCalibraFit::Instance()
1285   // corresponding to the given TObject
1286   //
1287
1288   // Some patterns
1289   const Char_t *patternz0    = "Nz0";
1290   const Char_t *patternz1    = "Nz1";
1291   const Char_t *patternz2    = "Nz2";
1292   const Char_t *patternz3    = "Nz3";
1293   const Char_t *patternz4    = "Nz4";
1294
1295   const Char_t *patternrphi10 = "Nrphi10";
1296   const Char_t *patternrphi100 = "Nrphi100";
1297   const Char_t *patternz10 = "Nz10";
1298   const Char_t *patternz100 = "Nz100";
1299
1300   if ((strstr(name,patternrphi100)) && (strstr(name,patternz100))) {
1301     calibMode->SetAllTogether(i);
1302     return kTRUE;
1303   }
1304   if ((strstr(name,patternrphi10)) && (strstr(name,patternz10))) {
1305     calibMode->SetPerSuperModule(i);
1306     return kTRUE;
1307   }
1308   if (strstr(name,patternz0)) {
1309     calibMode->SetNz(i, 0);
1310     return kTRUE;
1311   }
1312   if (strstr(name,patternz1)) {
1313     calibMode->SetNz(i ,1);
1314     return kTRUE;
1315   }
1316   if (strstr(name,patternz2)) {
1317     calibMode->SetNz(i ,2);
1318     return kTRUE;
1319   }
1320   if (strstr(name,patternz3)) {
1321     calibMode->SetNz(i ,3);
1322     return kTRUE;
1323   }
1324   if (strstr(name,patternz4)) {
1325     calibMode->SetNz(i ,4);
1326     return kTRUE;
1327   }
1328
1329   calibMode->SetNz(i ,0);
1330   return kFALSE;
1331 }
1332 //____________________________________________________________________________________________________
1333 Int_t AliTRDcalibration::GetNumberOfGroupsPRF(const char* nametitle) const
1334 {
1335   //
1336   // Get numberofgroupsprf
1337   //
1338
1339   // Some patterns
1340   const Char_t *pattern0 = "Ngp0";
1341   const Char_t *pattern1 = "Ngp1";
1342   const Char_t *pattern2 = "Ngp2";
1343   const Char_t *pattern3 = "Ngp3";
1344   const Char_t *pattern4 = "Ngp4";
1345   const Char_t *pattern5 = "Ngp5";
1346   const Char_t *pattern6 = "Ngp6";
1347
1348   // Nrphi mode
1349   if (strstr(nametitle,pattern0)) {
1350     return 0;
1351   }
1352   if (strstr(nametitle,pattern1)) {
1353     return 1;
1354   }
1355   if (strstr(nametitle,pattern2)) {
1356     return 2;
1357   }
1358   if (strstr(nametitle,pattern3)) {
1359     return 3;
1360   }
1361   if (strstr(nametitle,pattern4)) {
1362     return 4;
1363   }
1364   if (strstr(nametitle,pattern5)) {
1365     return 5;
1366   }
1367   if (strstr(nametitle,pattern6)){
1368     return 6;
1369   }
1370   else return -1;
1371
1372
1373 }
1374