cleanup of HLTTRDCalibration component, fixing bug in calibration: dQdl (Theodor)
[u/mrichter/AliRoot.git] / HLT / TRD / AliHLTTRDCalibrationComponent.cxx
1 // $Id$
2
3 /**************************************************************************
4  * Copyright(c) 1998-1999, ALICE Experiment at CERN, All rights reserved. *
5  *                                                                        *
6  * Authors: Matthias Richter <Matthias.Richter@ift.uib.no>                *
7  *          Timm Steinbeck <timm@kip.uni-heidelberg.de>                   *
8  *          for The ALICE Off-line Project.                               *
9  *                                                                        *
10  * Permission to use, copy, modify and distribute this software and its   *
11  * documentation strictly for non-commercial purposes is hereby granted   *
12  * without fee, provided that the above copyright notice appears in all   *
13  * copies and that both the copyright notice and this permission notice   *
14  * appear in the supporting documentation. The authors make no claims     *
15  * about the suitability of this software for any purpose. It is          *
16  * provided "as is" without express or implied warranty.                  *
17  **************************************************************************/
18
19 /** @file   AliHLTTRDCalibrationComponent.cxx
20     @author Timm Steinbeck, Matthias Richter
21     @date
22     @brief  A TRDCalibration processing component for the HLT. */
23
24 #if __GNUC__ >= 3
25 using namespace std;
26 #endif
27
28 #include "TTree.h"
29 #include "TFile.h"
30 #include "TBranch.h"
31 #include "TH2I.h"
32 #include "TH2.h"
33 #include "TProfile2D.h"
34
35 #include "AliHLTReadoutList.h"
36
37 #include "AliHLTTRDCalibrationComponent.h"
38 #include "AliHLTTRDDefinitions.h"
39 #include "AliHLTTRDUtils.h"
40
41 #include "AliCDBManager.h"
42 #include "AliCDBStorage.h"
43 #include "AliRawReaderMemory.h"
44
45 #include "AliTRDCalPad.h"
46 #include "AliTRDCalDet.h"
47
48 #include "AliTRDCalibraFillHisto.h"
49 #include "AliTRDtrackV1.h"
50
51 #include "AliTRDCalibraFit.h"
52 #include "AliTRDCalibraMode.h"
53 #include "AliTRDCalibraVector.h"
54 #include "AliTRDCalibraVdriftLinearFit.h"
55 #include "AliTRDReconstructor.h"
56 #include "AliTRDrecoParam.h"
57
58 #include <cstdlib>
59 #include <cerrno>
60 #include <string>
61
62 ClassImp(AliHLTTRDCalibrationComponent);
63
64 AliHLTTRDCalibrationComponent::AliHLTTRDCalibrationComponent()
65   : AliHLTCalibrationProcessor(),
66     fTRDCalibraFillHisto(NULL),
67     fOutputSize(500000),
68     fTracksArray(NULL),
69     fOutArray(NULL),
70     fAfterRunArray(NULL),
71     fDisplayArray(NULL),
72     fSavedTimeBins(kFALSE),
73     fTrgStrings(NULL),
74     fAccRejTrg(0)
75 {
76   // Default constructor
77 }
78
79 AliHLTTRDCalibrationComponent::~AliHLTTRDCalibrationComponent()
80 {
81   // Destructor
82 }
83
84 const char* AliHLTTRDCalibrationComponent::GetComponentID()
85 {
86   // Return the component ID const char *
87   return "TRDCalibration"; // The ID of this component
88 }
89
90 void AliHLTTRDCalibrationComponent::GetInputDataTypes( vector<AliHLTComponentDataType>& list)
91 {
92   // Get the list of input data
93   list.clear(); // We do not have any requirements for our input data type(s).
94   list.push_back(AliHLTTRDDefinitions::fgkTracksDataType);
95 }
96
97 AliHLTComponentDataType AliHLTTRDCalibrationComponent::GetOutputDataType()
98 {
99   // Get the output data type
100   return kAliHLTMultipleDataType;
101   //  return AliHLTTRDDefinitions::fgkCalibrationDataType;
102  
103 }
104
105 int AliHLTTRDCalibrationComponent::GetOutputDataTypes(AliHLTComponentDataTypeList& tgtList)
106 {
107   // Get the output data type
108   tgtList.clear();
109   tgtList.push_back(AliHLTTRDDefinitions::fgkCalibrationDataType);
110   tgtList.push_back(AliHLTTRDDefinitions::fgkEORCalibrationDataType);
111   return tgtList.size();
112 }
113
114 void AliHLTTRDCalibrationComponent::GetOutputDataSize( unsigned long& constBase, double& inputMultiplier )
115 {
116   // Get the output data size
117   constBase = fOutputSize;
118   inputMultiplier = 0;
119 }
120
121 AliHLTComponent* AliHLTTRDCalibrationComponent::Spawn()
122 {
123   // Spawn function, return new instance of this class
124   return new AliHLTTRDCalibrationComponent;
125 };
126
127 Int_t AliHLTTRDCalibrationComponent::ScanArgument( int argc, const char** argv )
128 {
129   // perform initialization. We check whether our relative output size is specified in the arguments.
130   int i = 0;
131   char* cpErr;
132   if(!fTrgStrings)
133     fTrgStrings = new TObjArray();
134
135   while ( i < argc )
136     {
137       HLTDebug("argv[%d] == %s", i, argv[i] );
138       if ( !strcmp( argv[i], "output_size" ) )
139         {
140           if ( i+1>=argc )
141             {
142               HLTError("Missing output_size parameter");
143               return ENOTSUP;
144             }
145           HLTDebug("argv[%d+1] == %s", i, argv[i+1] );
146           fOutputSize = strtoul( argv[i+1], &cpErr, 0 );
147           if ( *cpErr )
148             {
149               HLTError("Cannot convert output_size parameter '%s'", argv[i+1] );
150               return EINVAL;
151             }
152           HLTInfo("Output size set to %lu %%", fOutputSize );
153           i += 2;
154           continue;
155         }
156       if ( !strcmp( argv[i], "-everyNevent" ) )
157         {
158           if ( i+1>=argc )
159             {
160               HLTError("Missing everyNevent parameter");
161               return ENOTSUP;
162             }
163           HLTDebug("argv[%d+1] == %s", i, argv[i+1] );
164           HLTInfo("Option -everyNevent depreceated");
165           i += 2;
166           continue;
167         }
168       if ( !strcmp( argv[i], "-TrgStr" ) )
169         {
170           if ( i+1>=argc )
171             {
172               HLTError("Missing parameter for mbTriggerString");
173               return ENOTSUP;
174             }
175           HLTDebug("argv[%d+1] == %s", i, argv[i+1] );
176           fTrgStrings->Add(new TObjString(argv[i+1]));
177           i += 2;
178           continue;
179         }
180
181       if ( !strcmp( argv[i], "-acceptTrgStr" ) )
182         {
183           fAccRejTrg=1;
184           i += 1;
185           continue;
186         }
187       if ( !strcmp( argv[i], "-rejectTrgStr" ) )
188         {
189           fAccRejTrg=-1;
190           i += 1;
191           continue;
192         }
193
194       else {
195         HLTError("Unknown option '%s'", argv[i] );
196         return EINVAL;
197       }
198     }
199   return i;
200 }
201
202 Int_t AliHLTTRDCalibrationComponent::InitCalibration()
203 {
204
205   if(!fTrgStrings)
206     fTrgStrings = new TObjArray();
207
208   if(!AliCDBManager::Instance()->IsDefaultStorageSet()){
209     HLTError("DefaultStorage is not set in CDBManager");
210     return -EINVAL;
211   }
212   if(AliCDBManager::Instance()->GetRun()<0){
213     HLTError("Run Number is not set in CDBManager");
214     return -EINVAL;
215   }
216   HLTInfo("CDB default storage: %s; RunNo: %i", (AliCDBManager::Instance()->GetDefaultStorage()->GetBaseFolder()).Data(), AliCDBManager::Instance()->GetRun());
217
218   if(fTrgStrings->GetEntriesFast()>0 && !fAccRejTrg){
219     HLTError("Trigger string(s) given, but acceptTrgStr or rejectTrgStr not selected");
220     return -EINVAL;
221   }
222
223   fTRDCalibraFillHisto = AliTRDCalibraFillHisto::Instance();
224   fTRDCalibraFillHisto->SetIsHLT(kTRUE);
225   fTRDCalibraFillHisto->SetHisto2d(); // choose to use histograms
226   fTRDCalibraFillHisto->SetCH2dOn();  // choose to calibrate the gain
227   fTRDCalibraFillHisto->SetPH2dOn();  // choose to calibrate the drift velocity
228   fTRDCalibraFillHisto->SetPRF2dOn(); // choose to look at the PRF
229   fTRDCalibraFillHisto->SetIsHLT(); // per detector
230   //fTRDCalibraFillHisto->SetDebugLevel(1);// debug
231   fTRDCalibraFillHisto->SetFillWithZero(kTRUE);
232   fTRDCalibraFillHisto->SetLinearFitterOn(kTRUE);
233   fTRDCalibraFillHisto->SetNumberBinCharge(100);
234   
235   fTracksArray = new TClonesArray("AliTRDtrackV1");
236   fOutArray = new TObjArray(4);
237   fAfterRunArray=new TObjArray(5);
238   fDisplayArray=new TObjArray(4);
239
240   HLTDebug("run SetupCTPData");
241   SetupCTPData();
242
243   return 0;
244 }
245
246 Int_t AliHLTTRDCalibrationComponent::DeinitCalibration()
247 {
248   
249   // Deinitialization of the component
250   
251   HLTDebug("DeinitCalibration");
252   delete fTracksArray; fTracksArray=0;
253   //fTRDCalibraFillHisto->Destroy();
254   //fOutArray->Delete();
255   delete fOutArray; fOutArray=0;
256   fAfterRunArray->Delete();
257   delete fAfterRunArray; fAfterRunArray=0;
258   fDisplayArray->Delete();
259   delete fDisplayArray; fDisplayArray=0;
260   fTrgStrings->Delete();
261   delete fTrgStrings; fTrgStrings=0;
262   return 0;
263 }
264
265 Int_t AliHLTTRDCalibrationComponent::ProcessCalibration(const AliHLTComponent_EventData& evtData,
266                                                         const AliHLTComponent_BlockData* /*blocks*/,
267                                                         AliHLTComponent_TriggerData& trigData,
268                                                         AliHLTUInt8_t* /*outputPtr*/,
269                                                         AliHLTUInt32_t& /*size*/,
270                                                         vector<AliHLTComponent_BlockData>& /*outputBlocks*/)
271 {
272   HLTDebug("NofBlocks %lu", evtData.fBlockCnt );
273   // Process an event
274         
275   TClonesArray* TCAarray[18] = {0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0};
276   Int_t usedEntries = 0;
277   Int_t blockOrObject = 0;
278   Int_t nTimeBins = -1;
279
280   for (const AliHLTComponentBlockData* pBlock=GetFirstInputBlock(AliHLTTRDDefinitions::fgkTracksDataType); pBlock; pBlock=GetNextInputBlock()) 
281     {
282       TCAarray[0] = fTracksArray;
283       AliHLTTRDUtils::ReadTracks(TCAarray[0], pBlock->fPtr, pBlock->fSize, &nTimeBins);
284       usedEntries = 1;
285       blockOrObject = -1;
286     }  
287
288   for(const TObject *iter = GetFirstInputObject(AliHLTTRDDefinitions::fgkHiLvlTracksDataType); iter; iter = GetNextInputObject()) 
289     {
290       if(blockOrObject<0){
291         HLTError("You may not mix high level and low level!");
292         return -1;
293       }
294
295       TCAarray[usedEntries] = dynamic_cast<TClonesArray*>(const_cast<TObject*>(iter));
296       if(!TCAarray[usedEntries])continue;
297       TObjString* strg = dynamic_cast<TObjString*>(const_cast<TObject*>(GetNextInputObject()));
298       if(!strg)continue;
299
300       nTimeBins = strg->String().Atoi();
301       usedEntries++;
302       blockOrObject = 1;
303     }
304
305   if(!blockOrObject)
306     return 0;
307
308   if(!fSavedTimeBins){
309     if(nTimeBins<0){
310       HLTFatal("Number of timebins is negative!");
311       return -1;
312     }
313     HLTDebug("Saving number of time bins which was read from input block. Value is: %d", nTimeBins);
314     fTRDCalibraFillHisto->Init2Dhistos(nTimeBins); // initialise the histos
315     fTRDCalibraFillHisto->SetNumberClusters(0); // At least 1 clusters
316     fTRDCalibraFillHisto->SetNumberClustersf(nTimeBins); // Not more than %d  clusters
317     fSavedTimeBins=kTRUE;
318   }
319
320   Bool_t TriggerPassed=kFALSE;
321
322   if(fAccRejTrg){
323     if(fAccRejTrg>0){
324       TriggerPassed=kFALSE;
325       for(int i = 0; i < fTrgStrings->GetEntriesFast(); i++){
326         const TObjString *const obString=(TObjString*)fTrgStrings->At(i);
327         const TString tString=obString->GetString();
328         //printf("Trigger Output: %i\n",EvaluateCTPTriggerClass(tString.Data(),trigData));
329         if(EvaluateCTPTriggerClass(tString.Data(),trigData)){TriggerPassed=kTRUE; break;}
330       }
331     }
332     else{
333       TriggerPassed=kTRUE;
334       for(int i = 0; i < fTrgStrings->GetEntriesFast(); i++){
335         const TObjString *const obString=(TObjString*)fTrgStrings->At(i);
336         const TString tString=obString->GetString();
337         if(EvaluateCTPTriggerClass(tString.Data(),trigData)){TriggerPassed=kFALSE; break;}
338       }
339     }
340   }
341   
342   fTRDCalibraFillHisto->SetCH2dOn(TriggerPassed);
343   for(int i=0; i<usedEntries; i++){
344     const TClonesArray* inArr = TCAarray[i];
345     Int_t nbEntries = inArr->GetEntries();
346     HLTDebug(" %i TRDtracks in tracksArray", nbEntries);
347     AliTRDtrackV1* trdTrack = 0x0;
348     for (Int_t i = 0; i < nbEntries; i++){
349       HLTDebug("%i/%i: ", i+1, nbEntries);
350       trdTrack = (AliTRDtrackV1*)inArr->At(i);
351       // for(int i=0; i<7; i++)
352       //   if(trdTrack->GetTracklet(i))trdTrack->GetTracklet(i)->Bootstrap(fReconstructor);
353       fTRDCalibraFillHisto->UpdateHistogramsV1(trdTrack);
354     }
355   }
356
357   if(!fOutArray->At(0))FormOutput(0);
358   if(!fDisplayArray->At(0))FormOutput(1);
359   PushBack(fDisplayArray, AliHLTTRDDefinitions::fgkCalibrationDataType);
360
361   if(blockOrObject<0){
362     TCAarray[0]->Delete();
363   }
364
365   return 0;
366
367 }
368
369 /**
370  * Form output array of histrograms
371  */
372 //============================================================================
373 void AliHLTTRDCalibrationComponent::FormOutput(Int_t param)
374 {
375   // gain histo
376   TH2I *hCH2d = fTRDCalibraFillHisto->GetCH2d();
377   if(!param)fOutArray->Add(hCH2d);
378   else fDisplayArray->Add(hCH2d);
379
380   // drift velocity histo
381   TProfile2D *hPH2d = fTRDCalibraFillHisto->GetPH2d();
382   if(!param)fOutArray->Add(hPH2d);
383   else fDisplayArray->Add(hPH2d);
384
385   // PRF histo
386   TProfile2D *hPRF2d = fTRDCalibraFillHisto->GetPRF2d();
387   if(!param)fOutArray->Add(hPRF2d);
388   else fDisplayArray->Add(hPRF2d);
389
390   // Vdrift Linear Fit
391   if(!param){
392     AliTRDCalibraVdriftLinearFit *hVdriftLinearFitOne=(AliTRDCalibraVdriftLinearFit *)fTRDCalibraFillHisto->GetVdriftLinearFit();
393     fOutArray->Add(hVdriftLinearFitOne);
394   }
395   else{
396     TH2S *hVdriftLinearFitOne = (TH2S *)(((AliTRDCalibraVdriftLinearFit *)fTRDCalibraFillHisto->GetVdriftLinearFit())->GetLinearFitterHisto(10,kTRUE)); 
397     fDisplayArray->Add(hVdriftLinearFitOne);
398   }
399
400   HLTDebug("GetCH2d = 0x%x; NEntries = %i; size = %i", hCH2d, hCH2d->GetEntries(), sizeof(*hCH2d));
401   hCH2d->Print();
402   HLTDebug("GetPH2d = 0x%x; NEntries = %i; size = %i", hPH2d, hPH2d->GetEntries(), sizeof(*hPH2d));
403   hPH2d->Print();
404   HLTDebug("GetPRF2d = 0x%x; NEntries = %i; size = %i", hPRF2d, hPRF2d->GetEntries(), sizeof(*hPRF2d));
405   hPRF2d->Print();
406   //HLTDebug("GetVdriftLinearFit = 0x%x; size = %i", hVdriftLinearFitOne, sizeof(hVdriftLinearFitOne)); 
407   
408   HLTDebug("output Array: pointer = 0x%x; NEntries = %i; size = %i", fOutArray, fOutArray->GetEntries(), sizeof(fOutArray));
409    
410 }
411
412 Int_t AliHLTTRDCalibrationComponent::ShipDataToFXS(const AliHLTComponentEventData& /*evtData*/, AliHLTComponentTriggerData& /*trigData*/)
413 {
414   //fTRDCalibraFillHisto->DestroyDebugStreamer();
415
416   AliHLTReadoutList rdList(AliHLTReadoutList::kTRD);
417
418   EORCalibration();
419   
420   fOutArray->Remove(fOutArray->FindObject("AliTRDCalibraVdriftLinearFit"));
421   //fOutArray->Remove(fOutArray->FindObject("PRF2d"));
422   //fOutArray->Remove(fOutArray->FindObject("PH2d"));
423   //fOutArray->Remove(fOutArray->FindObject("CH2d"));
424
425   //if(!(fOutArray->FindObject("CH2d"))) {
426   //  TH2I * ch2d = new TH2I("CH2d","Nz0Nrphi0",100,0.0,300.0,540,0,540);
427   //  fOutArray->Add(ch2d);
428   //}
429
430   //if(!(fOutArray->FindObject("PH2d"))) {
431   //  TProfile2D * ph2d = new TProfile2D("PH2d","Nz0Nrphi0",30,-0.05,2.95,540,0,540);
432   //  fOutArray->Add(ph2d);
433   //}
434
435   //if(!(fOutArray->FindObject("PRF2d"))) {
436   //  TProfile2D * prf2d = new TProfile2D("PRF2d","Nz0Nrphi0Ngp3",60,-9.0,9.0,540,0,540);
437   //  fOutArray->Add(prf2d);
438   //}
439
440
441   HLTDebug("Size of the fOutArray is %d\n",fOutArray->GetEntriesFast());
442
443   /*
444   TString fileName="$ALIHLT_TOPDIR/build-debug/output/CalibHistoDump_run";
445   fileName+=".root";
446   HLTInfo("Dumping Histogram file to %s",fileName.Data());
447   TFile* file = TFile::Open(fileName, "RECREATE");
448   //fAfterRunArray->Write();
449   fOutArray->Write();
450   file->Close();
451   HLTInfo("Histogram file dumped");
452   */
453   
454   PushToFXS((TObject*)fOutArray, "TRD", "GAINDRIFTPRF", rdList.Buffer() );
455   //PushToFXS((TObject*)fOutArray->FindObject("CH2d"), "TRD", "GAINDRIFTPRF", rdList.Buffer() );
456
457         
458   return 0;
459 }
460 Int_t AliHLTTRDCalibrationComponent::EORCalibration()
461 {
462   //Also Fill histograms for the online display
463   TH2I *hCH2d=(TH2I*)fOutArray->FindObject("CH2d");
464   TProfile2D *hPH2d=(TProfile2D*)fOutArray->FindObject("PH2d");
465   TProfile2D *hPRF2d= (TProfile2D*)fOutArray->FindObject("PRF2d");
466   AliTRDCalibraVdriftLinearFit* hVdriftLinearFit = (AliTRDCalibraVdriftLinearFit*)fOutArray->FindObject("AliTRDCalibraVdriftLinearFit");
467  
468
469   if(!hCH2d || !hPH2d || !hPRF2d || !hVdriftLinearFit) return 0; 
470
471   //Fit
472   AliTRDCalibraFit *calibra = AliTRDCalibraFit::Instance();
473
474   //Gain
475   calibra->SetMinEntries(100);
476   calibra->AnalyseCH(hCH2d);
477   Int_t nbtg = 6*4*18*((Int_t) ((AliTRDCalibraMode *)calibra->GetCalibraMode())->GetDetChamb0(0))
478     + 6*  18*((Int_t) ((AliTRDCalibraMode *)calibra->GetCalibraMode())->GetDetChamb2(0));
479   Int_t nbfit       = calibra->GetNumberFit();
480   Int_t nbE         = calibra->GetNumberEnt();
481   TH1F *coefgain = 0x0;
482   // enough statistics
483   //if ((nbtg >                  0) && 
484   //   (nbfit        >= 0.2*nbE)) {
485   // create the cal objects
486   //calibra->PutMeanValueOtherVectorFit(1,kTRUE);
487   TObjArray object           = calibra->GetVectorFit();
488   AliTRDCalDet *objgaindet   = calibra->CreateDetObjectGain(&object,kFALSE);
489   coefgain                   = objgaindet->MakeHisto1DAsFunctionOfDet();
490   //}
491   calibra->ResetVectorFit();
492
493   // vdrift second method
494   calibra->SetMinEntries(100); // If there is less than 100
495   hVdriftLinearFit->FillPEArray();
496   calibra->AnalyseLinearFitters(hVdriftLinearFit);
497   nbtg = 540;
498   nbfit = calibra->GetNumberFit();
499   nbE   = calibra->GetNumberEnt();
500   TH1F *coefdriftsecond = 0x0;
501   // enough statistics
502   //if ((nbtg >                  0) && 
503   // (nbfit        >= 0.1*nbE)) {
504   // create the cal objects
505   //calibra->PutMeanValueOtherVectorFit(1,kTRUE);
506   object  = calibra->GetVectorFit();
507   AliTRDCalDet *objdriftvelocitydetsecond = calibra->CreateDetObjectVdrift(&object,kTRUE);
508   objdriftvelocitydetsecond->SetTitle("secondmethodvdrift");
509   coefdriftsecond  = objdriftvelocitydetsecond->MakeHisto1DAsFunctionOfDet();
510   //}
511   calibra->ResetVectorFit();
512   
513   // vdrift first method
514   calibra->SetMinEntries(100*20); // If there is less than 20000
515   calibra->AnalysePH(hPH2d);
516   nbtg = 6*4*18*((Int_t) ((AliTRDCalibraMode *)calibra->GetCalibraMode())->GetDetChamb0(1))
517     + 6*  18*((Int_t) ((AliTRDCalibraMode *)calibra->GetCalibraMode())->GetDetChamb2(1));
518   nbfit        = calibra->GetNumberFit();
519   nbE          = calibra->GetNumberEnt();
520   TH1F *coefdrift = 0x0;
521   TH1F *coeft0 = 0x0;
522   // enough statistics
523   //if ((nbtg >                  0) && 
524   // (nbfit        >= 0.2*nbE)) {
525   // create the cal objects
526   //calibra->PutMeanValueOtherVectorFit(1,kTRUE);
527   //calibra->PutMeanValueOtherVectorFit2(1,kTRUE);
528   object  = calibra->GetVectorFit();
529   AliTRDCalDet *objdriftvelocitydet = calibra->CreateDetObjectVdrift(&object,kTRUE);
530   coefdrift        = objdriftvelocitydet->MakeHisto1DAsFunctionOfDet();
531   object              = calibra->GetVectorFit2();
532   AliTRDCalDet *objtime0det  = calibra->CreateDetObjectT0(&object,kTRUE);
533   coeft0        = objtime0det->MakeHisto1DAsFunctionOfDet();
534   //}
535   calibra->ResetVectorFit();
536            
537
538   //PRF
539   calibra->SetMinEntries(200); 
540   calibra->AnalysePRFMarianFit(hPRF2d);
541   nbtg = 6*4*18*((Int_t) ((AliTRDCalibraMode *)calibra->GetCalibraMode())->GetDetChamb0(2))
542     + 6*  18*((Int_t) ((AliTRDCalibraMode *)calibra->GetCalibraMode())->GetDetChamb2(2));
543   nbfit        = calibra->GetNumberFit();
544   nbE          = calibra->GetNumberEnt();
545   TH1F *coefprf = 0x0;
546   // enough statistics
547   //if ((nbtg >                  0) && 
548   //  (nbfit        >= 0.95*nbE)) {
549   // create cal pad objects 
550   object            = calibra->GetVectorFit();
551   TObject *objPRFpad          = calibra->CreatePadObjectPRF(&object);
552   coefprf                     = ((AliTRDCalPad *) objPRFpad)->MakeHisto1D();
553   //}
554   calibra->ResetVectorFit();
555
556
557   coefgain->SetName("coefgain");
558   coefprf->SetName("coefprf");
559   coefdrift->SetName("coefdrift");
560   coefdriftsecond->SetName("coefdriftsecond");
561   coeft0->SetName("coeft0");
562   if(coefgain) fAfterRunArray->Add(coefgain);
563   if(coefprf) fAfterRunArray->Add(coefprf);
564   if(coefdrift) fAfterRunArray->Add(coefdrift);
565   if(coefdriftsecond) fAfterRunArray->Add(coefdriftsecond);
566   if(coeft0) fAfterRunArray->Add(coeft0);
567   
568
569   if(coefgain||coefprf||coefdrift||coeft0||coefdriftsecond) {
570     PushBack(fAfterRunArray, AliHLTTRDDefinitions::fgkEORCalibrationDataType);
571   }
572
573   /*
574   TString fileName="$ALIHLT_TOPDIR/build-debug/output/CalibHistoDump_run";
575   fileName+=".root";
576   HLTInfo("Dumping Histogram file to %s",fileName.Data());
577   TFile* file = TFile::Open(fileName, "RECREATE");
578   fAfterRunArray->Write();
579   fOutArray->Write();
580   file->Close();
581   HLTInfo("Histogram file dumped");
582   */
583
584   return 0;
585 }       
586