Changing once more (hopefully we get it correct this time...) the logic to trig the...
[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=18, 
31                            Int_t lastEv=20){ 
32
33
34   const Int_t kTotDDL=24;
35   const Int_t kModPerDDL=12;
36   const Int_t kSides=2;
37   Bool_t writtenoutput=kFALSE;
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
87     UChar_t cdhAttr=AliITSRawStreamSDD::ReadBlockAttributes(rd);
88     UInt_t amSamplFreq=AliITSRawStreamSDD::ReadAMSamplFreqFromCDH(cdhAttr);
89     AliITSRawStream* s=AliITSRawStreamSDD::CreateRawStreamSDD(rd,cdhAttr);
90     if(!writtenoutput){
91       printf("Use %s raw stream, sampling frequency %d MHz\n",s->ClassName(),amSamplFreq);
92       writtenoutput=kTRUE;
93     }
94
95     while(s->Next()){
96       Int_t iDDL=rd->GetDDLID();
97       Int_t iCarlos=s->GetCarlosId();
98       if(s->IsCompletedModule()) continue;
99       if(s->IsCompletedDDL()) continue;
100       if(iDDL>=0 && iDDL<kTotDDL){ 
101         Int_t index=kSides*(kModPerDDL*iDDL+iCarlos)+s->GetChannel(); 
102         histo[index]->Fill(s->GetCoord2(),s->GetCoord1(),s->GetSignal());
103       }
104     }
105     delete s;
106     iev++;
107     for(Int_t iddl=0; iddl<kTotDDL;iddl++){
108       for(Int_t imod=0; imod<kModPerDDL;imod++){
109         for(Int_t isid=0;isid<kSides;isid++){
110           Int_t index=kSides*(kModPerDDL*iddl+imod)+isid;
111           base[index]->AddEvent(histo[index]);
112           if(iddl==nDDL){
113             Int_t index2=kSides*imod+isid;
114             c0->cd(index2+1);
115             histo[index]->DrawCopy("colz");
116             sprintf(text,"DDL %d channel %d Side %d",nDDL,imod,isid);
117             t0->DrawLatex(0.15,0.92,text);
118             c0->Update();
119           }
120         }
121       }
122     }
123     printf(" --- OK\n");
124   }while(rd->NextEvent()&&iev<=lastEv);
125
126   for(Int_t iddl=0; iddl<kTotDDL;iddl++){
127     for(Int_t imod=0; imod<kModPerDDL;imod++){
128       for(Int_t isid=0;isid<kSides;isid++){
129         Int_t index=kSides*(kModPerDDL*iddl+imod)+isid;
130         base[index]->ValidateAnodes();
131         base[index]->WriteToASCII(); // fondamentale!!!!!!!!!
132         delete base[index];
133       }
134     }
135   }
136   delete rd;
137   delete [] base;
138   
139   printf("Start second analysis for Common Mode correction\n");
140   AliITSOnlineSDDCMN **corr=new AliITSOnlineSDDCMN*[kTotDDL*kModPerDDL*kSides];
141   Bool_t isFilled[kTotDDL*kModPerDDL*kSides];
142
143   for(Int_t iddl=0; iddl<kTotDDL;iddl++){
144     for(Int_t imod=0; imod<kModPerDDL;imod++){
145       for(Int_t isid=0;isid<kSides;isid++){
146         Int_t index=kSides*(kModPerDDL*iddl+imod)+isid;
147         corr[index]=new AliITSOnlineSDDCMN(iddl,imod,isid);
148         if(adcfreq==40) corr[index]->SetLastGoodTB(254);
149         else corr[index]->SetLastGoodTB(126);
150         isFilled[index]=0;
151       }
152     }
153   }
154
155   iev=firstEv;
156   AliRawReader *rd2; 
157   if(strstr(datafil,".root")!=0){
158     rd2=new AliRawReaderRoot(datafil,iev);
159   }else{
160     rd2=new AliRawReaderDate(datafil,iev);
161   }
162   do{
163     c0->Clear();
164     c0->Divide(4,6,0.001,0.001);
165     printf("Event # %d ",iev);
166     rd2->Reset();
167     for(Int_t iddl=0; iddl<kTotDDL;iddl++){
168       for(Int_t imod=0; imod<kModPerDDL;imod++){
169         for(Int_t isid=0;isid<kSides;isid++){
170           Int_t index=kSides*(kModPerDDL*iddl+imod)+isid;
171           histo[index]->Reset();
172         }
173       }
174     }
175     
176     UChar_t cdhAttr=AliITSRawStreamSDD::ReadBlockAttributes(rd2);
177     UInt_t amSamplFreq=AliITSRawStreamSDD::ReadAMSamplFreqFromCDH(cdhAttr);
178     AliITSRawStream* s=AliITSRawStreamSDD::CreateRawStreamSDD(rd2,cdhAttr);
179     if(!writtenoutput){
180       printf("Use %s raw stream, sampling frequency %d MHz\n",s->ClassName(),amSamplFreq);
181       writtenoutput=kTRUE;
182     }
183     while(s->Next()){
184       Int_t iDDL=rd2->GetDDLID();
185       Int_t iCarlos=s->GetCarlosId();
186       if(s->IsCompletedModule()) continue;
187       if(s->IsCompletedDDL()) continue;
188       if(iDDL>=0 && iDDL<kTotDDL){ 
189         Int_t index=kSides*(kModPerDDL*iDDL+iCarlos)+s->GetChannel(); 
190         histo[index]->Fill(s->GetCoord2(),s->GetCoord1(),s->GetSignal());
191         isFilled[index]=1;
192       }
193     }
194     delete s;
195     iev++;
196     for(Int_t iddl=0; iddl<kTotDDL;iddl++){
197       for(Int_t imod=0; imod<kModPerDDL;imod++){
198         for(Int_t isid=0;isid<kSides;isid++){
199           Int_t index=kSides*(kModPerDDL*iddl+imod)+isid;
200           if(isFilled[index]) corr[index]->AddEvent(histo[index]);
201           if(iddl==nDDL){
202             Int_t index2=kSides*imod+isid;
203             c0->cd(index2+1);
204             histo[index]->DrawCopy("colz");
205             sprintf(text,"DDL %d channel %d Side %d",nDDL,imod,isid);
206             t0->DrawLatex(0.15,0.92,text);
207             c0->Update();
208           }
209         }
210       }
211     }
212     printf(" --- OK\n");
213   }while(rd2->NextEvent()&&iev<=lastEv);
214
215   TH1F *htotbas=new TH1F("htotbas","",100,0.,150.);
216   TH1F *htotbaseq=new TH1F("htotbaseq","",100,0.,150.);
217   TH1F *htotnoise=new TH1F("htotnoise","",100,0.,10.);
218   TH1F *htotnoisecorr=new TH1F("htotnoisecorr","",100,0.,10.);
219   TH1F *hstatus=new TH1F("hstatus","",2,-0.5,1.5);
220
221   TFile *outfil=new TFile("SDDbase-results.root","recreate");
222   for(Int_t iddl=0; iddl<kTotDDL;iddl++){
223     for(Int_t imod=0; imod<kModPerDDL;imod++){
224       for(Int_t isid=0;isid<kSides;isid++){
225         Int_t index=kSides*(kModPerDDL*iddl+imod)+isid;
226         if(isFilled[index]){
227           corr[index]->ValidateAnodes();
228           corr[index]->WriteToASCII();
229           corr[index]->WriteToROOT(outfil);
230           for(Int_t ian=0; ian<256;ian++){
231             Float_t basl=corr[index]->GetAnodeBaseline(ian);
232             Float_t basleq=corr[index]->GetAnodeEqualizedBaseline(ian);
233             Float_t noi=corr[index]->GetAnodeRawNoise(ian);
234             Float_t cornoi=corr[index]->GetAnodeCorrNoise(ian);
235             Int_t anstatus=corr[index]->IsAnodeGood(ian);
236             hstatus->Fill(anstatus);
237             htotbas->Fill(basl);
238             htotbaseq->Fill(basleq);
239             htotnoise->Fill(noi);
240             htotnoisecorr->Fill(cornoi);
241           }
242         }
243       }
244     }
245   }
246   outfil->Close();
247   // Draw Statistics of baselines and noise
248   TCanvas *call=new TCanvas("call","General stats",700,700);
249   call->Divide(2,2);
250   call->cd(1);
251   htotbas->Draw();
252   htotbas->GetXaxis()->SetTitle("Baselines");
253   htotbas->GetXaxis()->SetTitleSize(0.07);
254   htotbas->GetXaxis()->SetTitleOffset(0.6);
255   call->cd(2);
256   htotbaseq->Draw();
257   htotbaseq->GetXaxis()->SetTitle("Baselines after equalization");
258   htotbaseq->GetXaxis()->SetTitleSize(0.07);
259   htotbaseq->GetXaxis()->SetTitleOffset(0.6);
260   call->cd(3);
261   htotnoisecorr->SetLineColor(2);
262   htotnoisecorr->Draw();
263   call->Update();
264   TPaveStats *st1=(TPaveStats*)htotnoisecorr->GetListOfFunctions()->FindObject("stats");
265   st1->SetY1NDC(0.51);
266   st1->SetY2NDC(0.7);
267   htotnoisecorr->GetXaxis()->SetTitle("Noise");
268   htotnoisecorr->GetXaxis()->SetTitleSize(0.07);
269   htotnoisecorr->GetXaxis()->SetTitleOffset(0.6);
270   htotnoise->Draw("SAMES");
271   call->Update();
272   TPaveStats *st2=(TPaveStats*)htotnoise->GetListOfFunctions()->FindObject("stats");
273   st2->SetY1NDC(0.71);
274   st2->SetY2NDC(0.9);
275   
276   call->cd(4);
277   hstatus->Draw();
278   hstatus->GetXaxis()->SetTitle("Anode Status (0=bad 1=good)");
279   hstatus->GetXaxis()->SetTitleSize(0.07);
280   hstatus->GetXaxis()->SetTitleOffset(0.6);
281   call->Update();
282   call->SaveAs("GenStatsPedestal.gif");
283
284   // Draw baselines and noise for all modules of the selected DDL
285
286   TH1F** hbas = new TH1F*[kSides*kModPerDDL];
287   TH1F** hrawn = new TH1F*[kSides*kModPerDDL];
288   TH1F** hcorrn = new TH1F*[kSides*kModPerDDL];
289   TH1F** hdbas = new TH1F*[kSides*kModPerDDL];
290   TH1F** hdrawn = new TH1F*[kSides*kModPerDDL];
291   TH1F** hdcorrn = new TH1F*[kSides*kModPerDDL];
292   TCanvas *c1=new TCanvas("c1","DDL: Baselines vs anode",900,900);
293   c1->SetBottomMargin(0.14);
294   c1->Divide(4,6,0.001,0.001);
295   TCanvas *c2=new TCanvas("c2","DDL: Noise vs anode",900,900);
296   c2->SetBottomMargin(0.14);
297   c2->Divide(4,6,0.001,0.001);
298   TCanvas *c3=new TCanvas("c3","DDL: Baselines distr",900,900);
299   c3->SetBottomMargin(0.14);
300   c3->Divide(4,6,0.001,0.001);
301   TCanvas *c4=new TCanvas("c4","DDL: Noise Distr",900,900);
302   c4->SetBottomMargin(0.14);
303   c4->Divide(4,6,0.001,0.001);
304   TLatex *t1=new TLatex(0.15,0.2,"Raw Noise");
305   t1->SetNDC();
306   t1->SetTextSize(0.05);
307   TLatex *t2=new TLatex(0.4,0.2,"Corrected Noise");
308   t2->SetNDC();
309   t2->SetTextSize(0.05);
310   t2->SetTextColor(2);
311   TLatex *t3=new TLatex();
312   t3->SetNDC();
313   t3->SetTextSize(0.06);
314   t3->SetTextColor(4);
315
316   for(Int_t imod=0; imod<kModPerDDL;imod++){
317     for(Int_t isid=0;isid<kSides;isid++){
318       Int_t index1=kSides*(kModPerDDL*nDDL+imod)+isid;
319       Int_t index2=kSides*imod+isid;
320       sprintf(text,"DDL %d channel %d Side %d",nDDL,imod,isid);
321       hbas[index2]=corr[index1]->GetBaselineAnodeHisto();
322       hrawn[index2]=corr[index1]->GetRawNoiseAnodeHisto();
323       hcorrn[index2]=corr[index1]->GetCorrNoiseAnodeHisto();
324       hdbas[index2]=corr[index1]->GetBaselineHisto();
325       hdrawn[index2]=corr[index1]->GetRawNoiseHisto();
326       hdcorrn[index2]=corr[index1]->GetCorrNoiseHisto();
327       c1->cd(index2+1);
328       hbas[index2]->Draw();
329       hbas[index2]->SetMinimum(0);
330       hbas[index2]->SetMaximum(75);
331       hbas[index2]->GetXaxis()->SetTitle("Anode");
332       hbas[index2]->GetYaxis()->SetTitle("Baseline");
333       hbas[index2]->GetXaxis()->SetTitleSize(0.07);
334       hbas[index2]->GetYaxis()->SetTitleSize(0.07);
335       hbas[index2]->GetXaxis()->SetTitleOffset(0.6);
336       hbas[index2]->GetYaxis()->SetTitleOffset(0.7);
337       t3->DrawLatex(0.15,0.92,text);
338       c1->Update();
339
340
341       c2->cd(index2+1); 
342       hrawn[index2]->SetMinimum(1.);
343       hrawn[index2]->SetMaximum(6.);
344       hrawn[index2]->Draw();
345       hrawn[index2]->GetXaxis()->SetTitle("Anode");
346       hrawn[index2]->GetYaxis()->SetTitle("Noise");
347       hrawn[index2]->GetXaxis()->SetTitleSize(0.07);
348       hrawn[index2]->GetYaxis()->SetTitleSize(0.07);
349       hrawn[index2]->GetXaxis()->SetTitleOffset(0.6);
350       hrawn[index2]->GetYaxis()->SetTitleOffset(0.7);
351       gStyle->SetOptStat(0);
352       hrawn[index2]->SetStats(0);
353       hcorrn[index2]->SetLineColor(2);
354       hcorrn[index2]->Draw("SAME");
355       t1->Draw();
356       t2->Draw();
357       t3->DrawLatex(0.15,0.92,text);
358       c2->Update();
359
360       c3->cd(index2+1);
361       hdbas[index2]->Draw();
362       hdbas[index2]->GetXaxis()->SetTitle("Baseline");
363       hdbas[index2]->GetXaxis()->SetTitleSize(0.07);
364       hdbas[index2]->GetXaxis()->SetTitleOffset(0.6);
365       t3->DrawLatex(0.15,0.92,text);
366       c3->Update();
367
368       c4->cd(index2+1); 
369       hdrawn[index2]->Draw();
370       hdrawn[index2]->GetXaxis()->SetTitle("Noise");
371       hdrawn[index2]->GetXaxis()->SetTitleSize(0.07);
372       hdrawn[index2]->GetXaxis()->SetTitleOffset(0.6);
373       hdcorrn[index2]->SetLineColor(2);
374       hdcorrn[index2]->Draw("SAME");
375       t1->Draw();
376       t2->Draw();
377       t3->DrawLatex(0.15,0.92,text);
378       c4->Update();
379     }
380   }
381
382   c1->SaveAs("Baselines.gif");
383   c2->SaveAs("Noise.gif");
384   c3->SaveAs("BaselinesDist.gif");
385   c4->SaveAs("NoiseDist.gif");
386   
387   Char_t delfil[100];
388   for(Int_t iddl=0; iddl<kTotDDL;iddl++){
389     for(Int_t imod=0; imod<kModPerDDL;imod++){
390       for(Int_t isid=0;isid<kSides;isid++){
391         sprintf(delfil,"rm SDDbase_step1_ddl%02dc%02d_sid%d.data",iddl,imod,isid);
392         gSystem->Exec(delfil);
393       }
394     }
395   }
396
397 }
398
399 void AnalyzeSDDNoiseAllMod(Int_t nrun, Int_t n2, Int_t year=2009, Char_t* dir="LHC09b_SDD",
400                            Int_t adcfreq=20, 
401                            Int_t nDDL=0, 
402                            Int_t firstEv=18, 
403                            Int_t lastEv=20){
404
405   TGrid::Connect("alien:",0,0,"t");
406   Char_t filnam[200];
407   sprintf(filnam,"alien:///alice/data/%d/%s/%09d/raw/%02d%09d%03d.10.root",year,dir,nrun,year-2000,nrun,n2);
408   printf("Open file %s\n",filnam);
409   AnalyzeSDDNoiseAllMod(filnam,adcfreq,nDDL,firstEv,lastEv);
410 }