Add offline preprocessor for pass0 calibration (Raphaelle)
authorcblume <cblume@f7af4fe6-9843-0410-8265-dc069ae4e863>
Fri, 28 May 2010 19:46:06 +0000 (19:46 +0000)
committercblume <cblume@f7af4fe6-9843-0410-8265-dc069ae4e863>
Fri, 28 May 2010 19:46:06 +0000 (19:46 +0000)
TRD/AliTRDPreprocessorOffline.cxx [new file with mode: 0644]
TRD/AliTRDPreprocessorOffline.h [new file with mode: 0644]
TRD/Cal/AliTRDCalROC.cxx
TRD/TRDcalibLinkDef.h
TRD/libTRDcalib.pkg

diff --git a/TRD/AliTRDPreprocessorOffline.cxx b/TRD/AliTRDPreprocessorOffline.cxx
new file mode 100644 (file)
index 0000000..91c622e
--- /dev/null
@@ -0,0 +1,575 @@
+/**************************************************************************
+ * Copyright(c) 1998-1999, ALICE Experiment at CERN, All rights reserved. *
+ *                                                                        *
+ * Author: The ALICE Off-line Project.                                    *
+ * Contributors are mentioned in the code where appropriate.              *
+ *                                                                        *
+ * Permission to use, copy, modify and distribute this software and its   *
+ * documentation strictly for non-commercial purposes is hereby granted   *
+ * without fee, provided that the above copyright notice appears in all   *
+ * copies and that both the copyright notice and this permission notice   *
+ * appear in the supporting documentation. The authors make no claims     *
+ * about the suitability of this software for any purpose. It is          *
+ * provided "as is" without express or implied warranty.                  *
+ **************************************************************************/
+
+
+
+/*
+  Responsible: marian.ivanov@cern.ch 
+  Code to analyze the TRD calibration and to produce OCDB entries  
+
+
+   .x ~/rootlogon.C
+   gSystem->Load("libANALYSIS");
+   gSystem->Load("libTRDcalib");
+
+   AliTRDPreprocessorOffline proces;
+   TString ocdbPath="local:////"
+   ocdbPath+=gSystem->GetFromPipe("pwd");
+
+   proces.CalibTimeGain("CalibObjects.root",run0,run1,ocdbPath);
+   proces.CalibTimeVdrift("CalibObjects.root",run0,run1,ocdbPath);
+  // take the raw calibration data from the file CalibObjects.root 
+  // and make a OCDB entry with run  validity run0-run1
+  // results are stored at the ocdbPath - local or alien ...
+  // default storage ""- data stored at current working directory 
+*/
+#include "Riostream.h"
+#include <fstream>
+#include "TFile.h"
+#include "TCanvas.h"
+#include "TLegend.h"
+#include "TH2I.h"
+#include "TH1I.h"
+#include "TH2F.h"
+#include "TH1F.h"
+#include "TProfile2D.h"
+#include "AliTRDCalDet.h"
+#include "AliTRDCalPad.h"
+#include "AliCDBMetaData.h"
+#include "AliCDBId.h"
+#include "AliCDBManager.h"
+#include "AliCDBStorage.h"
+#include "AliTRDCalibraMode.h"
+#include "AliTRDCalibraFit.h"
+#include "AliTRDCalibraVdriftLinearFit.h"
+#include "AliTRDPreprocessorOffline.h"
+
+
+ClassImp(AliTRDPreprocessorOffline)
+
+AliTRDPreprocessorOffline::AliTRDPreprocessorOffline():
+  TNamed("TRDPreprocessorOffline","TRDPreprocessorOffline"),
+  fMethodSecond(kTRUE),
+  fCH2d(0x0),
+  fPH2d(0x0),
+  fPRF2d(0x0),
+  fAliTRDCalibraVdriftLinearFit(0x0),
+  fNEvents(0x0),
+  fAbsoluteGain(0x0),
+  fPlots(new TObjArray(8)),
+  fCalibObjects(new TObjArray(8))
+{
+  //
+  // default constructor
+  //
+}
+
+AliTRDPreprocessorOffline::~AliTRDPreprocessorOffline() {
+  //
+  // Destructor
+  //
+
+  if(fCH2d) delete fCH2d;
+  if(fPH2d) delete fPH2d;
+  if(fPRF2d) delete fPRF2d;
+  if(fAliTRDCalibraVdriftLinearFit) delete fAliTRDCalibraVdriftLinearFit;
+  if(fNEvents) delete fNEvents;
+  if(fAbsoluteGain) delete fAbsoluteGain;
+  if(fPlots) delete fPlots;
+  if(fCalibObjects) delete fCalibObjects;
+  
+}
+
+void AliTRDPreprocessorOffline::CalibVdriftT0(const Char_t* file, Int_t startRunNumber, Int_t endRunNumber, TString ocdbStorage){
+  //
+  // make calibration of the drift velocity
+  // Input parameters:
+  //      file                             - the location of input file
+  //      startRunNumber, endRunNumber     - run validity period 
+  //      ocdbStorage                      - path to the OCDB storage
+  //                                       - if empty - local storage 'pwd' uesed
+  if (ocdbStorage.Length()<=0) ocdbStorage="local://"+gSystem->GetFromPipe("pwd")+"/OCDB";
+  //
+  // 1. Initialization 
+  //
+  //
+  // 2. extraction of the information
+  //
+  if(ReadVdriftT0Global(file)) AnalyzeVdriftT0();
+  if(ReadVdriftLinearFitGlobal(file)) AnalyzeVdriftLinearFit();
+  //
+  // 3. Append QA plots
+  //
+  //MakeDefaultPlots(fVdriftArray,fVdriftArray);
+  //
+  //
+  // 4. update of OCDB
+  //
+  //
+  UpdateOCDBVdrift(startRunNumber,endRunNumber,ocdbStorage);
+  UpdateOCDBT0(startRunNumber,endRunNumber,ocdbStorage);
+
+}
+
+
+void AliTRDPreprocessorOffline::CalibGain(const Char_t* file, Int_t startRunNumber, Int_t endRunNumber, TString ocdbStorage){
+  //
+  // make calibration of the drift velocity
+  // Input parameters:
+  //      file                             - the location of input file
+  //      startRunNumber, endRunNumber     - run validity period 
+  //      ocdbStorage                      - path to the OCDB storage
+  //                                       - if empty - local storage 'pwd' uesed
+  if (ocdbStorage.Length()<=0) ocdbStorage="local://"+gSystem->GetFromPipe("pwd")+"/OCDB";
+  //
+  // 1. Initialization 
+  if(!ReadGainGlobal(file)) return;
+  //
+  //
+  // 2. extraction of the information
+  //
+  AnalyzeGain();
+  //
+  // 3. Append QA plots
+  //
+  //MakeDefaultPlots(fVdriftArray,fVdriftArray);
+  //
+  //
+  // 4. update of OCDB
+  //
+  //
+  UpdateOCDBGain(startRunNumber,endRunNumber,ocdbStorage);
+  
+}
+
+
+void AliTRDPreprocessorOffline::CalibPRF(const Char_t* file, Int_t startRunNumber, Int_t endRunNumber, TString ocdbStorage){
+  //
+  // make calibration of the drift velocity
+  // Input parameters:
+  //      file                             - the location of input file
+  //      startRunNumber, endRunNumber     - run validity period 
+  //      ocdbStorage                      - path to the OCDB storage
+  //                                       - if empty - local storage 'pwd' uesed
+  if (ocdbStorage.Length()<=0) ocdbStorage="local://"+gSystem->GetFromPipe("pwd")+"/OCDB";
+  //
+  // 1. Initialization 
+  if(!ReadPRFGlobal(file)) return;
+  //
+  //
+  // 2. extraction of the information
+  //
+  AnalyzePRF();
+  //
+  // 3. Append QA plots
+  //
+  //MakeDefaultPlots(fVdriftArray,fVdriftArray);
+  //
+  //
+  // 4. update of OCDB
+  //
+  //
+  UpdateOCDBPRF(startRunNumber,endRunNumber,ocdbStorage);
+  
+}
+
+
+Bool_t AliTRDPreprocessorOffline::ReadGainGlobal(const Char_t* fileName){
+  //
+  // read calibration entries from file
+  // 
+  TFile fcalib(fileName);
+  TObjArray * array = (TObjArray*)fcalib.Get("TRDCalib");
+  if (array){
+    TH2I *ch2d = (TH2I *) array->FindObject("CH2d");
+    if(!ch2d) return kFALSE;
+    fCH2d = (TH2I*)ch2d->Clone();
+    //fNEvents = (TH1I *) array->FindObject("NEvents");
+    //fAbsoluteGain = (TH2F *) array->FindObject("AbsoluteGain");
+  }else{
+    TH2I *ch2d = (TH2I *) fcalib.Get("CH2d");
+    if(!ch2d) return kFALSE;
+    fCH2d = (TH2I*)ch2d->Clone();
+    //fNEvents = (TH1I *) fcalib.Get("NEvents");
+    //fAbsoluteGain = (TH2F *) fcalib.Get("AbsoluteGain");
+  }
+  fCH2d->SetDirectory(0);
+  //printf("title of CH2d %s\n",fCH2d->GetTitle());
+
+  return kTRUE;
+  
+}
+
+Bool_t AliTRDPreprocessorOffline::ReadVdriftT0Global(const Char_t* fileName){
+  //
+  // read calibration entries from file
+  // 
+  TFile fcalib(fileName);
+  TObjArray * array = (TObjArray*)fcalib.Get("TRDCalib");
+  if (array){
+    TProfile2D *ph2d = (TProfile2D *) array->FindObject("PH2d");
+    if(!ph2d) return kFALSE;
+    fPH2d = (TProfile2D*)ph2d->Clone();
+    //fNEvents = (TH1I *) array->FindObject("NEvents");
+  }else{
+    TProfile2D *ph2d = (TProfile2D *) fcalib.Get("PH2d");
+    if(!ph2d) return kFALSE;
+    fPH2d = (TProfile2D*)ph2d->Clone();
+    //fNEvents = (TH1I *) fcalib.Get("NEvents");
+  }
+  fPH2d->SetDirectory(0);
+  //printf("title of PH2d %s\n",fPH2d->GetTitle());
+  
+  return kTRUE;
+  
+}
+
+Bool_t AliTRDPreprocessorOffline::ReadVdriftLinearFitGlobal(const Char_t* fileName){
+  //
+  // read calibration entries from file
+  // 
+  TFile fcalib(fileName);
+  TObjArray * array = (TObjArray*)fcalib.Get("TRDCalib");
+  if (array){
+    fAliTRDCalibraVdriftLinearFit = (AliTRDCalibraVdriftLinearFit *) array->FindObject("AliTRDCalibraVdriftLinearFit");
+    //fNEvents = (TH1I *) array->FindObject("NEvents");
+  }else{
+    fAliTRDCalibraVdriftLinearFit = (AliTRDCalibraVdriftLinearFit *) fcalib.Get("AliTRDCalibraVdriftLinearFit");
+    //fNEvents = (TH1I *) fcalib.Get("NEvents");
+  }
+  if(!fAliTRDCalibraVdriftLinearFit) {
+    //printf("No AliTRDCalibraVdriftLinearFit\n");
+    return kFALSE;
+  }
+  return kTRUE;
+  
+}
+
+Bool_t AliTRDPreprocessorOffline::ReadPRFGlobal(const Char_t* fileName){
+  //
+  // read calibration entries from file
+  // 
+  TFile fcalib(fileName);
+  TObjArray * array = (TObjArray*)fcalib.Get("TRDCalib");
+  if (array){
+    TProfile2D *prf2d = (TProfile2D *) array->FindObject("PRF2d");
+    if(!prf2d) return kFALSE;
+    fPRF2d = (TProfile2D*)prf2d->Clone();
+    //fNEvents = (TH1I *) array->FindObject("NEvents");
+  }else{
+    TProfile2D *prf2d = (TProfile2D *) fcalib.Get("PRF2d");
+    if(!prf2d) return kFALSE;
+    fPRF2d = (TProfile2D*)prf2d->Clone();
+    //fNEvents = (TH1I *) fcalib.Get("NEvents");
+  }
+  fPRF2d->SetDirectory(0);
+  //printf("title of PRF2d %s\n",fPRF2d->GetTitle());
+  
+  return kTRUE;
+
+}
+
+
+
+Bool_t AliTRDPreprocessorOffline::AnalyzeGain(){
+  //
+  // Analyze gain - produce the calibration objects
+  //
+
+  AliTRDCalibraFit *calibra = AliTRDCalibraFit::Instance();
+  calibra->SetMinEntries(800); // If there is less than 1000 entries in the histo: no fit
+  calibra->AnalyseCH(fCH2d);
+
+  Int_t nbtg = 6*4*18*((Int_t) ((AliTRDCalibraMode *)calibra->GetCalibraMode())->GetDetChamb0(0))
+    + 6*  18*((Int_t) ((AliTRDCalibraMode *)calibra->GetCalibraMode())->GetDetChamb2(0));
+  Int_t nbfit       = calibra->GetNumberFit();
+  Int_t nbE         = calibra->GetNumberEnt();
+
+
+  Bool_t ok = kFALSE;
+  // enough statistics
+  if ((nbtg >                  0) && 
+      (nbfit        >= 0.5*nbE) && (nbE > 30)) {
+    // create the cal objects
+    calibra->PutMeanValueOtherVectorFit(1,kTRUE);
+    TObjArray object           = calibra->GetVectorFit();
+    AliTRDCalDet *calDetGain   = calibra->CreateDetObjectGain(&object);
+    TH1F *coefGain  = calDetGain->MakeHisto1DAsFunctionOfDet();
+    // Put them in the array
+    fCalibObjects->AddAt(calDetGain,kGain);
+    fPlots->AddAt(coefGain,kGain);
+    //
+    ok = kTRUE;
+  }
+  
+  calibra->ResetVectorFit();
+  
+  return ok;
+  
+}
+
+Bool_t AliTRDPreprocessorOffline::AnalyzeVdriftT0(){
+  //
+  // Analyze VdriftT0 - produce the calibration objects
+  //
+
+  AliTRDCalibraFit *calibra = AliTRDCalibraFit::Instance();
+  calibra->SetMinEntries(800*20); // If there is less than 1000 entries in the histo: no fit
+  calibra->AnalysePH(fPH2d);
+
+  Int_t nbtg = 6*4*18*((Int_t) ((AliTRDCalibraMode *)calibra->GetCalibraMode())->GetDetChamb0(1))
+    + 6*  18*((Int_t) ((AliTRDCalibraMode *)calibra->GetCalibraMode())->GetDetChamb2(1));
+  Int_t nbfit       = calibra->GetNumberFit();
+  Int_t nbfitSuccess = calibra->GetNumberFitSuccess();
+  Int_t nbE         = calibra->GetNumberEnt();
+
+  //printf("nbtg %d, nbfit %d, nbE %d, nbfitSuccess %d\n",nbtg,nbfit,nbE,nbfitSuccess);
+
+  Bool_t ok = kFALSE;
+  if ((nbtg >                  0) && 
+      (nbfit        >= 0.5*nbE) && (nbE > 30) && (nbfitSuccess > 30)) {
+    //printf("Pass the cut for VdriftT0\n");
+    // create the cal objects
+    calibra->RemoveOutliers(1,kTRUE);
+    calibra->PutMeanValueOtherVectorFit(1,kTRUE);
+    calibra->RemoveOutliers2(kTRUE);
+    calibra->PutMeanValueOtherVectorFit2(1,kTRUE);
+    TObjArray object = calibra->GetVectorFit();
+    AliTRDCalDet *calDetVdrift = calibra->CreateDetObjectVdrift(&object);
+    TH1F *coefVdriftPH  = calDetVdrift->MakeHisto1DAsFunctionOfDet();
+    AliTRDCalPad *calPadVdrift = (AliTRDCalPad *)calibra->CreatePadObjectVdrift(&object,calDetVdrift);
+    TH1F *coefPadVdrift   = calPadVdrift->MakeHisto1D();
+    object       = calibra->GetVectorFit2();
+    AliTRDCalDet *calDetT0  = calibra->CreateDetObjectT0(&object);
+    TH1F *coefT0  = calDetT0->MakeHisto1DAsFunctionOfDet();
+    AliTRDCalPad *calPadT0 = (AliTRDCalPad *)calibra->CreatePadObjectT0(&object,calDetT0);
+    TH1F *coefPadT0  = calPadT0->MakeHisto1D();
+    // Put them in the array
+    fCalibObjects->AddAt(calDetT0,kT0PHDet);
+    fCalibObjects->AddAt(calDetVdrift,kVdriftPHDet);
+    fCalibObjects->AddAt(calPadT0,kT0PHPad);
+    fCalibObjects->AddAt(calPadVdrift,kVdriftPHPad);
+    fPlots->AddAt(coefVdriftPH,kVdriftPHDet);
+    fPlots->AddAt(coefT0,kT0PHDet);
+    fPlots->AddAt(coefPadVdrift,kVdriftPHPad);
+    fPlots->AddAt(coefPadT0,kT0PHPad);
+    //
+    ok = kTRUE;
+  }
+  calibra->ResetVectorFit();
+  return ok;
+  
+}
+
+Bool_t AliTRDPreprocessorOffline::AnalyzeVdriftLinearFit(){
+  //
+  // Analyze vdrift linear fit - produce the calibration objects
+  //
+
+  AliTRDCalibraFit *calibra = AliTRDCalibraFit::Instance();
+  calibra->SetMinEntries(800); // If there is less than 1000 entries in the histo: no fit
+  calibra->AnalyseLinearFitters(fAliTRDCalibraVdriftLinearFit);
+
+  //Int_t nbtg        = 540;
+  Int_t nbfit       = calibra->GetNumberFit();
+  Int_t nbE         = calibra->GetNumberEnt();
+
+  
+  Bool_t ok = kFALSE;
+  // enough statistics
+  if ((nbfit        >= 0.5*nbE) && (nbE > 30)) {
+    // create the cal objects
+    TObjArray object  = calibra->GetVectorFit();
+    AliTRDCalDet *calDetVdrift = calibra->CreateDetObjectVdrift(&object,kTRUE);
+    TH1F *coefDriftLinear      = calDetVdrift->MakeHisto1DAsFunctionOfDet();
+    object                     = calibra->GetVectorFit2();
+    AliTRDCalDet *calDetLorentz = calibra->CreateDetObjectLorentzAngle(&object);
+    TH1F *coefLorentzAngle = calDetLorentz->MakeHisto1DAsFunctionOfDet();
+    // Put them in the array
+    fCalibObjects->AddAt(calDetVdrift,kVdriftLinear);
+    fCalibObjects->AddAt(calDetLorentz,kLorentzLinear);
+    fPlots->AddAt(coefDriftLinear,kVdriftLinear);
+    fPlots->AddAt(coefLorentzAngle,kLorentzLinear);
+    //
+    ok = kTRUE;
+  }
+  
+  calibra->ResetVectorFit();
+  
+  return ok;
+  
+}
+
+
+Bool_t AliTRDPreprocessorOffline::AnalyzePRF(){
+  //
+  // Analyze PRF - produce the calibration objects
+  //
+
+  AliTRDCalibraFit *calibra = AliTRDCalibraFit::Instance();
+  calibra->SetMinEntries(600); // If there is less than 1000 entries in the histo: no fit
+  calibra->AnalysePRFMarianFit(fPRF2d);
+
+  Int_t nbtg = 6*4*18*((Int_t) ((AliTRDCalibraMode *)calibra->GetCalibraMode())->GetDetChamb0(2))
+    + 6*  18*((Int_t) ((AliTRDCalibraMode *)calibra->GetCalibraMode())->GetDetChamb2(2));
+  Int_t nbfit       = calibra->GetNumberFit();
+  Int_t nbE         = calibra->GetNumberEnt();
+
+  
+  Bool_t ok = kFALSE;
+  // enough statistics
+  if ((nbtg >                  0) && 
+      (nbfit        >= 0.95*nbE) && (nbE > 30)) {
+    // create the cal objects
+    TObjArray object  = calibra->GetVectorFit();
+    AliTRDCalPad *calPadPRF = (AliTRDCalPad*) calibra->CreatePadObjectPRF(&object);
+    TH1F *coefPRF           = calPadPRF->MakeHisto1D();
+    // Put them in the array
+    fCalibObjects->AddAt(calPadPRF,kPRF);
+    fPlots->AddAt(coefPRF,kPRF);
+    //
+    ok = kTRUE;
+  }
+  
+  calibra->ResetVectorFit();
+  
+  return ok;
+  
+}
+
+
+void AliTRDPreprocessorOffline::UpdateOCDBGain(Int_t startRunNumber, Int_t endRunNumber, const Char_t *storagePath){
+  //
+  // Update OCDB entry
+  //
+  
+  AliCDBMetaData *metaData= new AliCDBMetaData();
+  metaData->SetObjectClassName("AliTRDCalDet");
+  metaData->SetResponsible("Raphaelle Bailhache");
+  metaData->SetBeamPeriod(1);
+  
+  AliCDBId id1("TRD/Calib/ChamberGainFactor", startRunNumber, endRunNumber);
+  AliCDBStorage * gStorage = AliCDBManager::Instance()->GetStorage(storagePath);
+  AliTRDCalDet *calDet = (AliTRDCalDet *) fCalibObjects->At(kGain);
+  if(calDet) gStorage->Put(calDet, id1, metaData);
+  //else {
+  //  printf("No calDet object for Gain\n");
+  //}
+
+    
+
+}
+
+void AliTRDPreprocessorOffline::UpdateOCDBVdrift(Int_t startRunNumber, Int_t endRunNumber, const Char_t *storagePath){
+  //
+  // Update OCDB entry
+  //
+
+  Int_t detVdrift = kVdriftPHDet;
+  
+  if(fMethodSecond) detVdrift = kVdriftLinear;
+  
+  AliCDBMetaData *metaData= new AliCDBMetaData();
+  metaData->SetObjectClassName("AliTRDCalDet");
+  metaData->SetResponsible("Raphaelle Bailhache");
+  metaData->SetBeamPeriod(1);
+  
+  AliCDBId id1("TRD/Calib/ChamberVdrift", startRunNumber, endRunNumber);
+  AliCDBStorage * gStorage = AliCDBManager::Instance()->GetStorage(storagePath);
+  AliTRDCalDet *calDet = (AliTRDCalDet *) fCalibObjects->At(detVdrift);
+  if(calDet) gStorage->Put(calDet, id1, metaData);
+  //else {
+  //  printf("No calDet object for Vdrift\n");
+  //}
+
+  //
+
+  if(!fMethodSecond) {
+    
+    AliCDBMetaData *metaDataPad= new AliCDBMetaData();
+    metaDataPad->SetObjectClassName("AliTRDCalPad");
+    metaDataPad->SetResponsible("Raphaelle Bailhache");
+    metaDataPad->SetBeamPeriod(1);
+    
+    AliCDBId id1Pad("TRD/Calib/LocalVdrift", startRunNumber, endRunNumber);
+    AliTRDCalPad *calPad = (AliTRDCalPad *) fCalibObjects->At(kVdriftPHPad);
+    if(calPad) gStorage->Put(calPad, id1Pad, metaDataPad);
+    //else {
+    //  printf("No calPad object for Vdrift\n");
+    //}
+
+  
+  }
+
+}
+
+void AliTRDPreprocessorOffline::UpdateOCDBT0(Int_t startRunNumber, Int_t endRunNumber, const Char_t *storagePath){
+  //
+  // Update OCDB entry
+  //
+  
+  AliCDBMetaData *metaData= new AliCDBMetaData();
+  metaData->SetObjectClassName("AliTRDCalDet");
+  metaData->SetResponsible("Raphaelle Bailhache");
+  metaData->SetBeamPeriod(1);
+  
+  AliCDBId id1("TRD/Calib/ChamberT0", startRunNumber, endRunNumber);
+  AliCDBStorage * gStorage = AliCDBManager::Instance()->GetStorage(storagePath);
+  AliTRDCalDet *calDet = (AliTRDCalDet *) fCalibObjects->At(kT0PHDet);
+  if(calDet) gStorage->Put(calDet, id1, metaData);
+  //else {
+  //  printf("No calDet object for T0\n");
+  //}
+
+  //
+
+  AliCDBMetaData *metaDataPad= new AliCDBMetaData();
+  metaDataPad->SetObjectClassName("AliTRDCalPad");
+  metaDataPad->SetResponsible("Raphaelle Bailhache");
+  metaDataPad->SetBeamPeriod(1);
+  
+  AliCDBId id1Pad("TRD/Calib/LocalT0", startRunNumber, endRunNumber);
+  AliTRDCalPad *calPad = (AliTRDCalPad *) fCalibObjects->At(kT0PHPad);
+  if(calPad) gStorage->Put(calPad, id1Pad, metaDataPad);
+  //else {
+  //  printf("No calPad object for T0\n");
+  //}
+
+
+}
+
+void AliTRDPreprocessorOffline::UpdateOCDBPRF(Int_t startRunNumber, Int_t endRunNumber, const Char_t *storagePath){
+  //
+  // Update OCDB entry
+  //
+  
+  AliCDBMetaData *metaData= new AliCDBMetaData();
+  metaData->SetObjectClassName("AliTRDCalPad");
+  metaData->SetResponsible("Raphaelle Bailhache");
+  metaData->SetBeamPeriod(1);
+  
+  AliCDBId id1("TRD/Calib/PRFWidth", startRunNumber, endRunNumber);
+  AliCDBStorage * gStorage = AliCDBManager::Instance()->GetStorage(storagePath);
+  AliTRDCalPad *calPad = (AliTRDCalPad *) fCalibObjects->At(kPRF);
+  if(calPad) gStorage->Put(calPad, id1, metaData);
+  //else {
+  //  printf("No calPad object for PRF\n");
+  //}
+
+}
diff --git a/TRD/AliTRDPreprocessorOffline.h b/TRD/AliTRDPreprocessorOffline.h
new file mode 100644 (file)
index 0000000..6f796e7
--- /dev/null
@@ -0,0 +1,81 @@
+#ifndef ALITRDPREPROCESSOROFFLINE_H
+#define ALITRDPREPROCESSOROFFLINE_H
+
+/* Copyright(c) 1998-1999, ALICE Experiment at CERN, All rights reserved. *
+ * See cxx source for full Copyright notice                               */
+
+//
+//
+//    Class to create OCDB entries - processing the results of the OFFLINE calibration
+//
+
+
+#include "TNamed.h"
+class TObjArray;
+class TH2I;
+class TProfile2D;
+class AliTRDCalibraVdriftLinearFit;
+class TH1I;
+class TH2F;
+
+
+class AliTRDPreprocessorOffline:public TNamed { 
+public:
+  enum{
+    kGain = 0,
+      kVdriftPHDet = 1,
+      kVdriftPHPad = 2,
+      kT0PHDet = 3,
+      kT0PHPad = 4,
+      kVdriftLinear = 5,
+      kLorentzLinear = 6,
+      kPRF = 7
+      };   
+  
+  AliTRDPreprocessorOffline();
+  virtual ~AliTRDPreprocessorOffline();
+
+  void SetLinearFitForVdrift(Bool_t methodsecond) { fMethodSecond = methodsecond;};
+  Bool_t GetLinearFitForVdrift() const { return fMethodSecond;};
+
+  void CalibVdriftT0(const Char_t* file, Int_t startRunNumber, Int_t endRunNumber, TString ocdbStorage="");
+  void CalibGain(const Char_t* file, Int_t startRunNumber, Int_t endRunNumber,  TString  ocdbStorage="");
+  void CalibPRF(const Char_t* file, Int_t startRunNumber, Int_t endRunNumber,  TString  ocdbStorage="");
+  
+  Bool_t ReadGainGlobal(const Char_t* fileName="CalibObjects.root");
+  Bool_t ReadVdriftT0Global(const Char_t* fileName="CalibObjects.root");
+  Bool_t ReadVdriftLinearFitGlobal(const Char_t* fileName="CalibObjects.root");
+  Bool_t ReadPRFGlobal(const Char_t* fileName="CalibObjects.root");
+
+  Bool_t AnalyzeGain(); 
+  Bool_t AnalyzeVdriftT0(); 
+  Bool_t AnalyzeVdriftLinearFit(); 
+  Bool_t AnalyzePRF(); 
+  
+  void UpdateOCDBT0(Int_t startRunNumber, Int_t endRunNumber, const char* storagePath);
+  void UpdateOCDBVdrift(Int_t startRunNumber, Int_t endRunNumber, const char* storagePath);
+  void UpdateOCDBGain(Int_t  startRunNumber, Int_t endRunNumber, const char* storagePath);
+  void UpdateOCDBPRF(Int_t  startRunNumber, Int_t endRunNumber, const char* storagePath);
+
+  
+private:
+  Bool_t fMethodSecond;                   // Second Method for drift velocity   
+  TH2I *fCH2d;                            // Gain
+  TProfile2D *fPH2d;                      // Drift velocity first method
+  TProfile2D *fPRF2d;                     // PRF
+  AliTRDCalibraVdriftLinearFit *fAliTRDCalibraVdriftLinearFit; // Drift velocity second method
+  TH1I *fNEvents;                         // Number of events 
+  TH2F *fAbsoluteGain;                    // Absolute Gain calibration
+  TObjArray * fPlots;                     // array with some plots to check
+  TObjArray * fCalibObjects;              // array with calibration objects 
+
+
+  
+
+private:
+  AliTRDPreprocessorOffline& operator=(const AliTRDPreprocessorOffline&); // not implemented
+  AliTRDPreprocessorOffline(const AliTRDPreprocessorOffline&); // not implemented
+  ClassDef(AliTRDPreprocessorOffline,1)
+};
+
+#endif
index 6e3514c..5acd431 100644 (file)
@@ -229,10 +229,10 @@ Double_t AliTRDCalROC::GetMean(AliTRDCalROC* const outlierROC) const
    Int_t nPoints = 0;
    for (Int_t i=0;i<fNchannels;i++) {
       if ((!outlierROC) || (!(outlierROC->GetValue(i)))) {
-       if(fData[i] > 0.000000000000001){
+       //if(fData[i] > 0.000000000000001){
          ddata[nPoints]= (Double_t) fData[i]/10000;
          nPoints++;
-       }
+        //}
       }
    }
    Double_t mean = TMath::Mean(nPoints,ddata);
@@ -273,10 +273,10 @@ Double_t AliTRDCalROC::GetRMS(AliTRDCalROC* const outlierROC) const
   Int_t nPoints = 0;
   for (Int_t i=0;i<fNchannels;i++) {
     if ((!outlierROC) || (!(outlierROC->GetValue(i)))) {
-       if(fData[i] > 0.000000000000001){
+      //if(fData[i] > 0.000000000000001){
          ddata[nPoints]= (Double_t)fData[i]/10000;
          nPoints++;
-       }
+        //}
     }
   }
   Double_t mean = TMath::RMS(nPoints,ddata);
index 6d32464..f95f3d3 100644 (file)
@@ -2,16 +2,15 @@
 /* Copyright(c) 1998-1999, ALICE Experiment at CERN, All rights reserved. *
  * See cxx source for full Copyright notice                               */
 
-/* $Id: TPCcalibLinkDef.h 34646 2009-09-05 13:12:53Z marian $ */
+/* $Id: TRDcalibLinkDef.h 34646 2009-09-05 13:12:53Z marian $ */
 
 #pragma link off all globals;
 #pragma link off all classes;
 #pragma link off all functions;
 
+#pragma link C++ class  AliTRDPreprocessorOffline+;
 #pragma link C++ class  AliTRDCalibTask+;
 
-
 #endif
 
 
index c0e7d22..27b5735 100644 (file)
@@ -1,6 +1,7 @@
 #-*- Mode: Makefile -*-
 
-SRCS = AliTRDCalibTask.cxx 
+SRCS = AliTRDCalibTask.cxx \
+       AliTRDPreprocessorOffline.cxx
 
 HDRS:= $(SRCS:.cxx=.h)