Bug fix in chamber status (Julian)
[u/mrichter/AliRoot.git] / TRD / AliTRDPreprocessorOffline.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
17
18 /*
19   Responsible: Raphaelle Bailhache (rbailhache@ikf.uni-frankfurt.de) 
20   Code to analyze the TRD calibration and to produce OCDB entries  
21
22
23    .x ~/rootlogon.C
24    gSystem->Load("libANALYSIS");
25    gSystem->Load("libTRDcalib");
26
27    AliTRDPreprocessorOffline proces;
28    TString ocdbPath="local:////"
29    ocdbPath+=gSystem->GetFromPipe("pwd");
30
31    proces.CalibTimeGain("CalibObjects.root",run0,run1,ocdbPath);
32    proces.CalibTimeVdrift("CalibObjects.root",run0,run1,ocdbPath);
33   // take the raw calibration data from the file CalibObjects.root 
34   // and make a OCDB entry with run  validity run0-run1
35   // results are stored at the ocdbPath - local or alien ...
36   // default storage ""- data stored at current working directory 
37  
38 */
39 #include "AliLog.h"
40 #include "Riostream.h"
41 #include <fstream>
42 #include "TFile.h"
43 #include "TCanvas.h"
44 #include "TLegend.h"
45 #include "TH2I.h"
46 #include "TH1I.h"
47 #include "TH2F.h"
48 #include "TH1F.h"
49 #include "TMath.h"
50 #include "THnSparse.h"
51 #include "TProfile2D.h"
52 #include "AliTRDCalDet.h"
53 #include "AliTRDCalPad.h"
54 #include "AliCDBMetaData.h"
55 #include "AliCDBId.h"
56 #include "AliCDBManager.h"
57 #include "AliCDBStorage.h"
58 #include "AliTRDCalibraMode.h"
59 #include "AliTRDCalibraFit.h"
60 #include "AliTRDCalibraVdriftLinearFit.h"
61 #include "AliTRDCalibraExbAltFit.h"
62 #include "AliTRDPreprocessorOffline.h"
63 #include "AliTRDCalChamberStatus.h"
64 #include "AliTRDCalibChamberStatus.h"
65 #include "AliTRDCommonParam.h"
66 #include "AliCDBManager.h"
67 #include "AliCDBEntry.h"
68
69
70 ClassImp(AliTRDPreprocessorOffline)
71
72   AliTRDPreprocessorOffline::AliTRDPreprocessorOffline():
73   TNamed("TPCPreprocessorOffline","TPCPreprocessorOffline"),
74   fMethodSecond(kTRUE),
75   fNameList("TRDCalib"),
76   fCalDetGainUsed(0x0),
77   fCalDetVdriftUsed(0x0),
78   fCalDetExBUsed(0x0),
79   fCH2d(0x0),
80   fPH2d(0x0),
81   fPRF2d(0x0),
82   fSparse(0x0),
83   fAliTRDCalibraVdriftLinearFit(0x0),
84   fAliTRDCalibraExbAltFit(0x0),
85   fNEvents(0x0),
86   fAbsoluteGain(0x0),
87   fPlots(new TObjArray(kNumCalibObjs)),
88   fCalibObjects(new TObjArray(kNumCalibObjs)),
89   fFirstRunGainUsed(0),
90   fVersionGainUsed(0),
91   fSubVersionGainUsed(0),
92   fFirstRunVdriftUsed(0),
93   fVersionVdriftUsed(0), 
94   fSubVersionVdriftUsed(0),
95   fFirstRunExBUsed(0),
96   fVersionExBUsed(0), 
97   fSubVersionExBUsed(0),
98   fNoExBUsedInReco(kFALSE),
99   fSwitchOnValidation(kTRUE),
100   fVdriftValidated(kFALSE),
101   fExBValidated(kFALSE),
102   fT0Validated(kFALSE),
103   fMinStatsVdriftT0PH(800*20),
104   fMinStatsVdriftLinear(800),
105   fMinStatsGain(800),
106   fMinStatsPRF(600),
107   fMinStatsChamberStatus(20),
108   fMinSingleStatsChamberStatus(0.05),
109   fBackCorrectGain(kFALSE),  
110   fBackCorrectVdrift(kTRUE),
111   fNotEnoughStatisticsForTheGain(kFALSE),
112   fNotEnoughStatisticsForTheVdriftLinear(kFALSE),
113   fStatusNeg(0),
114   fStatusPos(0),
115   fBadCalibValidate(40),
116   fNoDataValidate(40),
117   fRMSBadCalibratedGain(20.0),
118   fRMSBadCalibratedVdrift(20.0),
119   fRMSBadCalibratedExB(20.0),
120   fRobustFitDriftVelocity(kTRUE),
121   fRobustFitExbAlt(kFALSE),
122   fAlternativeVdrfitFit(kFALSE),
123   fAlternativeExbAltFit(kFALSE)
124 {
125   //
126   // default constructor
127   //
128
129   memset(fBadCalib, 0, sizeof(Int_t) * 18);
130   memset(fNoData, 0, sizeof(Int_t) * 18);
131 }
132 //_________________________________________________________________________________________________________________
133 AliTRDPreprocessorOffline::~AliTRDPreprocessorOffline() {
134   //
135   // Destructor
136   //
137
138   if(fCalDetGainUsed) delete fCalDetGainUsed;
139   if(fCalDetVdriftUsed) delete fCalDetVdriftUsed;
140   if(fCalDetExBUsed) delete fCalDetExBUsed;
141   if(fCH2d) delete fCH2d;
142   if(fPH2d) delete fPH2d;
143   if(fPRF2d) delete fPRF2d;
144   if(fSparse) delete fSparse;
145   if(fAliTRDCalibraVdriftLinearFit) delete fAliTRDCalibraVdriftLinearFit;
146   if(fAliTRDCalibraExbAltFit) delete fAliTRDCalibraExbAltFit;
147   if(fNEvents) delete fNEvents;
148   if(fAbsoluteGain) delete fAbsoluteGain;
149   if(fPlots) delete fPlots;
150   if(fCalibObjects) delete fCalibObjects;
151   
152 }
153 //___________________________________________________________________________________
154 void AliTRDPreprocessorOffline::Process(const Char_t* file, Int_t startRunNumber, Int_t endRunNumber, TString ocdbStorage) 
155 {
156   //
157   // Process to the gain, vdrift, timeoffset, exb and chamber status calibration
158   //
159   
160   if(SetCalDetGain(startRunNumber,fVersionGainUsed,fSubVersionGainUsed) && SetCalDetVdriftExB(startRunNumber,fVersionVdriftUsed,fSubVersionVdriftUsed,fVersionExBUsed,fSubVersionExBUsed)) {
161     
162     CalibVdriftT0(file,startRunNumber,endRunNumber,ocdbStorage);
163     CalibGain(file,startRunNumber,endRunNumber,ocdbStorage);
164     CalibChamberStatus(file,startRunNumber,endRunNumber,ocdbStorage);
165     CalibExbAlt(file,startRunNumber,endRunNumber,ocdbStorage);
166     
167   }
168   
169   PrintStatus();
170   
171 }
172 //___________________________________________________________________________________________________________________
173
174 void AliTRDPreprocessorOffline::CalibVdriftT0(const Char_t* file, Int_t startRunNumber, Int_t endRunNumber, TString ocdbStorage){
175   //
176   // make calibration of the drift velocity
177   // Input parameters:
178   //      file                             - the location of input file
179   //      startRunNumber, endRunNumber     - run validity period 
180   //      ocdbStorage                      - path to the OCDB storage
181   //                                       - if empty - local storage 'pwd' uesed
182   if (ocdbStorage.Length()<=0) ocdbStorage="local://"+gSystem->GetFromPipe("pwd")+"/OCDB";
183   //
184   // 1. Initialization 
185   //
186   fVdriftValidated = kTRUE;
187   fT0Validated = kTRUE;
188   fExBValidated = kTRUE;
189   fNotEnoughStatisticsForTheVdriftLinear = kFALSE;
190   //
191   // 2. extraction of the information
192   //
193   if(ReadVdriftLinearFitGlobal(file) && fCalDetVdriftUsed && fCalDetExBUsed) AnalyzeVdriftLinearFit();
194   if(ReadVdriftT0Global(file)) AnalyzeVdriftT0();
195   //
196   // 3. Append QA plots
197   //
198   //MakeDefaultPlots(fVdriftArray,fVdriftArray);
199   //
200   //
201   // 4. validate OCDB entries
202   //
203   if(fSwitchOnValidation==kTRUE && ValidateVdrift()==kFALSE) { 
204     //AliError("TRD vdrift OCDB parameters out of range!");
205     fVdriftValidated = kFALSE;
206   }
207   if(fSwitchOnValidation==kTRUE && ValidateT0()==kFALSE) { 
208     //AliError("TRD t0 OCDB parameters out of range!");
209     fT0Validated = kFALSE;
210   }
211   if(fSwitchOnValidation==kTRUE && ValidateExB()==kFALSE) { 
212     //AliError("TRD t0 OCDB parameters out of range!");
213     fExBValidated = kFALSE;
214   }
215   //
216   // 5. update of OCDB
217   //
218   //
219   if(fVdriftValidated) UpdateOCDBVdrift(startRunNumber,endRunNumber,ocdbStorage);
220   if(fT0Validated) UpdateOCDBT0(startRunNumber,endRunNumber,ocdbStorage);
221   if(fExBValidated) UpdateOCDBExB(startRunNumber,endRunNumber,ocdbStorage);
222   
223 }
224 //___________________________________________________________________________________________________________________
225
226 void AliTRDPreprocessorOffline::CalibExbAlt(const Char_t* file, Int_t startRunNumber, Int_t endRunNumber, TString ocdbStorage){
227   //
228   // make calibration of the drift velocity
229   // Input parameters:
230   //      file                             - the location of input file
231   //      startRunNumber, endRunNumber     - run validity period 
232   //      ocdbStorage                      - path to the OCDB storage
233   //                                       - if empty - local storage 'pwd' uesed
234   if (ocdbStorage.Length()<=0) ocdbStorage="local://"+gSystem->GetFromPipe("pwd")+"/OCDB";
235   //
236   // 1. Initialization 
237   //
238
239   //
240   // 2. extraction of the information
241   //
242   if(ReadExbAltFitGlobal(file)) AnalyzeExbAltFit();
243   //
244   // 3. Append QA plots
245   //
246   //MakeDefaultPlots(fVdriftArray,fVdriftArray);
247   //
248   //
249   // 4. validate OCDB entries
250   //
251   //
252   // 5. update of OCDB
253   //
254   //
255   UpdateOCDBExBAlt(startRunNumber,endRunNumber,ocdbStorage);
256   
257 }
258
259 //_________________________________________________________________________________________________________________
260
261 void AliTRDPreprocessorOffline::CalibGain(const Char_t* file, Int_t startRunNumber, Int_t endRunNumber, TString ocdbStorage){
262   //
263   // make calibration of the drift velocity
264   // Input parameters:
265   //      file                             - the location of input file
266   //      startRunNumber, endRunNumber     - run validity period 
267   //      ocdbStorage                      - path to the OCDB storage
268   //                                       - if empty - local storage 'pwd' uesed
269   if (ocdbStorage.Length()<=0) ocdbStorage="local://"+gSystem->GetFromPipe("pwd")+"/OCDB";
270   //
271   fNotEnoughStatisticsForTheGain = kFALSE;
272   //
273   // 1. Initialization 
274   if(!ReadGainGlobal(file)) return;
275   //
276   //
277   // 2. extraction of the information
278   //
279   AnalyzeGain();
280   if(fBackCorrectGain) CorrectFromDetGainUsed();
281   //if(fBackCorrectVdrift) CorrectFromDetVdriftUsed();
282   //
283   // 3. Append QA plots
284   //
285   //MakeDefaultPlots(fVdriftArray,fVdriftArray);
286   //
287   //
288   // 4. validate OCDB entries
289   //
290   if(fSwitchOnValidation==kTRUE && ValidateGain()==kFALSE) { 
291     //AliError("TRD gain OCDB parameters out of range!");
292     return;
293   }
294   //
295   // 5. update of OCDB
296   //
297   //
298   if((!fCalDetVdriftUsed) || (fCalDetVdriftUsed && fVdriftValidated)) UpdateOCDBGain(startRunNumber,endRunNumber,ocdbStorage);
299   
300   
301 }
302 //________________________________________________________________________________________________________________
303
304 void AliTRDPreprocessorOffline::CalibPRF(const Char_t* file, Int_t startRunNumber, Int_t endRunNumber, TString ocdbStorage){
305   //
306   // make calibration of the drift velocity
307   // Input parameters:
308   //      file                             - the location of input file
309   //      startRunNumber, endRunNumber     - run validity period 
310   //      ocdbStorage                      - path to the OCDB storage
311   //                                       - if empty - local storage 'pwd' uesed
312   if (ocdbStorage.Length()<=0) ocdbStorage="local://"+gSystem->GetFromPipe("pwd")+"/OCDB";
313   //
314   // 1. Initialization 
315   if(!ReadPRFGlobal(file)) return;
316   //
317   //
318   // 2. extraction of the information
319   //
320   AnalyzePRF();
321   //
322   // 3. Append QA plots
323   //
324   //MakeDefaultPlots(fVdriftArray,fVdriftArray);
325   //
326   //
327   //
328   // 4. validate OCDB entries
329   //
330   if(fSwitchOnValidation==kTRUE && ValidatePRF()==kFALSE) { 
331     //AliError("TRD prf OCDB parameters out of range!");
332     return;
333   }
334   //
335   // 5. update of OCDB
336   //
337   //
338   UpdateOCDBPRF(startRunNumber,endRunNumber,ocdbStorage);
339   
340 }
341 //________________________________________________________________________________________________________________
342
343 void AliTRDPreprocessorOffline::CalibChamberStatus(const Char_t* file, Int_t startRunNumber, Int_t endRunNumber, TString ocdbStorage){
344   //
345   // make calibration of the chamber status
346   // Input parameters:
347   //      startRunNumber, endRunNumber     - run validity period 
348   //      ocdbStorage                      - path to the OCDB storage
349   //                                       - if empty - local storage 'pwd' uesed
350   if (ocdbStorage.Length()<=0) ocdbStorage="local://"+gSystem->GetFromPipe("pwd")+"/OCDB";
351   //
352   //
353   // 1. Initialization  
354   if(!ReadStatusGlobal(file)) return;
355   //
356   //
357   //
358   // 2. extraction of the information
359   //
360   if(!AnalyzeChamberStatus()) return;
361   //
362   // 3. Append QA plots
363   //
364   //MakeDefaultPlots(fVdriftArray,fVdriftArray);
365   //
366   //
367   //
368   // 4. validate OCDB entries
369   //
370   if(fSwitchOnValidation==kTRUE && ValidateChamberStatus()==kFALSE) { 
371     //AliError("TRD Chamber status OCDB parameters not ok!");
372     return;
373   }
374   //
375   // 5. update of OCDB
376   //
377   //
378   if((!fNotEnoughStatisticsForTheVdriftLinear) && (!fNotEnoughStatisticsForTheGain)) UpdateOCDBChamberStatus(startRunNumber,endRunNumber,ocdbStorage);
379   //UpdateOCDBChamberStatus(startRunNumber,endRunNumber,ocdbStorage);
380   
381 }
382 //______________________________________________________________________________________________________
383 Bool_t AliTRDPreprocessorOffline::Init(const Char_t* fileName){
384   //
385   // read the calibration used during the reconstruction
386   // 
387
388   if(ReadVdriftT0Global(fileName)) {
389     
390     TString nameph = fPH2d->GetTitle();
391     fFirstRunVdriftUsed = GetFirstRun(nameph); 
392     fVersionVdriftUsed = GetVersion(nameph);  
393     fSubVersionVdriftUsed = GetSubVersion(nameph);    
394
395     //printf("Found Version %d, Subversion %d for vdrift\n",fVersionVdriftUsed,fSubVersionVdriftUsed);
396   
397   }
398
399   if(ReadGainGlobal(fileName)) {
400
401     TString namech = fCH2d->GetTitle();
402     fFirstRunGainUsed = GetFirstRun(namech); 
403     fVersionGainUsed = GetVersion(namech);  
404     fSubVersionGainUsed = GetSubVersion(namech);    
405
406     //printf("Found Version %d, Subversion %d for gain\n",fVersionGainUsed,fSubVersionGainUsed);
407
408   }
409   
410   if(ReadVdriftLinearFitGlobal(fileName)) {
411
412     TString namelinear = fAliTRDCalibraVdriftLinearFit->GetNameCalibUsed();
413     fFirstRunExBUsed = GetFirstRun(namelinear); 
414     fVersionExBUsed = GetVersion(namelinear);  
415     fSubVersionExBUsed = GetSubVersion(namelinear);   
416
417     //printf("Found Version %d, Subversion %d, run %d for ExB\n",fVersionExBUsed,fSubVersionExBUsed,fFirstRunExBUsed);
418     
419   }
420    
421   if(fVersionVdriftUsed == 0) fStatusPos = fStatusPos |kVdriftErrorOld;
422   if(fVersionGainUsed == 0) fStatusPos = fStatusPos | kGainErrorOld;
423  
424   return kTRUE;
425   
426 }
427 //___________________________________________________________________________________________________________________
428
429 Bool_t AliTRDPreprocessorOffline::ReadStatusGlobal(const Char_t* fileName){
430   //
431   // read calibration entries from file
432   // 
433   if(fSparse) return kTRUE;
434   TFile fcalib(fileName);
435   TObjArray * array = (TObjArray*)fcalib.Get(fNameList);
436   if (array){
437     fSparse = (THnSparseI *) array->FindObject("NumberOfEntries");
438     if(!fSparse) return kFALSE;
439   }
440   else 
441     return kFALSE;
442   
443   return kTRUE;
444   
445 }
446 //___________________________________________________________________________________________________________________
447
448 Bool_t AliTRDPreprocessorOffline::ReadGainGlobal(const Char_t* fileName){
449   //
450   // read calibration entries from file
451   // 
452   if(fCH2d) return kTRUE;
453   TFile fcalib(fileName);
454   TObjArray * array = (TObjArray*)fcalib.Get(fNameList);
455   if (array){
456     TH2I *ch2d = (TH2I *) array->FindObject("CH2d");
457     if(!ch2d) return kFALSE;
458     fCH2d = (TH2I*)ch2d->Clone();
459     //fNEvents = (TH1I *) array->FindObject("NEvents");
460     //fAbsoluteGain = (TH2F *) array->FindObject("AbsoluteGain");
461   }else{
462     TH2I *ch2d = (TH2I *) fcalib.Get("CH2d");
463     if(!ch2d) return kFALSE;
464     fCH2d = (TH2I*)ch2d->Clone();
465     //fNEvents = (TH1I *) fcalib.Get("NEvents");
466     //fAbsoluteGain = (TH2F *) fcalib.Get("AbsoluteGain");
467   }
468   fCH2d->SetDirectory(0);
469   //printf("title of CH2d %s\n",fCH2d->GetTitle());
470
471   return kTRUE;
472   
473 }
474 //_________________________________________________________________________________________________________________
475
476 Bool_t AliTRDPreprocessorOffline::ReadVdriftT0Global(const Char_t* fileName){
477   //
478   // read calibration entries from file
479   // 
480   if(fPH2d) return kTRUE;
481   TFile fcalib(fileName);
482   TObjArray * array = (TObjArray*)fcalib.Get(fNameList);
483   if (array){
484     TProfile2D *ph2d = (TProfile2D *) array->FindObject("PH2d");
485     if(!ph2d) return kFALSE;
486     fPH2d = (TProfile2D*)ph2d->Clone();
487     //fNEvents = (TH1I *) array->FindObject("NEvents");
488   }else{
489     TProfile2D *ph2d = (TProfile2D *) fcalib.Get("PH2d");
490     if(!ph2d) return kFALSE;
491     fPH2d = (TProfile2D*)ph2d->Clone();
492     //fNEvents = (TH1I *) fcalib.Get("NEvents");
493   }
494   fPH2d->SetDirectory(0);
495   //printf("title of PH2d %s\n",fPH2d->GetTitle());
496   
497   return kTRUE;
498   
499 }
500 //___________________________________________________________________________________________________________________
501
502 Bool_t AliTRDPreprocessorOffline::ReadVdriftLinearFitGlobal(const Char_t* fileName){
503   //
504   // read calibration entries from file
505   // 
506   if(fAliTRDCalibraVdriftLinearFit) return kTRUE;
507   TFile fcalib(fileName);
508   TObjArray * array = (TObjArray*)fcalib.Get(fNameList);
509   if (array){
510     fAliTRDCalibraVdriftLinearFit = (AliTRDCalibraVdriftLinearFit *) array->FindObject("AliTRDCalibraVdriftLinearFit");
511     //fNEvents = (TH1I *) array->FindObject("NEvents");
512   }else{
513     fAliTRDCalibraVdriftLinearFit = (AliTRDCalibraVdriftLinearFit *) fcalib.Get("AliTRDCalibraVdriftLinearFit");
514     //fNEvents = (TH1I *) fcalib.Get("NEvents");
515   }
516   if(!fAliTRDCalibraVdriftLinearFit) {
517     //printf("No AliTRDCalibraVdriftLinearFit\n");
518     return kFALSE;
519   }
520   return kTRUE;
521   
522 }
523 //_____________________________________________________________________________________________________________
524 Bool_t AliTRDPreprocessorOffline::ReadExbAltFitGlobal(const Char_t* fileName){
525   //
526   // read calibration entries from file
527   // 
528   if(fAliTRDCalibraExbAltFit) return kTRUE;
529   TFile fcalib(fileName);
530   TObjArray * array = (TObjArray*)fcalib.Get(fNameList);
531   if (array){
532     fAliTRDCalibraExbAltFit = (AliTRDCalibraExbAltFit *) array->FindObject("AliTRDCalibraExbAltFit");
533     //fNEvents = (TH1I *) array->FindObject("NEvents");
534   }else{
535     fAliTRDCalibraExbAltFit = (AliTRDCalibraExbAltFit *) fcalib.Get("AliTRDCalibraExbAltFit");
536     //fNEvents = (TH1I *) fcalib.Get("NEvents");
537   }
538   if(!fAliTRDCalibraExbAltFit) {
539     //printf("No AliTRDCalibraExbAltFit\n");
540     return kFALSE;
541   }
542   return kTRUE;
543   
544 }
545 //_____________________________________________________________________________________________________________
546
547 Bool_t AliTRDPreprocessorOffline::ReadPRFGlobal(const Char_t* fileName){
548   //
549   // read calibration entries from file
550   // 
551   if(fPRF2d) return kTRUE;
552   TFile fcalib(fileName);
553   TObjArray * array = (TObjArray*)fcalib.Get(fNameList);
554   if (array){
555     TProfile2D *prf2d = (TProfile2D *) array->FindObject("PRF2d");
556     if(!prf2d) return kFALSE;
557     fPRF2d = (TProfile2D*)prf2d->Clone();
558     //fNEvents = (TH1I *) array->FindObject("NEvents");
559   }else{
560     TProfile2D *prf2d = (TProfile2D *) fcalib.Get("PRF2d");
561     if(!prf2d) return kFALSE;
562     fPRF2d = (TProfile2D*)prf2d->Clone();
563     //fNEvents = (TH1I *) fcalib.Get("NEvents");
564   }
565   fPRF2d->SetDirectory(0);
566   //printf("title of PRF2d %s\n",fPRF2d->GetTitle());
567   
568   return kTRUE;
569
570 }
571 //__________________________________________________________________________________________________________
572
573 Bool_t AliTRDPreprocessorOffline::AnalyzeGain(){
574   //
575   // Analyze gain - produce the calibration objects
576   //
577
578   AliTRDCalibraFit *calibra = AliTRDCalibraFit::Instance();
579   calibra->SetMinEntries(fMinStatsGain); // If there is less than 1000 entries in the histo: no fit
580   calibra->AnalyseCH(fCH2d);
581
582   Int_t nbtg = 6*4*18*((Int_t) ((AliTRDCalibraMode *)calibra->GetCalibraMode())->GetDetChamb0(0))
583     + 6*  18*((Int_t) ((AliTRDCalibraMode *)calibra->GetCalibraMode())->GetDetChamb2(0));
584   Int_t nbfit       = calibra->GetNumberFit();
585   Int_t nbE         = calibra->GetNumberEnt();
586
587
588   Bool_t ok = kFALSE;
589   Bool_t meanother = kFALSE;
590   // enough statistics
591   if ((nbtg >                  0) && 
592       (nbfit        >= 0.5*nbE) && (nbE > 30)) {
593     // create the cal objects
594     if(!fBackCorrectGain) {
595       calibra->PutMeanValueOtherVectorFit(1,kTRUE);
596       meanother = kTRUE;
597     }
598     TObjArray object           = calibra->GetVectorFit();
599     AliTRDCalDet *calDetGain   = calibra->CreateDetObjectGain(&object,meanother);
600     TH1F *coefGain  = calDetGain->MakeHisto1DAsFunctionOfDet();
601     // Put them in the array
602     fCalibObjects->AddAt(calDetGain,kGain);
603     fPlots->AddAt(coefGain,kGain);
604     //
605     ok = kTRUE;
606   }
607   else {
608     fNotEnoughStatisticsForTheGain = kTRUE;
609     Int_t minStatsGain = fMinStatsGain*30;
610     calibra->SetMinEntries(minStatsGain); // Because we do it for all, we increase this
611     Double_t gainoverallnotnormalized =  calibra->AnalyseCHAllTogether(fCH2d);
612     if(fCalDetGainUsed && (gainoverallnotnormalized > 0.0)) {
613       AliTRDCalDet *calDetGain = new AliTRDCalDet(*fCalDetGainUsed);
614       Double_t oldmean = fCalDetGainUsed->CalcMean(kFALSE);
615       //printf("oldmean %f\n",oldmean);
616       if(oldmean > 0.0)  {
617         Double_t scalefactor = calibra->GetScaleFactorGain();
618         //printf("Correction factor %f\n",gainoverallnotnormalized*scalefactor);
619         calDetGain->Multiply(gainoverallnotnormalized*scalefactor/oldmean);
620         //printf("newmean %f\n",calDetGain->CalcMean(kFALSE));
621         TH1F *coefGain  = calDetGain->MakeHisto1DAsFunctionOfDet();
622         fCalibObjects->AddAt(calDetGain,kGain);
623         fPlots->AddAt(coefGain,kGain);
624         // 
625         ok = kTRUE;
626         fStatusNeg = fStatusNeg | kGainNotEnoughStatsButFill;
627       }
628       else {
629         fStatusPos = fStatusPos | kGainErrorOld;
630       }      
631     }
632     else {
633       if(gainoverallnotnormalized <= 0.0) fStatusNeg = fStatusNeg | kGainNotEnoughStatsNotFill;
634       if(!fCalDetGainUsed) fStatusPos = fStatusPos | kGainErrorOld;
635     }
636   }
637   
638   calibra->ResetVectorFit();
639   
640   return ok;
641   
642 }
643 //_____________________________________________________________________________________________________
644 Bool_t AliTRDPreprocessorOffline::AnalyzeVdriftT0(){
645   //
646   // Analyze VdriftT0 - produce the calibration objects
647   //
648
649   AliTRDCalibraFit *calibra = AliTRDCalibraFit::Instance();
650   calibra->SetMinEntries(fMinStatsVdriftT0PH); // If there is less than 1000 entries in the histo: no fit
651   calibra->AnalysePH(fPH2d);
652
653   Int_t nbtg = 6*4*18*((Int_t) ((AliTRDCalibraMode *)calibra->GetCalibraMode())->GetDetChamb0(1))
654     + 6*  18*((Int_t) ((AliTRDCalibraMode *)calibra->GetCalibraMode())->GetDetChamb2(1));
655   Int_t nbfit       = calibra->GetNumberFit();
656   Int_t nbfitSuccess = calibra->GetNumberFitSuccess();
657   Int_t nbE         = calibra->GetNumberEnt();
658
659   //printf("nbtg %d, nbfit %d, nbE %d, nbfitSuccess %d\n",nbtg,nbfit,nbE,nbfitSuccess);
660
661   Bool_t ok = kFALSE;
662   if ((nbtg >                  0) && 
663       (nbfit        >= 0.5*nbE) && (nbE > 30) && (nbfitSuccess > 30)) {
664     //printf("Pass the cut for VdriftT0\n");
665     // create the cal objects
666     calibra->RemoveOutliers(1,kFALSE);
667     calibra->PutMeanValueOtherVectorFit(1,kFALSE);
668     calibra->RemoveOutliers2(kFALSE);
669     calibra->PutMeanValueOtherVectorFit2(1,kFALSE);
670     //
671     TObjArray object = calibra->GetVectorFit();
672     AliTRDCalDet *calDetVdrift = calibra->CreateDetObjectVdrift(&object);
673     TH1F *coefVdriftPH  = calDetVdrift->MakeHisto1DAsFunctionOfDet();
674     AliTRDCalPad *calPadVdrift = (AliTRDCalPad *)calibra->CreatePadObjectVdrift(&object,calDetVdrift);
675     TH1F *coefPadVdrift   = calPadVdrift->MakeHisto1D();
676     object       = calibra->GetVectorFit2();
677     AliTRDCalDet *calDetT0  = calibra->CreateDetObjectT0(&object);
678     TH1F *coefT0  = calDetT0->MakeHisto1DAsFunctionOfDet();
679     AliTRDCalPad *calPadT0 = (AliTRDCalPad *)calibra->CreatePadObjectT0(&object,calDetT0);
680     TH1F *coefPadT0  = calPadT0->MakeHisto1D();
681     // Put them in the array
682     fCalibObjects->AddAt(calDetT0,kT0PHDet);
683     fCalibObjects->AddAt(calDetVdrift,kVdriftPHDet);
684     fCalibObjects->AddAt(calPadT0,kT0PHPad);
685     fCalibObjects->AddAt(calPadVdrift,kVdriftPHPad);
686     fPlots->AddAt(coefVdriftPH,kVdriftPHDet);
687     fPlots->AddAt(coefT0,kT0PHDet);
688     fPlots->AddAt(coefPadVdrift,kVdriftPHPad);
689     fPlots->AddAt(coefPadT0,kT0PHPad);
690     //
691     ok = kTRUE;
692   }
693   else {
694     //printf("Not enough stats timeoffset\n");
695     fStatusNeg = fStatusNeg | kTimeOffsetNotEnoughStatsNotFill;
696   }
697   calibra->ResetVectorFit();
698  
699   return ok;
700   
701 }
702 //____________________________________________________________________________________________________________________
703 Bool_t AliTRDPreprocessorOffline::AnalyzeVdriftLinearFit(){
704   //
705   // Analyze vdrift linear fit - produce the calibration objects
706   //
707
708   //printf("Analyse linear fit\n");
709
710   
711   AliTRDCalibraFit *calibra = AliTRDCalibraFit::Instance();
712   calibra->SetCalDetVdriftExB(fCalDetVdriftUsed,fCalDetExBUsed);
713   calibra->SetMinEntries(fMinStatsVdriftLinear); // If there is less than 1000 entries in the histo: no fit
714   //printf("Fill PE Array\n");
715   fAliTRDCalibraVdriftLinearFit->SetRobustFit(fRobustFitDriftVelocity);
716   if(!fAlternativeVdrfitFit)
717     fAliTRDCalibraVdriftLinearFit->FillPEArray();
718   else
719     fAliTRDCalibraVdriftLinearFit->FillPEArray2();
720   //printf("AliTRDCalibraFit\n");
721   calibra->AnalyseLinearFitters(fAliTRDCalibraVdriftLinearFit);
722   //printf("After\n");
723
724   //Int_t nbtg        = 540;
725   Int_t nbfit       = calibra->GetNumberFit();
726   Int_t nbE         = calibra->GetNumberEnt();
727
728   
729   Bool_t ok = kFALSE;
730   // enough statistics
731   if ((nbfit        >= 0.5*nbE) && (nbE > 30)) {
732     // create the cal objects
733     //calibra->RemoveOutliers(1,kTRUE);
734     calibra->PutMeanValueOtherVectorFit(1,kTRUE);
735     //calibra->RemoveOutliers2(kTRUE);
736     calibra->PutMeanValueOtherVectorFit2(1,kTRUE);
737     //
738     TObjArray object  = calibra->GetVectorFit();
739     AliTRDCalDet *calDetVdrift = calibra->CreateDetObjectVdrift(&object,kTRUE);
740     TH1F *coefDriftLinear      = calDetVdrift->MakeHisto1DAsFunctionOfDet();
741     object                     = calibra->GetVectorFit2();
742     AliTRDCalDet *calDetLorentz = calibra->CreateDetObjectLorentzAngle(&object);
743     TH1F *coefLorentzAngle = calDetLorentz->MakeHisto1DAsFunctionOfDet();
744     //if(!calDetLorentz) printf("No lorentz created\n");
745     // Put them in the array
746     fCalibObjects->AddAt(calDetVdrift,kVdriftLinear);
747     fCalibObjects->AddAt(calDetLorentz,kLorentzLinear);
748     fPlots->AddAt(coefDriftLinear,kVdriftLinear);
749     fPlots->AddAt(coefLorentzAngle,kLorentzLinear);
750     //
751     ok = kTRUE;
752   }
753   else {
754     fNotEnoughStatisticsForTheVdriftLinear = kTRUE;
755     Int_t minNumberOfEntriesForAll = fMinStatsVdriftLinear*30;
756     calibra->SetMinEntries(minNumberOfEntriesForAll); // Because we do it for all, we increase this
757     Double_t vdriftoverall = -100.0;
758     Double_t exboverall = 100.0;
759     calibra->AnalyseLinearFittersAllTogether(fAliTRDCalibraVdriftLinearFit,vdriftoverall,exboverall);
760     //printf("Found mean vdrift %f and exb %f\n",vdriftoverall,exboverall);
761     if(fCalDetVdriftUsed && (vdriftoverall > 0.0) && (exboverall < 70.0)) {
762       AliTRDCalDet *calDetVdrift = new AliTRDCalDet(*fCalDetVdriftUsed);
763       AliTRDCalDet *calDetLorentz = new AliTRDCalDet(*fCalDetExBUsed);
764       Double_t oldmeanvdrift = fCalDetVdriftUsed->CalcMean(kFALSE);
765       Double_t oldmeanexb = fCalDetExBUsed->CalcMean(kFALSE);
766       //printf("oldmean %f\n",oldmean);
767       if((oldmeanvdrift > 0.0) && (oldmeanexb < 70.0))  {
768         //printf("Correction factor %f\n",vdriftoverall);
769         calDetVdrift->Multiply(vdriftoverall/oldmeanvdrift);
770         calDetLorentz->Multiply(exboverall/oldmeanexb);
771         //printf("newmean %f\n",calDetVdrift->CalcMean(kFALSE));
772         TH1F *coefDriftLinear  = calDetVdrift->MakeHisto1DAsFunctionOfDet();
773         TH1F *coefLorentzAngle = calDetLorentz->MakeHisto1DAsFunctionOfDet();
774         // Put them in the array
775         fCalibObjects->AddAt(calDetVdrift,kVdriftLinear);
776         fCalibObjects->AddAt(calDetLorentz,kLorentzLinear);
777         fPlots->AddAt(coefDriftLinear,kVdriftLinear);
778         fPlots->AddAt(coefLorentzAngle,kLorentzLinear);
779         // 
780         ok = kTRUE;
781         fStatusNeg = fStatusNeg | kVdriftNotEnoughStatsButFill;
782       }
783       else {
784         if(oldmeanvdrift) fStatusPos = fStatusPos | kVdriftErrorOld;
785         if(oldmeanexb) fStatusPos = fStatusPos | kExBErrorOld;
786       }      
787     }
788     else {
789       if((vdriftoverall <= 0.0) && (exboverall > 70.0)) fStatusNeg = fStatusNeg | kVdriftNotEnoughStatsNotFill;
790       if(!fCalDetVdriftUsed) fStatusPos = fStatusPos | kVdriftErrorOld;
791       if(!fCalDetExBUsed) fStatusPos = fStatusPos | kExBErrorOld;
792     }
793   }
794   
795   calibra->ResetVectorFit();
796   
797   return ok;
798   
799 }
800 //________________________________________________________________________________________________________________
801
802 Bool_t AliTRDPreprocessorOffline::AnalyzeExbAltFit(){
803   //
804   // Analyze vdrift linear fit - produce the calibration objects
805   //
806
807   //printf("Analyse linear fit\n");
808
809   
810   AliTRDCalibraFit *calibra = AliTRDCalibraFit::Instance();
811   calibra->SetMinEntries(fMinStatsVdriftLinear); // If there is less than 1000 entries in the histo: no fit
812   //printf("Fill PE Array\n");
813   fAliTRDCalibraExbAltFit->SetRobustFit(fRobustFitExbAlt);
814   if(!fAlternativeExbAltFit)
815     fAliTRDCalibraExbAltFit->FillPEArray();
816   else
817     fAliTRDCalibraExbAltFit->FillPEArray2();
818   //printf("AliTRDCalibraFit\n");
819   calibra->AnalyseExbAltFit(fAliTRDCalibraExbAltFit);
820   //printf("After\n");
821
822   //Int_t nbtg        = 540;
823   Int_t nbfit       = calibra->GetNumberFit();
824   Int_t nbE         = calibra->GetNumberEnt();
825
826   
827   Bool_t ok = kFALSE;
828   // enough statistics
829   if ((nbfit        >= 0.5*nbE) && (nbE > 30)) {
830     // create the cal objects
831     //calibra->RemoveOutliers(1,kTRUE);
832     calibra->PutMeanValueOtherVectorFit2(1,kTRUE);
833     //
834     TObjArray object  = calibra->GetVectorFit2();
835     AliTRDCalDet *calDetLorentz = calibra->CreateDetObjectExbAlt(&object);
836     TH1F *coefLorentzAngle = calDetLorentz->MakeHisto1DAsFunctionOfDet();
837     //if(!calDetLorentz) printf("No lorentz created\n");
838     // Put them in the array
839     fCalibObjects->AddAt(calDetLorentz,kExbAlt);
840     fPlots->AddAt(coefLorentzAngle,kExbAlt);
841     //
842     ok = kTRUE;
843   }
844   
845   calibra->ResetVectorFit();
846   
847   return ok;
848   
849 }
850 //________________________________________________________________________________________________________________
851
852 Bool_t AliTRDPreprocessorOffline::AnalyzePRF(){
853   //
854   // Analyze PRF - produce the calibration objects
855   //
856
857   AliTRDCalibraFit *calibra = AliTRDCalibraFit::Instance();
858   calibra->SetMinEntries(fMinStatsPRF); // If there is less than 1000 entries in the histo: no fit
859   calibra->AnalysePRFMarianFit(fPRF2d);
860
861   Int_t nbtg = 6*4*18*((Int_t) ((AliTRDCalibraMode *)calibra->GetCalibraMode())->GetDetChamb0(2))
862     + 6*  18*((Int_t) ((AliTRDCalibraMode *)calibra->GetCalibraMode())->GetDetChamb2(2));
863   Int_t nbfit       = calibra->GetNumberFit();
864   Int_t nbE         = calibra->GetNumberEnt();
865
866   
867   Bool_t ok = kFALSE;
868   // enough statistics
869   if ((nbtg >                  0) && 
870       (nbfit        >= 0.95*nbE) && (nbE > 30)) {
871     // create the cal objects
872     TObjArray object  = calibra->GetVectorFit();
873     AliTRDCalPad *calPadPRF = (AliTRDCalPad*) calibra->CreatePadObjectPRF(&object);
874     TH1F *coefPRF           = calPadPRF->MakeHisto1D();
875     // Put them in the array
876     fCalibObjects->AddAt(calPadPRF,kPRF);
877     fPlots->AddAt(coefPRF,kPRF);
878     //
879     ok = kTRUE;
880   }
881   
882   calibra->ResetVectorFit();
883   
884   return ok;
885   
886 }
887
888 //_____________________________________________________________________________
889 Bool_t AliTRDPreprocessorOffline::AnalyzeChamberStatus()
890 {
891   //
892   // Produce AliTRDCalChamberStatus out of calibration results
893   //
894   
895   // set up AliTRDCalibChamberStatus
896   AliTRDCalibChamberStatus *ChamberStatus = new AliTRDCalibChamberStatus();
897   ChamberStatus->SetSparseI(fSparse);
898   ChamberStatus->AnalyseHisto(fMinStatsChamberStatus, fMinSingleStatsChamberStatus);
899   // get AliTRDCalChamberStatus
900   AliTRDCalChamberStatus *CalChamberStatus = ChamberStatus->GetCalChamberStatus();
901
902   // get calibration objects
903   AliTRDCalDet *calDetGain   = (AliTRDCalDet *) fCalibObjects->At(kGain);
904   AliTRDCalDet *calDetVDrift = (AliTRDCalDet *) fCalibObjects->At(kVdriftLinear);
905   AliTRDCalDet *calDetExB    = (AliTRDCalDet *) fCalibObjects->At(kLorentzLinear);
906
907   // Check
908   if((!calDetGain) || (!calDetVDrift) || (!fCH2d) || (!calDetExB) || (!CalChamberStatus)) return kFALSE;
909
910   // Gain
911   Double_t gainmean   = calDetGain->GetMean();
912   Double_t vdriftmean = calDetVDrift->GetMean();
913   Double_t exbmean    = calDetExB->GetMean();
914
915   Double_t gainrms    = calDetGain->GetRMSRobust();
916   Double_t vdriftrms  = calDetVDrift->GetRMSRobust();
917   Double_t exbrms     = calDetExB->GetRMSRobust();
918
919   //printf("Gain mean: %f, rms: %f\n",gainmean,gainrms);
920   //printf("Vdrift mean: %f, rms: %f\n",vdriftmean,vdriftrms);
921   //printf("ExB mean: %f, rms: %f\n",exbmean,exbrms);
922
923   // Check
924   if((TMath::Abs(gainrms) < 0.001) || (TMath::Abs(vdriftrms) < 0.001) || (TMath::Abs(exbrms) < 0.0000001)) return kFALSE;
925
926   // mask chambers with empty gain entries
927   //Int_t counter = 0;
928   for (Int_t idet = 0; idet < 540; idet++) {
929
930     // ch2d
931     TH1I *projch =  (TH1I *) fCH2d->ProjectionX("projch",idet+1,idet+1,(Option_t *)"e");
932     Double_t entries = projch->GetEntries();
933     //printf("Number of entries %f for det %d\n",entries,idet);
934
935     // gain
936     Double_t gain = calDetGain->GetValue(idet);
937
938     // vdrift
939     Double_t vdrift = calDetVDrift->GetValue(idet);
940
941     // exb
942     Double_t exb = calDetExB->GetValue(idet);
943
944
945     if( (entries<50 && !CalChamberStatus->IsNoData(idet))  ||
946         TMath::Abs(gainmean-gain) > (fRMSBadCalibratedGain*gainrms)          ||
947         TMath::Abs(vdriftmean-vdrift) > (fRMSBadCalibratedVdrift*vdriftrms)    ||
948         TMath::Abs(exbmean-exb) > (fRMSBadCalibratedExB*exbrms) ) {
949      
950       //printf(" chamber det %03d masked \n",idet);
951       //printf(" gainmean %f and gain %f, gainrms %f \n",gainmean,gain,gainrms);
952       //printf(" vdriftmean %f and vdrift %f, vdriftrms %f \n",vdriftmean,vdrift,vdriftrms);
953       //printf(" exbmean %f and exb %f, exbrms %f \n",exbmean,exb,exbrms);
954       
955       CalChamberStatus->SetStatus(idet,AliTRDCalChamberStatus::kBadCalibrated);
956       //counter++;
957     }
958
959     delete projch;
960     
961    }
962
963   // Security
964   for(Int_t sm=0; sm < 18; sm++) {
965     Int_t smnodata   = 0;
966     Int_t smbadcalib = 0;
967     for(Int_t det = 0; det < 30; det++){
968       Int_t detector = sm*30+det;
969       if(CalChamberStatus->IsNoData(detector)) smnodata++;
970       else {
971         if(CalChamberStatus->IsBadCalibrated(detector)) smbadcalib++;
972       }
973     }
974     fNoData[sm]  = smnodata;
975     fBadCalib[sm]= smbadcalib;
976     //printf("No Data %d, bad calibrated %d for %d\n",fNoData[sm],fBadCalib[sm],sm);
977   }
978   
979   // Security
980   //   for(Int_t sm=0; sm < 18; sm++) {
981   //     Int_t counter = 0;
982   //     for(Int_t det = 0; det < 30; det++){
983   //       Int_t detector = sm*30+det;
984   //       if(CalChamberStatus->IsBadCalibrated(detector)) counter++;
985   //     }
986   //     if(counter >= 20) {
987   //       for(Int_t det = 0; det < 30; det++){
988   //    Int_t detector = sm*30+det;
989   //    CalChamberStatus->SetStatus(detector,AliTRDCalChamberStatus::kGood);
990   //       }
991   //     }
992   //  }
993
994    fCalibObjects->AddAt(CalChamberStatus,kChamberStatus);
995    return kTRUE;
996
997  }
998
999
1000  //________________________________________________________________________________________________
1001  void AliTRDPreprocessorOffline::CorrectFromDetGainUsed() {
1002    //
1003    // Correct from the gas gain used afterwards
1004    //
1005    AliTRDCalDet *calDetGain = (AliTRDCalDet *) fCalibObjects->At(kGain);
1006    if(!calDetGain) return;
1007
1008    // Calculate mean
1009    Double_t mean = 0.0;
1010    Int_t nbdet = 0;
1011
1012    for(Int_t det = 0; det < 540; det++) {
1013
1014      Float_t gaininit = fCalDetGainUsed->GetValue(det);
1015      Float_t gainout = calDetGain->GetValue(det);
1016
1017
1018      if(TMath::Abs(gainout-1.0) > 0.000001) {
1019        mean += (gaininit*gainout);
1020        nbdet++;
1021      }  
1022    }
1023    if(nbdet > 0) mean = mean/nbdet;
1024
1025    for(Int_t det = 0; det < 540; det++) {
1026
1027      Float_t gaininit = fCalDetGainUsed->GetValue(det);
1028      Float_t gainout = calDetGain->GetValue(det);
1029
1030      if(TMath::Abs(gainout-1.0) > 0.000001) {
1031        Double_t newgain = gaininit*gainout;
1032        if(newgain < 0.1) newgain = 0.1;
1033        if(newgain > 1.9) newgain = 1.9;
1034        calDetGain->SetValue(det,newgain);
1035      }
1036      else {
1037        Double_t newgain = mean;
1038        if(newgain < 0.1) newgain = 0.1;
1039        if(newgain > 1.9) newgain = 1.9;
1040        calDetGain->SetValue(det,newgain);
1041      }
1042    }
1043
1044
1045  }
1046  //________________________________________________________________________________________________
1047  void AliTRDPreprocessorOffline::CorrectFromDetVdriftUsed() {
1048    //
1049    // Correct from the drift velocity
1050    //
1051
1052    //printf("Correct for vdrift\n");
1053
1054    AliTRDCalDet *calDetGain = (AliTRDCalDet *) fCalibObjects->At(kGain);
1055    if(!calDetGain) return;
1056
1057    Int_t detVdrift = kVdriftPHDet;
1058    if(fMethodSecond) detVdrift = kVdriftLinear;
1059    AliTRDCalDet *calDetVdrift = (AliTRDCalDet *) fCalibObjects->At(detVdrift);
1060    if(!calDetVdrift) return;
1061
1062    // Calculate mean
1063    if(!fNotEnoughStatisticsForTheVdriftLinear) {
1064      for(Int_t det = 0; det < 540; det++) {
1065        
1066        Float_t vdriftinit = fCalDetVdriftUsed->GetValue(det);
1067        Float_t vdriftout = calDetVdrift->GetValue(det);
1068        
1069        Float_t gain = calDetGain->GetValue(det);
1070        if(vdriftout > 0.0) gain = gain*vdriftinit/vdriftout;
1071        if(gain < 0.1) gain = 0.1;
1072        if(gain > 1.9) gain = 1.9;
1073        calDetGain->SetValue(det,gain);
1074      }
1075    }
1076    else {
1077      
1078      Float_t vdriftinit = fCalDetVdriftUsed->CalcMean(kFALSE);
1079      Float_t vdriftout = calDetVdrift->CalcMean(kFALSE);
1080      Float_t factorcorrectif = 1.0;
1081      if(vdriftout > 0.0) factorcorrectif = vdriftinit/vdriftout;
1082      for(Int_t det = 0; det < 540; det++) {
1083        Float_t gain = calDetGain->GetValue(det);
1084        gain = gain*factorcorrectif;
1085        if(gain < 0.1) gain = 0.1;
1086        if(gain > 1.9) gain = 1.9;
1087        calDetGain->SetValue(det,gain);
1088      }
1089      
1090    }
1091    
1092  }
1093 //_________________________________________________________________________________________________________________
1094  void AliTRDPreprocessorOffline::UpdateOCDBGain(Int_t startRunNumber, Int_t endRunNumber, const Char_t *storagePath){
1095    //
1096    // Update OCDB entry
1097    //
1098
1099    AliCDBMetaData *metaData= new AliCDBMetaData();
1100    metaData->SetObjectClassName("AliTRDCalDet");
1101    metaData->SetResponsible("Raphaelle Bailhache");
1102    metaData->AddDateToComment();
1103    metaData->SetBeamPeriod(1);
1104
1105    AliCDBId id1("TRD/Calib/ChamberGainFactor", startRunNumber, endRunNumber);
1106    AliCDBStorage * gStorage = AliCDBManager::Instance()->GetStorage(storagePath);
1107    AliTRDCalDet *calDet = (AliTRDCalDet *) fCalibObjects->At(kGain);
1108    if(calDet) gStorage->Put(calDet, id1, metaData);
1109
1110
1111  }
1112  //___________________________________________________________________________________________________________________
1113  void AliTRDPreprocessorOffline::UpdateOCDBExB(Int_t startRunNumber, Int_t endRunNumber, const Char_t *storagePath){
1114    //
1115    // Update OCDB entry
1116    //
1117
1118    Int_t detExB = kLorentzLinear;
1119    if(!fMethodSecond) return;
1120
1121    //printf("Pass\n");
1122
1123    AliCDBMetaData *metaData= new AliCDBMetaData();
1124    metaData->SetObjectClassName("AliTRDCalDet");
1125    metaData->SetResponsible("Raphaelle Bailhache");
1126    metaData->AddDateToComment();
1127    metaData->SetBeamPeriod(1);
1128
1129    AliCDBId id1("TRD/Calib/ChamberExB", startRunNumber, endRunNumber);
1130    AliCDBStorage * gStorage = AliCDBManager::Instance()->GetStorage(storagePath);
1131    AliTRDCalDet *calDet = (AliTRDCalDet *) fCalibObjects->At(detExB);
1132    if(calDet) gStorage->Put(calDet, id1, metaData);
1133    //if(!calDet) printf("No caldet\n");
1134
1135  }
1136 //___________________________________________________________________________________________________________________
1137 void AliTRDPreprocessorOffline::UpdateOCDBExBAlt(Int_t startRunNumber, Int_t endRunNumber, const Char_t *storagePath){
1138   //
1139   // Update OCDB entry
1140   //
1141
1142   Int_t detExB = kExbAlt;
1143   if(!fMethodSecond) return;
1144
1145   //printf("Pass\n");
1146
1147   AliCDBMetaData *metaData= new AliCDBMetaData();
1148   metaData->SetObjectClassName("AliTRDCalDet");
1149   metaData->SetResponsible("Theo Rascanu");
1150   metaData->AddDateToComment();
1151   metaData->SetBeamPeriod(1);
1152
1153   AliCDBId id1("TRD/Calib/ChamberExBAlt", startRunNumber, endRunNumber);
1154   AliCDBStorage * gStorage = AliCDBManager::Instance()->GetStorage(storagePath);
1155   AliTRDCalDet *calDet = (AliTRDCalDet *) fCalibObjects->At(detExB);
1156   if(calDet) gStorage->Put(calDet, id1, metaData);
1157   //if(!calDet) printf("No caldet\n");
1158
1159 }
1160  //___________________________________________________________________________________________________________________
1161  void AliTRDPreprocessorOffline::UpdateOCDBVdrift(Int_t startRunNumber, Int_t endRunNumber, const Char_t *storagePath){
1162    //
1163    // Update OCDB entry
1164    //
1165
1166    Int_t detVdrift = kVdriftPHDet;
1167
1168    if(fMethodSecond) detVdrift = kVdriftLinear;
1169
1170    AliCDBMetaData *metaData= new AliCDBMetaData();
1171    metaData->SetObjectClassName("AliTRDCalDet");
1172    metaData->SetResponsible("Raphaelle Bailhache");
1173    metaData->AddDateToComment();
1174    metaData->SetBeamPeriod(1);
1175
1176    AliCDBId id1("TRD/Calib/ChamberVdrift", startRunNumber, endRunNumber);
1177    AliCDBStorage * gStorage = AliCDBManager::Instance()->GetStorage(storagePath);
1178    AliTRDCalDet *calDet = (AliTRDCalDet *) fCalibObjects->At(detVdrift);
1179    if(calDet) gStorage->Put(calDet, id1, metaData);
1180
1181    //
1182
1183    if(!fMethodSecond) {
1184
1185      AliCDBMetaData *metaDataPad= new AliCDBMetaData();
1186      metaDataPad->SetObjectClassName("AliTRDCalPad");
1187      metaDataPad->SetResponsible("Raphaelle Bailhache");
1188      metaDataPad->AddDateToComment();
1189      metaDataPad->SetBeamPeriod(1);
1190
1191      AliCDBId id1Pad("TRD/Calib/LocalVdrift", startRunNumber, endRunNumber);
1192      AliTRDCalPad *calPad = (AliTRDCalPad *) fCalibObjects->At(kVdriftPHPad);
1193      if(calPad) gStorage->Put(calPad, id1Pad, metaDataPad);
1194
1195    }
1196
1197  }
1198  //________________________________________________________________________________________________________________________
1199  void AliTRDPreprocessorOffline::UpdateOCDBT0(Int_t startRunNumber, Int_t endRunNumber, const Char_t *storagePath){
1200    //
1201    // Update OCDB entry
1202    //
1203
1204    AliCDBMetaData *metaData= new AliCDBMetaData();
1205    metaData->SetObjectClassName("AliTRDCalDet");
1206    metaData->SetResponsible("Raphaelle Bailhache");
1207    metaData->AddDateToComment();
1208    metaData->SetBeamPeriod(1);
1209
1210    AliCDBId id1("TRD/Calib/ChamberT0", startRunNumber, endRunNumber);
1211    AliCDBStorage * gStorage = AliCDBManager::Instance()->GetStorage(storagePath);
1212    AliTRDCalDet *calDet = (AliTRDCalDet *) fCalibObjects->At(kT0PHDet);
1213    if(calDet) gStorage->Put(calDet, id1, metaData);
1214
1215    //
1216
1217    AliCDBMetaData *metaDataPad= new AliCDBMetaData();
1218    metaDataPad->SetObjectClassName("AliTRDCalPad");
1219    metaDataPad->SetResponsible("Raphaelle Bailhache");
1220    metaDataPad->AddDateToComment();
1221    metaDataPad->SetBeamPeriod(1);
1222
1223    AliCDBId id1Pad("TRD/Calib/LocalT0", startRunNumber, endRunNumber);
1224    AliTRDCalPad *calPad = (AliTRDCalPad *) fCalibObjects->At(kT0PHPad);
1225    if(calPad) gStorage->Put(calPad, id1Pad, metaDataPad);
1226
1227
1228
1229  }
1230  //_________________________________________________________________________________________________________________
1231  void AliTRDPreprocessorOffline::UpdateOCDBPRF(Int_t startRunNumber, Int_t endRunNumber, const Char_t *storagePath){
1232    //
1233    // Update OCDB entry
1234    //
1235
1236    AliCDBMetaData *metaData= new AliCDBMetaData();
1237    metaData->SetObjectClassName("AliTRDCalPad");
1238    metaData->SetResponsible("Raphaelle Bailhache");
1239    metaData->AddDateToComment();
1240    metaData->SetBeamPeriod(1);
1241
1242
1243    AliCDBId id1("TRD/Calib/PRFWidth", startRunNumber, endRunNumber);
1244    AliCDBStorage * gStorage = AliCDBManager::Instance()->GetStorage(storagePath);
1245    AliTRDCalPad *calPad = (AliTRDCalPad *) fCalibObjects->At(kPRF);
1246    if(calPad) gStorage->Put(calPad, id1, metaData);
1247
1248
1249  }
1250  //_________________________________________________________________________________________________________________
1251  void AliTRDPreprocessorOffline::UpdateOCDBChamberStatus(Int_t startRunNumber, Int_t endRunNumber, const Char_t *storagePath){
1252    //
1253    // Update OCDB entry
1254    //
1255
1256    AliCDBMetaData *metaData= new AliCDBMetaData();
1257    metaData->SetObjectClassName("AliTRDCalChamberStatus");
1258    metaData->SetResponsible("Raphaelle Bailhache and Julian Book");
1259    metaData->AddDateToComment();
1260    metaData->SetBeamPeriod(1);
1261
1262    AliCDBId id1("TRD/Calib/ChamberStatus", startRunNumber, endRunNumber);
1263    AliCDBStorage * gStorage = AliCDBManager::Instance()->GetStorage(storagePath);
1264    AliTRDCalChamberStatus *calChamberStatus = (AliTRDCalChamberStatus *) fCalibObjects->At(kChamberStatus);
1265    if(calChamberStatus) gStorage->Put(calChamberStatus, id1, metaData);
1266
1267
1268  }
1269  //__________________________________________________________________________________________________________________________
1270  Bool_t AliTRDPreprocessorOffline::ValidateGain() {
1271    //
1272    // Validate OCDB entry
1273    //
1274
1275    AliTRDCalDet *calDet = (AliTRDCalDet *) fCalibObjects->At(kGain);
1276    if(calDet) {
1277      Double_t mean = calDet->GetMean();
1278      Double_t rms = calDet->GetRMSRobust();
1279      if((mean > 0.2) && (mean < 1.4) && (rms < 0.5)) return kTRUE;
1280      //if((mean > 0.2) && (mean < 1.4)) return kTRUE;
1281      else {
1282        fStatusPos = fStatusPos | kGainErrorRange;
1283        return kFALSE;
1284      }
1285    }
1286    else return kFALSE;
1287    
1288
1289
1290  }
1291  //__________________________________________________________________________________________________________________________
1292  Bool_t AliTRDPreprocessorOffline::ValidateVdrift(){
1293    //
1294    // Update OCDB entry
1295    //
1296
1297    Int_t detVdrift = kVdriftPHDet;
1298    Bool_t ok = kTRUE;
1299
1300    if(fMethodSecond) detVdrift = kVdriftLinear;
1301
1302    AliTRDCalDet *calDet = (AliTRDCalDet *) fCalibObjects->At(detVdrift);
1303    if(calDet) {
1304      Double_t mean = calDet->GetMean();
1305      Double_t rms = calDet->GetRMSRobust();
1306      //printf("Vdrift::mean %f, rms %f\n",mean,rms);
1307      if(!((mean > 1.0) && (mean < 2.0) && (rms < 0.5))) {
1308        fStatusPos = fStatusPos | kVdriftErrorRange;
1309        ok = kFALSE;
1310      }
1311    }
1312    else return kFALSE; 
1313
1314    if(!fMethodSecond) {
1315      AliTRDCalPad *calPad = (AliTRDCalPad *) fCalibObjects->At(kVdriftPHPad);
1316      if(calPad) {
1317        Double_t mean = calPad->GetMean();
1318        Double_t rms = calPad->GetRMS();
1319        //printf("Vdrift::meanpad %f, rmspad %f\n",mean,rms);
1320        if(!((mean > 0.9) && (mean < 1.1) && (rms < 0.6))) {
1321          fStatusPos = fStatusPos | kVdriftErrorRange;
1322          ok = kFALSE;
1323        }
1324      }
1325      else return kFALSE;
1326    }
1327
1328    return ok;
1329
1330  }
1331  //__________________________________________________________________________________________________________________________
1332  Bool_t AliTRDPreprocessorOffline::ValidateExB(){
1333    //
1334    // Update OCDB entry
1335    //
1336
1337    AliTRDCalDet *calDet = (AliTRDCalDet *) fCalibObjects->At(kLorentzLinear);
1338    if(calDet) {
1339      Double_t mean = calDet->GetMean();
1340      Double_t rms = calDet->GetRMSRobust();
1341      //printf("Vdrift::mean %f, rms %f\n",mean,rms);
1342      if(!((mean > -1.0) && (mean < 1.0) && (rms < 0.5))) {
1343        fStatusNeg = fStatusNeg | kExBErrorRange;
1344        return kFALSE;
1345      }
1346      else return kTRUE;
1347    }
1348    else return kFALSE; 
1349    
1350  }
1351  //__________________________________________________________________________________________________________________________
1352  Bool_t AliTRDPreprocessorOffline::ValidateT0(){
1353    //
1354    // Update OCDB entry
1355    //
1356
1357    AliTRDCalDet *calDet = (AliTRDCalDet *) fCalibObjects->At(kT0PHDet);
1358    AliTRDCalPad *calPad = (AliTRDCalPad *) fCalibObjects->At(kT0PHPad);
1359    if(calDet && calPad) {
1360      Double_t meandet = calDet->GetMean();
1361      Double_t rmsdet = calDet->GetRMSRobust();
1362      Double_t meanpad = calPad->GetMean();
1363      //Double_t rmspad = calPad->GetRMS();
1364      //printf("T0::minimum %f, rmsdet %f,meanpad %f, rmspad %f\n",meandet,rmsdet,meanpad,rmspad);
1365      if((meandet > -1.5) && (meandet < 5.0) && (rmsdet < 4.0) && (meanpad < 5.0) && (meanpad > -0.5)) return kTRUE;
1366      else {
1367        fStatusPos = fStatusPos | kTimeOffsetErrorRange;
1368        return kFALSE;
1369      }
1370    }
1371    else return kFALSE;
1372
1373  }
1374  //__________________________________________________________________________________________________________________________
1375  Bool_t AliTRDPreprocessorOffline::ValidatePRF() const{
1376    //
1377    // Update OCDB entry
1378    //
1379
1380    AliTRDCalPad *calPad = (AliTRDCalPad *) fCalibObjects->At(kPRF);
1381    if(calPad) {
1382      Double_t meanpad = calPad->GetMean();
1383      Double_t rmspad = calPad->GetRMS();
1384      //printf("PRF::meanpad %f, rmspad %f\n",meanpad,rmspad);
1385      if((meanpad < 1.0) && (rmspad < 0.8)) return kTRUE;
1386      else return kFALSE;
1387    }
1388    else return kFALSE;
1389
1390
1391  }
1392  //__________________________________________________________________________________________________________________________
1393 Bool_t AliTRDPreprocessorOffline::ValidateChamberStatus(){
1394   //
1395   // Update OCDB entry
1396   //
1397   
1398   AliTRDCalChamberStatus *calChamberStatus = (AliTRDCalChamberStatus *) fCalibObjects->At(kChamberStatus);
1399   if(calChamberStatus) {
1400
1401     Int_t detectornodata   = 0;
1402     Int_t detectorbadcalib = 0;
1403     
1404     for(Int_t sm=0; sm < 18; sm++) {
1405       //printf("%d chambers w/o data in sm %d\n",fNoData[sm],sm);
1406       //printf("%d bad calibrated chambers in sm %d\n",fBadCalib[sm],sm);
1407       if(fNoData[sm] != 30) detectornodata += fNoData[sm];
1408       detectorbadcalib+=fBadCalib[sm];
1409     }
1410     //printf("Number of chambers w/o data %d\n",detectornodata);
1411     //printf("Number of chambers bad calibrated %d\n",detectorbadcalib);
1412
1413     if((detectornodata > fNoDataValidate) ||
1414        (detectorbadcalib > fBadCalibValidate)){
1415       fStatusPos = fStatusPos | kChamberStatusErrorRange;
1416       return kFALSE;
1417     }
1418     return kTRUE;
1419   }
1420   else return kFALSE;
1421   
1422 }
1423 //_____________________________________________________________________________
1424 Int_t AliTRDPreprocessorOffline::GetVersion(TString name) const
1425 {
1426   //
1427   // Get version from the title
1428   //
1429   
1430   // Some patterns
1431   const Char_t *version = "Ver";
1432   if(!strstr(name.Data(),version)) return -1;
1433   const Char_t *after = "Subver";  
1434   if(!strstr(name.Data(),after)) return -1;
1435
1436   for(Int_t ver = 0; ver < 999999999; ver++) {
1437
1438     TString vertry(version);
1439     vertry += ver;
1440     vertry += after;
1441
1442     //printf("vertry %s and name %s\n",vertry.Data(),name.Data());
1443
1444     if(strstr(name.Data(),vertry.Data())) return ver;
1445     
1446   }
1447   
1448   return -1;
1449
1450 }
1451
1452 //_____________________________________________________________________________
1453 Int_t AliTRDPreprocessorOffline::GetSubVersion(TString name) const
1454 {
1455   //
1456   // Get subversion from the title
1457   //
1458   
1459   // Some patterns
1460   const Char_t *subversion = "Subver";
1461   if(!strstr(name.Data(),subversion)) return -1;
1462   const Char_t *after = "FirstRun";
1463   if(!strstr(name.Data(),after)) {
1464     after = "Nz";
1465   }
1466   if(!strstr(name.Data(),after)) return -1;
1467
1468   
1469   for(Int_t ver = 0; ver < 999999999; ver++) {
1470     
1471     TString vertry(subversion);
1472     vertry += ver;
1473     vertry += after;
1474
1475     //printf("vertry %s and name %s\n",vertry.Data(),name.Data());
1476
1477     if(strstr(name.Data(),vertry.Data())) return ver;
1478     
1479   }
1480   
1481   return -1;
1482
1483 }
1484
1485 //_____________________________________________________________________________
1486 Int_t AliTRDPreprocessorOffline::GetFirstRun(TString name) const
1487 {
1488   //
1489   // Get first run from the title
1490   //
1491   
1492   // Some patterns
1493   const Char_t *firstrun = "FirstRun";
1494   if(!strstr(name.Data(),firstrun)) return -1;
1495   const Char_t *after = "Nz";  
1496   if(!strstr(name.Data(),after)) return -1;
1497   
1498   
1499   for(Int_t ver = 0; ver < 999999999; ver++) {
1500
1501     TString vertry(firstrun);
1502     vertry += ver;
1503     vertry += after;
1504
1505     //printf("vertry %s and name %s\n",vertry.Data(),name.Data());
1506
1507     if(strstr(name.Data(),vertry.Data())) return ver;
1508     
1509   }
1510   
1511   return -1;
1512
1513 }
1514 //_____________________________________________________________________________
1515 Bool_t AliTRDPreprocessorOffline::CheckStatus(Int_t status, Int_t bitMask) const
1516 {
1517   //
1518   // Checks the status
1519   //
1520
1521   return (status & bitMask) ? kTRUE : kFALSE;
1522   
1523 }
1524 //_____________________________________________________________________________
1525 Int_t AliTRDPreprocessorOffline::GetStatus() const
1526 {
1527   //
1528   // Checks the status
1529   // fStatusPos: errors
1530   // fStatusNeg: only info
1531   //
1532
1533   if(fStatusPos > 0) return fStatusPos;
1534   else return (-TMath::Abs(fStatusNeg));
1535   
1536 }
1537 //_____________________________________________________________________________
1538 void AliTRDPreprocessorOffline::PrintStatus() const
1539 {
1540   //
1541   // Do Summary
1542   //
1543
1544   AliInfo(Form("The error status is %d",fStatusPos));
1545   AliInfo(Form("IsGainErrorOld? %d",(Int_t)IsGainErrorOld()));
1546   AliInfo(Form("IsVdriftErrorOld? %d",(Int_t)IsVdriftErrorOld()));
1547   AliInfo(Form("IsGainErrorRange? %d",(Int_t)IsGainErrorRange()));
1548   AliInfo(Form("IsVdriftErrorRange? %d",(Int_t)IsVdriftErrorRange()));
1549   AliInfo(Form("IsTimeOffsetErrorRange? %d",(Int_t)IsTimeOffsetErrorRange()));
1550   AliInfo(Form("IsChamberStatusErrorRange? %d",(Int_t)IsChamberStatusErrorRange()));
1551
1552  
1553   AliInfo(Form("The info status is %d",fStatusNeg));
1554   AliInfo(Form("IsGainNotEnoughStatsButFill? %d",(Int_t)IsGainNotEnoughStatsButFill()));
1555   AliInfo(Form("IsVdriftNotEnoughStatsButFill? %d",(Int_t)IsVdriftNotEnoughStatsButFill()));
1556   AliInfo(Form("IsGainNotEnoughStatsNotFill? %d",(Int_t)IsGainNotEnoughStatsNotFill()));
1557   AliInfo(Form("IsVdriftNotEnoughStatsNotFill? %d",(Int_t)IsVdriftNotEnoughStatsNotFill()));
1558   AliInfo(Form("IsTimeOffsetNotEnoughStatsNotFill? %d",(Int_t)IsTimeOffsetNotEnoughStatsNotFill()));
1559
1560   AliInfo(Form("IsExBErrorRange? %d",(Int_t)IsExBErrorRange()));
1561   AliInfo(Form("IsExBErrorOld? %d",(Int_t)IsExBErrorOld()));
1562   
1563 }
1564 //___________________________________________________________________________________
1565 void AliTRDPreprocessorOffline::SetCalDetVdrift(AliTRDCalDet *calDetVdriftUsed) 
1566 {
1567
1568   fCalDetVdriftUsed = calDetVdriftUsed;
1569
1570   fCalDetExBUsed = new AliTRDCalDet("lorentz angle tan","lorentz angle tan (detector value)");
1571   for(Int_t k = 0; k < 540; k++){
1572     fCalDetExBUsed->SetValue(k,AliTRDCommonParam::Instance()->GetOmegaTau(fCalDetVdriftUsed->GetValue(k)));
1573     //printf("Set the exb object for detector %d, vdrift %f and exb %f\n",k,fCalDetVdriftUsed->GetValue(k),fCalDetExBUsed->GetValue(k));
1574   }
1575   
1576 };
1577 //___________________________________________________________________________________
1578 Bool_t AliTRDPreprocessorOffline::SetCalDetGain(Int_t runNumber, Int_t version, Int_t subversion) 
1579 {
1580   //
1581   // Set the fCalDetGainUsed
1582   //
1583
1584   if((version == 0) && (subversion == 0)) return kFALSE;
1585
1586   AliCDBEntry *entry = AliCDBManager::Instance()->Get("TRD/Calib/ChamberGainFactor",runNumber, version, subversion);
1587   if(!entry) {
1588     AliError("Found no entry\n");
1589     fStatusPos = fStatusPos | kGainErrorOld;
1590     return kFALSE;
1591   }
1592   //const AliCDBId id = entry->GetId();
1593   //version = id.GetVersion();
1594   //subversion = id.GetSubVersion();
1595   //printf("Found version %d and subversion %d for vdrift\n",version,subversion);
1596   AliTRDCalDet* calDet = (AliTRDCalDet *)entry->GetObject();
1597   if(calDet) fCalDetGainUsed = calDet;
1598   else {
1599     fStatusPos = fStatusPos | kGainErrorOld;
1600     return kFALSE;
1601   }
1602   
1603   return kTRUE;
1604
1605 }
1606 //___________________________________________________________________________________
1607 Bool_t AliTRDPreprocessorOffline::SetCalDetVdriftExB(Int_t runNumber, Int_t versionv, Int_t subversionv, Int_t versionexb, Int_t subversionexb) 
1608 {
1609   //
1610   // Set the fCalDetVdriftUsed and fCalDetExBUsed
1611   //
1612
1613   if((versionv == 0) && (subversionv == 0)) return kFALSE;
1614
1615   AliCDBEntry *entry = AliCDBManager::Instance()->Get("TRD/Calib/ChamberVdrift",runNumber, versionv, subversionv);
1616   if(!entry) {
1617     AliError("Found no entry\n");
1618     fStatusPos = fStatusPos | kVdriftErrorOld;
1619     return kFALSE;
1620   }
1621   AliTRDCalDet* calDet = (AliTRDCalDet *)entry->GetObject();
1622   if(calDet) fCalDetVdriftUsed = calDet;
1623   else {
1624     fStatusPos = fStatusPos | kVdriftErrorOld;
1625     return kFALSE;
1626   }
1627
1628   // ExB object
1629
1630   if((versionexb == 0) && (subversionexb == 0)) {
1631     
1632     fCalDetExBUsed = new AliTRDCalDet("lorentz angle tan","lorentz angle tan (detector value)");
1633     for(Int_t k = 0; k < 540; k++){
1634       fCalDetExBUsed->SetValue(k,AliTRDCommonParam::Instance()->GetOmegaTau(fCalDetVdriftUsed->GetValue(k)));
1635       //printf("Nothing found: set the exb object for detector %d, vdrift %f and exb %f\n",k,fCalDetVdriftUsed->GetValue(k),fCalDetExBUsed->GetValue(k));
1636     }
1637   }
1638   else {
1639
1640     entry = 0x0;
1641     entry = AliCDBManager::Instance()->Get("TRD/Calib/ChamberExB",runNumber, versionexb, subversionexb);
1642     if(!entry) {
1643       //printf("Found no entry\n");
1644       fStatusPos = fStatusPos | kExBErrorOld;   
1645       return kFALSE;
1646     }
1647     AliTRDCalDet* calDetexb = (AliTRDCalDet *)entry->GetObject();
1648     if(!calDetexb) {
1649       fStatusPos = fStatusPos | kExBErrorOld;   
1650       return kFALSE;
1651     }
1652     
1653     Double_t meanexb = calDetexb->GetMean();
1654     //printf("Mean value %f\n",meanexb);
1655     if((meanexb > 70) || (fNoExBUsedInReco)) {
1656       fCalDetExBUsed = new AliTRDCalDet("lorentz angle tan","lorentz angle tan (detector value)");
1657       for(Int_t k = 0; k < 540; k++){
1658         fCalDetExBUsed->SetValue(k,AliTRDCommonParam::Instance()->GetOmegaTau(fCalDetVdriftUsed->GetValue(k)));
1659         //printf("Found but: set the exb object for detector %d, vdrift %f and exb %f\n",k,fCalDetVdriftUsed->GetValue(k),fCalDetExBUsed->GetValue(k));
1660       }
1661     }
1662     else {
1663       fCalDetExBUsed = calDetexb;
1664     }
1665     
1666   }
1667
1668   
1669   return kTRUE;
1670
1671 }
1672