Update of calibration procedures
authorcoppedis <coppedis@f7af4fe6-9843-0410-8265-dc069ae4e863>
Fri, 20 Jan 2006 15:18:44 +0000 (15:18 +0000)
committercoppedis <coppedis@f7af4fe6-9843-0410-8265-dc069ae4e863>
Fri, 20 Jan 2006 15:18:44 +0000 (15:18 +0000)
ZDC/AliZDCCalibData.cxx
ZDC/AliZDCCalibData.h
ZDC/AliZDCDigitizer.cxx
ZDC/AliZDCDigitizer.h
ZDC/AliZDCReconstructor.cxx
ZDC/AliZDCReconstructor.h
ZDC/Calib/Data/Run0_100_v0_s0.root

index b12463d..7535cf3 100644 (file)
@@ -95,10 +95,8 @@ void AliZDCCalibData::PrepHistos()
   Int_t   kNChannels = 47;
   Float_t kMaxPedVal = 47.;
   TString hname = GetName();  hname += "_Pedestals";
-  fHistMeanPed   = new TH1F(hname.Data(),hname.Data(),kNChannels,0.,kMaxPedVal);
-  for (int i=0; i<47; i++) {
-    fHistMeanPed->SetBinContent(i+1,GetMeanPed(i));
-  }
+  fHistMeanPed = new TH1F(hname.Data(),hname.Data(),kNChannels,0.,kMaxPedVal);
+  for(int i=0; i<47; i++)  fHistMeanPed->SetBinContent(i+1,GetMeanPed(i));
 }
 
 //________________________________________________________________
@@ -118,7 +116,7 @@ void  AliZDCCalibData::Print(Option_t *) const
    }
  
    printf("\n  ----    Energy calibration coefficients         ----\n\n");
-   for(int t=0; t<4; t++){
+   for(int t=0; t<3; t++){
       printf(" En Calib Coeff. [ZDC%d] = %f\n",t,fEnCalibration[t]);
    }
 } 
index 89ca614..1ed3bef 100644 (file)
@@ -29,9 +29,9 @@ class AliZDCCalibData: public TNamed {
   Float_t  GetEnCalib(Int_t channel)   const {return fEnCalibration[channel];}
   Float_t* GetEnCalib()   const {return (float*)fEnCalibration;}
   //
-  void     SetMeanPed(Float_t val, Int_t channel) {fMeanPedestal[channel]=val;}
+  void     SetMeanPed(Int_t channel, Float_t val) {fMeanPedestal[channel]=val;}
   void     SetMeanPed(Float_t* MeanPed);
-  void            SetEnCalib(Float_t val, Int_t channel) {fEnCalibration[channel]=val;}
+  void            SetEnCalib(Int_t channel, Float_t val) {fEnCalibration[channel]=val;}
   void            SetEnCalib(Float_t* EnCalib);
   //
   void     PrepHistos();
@@ -40,7 +40,7 @@ class AliZDCCalibData: public TNamed {
 
  protected:
   Float_t  fMeanPedestal[47];  // Mean pedestal values
-  Float_t  fEnCalibration[4];  // Coeff. for energy calibration (4 different ZDC's?)
+  Float_t  fEnCalibration[3];  // Coeff. for energy calibration
   TH1F*    fHistMeanPed;        //! histos for drawing mean pedestals
   //
   ClassDef(AliZDCCalibData,1)    // ZDC Sensor Calibration data
index e010758..cc60253 100644 (file)
@@ -36,6 +36,7 @@
 #include "AliGenHijingEventHeader.h"
 #include "AliRunDigitizer.h"
 #include "AliRunLoader.h"
+#include "AliCDBManager.h"
 #include "AliCDBStorage.h"
 #include "AliCDBEntry.h"
 #include "AliZDCSDigit.h"
@@ -58,8 +59,14 @@ AliZDCDigitizer::AliZDCDigitizer()
 AliZDCDigitizer::AliZDCDigitizer(AliRunDigitizer* manager):
   AliDigitizer(manager)
 {
+
   // Constructor    
-  // if(!fStorage) fStorage =  AliCDBManager::Instance()->GetStorage("local://DBlocal");
+  fStorage =  SetStorage("local://$ALICE_ROOT");
+
+  // Get calibration data
+  int runNumber = 0;
+  fCalibData = GetCalibData(runNumber); 
+
 }
 
 //____________________________________________________________________________
@@ -153,9 +160,9 @@ void AliZDCDigitizer::Exec(Option_t* /*option*/)
     // 
     specN = ((AliGenHijingEventHeader*) genHeader)->ProjSpectatorsn();
     specP = ((AliGenHijingEventHeader*) genHeader)->ProjSpectatorsp();
-    AliDebug(2, Form("\n b = %f fm, Nspecn = %d, Nspecp = %d\n",
+    AliDebug(2, Form("\n AliZDCDigitizer -> b = %f fm, Nspecn = %d, Nspecp = %d\n",
                      impPar, specN, specP));
-   printf("\n\t **** b = %f fm, Nspecn = %d, Nspecp = %d\n",
+   printf("\n\t AliZDCDigitizer -> b = %f fm, Nspecn = %d, Nspecp = %d\n",
                      impPar, specN, specP);
   }
 
@@ -163,9 +170,9 @@ void AliZDCDigitizer::Exec(Option_t* /*option*/)
   if (impPar >= 0) {
     Int_t freeSpecN, freeSpecP;
     Fragmentation(impPar, specN, specP, freeSpecN, freeSpecP);
-    printf("\n\t ---- Adding signal for %d free spectator n\n",freeSpecN);
+    printf("\n\t AliZDCDigitizer ---- Adding signal for %d free spectator n\n",freeSpecN);
     SpectatorSignal(1, freeSpecN, pm);
-    printf("\t ---- Adding signal for %d free spectator p\n\n",freeSpecP);
+    printf("\t AliZDCDigitizer ---- Adding signal for %d free spectator p\n\n",freeSpecP);
     SpectatorSignal(2, freeSpecP, pm);
   }
 
@@ -198,8 +205,8 @@ void AliZDCDigitizer::Exec(Option_t* /*option*/)
         for (Int_t res = 0; res < 2; res++){
           digi[res] = Phe2ADCch(sector[0], sector[1], pm[sector[0]-1][sector[1]], res) 
                    + Pedestal(sector[0], sector[1], res);
-         printf("\t DIGIT added -> det = %d, quad = %d - digi[%d] = %d\n\n",
-             sector[0], sector[1], res, digi[res]); // Chiara debugging!
+         //printf("\t DIGIT added -> det = %d, quad = %d - digi[%d] = %d\n\n",
+         //    sector[0], sector[1], res, digi[res]); // Chiara debugging!
        }
         new(pdigit) AliZDCDigit(sector, digi);
         treeD->Fill();
@@ -283,7 +290,7 @@ void AliZDCDigitizer::SpectatorSignal(Int_t SpecType, Int_t numEvents,
          //
          Float_t lightQ = hitsSpec[4];
          Float_t lightC = hitsSpec[5];
-         AliDebug(3, Form("Volume = (%d, %d), lightQ = %.0f, lightC = %.0f",
+         AliDebug(3, Form("SpectatorSignal -> vol = (%d, %d), lightQ = %.0f, lightC = %.0f",
                            volume[0], volume[1], lightQ, lightC));
          //printf("\n   Volume = (%d, %d), lightQ = %.0f, lightC = %.0f",
           //                 volume[0], volume[1], lightQ, lightC);
@@ -317,45 +324,60 @@ Int_t AliZDCDigitizer::Phe2ADCch(Int_t Det, Int_t Quad, Float_t Light,
 {
   // Evaluation of the ADC channel corresponding to the light yield Light
   Int_t ADCch = (Int_t) (Light * fPMGain[Det-1][Quad] * fADCRes[Res]);
-  printf("\t Phe2ADCch -> det %d quad %d - phe %.0f  ADC %d\n", Det,Quad,Light,ADCch);
+  //printf("\t Phe2ADCch -> det %d quad %d - phe %.0f  ADC %d\n", Det,Quad,Light,ADCch);
   return ADCch;
 }
 
 //_____________________________________________________________________________
 Int_t AliZDCDigitizer::Pedestal(Int_t Det, Int_t Quad, Int_t Res) const
 {
-  // Get calibration data
-  int runNumber = 0;
-  AliZDCCalibData *calibda = GetCalibData(runNumber);  
-  //calibda->Print("");
   
   Float_t meanPed;
-  if(Det != 3) meanPed = calibda->GetMeanPed(10*(Det-1)+Quad+5*Res);
-  else         meanPed = calibda->GetMeanPed(10*(Det-1)+Quad+1*Res);
+  if(Det != 3) meanPed = fCalibData->GetMeanPed(10*(Det-1)+Quad+5*Res);
+  else         meanPed = fCalibData->GetMeanPed(10*(Det-1)+Quad+1*Res);
   
-  printf("\t Pedestal -> det = %d, quad = %d, res = %d - Ped[%d] = %d\n",
-       Det, Quad, Res,10*(Det-1)+Quad+5*Res,(Int_t) meanPed); // Chiara debugging!
+  //printf("\t Pedestal -> det = %d, quad = %d, res = %d - Ped[%d] = %d\n",
+  //   Det, Quad, Res,10*(Det-1)+Quad+5*Res,(Int_t) meanPed); // Chiara debugging!
   
   return (Int_t) meanPed;
 }
 
 //_____________________________________________________________________________
-AliZDCCalibData* AliZDCDigitizer::GetCalibData(int runNumber) const
+AliCDBStorage* AliZDCDigitizer::SetStorage(const char *uri) 
 {
+  //printf("\n\t AliZDCDigitizer::SetStorage \n");
 
-  printf("\t AliZDCReconstructor::GetCalibData for RUN #%d\n",runNumber);
-  //fStorage->PrintSelectionList();
-  //AliCDBEntry *entry = fStorage->Get("ZDC/Calib/Data",runNumber);
+  Bool_t deleteManager = kFALSE;
   
-  AliCDBStorage *fStorage = AliCDBManager::Instance()->GetStorage("local://$ALICE_ROOT");
-  AliCDBEntry  *entry = fStorage->Get("ZDC/Calib/Data",0);
+  AliCDBManager *manager = AliCDBManager::Instance();
+  AliCDBStorage *defstorage = manager->GetDefaultStorage();
   
-  AliZDCCalibData *calibda = 0;
-  if (entry) calibda = (AliZDCCalibData*) entry->GetObject();
-  //calibda->Print("");
+  if(!defstorage || !(defstorage->Contains("ZDC"))){ 
+     AliWarning("No default storage set or default storage doesn't contain ZDC!");
+     manager->SetDefaultStorage(uri);
+     deleteManager = kTRUE;
+  }
+  AliCDBStorage *storage = manager->GetDefaultStorage();
+
+  if(deleteManager){
+    AliCDBManager::Instance()->UnsetDefaultStorage();
+    defstorage = 0;   // the storage is killed by AliCDBManager::Instance()->Destroy()
+  }
+
+  return storage; 
+}
 
- // AliCDBManager::Instance()->Destroy();
+//_____________________________________________________________________________
+AliZDCCalibData* AliZDCDigitizer::GetCalibData(int runNumber) const
+{
 
-  return calibda;
+  //printf("\n\t AliZDCDigitizer::GetCalibData \n");
+      
+  AliCDBEntry  *entry = fStorage->Get("ZDC/Calib/Data",runNumber);  
+  AliZDCCalibData *calibdata = (AliZDCCalibData*) entry->GetObject();
+    
+  if (!calibdata)  AliWarning("No calibration data from calibration database !");
 
+  return calibdata;
 }
index c2a2dc9..ba786d0 100644 (file)
@@ -44,6 +44,7 @@ public:
   Float_t GetADCRes(Int_t i) const {return fADCRes[i];}
   
   void   GetStorage(const char* uri) {fStorage = AliCDBManager::Instance()->GetStorage(uri);}
+  AliCDBStorage   *SetStorage(const char* uri);
   AliZDCCalibData *GetCalibData(int runNumber) const; 
 
 private:
@@ -56,11 +57,11 @@ private:
                     Int_t Res) const;
   Int_t   Pedestal(Int_t Detector, Int_t Quadrant, Int_t Res) const;
 
-  Float_t fPMGain[3][5];      // PM gain
-  Float_t fADCRes[2];        // ADC conversion factors
+  Float_t fPMGain[3][5];       // PM gain
+  Float_t fADCRes[2];          // ADC conversion factors
   
-  AliCDBStorage *fStorage; //! storage
-
+  AliCDBStorage *fStorage;     //! storage
+  AliZDCCalibData *fCalibData;         //! calibration data
        
   ClassDef(AliZDCDigitizer, 3)     // digitizer for ZDC
 };    
index 6bc6675..331cf51 100644 (file)
@@ -41,7 +41,7 @@ ClassImp(AliZDCReconstructor)
 AliZDCReconstructor:: AliZDCReconstructor()
 {
   // **** Default constructor
-  // if(!fStorage) fStorage =  AliCDBManager::Instance()->GetStorage("local://DBlocal");
+  fStorage =  0;
 
   //  ---      Number of generated spectator nucleons and impact parameter
   // --------------------------------------------------------------------------------------------------
@@ -73,6 +73,13 @@ AliZDCReconstructor:: AliZDCReconstructor()
   fZEMp  = new TF1("fZEMp","82.49-0.03611*x+0.00000385*x*x",0.,4000.);
   fZEMsp = new TF1("fZEMsp","208.7-0.09006*x+0.000009526*x*x",0.,4000.);
   fZEMb  = new TF1("fZEMb","16.06-0.01633*x+1.44e-5*x*x-6.778e-9*x*x*x+1.438e-12*x*x*x*x-1.112e-16*x*x*x*x*x",0.,4000.);
+  
+  // Setting storage
+  fStorage =  SetStorage("local://$ALICE_ROOT");
+
+  // Get calibration data
+  int runNumber = 0;
+  fCalibData = GetCalibData(runNumber); 
 }
 
 //_____________________________________________________________________________
@@ -119,13 +126,9 @@ AliZDCReconstructor::~AliZDCReconstructor()
 void AliZDCReconstructor::Reconstruct(AliRunLoader* runLoader) const
 {
   // *** Local ZDC reconstruction for digits
-  
-  // Get calibration data
-  int runNumber = 0;
-  AliZDCCalibData *calibda = GetCalibData(runNumber);  
-  
+    
   Float_t meanPed[47];
-  for(Int_t jj=0; jj<47; jj++) meanPed[jj] = calibda->GetMeanPed(jj);
+  for(Int_t jj=0; jj<47; jj++) meanPed[jj] = fCalibData->GetMeanPed(jj);
 
   AliLoader* loader = runLoader->GetLoader("ZDCLoader");
   if (!loader) return;
@@ -151,12 +154,14 @@ void AliZDCReconstructor::Reconstruct(AliRunLoader* runLoader) const
       if (!pdigit) continue;
 
       if(digit.GetSector(0) == 1)
-                zncorr  += (Float_t) (digit.GetADCValue(0)-meanPed[digit.GetSector(1)]);          // ped 4 high gain ZN ADCs
+                zncorr  += (Float_t) (digit.GetADCValue(0)-meanPed[digit.GetSector(1)]);     // ped 4 high gain ZN ADCs
       else if(digit.GetSector(0) == 2)
-        zpcorr  += (Float_t) (digit.GetADCValue(0)-meanPed[digit.GetSector(1)+10]); // ped 4 high gain ZP ADCs
+        zpcorr  += (Float_t) (digit.GetADCValue(0)-meanPed[digit.GetSector(1)+10]);  // ped 4 high gain ZP ADCs
       else if(digit.GetSector(0) == 3){
-        if(digit.GetSector(1)==1)      zemcorr += (Float_t) (digit.GetADCValue(0)-meanPed[digit.GetSector(1)+20]); // ped 4 high gain ZEM1 ADCs
-        else if(digit.GetSector(1)==2) zemcorr += (Float_t) (digit.GetADCValue(0)-meanPed[digit.GetSector(1)+22]); // ped 4 high gain ZEM2 ADCs
+        if(digit.GetSector(1)==1)      
+          zemcorr += (Float_t) (digit.GetADCValue(0)-meanPed[digit.GetSector(1)+20]); // ped 4 high gain ZEM1 ADCs
+        else if(digit.GetSector(1)==2) 
+          zemcorr += (Float_t) (digit.GetADCValue(0)-meanPed[digit.GetSector(1)+22]); // ped 4 high gain ZEM2 ADCs
       }
     }
     if(zncorr<0)  zncorr=0;
@@ -164,7 +169,7 @@ void AliZDCReconstructor::Reconstruct(AliRunLoader* runLoader) const
     if(zemcorr<0) zemcorr=0;
 
     // reconstruct the event
-    printf("\n \t ZDCReco from digits-> Ev.#%d ZN = %.0f, ZP = %.0f, ZEM = %.0f\n",iEvent,zncorr,zpcorr,zemcorr);
+    //printf("\n \t ZDCReco from digits-> Ev.#%d ZN = %.0f, ZP = %.0f, ZEM = %.0f\n",iEvent,zncorr,zpcorr,zemcorr);
     ReconstructEvent(loader, zncorr, zpcorr, zemcorr);
   }
 
@@ -178,12 +183,8 @@ void AliZDCReconstructor::Reconstruct(AliRunLoader* runLoader,
 {
   // *** Local ZDC reconstruction for raw data
   
-  // Calibration data
-  int runNumber = 0;
-  AliZDCCalibData *calibda = GetCalibData(runNumber);  
-  
   Float_t meanPed[47];
-  for(Int_t jj=0; jj<47; jj++) meanPed[jj] = calibda->GetMeanPed(jj);
+  for(Int_t jj=0; jj<47; jj++) meanPed[jj] = fCalibData->GetMeanPed(jj);
 
   AliLoader* loader = runLoader->GetLoader("ZDCLoader");
   if (!loader) return;
@@ -210,7 +211,7 @@ void AliZDCReconstructor::Reconstruct(AliRunLoader* runLoader,
     if(zemcorr<0) zemcorr=0;
     
     // reconstruct the event
-    printf("\n\t ZDCReco from raw-> Ev.#%d ZN = %.0f, ZP = %.0f, ZEM = %.0f\n",iEvent,zncorr,zpcorr,zemcorr);
+    //printf("\n\t ZDCReco from raw-> Ev.#%d ZN = %.0f, ZP = %.0f, ZEM = %.0f\n",iEvent,zncorr,zpcorr,zemcorr);
     ReconstructEvent(loader, zncorr, zpcorr, zemcorr);
   }
 
@@ -250,7 +251,6 @@ void AliZDCReconstructor::ReconstructEvent(AliLoader* loader, Float_t zncorr,
   //  if AliDebug(1,Form("    znenergy = %f TeV, zpenergy = %f TeV, zdcenergy = %f GeV, "
   //                      "\n          zemenergy = %f TeV\n", znenergy, zpenergy, 
   //                      zdcenergy, zemenergy);
-  
   //  if(zdcenergy==0)
   //    if AliDebug(1,Form("\n\n       ###     ATTENZIONE!!! -> ev# %d: znenergy = %f TeV, zpenergy = %f TeV, zdcenergy = %f GeV, "
   //                        " zemenergy = %f TeV\n\n", fMerger->EvNum(), znenergy, zpenergy, zdcenergy, zemenergy); 
@@ -303,8 +303,8 @@ void AliZDCReconstructor::ReconstructEvent(AliLoader* loader, Float_t zncorr,
   Int_t nPart, nPartTot;
   nPart = 207-nGenSpecN-nGenSpecP;
   nPartTot = 207-nGenSpec;
-  printf("\t  ZDCeventReco-> ZNEn = %.0f GeV, ZPEn = %.0f GeV, ZEMEn = %.0f GeV\n",
-       znenergy, zpenergy, zemenergy);
+  //printf("\t  ZDCeventReco-> ZNEn = %.0f GeV, ZPEn = %.0f GeV, ZEMEn = %.0f GeV\n",
+  //   znenergy, zpenergy, zemenergy);
 
   // create the output tree
   loader->MakeTree("R");
@@ -345,21 +345,41 @@ void AliZDCReconstructor::FillESD(AliRunLoader* runLoader,
 }
 
 //_____________________________________________________________________________
-AliZDCCalibData* AliZDCReconstructor::GetCalibData(int runNumber) const
+AliCDBStorage* AliZDCReconstructor::SetStorage(const char *uri) 
 {
+  //printf("\n\t AliZDCReconstructor::SetStorage \n");
 
-  //printf("\n\t AliZDCReconstructor::GetCalibData \n");
-  //fStorage->PrintSelectionList();
-  //AliCDBEntry *entry = fStorage->Get("DBlocal/ZDC/Calib/Data",runNumber);
-
+  Bool_t deleteManager = kFALSE;
   
-  AliCDBStorage *fStorage = AliCDBManager::Instance()->GetStorage("local://$ALICE_ROOT");
-  AliCDBEntry  *entry = fStorage->Get("ZDC/Calib/Data",0);
+  AliCDBManager *manager = AliCDBManager::Instance();
+  AliCDBStorage *defstorage = manager->GetDefaultStorage();
   
-  AliZDCCalibData *calibda = (AliZDCCalibData*) entry->GetObject();
+  if(!defstorage || !(defstorage->Contains("ZDC"))){ 
+     AliWarning("No default storage set or default storage doesn't contain ZDC!");
+     manager->SetDefaultStorage(uri);
+     deleteManager = kTRUE;
+  }
+  AliCDBStorage *storage = manager->GetDefaultStorage();
+
+  if(deleteManager){
+    AliCDBManager::Instance()->UnsetDefaultStorage();
+    defstorage = 0;   // the storage is killed by AliCDBManager::Instance()->Destroy()
+  }
+
+  return storage; 
+}
 
-  //AliCDBManager::Instance()->Destroy();
+//_____________________________________________________________________________
+AliZDCCalibData* AliZDCReconstructor::GetCalibData(int runNumber) const
+{
 
-  return calibda;
+  //printf("\n\t AliZDCReconstructor::GetCalibData \n");
+      
+  AliCDBEntry  *entry = fStorage->Get("ZDC/Calib/Data",runNumber);  
+  AliZDCCalibData *calibdata = (AliZDCCalibData*) entry->GetObject();
+    
+  if (!calibdata)  AliWarning("No calibration data from calibration database !");
 
+  return calibdata;
 }
index 79e57ad..7d42e26 100644 (file)
@@ -40,6 +40,7 @@ public:
                  {AliReconstructor::FillESD(runLoader,rawReader,esd);}
   
   void   GetStorage(const char* uri) {fStorage = AliCDBManager::Instance()->GetStorage(uri);}
+  AliCDBStorage   *SetStorage(const char* uri);
   AliZDCCalibData *GetCalibData(int runNumber) const; 
 
 private:
@@ -61,7 +62,8 @@ private:
   TF1*   fZEMsp;     //! Nspectators from ZEM energy
   TF1*   fZEMb;      //! b from ZEM energy
   
-  AliCDBStorage *fStorage; //! storage
+  AliCDBStorage *fStorage;     //! storage
+  AliZDCCalibData *fCalibData;         //! calibration data
 
   ClassDef(AliZDCReconstructor, 1)   // class for the ZDC reconstruction
 };
index 61525ca..fb0d2f7 100644 (file)
Binary files a/ZDC/Calib/Data/Run0_100_v0_s0.root and b/ZDC/Calib/Data/Run0_100_v0_s0.root differ