First set of modifications to handle in SDD DAs the information about ADC sampling...
authorprino <prino@f7af4fe6-9843-0410-8265-dc069ae4e863>
Mon, 8 Jun 2009 16:47:38 +0000 (16:47 +0000)
committerprino <prino@f7af4fe6-9843-0410-8265-dc069ae4e863>
Mon, 8 Jun 2009 16:47:38 +0000 (16:47 +0000)
ITS/AliITSOnlineSDDCMN.cxx
ITS/AliITSOnlineSDDCMN.h
ITS/AliITSOnlineSDDInjectors.h
ITS/AnalyzeSDDGainAllMod.C
ITS/AnalyzeSDDInjectorsAllMod.C
ITS/AnalyzeSDDNoiseAllMod.C

index 83450fd..c1999bf 100644 (file)
@@ -106,10 +106,8 @@ void  AliITSOnlineSDDCMN::ValidateAnodes(){
 }
 
 //______________________________________________________________________
-void AliITSOnlineSDDCMN::AddEvent(TH2F* hrawd){
+TH2F* AliITSOnlineSDDCMN::GetCleanEvent(TH2F* hrawd) const {
   // 
-  fNEvents++;
-  const Int_t kTimeBins=fLastGoodTB+1;
   TH2F* hcorrd=new TH2F("hcorrd","",hrawd->GetNbinsX(),hrawd->GetXaxis()->GetXmin(),hrawd->GetXaxis()->GetXmax(),hrawd->GetNbinsY(),hrawd->GetYaxis()->GetXmin(),hrawd->GetYaxis()->GetXmax());
   for(Int_t itb=fFirstGoodTB;itb<=fLastGoodTB;itb++){
     Float_t sumEven=0., sumOdd=0.;
@@ -133,14 +131,23 @@ void AliITSOnlineSDDCMN::AddEvent(TH2F* hrawd){
       hcorrd->SetBinContent(itb+1,ian+1,cntCorr);
     }
   }
+  return hcorrd;
+}
+//______________________________________________________________________
+void AliITSOnlineSDDCMN::AddEvent(TH2F* hrawd){
+  // 
+  fNEvents++;
+  TH2F* hcorrd=GetCleanEvent(hrawd);
 
   for(Int_t ian=0;ian<fgkNAnodes;ian++){
     if(!fGoodAnode[ian]) continue;
     Float_t sumQ=0.;
+    Int_t cnt=0;
     for(Int_t itb=fFirstGoodTB;itb<=fLastGoodTB;itb++){
-      sumQ+=TMath::Power(hcorrd->GetBinContent(itb+1,ian+1)-fBaseline[ian],2);      
+      sumQ+=TMath::Power(hcorrd->GetBinContent(itb+1,ian+1)-fBaseline[ian],2); 
+      cnt++;    
     }
-    fSumCorrNoise[ian]+=TMath::Sqrt(sumQ/(Float_t)kTimeBins);
+    fSumCorrNoise[ian]+=TMath::Sqrt(sumQ/(Float_t)cnt);
   }
   delete hcorrd;
 }
index 1e1024e..d6de843 100644 (file)
@@ -19,6 +19,7 @@ class AliITSOnlineSDDCMN : public AliITSOnlineSDD {
   AliITSOnlineSDDCMN(Int_t nddl, Int_t ncarlos, Int_t sid);
   virtual ~AliITSOnlineSDDCMN();
   void Reset();
+  TH2F* GetCleanEvent(TH2F* hrawd) const;
   void AddEvent(TH2F* hrawd);
   void ValidateAnodes();
   void ReadBaselines();
index 6490198..8682404 100644 (file)
@@ -23,6 +23,7 @@ class AliITSOnlineSDDInjectors : public AliITSOnlineSDD {
   AliITSOnlineSDDInjectors(Int_t nddl, Int_t ncarlos, Int_t sid);
   virtual ~AliITSOnlineSDDInjectors();
 
+
   void SetThresholds(Float_t tl, Float_t th){
     fLowThreshold=tl;
     fHighThreshold=th;
@@ -31,6 +32,18 @@ class AliITSOnlineSDDInjectors : public AliITSOnlineSDD {
     fTbMin[jlin]=tbmin;
     fTbMax[jlin]=tbmax;
   }
+  void Set20MHzConfig(){
+    SetInjLineRange(0,10,20);
+    SetInjLineRange(1,50,70);
+    SetInjLineRange(2,100,120);
+    SetTimeStep(50.);
+  }
+  void Set40MHzConfig(){
+    SetInjLineRange(0,20,50);
+    SetInjLineRange(1,90,160);
+    SetInjLineRange(2,170,240);
+    SetTimeStep(25.);
+  }
   void SetPolOrder(Int_t n){fPolOrder=n;}
   void SetMinDriftSpeed(Float_t vmin){fMinDriftSpeed=vmin;}
   void SetMaxDriftSpeed(Float_t vmax){fMaxDriftSpeed=vmax;}
index 2ed1c2f..c0f388f 100644 (file)
@@ -23,7 +23,7 @@
 // Origin: F. Prino (prino@to.infn.it)
 
 
-void AnalyzeSDDGainAllMod(Char_t *datafil, Int_t nDDL, Int_t firstEv=10, Int_t lastEv=16, Float_t pascalDAC=100){
+void AnalyzeSDDGainAllMod(Char_t *datafil, Int_t adcfreq=20, Int_t nDDL=0, Int_t firstEv=10, Int_t lastEv=16, Float_t pascalDAC=100){
 
   const Int_t kTotDDL=24;
   const Int_t kModPerDDL=12;
@@ -40,6 +40,8 @@ void AnalyzeSDDGainAllMod(Char_t *datafil, Int_t nDDL, Int_t firstEv=10, Int_t l
       for(Int_t isid=0;isid<kSides;isid++){
        Int_t index=kSides*(kModPerDDL*iddl+imod)+isid;
        anal[index]=new AliITSOnlineSDDTP(iddl,imod,isid,pascalDAC);
+       if(adcfreq==40) anal[index]->SetLastGoodTB(254);
+       else anal[index]->SetLastGoodTB(126);
        sprintf(hisnam,"h%02dc%02ds%d",iddl,imod,isid);
        histo[index]=new TH2F(hisnam,"",256,-0.5,255.5,256,-0.5,255.5);
        isFilled[index]=0;
@@ -222,10 +224,10 @@ void AnalyzeSDDGainAllMod(Char_t *datafil, Int_t nDDL, Int_t firstEv=10, Int_t l
 
 }
 
-void AnalyzeSDDGainAllMod(Int_t nrun, Int_t n2, Char_t* dir="LHC08d_SDD", Int_t nDDL=0, Int_t firstEv=15, Int_t lastEv=20, Float_t pascalDAC=100){
+void AnalyzeSDDGainAllMod(Int_t nrun, Int_t n2, Char_t* dir="LHC08d_SDD",Int_t adcfreq=20, Int_t nDDL=0, Int_t firstEv=15, Int_t lastEv=20, Float_t pascalDAC=100){
   TGrid::Connect("alien:",0,0,"t");
   Char_t filnam[200];
   sprintf(filnam,"alien:///alice/data/2008/%s/%09d/raw/08%09d%03d.10.root",dir,nrun,nrun,n2);
   printf("Open file %s\n",filnam);
-  AnalyzeSDDGainAllMod(filnam,nDDL,firstEv,lastEv,pascalDAC);
+  AnalyzeSDDGainAllMod(filnam,adcfreq,nDDL,firstEv,lastEv,pascalDAC);
 }
index f4ad174..fb67cf1 100644 (file)
@@ -26,7 +26,7 @@
 // Origin: F. Prino (prino@to.infn.it)
 
 
-void AnalyzeSDDInjectorsAllMod(Char_t *datafil, Int_t nDDL, Int_t firstEv=10, Int_t lastEv=15,Int_t jpad=16, Int_t statuscut=7){
+void AnalyzeSDDInjectorsAllMod(Char_t *datafil, Int_t adcfreq=20, Int_t nDDL=0, Int_t firstEv=10, Int_t lastEv=15,Int_t jpad=16, Int_t statuscut=7){
 
   const Int_t kTotDDL=24;
   const Int_t kModPerDDL=12;
@@ -55,12 +55,8 @@ void AnalyzeSDDInjectorsAllMod(Char_t *datafil, Int_t nDDL, Int_t firstEv=10, In
        sprintf(hisnam,"h%02dc%02ds%d",iddl,imod,isid);
        histo[index]=new TH2F(hisnam,"",256,-0.5,255.5,256,-0.5,255.5);
        anal[index]=new AliITSOnlineSDDInjectors(iddl,imod,isid);
-/* Uncomment these lines for analysis of runs with 40 MHz sapling */
-//     anal[index]->SetInjLineRange(0,20,50);
-//     anal[index]->SetInjLineRange(1,90,160);
-//     anal[index]->SetInjLineRange(2,170,240);
-//     anal[index]->SetTimeStep(25.);
-/* END of lines to be uncommented */
+       if(adcfreq==40) anal[index]->Set40MHzConfig();
+       else anal[index]->Set20MHzConfig();
        nWrittenEv[index]=0;
       }
     }
@@ -172,7 +168,7 @@ void AnalyzeSDDInjectorsAllMod(Char_t *datafil, Int_t nDDL, Int_t firstEv=10, In
            gvel[index]->GetXaxis()->SetTitleOffset(0.6);
            gvel[index]->GetYaxis()->SetTitleOffset(0.6);
            if(gvel[index]->GetN()>0) gvel[index]->Draw("AP");
-           Float_t *param=anal[index]->GetDriftSpeedFitParam();
+           Double_t *param=anal[index]->GetDriftSpeedFitParam();
            funz->SetParameters(param[0],param[1],param[2],param[3]);
            funz->SetLineColor(2);
            funz->DrawCopy("LSAME");
@@ -314,12 +310,12 @@ void AnalyzeSDDInjectorsAllMod(Char_t *datafil, Int_t nDDL, Int_t firstEv=10, In
 
 }
 
-void AnalyzeSDDInjectorsAllMod(Int_t nrun, Int_t n2, Char_t* dir="LHC08d_SDD", Int_t nDDL=0, Int_t firstEv=15, Int_t lastEv=15){
+void AnalyzeSDDInjectorsAllMod(Int_t nrun, Int_t n2, Char_t* dir="LHC08d_SDD", Int_t adcfreq=20, Int_t nDDL=0, Int_t firstEv=15, Int_t lastEv=15){
   TGrid::Connect("alien:",0,0,"t");
   Char_t filnam[200];
   sprintf(filnam,"alien:///alice/data/2008/%s/%09d/raw/08%09d%03d.10.root",dir,nrun,nrun,n2);
   printf("Open file %s\n",filnam);
-  AnalyzeSDDInjectorsAllMod(filnam,nDDL,firstEv,lastEv);
+  AnalyzeSDDInjectorsAllMod(filnam,adcfreq,nDDL,firstEv,lastEv);
 }
 
 
index 7ea1aac..b480d20 100644 (file)
@@ -23,7 +23,7 @@
 // All DDLs are analyzed, the argument nDDL selects the DDL to be plotted
 // Origin: F. Prino (prino@to.infn.it)
 
-void AnalyzeSDDNoiseAllMod(Char_t *datafil, Int_t nDDL, Int_t firstEv=10, Int_t lastEv=12){
+void AnalyzeSDDNoiseAllMod(Char_t *datafil, Int_t adcfreq=20, Int_t nDDL=0, Int_t firstEv=10, Int_t lastEv=12){
 
   const Int_t kTotDDL=24;
   const Int_t kModPerDDL=12;
@@ -37,7 +37,9 @@ void AnalyzeSDDNoiseAllMod(Char_t *datafil, Int_t nDDL, Int_t firstEv=10, Int_t
     for(Int_t imod=0; imod<kModPerDDL;imod++){
       for(Int_t isid=0;isid<kSides;isid++){
        Int_t index=kSides*(kModPerDDL*iddl+imod)+isid;
-       base[index]=new AliITSOnlineSDDBase(iddl,imod,isid);
+       base[index]=new AliITSOnlineSDDBase(iddl,imod,isid);    
+       if(adcfreq==40) base[index]->SetLastGoodTB(254);
+       else base[index]->SetLastGoodTB(126);
        sprintf(hisnam,"h%02dc%02ds%d",iddl,imod,isid);
        histo[index]=new TH2F(hisnam,"",128,-0.5,127.5,256,-0.5,255.5);
       }
@@ -127,6 +129,8 @@ void AnalyzeSDDNoiseAllMod(Char_t *datafil, Int_t nDDL, Int_t firstEv=10, Int_t
       for(Int_t isid=0;isid<kSides;isid++){
        Int_t index=kSides*(kModPerDDL*iddl+imod)+isid;
        corr[index]=new AliITSOnlineSDDCMN(iddl,imod,isid);
+       if(adcfreq==40) corr[index]->SetLastGoodTB(254);
+       else corr[index]->SetLastGoodTB(126);
        isFilled[index]=0;
       }
     }
@@ -369,10 +373,10 @@ void AnalyzeSDDNoiseAllMod(Char_t *datafil, Int_t nDDL, Int_t firstEv=10, Int_t
 
 }
 
-void AnalyzeSDDNoiseAllMod(Int_t nrun, Int_t n2, Char_t* dir="LHC08d_SDD",Int_t nDDL=0, Int_t firstEv=15, Int_t lastEv=18){
+void AnalyzeSDDNoiseAllMod(Int_t nrun, Int_t n2, Char_t* dir="LHC08d_SDD",Int_t adcfreq=20, Int_t nDDL=0, Int_t firstEv=15, Int_t lastEv=18){
   TGrid::Connect("alien:",0,0,"t");
   Char_t filnam[200];
   sprintf(filnam,"alien:///alice/data/2008/%s/%09d/raw/08%09d%03d.10.root",dir,nrun,nrun,n2);
   printf("Open file %s\n",filnam);
-  AnalyzeSDDNoiseAllMod(filnam,nDDL,firstEv,lastEv);
+  AnalyzeSDDNoiseAllMod(filnam,adcfreq,nDDL,firstEv,lastEv);
 }