]> git.uio.no Git - u/mrichter/AliRoot.git/blob - TPC/AliTPCcalibDB.cxx
Adding interface for Run info
[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 // Calibration data:
27 // 0.)  Altro mapping
28 //          Simulation      - not yet 
29 //          Reconstruction  - AliTPCclustererMI::Digits2Clusters(AliRawReader* rawReader)
30 //
31 // 1.)  pad by pad calibration -  AliTPCCalPad
32 //      
33 //      a.) fPadGainFactor
34 //          Simulation: AliTPCDigitizer::ExecFast - Multiply by gain
35 //          Reconstruction : AliTPCclustererMI::Digits2Clusters - Divide by gain  
36 //
37 //      b.) fPadNoise -
38 //          Simulation:        AliTPCDigitizer::ExecFast
39 //          Reconstruction:    AliTPCclustererMI::FindClusters(AliTPCCalROC * noiseROC)
40 //                             Noise depending cut on clusters charge (n sigma)
41 //      c.) fPedestal:
42 //          Simulation:     Not used yet - To be impleneted - Rounding to the nearest integer
43 //          Reconstruction: Used in AliTPCclustererMI::Digits2Clusters(AliRawReader* rawReader) 
44 //                          if data taken without zero suppression  
45 //                          Currently switch in  fRecoParam->GetCalcPedestal();
46 //      
47 //      d.) fPadTime0
48 //          Simulation:      applied in the AliTPC::MakeSector - adding offset
49 //          Reconstruction:  AliTPCTransform::Transform() - remove offset
50 //                           AliTPCTransform::Transform() - to be called
51 //                           in AliTPCtracker::Transform()      
52 //
53 // 
54 // 2.)  Space points transformation:
55 //
56 //      a.) General coordinate tranformation - AliTPCtransform (see $ALICE_ROOT/TPC/AliTPCtransform.cxx)
57 //          Created on fly - use the other calibration components
58 //                 Unisochronity  - (substract time0 - pad by pad)
59 //                 Drift velocity - Currently common drift velocity - functionality of AliTPCParam
60 //                 ExB effect    
61 //          Simulation     - Not used directly (the effects are applied one by one (see AliTPC::MakeSector)
62 //          Reconstruction - 
63 //                           AliTPCclustererMI::AddCluster
64 //                           AliTPCtrackerMI::Transform
65 //      b.) ExB effect calibration - 
66 //             classes (base class AliTPCExB, implementation- AliTPCExBExact.h  AliTPCExBFirst.h)
67 //             a.a) Simulation:   applied in the AliTPC::MakeSector - 
68 //                                calib->GetExB()->CorrectInverse(dxyz0,dxyz1);
69 //             a.b) Reconstruction -  
70 //                  
71 //                  in AliTPCtransform::Correct() - called calib->GetExB()->Correct(dxyz0,dxyz1)
72 //
73 //  3.)   cluster error, shape and Q parameterization
74 //
75 //
76 //
77 ///////////////////////////////////////////////////////////////////////////////
78
79 #include <iostream>
80 #include <fstream>
81
82
83 #include <AliCDBManager.h>
84 #include <AliCDBEntry.h>
85 #include <AliLog.h>
86 #include <AliMagF.h>
87 #include <AliMagWrapCheb.h>
88
89 #include "AliTPCcalibDB.h"
90 #include "AliTPCAltroMapping.h"
91 #include "AliTPCExB.h"
92
93 #include "AliTPCCalROC.h"
94 #include "AliTPCCalPad.h"
95 #include "AliTPCSensorTempArray.h"
96 #include "AliGRPObject.h"
97 #include "AliTPCTransform.h"
98
99 class AliCDBStorage;
100 class AliTPCCalDet;
101 //
102 //
103
104 #include "TFile.h"
105 #include "TKey.h"
106
107 #include "TObjArray.h"
108 #include "TObjString.h"
109 #include "TString.h"
110 #include "AliTPCCalPad.h"
111 #include "AliTPCCalibPulser.h"
112 #include "AliTPCCalibPedestal.h"
113 #include "AliTPCCalibCE.h"
114 #include "AliTPCExBFirst.h"
115 #include "AliTPCTempMap.h"
116 #include "AliTPCCalibVdrift.h"
117
118
119
120
121 ClassImp(AliTPCcalibDB)
122
123 AliTPCcalibDB* AliTPCcalibDB::fgInstance = 0;
124 Bool_t AliTPCcalibDB::fgTerminated = kFALSE;
125 TObjArray    AliTPCcalibDB::fgExBArray;    // array of ExB corrections
126
127
128 //_ singleton implementation __________________________________________________
129 AliTPCcalibDB* AliTPCcalibDB::Instance()
130 {
131   //
132   // Singleton implementation
133   // Returns an instance of this class, it is created if neccessary
134   //
135   
136   if (fgTerminated != kFALSE)
137     return 0;
138
139   if (fgInstance == 0)
140     fgInstance = new AliTPCcalibDB();
141   
142   return fgInstance;
143 }
144
145 void AliTPCcalibDB::Terminate()
146 {
147   //
148   // Singleton implementation
149   // Deletes the instance of this class and sets the terminated flag, instances cannot be requested anymore
150   // This function can be called several times.
151   //
152   
153   fgTerminated = kTRUE;
154   
155   if (fgInstance != 0)
156   {
157     delete fgInstance;
158     fgInstance = 0;
159   }
160 }
161
162 //_____________________________________________________________________________
163 AliTPCcalibDB::AliTPCcalibDB():
164   TObject(),
165   fRun(-1),
166   fTransform(0),
167   fExB(0),
168   fPadGainFactor(0),
169   fDedxGainFactor(0),
170   fPadTime0(0),
171   fPadNoise(0),
172   fPedestals(0),
173   fTemperature(0),
174   fMapping(0),
175   fParam(0),
176   fClusterParam(0),  
177   fGRPArray(100000),            //! array of GRPs  -  per run  - JUST for calibration studies
178   fGRPMaps(100000),            //! array of GRPs  -  per run  - JUST for calibration studies
179   fGoofieArray(100000),         //! array of GOOFIE values -per run - Just for calibration studies
180   fTemperatureArray(100000),    //! array of temperature sensors - per run - Just for calibration studies
181   fVdriftArray(100000),                 //! array of v drift interfaces
182   fRunList(100000)              //! run list - indicates try to get the run param 
183
184 {
185   //
186   // constructor
187   //  
188   //
189   Update();    // temporary
190 }
191
192 AliTPCcalibDB::AliTPCcalibDB(const AliTPCcalibDB& ):
193   TObject(),
194   fRun(-1),
195   fTransform(0),
196   fExB(0),
197   fPadGainFactor(0),
198   fDedxGainFactor(0),
199   fPadTime0(0),
200   fPadNoise(0),
201   fPedestals(0),
202   fTemperature(0),
203   fMapping(0),
204   fParam(0),
205   fClusterParam(0),
206   fGRPArray(0),          //! array of GRPs  -  per run  - JUST for calibration studies
207   fGRPMaps(0),          //! array of GRPs  -  per run  - JUST for calibration studies
208   fGoofieArray(0),        //! array of GOOFIE values -per run - Just for calibration studies
209   fTemperatureArray(0),   //! array of temperature sensors - per run - Just for calibration studies
210   fVdriftArray(0),         //! array of v drift interfaces
211   fRunList(0)              //! run list - indicates try to get the run param 
212 {
213   //
214   // Copy constructor invalid -- singleton implementation
215   //
216    Error("copy constructor","invalid -- singleton implementation");
217 }
218
219 AliTPCcalibDB& AliTPCcalibDB::operator= (const AliTPCcalibDB& )
220 {
221 //
222 // Singleton implementation - no assignment operator
223 //
224   Error("operator =", "assignment operator not implemented");
225   return *this;
226 }
227
228
229
230 //_____________________________________________________________________________
231 AliTPCcalibDB::~AliTPCcalibDB() 
232 {
233   //
234   // destructor
235   //
236   
237   // don't delete anything, CDB cache is active!
238   //if (fPadGainFactor) delete fPadGainFactor;
239   //if (fPadTime0) delete fPadTime0;
240   //if (fPadNoise) delete fPadNoise;
241 }
242
243
244 //_____________________________________________________________________________
245 AliCDBEntry* AliTPCcalibDB::GetCDBEntry(const char* cdbPath)
246 {
247   // 
248   // Retrieves an entry with path <cdbPath> from the CDB.
249   //
250   char chinfo[1000];
251     
252   AliCDBEntry* entry = AliCDBManager::Instance()->Get(cdbPath, fRun); 
253   if (!entry) 
254   { 
255     sprintf(chinfo,"AliTPCcalibDB: Failed to get entry:\t%s ", cdbPath);
256     AliError(chinfo); 
257     return 0; 
258   }
259   return entry;
260 }
261
262
263 //_____________________________________________________________________________
264 void AliTPCcalibDB::SetRun(Long64_t run)
265 {
266   //
267   // Sets current run number. Calibration data is read from the corresponding file. 
268   //  
269   if (fRun == run)
270     return;  
271   fRun = run;
272   Update();
273 }
274   
275
276
277 void AliTPCcalibDB::Update(){
278   //
279   AliCDBEntry * entry=0;
280   
281   Bool_t cdbCache = AliCDBManager::Instance()->GetCacheFlag(); // save cache status
282   AliCDBManager::Instance()->SetCacheFlag(kTRUE); // activate CDB cache
283   
284   //
285   entry          = GetCDBEntry("TPC/Calib/PadGainFactor");
286   if (entry){
287     //if (fPadGainFactor) delete fPadGainFactor;
288     entry->SetOwner(kTRUE);
289     fPadGainFactor = (AliTPCCalPad*)entry->GetObject();
290   }
291   //
292   entry          = GetCDBEntry("TPC/Calib/GainFactorDedx");
293   if (entry){
294     entry->SetOwner(kTRUE);
295     fDedxGainFactor = (AliTPCCalPad*)entry->GetObject();
296   }
297   //
298   entry          = GetCDBEntry("TPC/Calib/PadTime0");
299   if (entry){
300     //if (fPadTime0) delete fPadTime0;
301     entry->SetOwner(kTRUE);
302     fPadTime0 = (AliTPCCalPad*)entry->GetObject();
303   }
304   //
305   //
306   entry          = GetCDBEntry("TPC/Calib/PadNoise");
307   if (entry){
308     //if (fPadNoise) delete fPadNoise;
309     entry->SetOwner(kTRUE);
310     fPadNoise = (AliTPCCalPad*)entry->GetObject();
311   }
312
313   entry          = GetCDBEntry("TPC/Calib/Pedestals");
314   if (entry){
315     //if (fPedestals) delete fPedestals;
316     entry->SetOwner(kTRUE);
317     fPedestals = (AliTPCCalPad*)entry->GetObject();
318   }
319
320   entry          = GetCDBEntry("TPC/Calib/Temperature");
321   if (entry){
322     //if (fTemperature) delete fTemperature;
323     entry->SetOwner(kTRUE);
324     fTemperature = (AliTPCSensorTempArray*)entry->GetObject();
325   }
326
327   entry          = GetCDBEntry("TPC/Calib/Parameters");
328   if (entry){
329     //if (fPadNoise) delete fPadNoise;
330     entry->SetOwner(kTRUE);
331     fParam = (AliTPCParam*)(entry->GetObject()->Clone());
332   }
333
334   entry          = GetCDBEntry("TPC/Calib/ClusterParam");
335   if (entry){
336     //if (fPadNoise) delete fPadNoise;
337     entry->SetOwner(kTRUE);
338     fClusterParam = (AliTPCClusterParam*)(entry->GetObject()->Clone());
339   }
340
341   entry          = GetCDBEntry("TPC/Calib/Mapping");
342   if (entry){
343     //if (fPadNoise) delete fPadNoise;
344     entry->SetOwner(kTRUE);
345     TObjArray * array = dynamic_cast<TObjArray*>(entry->GetObject());
346     if (array && array->GetEntriesFast()==6){
347       fMapping = new AliTPCAltroMapping*[6];
348       for (Int_t i=0; i<6; i++){
349         fMapping[i] =  dynamic_cast<AliTPCAltroMapping*>(array->At(i));
350       }
351     }
352   }
353
354
355
356   //entry          = GetCDBEntry("TPC/Calib/ExB");
357   //if (entry) {
358   //  entry->SetOwner(kTRUE);
359   //  fExB=dynamic_cast<AliTPCExB*>(entry->GetObject()->Clone());
360   //}
361   //
362   // ExB  - calculate during initialization 
363   //      - 
364   fExB =  GetExB(-5,kTRUE);
365      //
366   if (!fTransform) {
367     fTransform=new AliTPCTransform(); 
368   }
369
370   //
371   AliCDBManager::Instance()->SetCacheFlag(cdbCache); // reset original CDB cache
372   
373 }
374
375
376
377 void AliTPCcalibDB::CreateObjectList(const Char_t *filename, TObjArray *calibObjects)
378 {
379 //
380 // Create calibration objects and read contents from OCDB
381 //
382    if ( calibObjects == 0x0 ) return;
383    ifstream in;
384    in.open(filename);
385    if ( !in.is_open() ){
386       fprintf(stderr,"Error: cannot open list file '%s'", filename);
387       return;
388    }
389    
390    AliTPCCalPad *calPad=0x0;
391    
392    TString sFile;
393    sFile.ReadFile(in);
394    in.close();
395    
396    TObjArray *arrFileLine = sFile.Tokenize("\n");
397    
398    TIter nextLine(arrFileLine);
399    
400    TObjString *sObjLine=0x0;
401    while ( (sObjLine = (TObjString*)nextLine()) ){
402       TString sLine(sObjLine->GetString());
403       
404       TObjArray *arrNextCol = sLine.Tokenize("\t");
405       
406       TObjString *sObjType     = (TObjString*)(arrNextCol->At(0));
407       TObjString *sObjFileName = (TObjString*)(arrNextCol->At(1));
408       
409       if ( !sObjType || ! sObjFileName ) continue;
410       TString sType(sObjType->GetString());
411       TString sFileName(sObjFileName->GetString());
412       printf("%s\t%s\n",sType.Data(),sFileName.Data());
413       
414       TFile *fIn = TFile::Open(sFileName);
415       if ( !fIn ){
416          fprintf(stderr,"File not found: '%s'", sFileName.Data());
417          continue;
418       }
419       
420       if ( sType == "CE" ){
421          AliTPCCalibCE *ce = (AliTPCCalibCE*)fIn->Get("AliTPCCalibCE");
422          
423          calPad = new AliTPCCalPad((TObjArray*)ce->GetCalPadT0());         
424          calPad->SetNameTitle("CETmean","CETmean");
425          calibObjects->Add(calPad);
426          
427          calPad = new AliTPCCalPad((TObjArray*)ce->GetCalPadQ());         
428          calPad->SetNameTitle("CEQmean","CEQmean");
429          calibObjects->Add(calPad);        
430          
431          calPad = new AliTPCCalPad((TObjArray*)ce->GetCalPadRMS());
432          calPad->SetNameTitle("CETrms","CETrms");
433          calibObjects->Add(calPad);         
434                   
435       } else if ( sType == "Pulser") {
436          AliTPCCalibPulser *sig = (AliTPCCalibPulser*)fIn->Get("AliTPCCalibPulser");
437          
438          calPad = new AliTPCCalPad((TObjArray*)sig->GetCalPadT0());         
439          calPad->SetNameTitle("PulserTmean","PulserTmean");
440          calibObjects->Add(calPad);
441          
442          calPad = new AliTPCCalPad((TObjArray*)sig->GetCalPadQ());         
443          calPad->SetNameTitle("PulserQmean","PulserQmean");
444          calibObjects->Add(calPad);        
445          
446          calPad = new AliTPCCalPad((TObjArray*)sig->GetCalPadRMS());
447          calPad->SetNameTitle("PulserTrms","PulserTrms");
448          calibObjects->Add(calPad);         
449       
450       } else if ( sType == "Pedestals") {
451          AliTPCCalibPedestal *ped = (AliTPCCalibPedestal*)fIn->Get("AliTPCCalibPedestal");
452          
453          calPad = new AliTPCCalPad((TObjArray*)ped->GetCalPadPedestal());         
454          calPad->SetNameTitle("Pedestals","Pedestals");
455          calibObjects->Add(calPad);
456          
457          calPad = new AliTPCCalPad((TObjArray*)ped->GetCalPadRMS());         
458          calPad->SetNameTitle("Noise","Noise");
459          calibObjects->Add(calPad);        
460      
461       } else {
462          fprintf(stderr,"Undefined Type: '%s'",sType.Data());
463          
464       }
465       delete fIn;
466    }
467 }
468
469
470
471 void AliTPCcalibDB::MakeTree(const char * fileName, TObjArray * array, const char * mapFileName, AliTPCCalPad* outlierPad, Float_t ltmFraction) {
472   //
473   // Write a tree with all available information
474   // if mapFileName is specified, the Map information are also written to the tree
475   // pads specified in outlierPad are not used for calculating statistics
476   //  - the same function as AliTPCCalPad::MakeTree - 
477   //
478    AliTPCROC* tpcROCinstance = AliTPCROC::Instance();
479
480    TObjArray* mapIROCs = 0;
481    TObjArray* mapOROCs = 0;
482    TVectorF *mapIROCArray = 0;
483    TVectorF *mapOROCArray = 0;
484    Int_t mapEntries = 0;
485    TString* mapNames = 0;
486    
487    if (mapFileName) {
488       TFile mapFile(mapFileName, "read");
489       
490       TList* listOfROCs = mapFile.GetListOfKeys();
491       mapEntries = listOfROCs->GetEntries()/2;
492       mapIROCs = new TObjArray(mapEntries*2);
493       mapOROCs = new TObjArray(mapEntries*2);
494       mapIROCArray = new TVectorF[mapEntries];
495       mapOROCArray = new TVectorF[mapEntries];
496       
497       mapNames = new TString[mapEntries];
498       for (Int_t ivalue = 0; ivalue < mapEntries; ivalue++) {
499         TString nameROC(((TKey*)(listOfROCs->At(ivalue*2)))->GetName());
500          nameROC.Remove(nameROC.Length()-4, 4);
501          mapIROCs->AddAt((AliTPCCalROC*)mapFile.Get((nameROC + "IROC").Data()), ivalue);
502          mapOROCs->AddAt((AliTPCCalROC*)mapFile.Get((nameROC + "OROC").Data()), ivalue);
503          mapNames[ivalue].Append(nameROC);
504       }
505       
506       for (Int_t ivalue = 0; ivalue < mapEntries; ivalue++) {
507          mapIROCArray[ivalue].ResizeTo(tpcROCinstance->GetNChannels(0));
508          mapOROCArray[ivalue].ResizeTo(tpcROCinstance->GetNChannels(36));
509       
510          for (UInt_t ichannel = 0; ichannel < tpcROCinstance->GetNChannels(0); ichannel++)
511             (mapIROCArray[ivalue])[ichannel] = ((AliTPCCalROC*)(mapIROCs->At(ivalue)))->GetValue(ichannel);
512          for (UInt_t ichannel = 0; ichannel < tpcROCinstance->GetNChannels(36); ichannel++)
513             (mapOROCArray[ivalue])[ichannel] = ((AliTPCCalROC*)(mapOROCs->At(ivalue)))->GetValue(ichannel);
514       }
515
516    } //  if (mapFileName)
517   
518    TTreeSRedirector cstream(fileName);
519    Int_t arrayEntries = array->GetEntries();
520    
521    TString* names = new TString[arrayEntries];
522    for (Int_t ivalue = 0; ivalue < arrayEntries; ivalue++)
523       names[ivalue].Append(((AliTPCCalPad*)array->At(ivalue))->GetName());
524
525    for (UInt_t isector = 0; isector < tpcROCinstance->GetNSectors(); isector++) {
526       //
527       // get statistic for given sector
528       //
529       TVectorF median(arrayEntries);
530       TVectorF mean(arrayEntries);
531       TVectorF rms(arrayEntries);
532       TVectorF ltm(arrayEntries);
533       TVectorF ltmrms(arrayEntries);
534       TVectorF medianWithOut(arrayEntries);
535       TVectorF meanWithOut(arrayEntries);
536       TVectorF rmsWithOut(arrayEntries);
537       TVectorF ltmWithOut(arrayEntries);
538       TVectorF ltmrmsWithOut(arrayEntries);
539       
540       TVectorF *vectorArray = new TVectorF[arrayEntries];
541       for (Int_t ivalue = 0; ivalue < arrayEntries; ivalue++)
542          vectorArray[ivalue].ResizeTo(tpcROCinstance->GetNChannels(isector));
543       
544       for (Int_t ivalue = 0; ivalue < arrayEntries; ivalue++) {
545          AliTPCCalPad* calPad = (AliTPCCalPad*) array->At(ivalue);
546          AliTPCCalROC* calROC = calPad->GetCalROC(isector);
547          AliTPCCalROC* outlierROC = 0;
548          if (outlierPad) outlierROC = outlierPad->GetCalROC(isector);
549          if (calROC) {
550             median[ivalue] = calROC->GetMedian();
551             mean[ivalue] = calROC->GetMean();
552             rms[ivalue] = calROC->GetRMS();
553             Double_t ltmrmsValue = 0;
554             ltm[ivalue] = calROC->GetLTM(&ltmrmsValue, ltmFraction);
555             ltmrms[ivalue] = ltmrmsValue;
556             if (outlierROC) {
557                medianWithOut[ivalue] = calROC->GetMedian(outlierROC);
558                meanWithOut[ivalue] = calROC->GetMean(outlierROC);
559                rmsWithOut[ivalue] = calROC->GetRMS(outlierROC);
560                ltmrmsValue = 0;
561                ltmWithOut[ivalue] = calROC->GetLTM(&ltmrmsValue, ltmFraction, outlierROC);
562                ltmrmsWithOut[ivalue] = ltmrmsValue;
563             }
564          }
565          else {
566             median[ivalue] = 0.;
567             mean[ivalue] = 0.;
568             rms[ivalue] = 0.;
569             ltm[ivalue] = 0.;
570             ltmrms[ivalue] = 0.;
571             medianWithOut[ivalue] = 0.;
572             meanWithOut[ivalue] = 0.;
573             rmsWithOut[ivalue] = 0.;
574             ltmWithOut[ivalue] = 0.;
575             ltmrmsWithOut[ivalue] = 0.;
576          }
577       }
578       
579       //
580       // fill vectors of variable per pad
581       //
582       TVectorF *posArray = new TVectorF[8];
583       for (Int_t ivalue = 0; ivalue < 8; ivalue++)
584          posArray[ivalue].ResizeTo(tpcROCinstance->GetNChannels(isector));
585
586       Float_t posG[3] = {0};
587       Float_t posL[3] = {0};
588       Int_t ichannel = 0;
589       for (UInt_t irow = 0; irow < tpcROCinstance->GetNRows(isector); irow++) {
590          for (UInt_t ipad = 0; ipad < tpcROCinstance->GetNPads(isector, irow); ipad++) {
591             tpcROCinstance->GetPositionLocal(isector, irow, ipad, posL);
592             tpcROCinstance->GetPositionGlobal(isector, irow, ipad, posG);
593             posArray[0][ichannel] = irow;
594             posArray[1][ichannel] = ipad;
595             posArray[2][ichannel] = posL[0];
596             posArray[3][ichannel] = posL[1];
597             posArray[4][ichannel] = posG[0];
598             posArray[5][ichannel] = posG[1];
599             posArray[6][ichannel] = (Int_t)(ipad - (Double_t)(tpcROCinstance->GetNPads(isector, irow))/2);
600             posArray[7][ichannel] = ichannel;
601             
602             // loop over array containing AliTPCCalPads
603             for (Int_t ivalue = 0; ivalue < arrayEntries; ivalue++) {
604                AliTPCCalPad* calPad = (AliTPCCalPad*) array->At(ivalue);
605                AliTPCCalROC* calROC = calPad->GetCalROC(isector);
606                if (calROC)
607                   (vectorArray[ivalue])[ichannel] = calROC->GetValue(irow, ipad);
608                else
609                   (vectorArray[ivalue])[ichannel] = 0;
610             }
611             ichannel++;
612          }
613       }
614       
615       cstream << "calPads" <<
616          "sector=" << isector;
617       
618       for  (Int_t ivalue = 0; ivalue < arrayEntries; ivalue++) {
619          cstream << "calPads" <<
620             (Char_t*)((names[ivalue] + "_Median=").Data()) << median[ivalue] <<
621             (Char_t*)((names[ivalue] + "_Mean=").Data()) << mean[ivalue] <<
622             (Char_t*)((names[ivalue] + "_RMS=").Data()) << rms[ivalue] <<
623             (Char_t*)((names[ivalue] + "_LTM=").Data()) << ltm[ivalue] <<
624             (Char_t*)((names[ivalue] + "_RMS_LTM=").Data()) << ltmrms[ivalue];
625          if (outlierPad) {
626             cstream << "calPads" <<
627                (Char_t*)((names[ivalue] + "_Median_OutlierCutted=").Data()) << medianWithOut[ivalue] <<
628                (Char_t*)((names[ivalue] + "_Mean_OutlierCutted=").Data()) << meanWithOut[ivalue] <<
629                (Char_t*)((names[ivalue] + "_RMS_OutlierCutted=").Data()) << rmsWithOut[ivalue] <<
630                (Char_t*)((names[ivalue] + "_LTM_OutlierCutted=").Data()) << ltmWithOut[ivalue] <<
631                (Char_t*)((names[ivalue] + "_RMS_LTM_OutlierCutted=").Data()) << ltmrmsWithOut[ivalue];
632          }
633       }
634
635       for  (Int_t ivalue = 0; ivalue < arrayEntries; ivalue++) {
636          cstream << "calPads" <<
637             (Char_t*)((names[ivalue] + ".=").Data()) << &vectorArray[ivalue];
638       }
639
640       if (mapFileName) {
641          for  (Int_t ivalue = 0; ivalue < mapEntries; ivalue++) {
642             if (isector < 36)
643                cstream << "calPads" <<
644                   (Char_t*)((mapNames[ivalue] + ".=").Data()) << &mapIROCArray[ivalue];
645             else
646                cstream << "calPads" <<
647                   (Char_t*)((mapNames[ivalue] + ".=").Data()) << &mapOROCArray[ivalue];
648          }
649       }
650
651       cstream << "calPads" <<
652          "row.=" << &posArray[0] <<
653          "pad.=" << &posArray[1] <<
654          "lx.=" << &posArray[2] <<
655          "ly.=" << &posArray[3] <<
656          "gx.=" << &posArray[4] <<
657          "gy.=" << &posArray[5] <<
658          "rpad.=" << &posArray[6] <<
659          "channel.=" << &posArray[7];
660          
661       cstream << "calPads" <<
662          "\n";
663
664       delete[] posArray;
665       delete[] vectorArray;
666    }
667    
668
669    delete[] names;
670    if (mapFileName) {
671       delete mapIROCs;
672       delete mapOROCs;
673       delete[] mapIROCArray;
674       delete[] mapOROCArray;
675       delete[] mapNames;
676    }
677 }
678
679
680
681 void AliTPCcalibDB::RegisterExB(Int_t index, Float_t bz, Bool_t bdelete){
682   //
683   // Register static ExB correction map
684   // index - registration index - used for visualization
685   // bz    - bz field in kGaus
686
687   Float_t factor =  bz/(-5.);  // default b filed in Cheb with minus sign
688   
689   AliMagF*   bmap = new AliMagWrapCheb("MapsExB","MapsExB", 2, factor, 10., AliMagWrapCheb::k5kG,kTRUE,"$(ALICE_ROOT)/data/maps/mfchebKGI_sym.root");
690   
691   AliTPCExBFirst *exb  = new  AliTPCExBFirst(bmap,0.88*2.6400e+04,50,50,50);
692   AliTPCExB::SetInstance(exb);
693   
694   if (bdelete){
695     delete bmap;
696   }else{
697     AliTPCExB::RegisterField(index,bmap);
698   }
699   if (index>=fgExBArray.GetEntries()) fgExBArray.Expand((index+1)*2+11);
700   fgExBArray.AddAt(exb,index);
701 }
702
703
704 AliTPCExB*    AliTPCcalibDB::GetExB(Float_t bz, Bool_t deleteB) {
705   //
706   // bz filed in KGaus not in tesla
707   // Get ExB correction map
708   // if doesn't exist - create it
709   //
710   Int_t index = TMath::Nint(5+bz);
711   if (index>fgExBArray.GetEntries()) fgExBArray.Expand((index+1)*2+11);
712   if (!fgExBArray.At(index)) AliTPCcalibDB::RegisterExB(index,bz,deleteB);
713   return (AliTPCExB*)fgExBArray.At(index);
714 }
715
716
717 void  AliTPCcalibDB::SetExBField(Float_t bz){
718   //
719   // Set magnetic filed for ExB correction
720   //
721   fExB = GetExB(bz,kFALSE);
722 }
723
724
725
726 void AliTPCcalibDB::GetRunInformations( Int_t run){
727   //
728   // - > Don't use it for reconstruction - Only for Calibration studies
729   //
730   AliCDBEntry * entry = 0;
731   if (run>= fRunList.GetSize()){
732     fRunList.Set(run*2+1);
733     fGRPArray.Expand(run*2+1);
734     fGRPMaps.Expand(run*2+1);
735     fGoofieArray.Expand(run*2+1); 
736     fTemperatureArray.Expand(run*2+1);
737     fVdriftArray.Expand(run*2+1);
738   }
739   if (fRunList[run]>0) return;
740   entry = AliCDBManager::Instance()->Get("GRP/GRP/Data",run);
741   if (entry)  {
742     AliGRPObject * grpRun = dynamic_cast<AliGRPObject*>(entry->GetObject());
743     if (!grpRun){
744       TMap* map = dynamic_cast<TMap*>(entry->GetObject());
745       if (map){
746         grpRun = new AliGRPObject; 
747         grpRun->ReadValuesFromMap(map);
748         fGRPMaps.AddAt(map,run);
749       }
750     }
751     fGRPArray.AddAt(grpRun,run);
752   }
753   entry = AliCDBManager::Instance()->Get("TPC/Calib/Goofie",run);
754   if (entry)  fGoofieArray.AddAt(entry->GetObject(),run);
755   entry = AliCDBManager::Instance()->Get("TPC/Calib/Temperature",run);
756   if (entry)  fTemperatureArray.AddAt(entry->GetObject(),run);
757   fRunList[run]=1;  // sign as used
758
759   AliDCSSensor * press = GetPressureSensor(run);
760   AliTPCSensorTempArray * temp = GetTemperatureSensor(run);
761   if (press && temp){
762     AliTPCCalibVdrift * vdrift = new AliTPCCalibVdrift(temp, press,0);
763     fVdriftArray.AddAt(vdrift,run);
764   }
765 }
766
767
768 Float_t AliTPCcalibDB::GetGain(Int_t sector, Int_t row, Int_t pad){
769   //
770   //
771   AliTPCCalPad *calPad = Instance()->fDedxGainFactor;;
772   if (!calPad) return 0;
773   return calPad->GetCalROC(sector)->GetValue(row,pad);
774 }
775
776
777 AliGRPObject *AliTPCcalibDB::GetGRP(Int_t run){
778   //
779   // Get GRP object for given run 
780   //
781   AliGRPObject * grpRun = dynamic_cast<AliGRPObject *>((Instance()->fGRPArray).At(run));
782   if (!grpRun) {
783     Instance()->GetRunInformations(run);
784     grpRun = dynamic_cast<AliGRPObject *>(Instance()->fGRPArray.At(run));
785     if (!grpRun) return 0; 
786   }
787   return grpRun;
788 }
789
790 TMap *  AliTPCcalibDB::GetGRPMap(Int_t run){
791   //
792   //
793   //
794   TMap * grpRun = dynamic_cast<TMap *>((Instance()->fGRPMaps).At(run));
795   if (!grpRun) {
796     Instance()->GetRunInformations(run);
797     grpRun = dynamic_cast<TMap *>(Instance()->fGRPMaps.At(run));
798     if (!grpRun) return 0; 
799   }
800   return grpRun;
801 }
802
803
804 AliDCSSensor * AliTPCcalibDB::GetPressureSensor(Int_t run, Int_t type){
805   //
806   // Get Pressure sensor
807   //
808   //
809   // First try to get if trom map - if existing  (Old format of data storing)
810   //
811   TMap *map = GetGRPMap(run);  
812   if (map){
813     AliDCSSensor * sensor = 0;
814     TObject *osensor=0;
815     if (type==0) osensor = ((*map)("fCavernPressure"));
816     if (type==1) osensor = ((*map)("fP2Pressure"));
817     sensor =dynamic_cast<AliDCSSensor *>(osensor); 
818     if (sensor) return sensor;
819   }
820   //
821   // If not map try to get it from the GRPObject
822   //
823   AliGRPObject * grpRun = dynamic_cast<AliGRPObject *>(fGRPArray.At(run)); 
824   if (!grpRun) {
825     GetRunInformations(run);
826     grpRun = dynamic_cast<AliGRPObject *>(fGRPArray.At(run));
827     if (!grpRun) return 0; 
828   }
829   AliDCSSensor * sensor = grpRun->GetCavernAtmosPressure();
830   if (type==1) sensor = grpRun->GetSurfaceAtmosPressure();
831   return sensor; 
832 }
833
834 AliTPCSensorTempArray * AliTPCcalibDB::GetTemperatureSensor(Int_t run){
835   //
836   // Get temperature sensor array
837   //
838   AliTPCSensorTempArray * tempArray = (AliTPCSensorTempArray *)fTemperatureArray.At(run);
839   if (!tempArray) {
840     GetRunInformations(run);
841     tempArray = (AliTPCSensorTempArray *)fTemperatureArray.At(run);
842   }
843   return tempArray;
844 }
845
846 AliDCSSensorArray * AliTPCcalibDB::GetGoofieSensors(Int_t run){
847   //
848   // Get temperature sensor array
849   //
850   AliDCSSensorArray * goofieArray = (AliDCSSensorArray *)fGoofieArray.At(run);
851   if (!goofieArray) {
852     GetRunInformations(run);
853     goofieArray = (AliDCSSensorArray *)fGoofieArray.At(run);
854   }
855   return goofieArray;
856 }
857
858 AliTPCCalibVdrift *     AliTPCcalibDB::GetVdrift(Int_t run){
859   //
860   // Get the interface to the the vdrift 
861   //
862   AliTPCCalibVdrift  * vdrift = (AliTPCCalibVdrift*)fVdriftArray.At(run);
863   if (!vdrift) {
864     GetRunInformations(run);
865     vdrift= (AliTPCCalibVdrift*)fVdriftArray.At(run);
866   }
867   return vdrift;
868 }
869
870
871
872
873 Float_t AliTPCcalibDB::GetPressure(Int_t timeStamp, Int_t run, Int_t type){
874   //
875   // GetPressure for given time stamp and runt
876   //
877   TTimeStamp stamp(timeStamp);
878   AliDCSSensor * sensor = Instance()->GetPressureSensor(run,type);
879   if (!sensor) return 0;
880   
881   if (!sensor->GetFit()) return 0;
882   return sensor->GetValue(stamp);
883 }
884
885 Float_t AliTPCcalibDB::GetValueGoofie(Int_t timeStamp, Int_t run, Int_t type){
886   //
887   // GetPressure for given time stamp and runt
888   //
889   TTimeStamp stamp(timeStamp);
890   AliDCSSensorArray* goofieArray = AliTPCcalibDB::Instance()->GetGoofieSensors(run);
891   if (!goofieArray) return 0;
892   AliDCSSensor *sensor = goofieArray->GetSensor(type);
893   return sensor->GetValue(stamp);
894 }
895
896
897
898
899
900
901 Bool_t  AliTPCcalibDB::GetTemperatureFit(Int_t timeStamp, Int_t run, Int_t side,TVectorD& fit){
902   //
903   //
904   //
905   TTimeStamp tstamp(timeStamp);
906   AliTPCSensorTempArray* tempArray  = Instance()->GetTemperatureSensor(run);
907   if (! tempArray) return kFALSE;
908   AliTPCTempMap * tempMap = new AliTPCTempMap(tempArray);
909   TLinearFitter * fitter = tempMap->GetLinearFitter(3,side,tstamp);
910   if (fitter){
911     fitter->Eval(); 
912     fitter->GetParameters(fit);
913   }
914   delete fitter;
915   delete tempMap;
916   if (!fitter) return kFALSE;
917   return kTRUE;
918 }
919
920 Float_t AliTPCcalibDB::GetTemperature(Int_t timeStamp, Int_t run, Int_t side){
921   //
922   //
923   //
924   TVectorD vec(5);
925   if (side==0) {
926     GetTemperatureFit(timeStamp,run,0,vec);
927     return vec[0];
928   }
929   if (side==1){
930     GetTemperatureFit(timeStamp,run,0,vec);
931     return vec[0];
932   }
933   return 0;
934 }
935
936
937 Double_t AliTPCcalibDB::GetPTRelative(UInt_t timeSec, Int_t run, Int_t side){
938   //
939   // Get relative P/T 
940   // time - absolute time
941   // run  - run number
942   // side - 0 - A side   1-C side
943   AliTPCCalibVdrift * vdrift =  Instance()->GetVdrift(run);
944   if (!vdrift) return 0;
945   return vdrift->GetPTRelative(timeSec,side);
946 }
947
948
949 void AliTPCcalibDB::ProcessEnv(const char * runList){
950   //
951   // Example test function  - how to use the environment variables
952   // runList  -  ascii file with run numbers
953   // output   -  dcsTime.root file with tree 
954
955   ifstream in;
956   in.open(runList);
957   Int_t irun=0;
958   TTreeSRedirector *pcstream = new TTreeSRedirector("dcsTime.root");
959   while(in.good()) {
960     in >> irun;
961     if (irun==0) continue;
962     printf("Processing run %d\n",irun);
963     AliDCSSensor * sensorPressure = AliTPCcalibDB::Instance()->GetPressureSensor(irun);
964     if (!sensorPressure) continue;
965     AliTPCSensorTempArray * tempArray = AliTPCcalibDB::Instance()->GetTemperatureSensor(irun);
966     AliTPCTempMap * tempMap = new AliTPCTempMap(tempArray);
967     AliDCSSensorArray* goofieArray = AliTPCcalibDB::Instance()->GetGoofieSensors(irun);
968     //
969     Int_t startTime = sensorPressure->GetStartTime();
970     Int_t endTime = sensorPressure->GetEndTime();
971     Int_t dtime = TMath::Max((endTime-startTime)/20,10*60);
972     for (Int_t itime=startTime; itime<endTime; itime+=dtime){
973       //
974       TTimeStamp tstamp(itime);
975       Float_t valuePressure = sensorPressure->GetValue(tstamp);
976
977       TLinearFitter * fitter = 0;
978       TVectorD vecTemp[10];
979       if (itime<tempArray->GetStartTime().GetSec() || itime>tempArray->GetEndTime().GetSec()){  
980       }else{
981         for (Int_t itype=0; itype<5; itype++)
982           for (Int_t iside=0; iside<2; iside++){
983             fitter= tempMap->GetLinearFitter(itype,iside,tstamp);
984             if (!fitter) continue;
985             fitter->Eval(); fitter->GetParameters(vecTemp[itype+iside*5]);
986             delete fitter;
987           }
988       }
989       
990       TVectorD vecGoofie, vecEntries, vecMean, vecMedian,vecRMS;
991       if (goofieArray){ 
992         vecGoofie.ResizeTo(goofieArray->NumSensors());
993         ProcessGoofie(goofieArray, vecEntries ,vecMedian, vecMean, vecRMS);
994   //
995         for (Int_t isensor=0; isensor<goofieArray->NumSensors();isensor++){
996           AliDCSSensor *gsensor = goofieArray->GetSensor(isensor);
997           if (gsensor){
998             vecGoofie[isensor] = gsensor->GetValue(tstamp);
999           }
1000         }
1001       }
1002
1003
1004       //tempMap->GetLinearFitter(0,0,itime);
1005       (*pcstream)<<"dcs"<<
1006         "run="<<irun<<
1007         "time="<<itime<<
1008         "goofie.="<<&vecGoofie<<
1009         "goofieE.="<<&vecEntries<<
1010         "goofieMean.="<<&vecMean<<
1011         "goofieMedian.="<<&vecMedian<<
1012         "goofieRMS.="<<&vecRMS<<
1013         "press="<<valuePressure<<
1014         "temp00.="<<&vecTemp[0]<<
1015         "temp10.="<<&vecTemp[1]<<
1016         "temp20.="<<&vecTemp[2]<<
1017         "temp30.="<<&vecTemp[3]<<
1018         "temp40.="<<&vecTemp[4]<<
1019         "temp01.="<<&vecTemp[5]<<
1020         "temp11.="<<&vecTemp[6]<<
1021         "temp21.="<<&vecTemp[7]<<
1022         "temp31.="<<&vecTemp[8]<<
1023         "temp41.="<<&vecTemp[9]<<
1024         "\n";
1025     }
1026   }
1027   delete pcstream;
1028 }
1029
1030
1031 void AliTPCcalibDB::ProcessGoofie( AliDCSSensorArray* goofieArray, TVectorD & vecEntries, TVectorD & vecMedian, TVectorD &vecMean, TVectorD &vecRMS){
1032   /*
1033     
1034   1       TPC_ANODE_I_A00_STAT
1035   2       TPC_DVM_CO2
1036   3       TPC_DVM_DriftVelocity
1037   4       TPC_DVM_FCageHV
1038   5       TPC_DVM_GainFar
1039   6       TPC_DVM_GainNear
1040   7       TPC_DVM_N2
1041   8       TPC_DVM_NumberOfSparks
1042   9       TPC_DVM_PeakAreaFar
1043   10      TPC_DVM_PeakAreaNear
1044   11      TPC_DVM_PeakPosFar
1045   12      TPC_DVM_PeakPosNear
1046   13      TPC_DVM_PickupHV
1047   14      TPC_DVM_Pressure
1048   15      TPC_DVM_T1_Over_P
1049   16      TPC_DVM_T2_Over_P
1050   17      TPC_DVM_T_Over_P
1051   18      TPC_DVM_TemperatureS1
1052    */
1053   //
1054   //
1055   //  TVectorD  vecMedian; TVectorD  vecEntries; TVectorD  vecMean; TVectorD  vecRMS;
1056   Double_t kEpsilon=0.0000000001;
1057   Double_t kBig=100000000000.;
1058   Int_t nsensors = goofieArray->NumSensors();
1059   vecEntries.ResizeTo(nsensors);
1060   vecMedian.ResizeTo(nsensors);
1061   vecMean.ResizeTo(nsensors);
1062   vecRMS.ResizeTo(nsensors);
1063   TVectorF values;
1064   for (Int_t isensor=0; isensor<goofieArray->NumSensors();isensor++){
1065     AliDCSSensor *gsensor = goofieArray->GetSensor(isensor);
1066     if (gsensor &&  gsensor->GetGraph()){
1067       Int_t npoints = gsensor->GetGraph()->GetN();
1068       // filter zeroes
1069       values.ResizeTo(npoints);
1070       Int_t nused =0;
1071       for (Int_t ipoint=0; ipoint<npoints; ipoint++){
1072         if (TMath::Abs(gsensor->GetGraph()->GetY()[ipoint])>kEpsilon && 
1073            TMath::Abs(gsensor->GetGraph()->GetY()[ipoint])<kBig ){
1074           values[nused]=gsensor->GetGraph()->GetY()[ipoint];
1075           nused++;
1076         }
1077       }
1078       //
1079       vecEntries[isensor]= nused;      
1080       if (nused>1){
1081         vecMedian[isensor] = TMath::Median(nused,values.GetMatrixArray());
1082         vecMean[isensor]   = TMath::Mean(nused,values.GetMatrixArray());
1083         vecRMS[isensor]    = TMath::RMS(nused,values.GetMatrixArray());
1084       }
1085     }
1086   }
1087 }
1088