AliTPCcalibDButil.cxx.diff Add Time0 creation using a combination of CE and...
[u/mrichter/AliRoot.git] / TPC / CalibMacros / CalibEnv.C
1 /*
2 // .x ~/NimStyle.C
3 // .x ~/rootlogon.C
4
5 gSystem->AddIncludePath("-I$ALICE_ROOT/TPC");
6
7
8
9 .L $ALICE_ROOT/TPC/CalibMacros/CalibEnv.C+
10 Init();
11 CalibEnv("run.list");
12 GetTree();
13
14 TFile f("dcsTime.root")
15
16 //
17 // if you want to use alien OCDB
18 // 
19 gSystem->Load("libXrdClient.so");
20 gSystem->Load("libNetx.so");
21 if (!gGrid) TGrid::Connect("alien://",0,0,"t");
22
23
24
25 */
26
27 #include <iostream>
28 #include <fstream>
29 #include <stdio.h>
30 #include <AliCDBManager.h>
31 #include <AliCDBEntry.h>
32 #include <AliLog.h>
33 #include <AliMagF.h>
34 #include "AliTPCcalibDB.h"
35 #include "AliTPCcalibDButil.h"
36 #include "AliTPCAltroMapping.h"
37 #include "AliTPCExB.h"
38 #include "AliTPCCalROC.h"
39 #include "AliTPCCalPad.h"
40 #include "AliTPCSensorTempArray.h"
41 #include "AliGRPObject.h"
42 #include "AliTPCTransform.h"
43 #include "TFile.h"
44 #include "TKey.h"
45 #include "TObjArray.h"
46 #include "TObjString.h"
47 #include "TString.h"
48 #include "AliTPCCalPad.h"
49 #include "AliTPCROC.h"
50 #include "AliTPCParam.h"
51 #include "AliTPCCalibPulser.h"
52 #include "AliTPCCalibPedestal.h"
53 #include "AliTPCCalibCE.h"
54 #include "AliTPCExBFirst.h"
55 #include "TTreeStream.h"
56 #include "AliTPCTempMap.h"
57 #include "TVectorD.h"
58 #include "TMatrixD.h"
59 #include "AliTPCCalibRaw.h"
60 #include "AliSplineFit.h"
61
62
63 TTree * dcsTree=0;
64 TString refFile="/data/Work/data/calib/guiTrees/RefCalPads_83680.root";
65 // TString refFile="/lustre/alice/wiechula/calib/guiTrees/RefCalPads_83680.root"
66 void GetProductionInfo(Int_t run, Int_t &nalien, Int_t &nRawAlien, Int_t &nlocal, Int_t &nRawLocal);
67   
68 void Init2008(){
69   //
70   //
71   //
72   AliCDBManager::Instance()->SetDefaultStorage("local://$ALICE_ROOT/OCDB");
73   AliCDBManager::Instance()->SetSpecificStorage("TPC/Calib/Parameters","local://$ALICE_ROOT/OCDB");
74   AliCDBManager::Instance()->SetSpecificStorage("GRP/GRP/Data","local:///lustre/alice/alien/alice/data/2008/LHC08d/OCDB/");
75   AliCDBManager::Instance()->SetSpecificStorage("TPC/Calib/Temperature","local:///lustre/alice/alien/alice/data/2008/LHC08d/OCDB/");
76   AliCDBManager::Instance()->SetSpecificStorage("TPC/Calib/HighVoltage","local:///lustre/alice/alien/alice/data/2008/LHC08d/OCDB/");
77   AliCDBManager::Instance()->SetSpecificStorage("TPC/Calib/Goofie","local:///lustre/alice/alien/alice/data/2008/LHC08d/OCDB/");
78   AliCDBManager::Instance()->SetRun(1);
79 }
80
81 void Init(){  
82   AliCDBManager::Instance()->SetDefaultStorage("local:///lustre/alice/alien/alice/data/2009/OCDB/");
83   AliCDBManager::Instance()->SetRun(1);
84 }
85
86 void InitAlien(const char *path="LHC08b"){
87   //
88   //
89   //
90   TString alpath="alien://folder=/alice/data/2008/";
91   alpath+=path;
92   alpath+="/OCDB";
93   
94   AliCDBManager::Instance()->SetDefaultStorage("local://$ALICE_ROOT/OCDB");
95   AliCDBManager::Instance()->SetSpecificStorage("TPC/Calib/Parameters","local://$ALICE_ROOT/OCDB");
96   AliCDBManager::Instance()->SetSpecificStorage("GRP/GRP/Data",alpath.Data());
97   AliCDBManager::Instance()->SetSpecificStorage("TPC/Calib/Temperature",alpath.Data());
98   AliCDBManager::Instance()->SetSpecificStorage("TPC/Calib/Goofie",alpath.Data());
99   AliCDBManager::Instance()->SetRun(1);
100 }
101
102
103 void CalibEnv(const char * runList){
104   //
105   //
106   //
107   AliTPCcalibDB * calibDB = AliTPCcalibDB::Instance();
108   ifstream in;
109   in.open(runList);
110   Int_t irun=0;
111   TTreeSRedirector *pcstream = new TTreeSRedirector("dcsTime.root");
112   AliTPCcalibDButil dbutil;
113   dbutil.SetRefFile(refFile.Data());
114   Int_t startTime = 0;
115   Int_t endTime = 0;
116   Int_t startTimeGRP = 0;
117   Int_t stopTimeGRP  = 0;
118
119   AliSplineFit *fitVdrift=0x0;
120   //  for (Int_t irun=startRun; irun<stopRun; irun++){
121   while(in.good()) {
122     in >> irun;
123     if (in.eof()) break;
124     if (irun==0) continue;
125     printf("Processing run %d ...\n",irun);
126     AliTPCcalibDB::Instance()->SetRun(irun);
127     dbutil.UpdateFromCalibDB();
128     AliDCSSensorArray *arrHV=calibDB->GetVoltageSensors(irun);
129     if (!arrHV) continue;
130     for  (Int_t isenHV=0; isenHV<arrHV->NumSensors(); ++isenHV){
131       AliDCSSensor *senHV=arrHV->GetSensorNum(isenHV);
132       if (!senHV) continue;
133       startTime=senHV->GetStartTime();
134       endTime  =senHV->GetEndTime();
135       if (startTime>0&&endTime>0) break;
136     }
137     if (calibDB->GetGRP(irun)){
138       startTimeGRP = AliTPCcalibDB::GetGRP(irun)->GetTimeStart();
139       stopTimeGRP  = AliTPCcalibDB::GetGRP(irun)->GetTimeEnd();
140     }
141     //    AliDCSSensor * sensorPressure = AliTPCcalibDB::Instance()->GetPressureSensor(irun);
142 //     if (!sensorPressure) continue;
143 //     Int_t startTime = sensorPressure->GetStartTime();
144 //     Int_t endTime = sensorPressure->GetEndTime();
145 //     Int_t startTimeGRP = AliTPCcalibDB::GetGRP(irun)->GetTimeStart();
146 //     Int_t stopTimeGRP  = AliTPCcalibDB::GetGRP(irun)->GetTimeEnd();
147     AliTPCSensorTempArray * tempArray = AliTPCcalibDB::Instance()->GetTemperatureSensor(irun);
148     AliTPCTempMap * tempMap = new AliTPCTempMap(tempArray);
149     AliDCSSensorArray* goofieArray = AliTPCcalibDB::Instance()->GetGoofieSensors(irun);
150     //
151     Int_t dtime = TMath::Max((endTime-startTime)/20,10*60);
152     //Goofie statistical data
153     TVectorD vecEntries, vecMean, vecMedian,vecRMS;
154     dbutil.ProcessGoofie(vecEntries ,vecMedian, vecMean, vecRMS);
155     //CE data processing - see ProcessCEdata function for description of the results
156     TVectorD fitResultsA, fitResultsC;
157     Int_t nmaskedCE;
158 //     dbutil.ProcessCEdata("(sector<36)++gx++gy++lx++lx**2",fitResultsA,fitResultsC,nmaskedCE);
159     dbutil.ProcessCEdata("(sector<36)++gy++gx++(lx-134)++(sector<36)*(lx-134)",fitResultsA,fitResultsC,nmaskedCE);
160     TVectorD vecTEntries, vecTMean, vecTRMS, vecTMedian, vecQEntries, vecQMean, vecQRMS, vecQMedian;
161     Float_t driftTimeA, driftTimeC;
162     dbutil.ProcessCEgraphs(vecTEntries, vecTMean, vecTRMS, vecTMedian,
163                            vecQEntries, vecQMean, vecQRMS, vecQMedian,
164                            driftTimeA, driftTimeC );
165     //drift velocity using tracks
166     fitVdrift=calibDB->GetVdriftSplineFit("ALISPLINEFIT_MEAN_VDRIFT_COSMICS_ALL",irun);
167     //noise data Processing - see ProcessNoiseData function for description of the results
168     TVectorD vNoiseMean, vNoiseMeanSenRegions, vNoiseRMS, vNoiseRMSSenRegions;
169     Int_t nonMaskedZero=0;
170     dbutil.ProcessNoiseData(vNoiseMean, vNoiseMeanSenRegions, vNoiseRMS, vNoiseRMSSenRegions, nonMaskedZero);
171     // comparisons
172     TVectorF pedestalDeviations;
173     TVectorF noiseDeviations;
174     TVectorF pulserQdeviations;
175     Float_t varQMean;
176     Int_t npadsOutOneTB;
177     Int_t npadsOffAdd;
178     dbutil.ProcessPedestalVariations(pedestalDeviations);
179     dbutil.ProcessNoiseVariations(noiseDeviations);
180     dbutil.ProcessPulserVariations(pulserQdeviations,varQMean,npadsOutOneTB,npadsOffAdd);
181     
182     //L3 data 
183     Float_t bz=AliTPCcalibDB::GetBz(irun);
184     Char_t  l3pol=AliTPCcalibDB::GetL3Polarity(irun);
185     //calibration Pulser data processing
186     Int_t nOffChannels=0;
187     TVectorD vTimePulser;
188     nOffChannels=dbutil.GetNPulserOutliers();
189     dbutil.ProcessPulser(vTimePulser);
190     //ALTRO data
191     Int_t nMasked=0;
192     dbutil.ProcessALTROConfig(nMasked);
193     //Calib RAW data
194     Int_t nFailL1=-1;
195     if (calibDB->GetCalibRaw())
196       nFailL1=calibDB->GetCalibRaw()->GetNFailL1Phase();
197     //production information
198     Int_t nalien=0,nRawAlien=0,nlocal=0,nRawLocal=0;
199 //     GetProductionInfo(irun, nalien, nRawAlien, nlocal,nRawLocal);
200     //run type
201     TObjString runType(AliTPCcalibDB::GetRunType(irun).Data());
202     
203     for (Int_t itime=startTime; itime<endTime; itime+=dtime){
204       //
205       TTimeStamp tstamp(itime);
206       Float_t valuePressure  = calibDB->GetPressure(tstamp,irun,0);
207       Float_t valuePressure2 = calibDB->GetPressure(tstamp,irun,1);
208       //temperature fits
209       TLinearFitter * fitter = 0;
210       TVectorD vecTemp[10];
211       for (Int_t itype=0; itype<5; itype++)
212         for (Int_t iside=0; iside<2; iside++){
213           fitter= tempMap->GetLinearFitter(itype,iside,tstamp);   
214           if (!fitter) continue;
215           fitter->Eval();
216           fitter->GetParameters(vecTemp[itype+iside*5]);
217           delete fitter;
218         }
219       //measured skirt temperatures
220       TVectorD vecSkirtTempA(18);
221       TVectorD vecSkirtTempC(18);
222       Int_t nsenTemp=tempArray->NumSensors();
223       for (Int_t isenTemp=0;isenTemp<nsenTemp;++isenTemp){
224         AliTPCSensorTemp *senTemp=(AliTPCSensorTemp*)tempArray->GetSensorNum(isenTemp);
225         if (senTemp->GetType()!=3) continue;
226         if (TMath::Sqrt(senTemp->GetX()*senTemp->GetX()+senTemp->GetY()*senTemp->GetY())<100) continue; //only skirt, outer FC vessel
227         Double_t val=senTemp->GetValue(tstamp);
228         if (senTemp->GetSide()==0)
229           vecSkirtTempA[senTemp->GetSector()]=val;
230         else
231           vecSkirtTempC[senTemp->GetSector()]=val;
232       }
233       //goofie data
234       TVectorD vecGoofie;
235       if (goofieArray){
236         vecGoofie.ResizeTo(goofieArray->NumSensors());
237         for (Int_t isensor=0; isensor<goofieArray->NumSensors();isensor++){
238           AliDCSSensor *gsensor = goofieArray->GetSensor(isensor);
239           if (gsensor){
240             vecGoofie[isensor] = gsensor->GetValue(tstamp);
241           }
242         }
243       } else {
244         vecGoofie.ResizeTo(19);
245       }
246       Double_t ptrelative0 = AliTPCcalibDB::GetPTRelative((UInt_t)itime,irun,0);
247       Double_t ptrelative1 = AliTPCcalibDB::GetPTRelative((UInt_t)itime,irun,1);
248       //
249       TVectorD voltagesIROC(36);
250       TVectorD voltagesOROC(36);
251       for(Int_t j=1; j<36; j++) voltagesIROC[j-1] = AliTPCcalibDB::Instance()->GetChamberHighVoltage(irun, j,itime);
252       for(Int_t j=36; j<72; j++) voltagesOROC[j-36] = AliTPCcalibDB::Instance()->GetChamberHighVoltage(irun, j,itime);
253       Double_t voltIROC = TMath::Median(36, voltagesIROC.GetMatrixArray());
254       Double_t voltOROC = TMath::Median(36, voltagesOROC.GetMatrixArray());
255       //
256       Float_t  coverIA=AliTPCcalibDB::GetCoverVoltage(irun,0,itime);
257       Float_t  coverIC=AliTPCcalibDB::GetCoverVoltage(irun,18,itime);
258       Float_t  coverOA=AliTPCcalibDB::GetCoverVoltage(irun,36,itime);
259       Float_t  coverOC=AliTPCcalibDB::GetCoverVoltage(irun,54,itime);
260       Float_t  skirtA=AliTPCcalibDB::GetSkirtVoltage(irun,0,itime);
261       Float_t  skirtC=AliTPCcalibDB::GetSkirtVoltage(irun,18,itime);
262       Float_t  ggOffA=AliTPCcalibDB::GetGGoffsetVoltage(irun,0,itime);
263       Float_t  ggOffC=AliTPCcalibDB::GetGGoffsetVoltage(irun,18,itime);
264       //drift velocity
265       Float_t dvCorr=-5;
266       if (fitVdrift) dvCorr=fitVdrift->Eval(itime);
267       
268       
269       //tempMap->GetLinearFitter(0,0,itime);
270       (*pcstream)<<"dcs"<<
271         "run="<<irun<<
272         "time="<<itime<<
273         "startTimeGRP="<<startTimeGRP<<
274         "stopTimeGRP="<<stopTimeGRP<<
275         //run type
276         "runType.="<<&runType<<
277         // voltage setting
278         "VIROC.="<<&voltagesIROC<<
279         "VOROC.="<<&voltagesOROC<<
280         "medianVIROC="<<voltIROC<<
281         "medianVOROC="<<voltOROC<<
282         "coverIA=" << coverIA <<
283         "coverIC=" << coverIC <<
284         "coverOA=" << coverOA <<
285         "coverOC=" << coverOC <<
286         "skirtA=" << skirtA <<
287         "skirtC=" << skirtC <<
288         "ggOffA=" << ggOffA <<
289         "ggOffC=" << ggOffC <<
290         //
291         "ptrel0="<<ptrelative0<<  // deltaTP/TP  - A side
292         "ptrel1="<<ptrelative1<<  // deltaTP/TPC - C side
293         "goofie.="<<&vecGoofie<<
294         "goofieE.="<<&vecEntries<<
295         "goofieMean.="<<&vecMean<<
296         "goofieMedian.="<<&vecMedian<<
297         "goofieRMS.="<<&vecRMS<<
298         //
299         "press="<<valuePressure<<
300         "press2="<<valuePressure2<<
301         "temp00.="<<&vecTemp[0]<<
302         "temp10.="<<&vecTemp[1]<<
303         "temp20.="<<&vecTemp[2]<<
304         "temp30.="<<&vecTemp[3]<<
305         "temp40.="<<&vecTemp[4]<<
306         "temp01.="<<&vecTemp[5]<<
307         "temp11.="<<&vecTemp[6]<<
308         "temp21.="<<&vecTemp[7]<<
309         "temp31.="<<&vecTemp[8]<<
310         "temp41.="<<&vecTemp[9]<<
311         "tempSkirtA.="<<&vecSkirtTempA<<
312         "tempSkirtC.="<<&vecSkirtTempC;
313       
314       (*pcstream)<<"dcs"<<      
315         //noise data
316         "meanNoise.="<<&vNoiseMean<<
317         "meanNoiseSen.="<<&vNoiseMeanSenRegions<<
318         "rmsNoise.="<<&vNoiseRMS<<
319         "rmsNoiseSen.="<<&vNoiseRMSSenRegions<<
320         "zeroNoise="<<nonMaskedZero<<
321         //pulser data
322         "timePulser.=" << &vTimePulser <<
323         "nOffPulser="<<nOffChannels<<
324         //altro data
325         "nMasked="<<nMasked<<
326         //ce data
327         "CEfitA.="<<&fitResultsA<<
328         "CEfitC.="<<&fitResultsC<<
329         "nmaskedCE="<<nmaskedCE<<
330         //ce graph data
331         "CEgrTEntries.="<<&vecTEntries<<
332         "CEgrTMean.="<<&vecTMean<<
333         "CEgrTRMS.="<<&vecTRMS<<
334         "CEgrTMedian.="<<&vecTMedian<<
335         "CEgrQEntries.="<<&vecQEntries<<
336         "CEgrQMean.="<<&vecQMean<<
337         "CEgrQRMS.="<<&vecQRMS<<
338         "CEgrQMedian.="<<&vecQMedian<<
339         "CEgrDriftA="<<driftTimeA<<
340         "CEgrDriftC="<<driftTimeC<<
341         //calib raw data
342         "nFailL1="<<nFailL1<<
343         // b field
344         "Bz="<< bz <<
345         "L3polarity="<<l3pol<<
346         // production information
347         "nalien="<<nalien<<
348         "nRawAlien="<<nRawAlien<<
349         "nlocal="<<nlocal<<
350         "nRawLocal="<<nRawLocal<<
351         //comparisons with ref data
352         "pedestalDeviations.="<<&pedestalDeviations<<
353         "noiseDeviations.="<<&noiseDeviations<<
354         "pulserQdeviations.="<<&pulserQdeviations<<
355 //         "pulserVarQMean="<<varQMean<<
356         "pulserNpadsOutOneTB="<<npadsOutOneTB<<
357         "pulserNpadsOffAdd="<<npadsOffAdd<<
358         "driftCorrCosmAll="<<dvCorr<<
359       
360         "\n";
361     }
362   }
363   delete pcstream;
364 }
365
366
367
368
369 void GetProductionInfo(Int_t run, Int_t &nalien, Int_t &nRawAlien, Int_t &nlocal, Int_t &nRawLocal){
370   //
371   // find number of ESDs in central and local production for run
372   //
373
374   nalien=0;
375   nRawAlien=0;
376   nlocal=0;
377   nRawLocal=0;
378   TString sNlines;
379   //find number of ESDs in alien
380   TString command="alien_find /alice/data/2009 ";
381   command += Form("%09d",run);
382   command += " | grep AliESDs.root | wc -l";
383   sNlines = gSystem->GetFromPipe(command.Data());
384   nalien=sNlines.Atoi();
385   //find number of raw files on alien
386   command="alien_find /alice/data/2009 ";
387   command += Form("%09d",run);
388   command += " | grep raw | grep -v tag | wc -l";
389   sNlines = gSystem->GetFromPipe(command.Data());
390   nRawAlien=sNlines.Atoi();
391   //find number of ESDs local
392   command="find /lustre/alice/alien/alice/data/2009 -name AliESDs.root | grep ";
393   command += Form("%09d",run);
394   command += " | wc -l";
395   sNlines = gSystem->GetFromPipe(command.Data());
396   nlocal=sNlines.Atoi();
397   //find number of local raw data files
398   command="find /lustre/alice/alien/alice/data/2009 -name \"*.root\" | grep ";
399   command += Form("%09d",run);
400   command += " | grep raw | grep -v tag | wc -l";
401   sNlines = gSystem->GetFromPipe(command.Data());
402   nRawLocal=sNlines.Atoi();
403 }
404
405 void FilterMag(const char * runList){
406   //
407   //
408   //
409   //  AliTPCcalibDB * calibDB = AliTPCcalibDB::Instance();
410   ifstream in;
411   in.open(runList);
412   Int_t irun=0;
413   while(in.good()) {
414     in >> irun;
415     if (irun==0) continue;
416     AliGRPObject *grp = AliTPCcalibDB::GetGRP(irun);
417     Float_t current = -1;
418     Float_t bz      = -1;
419 //     Float_t press   =  0;
420     if (grp){
421       current = grp->GetL3Current((AliGRPObject::Stats)0);
422       bz = 5*current/30000.;
423       printf("Run%d\tL3 current%f\tBz\t%f\n",irun,current,bz);
424     }
425     else{
426       printf("Run%d\tL3 current%f\tBz\t%f\n",irun,current,bz);
427     }
428   }
429   
430 }
431
432
433 void GetTree(){
434   TFile *fdcs = new TFile("dcsTime.root");
435   dcsTree  = (TTree*)fdcs->Get("dcs");
436   //
437   // mean temp A
438   
439   dcsTree->Draw("temp30.fElements[0]");
440   
441 }
442
443 void GetNominalValues(){
444   //
445   if (!dcsTree) return;
446 }
447
448
449
450
451 /*
452
453 AliDCSSensor * sensorPressure = AliTPCcalibDB::Instance()->GetPressureSensor(62084);
454 entry = AliCDBManager::Instance()->Get("TPC/Calib/Temperature",run);
455 AliTPCSensorTempArray * tempArray = (AliTPCSensorTempArray *)entry->GetObject();
456 AliTPCSensorTempArray * tempArray = (AliTPCSensorTempArray *)AliTPCcalibDB::Instance()->GetTemperatureSensor(62084)
457 AliTPCTempMap * tempMap = new AliTPCTempMap(tempArray);
458 TLinearFitter * fitter = tempMap->GetLinearFitter(0,0,tempArray->GetStartTime());
459
460 AliDCSSensorArray* goofieArray = AliTPCcalibDB::Instance()->GetGoofieSensors(62084);
461
462 */
463
464 /*
465
466 void PlotPressureResol(){
467   //
468   // Example
469   //
470   dcs->Draw("100000*(press-press2-4.782)/press/sqrt(2.)>>his(100,-50,50)","run>61400","")
471   his->SetXTitle("#sigma_{P/P_{0}}(10^{-5})");
472   gPad->SaveAs("picDCS/deltaPoverP.eps");
473   gPad->SaveAs("picDCS/deltaPoverP.gif");
474
475 }
476 void PlotTresol(){
477   //
478   // T resolution example
479   // plot difference of the temperature from A and C side
480   // Supposing the error is independent - (division by sqrt(2))
481   dcs->Draw("100000*(temp30.fElements[0]-temp31.fElements[0]+0.00509)/(temp31.fElements[0]+273.15)/sqrt(2.)>>his(100,-5,5)","run>61400","");
482   his->SetXTitle("#sigma_{T/T_{0}}(10^{-5})");
483   gPad->SaveAs("picDCS/deltaToverT.eps");
484   gPad->SaveAs("picDCS/deltaToverT.gif");
485 }
486 */