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 }