Replace AliTRDCalibra
[u/mrichter/AliRoot.git] / 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
26 // calibration (see AliTRDCalibraMode). 
27 // 2D Histograms (Histo2d) or vectors (Vector2d), then converted in Trees, will be filled
28 // from RAW DATA in a run or from reconstructed TRD tracks during the offline tracking 
29 // in the function "FollowBackProlongation" (AliTRDtracker)
30 // Per default the functions to fill are off.                                   
31 //                        
32 // Author:
33 //   R. Bailhache (R.Bailhache@gsi.de)
34 //                            
35 //////////////////////////////////////////////////////////////////////////////////////
36
37 #include <TTree.h>
38 #include <TProfile2D.h>
39 #include <TProfile.h>
40 #include <TFile.h>
41 #include <TChain.h>
42 #include <TStyle.h>
43 #include <TCanvas.h>
44 #include <TGraphErrors.h>
45 #include <TObjArray.h>
46 #include <TH1F.h>
47 #include <TH2I.h>
48 #include <TH2.h>
49 #include <TStopwatch.h>
50 #include <TMath.h>
51 #include <TDirectory.h>
52 #include <TROOT.h>
53
54 #include "AliLog.h"
55 #include "AliCDBManager.h"
56
57 #include "AliTRDCalibraFillHisto.h"
58 #include "AliTRDCalibraMode.h"
59 #include "AliTRDCalibraVector.h"
60 #include "AliTRDcalibDB.h"
61 #include "AliTRDCommonParam.h"
62 #include "AliTRDmcmTracklet.h"
63 #include "AliTRDpadPlane.h"
64 #include "AliTRDcluster.h"
65 #include "AliTRDtrack.h"
66
67
68 ClassImp(AliTRDCalibraFillHisto)
69
70 AliTRDCalibraFillHisto* AliTRDCalibraFillHisto::fgInstance = 0;
71 Bool_t AliTRDCalibraFillHisto::fgTerminated = kFALSE;
72
73 //_____________singleton implementation_________________________________________________
74 AliTRDCalibraFillHisto *AliTRDCalibraFillHisto::Instance()
75 {
76   //
77   // Singleton implementation
78   //
79
80   if (fgTerminated != kFALSE) {
81     return 0;
82   }
83
84   if (fgInstance == 0) {
85     fgInstance = new AliTRDCalibraFillHisto();
86   }
87
88   return fgInstance;
89
90 }
91
92 //______________________________________________________________________________________
93 void AliTRDCalibraFillHisto::Terminate()
94 {
95   //
96   // Singleton implementation
97   // Deletes the instance of this class
98   //
99
100   fgTerminated = kTRUE;
101
102   if (fgInstance != 0) {
103     delete fgInstance;
104     fgInstance = 0;
105   }
106
107 }
108
109 //______________________________________________________________________________________
110 AliTRDCalibraFillHisto::AliTRDCalibraFillHisto()
111   :TObject()
112   ,fMITracking(kFALSE)
113   ,fMcmTracking(kFALSE)
114   ,fMcmCorrectAngle(kFALSE)
115   ,fCH2dOn(kFALSE)
116   ,fPH2dOn(kFALSE)
117   ,fPRF2dOn(kFALSE)
118   ,fHisto2d(kFALSE)
119   ,fVector2d(kFALSE)
120   ,fRelativeScale(0)
121   ,fCountRelativeScale(0)
122   ,fRelativeScaleAuto(kFALSE)
123   ,fThresholdClusterPRF1(0.0)
124   ,fThresholdClusterPRF2(0.0)
125   ,fCenterOfflineCluster(kFALSE)
126   ,fWriteName(0)
127   ,fDetectorAliTRDtrack(kFALSE)
128   ,fChamberAliTRDtrack(-1)
129   ,fDetectorPreviousTrack(-1)
130   ,fGoodTrack(kTRUE)
131   ,fAmpTotal(0x0)
132   ,fPHPlace(0x0)
133   ,fPHValue(0x0)
134   ,fNumberClusters(0)
135   ,fProcent(0.0)
136   ,fDifference(0)
137   ,fNumberTrack(0)
138   ,fTimeMax(0)
139   ,fSf(0.0)
140   ,fNumberBinCharge(0)
141   ,fNumberBinPRF(0)
142   ,fCalibraVector(0)
143   ,fPH2d(0x0)
144   ,fPRF2d(0x0)
145   ,fCH2d(0x0)
146 {
147   //
148   // Default constructor
149   //
150
151   fCalibraMode = new AliTRDCalibraMode();
152
153   // Write
154   for (Int_t i = 0; i < 3; i++) {
155     fWrite[i]     = kFALSE;
156   }
157
158   // Init
159   Init();
160   
161 }
162
163 //______________________________________________________________________________________
164 AliTRDCalibraFillHisto::AliTRDCalibraFillHisto(const AliTRDCalibraFillHisto &c)
165   :TObject(c)
166   ,fMITracking(kFALSE)
167   ,fMcmTracking(kFALSE)
168   ,fMcmCorrectAngle(kFALSE)
169   ,fCH2dOn(kFALSE)
170   ,fPH2dOn(kFALSE)
171   ,fPRF2dOn(kFALSE)
172   ,fHisto2d(kFALSE)
173   ,fVector2d(kFALSE)
174   ,fRelativeScale(0)
175   ,fCountRelativeScale(0)
176   ,fRelativeScaleAuto(kFALSE)
177   ,fThresholdClusterPRF1(0.0)
178   ,fThresholdClusterPRF2(0.0)
179   ,fCenterOfflineCluster(kFALSE)
180   ,fWriteName(0)
181   ,fDetectorAliTRDtrack(kFALSE)
182   ,fChamberAliTRDtrack(-1)
183   ,fDetectorPreviousTrack(-1)
184   ,fGoodTrack(kTRUE)
185   ,fAmpTotal(0x0)
186   ,fPHPlace(0x0)
187   ,fPHValue(0x0)
188   ,fNumberClusters(0)
189   ,fProcent(0.0)
190   ,fDifference(0)
191   ,fNumberTrack(0)
192   ,fTimeMax(0)
193   ,fSf(0.0)
194   ,fNumberBinCharge(0)
195   ,fNumberBinPRF(0)
196   ,fCalibraVector(0)
197   ,fPH2d(0x0)
198   ,fPRF2d(0x0)
199   ,fCH2d(0x0) 
200 {
201   //
202   // Copy constructor
203   //
204
205 }
206
207 //____________________________________________________________________________________
208 AliTRDCalibraFillHisto::~AliTRDCalibraFillHisto()
209 {
210   //
211   // AliTRDCalibraFillHisto destructor
212   //
213
214   ClearHistos();
215   
216 }
217
218 //_____________________________________________________________________________
219 void AliTRDCalibraFillHisto::Destroy() 
220 {
221   //
222   // Delete instance 
223   //
224
225   if (fgInstance) {
226     delete fgInstance;
227     fgInstance = 0x0;
228   }
229
230 }
231
232 //_____________________________________________________________________________
233 void AliTRDCalibraFillHisto::ClearHistos() 
234 {
235   //
236   // Delete the histos
237   //
238
239   if (fPH2d) {
240     delete fPH2d;
241     fPH2d  = 0x0;
242   }
243   if (fCH2d) {
244     delete fCH2d;
245     fCH2d  = 0x0;
246   }
247   if (fPRF2d) {
248     delete fPRF2d;
249     fPRF2d = 0x0;
250   }
251
252 }
253
254 //_____________________________________________________________________________
255 void AliTRDCalibraFillHisto::Init() 
256 {
257   //
258   // Init some default values
259   //
260
261   // How to fill the 2D
262   fThresholdClusterPRF1 = 2.0;
263   fThresholdClusterPRF2 = 15.0;
264   
265   // Store the Info
266   fNumberBinCharge      = 100;
267   fNumberBinPRF         = 40;
268   
269   // Write
270   fWriteName            = "TRD.calibration.root";
271   
272   // Internal variables
273   
274   // Fill the 2D histos in the offline tracking
275   fDetectorPreviousTrack = -1;
276   fChamberAliTRDtrack    = -1;
277   fGoodTrack             = kTRUE;
278
279   fProcent               = 6.0;
280   fDifference            = 17;
281   fNumberClusters        = 18;
282   fNumberTrack           = 0;
283   fNumberUsedCh[0]       = 0;
284   fNumberUsedCh[1]       = 0;
285   fNumberUsedPh[0]       = 0;
286   fNumberUsedPh[1]       = 0;
287  
288 }
289
290 //____________Functions for initialising the AliTRDCalibraFillHisto in the code_________
291 Bool_t AliTRDCalibraFillHisto::Init2Dhistos()
292 {
293   //
294   // For the offline tracking
295   // This function will be called in the function AliReconstruction::Run() 
296   // Init the calibration mode (Nz, Nrphi), the 2D histograms if fHisto2d = kTRUE, 
297   //
298
299   // DB Setting
300   // Get cal
301   AliTRDcalibDB *cal = AliTRDcalibDB::Instance();
302   if (!cal) {
303     AliInfo("Could not get calibDB");
304     return kFALSE;
305   }
306   AliTRDCommonParam *parCom = AliTRDCommonParam::Instance();
307   if (!parCom) {
308     AliInfo("Could not get CommonParam");
309     return kFALSE;
310   }
311
312   // Some parameters
313   fTimeMax = cal->GetNumberOfTimeBins();
314   fSf      = parCom->GetSamplingFrequency();
315   if (fRelativeScaleAuto) {
316     fRelativeScale = 0;
317   }
318   else {
319     fRelativeScale = 20;
320   }
321
322   //If vector method On initialised all the stuff
323   if(fVector2d){
324     fCalibraVector = new AliTRDCalibraVector();
325   }
326
327
328   // Create the 2D histos corresponding to the pad groupCalibration mode
329   if (fCH2dOn) {
330
331     AliInfo(Form("The pad calibration mode for the relative gain calibration: Nz %d, and Nrphi %d"
332                 ,fCalibraMode->GetNz(0)
333                 ,fCalibraMode->GetNrphi(0)));
334     
335     // Calcul the number of Xbins
336     Int_t Ntotal0 = 0;
337     fCalibraMode->ModePadCalibration(2,0);
338     fCalibraMode->ModePadFragmentation(0,2,0,0);
339     fCalibraMode->SetDetChamb2(0);
340     Ntotal0 += 6 * 18 * fCalibraMode->GetDetChamb2(0);
341     fCalibraMode->ModePadCalibration(0,0);
342     fCalibraMode->ModePadFragmentation(0,0,0,0);
343     fCalibraMode->SetDetChamb0(0);
344     Ntotal0 += 6 * 4 * 18 * fCalibraMode->GetDetChamb0(0);
345     AliInfo(Form("Total number of Xbins: %d",Ntotal0));
346
347     // Create the 2D histo
348     if (fHisto2d) {
349       CreateCH2d(Ntotal0);
350     }
351     if (fVector2d) {
352       fCalibraVector->SetNumberBinCharge(fNumberBinCharge);
353     }
354
355     // Variable
356     fAmpTotal = new Float_t[TMath::Max(fCalibraMode->GetDetChamb2(0),fCalibraMode->GetDetChamb0(0))];
357     for (Int_t k = 0; k < TMath::Max(fCalibraMode->GetDetChamb2(0),fCalibraMode->GetDetChamb0(0)); k++) {
358       fAmpTotal[k] = 0.0;
359     } 
360
361   }
362
363   if (fPH2dOn) {
364
365     AliInfo(Form("The pad calibration mode for the drift velocity calibration: Nz %d, and Nrphi %d"
366                 ,fCalibraMode->GetNz(1)
367                 ,fCalibraMode->GetNrphi(1)));
368     
369     // Calcul the number of Xbins
370     Int_t Ntotal1 = 0;
371     fCalibraMode->ModePadCalibration(2,1);
372     fCalibraMode->ModePadFragmentation(0,2,0,1);
373     fCalibraMode->SetDetChamb2(1);
374     Ntotal1 += 6 * 18 * fCalibraMode->GetDetChamb2(1);
375     fCalibraMode->ModePadCalibration(0,1);
376     fCalibraMode->ModePadFragmentation(0,0,0,1);
377     fCalibraMode->SetDetChamb0(1);
378     Ntotal1 += 6 * 4 * 18 * fCalibraMode->GetDetChamb0(1);
379     AliInfo(Form("Total number of Xbins: %d",Ntotal1));
380
381     // Create the 2D histo
382     if (fHisto2d) {
383       CreatePH2d(Ntotal1);
384     }
385     if (fVector2d) {
386       fCalibraVector->SetTimeMax(fTimeMax);
387     }
388    
389     // Variable
390     fPHPlace = new Short_t[fTimeMax];
391     for (Int_t k = 0; k < fTimeMax; k++) {
392       fPHPlace[k] = -1;
393     } 
394     fPHValue = new Float_t[fTimeMax];
395     for (Int_t k = 0; k < fTimeMax; k++) {
396       fPHValue[k] = 0.0;
397     }
398
399   }
400
401   if (fPRF2dOn) {
402
403     AliInfo(Form("The pad calibration mode for the PRF calibration: Nz %d, and Nrphi %d"
404                 ,fCalibraMode->GetNz(2)
405                 ,fCalibraMode->GetNrphi(2)));
406     
407     // Calcul the number of Xbins
408     Int_t Ntotal2 = 0;
409     fCalibraMode->ModePadCalibration(2,2);
410     fCalibraMode->ModePadFragmentation(0,2,0,2);
411     fCalibraMode->SetDetChamb2(2);
412     Ntotal2 += 6 * 18 * fCalibraMode->GetDetChamb2(2);
413     fCalibraMode->ModePadCalibration(0,2);
414     fCalibraMode->ModePadFragmentation(0,0,0,2);
415     fCalibraMode->SetDetChamb0(2);
416     Ntotal2 += 6 * 4 * 18 * fCalibraMode->GetDetChamb0(2);
417     AliInfo(Form("Total number of Xbins: %d",Ntotal2));
418
419     // Create the 2D histo
420     if (fHisto2d) {
421       CreatePRF2d(Ntotal2);
422     }
423     if (fVector2d) {
424       fCalibraVector->SetNumberBinPRF(fNumberBinPRF);
425     }
426   
427   }
428
429   return kTRUE;
430
431 }
432
433 //____________Functions for filling the histos in the code_____________________
434
435 //____________Offine tracking in the AliTRDtracker_____________________________
436 Bool_t AliTRDCalibraFillHisto::ResetTrack()
437 {
438   //
439   // For the offline tracking
440   // This function will be called in the function
441   // AliTRDtracker::FollowBackPropagation() at the beginning 
442   // Reset the parameter to know we have a new TRD track
443   //
444   
445   fDetectorAliTRDtrack = kFALSE;
446   return kTRUE;
447
448 }
449
450 //____________Offline tracking in the AliTRDtracker____________________________
451 Bool_t AliTRDCalibraFillHisto::UpdateHistograms(AliTRDcluster *cl, AliTRDtrack *t)
452 {
453   //
454   // For the offline tracking
455   // This function will be called in the function
456   // AliTRDtracker::FollowBackPropagation() in the loop over the clusters
457   // of TRD tracks 
458   // Fill the 2D histos or the vectors with the info of the clusters at
459   // the end of a detectors if the track is "good"
460   //
461
462   // Get the parameter object
463   AliTRDCommonParam *parCom = AliTRDCommonParam::Instance();
464   if (!parCom) {
465     AliInfo("Could not get CommonParam");
466     return kFALSE;
467   }
468
469   // Get the parameter object
470   AliTRDcalibDB     *cal    = AliTRDcalibDB::Instance();
471   if (!cal) {
472     AliInfo("Could not get calibDB");
473     return kFALSE;
474   }
475  
476   // Localisation of the detector
477   Int_t detector = cl->GetDetector();
478   Int_t chamber  = GetChamber(detector);
479   Int_t plane    = GetPlane(detector);
480
481   // Fill the infos for the previous clusters if not the same
482   // detector anymore or if not the same track
483   if (((detector != fDetectorPreviousTrack) || (!fDetectorAliTRDtrack)) && 
484       (fDetectorPreviousTrack != -1)) {
485     
486     fNumberTrack++;   
487     
488     // If the same track, then look if the previous detector is in
489     // the same plane, if yes: not a good track
490     if (fDetectorAliTRDtrack && 
491         (GetPlane(detector) <= GetPlane(fDetectorPreviousTrack))) {
492       fGoodTrack = kFALSE;
493     }
494
495     // Fill only if the track doesn't touch a masked pad or doesn't
496     // appear in the middle (fGoodTrack)
497     if (fGoodTrack) {
498
499       // Gain calibration
500       if (fCH2dOn) {
501         FillTheInfoOfTheTrackCH();
502       }
503       
504       // PH calibration
505       if (fPH2dOn) {
506         FillTheInfoOfTheTrackPH();    
507       }
508       
509     } // if a good track
510     
511     ResetfVariables();
512     
513   } // Fill at the end the charge
514   
515   // Calcul the position of the detector
516   if (detector != fDetectorPreviousTrack) {
517     LocalisationDetectorXbins(detector);
518   }
519
520   // Reset the good track for the PRF
521   Bool_t good = kTRUE;
522   
523   // Localisation of the cluster
524   Double_t pos[3] = { 0.0, 0.0, 0.0 };
525   pos[0] = cl->GetX();
526   pos[1] = cl->GetY();
527   pos[2] = cl->GetZ();
528   Int_t    time   = cl->GetLocalTimeBin();
529   
530   // Reset the detector
531   fDetectorPreviousTrack = detector;
532   fDetectorAliTRDtrack   = kTRUE;
533   
534   // Position of the cluster
535   AliTRDpadPlane *padplane = parCom->GetPadPlane(plane,chamber);
536   Int_t    row        = padplane->GetPadRowNumber(pos[2]);
537   Double_t offsetz    = padplane->GetPadRowOffset(row,pos[2]);
538   Double_t offsettilt = padplane->GetTiltOffset(offsetz);
539   //Int_t    col        = padplane->GetPadColNumber(pos[1] + offsettilt,offsetz);
540   Int_t    col        = padplane->GetPadColNumber(pos[1]+offsettilt);
541   
542   // See if we are not near a masked pad
543   if (!IsPadOn(detector,col,row)) {
544     good       = kFALSE;
545     fGoodTrack = kFALSE;
546   }
547
548   if (col > 0) {
549     if (!IsPadOn(detector,col-1,row)) {
550       fGoodTrack = kFALSE;
551       good       = kFALSE;
552     }
553   }
554
555   if (col < 143) {
556     if (!IsPadOn(detector,col+1,row)) {
557       fGoodTrack = kFALSE;
558       good       = kFALSE;
559     }
560   }
561
562   // Row of the cluster and position in the pad groups
563   Int_t posr[3] = { 0, 0, 0 };
564   if ((fCH2dOn)  && (fCalibraMode->GetNnZ(0) != 0)) {
565     posr[0] = (Int_t) row / fCalibraMode->GetNnZ(0);
566   }
567   if ((fPH2dOn)  && (fCalibraMode->GetNnZ(1) != 0)) {
568     posr[1] = (Int_t) row / fCalibraMode->GetNnZ(1);
569   }
570   if ((fPRF2dOn) && (fCalibraMode->GetNnZ(2) != 0)) {
571     posr[2] = (Int_t) row / fCalibraMode->GetNnZ(2);
572   }  
573       
574   // Col of the cluster and position in the pad groups
575   Int_t posc[3] = { 0, 0, 0 };
576   if ((fCH2dOn)  && (fCalibraMode->GetNnRphi(0) != 0)) {
577     posc[0] = (Int_t) col / fCalibraMode->GetNnRphi(0);
578   }
579   if ((fPH2dOn)  && (fCalibraMode->GetNnRphi(1) != 0)) {
580     posc[1] = (Int_t) col / fCalibraMode->GetNnRphi(1);
581   }
582   if ((fPRF2dOn) && (fCalibraMode->GetNnRphi(2) != 0)) {
583     posc[2] = (Int_t) col / fCalibraMode->GetNnRphi(2);
584   }
585
586   // Charge in the cluster
587   // For the moment take the abs
588   Float_t  q        = TMath::Abs(cl->GetQ());
589   Short_t  *signals = cl->GetSignals();
590  
591   // Correction due to the track angle
592   Float_t correction    = 1.0;
593   Float_t normalisation = 6.67;
594   if ((q >0) && (t->GetNdedx() > 0)) {
595     correction = t->GetClusterdQdl((t->GetNdedx() - 1)) / (normalisation);
596   }
597
598   // Fill the fAmpTotal with the charge
599   if (fCH2dOn) {
600     fAmpTotal[(Int_t) (posc[0]*fCalibraMode->GetNfragZ(0)+posr[0])] += correction;
601   }
602
603   // Fill the fPHPlace and value
604   if (fPH2dOn) {
605     fPHPlace[time] = posc[1]*fCalibraMode->GetNfragZ(1)+posr[1];
606     fPHValue[time] = correction;
607   }
608
609   // Fill direct the PRF
610   if ((fPRF2dOn) && (good)) {
611
612     Float_t yminus  = 0.0;
613     Float_t xcenter = 0.0;
614     Float_t ycenter = 0.0;
615     Float_t ymax    = 0.0;
616     Bool_t  echec   = kFALSE;
617     
618     if ((cl->From3pad()) && (!cl->IsUsed())) { 
619          
620       // Center 3 balanced and cut on the cluster shape
621       if ((((Float_t) signals[3]) > fThresholdClusterPRF2) && 
622           (((Float_t) signals[2]) > fThresholdClusterPRF2) && 
623           (((Float_t) signals[4]) > fThresholdClusterPRF2) && 
624           (((Float_t) signals[1]) < fThresholdClusterPRF1) && 
625           (((Float_t) signals[5]) < fThresholdClusterPRF1) && 
626           ((((Float_t) signals[2])*((Float_t) signals[4])/(((Float_t) signals[3])*((Float_t) signals[3]))) < 0.06)) {
627
628
629         //First calculate the position of the cluster and the y
630         //echec enables to repair cases where it fails
631         //
632         // Col correspond to signals[3]
633         if (fCenterOfflineCluster) {
634           xcenter = cl->GetCenter();
635         }
636         else {
637           // Security if the denomiateur is 0 
638           if ((((Float_t) (((Float_t) signals[3]) * ((Float_t) signals[3]))) / 
639                            ((Float_t) (((Float_t) signals[2]) * ((Float_t) signals[4])))) != 1.0) {
640             xcenter = 0.5 * (TMath::Log((Float_t) (((Float_t) signals[4]) / ((Float_t) signals[2]))))
641                           / (TMath::Log(((Float_t) (((Float_t) signals[3]) * ((Float_t) signals[3]))) 
642                                       / ((Float_t) (((Float_t) signals[2]) * ((Float_t) signals[4])))));
643           }
644           else {
645             echec = kTRUE;
646           }
647         }
648         //after having calculating the position calculate the y
649         if (TMath::Abs(xcenter) < 0.5) {
650           ycenter = (Float_t) (((Float_t) signals[3]) 
651                             / (((Float_t) signals[2]) + ((Float_t) signals[3]) + (((Float_t) signals[4]))));
652           yminus  = (Float_t) (((Float_t) signals[2]) 
653                             / (((Float_t) signals[2]) + ((Float_t) signals[3]) + (((Float_t) signals[4]))));
654           ymax    = (Float_t) (((Float_t) signals[4]) 
655                             / (((Float_t) signals[2]) + ((Float_t) signals[3]) + (((Float_t) signals[4]))));
656           //If the charge of the cluster is too far away from the corrected one cut
657           if ((TMath::Abs(((Float_t) signals[2]) + ((Float_t) signals[3]) + (((Float_t) signals[4])) - q) > 10.0)) {
658             echec = kTRUE;
659           }
660         }
661         else {
662           echec = kTRUE;
663         }
664       
665         //Then Fill the histo if no echec
666         //
667         // Fill only if it is in the drift region!
668         if ((((Float_t) (((Float_t) time) / fSf)) > 0.3) && (!echec)) {
669           if (fHisto2d) {
670             fPRF2d->Fill((fCalibraMode->GetXbins(2)+posc[2]*fCalibraMode->GetNfragZ(2)+posr[2]+0.5),xcenter,ycenter);
671             fPRF2d->Fill((fCalibraMode->GetXbins(2)+posc[2]*fCalibraMode->GetNfragZ(2)+posr[2]+0.5),-(xcenter+1.0),yminus);
672             fPRF2d->Fill((fCalibraMode->GetXbins(2)+posc[2]*fCalibraMode->GetNfragZ(2)+posr[2]+0.5),1.0-xcenter,ymax);
673           }
674           if (fVector2d) {
675             fCalibraVector->UpdateVectorPRF(fCalibraMode->GetXbins(2)+posc[2]*fCalibraMode->GetNfragZ(2)+posr[2],xcenter,ycenter);
676             fCalibraVector->UpdateVectorPRF(fCalibraMode->GetXbins(2)+posc[2]*fCalibraMode->GetNfragZ(2)+posr[2],-(xcenter+1.0),yminus);
677             fCalibraVector->UpdateVectorPRF(fCalibraMode->GetXbins(2)+posc[2]*fCalibraMode->GetNfragZ(2)+posr[2],1.0-xcenter,ymax);
678           }
679         } // If in the drift region
680       } // center 3 balanced and cut on the cluster shape
681     } // Cluster isole
682   } // PRF2dOn  
683   
684   return kTRUE;
685   
686 }
687
688 //____________Online trackling in AliTRDtrigger________________________________
689 Bool_t AliTRDCalibraFillHisto::UpdateHistogramcm(AliTRDmcmTracklet *trk)
690 {
691   //
692   // For the tracking
693   // This function will be called in the function AliTRDtrigger::TestTracklet
694   // before applying the pt cut on the tracklets 
695   // Fill the infos for the tracklets fTrkTest if the tracklets is "good"
696   //
697   
698   // Localisation of the Xbins involved
699   Int_t idect = trk->GetDetector();
700   LocalisationDetectorXbins(idect);
701
702   // Get the parameter object
703   AliTRDcalibDB *cal = AliTRDcalibDB::Instance();
704   if (!cal) {
705     AliInfo("Could not get calibDB");
706     return kFALSE;
707   }
708    
709   // Reset
710   ResetfVariables();
711  
712   // Row of the tracklet and position in the pad groups
713   Int_t row     = trk->GetRow();
714   Int_t posr[3] = { 0, 0, 0 };
715   if ((fCH2dOn)  && (fCalibraMode->GetNnZ(0) != 0)) {
716     posr[0] = (Int_t) row / fCalibraMode->GetNnZ(0);
717   }
718   if ((fPH2dOn)  && (fCalibraMode->GetNnZ(1) != 0)) {
719     posr[1] = (Int_t) row / fCalibraMode->GetNnZ(1);
720   }
721   if ((fPRF2dOn) && (fCalibraMode->GetNnZ(2) != 0)) {
722     posr[2] = (Int_t) row / fCalibraMode->GetNnZ(2);
723   }
724  
725   // Eventuelle correction due to track angle in z direction
726   Float_t correction = 1.0;
727   if (fMcmCorrectAngle) {
728     Float_t z = trk->GetRowz();
729     Float_t r = trk->GetTime0();
730     correction = r / TMath::Sqrt((r*r+z*z));
731   }
732
733   // Boucle sur les clusters
734   // Condition on number of cluster: don't come from the middle of the detector
735   if (trk->GetNclusters() >= fNumberClusters) {
736
737     for (Int_t icl = 0; icl < trk->GetNclusters(); icl++) {
738
739       Float_t amp[3] = { 0.0, 0.0, 0.0 };
740       Int_t   time   = trk->GetClusterTime(icl);
741       Int_t   col    = trk->GetClusterCol(icl);
742             
743       amp[0] = trk->GetClusterADC(icl)[0] * correction;
744       amp[1] = trk->GetClusterADC(icl)[1] * correction;
745       amp[2] = trk->GetClusterADC(icl)[2] * correction;
746
747       
748       if ((amp[0] < 0.0) || 
749           (amp[1] < 0.0) || 
750           (amp[2] < 0.0)) {
751         continue;
752       }
753
754       // Col of cluster and position in the pad groups
755       Int_t posc[3] = { 0, 0, 0 };
756       if ((fCH2dOn)  && (fCalibraMode->GetNnRphi(0) != 0)) {
757         posc[0] = (Int_t) col / fCalibraMode->GetNnRphi(0);
758       }
759       if ((fPH2dOn)  && (fCalibraMode->GetNnRphi(1) != 0)) {
760         posc[1] = (Int_t) col / fCalibraMode->GetNnRphi(1);
761       }
762       if ((fPRF2dOn) && (fCalibraMode->GetNnRphi(2) != 0)) {
763         posc[2] = (Int_t) col / fCalibraMode->GetNnRphi(2);
764       }
765
766       // See if we are not near a masked pad
767       Bool_t good = kTRUE;
768       if (!IsPadOn(idect,col,row)) {
769         fGoodTrack = kFALSE;
770         good       = kFALSE;
771       }
772
773       if (col >   0) {
774         if (!IsPadOn(idect,col-1,row)) {
775           fGoodTrack = kFALSE;
776           good       = kFALSE;
777         }
778       }
779       
780       if (col < 143) {
781         if (!IsPadOn(idect,col+1,row)) {
782           fGoodTrack = kFALSE;
783           good       = kFALSE;
784         }
785       }
786
787       // Total spectrum
788       if (fPH2dOn) {
789         fPHPlace[time] = posc[1] * fCalibraMode->GetNfragZ(1) + posr[1];
790       }
791
792       if (fCH2dOn) {
793         fAmpTotal[(Int_t) (posc[0]*fCalibraMode->GetNfragZ(0)+posr[0])] += (Float_t) (amp[0]+amp[1]+amp[2]);
794       }
795       if (fPH2dOn) {
796         fPHValue[time] = (Float_t) (amp[0]+amp[1]+amp[2]);
797       }
798       
799             
800       // Fill PRF direct
801       if (fPRF2dOn && good) {
802         
803         if ((amp[0] > fThresholdClusterPRF2) && 
804             (amp[1] > fThresholdClusterPRF2) && 
805             (amp[2] > fThresholdClusterPRF2) && 
806             ((amp[0]*amp[2]/(amp[1]*amp[1])) < 0.06)) {
807         
808           // Security of the denomiateur is 0
809           if ((((Float_t) (((Float_t) amp[1]) * ((Float_t) amp[1]))) 
810              / ((Float_t) (((Float_t) amp[0]) * ((Float_t) amp[2])))) != 1.0) {
811             Float_t xcenter = 0.5 * (TMath::Log(amp[2] / amp[0]))
812                                   / (TMath::Log((amp[1]*amp[1]) / (amp[0]*amp[2])));
813             Float_t ycenter = amp[1] / (amp[0] + amp[1] + amp[2]);
814
815             if (TMath::Abs(xcenter) < 0.5) {
816               Float_t yminus = amp[0] / (amp[0]+amp[1]+amp[2]);
817               Float_t ymax   = amp[2] / (amp[0]+amp[1]+amp[2]);
818               // Fill only if it is in the drift region!
819               if (((Float_t) time / fSf) > 0.3) {
820                 if (fHisto2d) {
821                   fPRF2d->Fill((fCalibraMode->GetXbins(2)+posc[2]*fCalibraMode->GetNfragZ(2)+posr[2]+0.5),xcenter,ycenter);
822                   fPRF2d->Fill((fCalibraMode->GetXbins(2)+posc[2]*fCalibraMode->GetNfragZ(2)+posr[2]+0.5),-(xcenter+1.0),yminus);
823                   fPRF2d->Fill((fCalibraMode->GetXbins(2)+posc[2]*fCalibraMode->GetNfragZ(2)+posr[2]+0.5),(1.0-xcenter),ymax);
824                 }
825                 if (fVector2d) {
826                   fCalibraVector->UpdateVectorPRF((fCalibraMode->GetXbins(2)+posc[2]*fCalibraMode->GetNfragZ(2)+posr[2]),xcenter,ycenter);
827                   fCalibraVector->UpdateVectorPRF(fCalibraMode->GetXbins(2)+posc[2]*fCalibraMode->GetNfragZ(2)+posr[2],-(xcenter+1.0),yminus);
828                   fCalibraVector->UpdateVectorPRF(fCalibraMode->GetXbins(2)+posc[2]*fCalibraMode->GetNfragZ(2)+posr[2],(1.0-xcenter),ymax);
829                 }
830               }//in the drift region 
831             }//in the middle
832           }//denominateur security
833         }//cluster shape and thresholds
834       }//good and PRF On
835       
836     } // Boucle clusters
837     
838     // Fill the charge
839     if (fCH2dOn && fGoodTrack) {
840       FillTheInfoOfTheTrackCH();
841     }
842
843     // PH calibration
844     if (fPH2dOn && fGoodTrack) {
845       FillTheInfoOfTheTrackPH();        
846     }
847
848     fNumberTrack++;
849         
850   } // Condition on number of clusters
851
852   return kTRUE;
853   
854 }
855
856 //_____________________________________________________________________________
857 Bool_t AliTRDCalibraFillHisto::IsPadOn(Int_t detector, Int_t col, Int_t row) const
858 {
859   //
860   // Look in the choosen database if the pad is On.
861   // If no the track will be "not good"
862   //
863
864   // Get the parameter object
865   AliTRDcalibDB *cal = AliTRDcalibDB::Instance();
866   if (!cal) {
867     AliInfo("Could not get calibDB");
868     return kFALSE;
869   }
870   
871   if (!cal->IsChamberInstalled(detector)     || 
872        cal->IsChamberMasked(detector)        ||
873        cal->IsPadMasked(detector,col,row)) {
874     return kFALSE;
875   }
876   else {
877     return kTRUE;
878   }
879   
880 }
881
882 //____________Functions for plotting the 2D____________________________________
883
884 //_____________________________________________________________________________
885 void AliTRDCalibraFillHisto::Plot2d()
886 {
887   //
888   // Plot the 2D histos 
889   //
890  
891   if (fPH2dOn) {
892     TCanvas *cph2d = new TCanvas("cph2d","",50,50,600,800);
893     cph2d->cd();
894     fPH2d->Draw("LEGO");
895   }
896   if (fCH2dOn) {
897     TCanvas *cch2d = new TCanvas("cch2d","",50,50,600,800);
898     cch2d->cd();
899     fCH2d->Draw("LEGO");
900   }
901   if (fPRF2dOn) {
902     TCanvas *cPRF2d = new TCanvas("cPRF2d","",50,50,600,800);
903     cPRF2d->cd();
904     fPRF2d->Draw("LEGO");
905   }
906
907 }
908
909 //____________Writing the 2D___________________________________________________
910
911 //_____________________________________________________________________________
912 Bool_t AliTRDCalibraFillHisto::Write2d()
913 {
914   //
915   // Write the 2D histograms or the vectors converted in trees in the file
916   // "TRD.calibration.root" 
917   //
918   
919   TFile *fout = TFile::Open(fWriteName,"RECREATE");
920   // Check if the file could be opened
921   if (!fout || !fout->IsOpen()) {
922     AliInfo("No File found!");
923     return kFALSE;
924   }
925   AliInfo(Form("Numbertrack: %d Numberusedch[0]: %d, Numberusedch[1]: %d Numberusedph[0]: %d, Numberusedph[1]: %d"
926               ,fNumberTrack
927               ,fNumberUsedCh[0]
928               ,fNumberUsedCh[1]
929               ,fNumberUsedPh[0]
930               ,fNumberUsedPh[1]));
931   
932   TStopwatch stopwatch;
933   stopwatch.Start();
934   AliInfo("Write2d");
935
936   if ((fCH2dOn ) && (fWrite[0])) {
937     if (fHisto2d) {
938       fout->WriteTObject(fCH2d);
939     }
940     if (fVector2d) {
941       TString name("Nz");
942       name += fCalibraMode->GetNz(0);
943       name += "Nrphi";
944       name += fCalibraMode->GetNrphi(0);
945       TTree *treeCH2d = fCalibraVector->ConvertVectorCTTreeHisto(fCalibraVector->GetVectorCH(),fCalibraVector->GetPlaCH(),"treeCH2d",(const char *) name);
946       fout->WriteTObject(treeCH2d);
947     }
948   }
949   if ((fPH2dOn ) && (fWrite[1])) {
950     if (fHisto2d) {
951       fout->WriteTObject(fPH2d);
952     }
953     if (fVector2d) {
954       TString name("Nz");
955       name += fCalibraMode->GetNz(1);
956       name += "Nrphi";
957       name += fCalibraMode->GetNrphi(1);
958       TTree *treePH2d = fCalibraVector->ConvertVectorPTreeHisto(fCalibraVector->GetVectorPH(),fCalibraVector->GetPlaPH(),"treePH2d",(const char *) name);
959       fout->WriteTObject(treePH2d);
960     }
961   }
962   if ((fPRF2dOn ) && (fWrite[2])) {
963     if (fHisto2d) {
964       fout->WriteTObject(fPRF2d);
965     }
966     if (fVector2d) {
967       TString name("Nz");
968       name += fCalibraMode->GetNz(2);
969       name += "Nrphi";
970       name += fCalibraMode->GetNrphi(2);
971       TTree *treePRF2d = fCalibraVector->ConvertVectorPTreeHisto(fCalibraVector->GetVectorPRF(),fCalibraVector->GetPlaPRF(),"treePRF2d",(const char *) name);
972       fout->WriteTObject(treePRF2d);
973     }
974   }
975   
976   fout->Close();
977   
978   AliInfo(Form("Execution time Write2d: R:%.2fs C:%.2fs"
979               ,stopwatch.RealTime(),stopwatch.CpuTime()));
980
981   return kTRUE;
982   
983 }
984
985 //____________Probe the histos_________________________________________________
986 Double_t *AliTRDCalibraFillHisto::StatH(TH2 *h, Int_t i)
987 {
988   //
989   // Check the number of stats in h, 0 is TH2I 1 is TProfile2D
990   // debug mode with 2 for TH2I and 3 for TProfile2D
991   // It gives a pointer to a Double_t[7] with the info following...
992   // [0] : number of calibration groups with entries
993   // [1] : minimal number of entries found
994   // [2] : calibration group number of the min
995   // [3] : maximal number of entries found
996   // [4] : calibration group number of the max
997   // [5] : mean number of entries found
998   // [6] : mean relativ error
999   //
1000
1001   Double_t *info = new Double_t[7];
1002    
1003   // Number of Xbins (detectors or groups of pads)
1004   Int_t    nbins   = h->GetNbinsX(); //number of calibration groups
1005   Int_t    nybins  = h->GetNbinsY(); //number of bins per histo
1006
1007   // Initialise
1008   Double_t nbwe = 0; //number of calibration groups with entries
1009   Double_t minentries = 0; //minimal number of entries found
1010   Double_t maxentries = 0; //maximal number of entries found
1011   Double_t placemin = 0; //calibration group number of the min
1012   Double_t placemax = -1; //calibration group number of the max
1013   Double_t meanstats = 0.0; //mean number of entries over the calibration group with at least ome entry
1014   Double_t meanrelativerror = 0.0; //mean relativ error in the TProfile2D
1015
1016   Double_t counter = 0;
1017
1018   //Debug
1019   TH1F *NbEntries = 0x0;//distribution of the number of entries
1020   TH1F *NbEntriesPerGroup = 0x0;//Number of entries per group
1021   TProfile *NbEntriesPerSp = 0x0;//Number of entries for one supermodule
1022     
1023   // Beginning of the loop over the calibration groups 
1024   for (Int_t idect = 0; idect < nbins; idect++) {
1025
1026     TH1I *projch = (TH1I *) h->ProjectionY("projch",idect+1,idect+1,(Option_t *)"e");
1027     projch->SetDirectory(0);
1028     
1029     // Number of entries for this calibration group
1030     Double_t nentries = 0.0;
1031     if((i%2) == 0){
1032       for (Int_t k = 0; k < nybins; k++) {
1033         nentries += h->GetBinContent(h->GetBin(idect+1,k+1));
1034       }
1035     }
1036     else{
1037       for (Int_t k = 0; k < nybins; k++) {
1038         nentries += ((TProfile2D *)h)->GetBinEntries(h->GetBin(idect+1,k+1));
1039         if(h->GetBinContent(h->GetBin(idect+1,k+1)) != 0) {
1040           meanrelativerror += (h->GetBinError(h->GetBin(idect+1,k+1))
1041                             / (TMath::Abs(h->GetBinContent(h->GetBin(idect+1,k+1)))));
1042           counter++;
1043         } 
1044       }
1045     }
1046
1047     //Debug
1048     if(i > 1){
1049       if((!((Bool_t)NbEntries)) && (nentries > 0)){
1050         NbEntries = new TH1F("Number of entries","Number of entries"
1051                                ,100,(Int_t)nentries/2,nentries*2);
1052         NbEntries->SetDirectory(0);
1053         NbEntriesPerGroup = new TH1F("Number of entries per group","Number of entries per group"
1054                                ,nbins,0,nbins);
1055         NbEntriesPerGroup->SetDirectory(0);
1056         NbEntriesPerSp = new TProfile("Number of entries per supermodule","Number of entries per supermodule"
1057                                ,(Int_t)(nbins/18),0,(Int_t)(nbins/18));
1058         NbEntriesPerSp->SetDirectory(0);
1059       }
1060       if(NbEntries){
1061         if(nentries > 0) NbEntries->Fill(nentries);
1062         NbEntriesPerGroup->Fill(idect+0.5,nentries);
1063         NbEntriesPerSp->Fill((idect%((Int_t)(nbins/18)))+0.5,nentries);
1064       }
1065     }
1066
1067     //min amd max
1068     if(nentries > maxentries){
1069       maxentries = nentries;
1070       placemax = idect;
1071     }
1072     if(idect == 0) {
1073       minentries = nentries;
1074     }
1075     if(nentries < minentries){
1076       minentries = nentries;
1077       placemin = idect;
1078     }
1079     //nbwe
1080     if(nentries > 0) {
1081       nbwe++;
1082       meanstats += nentries;
1083     }
1084
1085   }//calibration groups loop
1086   
1087   if(nbwe > 0) meanstats /= nbwe;
1088   if(counter > 0) meanrelativerror /= counter;
1089
1090   AliInfo(Form("There are %f calibration groups with entries",nbwe));
1091   AliInfo(Form("The minimum number of entries is %f for the group %f",minentries,placemin));
1092   AliInfo(Form("The maximum number of entries is %f for the group %f",maxentries,placemax));
1093   AliInfo(Form("The mean number of entries is %f",meanstats));
1094   if((i%2) == 1) AliInfo(Form("The mean relative error is %f",meanrelativerror));
1095
1096   info[0] = nbwe;
1097   info[1] = minentries;
1098   info[2] = placemin;
1099   info[3] = maxentries;
1100   info[4] = placemax;
1101   info[5] = meanstats;
1102   info[6] = meanrelativerror;
1103
1104   if(i > 1){
1105     gStyle->SetPalette(1);
1106     gStyle->SetOptStat(1111);
1107     gStyle->SetPadBorderMode(0);
1108     gStyle->SetCanvasColor(10);
1109     gStyle->SetPadLeftMargin(0.13);
1110     gStyle->SetPadRightMargin(0.01);
1111     TCanvas *stat = new TCanvas("stat","",50,50,600,800);
1112     stat->Divide(2,1);
1113     stat->cd(1);
1114     NbEntries->Draw("");
1115     stat->cd(2);
1116     NbEntriesPerSp->SetStats(0);
1117     NbEntriesPerSp->Draw("");
1118     TCanvas *stat1 = new TCanvas("stat1","",50,50,600,800);
1119     stat1->cd();
1120     NbEntriesPerGroup->SetStats(0);
1121     NbEntriesPerGroup->Draw("");
1122   }
1123
1124   return info;
1125
1126 }
1127
1128 //_____________________________________________________________________________
1129 void AliTRDCalibraFillHisto::SetRelativeScale(Float_t RelativeScale)
1130 {
1131   //
1132   // Set the factor that will divide the deposited charge
1133   // to fit in the histo range [0,300]
1134   //
1135  
1136   if (RelativeScale > 0.0) {
1137     fRelativeScale = RelativeScale;
1138   } 
1139   else {
1140     AliInfo("RelativeScale must be strict positif!");
1141   }
1142
1143
1144
1145 //_____________________________________________________________________________
1146 void AliTRDCalibraFillHisto::SetNz(Int_t i, Short_t Nz)
1147 {
1148   //
1149   // Set the mode of calibration group in the z direction for the parameter i
1150   // 
1151
1152   if ((Nz >= 0) && 
1153       (Nz <  5)) {
1154     fCalibraMode->SetNz(i, Nz); 
1155   }
1156   else { 
1157     AliInfo("You have to choose between 0 and 4");
1158   }
1159
1160 }
1161
1162 //_____________________________________________________________________________
1163 void AliTRDCalibraFillHisto::SetNrphi(Int_t i, Short_t Nrphi)
1164 {
1165   //
1166   // Set the mode of calibration group in the rphi direction for the parameter i
1167   //
1168  
1169   if ((Nrphi >= 0) && 
1170       (Nrphi <  7)) {
1171     fCalibraMode->SetNrphi(i ,Nrphi); 
1172   }
1173   else {
1174     AliInfo("You have to choose between 0 and 6");
1175   }
1176
1177 }
1178
1179 //____________Protected Functions______________________________________________
1180 //____________Create the 2D histo to be filled online__________________________
1181 //
1182
1183 //_____________________________________________________________________________
1184 void AliTRDCalibraFillHisto::CreatePRF2d(Int_t nn)
1185 {
1186   //
1187   // Create the 2D histos
1188   //
1189
1190   TString name("Nz");
1191   name += fCalibraMode->GetNz(2);
1192   name += "Nrphi";
1193   name += fCalibraMode->GetNrphi(2);
1194
1195   fPRF2d = new TProfile2D("PRF2d",(const Char_t *) name
1196                                  ,nn,0,nn,fNumberBinPRF,-1.5,1.5);
1197   fPRF2d->SetXTitle("Det/pad groups");
1198   fPRF2d->SetYTitle("Position x/W [pad width units]");
1199   fPRF2d->SetZTitle("Q_{i}/Q_{total}");
1200   fPRF2d->SetStats(0);
1201
1202 }
1203
1204 //_____________________________________________________________________________
1205 void AliTRDCalibraFillHisto::CreatePH2d(Int_t nn)
1206 {
1207   //
1208   // Create the 2D histos
1209   //
1210
1211   TString name("Nz");
1212   name += fCalibraMode->GetNz(1);
1213   name += "Nrphi";
1214   name += fCalibraMode->GetNrphi(1);
1215
1216   fPH2d = new TProfile2D("PH2d",(const Char_t *) name
1217                                ,nn,0,nn,fTimeMax
1218                                ,-0.5/fSf,(Float_t) (fTimeMax-0.5)/fSf);
1219   fPH2d->SetXTitle("Det/pad groups");
1220   fPH2d->SetYTitle("time [#mus]");
1221   fPH2d->SetZTitle("<PH> [a.u.]");
1222   fPH2d->SetStats(0);
1223
1224 }
1225
1226 //_____________________________________________________________________________
1227 void AliTRDCalibraFillHisto::CreateCH2d(Int_t nn)
1228 {
1229   //
1230   // Create the 2D histos
1231   //
1232
1233   TString name("Nz");
1234   name += fCalibraMode->GetNz(0);
1235   name += "Nrphi";
1236   name += fCalibraMode->GetNrphi(0);
1237
1238   fCH2d = new TH2I("CH2d",(const Char_t *) name
1239                          ,nn,0,nn,fNumberBinCharge,0,300);
1240   fCH2d->SetXTitle("Det/pad groups");
1241   fCH2d->SetYTitle("charge deposit [a.u]");
1242   fCH2d->SetZTitle("counts");
1243   fCH2d->SetStats(0);
1244   fCH2d->Sumw2();
1245
1246 }
1247
1248 //____________Offine tracking in the AliTRDtracker_____________________________
1249 void AliTRDCalibraFillHisto::FillTheInfoOfTheTrackCH()
1250 {
1251   //
1252   // For the offline tracking or mcm tracklets
1253   // This function will be called in the functions UpdateHistogram... 
1254   // to fill the info of a track for the relativ gain calibration
1255   //
1256         
1257   Int_t nb =  0; // Nombre de zones traversees
1258   Int_t fd = -1; // Premiere zone non nulle
1259   
1260   
1261   // See if the track goes through different zones
1262   for (Int_t k = 0; k < fCalibraMode->GetNfragZ(0)*fCalibraMode->GetNfragRphi(0); k++) {
1263     if (fAmpTotal[k] > 0.0) {
1264       nb++;
1265       if (nb == 1) {
1266         fd = k;
1267       }
1268     }
1269   }
1270  
1271   // If automatic scale
1272   if ((fCountRelativeScale < 100) && (fRelativeScaleAuto)) {
1273     // Take only the one zone track
1274     if (nb == 1) {
1275       fRelativeScale += fAmpTotal[fd] * 0.014 * 0.01;
1276       fCountRelativeScale++;
1277     }
1278   }
1279
1280   // We fill the CH2d after having scale with the first 100
1281   if ((fCountRelativeScale >= 100) && (fRelativeScaleAuto)) {
1282     // Case of track with only one zone
1283     if (nb == 1) {
1284       if (fHisto2d) {
1285         fCH2d->Fill(fCalibraMode->GetXbins(0)+fd+0.5,fAmpTotal[fd]/fRelativeScale);
1286       }
1287       if (fVector2d) {
1288         fCalibraVector->UpdateVectorCH(fCalibraMode->GetXbins(0)+fd,fAmpTotal[fd]/fRelativeScale);
1289       }
1290     } // Case 1 zone
1291     // Case of track with two zones
1292     if (nb == 2) {
1293       // Two zones voisines sinon rien!
1294       if ((fAmpTotal[fd]   > 0.0) && 
1295           (fAmpTotal[fd+1] > 0.0)) {
1296         // One of the two very big
1297         if (fAmpTotal[fd] > fProcent*fAmpTotal[fd+1]) {
1298           if (fHisto2d) {
1299             fCH2d->Fill(fCalibraMode->GetXbins(0)+fd+0.5,fAmpTotal[fd]/fRelativeScale);
1300           }
1301           if (fVector2d) {
1302             fCalibraVector->UpdateVectorCH(fCalibraMode->GetXbins(0)+fd,fAmpTotal[fd]/fRelativeScale);
1303           }
1304         }
1305         if (fAmpTotal[fd+1] > fProcent*fAmpTotal[fd])  {
1306           if (fHisto2d) {
1307             fCH2d->Fill(fCalibraMode->GetXbins(0)+fd+1.5,fAmpTotal[fd+1]/fRelativeScale);
1308           }
1309           if (fVector2d) {
1310             fCalibraVector->UpdateVectorCH(fCalibraMode->GetXbins(0)+fd,fAmpTotal[fd+1]/fRelativeScale);
1311           }
1312         }
1313       }
1314     } // Case 2 zones
1315   }
1316
1317   // Fill with no automatic scale
1318   if (!fRelativeScaleAuto) {
1319     // Case of track with only one zone
1320     if (nb == 1) {
1321       fNumberUsedCh[0]++;
1322       if (fHisto2d) {
1323         fCH2d->Fill(fCalibraMode->GetXbins(0)+fd+0.5,fAmpTotal[fd]/fRelativeScale);
1324       }
1325       if (fVector2d) {
1326         fCalibraVector->UpdateVectorCH(fCalibraMode->GetXbins(0)+fd,fAmpTotal[fd]/fRelativeScale);
1327       }
1328     } // Case 1 zone
1329     // Case of track with two zones
1330     if (nb == 2) {
1331       // Two zones voisines sinon rien!
1332       // Case 1
1333       if ((fAmpTotal[fd]   > 0.0) && 
1334           (fAmpTotal[fd+1] > 0.0)) {
1335         // One of the two very big
1336         if (fAmpTotal[fd] > fProcent*fAmpTotal[fd+1]) {
1337           if (fHisto2d) {
1338             fCH2d->Fill(fCalibraMode->GetXbins(0)+fd+0.5,fAmpTotal[fd]/fRelativeScale);
1339           }
1340           if (fVector2d) {
1341             fCalibraVector->UpdateVectorCH(fCalibraMode->GetXbins(0)+fd,fAmpTotal[fd]/fRelativeScale);
1342           }
1343           fNumberUsedCh[1]++;
1344         }
1345         if (fAmpTotal[fd+1] > fProcent*fAmpTotal[fd]) {
1346           if (fHisto2d) {
1347             fCH2d->Fill(fCalibraMode->GetXbins(0)+fd+1.5,fAmpTotal[fd+1]/fRelativeScale);
1348           }
1349           if (fVector2d) {
1350             fCalibraVector->UpdateVectorCH(fCalibraMode->GetXbins(0)+fd+1,fAmpTotal[fd+1]/fRelativeScale);
1351           }
1352           fNumberUsedCh[1]++;
1353         }
1354       }
1355       // Case 2
1356       if (fCalibraMode->GetNfragZ(0) > 1) {
1357         if (fAmpTotal[fd] > 0.0) {
1358           if ((fd+fCalibraMode->GetNfragZ(0)) < (fCalibraMode->GetNfragZ(0)*fCalibraMode->GetNfragRphi(0))) {
1359             if (fAmpTotal[fd+fCalibraMode->GetNfragZ(0)] > 0.0) {
1360               // One of the two very big
1361               if (fAmpTotal[fd] > fProcent*fAmpTotal[fd+fCalibraMode->GetNfragZ(0)]) {
1362                 if (fHisto2d) {
1363                   fCH2d->Fill(fCalibraMode->GetXbins(0)+fd+0.5,fAmpTotal[fd]/fRelativeScale);
1364                 }
1365                 if (fVector2d) {
1366                   fCalibraVector->UpdateVectorCH(fCalibraMode->GetXbins(0)+fd,fAmpTotal[fd]/fRelativeScale);
1367                 }
1368                 fNumberUsedCh[1]++;
1369               }
1370               if (fAmpTotal[fd+fCalibraMode->GetNfragZ(0)] > fProcent*fAmpTotal[fd]) {
1371                 if (fHisto2d) {
1372                   fCH2d->Fill(fCalibraMode->GetXbins(0)+fd+fCalibraMode->GetNfragZ(0)
1373                             + 0.5,fAmpTotal[fd+fCalibraMode->GetNfragZ(0)]/fRelativeScale);
1374                 }
1375                 fNumberUsedCh[1]++;
1376                 if (fVector2d) {
1377                   fCalibraVector->UpdateVectorCH(fCalibraMode->GetXbins(0)+fd+fCalibraMode->GetNfragZ(0)
1378                                                 ,fAmpTotal[fd+fCalibraMode->GetNfragZ(0)]/fRelativeScale);
1379                 }
1380               }
1381             }
1382           }
1383         }
1384       }
1385     } // Case 2 zones
1386
1387   }
1388
1389 }
1390
1391 //____________Offine tracking in the AliTRDtracker_____________________________
1392 void AliTRDCalibraFillHisto::ResetfVariables()
1393 {
1394   //
1395   // Reset values of fAmpTotal, fPHValue and fPHPlace for
1396   // the updateHistogram... functions
1397   //
1398
1399   // Reset the good track
1400   fGoodTrack = kTRUE;
1401   
1402   // Reset the fAmpTotal where we put value
1403   if (fCH2dOn) {
1404     for (Int_t k = 0; k < fCalibraMode->GetNfragZ(0)*fCalibraMode->GetNfragRphi(0); k++) {
1405       fAmpTotal[k] = 0.0;
1406     }
1407   }
1408   
1409   // Reset the fPHValue
1410   if (fPH2dOn) {
1411     for (Int_t k = 0; k < fTimeMax; k++) {
1412       fPHValue[k] = 0.0;
1413       fPHPlace[k] = -1;
1414     }
1415   }
1416
1417 }
1418
1419 //____________Offine tracking in the AliTRDtracker_____________________________
1420 void AliTRDCalibraFillHisto::FillTheInfoOfTheTrackPH()
1421 {
1422   //
1423   // For the offline tracking or mcm tracklets
1424   // This function will be called in the functions UpdateHistogram... 
1425   // to fill the info of a track for the drift velocity  calibration
1426   //
1427     
1428   Int_t nb  =  1; // Nombre de zones traversees 1, 2 ou plus de 3
1429   Int_t fd1 = -1; // Premiere zone non nulle
1430   Int_t fd2 = -1; // Deuxieme zone non nulle
1431   Int_t k1  = -1; // Debut de la premiere zone
1432   Int_t k2  = -1; // Debut de la seconde zone
1433
1434   // See if the track goes through different zones
1435   for (Int_t k = 0; k < fTimeMax; k++) {
1436     if (fPHValue[k] > 0.0) {
1437       if (fd1 == -1) {
1438         fd1 = fPHPlace[k];
1439         k1  = k;              
1440       }
1441       if (fPHPlace[k] != fd1) {
1442         if (fd2 == -1) {
1443           k2  = k;
1444           fd2 = fPHPlace[k];
1445           nb  = 2;
1446         }
1447         if (fPHPlace[k] != fd2) {
1448           nb = 3;
1449         }
1450       }
1451     }
1452   }
1453   
1454   // Fill 
1455   // Case of track with only one zone
1456   if (nb == 1) {
1457     fNumberUsedPh[0]++;
1458     //fd1 is the only zone
1459     for (Int_t i = 0; i < fTimeMax; i++) {
1460       if (fHisto2d) {
1461         fPH2d->Fill((fCalibraMode->GetXbins(1)+fd1)+0.5,(Float_t) i/fSf,(Float_t) fPHValue[i]);
1462       }
1463       if (fVector2d) {
1464         fCalibraVector->UpdateVectorPH((fCalibraMode->GetXbins(1)+fd1),i,fPHValue[i]);
1465       }
1466     }
1467   } // Case 1 zone
1468   // Case of track with two zones
1469   if (nb == 2) {
1470     // Two zones voisines sinon rien!
1471     // Case 1
1472     if ((fd1 == fd2+1) || 
1473         (fd2 == fd1+1)) {
1474       // One of the two fast all the think
1475       if (k2 > (k1+fDifference)) {
1476         //we choose to fill the fd1 with all the values
1477         fNumberUsedPh[1]++;
1478         for (Int_t i = 0; i < fTimeMax; i++) {
1479           if (fHisto2d) {
1480             fPH2d->Fill((fCalibraMode->GetXbins(1)+fd1)+0.5,(Float_t) i/fSf,(Float_t) fPHValue[i]);
1481           }
1482           if (fVector2d) {
1483             fCalibraVector->UpdateVectorPH((fCalibraMode->GetXbins(1)+fd1),i,fPHValue[i]);
1484           }
1485         }
1486       }
1487       if ((k2+fDifference) < fTimeMax) {
1488         //we choose to fill the fd2 with all the values
1489         fNumberUsedPh[1]++;
1490         for (Int_t i = 0; i < fTimeMax; i++) {
1491           if (fHisto2d) {
1492             fPH2d->Fill((fCalibraMode->GetXbins(1)+fd2)+0.5,(Float_t) i/fSf,(Float_t) fPHValue[i]);
1493           }
1494           if (fVector2d) {
1495             fCalibraVector->UpdateVectorPH((fCalibraMode->GetXbins(1)+fd2),i,fPHValue[i]);
1496           }
1497         }
1498       }
1499     }
1500     // Two zones voisines sinon rien!
1501     if (fCalibraMode->GetNfragZ(1) > 1) {
1502       // Case 2
1503       if ((fd1+fCalibraMode->GetNfragZ(1)) < (fCalibraMode->GetNfragZ(1)*fCalibraMode->GetNfragRphi(1))) {
1504         if (fd2 == (fd1+fCalibraMode->GetNfragZ(1))) {
1505           // One of the two fast all the think
1506           if (k2 > (k1+fDifference)) {
1507             //we choose to fill the fd1 with all the values
1508             fNumberUsedPh[1]++;
1509             for (Int_t i = 0; i < fTimeMax; i++) {
1510               if (fHisto2d) {
1511                 fPH2d->Fill((fCalibraMode->GetXbins(1)+fd1)+0.5,(Float_t) i/fSf,(Float_t) fPHValue[i]);
1512               }
1513               if (fVector2d) {
1514                 fCalibraVector->UpdateVectorPH((fCalibraMode->GetXbins(1)+fd1),i,fPHValue[i]);
1515               }
1516             }
1517           }
1518           if ((k2+fDifference) < fTimeMax) {
1519             //we choose to fill the fd2 with all the values
1520             fNumberUsedPh[1]++;
1521             for (Int_t i = 0; i < fTimeMax; i++) {
1522               if (fHisto2d) {
1523                 fPH2d->Fill((fCalibraMode->GetXbins(1)+fd2)+0.5,(Float_t) i/fSf,(Float_t) fPHValue[i]);
1524               }
1525               if (fVector2d) {
1526                 fCalibraVector->UpdateVectorPH((fCalibraMode->GetXbins(1)+fd2),i,fPHValue[i]);
1527               }
1528             }
1529           }
1530         }
1531       }
1532       // Two zones voisines sinon rien!
1533       // Case 3
1534       if ((fd1 - fCalibraMode->GetNfragZ(1)) >= 0) {
1535         if (fd2 == (fd1 - fCalibraMode->GetNfragZ(1))) {
1536           // One of the two fast all the think
1537           if (k2 > (k1 + fDifference)) {
1538             //we choose to fill the fd1 with all the values
1539             fNumberUsedPh[1]++;
1540             for (Int_t i = 0; i < fTimeMax; i++) {
1541               if (fHisto2d) {
1542                 fPH2d->Fill((fCalibraMode->GetXbins(1)+fd1)+0.5,(Float_t) i/fSf,(Float_t) fPHValue[i]);
1543               }
1544               if (fVector2d) {
1545                 fCalibraVector->UpdateVectorPH((fCalibraMode->GetXbins(1)+fd1),i,fPHValue[i]);
1546               }
1547             }
1548           }
1549           if ((k2+fDifference) < fTimeMax) {
1550             //we choose to fill the fd2 with all the values
1551             fNumberUsedPh[1]++;
1552             for (Int_t i = 0; i < fTimeMax; i++) {
1553               if (fHisto2d) {
1554                 fPH2d->Fill((fCalibraMode->GetXbins(1)+fd2)+0.5,(Float_t) i/fSf,(Float_t) fPHValue[i]);
1555               }
1556               if (fVector2d) {
1557                 fCalibraVector->UpdateVectorPH((fCalibraMode->GetXbins(1)+fd2),i,fPHValue[i]);
1558               }
1559             }
1560           }
1561         }
1562       }
1563     }
1564
1565   } // case 2 zones
1566
1567 }
1568
1569 //____________Set the pad calibration variables for the detector_______________
1570 Bool_t AliTRDCalibraFillHisto::LocalisationDetectorXbins(Int_t detector)
1571 {
1572   //
1573   // For the detector calcul the first Xbins and set the number of row
1574   // and col pads per calibration groups, the number of calibration
1575   // groups in the detector.
1576   //
1577   
1578   // first Xbins of the detector
1579   if (fCH2dOn) {
1580     fCalibraMode->CalculXBins(detector,0);
1581   }
1582   if (fPH2dOn) {
1583     fCalibraMode->CalculXBins(detector,1);
1584   }
1585   if (fPRF2dOn) {
1586     fCalibraMode->CalculXBins(detector,2);
1587   }
1588
1589   // fragmentation of idect
1590   for (Int_t i = 0; i < 3; i++) {
1591     fCalibraMode->ModePadCalibration((Int_t) GetChamber(detector),i);
1592     fCalibraMode->ModePadFragmentation((Int_t) GetPlane(detector)
1593                        , (Int_t) GetChamber(detector)
1594                        , (Int_t) GetSector(detector),i);
1595   }
1596   
1597   return kTRUE;
1598
1599 }
1600
1601 //
1602 //____________Some basic geometry function_____________________________________
1603 //
1604
1605 //_____________________________________________________________________________
1606 Int_t AliTRDCalibraFillHisto::GetPlane(Int_t d) const
1607 {
1608   //
1609   // Reconstruct the plane number from the detector number
1610   //
1611
1612   return ((Int_t) (d % 6));
1613
1614 }
1615
1616 //_____________________________________________________________________________
1617 Int_t AliTRDCalibraFillHisto::GetChamber(Int_t d) const
1618 {
1619   //
1620   // Reconstruct the chamber number from the detector number
1621   //
1622   Int_t fgkNplan = 6;
1623
1624   return ((Int_t) (d % 30) / fgkNplan);
1625
1626 }
1627
1628 //_____________________________________________________________________________
1629 Int_t AliTRDCalibraFillHisto::GetSector(Int_t d) const
1630 {
1631   //
1632   // Reconstruct the sector number from the detector number
1633   //
1634   Int_t fg = 30;
1635
1636   return ((Int_t) (d / fg));
1637
1638 }