From 0ab3a530ed036a9f0cf64444d78f5f8da1df2f58 Mon Sep 17 00:00:00 2001 From: cvetan Date: Fri, 7 Dec 2007 14:06:31 +0000 Subject: [PATCH] New version of the PMD DA avoiding the creation of huge amount of histograms. Now the algorithm should be much faster (Basanta) --- PMD/AliPMDCalibGain.cxx | 93 +++++++++++++++++--------------- PMD/AliPMDCalibGain.h | 11 ++-- PMD/AliPMDCalibPedestal.cxx | 104 +++++++++++++++++++----------------- PMD/AliPMDCalibPedestal.h | 18 +++++-- PMD/PMDda.cxx | 16 +++++- 5 files changed, 139 insertions(+), 103 deletions(-) diff --git a/PMD/AliPMDCalibGain.cxx b/PMD/AliPMDCalibGain.cxx index 184298a01fa..6fc0b6d49ec 100644 --- a/PMD/AliPMDCalibGain.cxx +++ b/PMD/AliPMDCalibGain.cxx @@ -44,16 +44,19 @@ AliPMDCalibGain::AliPMDCalibGain(): TObject() { for(Int_t ismn = 0; ismn < kMaxSMN; ismn++) { - fHsmIso[idet][ismn] = new TH1F(Form("HmsIso_%d_%d",idet,ismn),"",100,0.,1000.); + fSMIso[idet][ismn] = 0.; + fSMCount[idet][ismn] = 0.; for(Int_t jrow = 0; jrow < kMaxRow; jrow++) { for(Int_t kcol = 0; kcol < kMaxCol; kcol++) { - fHadcIso[idet][ismn][jrow][kcol] = new TH1F(Form("HadcIso_%d_%d_%d_%d",idet,ismn,jrow,kcol),"",100,0.,1000.); + fCellIso[idet][ismn][jrow][kcol] = 0.; + fCellCount[idet][ismn][jrow][kcol] = 0.; } } } } + } // ------------------------------------------------------------------------ // AliPMDCalibGain::AliPMDCalibGain(const AliPMDCalibGain &pmdcalibgain): @@ -63,12 +66,14 @@ AliPMDCalibGain::AliPMDCalibGain(const AliPMDCalibGain &pmdcalibgain): { for(Int_t ismn = 0; ismn < kMaxSMN; ismn++) { - fHsmIso[idet][ismn] = pmdcalibgain.fHsmIso[idet][ismn] ; + fSMIso[idet][ismn] = pmdcalibgain.fSMIso[idet][ismn] ; + fSMCount[idet][ismn] = pmdcalibgain.fSMCount[idet][ismn] ; for(Int_t jrow = 0; jrow < kMaxRow; jrow++) { for(Int_t kcol = 0; kcol < kMaxCol; kcol++) { - fHadcIso[idet][ismn][jrow][kcol] = pmdcalibgain.fHadcIso[idet][ismn][jrow][kcol]; + fCellIso[idet][ismn][jrow][kcol] = pmdcalibgain.fCellIso[idet][ismn][jrow][kcol]; + fCellCount[idet][ismn][jrow][kcol] = pmdcalibgain.fCellCount[idet][ismn][jrow][kcol]; } } } @@ -84,13 +89,16 @@ AliPMDCalibGain &AliPMDCalibGain::operator=(const AliPMDCalibGain &pmdcalibgain) { for(Int_t ismn = 0; ismn < kMaxSMN; ismn++) { - fHsmIso[idet][ismn] = pmdcalibgain.fHsmIso[idet][ismn] ; + fSMIso[idet][ismn] = pmdcalibgain.fSMIso[idet][ismn]; + fSMCount[idet][ismn] = pmdcalibgain.fSMCount[idet][ismn]; for(Int_t jrow = 0; jrow < kMaxRow;jrow++) { for(Int_t kcol = 0; kcol < kMaxCol; kcol++) { - fHadcIso[idet][ismn][jrow][kcol] = - pmdcalibgain.fHadcIso[idet][ismn][jrow][kcol]; + fCellIso[idet][ismn][jrow][kcol] = + pmdcalibgain.fCellIso[idet][ismn][jrow][kcol]; + fCellCount[idet][ismn][jrow][kcol] = + pmdcalibgain.fCellCount[idet][ismn][jrow][kcol]; } } } @@ -102,20 +110,7 @@ AliPMDCalibGain &AliPMDCalibGain::operator=(const AliPMDCalibGain &pmdcalibgain) AliPMDCalibGain::~AliPMDCalibGain() { // dtor - for(Int_t idet = 0; idet < kDet; idet++) - { - for(Int_t ismn = 0; ismn < kMaxSMN; ismn++) - { - delete fHsmIso[idet][ismn]; - for(Int_t jrow = 0; jrow < kMaxRow; jrow++) - { - for(Int_t kcol = 0; kcol < kMaxCol; kcol++) - { - delete fHadcIso[idet][ismn][jrow][kcol]; - } - } - } - } + } // ------------------------------------------------------------------------ // Bool_t AliPMDCalibGain::ProcessEvent(AliRawReader *rawReader) @@ -133,8 +128,8 @@ Bool_t AliPMDCalibGain::ProcessEvent(AliRawReader *rawReader) Int_t id1,jd1; //neighbour row/col Int_t isocount; //number of neighbours with 0 signal - Int_t d1[kDet][kMaxSMN][kMaxRow][kMaxCol]; - Bool_t streamout = kFALSE; + Float_t d1[kDet][kMaxSMN][kMaxRow][kMaxCol]; + Bool_t streamout = kFALSE; for(Int_t idet = 0; idet < kDet; idet++) { @@ -144,7 +139,7 @@ Bool_t AliPMDCalibGain::ProcessEvent(AliRawReader *rawReader) { for(Int_t icol = 0; icol < kMaxCol; icol++) { - d1[idet][ismn][irow][icol] = 0; + d1[idet][ismn][irow][icol] = 0.; } } } @@ -172,7 +167,7 @@ Bool_t AliPMDCalibGain::ProcessEvent(AliRawReader *rawReader) if (isig>0) { - d1[idet][ismn][irow][icol] = isig; + d1[idet][ismn][irow][icol] = (Float_t) isig; } } pmdddlcont.Clear(); @@ -198,8 +193,13 @@ Bool_t AliPMDCalibGain::ProcessEvent(AliRawReader *rawReader) isocount++; if(isocount == kCellNeighbour) { - fHsmIso[idet][ismn]->Fill(d1[idet][ismn][irow][icol]); - fHadcIso[idet][ismn][irow][icol]->Fill(d1[idet][ismn][irow][icol]); + + fSMIso[idet][ismn] += d1[idet][ismn][irow][icol]; + fCellIso[idet][ismn][irow][icol] += d1[idet][ismn][irow][icol]; + + fSMCount[idet][ismn]++; + fCellCount[idet][ismn][irow][icol]++; + } } } // neigh cell cond. @@ -215,35 +215,42 @@ Bool_t AliPMDCalibGain::ProcessEvent(AliRawReader *rawReader) void AliPMDCalibGain::Analyse(TTree *gaintree) { // Calculates the mean - Int_t DET, SM, ROW, COL; - Float_t GAIN; - Float_t modmean, cellmean; + Int_t det, sm, row, col; + Float_t gain; + Float_t modmean = 0.; + Float_t cellmean = 0.; - gaintree->Branch("DET",&DET,"DET/I"); - gaintree->Branch("SM",&SM,"SM/I"); - gaintree->Branch("ROW",&ROW,"ROW/I"); - gaintree->Branch("COL",&COL,"COL/I"); - gaintree->Branch("GAIN",&GAIN,"GAIN/F"); + gaintree->Branch("det",&det,"det/I"); + gaintree->Branch("sm",&sm,"sm/I"); + gaintree->Branch("row",&row,"row/I"); + gaintree->Branch("col",&col,"col/I"); + gaintree->Branch("gain",&gain,"gain/F"); for(Int_t idet = 0; idet < kDet; idet++) { for(Int_t ism = 0; ism < kMaxSMN; ism++) { - modmean = fHsmIso[idet][ism]->GetMean(); + + if (fSMCount[idet][ism] > 0) + modmean = fSMIso[idet][ism]/fSMCount[idet][ism]; + for(Int_t irow = 0; irow < kMaxRow; irow++) { for(Int_t icol = 0; icol < kMaxCol; icol++) { - cellmean = fHadcIso[idet][ism][irow][icol]->GetMean(); - DET = idet; - SM = ism; - ROW = irow; - COL = icol; - GAIN = 0.; + if (fCellCount[idet][ism][irow][icol] > 0) + cellmean = fCellIso[idet][ism][irow][icol]/fCellCount[idet][ism][irow][icol]; + + + det = idet; + sm = ism; + row = irow; + col = icol; + gain = 1.; if(modmean > 0.0) { - GAIN = cellmean/modmean; + gain = cellmean/modmean; } gaintree->Fill(); } diff --git a/PMD/AliPMDCalibGain.h b/PMD/AliPMDCalibGain.h index cf7a291ceb6..a27c08c593f 100644 --- a/PMD/AliPMDCalibGain.h +++ b/PMD/AliPMDCalibGain.h @@ -25,16 +25,17 @@ class AliPMDCalibGain : public TObject enum { - kDet = 2, // Number of Planes + kDet = 2, // Number of Planes kMaxSMN = 24, // Number of Modules kMaxRow = 48, // Number of Rows kMaxCol = 96 // Number of Columns }; - TH1F *fHsmIso[kDet][kMaxSMN]; //histos of isolated cells modulewise - TH1F *fHadcIso[kDet][kMaxSMN][kMaxRow][kMaxCol]; // histos of iso cells cellwise + Float_t fSMIso[kDet][kMaxSMN]; + Float_t fSMCount[kDet][kMaxSMN]; // counter + Float_t fCellIso[kDet][kMaxSMN][kMaxRow][kMaxCol]; // adc of iso cells + Float_t fCellCount[kDet][kMaxSMN][kMaxRow][kMaxCol]; // counter - -ClassDef(AliPMDCalibGain,1) // description +ClassDef(AliPMDCalibGain,2) // description }; #endif // ALIPMDCALIBGAIN_H diff --git a/PMD/AliPMDCalibPedestal.cxx b/PMD/AliPMDCalibPedestal.cxx index 4d2e690bb62..e7d6239ee16 100644 --- a/PMD/AliPMDCalibPedestal.cxx +++ b/PMD/AliPMDCalibPedestal.cxx @@ -46,16 +46,17 @@ AliPMDCalibPedestal::AliPMDCalibPedestal() : // default constructor // - for (int i = 0; i < 2; i++) + for (int i = 0; i < kDet; i++) { - for (int j = 0; j < 24; j++) + for (int j = 0; j < kMaxSMN; j++) { - for (int k = 0; k < 96; k++) + for (int k = 0; k < kMaxRow; k++) { - for (int l = 0; l < 96; l++) + for (int l = 0; l < kMaxCol; l++) { - - fPedHisto[i][j][k][l] = new TH1F(Form("PedHisto_%d_%d_%d_%d",i,j,k,l),"",300,0.,300.); + fPedVal[i][j][k][l] = 0.; + fPedValSq[i][j][k][l] = 0.; + fPedCount[i][j][k][l] = 0.; } } } @@ -70,16 +71,17 @@ AliPMDCalibPedestal::AliPMDCalibPedestal(const AliPMDCalibPedestal &ped) : // // copy constructor // - for (int i = 0; i < 2; i++) + for (int i = 0; i < kDet; i++) { - for (int j = 0; j < 24; j++) + for (int j = 0; j < kMaxSMN; j++) { - for (int k = 0; k < 48; k++) + for (int k = 0; k < kMaxRow; k++) { - for (int l = 0; l < 96; l++) + for (int l = 0; l < kMaxCol; l++) { - - fPedHisto[i][j][k][l] = ped.fPedHisto[i][j][k][l]; + fPedVal[i][j][k][l] = ped.fPedVal[i][j][k][l]; + fPedValSq[i][j][k][l] = ped.fPedValSq[i][j][k][l]; + fPedCount[i][j][k][l] = ped.fPedCount[i][j][k][l]; } } } @@ -103,20 +105,6 @@ AliPMDCalibPedestal::~AliPMDCalibPedestal() // // destructor // - for (int i = 0; i < 2; i++) - { - for (int j = 0; j < 24; j++) - { - for (int k = 0; k < 96; k++) - { - for (int l = 0; l < 96; l++) - { - - delete fPedHisto[i][j][k][l]; - } - } - } - } } //_____________________________________________________________________ Bool_t AliPMDCalibPedestal::ProcessEvent(AliRawReader *rawReader) @@ -134,7 +122,6 @@ Bool_t AliPMDCalibPedestal::ProcessEvent(AliRawReader *rawReader) for (Int_t iddl = 0; iddl < kDDL; iddl++) { - rawReader->Select("PMD", iddl, iddl); //cout << reader.GetDataSize() << endl; streamout = rawStream.DdlData(iddl, &pmdddlcont); @@ -149,10 +136,11 @@ Bool_t AliPMDCalibPedestal::ProcessEvent(AliRawReader *rawReader) //Int_t chno = pmdddl->GetChannel(); Int_t row = pmdddl->GetRow(); Int_t col = pmdddl->GetColumn(); - Int_t sig = pmdddl->GetSignal(); + Float_t sig = (Float_t) pmdddl->GetSignal(); - fPedHisto[det][smn][row][col]->Fill((Float_t) sig); - + fPedVal[det][smn][row][col] += sig; + fPedValSq[det][smn][row][col] += sig*sig; + fPedCount[det][smn][row][col]++; } pmdddlcont.Clear(); } @@ -165,29 +153,47 @@ void AliPMDCalibPedestal::Analyse(TTree *pedtree) // // Calculate pedestal Mean and RMS // - Int_t DET, SM, ROW, COL; - Float_t MEAN, RMS; - pedtree->Branch("DET",&DET,"DET/I"); - pedtree->Branch("SM",&SM,"SM/I"); - pedtree->Branch("ROW",&ROW,"ROW/I"); - pedtree->Branch("COL",&COL,"COL/I"); - pedtree->Branch("MEAN",&MEAN,"MEAN/F"); - pedtree->Branch("RMS",&RMS,"RMS/F"); - - for (int idet = 0; idet < 2; idet++) + Int_t det, sm, row, col; + Float_t mean, rms; + Float_t meansq, diff; + + + pedtree->Branch("det",&det,"det/I"); + pedtree->Branch("sm",&sm,"sm/I"); + pedtree->Branch("row",&row,"row/I"); + pedtree->Branch("col",&col,"col/I"); + pedtree->Branch("mean",&mean,"mean/F"); + pedtree->Branch("rms",&rms,"rms/F"); + + for (int idet = 0; idet < kDet; idet++) { - for (int ism = 0; ism < 24; ism++) + for (int ism = 0; ism < kMaxSMN; ism++) { - for (int irow = 0; irow < 48; irow++) + for (int irow = 0; irow < kMaxRow; irow++) { - for (int icol = 0; icol < 96; icol++) + for (int icol = 0; icol < kMaxCol; icol++) { - DET = idet; - SM = ism; - ROW = irow; - COL = icol; - MEAN = fPedHisto[idet][ism][irow][icol]->GetMean(); - RMS = fPedHisto[idet][ism][irow][icol]->GetRMS(); + det = idet; + sm = ism; + row = irow; + col = icol; + if (fPedCount[idet][ism][irow][icol] > 0) + { + mean = fPedVal[idet][ism][irow][icol]/fPedCount[idet][ism][irow][icol]; + + meansq = fPedValSq[idet][ism][irow][icol]/fPedCount[idet][ism][irow][icol]; + + diff = meansq - mean*mean; + if (diff > 0.) + { + rms = sqrt(diff); + } + else + { + rms = 0.; + } + } + pedtree->Fill(); } } diff --git a/PMD/AliPMDCalibPedestal.h b/PMD/AliPMDCalibPedestal.h index 172b144e0e7..861c557dad0 100644 --- a/PMD/AliPMDCalibPedestal.h +++ b/PMD/AliPMDCalibPedestal.h @@ -22,10 +22,20 @@ public: private: - TH1F *fPedHisto[2][24][48][96]; - - - ClassDef(AliPMDCalibPedestal,2) + enum + { + kDet = 2, // Number of Planes + kMaxSMN = 24, // Number of Modules + kMaxRow = 48, // Number of Rows + kMaxCol = 96 // Number of Columns + }; + + Float_t fPedVal[kDet][kMaxSMN][kMaxRow][kMaxCol]; + Float_t fPedValSq[kDet][kMaxSMN][kMaxRow][kMaxCol]; + Float_t fPedCount[kDet][kMaxSMN][kMaxRow][kMaxCol]; + + + ClassDef(AliPMDCalibPedestal,3) }; diff --git a/PMD/PMDda.cxx b/PMD/PMDda.cxx index 970880b756b..f1c488c1dc0 100644 --- a/PMD/PMDda.cxx +++ b/PMD/PMDda.cxx @@ -36,6 +36,10 @@ extern "C" { #include "TH1F.h" #include "TBenchmark.h" #include "TTree.h" +#include "TROOT.h" +#include "TPluginManager.h" + + /* Main routine Arguments: @@ -43,8 +47,16 @@ extern "C" { */ int main(int argc, char **argv) { - AliPMDCalibPedestal calibped; - AliPMDCalibGain calibgain; + /* magic line from Rene */ + gROOT->GetPluginManager()->AddHandler("TVirtualStreamerInfo", + "*", + "TStreamerInfo", + "RIO", + "TStreamerInfo()"); + + + AliPMDCalibPedestal calibped; + AliPMDCalibGain calibgain; TTree *ped = new TTree("ped","PMD Pedestal tree"); TTree *gain = new TTree("gain","PMD Gain tree"); -- 2.31.1