Using the Altro mapping from calib DB
[u/mrichter/AliRoot.git] / TPC / AliTPCcalibDB.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 // Class providing the calibration parameters by accessing the CDB           //
20 //                                                                           //
21 // Request an instance with AliTPCcalibDB::Instance()                        //
22 // If a new event is processed set the event number with SetRun              //
23 // Then request the calibration data                                         ////
24
25 //
26 //
27 // Calibration data:
28 // 1.)  pad by pad calibration -  AliTPCCalPad
29 //      
30 //      a.) fPadGainFactor
31 //          Simulation: AliTPCDigitizer::ExecFast - Multiply by gain
32 //          Reconstruction : AliTPCclustererMI::Digits2Clusters - Divide by gain  
33 //
34 //      b.) fPadNoise -
35 //          Simulation:        AliTPCDigitizer::ExecFast
36 //          Reconstruction:    AliTPCclustererMI::FindClusters(AliTPCCalROC * noiseROC)
37 //                             Noise depending cut on clusters (n sigma)
38 //      c.) fPedestal:
39 //          Simulation:     Not used yet - To be impleneted - Rounding to the nearest integer
40 //          Reconstruction: Used in AliTPCclustererMI::Digits2Clusters(AliRawReader* rawReader) 
41 //                          if data taken without zero suppression  
42 //                          Currently switch in  fRecoParam->GetCalcPedestal();
43 //      
44 //      d.) fPadTime0
45 //          Simulation:      applied in the AliTPC::MakeSector - adding offset
46 //          Reconstruction:  AliTPCTransform::Transform() - remove offset
47 //                           AliTPCTransform::Transform() - to be called
48 //                           in AliTPCtracker::Transform()      
49 ///////////////////////////////////////////////////////////////////////////////
50
51 #include <iostream>
52 #include <fstream>
53
54
55 #include <AliCDBManager.h>
56 #include <AliCDBEntry.h>
57 #include <AliLog.h>
58
59 #include "AliTPCcalibDB.h"
60 #include "AliTPCAltroMapping.h"
61 #include "AliTPCExB.h"
62
63 #include "AliTPCCalROC.h"
64 #include "AliTPCCalPad.h"
65 #include "AliTPCSensorTempArray.h"
66 #include "AliTPCTransform.h"
67
68 class AliCDBStorage;
69 class AliTPCCalDet;
70 //
71 //
72
73 #include "TFile.h"
74 #include "TKey.h"
75
76 #include "TObjArray.h"
77 #include "TObjString.h"
78 #include "TString.h"
79 #include "AliTPCCalPad.h"
80 #include "AliTPCCalibPulser.h"
81 #include "AliTPCCalibPedestal.h"
82 #include "AliTPCCalibCE.h"
83
84
85
86
87
88 ClassImp(AliTPCcalibDB)
89
90 AliTPCcalibDB* AliTPCcalibDB::fgInstance = 0;
91 Bool_t AliTPCcalibDB::fgTerminated = kFALSE;
92
93
94 //_ singleton implementation __________________________________________________
95 AliTPCcalibDB* AliTPCcalibDB::Instance()
96 {
97   //
98   // Singleton implementation
99   // Returns an instance of this class, it is created if neccessary
100   //
101   
102   if (fgTerminated != kFALSE)
103     return 0;
104
105   if (fgInstance == 0)
106     fgInstance = new AliTPCcalibDB();
107   
108   return fgInstance;
109 }
110
111 void AliTPCcalibDB::Terminate()
112 {
113   //
114   // Singleton implementation
115   // Deletes the instance of this class and sets the terminated flag, instances cannot be requested anymore
116   // This function can be called several times.
117   //
118   
119   fgTerminated = kTRUE;
120   
121   if (fgInstance != 0)
122   {
123     delete fgInstance;
124     fgInstance = 0;
125   }
126 }
127
128 //_____________________________________________________________________________
129 AliTPCcalibDB::AliTPCcalibDB():
130   fRun(-1),
131   fTransform(0),
132   fExB(0),
133   fPadGainFactor(0),
134   fPadTime0(0),
135   fPadNoise(0),
136   fPedestals(0),
137   fTemperature(0),
138   fMapping(0),
139   fParam(0)
140 {
141   //
142   // constructor
143   //  
144   //
145   Update();    // temporary
146 }
147
148 //_____________________________________________________________________________
149 AliTPCcalibDB::~AliTPCcalibDB() 
150 {
151   //
152   // destructor
153   //
154   
155   // don't delete anything, CDB cache is active!
156   //if (fPadGainFactor) delete fPadGainFactor;
157   //if (fPadTime0) delete fPadTime0;
158   //if (fPadNoise) delete fPadNoise;
159 }
160
161
162 //_____________________________________________________________________________
163 AliCDBEntry* AliTPCcalibDB::GetCDBEntry(const char* cdbPath)
164 {
165   // 
166   // Retrieves an entry with path <cdbPath> from the CDB.
167   //
168   char chinfo[1000];
169     
170   AliCDBEntry* entry = AliCDBManager::Instance()->Get(cdbPath, fRun); 
171   if (!entry) 
172   { 
173     sprintf(chinfo,"AliTPCcalibDB: Failed to get entry:\t%s ", cdbPath);
174     AliError(chinfo); 
175     return 0; 
176   }
177   return entry;
178 }
179
180
181 //_____________________________________________________________________________
182 void AliTPCcalibDB::SetRun(Long64_t run)
183 {
184   //
185   // Sets current run number. Calibration data is read from the corresponding file. 
186   //  
187   if (fRun == run)
188     return;  
189   fRun = run;
190   Update();
191 }
192   
193
194
195 void AliTPCcalibDB::Update(){
196   //
197   AliCDBEntry * entry=0;
198   
199   Bool_t cdbCache = AliCDBManager::Instance()->GetCacheFlag(); // save cache status
200   AliCDBManager::Instance()->SetCacheFlag(kTRUE); // activate CDB cache
201   
202   //
203   entry          = GetCDBEntry("TPC/Calib/PadGainFactor");
204   if (entry){
205     //if (fPadGainFactor) delete fPadGainFactor;
206     entry->SetOwner(kTRUE);
207     fPadGainFactor = (AliTPCCalPad*)entry->GetObject();
208   }
209   //
210   entry          = GetCDBEntry("TPC/Calib/PadTime0");
211   if (entry){
212     //if (fPadTime0) delete fPadTime0;
213     entry->SetOwner(kTRUE);
214     fPadTime0 = (AliTPCCalPad*)entry->GetObject();
215   }
216   //
217   //
218   entry          = GetCDBEntry("TPC/Calib/PadNoise");
219   if (entry){
220     //if (fPadNoise) delete fPadNoise;
221     entry->SetOwner(kTRUE);
222     fPadNoise = (AliTPCCalPad*)entry->GetObject();
223   }
224
225   entry          = GetCDBEntry("TPC/Calib/Pedestals");
226   if (entry){
227     //if (fPedestals) delete fPedestals;
228     entry->SetOwner(kTRUE);
229     fPedestals = (AliTPCCalPad*)entry->GetObject();
230   }
231
232   entry          = GetCDBEntry("TPC/Calib/Temperature");
233   if (entry){
234     //if (fTemperature) delete fTemperature;
235     entry->SetOwner(kTRUE);
236     fTemperature = (AliTPCSensorTempArray*)entry->GetObject();
237   }
238
239   entry          = GetCDBEntry("TPC/Calib/Parameters");
240   if (entry){
241     //if (fPadNoise) delete fPadNoise;
242     entry->SetOwner(kTRUE);
243     fParam = (AliTPCParam*)(entry->GetObject()->Clone());
244   }
245
246   entry          = GetCDBEntry("TPC/Calib/Mapping");
247   if (entry){
248     //if (fPadNoise) delete fPadNoise;
249     entry->SetOwner(kTRUE);
250     TObjArray * array = dynamic_cast<TObjArray*>(entry->GetObject());
251     if (array && array->GetEntriesFast()==6){
252       fMapping = new AliTPCAltroMapping*[6];
253       for (Int_t i=0; i<6; i++){
254         fMapping[i] =  dynamic_cast<AliTPCAltroMapping*>(array->At(i));
255       }
256     }
257   }
258
259
260
261   entry          = GetCDBEntry("TPC/Calib/ExB");
262   if (entry) {
263     entry->SetOwner(kTRUE);
264     fExB=dynamic_cast<AliTPCExB*>(entry->GetObject()->Clone());
265   }
266
267   if (!fTransform) {
268     fTransform=new AliTPCTransform(); 
269   }
270
271   //
272   AliCDBManager::Instance()->SetCacheFlag(cdbCache); // reset original CDB cache
273   
274 }
275 AliTPCcalibDB::AliTPCcalibDB(const AliTPCcalibDB& org)
276 {
277   //
278   // Copy constructor invalid -- singleton implementation
279   //
280    Error("copy constructor","invalid -- singleton implementation");
281 }
282
283 AliTPCcalibDB& AliTPCcalibDB::operator= (const AliTPCcalibDB& rhs)
284 {
285 //
286 // Singleton implementation - no assignment operator
287 //
288   Error("operator =", "assignment operator not implemented");
289   return *this;
290 }
291
292
293
294 void AliTPCcalibDB::CreateObjectList(const Char_t *filename, TObjArray *calibObjects)
295 {
296 //
297 // Create calibration objects and read contents from OCDB
298 //
299    if ( calibObjects == 0x0 ) return;
300    ifstream in;
301    in.open(filename);
302    if ( !in.is_open() ){
303       fprintf(stderr,"Error: cannot open list file '%s'", filename);
304       return;
305    }
306    
307    AliTPCCalPad *calPad=0x0;
308    
309    TString sFile;
310    sFile.ReadFile(in);
311    in.close();
312    
313    TObjArray *arrFileLine = sFile.Tokenize("\n");
314    
315    TIter nextLine(arrFileLine);
316    
317    TObjString *sObjLine=0x0;
318    while ( sObjLine = (TObjString*)nextLine() ){
319       TString sLine(sObjLine->GetString());
320       
321       TObjArray *arrNextCol = sLine.Tokenize("\t");
322       
323       TObjString *sObjType     = (TObjString*)(arrNextCol->At(0));
324       TObjString *sObjFileName = (TObjString*)(arrNextCol->At(1));
325       
326       if ( !sObjType || ! sObjFileName ) continue;
327       TString sType(sObjType->GetString());
328       TString sFileName(sObjFileName->GetString());
329       printf("%s\t%s\n",sType.Data(),sFileName.Data());
330       
331       TFile *fIn = TFile::Open(sFileName);
332       if ( !fIn ){
333          fprintf(stderr,"File not found: '%s'", sFileName.Data());
334          continue;
335       }
336       
337       if ( sType == "CE" ){
338          AliTPCCalibCE *ce = (AliTPCCalibCE*)fIn->Get("AliTPCCalibCE");
339          
340          calPad = new AliTPCCalPad((TObjArray*)ce->GetCalPadT0());         
341          calPad->SetNameTitle("CETmean","CETmean");
342          calibObjects->Add(calPad);
343          
344          calPad = new AliTPCCalPad((TObjArray*)ce->GetCalPadQ());         
345          calPad->SetNameTitle("CEQmean","CEQmean");
346          calibObjects->Add(calPad);        
347          
348          calPad = new AliTPCCalPad((TObjArray*)ce->GetCalPadRMS());
349          calPad->SetNameTitle("CETrms","CETrms");
350          calibObjects->Add(calPad);         
351                   
352       } else if ( sType == "Pulser") {
353          AliTPCCalibPulser *sig = (AliTPCCalibPulser*)fIn->Get("AliTPCCalibPulser");
354          
355          calPad = new AliTPCCalPad((TObjArray*)sig->GetCalPadT0());         
356          calPad->SetNameTitle("PulserTmean","PulserTmean");
357          calibObjects->Add(calPad);
358          
359          calPad = new AliTPCCalPad((TObjArray*)sig->GetCalPadQ());         
360          calPad->SetNameTitle("PulserQmean","PulserQmean");
361          calibObjects->Add(calPad);        
362          
363          calPad = new AliTPCCalPad((TObjArray*)sig->GetCalPadRMS());
364          calPad->SetNameTitle("PulserTrms","PulserTrms");
365          calibObjects->Add(calPad);         
366       
367       } else if ( sType == "Pedestals") {
368          AliTPCCalibPedestal *ped = (AliTPCCalibPedestal*)fIn->Get("AliTPCCalibPedestal");
369          
370          calPad = new AliTPCCalPad((TObjArray*)ped->GetCalPadPedestal());         
371          calPad->SetNameTitle("Pedestals","Pedestals");
372          calibObjects->Add(calPad);
373          
374          calPad = new AliTPCCalPad((TObjArray*)ped->GetCalPadRMS());         
375          calPad->SetNameTitle("Noise","Noise");
376          calibObjects->Add(calPad);        
377      
378       } else {
379          fprintf(stderr,"Undefined Type: '%s'",sType.Data());
380          
381       }
382       delete fIn;
383    }
384 }
385
386
387
388 void AliTPCcalibDB::MakeTree(const char * fileName, TObjArray * array, const char * mapFileName, AliTPCCalPad* outlierPad, Float_t ltmFraction) {
389   //
390   // Write a tree with all available information
391   // if mapFileName is specified, the Map information are also written to the tree
392   // pads specified in outlierPad are not used for calculating statistics
393   //  - the same function as AliTPCCalPad::MakeTree - 
394   //
395    AliTPCROC* tpcROCinstance = AliTPCROC::Instance();
396
397    TObjArray* mapIROCs = 0;
398    TObjArray* mapOROCs = 0;
399    TVectorF *mapIROCArray = 0;
400    TVectorF *mapOROCArray = 0;
401    Int_t mapEntries = 0;
402    TString* mapNames = 0;
403    
404    if (mapFileName) {
405       TFile mapFile(mapFileName, "read");
406       
407       TList* listOfROCs = mapFile.GetListOfKeys();
408       mapEntries = listOfROCs->GetEntries()/2;
409       mapIROCs = new TObjArray(mapEntries*2);
410       mapOROCs = new TObjArray(mapEntries*2);
411       mapIROCArray = new TVectorF[mapEntries];
412       mapOROCArray = new TVectorF[mapEntries];
413       
414       mapNames = new TString[mapEntries];
415       for (Int_t ivalue = 0; ivalue < mapEntries; ivalue++) {
416         TString nameROC(((TKey*)(listOfROCs->At(ivalue*2)))->GetName());
417          nameROC.Remove(nameROC.Length()-4, 4);
418          mapIROCs->AddAt((AliTPCCalROC*)mapFile.Get((nameROC + "IROC").Data()), ivalue);
419          mapOROCs->AddAt((AliTPCCalROC*)mapFile.Get((nameROC + "OROC").Data()), ivalue);
420          mapNames[ivalue].Append(nameROC);
421       }
422       
423       for (Int_t ivalue = 0; ivalue < mapEntries; ivalue++) {
424          mapIROCArray[ivalue].ResizeTo(tpcROCinstance->GetNChannels(0));
425          mapOROCArray[ivalue].ResizeTo(tpcROCinstance->GetNChannels(36));
426       
427          for (UInt_t ichannel = 0; ichannel < tpcROCinstance->GetNChannels(0); ichannel++)
428             (mapIROCArray[ivalue])[ichannel] = ((AliTPCCalROC*)(mapIROCs->At(ivalue)))->GetValue(ichannel);
429          for (UInt_t ichannel = 0; ichannel < tpcROCinstance->GetNChannels(36); ichannel++)
430             (mapOROCArray[ivalue])[ichannel] = ((AliTPCCalROC*)(mapOROCs->At(ivalue)))->GetValue(ichannel);
431       }
432
433    } //  if (mapFileName)
434   
435    TTreeSRedirector cstream(fileName);
436    Int_t arrayEntries = array->GetEntries();
437    
438    TString* names = new TString[arrayEntries];
439    for (Int_t ivalue = 0; ivalue < arrayEntries; ivalue++)
440       names[ivalue].Append(((AliTPCCalPad*)array->At(ivalue))->GetName());
441
442    for (UInt_t isector = 0; isector < tpcROCinstance->GetNSectors(); isector++) {
443       //
444       // get statistic for given sector
445       //
446       TVectorF median(arrayEntries);
447       TVectorF mean(arrayEntries);
448       TVectorF rms(arrayEntries);
449       TVectorF ltm(arrayEntries);
450       TVectorF ltmrms(arrayEntries);
451       TVectorF medianWithOut(arrayEntries);
452       TVectorF meanWithOut(arrayEntries);
453       TVectorF rmsWithOut(arrayEntries);
454       TVectorF ltmWithOut(arrayEntries);
455       TVectorF ltmrmsWithOut(arrayEntries);
456       
457       TVectorF *vectorArray = new TVectorF[arrayEntries];
458       for (Int_t ivalue = 0; ivalue < arrayEntries; ivalue++)
459          vectorArray[ivalue].ResizeTo(tpcROCinstance->GetNChannels(isector));
460       
461       for (Int_t ivalue = 0; ivalue < arrayEntries; ivalue++) {
462          AliTPCCalPad* calPad = (AliTPCCalPad*) array->At(ivalue);
463          AliTPCCalROC* calROC = calPad->GetCalROC(isector);
464          AliTPCCalROC* outlierROC = 0;
465          if (outlierPad) outlierROC = outlierPad->GetCalROC(isector);
466          if (calROC) {
467             median[ivalue] = calROC->GetMedian();
468             mean[ivalue] = calROC->GetMean();
469             rms[ivalue] = calROC->GetRMS();
470             Double_t ltmrmsValue = 0;
471             ltm[ivalue] = calROC->GetLTM(&ltmrmsValue, ltmFraction);
472             ltmrms[ivalue] = ltmrmsValue;
473             if (outlierROC) {
474                medianWithOut[ivalue] = calROC->GetMedian(outlierROC);
475                meanWithOut[ivalue] = calROC->GetMean(outlierROC);
476                rmsWithOut[ivalue] = calROC->GetRMS(outlierROC);
477                ltmrmsValue = 0;
478                ltmWithOut[ivalue] = calROC->GetLTM(&ltmrmsValue, ltmFraction, outlierROC);
479                ltmrmsWithOut[ivalue] = ltmrmsValue;
480             }
481          }
482          else {
483             median[ivalue] = 0.;
484             mean[ivalue] = 0.;
485             rms[ivalue] = 0.;
486             ltm[ivalue] = 0.;
487             ltmrms[ivalue] = 0.;
488             medianWithOut[ivalue] = 0.;
489             meanWithOut[ivalue] = 0.;
490             rmsWithOut[ivalue] = 0.;
491             ltmWithOut[ivalue] = 0.;
492             ltmrmsWithOut[ivalue] = 0.;
493          }
494       }
495       
496       //
497       // fill vectors of variable per pad
498       //
499       TVectorF *posArray = new TVectorF[8];
500       for (Int_t ivalue = 0; ivalue < 8; ivalue++)
501          posArray[ivalue].ResizeTo(tpcROCinstance->GetNChannels(isector));
502
503       Float_t posG[3] = {0};
504       Float_t posL[3] = {0};
505       Int_t ichannel = 0;
506       for (UInt_t irow = 0; irow < tpcROCinstance->GetNRows(isector); irow++) {
507          for (UInt_t ipad = 0; ipad < tpcROCinstance->GetNPads(isector, irow); ipad++) {
508             tpcROCinstance->GetPositionLocal(isector, irow, ipad, posL);
509             tpcROCinstance->GetPositionGlobal(isector, irow, ipad, posG);
510             posArray[0][ichannel] = irow;
511             posArray[1][ichannel] = ipad;
512             posArray[2][ichannel] = posL[0];
513             posArray[3][ichannel] = posL[1];
514             posArray[4][ichannel] = posG[0];
515             posArray[5][ichannel] = posG[1];
516             posArray[6][ichannel] = (Int_t)(ipad - (Double_t)(tpcROCinstance->GetNPads(isector, irow))/2);
517             posArray[7][ichannel] = ichannel;
518             
519             // loop over array containing AliTPCCalPads
520             for (Int_t ivalue = 0; ivalue < arrayEntries; ivalue++) {
521                AliTPCCalPad* calPad = (AliTPCCalPad*) array->At(ivalue);
522                AliTPCCalROC* calROC = calPad->GetCalROC(isector);
523                if (calROC)
524                   (vectorArray[ivalue])[ichannel] = calROC->GetValue(irow, ipad);
525                else
526                   (vectorArray[ivalue])[ichannel] = 0;
527             }
528             ichannel++;
529          }
530       }
531       
532       cstream << "calPads" <<
533          "sector=" << isector;
534       
535       for  (Int_t ivalue = 0; ivalue < arrayEntries; ivalue++) {
536          cstream << "calPads" <<
537             (Char_t*)((names[ivalue] + "_Median=").Data()) << median[ivalue] <<
538             (Char_t*)((names[ivalue] + "_Mean=").Data()) << mean[ivalue] <<
539             (Char_t*)((names[ivalue] + "_RMS=").Data()) << rms[ivalue] <<
540             (Char_t*)((names[ivalue] + "_LTM=").Data()) << ltm[ivalue] <<
541             (Char_t*)((names[ivalue] + "_RMS_LTM=").Data()) << ltmrms[ivalue];
542          if (outlierPad) {
543             cstream << "calPads" <<
544                (Char_t*)((names[ivalue] + "_Median_OutlierCutted=").Data()) << medianWithOut[ivalue] <<
545                (Char_t*)((names[ivalue] + "_Mean_OutlierCutted=").Data()) << meanWithOut[ivalue] <<
546                (Char_t*)((names[ivalue] + "_RMS_OutlierCutted=").Data()) << rmsWithOut[ivalue] <<
547                (Char_t*)((names[ivalue] + "_LTM_OutlierCutted=").Data()) << ltmWithOut[ivalue] <<
548                (Char_t*)((names[ivalue] + "_RMS_LTM_OutlierCutted=").Data()) << ltmrmsWithOut[ivalue];
549          }
550       }
551
552       for  (Int_t ivalue = 0; ivalue < arrayEntries; ivalue++) {
553          cstream << "calPads" <<
554             (Char_t*)((names[ivalue] + ".=").Data()) << &vectorArray[ivalue];
555       }
556
557       if (mapFileName) {
558          for  (Int_t ivalue = 0; ivalue < mapEntries; ivalue++) {
559             if (isector < 36)
560                cstream << "calPads" <<
561                   (Char_t*)((mapNames[ivalue] + ".=").Data()) << &mapIROCArray[ivalue];
562             else
563                cstream << "calPads" <<
564                   (Char_t*)((mapNames[ivalue] + ".=").Data()) << &mapOROCArray[ivalue];
565          }
566       }
567
568       cstream << "calPads" <<
569          "row.=" << &posArray[0] <<
570          "pad.=" << &posArray[1] <<
571          "lx.=" << &posArray[2] <<
572          "ly.=" << &posArray[3] <<
573          "gx.=" << &posArray[4] <<
574          "gy.=" << &posArray[5] <<
575          "rpad.=" << &posArray[6] <<
576          "channel.=" << &posArray[7];
577          
578       cstream << "calPads" <<
579          "\n";
580
581       delete[] posArray;
582       delete[] vectorArray;
583    }
584    
585
586    delete[] names;
587    if (mapFileName) {
588       delete mapIROCs;
589       delete mapOROCs;
590       delete[] mapIROCArray;
591       delete[] mapOROCArray;
592       delete[] mapNames;
593    }
594 }
595