Changes related to the extraction of the V0 finder into a separate class (A. Dainese...
[u/mrichter/AliRoot.git] / ITS / AnalyzeSDDNoiseAllMod.C
1 #if !defined(__CINT__) || defined(__MAKECINT__)
2 #include <TH2F.h>
3 #include <TCanvas.h>
4 #include <TStopwatch.h>
5 #include <TStyle.h>
6 #include <TGrid.h>
7 #include <TSystem.h>
8 #include <TLatex.h>
9 #include <TFile.h>
10 #include "AliRawReader.h"
11 #include "AliRawReaderDate.h"
12 #include "AliRawReaderRoot.h"
13 #include "AliITSOnlineSDDBase.h"
14 #include "AliITSOnlineSDDCMN.h"
15 #include "AliITSRawStreamSDD.h"
16 #include "AliITSRawStreamSDDCompressed.h"
17 #include "TPaveStats.h"
18 #endif
19
20 // Macro for the analysis of PEDESTAL runs (equivalent to ITSSDDBASda.cxx)
21 // Two functions named AnalyzeSDDNoiseAllModules: 
22 // The first is for analyzing a local raw data file and takes as agrument the file name.
23 // The second is for running on ALIEN
24 // All DDLs are analyzed, the argument nDDL selects the DDL to be plotted
25 // Origin: F. Prino (prino@to.infn.it)
26
27 void AnalyzeSDDNoiseAllMod(Char_t *datafil, 
28                            Int_t adcfreq=20, 
29                            Int_t nDDL=0, 
30                            Int_t firstEv=10, 
31                            Int_t lastEv=12, 
32                            Int_t dataformat=1){
33
34
35   const Int_t kTotDDL=24;
36   const Int_t kModPerDDL=12;
37   const Int_t kSides=2;
38
39   AliITSOnlineSDDBase **base=new AliITSOnlineSDDBase*[kTotDDL*kModPerDDL*kSides];
40   TH2F **histo=new TH2F*[kTotDDL*kModPerDDL*kSides];
41
42   Char_t hisnam[20];
43   for(Int_t iddl=0; iddl<kTotDDL;iddl++){
44     for(Int_t imod=0; imod<kModPerDDL;imod++){
45       for(Int_t isid=0;isid<kSides;isid++){
46         Int_t index=kSides*(kModPerDDL*iddl+imod)+isid;
47         base[index]=new AliITSOnlineSDDBase(iddl,imod,isid);    
48         if(adcfreq==40) base[index]->SetLastGoodTB(254);
49         else base[index]->SetLastGoodTB(126);
50         sprintf(hisnam,"h%02dc%02ds%d",iddl,imod,isid);
51         histo[index]=new TH2F(hisnam,"",256,-0.5,255.5,256,-0.5,255.5);
52       }
53     }
54   }
55
56
57   TCanvas* c0 = new TCanvas("c0","Ev Display",900,900);
58   gStyle->SetPalette(1);
59   Char_t text[50];
60
61   Int_t iev=firstEv;
62   AliRawReader *rd; 
63   if(strstr(datafil,".root")!=0){
64     rd=new AliRawReaderRoot(datafil,iev);
65   }else{
66     rd=new AliRawReaderDate(datafil,iev);
67   }
68   TLatex *t0=new TLatex();
69   t0->SetNDC();
70   t0->SetTextSize(0.06);
71   t0->SetTextColor(4);
72
73   do{
74     c0->Clear();
75     c0->Divide(4,6,0.001,0.001);
76     printf("Event # %d ",iev);
77     rd->Reset();
78     for(Int_t iddl=0; iddl<kTotDDL;iddl++){
79       for(Int_t imod=0; imod<kModPerDDL;imod++){
80         for(Int_t isid=0;isid<kSides;isid++){
81           Int_t index=kSides*(kModPerDDL*iddl+imod)+isid;
82           histo[index]->Reset();
83         }
84       }
85     }
86     AliITSRawStream* s;
87     if(dataformat==0){
88       s=new AliITSRawStreamSDD(rd);
89     }else{
90       s=new AliITSRawStreamSDDCompressed(rd);
91       if(dataformat==1) s->SetADCEncoded(kTRUE);
92     }
93     while(s->Next()){
94       Int_t iDDL=rd->GetDDLID();
95       Int_t iCarlos=s->GetCarlosId();
96       if(s->IsCompletedModule()) continue;
97       if(s->IsCompletedDDL()) continue;
98       if(iDDL>=0 && iDDL<kTotDDL){ 
99         Int_t index=kSides*(kModPerDDL*iDDL+iCarlos)+s->GetChannel(); 
100         histo[index]->Fill(s->GetCoord2(),s->GetCoord1(),s->GetSignal());
101       }
102     }
103     delete s;
104     iev++;
105     for(Int_t iddl=0; iddl<kTotDDL;iddl++){
106       for(Int_t imod=0; imod<kModPerDDL;imod++){
107         for(Int_t isid=0;isid<kSides;isid++){
108           Int_t index=kSides*(kModPerDDL*iddl+imod)+isid;
109           base[index]->AddEvent(histo[index]);
110           if(iddl==nDDL){
111             Int_t index2=kSides*imod+isid;
112             c0->cd(index2+1);
113             histo[index]->DrawCopy("colz");
114             sprintf(text,"DDL %d channel %d Side %d",nDDL,imod,isid);
115             t0->DrawLatex(0.15,0.92,text);
116             c0->Update();
117           }
118         }
119       }
120     }
121     printf(" --- OK\n");
122   }while(rd->NextEvent()&&iev<=lastEv);
123
124   for(Int_t iddl=0; iddl<kTotDDL;iddl++){
125     for(Int_t imod=0; imod<kModPerDDL;imod++){
126       for(Int_t isid=0;isid<kSides;isid++){
127         Int_t index=kSides*(kModPerDDL*iddl+imod)+isid;
128         base[index]->ValidateAnodes();
129         base[index]->WriteToASCII(); // fondamentale!!!!!!!!!
130         delete base[index];
131       }
132     }
133   }
134   delete rd;
135   delete [] base;
136   
137   printf("Start second analysis for Common Mode correction\n");
138   AliITSOnlineSDDCMN **corr=new AliITSOnlineSDDCMN*[kTotDDL*kModPerDDL*kSides];
139   Bool_t isFilled[kTotDDL*kModPerDDL*kSides];
140
141   for(Int_t iddl=0; iddl<kTotDDL;iddl++){
142     for(Int_t imod=0; imod<kModPerDDL;imod++){
143       for(Int_t isid=0;isid<kSides;isid++){
144         Int_t index=kSides*(kModPerDDL*iddl+imod)+isid;
145         corr[index]=new AliITSOnlineSDDCMN(iddl,imod,isid);
146         if(adcfreq==40) corr[index]->SetLastGoodTB(254);
147         else corr[index]->SetLastGoodTB(126);
148         isFilled[index]=0;
149       }
150     }
151   }
152
153   iev=firstEv;
154   AliRawReader *rd2; 
155   if(strstr(datafil,".root")!=0){
156     rd2=new AliRawReaderRoot(datafil,iev);
157   }else{
158     rd2=new AliRawReaderDate(datafil,iev);
159   }
160   do{
161     c0->Clear();
162     c0->Divide(4,6,0.001,0.001);
163     printf("Event # %d ",iev);
164     rd2->Reset();
165     for(Int_t iddl=0; iddl<kTotDDL;iddl++){
166       for(Int_t imod=0; imod<kModPerDDL;imod++){
167         for(Int_t isid=0;isid<kSides;isid++){
168           Int_t index=kSides*(kModPerDDL*iddl+imod)+isid;
169           histo[index]->Reset();
170         }
171       }
172     }
173     
174     AliITSRawStream* s;
175     if(dataformat==0){
176       s=new AliITSRawStreamSDD(rd2);
177     }else{
178       s=new AliITSRawStreamSDDCompressed(rd2);
179       if(dataformat==1) s->SetADCEncoded(kTRUE);
180     }
181     while(s->Next()){
182       Int_t iDDL=rd2->GetDDLID();
183       Int_t iCarlos=s->GetCarlosId();
184       if(s->IsCompletedModule()) continue;
185       if(s->IsCompletedDDL()) continue;
186       if(iDDL>=0 && iDDL<kTotDDL){ 
187         Int_t index=kSides*(kModPerDDL*iDDL+iCarlos)+s->GetChannel(); 
188         histo[index]->Fill(s->GetCoord2(),s->GetCoord1(),s->GetSignal());
189         isFilled[index]=1;
190       }
191     }
192     delete s;
193     iev++;
194     for(Int_t iddl=0; iddl<kTotDDL;iddl++){
195       for(Int_t imod=0; imod<kModPerDDL;imod++){
196         for(Int_t isid=0;isid<kSides;isid++){
197           Int_t index=kSides*(kModPerDDL*iddl+imod)+isid;
198           if(isFilled[index]) corr[index]->AddEvent(histo[index]);
199           if(iddl==nDDL){
200             Int_t index2=kSides*imod+isid;
201             c0->cd(index2+1);
202             histo[index]->DrawCopy("colz");
203             sprintf(text,"DDL %d channel %d Side %d",nDDL,imod,isid);
204             t0->DrawLatex(0.15,0.92,text);
205             c0->Update();
206           }
207         }
208       }
209     }
210     printf(" --- OK\n");
211   }while(rd2->NextEvent()&&iev<=lastEv);
212
213   TH1F *htotbas=new TH1F("htotbas","",100,0.,150.);
214   TH1F *htotbaseq=new TH1F("htotbaseq","",100,0.,150.);
215   TH1F *htotnoise=new TH1F("htotnoise","",100,0.,10.);
216   TH1F *htotnoisecorr=new TH1F("htotnoisecorr","",100,0.,10.);
217   TH1F *hstatus=new TH1F("hstatus","",2,-0.5,1.5);
218
219   TFile *outfil=new TFile("SDDbase-results.root","recreate");
220   for(Int_t iddl=0; iddl<kTotDDL;iddl++){
221     for(Int_t imod=0; imod<kModPerDDL;imod++){
222       for(Int_t isid=0;isid<kSides;isid++){
223         Int_t index=kSides*(kModPerDDL*iddl+imod)+isid;
224         if(isFilled[index]){
225           corr[index]->ValidateAnodes();
226           corr[index]->WriteToASCII();
227           corr[index]->WriteToROOT(outfil);
228           for(Int_t ian=0; ian<256;ian++){
229             Float_t basl=corr[index]->GetAnodeBaseline(ian);
230             Float_t basleq=corr[index]->GetAnodeEqualizedBaseline(ian);
231             Float_t noi=corr[index]->GetAnodeRawNoise(ian);
232             Float_t cornoi=corr[index]->GetAnodeCorrNoise(ian);
233             Int_t anstatus=corr[index]->IsAnodeGood(ian);
234             hstatus->Fill(anstatus);
235             htotbas->Fill(basl);
236             htotbaseq->Fill(basleq);
237             htotnoise->Fill(noi);
238             htotnoisecorr->Fill(cornoi);
239           }
240         }
241       }
242     }
243   }
244   outfil->Close();
245   // Draw Statistics of baselines and noise
246   TCanvas *call=new TCanvas("call","General stats",700,700);
247   call->Divide(2,2);
248   call->cd(1);
249   htotbas->Draw();
250   htotbas->GetXaxis()->SetTitle("Baselines");
251   htotbas->GetXaxis()->SetTitleSize(0.07);
252   htotbas->GetXaxis()->SetTitleOffset(0.6);
253   call->cd(2);
254   htotbaseq->Draw();
255   htotbaseq->GetXaxis()->SetTitle("Baselines after equalization");
256   htotbaseq->GetXaxis()->SetTitleSize(0.07);
257   htotbaseq->GetXaxis()->SetTitleOffset(0.6);
258   call->cd(3);
259   htotnoisecorr->SetLineColor(2);
260   htotnoisecorr->Draw();
261   call->Update();
262   TPaveStats *st1=(TPaveStats*)htotnoisecorr->GetListOfFunctions()->FindObject("stats");
263   st1->SetY1NDC(0.51);
264   st1->SetY2NDC(0.7);
265   htotnoisecorr->GetXaxis()->SetTitle("Noise");
266   htotnoisecorr->GetXaxis()->SetTitleSize(0.07);
267   htotnoisecorr->GetXaxis()->SetTitleOffset(0.6);
268   htotnoise->Draw("SAMES");
269   call->Update();
270   TPaveStats *st2=(TPaveStats*)htotnoise->GetListOfFunctions()->FindObject("stats");
271   st2->SetY1NDC(0.71);
272   st2->SetY2NDC(0.9);
273   
274   call->cd(4);
275   hstatus->Draw();
276   hstatus->GetXaxis()->SetTitle("Anode Status (0=bad 1=good)");
277   hstatus->GetXaxis()->SetTitleSize(0.07);
278   hstatus->GetXaxis()->SetTitleOffset(0.6);
279   call->Update();
280   call->SaveAs("GenStatsPedestal.gif");
281
282   // Draw baselines and noise for all modules of the selected DDL
283
284   TH1F** hbas = new TH1F*[kSides*kModPerDDL];
285   TH1F** hrawn = new TH1F*[kSides*kModPerDDL];
286   TH1F** hcorrn = new TH1F*[kSides*kModPerDDL];
287   TH1F** hdbas = new TH1F*[kSides*kModPerDDL];
288   TH1F** hdrawn = new TH1F*[kSides*kModPerDDL];
289   TH1F** hdcorrn = new TH1F*[kSides*kModPerDDL];
290   TCanvas *c1=new TCanvas("c1","DDL: Baselines vs anode",900,900);
291   c1->SetBottomMargin(0.14);
292   c1->Divide(4,6,0.001,0.001);
293   TCanvas *c2=new TCanvas("c2","DDL: Noise vs anode",900,900);
294   c2->SetBottomMargin(0.14);
295   c2->Divide(4,6,0.001,0.001);
296   TCanvas *c3=new TCanvas("c3","DDL: Baselines distr",900,900);
297   c3->SetBottomMargin(0.14);
298   c3->Divide(4,6,0.001,0.001);
299   TCanvas *c4=new TCanvas("c4","DDL: Noise Distr",900,900);
300   c4->SetBottomMargin(0.14);
301   c4->Divide(4,6,0.001,0.001);
302   TLatex *t1=new TLatex(0.15,0.2,"Raw Noise");
303   t1->SetNDC();
304   t1->SetTextSize(0.05);
305   TLatex *t2=new TLatex(0.4,0.2,"Corrected Noise");
306   t2->SetNDC();
307   t2->SetTextSize(0.05);
308   t2->SetTextColor(2);
309   TLatex *t3=new TLatex();
310   t3->SetNDC();
311   t3->SetTextSize(0.06);
312   t3->SetTextColor(4);
313
314   for(Int_t imod=0; imod<kModPerDDL;imod++){
315     for(Int_t isid=0;isid<kSides;isid++){
316       Int_t index1=kSides*(kModPerDDL*nDDL+imod)+isid;
317       Int_t index2=kSides*imod+isid;
318       sprintf(text,"DDL %d channel %d Side %d",nDDL,imod,isid);
319       hbas[index2]=corr[index1]->GetBaselineAnodeHisto();
320       hrawn[index2]=corr[index1]->GetRawNoiseAnodeHisto();
321       hcorrn[index2]=corr[index1]->GetCorrNoiseAnodeHisto();
322       hdbas[index2]=corr[index1]->GetBaselineHisto();
323       hdrawn[index2]=corr[index1]->GetRawNoiseHisto();
324       hdcorrn[index2]=corr[index1]->GetCorrNoiseHisto();
325       c1->cd(index2+1);
326       hbas[index2]->Draw();
327       hbas[index2]->SetMinimum(0);
328       hbas[index2]->SetMaximum(75);
329       hbas[index2]->GetXaxis()->SetTitle("Anode");
330       hbas[index2]->GetYaxis()->SetTitle("Baseline");
331       hbas[index2]->GetXaxis()->SetTitleSize(0.07);
332       hbas[index2]->GetYaxis()->SetTitleSize(0.07);
333       hbas[index2]->GetXaxis()->SetTitleOffset(0.6);
334       hbas[index2]->GetYaxis()->SetTitleOffset(0.7);
335       t3->DrawLatex(0.15,0.92,text);
336       c1->Update();
337
338
339       c2->cd(index2+1); 
340       hrawn[index2]->SetMinimum(1.);
341       hrawn[index2]->SetMaximum(6.);
342       hrawn[index2]->Draw();
343       hrawn[index2]->GetXaxis()->SetTitle("Anode");
344       hrawn[index2]->GetYaxis()->SetTitle("Noise");
345       hrawn[index2]->GetXaxis()->SetTitleSize(0.07);
346       hrawn[index2]->GetYaxis()->SetTitleSize(0.07);
347       hrawn[index2]->GetXaxis()->SetTitleOffset(0.6);
348       hrawn[index2]->GetYaxis()->SetTitleOffset(0.7);
349       gStyle->SetOptStat(0);
350       hrawn[index2]->SetStats(0);
351       hcorrn[index2]->SetLineColor(2);
352       hcorrn[index2]->Draw("SAME");
353       t1->Draw();
354       t2->Draw();
355       t3->DrawLatex(0.15,0.92,text);
356       c2->Update();
357
358       c3->cd(index2+1);
359       hdbas[index2]->Draw();
360       hdbas[index2]->GetXaxis()->SetTitle("Baseline");
361       hdbas[index2]->GetXaxis()->SetTitleSize(0.07);
362       hdbas[index2]->GetXaxis()->SetTitleOffset(0.6);
363       t3->DrawLatex(0.15,0.92,text);
364       c3->Update();
365
366       c4->cd(index2+1); 
367       hdrawn[index2]->Draw();
368       hdrawn[index2]->GetXaxis()->SetTitle("Noise");
369       hdrawn[index2]->GetXaxis()->SetTitleSize(0.07);
370       hdrawn[index2]->GetXaxis()->SetTitleOffset(0.6);
371       hdcorrn[index2]->SetLineColor(2);
372       hdcorrn[index2]->Draw("SAME");
373       t1->Draw();
374       t2->Draw();
375       t3->DrawLatex(0.15,0.92,text);
376       c4->Update();
377     }
378   }
379
380   c1->SaveAs("Baselines.gif");
381   c2->SaveAs("Noise.gif");
382   c3->SaveAs("BaselinesDist.gif");
383   c4->SaveAs("NoiseDist.gif");
384   
385   Char_t delfil[100];
386   for(Int_t iddl=0; iddl<kTotDDL;iddl++){
387     for(Int_t imod=0; imod<kModPerDDL;imod++){
388       for(Int_t isid=0;isid<kSides;isid++){
389         sprintf(delfil,"rm SDDbase_step1_ddl%02dc%02d_sid%d.data",iddl,imod,isid);
390         gSystem->Exec(delfil);
391       }
392     }
393   }
394
395 }
396
397 void AnalyzeSDDNoiseAllMod(Int_t nrun, Int_t n2, Char_t* dir="LHC08d_SDD",
398                            Int_t adcfreq=20, 
399                            Int_t nDDL=0, 
400                            Int_t firstEv=15, 
401                            Int_t lastEv=18, 
402                            Int_t dataformat=1){
403
404   TGrid::Connect("alien:",0,0,"t");
405   Char_t filnam[200];
406   sprintf(filnam,"alien:///alice/data/2008/%s/%09d/raw/08%09d%03d.10.root",dir,nrun,nrun,n2);
407   printf("Open file %s\n",filnam);
408   AnalyzeSDDNoiseAllMod(filnam,adcfreq,nDDL,firstEv,lastEv,dataformat);
409 }