Changes related to the extraction of the V0 finder into a separate class (A. Dainese...
[u/mrichter/AliRoot.git] / ITS / AnalyzeSDDInjectorsAllMod.C
1 #if !defined(__CINT__) || defined(__MAKECINT__)
2 #include <TH2F.h>
3 #include <TCanvas.h>
4 #include <TGraphErrors.h>
5 #include <TStopwatch.h>
6 #include <TStyle.h>
7 #include <TLatex.h>
8 #include <TFile.h>
9 #include <TMath.h>
10 #include <TGrid.h>
11 #include <TF1.h>
12 #include <TLine.h>
13 #include "AliRawReader.h"
14 #include "AliRawReaderDate.h"
15 #include "AliRawReaderRoot.h"
16 #include "AliITSOnlineSDDInjectors.h"
17 #include "AliITSRawStreamSDD.h"
18 #include "AliITSRawStreamSDDCompressed.h"
19 #include "AliITSDDLModuleMapSDD.h"
20 #endif
21
22 // Macro for the analysis of PULSER runs (equivalent to ITSSDDINJda.cxx)
23 // Two functions named AnalyzeSDDInjectorsAllModules: 
24 // The first is for analyzing a local raw data file and takes as agrument the file name.
25 // The second is for running on ALIEN
26 // All DDLs are analyzed, the argument nDDL selects the DDL to be plotted
27 // Origin: F. Prino (prino@to.infn.it)
28
29
30 void AnalyzeSDDInjectorsAllMod(Char_t *datafil, 
31                                Int_t adcfreq=20, 
32                                Int_t nDDL=0, 
33                                Int_t firstEv=10, 
34                                Int_t lastEv=15,
35                                Int_t jpad=20, 
36                                Int_t statuscut=7, 
37                                Int_t dataformat=1){
38
39   const Int_t kTotDDL=24;
40   const Int_t kModPerDDL=12;
41   const Int_t kSides=2;
42
43   AliITSDDLModuleMapSDD* dmap=new AliITSDDLModuleMapSDD();
44   dmap->SetJun08Map();
45
46   TH2F** histo = new TH2F*[kTotDDL*kModPerDDL*kSides];
47   Int_t nWrittenEv[kTotDDL*kModPerDDL*kSides];
48   TGraphErrors** gvel = new TGraphErrors*[kTotDDL*kModPerDDL*kSides];
49   AliITSOnlineSDDInjectors **anal=new AliITSOnlineSDDInjectors*[kTotDDL*kModPerDDL*kSides];
50   TH1F** hvdriftl=new TH1F*[260];  
51   TH1F** hvdriftr=new TH1F*[260];  
52   Char_t hisnam[20];
53   for(Int_t idet=0; idet<260;idet++){
54     sprintf(hisnam,"vdriftl%03d",idet);
55     hvdriftl[idet]=new TH1F(hisnam,"",500,5.5,8.0);
56     sprintf(hisnam,"vdriftr%03d",idet);
57     hvdriftr[idet]=new TH1F(hisnam,"",500,5.5,8.0);
58   }
59   for(Int_t iddl=0; iddl<kTotDDL;iddl++){
60     for(Int_t imod=0; imod<kModPerDDL;imod++){
61       for(Int_t isid=0;isid<kSides;isid++){
62         Int_t index=kSides*(kModPerDDL*iddl+imod)+isid;
63         sprintf(hisnam,"h%02dc%02ds%d",iddl,imod,isid);
64         histo[index]=new TH2F(hisnam,"",256,-0.5,255.5,256,-0.5,255.5);
65         anal[index]=new AliITSOnlineSDDInjectors(iddl,imod,isid);
66         if(adcfreq==40) anal[index]->Set40MHzConfig();
67         else anal[index]->Set20MHzConfig();
68         nWrittenEv[index]=0;
69       }
70     }
71   }
72   TGraphErrors *gvvsmod0=new TGraphErrors(0);
73   TGraphErrors *gvvsmod1=new TGraphErrors(0);
74   TGraphErrors *gtvsmod0=new TGraphErrors(0);
75   TGraphErrors *gtvsmod1=new TGraphErrors(0);
76   Float_t gvmin=6.0, gvmax=7.5;
77   Float_t gtmin=288., gtmax=308.;
78   TH1F* hanst=new TH1F("hanst","",8,-0.5,7.5);
79   TH1F* hpad7l=new TH1F("hpad7l","",33,-0.5,32.5);
80   TH1F* hpad7r=new TH1F("hpad7r","",33,-0.5,32.5);
81
82   TCanvas* c0 = new TCanvas("c0","Event display",900,900);
83   gStyle->SetPalette(1);
84   TCanvas* c1 = new TCanvas("c1","Drift Speed vs. anode",900,900);
85   Char_t text[50];
86
87   Int_t iev=firstEv;
88   AliRawReader *rd; 
89   if(strstr(datafil,".root")!=0){
90     rd=new AliRawReaderRoot(datafil,iev);
91   }else{
92     rd=new AliRawReaderDate(datafil,iev);
93   }
94
95   Char_t gname[15];
96   TF1 *funz=new TF1("funz","[0]+[1]*x+[2]*x*x+[3]*x*x*x",0.,255.);
97   TLatex *t0=new TLatex();
98   t0->SetNDC();
99   t0->SetTextSize(0.06);
100   t0->SetTextColor(4);
101   Int_t readEv=0;
102   do{
103     c0->Clear();
104     c0->Divide(4,6,0.001,0.001);
105     c1->Clear();
106     c1->Divide(4,6,0.001,0.001);
107     printf("Event # %d\n",iev);
108     UInt_t timeSt=rd->GetTimestamp();
109     rd->Reset();
110     for(Int_t iddl=0; iddl<kTotDDL;iddl++){
111       for(Int_t imod=0; imod<kModPerDDL;imod++){
112         for(Int_t isid=0;isid<kSides;isid++){
113           Int_t index=kSides*(kModPerDDL*iddl+imod)+isid;
114           histo[index]->Reset();
115         }
116       }
117     }
118
119     AliITSRawStream* s;
120     if(dataformat==0){
121       s=new AliITSRawStreamSDD(rd);
122     }else{
123       s=new AliITSRawStreamSDDCompressed(rd);
124       if(dataformat==1) s->SetADCEncoded(kTRUE);
125     }
126     while(s->Next()){
127       Int_t iDDL=rd->GetDDLID();
128       Int_t iCarlos=s->GetCarlosId();
129       if(s->IsCompletedModule()) continue;
130       if(s->IsCompletedDDL()) continue;
131       if(iDDL>=0 && iDDL<kTotDDL){ 
132         Int_t index=kSides*(kModPerDDL*iDDL+iCarlos)+s->GetChannel(); 
133         histo[index]->Fill(s->GetCoord2(),s->GetCoord1(),s->GetSignal());
134       }
135     }
136     
137     for(Int_t iddl=0; iddl<kTotDDL;iddl++){
138       for(Int_t imod=0; imod<kModPerDDL;imod++){
139         for(Int_t isid=0;isid<kSides;isid++){
140           Int_t index=kSides*(kModPerDDL*iddl+imod)+isid;
141           anal[index]->AnalyzeEvent(histo[index]); 
142           anal[index]->WriteToASCII(iev,timeSt,nWrittenEv[index]);
143           nWrittenEv[index]++;
144           Int_t iMod=dmap->GetModuleNumber(iddl,imod);
145           if(iMod!=-1){
146             for(Int_t ipad=0;ipad<33;ipad++){
147               Int_t st=anal[index]->GetInjPadStatus(ipad);
148               hanst->Fill(st);
149               if(anal[index]->GetInjPadStatus(ipad)>=statuscut){
150                 if(isid==0) hpad7l->Fill(ipad);
151                 if(isid==1) hpad7r->Fill(ipad);
152               }
153             }
154             if(anal[index]->GetInjPadStatus(jpad)>=statuscut){
155               Float_t vel=anal[index]->GetDriftSpeed(jpad);
156               if(isid==0) hvdriftl[iMod-240]->Fill(vel);
157               if(isid==1) hvdriftr[iMod-240]->Fill(vel);
158             }
159           }
160           if(iddl==nDDL){
161             Int_t index2=kSides*imod+isid;
162             c0->cd(index2+1);
163             histo[index]->SetMaximum(100.);
164             histo[index]->DrawCopy("colz");
165             sprintf(text,"DDL %d channel %d Side %d",nDDL,imod,isid);
166             t0->DrawLatex(0.15,0.92,text);
167             c0->Update();
168             c1->cd(index2+1);
169             gvel[index]=anal[index]->GetDriftSpeedGraph();
170             gvel[index]->SetMarkerStyle(20);
171             gvel[index]->SetTitle("");
172             sprintf(gname,"gvel%dev%d",index,iev);
173             gvel[index]->SetName(gname);
174             //      gvel[index]->SetMinimum(0);
175             //gvel[index]->SetMaximum(200);
176
177             gvel[index]->GetXaxis()->SetLimits(0,256);
178             gvel[index]->GetXaxis()->SetTitle("Anode");
179             gvel[index]->GetYaxis()->SetTitle("Drift vel.");
180             gvel[index]->GetXaxis()->SetTitleSize(0.07);
181             gvel[index]->GetYaxis()->SetTitleSize(0.07);
182             gvel[index]->GetXaxis()->SetTitleOffset(0.6);
183             gvel[index]->GetYaxis()->SetTitleOffset(0.6);
184             if(gvel[index]->GetN()>0) gvel[index]->Draw("AP");
185             Double_t *param=anal[index]->GetDriftSpeedFitParam();
186             funz->SetParameters(param[0],param[1],param[2],param[3]);
187             funz->SetLineColor(2);
188             funz->DrawCopy("LSAME");
189             t0->DrawLatex(0.15,0.92,text);
190             c1->Update();
191           }
192         }
193       }
194     }
195     delete s;
196     iev++;
197     readEv++;
198     printf(" --- OK\n");
199   }while(rd->NextEvent()&&iev<=lastEv);
200   printf("Total number of events = %d\n",readEv);
201   Float_t nfac=1./(Float_t)readEv/33./520.;
202   hanst->Scale(nfac);
203   nfac=1./(Float_t)readEv;
204   hpad7l->Scale(nfac);
205   hpad7r->Scale(nfac);
206
207   Int_t ipt0=0, ipt1=0;
208   Float_t Edrift=(1800-45)/291/0.012;  
209   TFile *outfil=new TFile("DriftSpeedHistos.root","recreate");
210   for(Int_t iMod=0; iMod<260; iMod++){
211     outfil->cd();
212     hvdriftl[iMod]->Write();    
213     hvdriftr[iMod]->Write();
214     Float_t modid=iMod+240;
215     if(hvdriftl[iMod]->GetEntries()>0){
216       Float_t avevell=hvdriftl[iMod]->GetMean();
217       Float_t rmsvell=hvdriftl[iMod]->GetRMS();
218       if(avevell > 5.5 && avevell < 8.5){
219         gvvsmod0->SetPoint(ipt0,modid,avevell);
220         gvvsmod0->SetPointError(ipt0,0,rmsvell);
221         Float_t mob=avevell*1.E5/Edrift;  
222         Float_t temper=293.15*TMath::Power((mob/1350.),-1/2.4); 
223         gtvsmod0->SetPoint(ipt0,modid,temper);
224         ++ipt0;
225       }
226     }
227     if(hvdriftr[iMod]->GetEntries()>0){
228       Float_t avevelr=hvdriftr[iMod]->GetMean();
229       Float_t rmsvelr=hvdriftr[iMod]->GetRMS();
230       if(avevelr > 5.5 && avevelr < 8.5){
231         gvvsmod1->SetPoint(ipt1,modid,avevelr);
232         gvvsmod1->SetPointError(ipt1,0,rmsvelr);
233         Float_t mob=avevelr*1.E5/Edrift;
234         Float_t temper=293.15*TMath::Power((mob/1350.),-1./2.4); 
235         gtvsmod1->SetPoint(ipt1,modid,temper);
236         ++ipt1;
237       }
238     }
239   }
240   gvvsmod0->SetName("gvvsmod0");
241   gvvsmod1->SetName("gvvsmod1");
242   gtvsmod0->SetName("gtvsmod0");
243   gtvsmod1->SetName("gtvsmod1");
244   outfil->cd();
245   gvvsmod0->Write();
246   gvvsmod1->Write();
247   gtvsmod0->Write();
248   gtvsmod1->Write();
249   outfil->Close();
250
251   TCanvas* c8=new TCanvas("c8","Drift Speed vs. mod");
252   gvvsmod0->SetTitle("");
253   gvvsmod1->SetTitle("");
254   gvvsmod0->SetMarkerStyle(20);
255   gvvsmod1->SetMarkerStyle(21);
256   gvvsmod1->SetMarkerColor(2);
257   gvvsmod0->Draw("AP");
258   gvvsmod0->SetMinimum(gvmin);
259   gvvsmod0->SetMaximum(gvmax);
260   gvvsmod0->GetXaxis()->SetTitle("Module Number");
261   Char_t title[50];
262   sprintf(title,"Vdrift at injector pad %d",jpad);
263   gvvsmod0->GetYaxis()->SetTitle(title);  
264   gvvsmod1->Draw("PSAME");
265   TLatex* tleft=new TLatex(0.7,0.82,"Side 0");
266   tleft->SetNDC();
267   tleft->SetTextColor(1);
268   tleft->Draw();
269   TLatex* tright=new TLatex(0.7,0.75,"Side 1");
270   tright->SetNDC();
271   tright->SetTextColor(2);
272   tright->Draw();
273
274   TLine *lin=new TLine(323,gvmin,323,gvmax);
275   lin->SetLineColor(4);
276   lin->Draw();
277   c8->Update();
278   c8->SaveAs("VdriftVsMod.gif");
279
280   TCanvas* c8t=new TCanvas("c8t","Temeprature vs. mod");
281   gtvsmod0->SetTitle("");
282   gtvsmod1->SetTitle("");
283   gtvsmod0->SetMarkerStyle(20);
284   gtvsmod1->SetMarkerStyle(21);
285   gtvsmod1->SetMarkerColor(2);
286   gtvsmod0->Draw("AP");
287   gtvsmod0->SetMinimum(gtmin);
288   gtvsmod0->SetMaximum(gtmax);
289   gtvsmod0->GetXaxis()->SetTitle("Module Number");
290   sprintf(title,"Estimated Temperature (K)");
291   gtvsmod0->GetYaxis()->SetTitle(title);  
292   gtvsmod1->Draw("PSAME");
293   tleft->Draw();
294   tright->Draw();
295   TLine *lint=new TLine(323,gtmin,323,gtmax);
296   lint->SetLineColor(4);
297   lint->Draw();
298   c8t->Update();
299   c8t->SaveAs("TempVsMod.gif");
300
301   TCanvas* c9=new TCanvas("c9","Injector status");
302   hanst->SetStats(0);
303   hanst->Draw();
304   hanst->GetXaxis()->SetTitle("Injector pad status");
305   hanst->GetXaxis()->CenterTitle();
306   c9->SaveAs("InjStatus.gif");
307
308 //   TCanvas* c10=new TCanvas("c10","Pad status 7",1200,600);
309 //   hpad7l->SetStats(0);
310 //   hpad7r->SetStats(0);
311 //   c10->Divide(2,1);
312 //   c10->cd(1);
313 //   hpad7l->Draw(); 
314 //   hpad7l->GetXaxis()->SetTitle("Side Left -- Pad number");
315 //   hpad7l->GetXaxis()->CenterTitle();
316 //   hpad7l->GetYaxis()->SetTitle("Number of status 7");
317 //   c10->cd(2);
318 //   hpad7r->Draw();
319 //   hpad7r->GetXaxis()->SetTitle("Side Right -- Pad number");
320 //   hpad7r->GetXaxis()->CenterTitle();
321 //   hpad7r->GetYaxis()->SetTitle("Number of status 7");
322 //   printf("Side 0, maximum pad=%d\n",hpad7l->GetMaximumBin());
323 //   printf("Side 1, maximum pad=%d\n",hpad7r->GetMaximumBin());
324   
325
326 }
327
328 void AnalyzeSDDInjectorsAllMod(Int_t nrun, Int_t n2, Char_t* dir="LHC08d_SDD", 
329                                Int_t adcfreq=20, 
330                                Int_t nDDL=0, 
331                                Int_t firstEv=15, 
332                                Int_t lastEv=15,
333                                Int_t jpad=20, 
334                                Int_t statuscut=7, 
335                                Int_t dataformat=1){
336
337   TGrid::Connect("alien:",0,0,"t");
338   Char_t filnam[200];
339   sprintf(filnam,"alien:///alice/data/2008/%s/%09d/raw/08%09d%03d.10.root",dir,nrun,nrun,n2);
340   printf("Open file %s\n",filnam);
341   AnalyzeSDDInjectorsAllMod(filnam,adcfreq,nDDL,firstEv,lastEv,jpad,statuscut,dataformat);
342 }
343
344
345