8601b40ee1c4f5abb2896fcd2eb76070800aa85c
[u/mrichter/AliRoot.git] / ITS / macrosSDD / ShowResponseSDD.C
1 #if !defined(__CINT__) || defined(__MAKECINT__)
2 #include <TFile.h>
3 #include <TH1F.h>
4 #include <TH2I.h>
5 #include <TGraph.h>
6 #include <TExec.h>
7 #include <TStyle.h>
8 #include <TString.h>
9 #include <TSystem.h>
10 #include <TGrid.h>
11 #include <TLatex.h>
12 #include <TCanvas.h>
13 #include <TObjArray.h>
14 #include "AliCDBEntry.h"
15 #include "AliITSresponseSDD.h"
16 #endif
17
18 /*  $Id: ShowResponseSDD.C 44072 2010-10-04 15:48:32Z prino $    */
19
20 // Macro to plot the calibration parameters 
21 // (time zero, ADCtokeV and Vdrift correction)
22 // from the OCDB file created from SDD offline calibration 
23 // (OCDB/ITS/Calib/RespSDD)
24
25 void ShowResponseSDD(TString filename="$ALICE_ROOT/OCDB/ITS/Calib/RespSDD/Run0_999999999_v0_s0.root"){
26   if(filename.Contains("alien")){
27     TGrid::Connect("alien:");
28   }
29   TFile* fr=TFile::Open(filename.Data());
30   AliCDBEntry* e=(AliCDBEntry*)fr->Get("AliCDBEntry");
31   AliITSresponseSDD* r=(AliITSresponseSDD*)e->GetObject();
32   TH1F* hTimeZero=new TH1F("hTimeZero","",260,239.5,499.5);
33   TH1F* hVdriftCorrLeft=new TH1F("hVdriftCorrLeft","",260,239.5,499.5);
34   TH1F* hVdriftCorrRight=new TH1F("hVdriftCorrRight","",260,239.5,499.5);
35   TH1F* hADCtokeV=new TH1F("hADCtokeV","",260,239.5,499.5);
36   Float_t averTz=0.;
37   Float_t averCv0=0.;
38   Float_t averCv1=0.;
39   Float_t averAk=0.;
40   for(Int_t iMod=240; iMod<500; iMod++){
41     Float_t tz=r->GetTimeZero(iMod);
42     Float_t cv0=r->GetDeltaVDrift(iMod,kFALSE);
43     Float_t cv1=r->GetDeltaVDrift(iMod,kTRUE);
44     Float_t ak=r->GetADCtokeV(iMod);
45     hTimeZero->SetBinContent(iMod-240+1,tz);
46     hVdriftCorrLeft->SetBinContent(iMod-240+1,cv0);
47     hVdriftCorrRight->SetBinContent(iMod-240+1,cv1);
48     hADCtokeV->SetBinContent(iMod-240+1,ak);
49     averTz+=tz;
50     averCv0+=cv0;
51     averCv1+=cv1;
52     averAk+=ak;
53   }
54   averTz/=260.;
55   averCv0/=260.;
56   averCv1/=260.;
57   averAk/=260.;
58
59   hTimeZero->SetMarkerStyle(20);
60   hADCtokeV->SetMarkerStyle(20);
61   hVdriftCorrLeft->SetMarkerStyle(22);
62   hVdriftCorrLeft->SetMarkerColor(2);
63   hVdriftCorrRight->SetMarkerStyle(24);
64   hVdriftCorrRight->SetMarkerColor(4);
65   hTimeZero->SetMaximum(hTimeZero->GetMaximum()*1.2);
66   hADCtokeV->SetMaximum(hADCtokeV->GetMaximum()*1.2);
67   Float_t maxV=TMath::Max(hVdriftCorrLeft->GetMaximum(),hVdriftCorrRight->GetMaximum());
68   Float_t minV=TMath::Min(hVdriftCorrLeft->GetMinimum(),hVdriftCorrRight->GetMinimum());
69   Float_t scale=TMath::Max(TMath::Abs(maxV),TMath::Abs(minV));
70   if(scale<0.02){
71     hVdriftCorrLeft->SetMinimum(-0.05);
72     hVdriftCorrLeft->SetMaximum(0.05);
73   }else{
74     hVdriftCorrLeft->SetMaximum(1.2*scale);
75     hVdriftCorrLeft->SetMinimum(-1.2*scale);
76   }
77
78   gStyle->SetOptStat(0);
79
80   printf("Charge vs. Time correction factor = %f\n",r->GetChargevsTime());
81
82   TCanvas* c1=new TCanvas("c1","Time Zero");
83   hTimeZero->Draw("P");
84   hTimeZero->GetXaxis()->SetTitle("Module Id");
85   hTimeZero->GetYaxis()->SetTitle("Time Zero (ns)");
86   TLatex* tt=new TLatex(0.2,0.8,Form("Average Tzero = %.2f",averTz));
87   tt->SetNDC();
88   tt->Draw();
89   c1->Modified();
90
91   TCanvas* c2=new TCanvas("c2","Vdrift Corr");
92   hVdriftCorrLeft->Draw("P");
93   hVdriftCorrRight->Draw("PSAME");
94   hVdriftCorrLeft->GetXaxis()->SetTitle("Module Id");
95   hVdriftCorrLeft->GetYaxis()->SetTitle("Vdrift correction (#mum/ns)");
96   TLatex* tc0=new TLatex(0.2,0.8,Form("Average Vdrift corr Left = %.3f",averCv0));
97   tc0->SetNDC();
98   tc0->SetTextColor(2);
99   tc0->Draw();
100   TLatex* tc1=new TLatex(0.2,0.72,Form("Average Vdrift corr Right = %.3f",averCv1));
101   tc1->SetNDC();
102   tc1->SetTextColor(4);
103   tc1->Draw();
104   c2->Modified();
105
106   TCanvas* c3=new TCanvas("c3","ADC calib");
107   hADCtokeV->Draw("P");
108   hADCtokeV->GetXaxis()->SetTitle("Module Id");
109   hADCtokeV->GetYaxis()->SetTitle("ADC to keV conv. factor");
110   TLatex* ta=new TLatex(0.2,0.8,Form("Average ADCtokeV = %.3f",averAk));
111   ta->SetNDC();
112   ta->Draw(); 
113   TLatex* tb=new TLatex(0.2,0.72,Form("Charge vs. Time correction factor = %f\n",r->GetChargevsTime()));
114   tb->SetNDC();
115   tb->SetTextColor(kGreen+1);
116   tb->Draw();
117   c3->Modified();
118
119 }