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