From 1e7108ce91875f8261cd50a506dc52bf1d16c40f Mon Sep 17 00:00:00 2001 From: alla Date: Sun, 27 Feb 2011 22:11:03 +0000 Subject: [PATCH] 1st draft for pass0 --- T0/AliT0CalibTimeEq.cxx | 160 +++++++++++++++++++++++++++++++++------- T0/AliT0CalibTimeEq.h | 5 ++ T0/T0Physda.cxx | 10 ++- 3 files changed, 145 insertions(+), 30 deletions(-) diff --git a/T0/AliT0CalibTimeEq.cxx b/T0/AliT0CalibTimeEq.cxx index bbc77285e75..f415d0d3cfa 100644 --- a/T0/AliT0CalibTimeEq.cxx +++ b/T0/AliT0CalibTimeEq.cxx @@ -109,15 +109,10 @@ void AliT0CalibTimeEq::Print(Option_t*) const Bool_t AliT0CalibTimeEq::ComputeOnlineParams(const char* filePhys) { // compute online equalized time - Double_t mean=0, meanver=0; - Double_t rms=0, rmsver=0; + Float_t meandiff, sigmadiff, meanver, meancfdtime, sigmacfdtime; + Double_t rms=0, rmsver=0, rmscfd=0; Int_t nent=0; Bool_t ok=false; - Int_t npeaks = 20; - Int_t sigma=3; - Bool_t down=false; - Int_t index[20]; - Int_t nfound=0; gFile = TFile::Open(filePhys); if(!gFile) { AliError("No input PHYS data found "); @@ -126,43 +121,131 @@ Bool_t AliT0CalibTimeEq::ComputeOnlineParams(const char* filePhys) { // gFile->ls(); ok=true; - Char_t buf1[30]; for (Int_t i=0; i<24; i++) { - sprintf(buf1,"CFD1minCFD%d",i+1); - TH1F *cfd = (TH1F*) gFile->Get(buf1); + TH1F *cfd = (TH1F*) gFile->Get(Form("CFD1minCFD%d",i+1)); + TH1F *cfdtime = (TH1F*) gFile->Get(Form("CFD%d",i+1)); if(!cfd) AliWarning(Form("no histograms collected by PHYS DA for channel %i", i)); // printf(" i = %d buf1 = %s\n", i, buf1); if(cfd) { - TSpectrum *s = new TSpectrum(2*npeaks,1); - nfound = s->Search(cfd,sigma," ",0.1); - if(nfound!=0){ - Float_t *xpeak = s->GetPositionX(); - TMath::Sort(nfound, xpeak, index,down); - Float_t xp = xpeak[index[0]]; - Double_t hmax = xp+3*sigma; - Double_t hmin = xp-3*sigma; - cfd->GetXaxis()->SetRangeUser(hmin-20,hmax+20); - mean=cfd->GetMean(); - rms=cfd->GetRMS(); - nent=(Int_t)cfd->GetEntries(); - if(nent<500 || rms>20. ) { + GetMeanAndSigma(cfd, meandiff, sigmadiff); + nent=cfd->GetEntries(); + if(nent<500 || cfd->GetRMS()>20. ) { + ok=false; + AliWarning(Form("Data is not good enouph in PMT %i - mean %f rsm %f nentries %i", i,meandiff,sigmadiff , nent)); + + } + else + { ok=false; - AliWarning(Form("Data is not good enouph in PMT %i - mean %f rsm %f nentries %i", i,mean, rms, nent)); + AliWarning(Form("Data is not good enouph in PMT %i , no clean peak", i)); + } + if(!cfd) AliWarning(Form("no histograms collected by PHYS DA for channel %i", i)); + } + // printf(" i = %d buf1 = %s\n", i, buf1); + if(cfdtime) { + GetMeanAndSigma(cfdtime,meancfdtime, sigmacfdtime); + nent=cfdtime->GetEntries(); + if(nent<500 || sigmacfdtime>30. ) { + ok=false; + AliWarning(Form("Data is not good enouph in PMT %i CFD data - meancfdtime %f rsm %f nentries %i", i,meancfdtime, sigmacfdtime, nent)); + } + } + else + { + ok=false; + AliWarning(Form("Data is not good enouph in PMT %i , no clean peak", i)); } + + + + SetTimeEq(i,meandiff); + SetTimeEqRms(i,sigmadiff); + SetCFDvalue(i,0,meancfdtime); + SetCFDvalue(i,0,sigmacfdtime); + if (cfd) delete cfd; + if (cfdtime) delete cfdtime; + + } + TH1F *ver = (TH1F*) gFile->Get("hVertex"); + if(!ver) AliWarning("no T0 histogram collected by PHYS DA "); + if(ver) { + meanver = ver->GetMean(); + rmsver = ver->GetRMS(); + } + SetMeanVertex(meanver); + SetRmsVertex(rmsver); + + gFile->Close(); + delete gFile; + + } + return ok; +} + +//________________________________________________________________ +Bool_t AliT0CalibTimeEq::ComputeOfflineParams(const char* filePhys) +{ + // compute online equalized time + Float_t meandiff, sigmadiff, meanver, meancfdtime, sigmacfdtime; + Double_t rms=0, rmsver=0, rmscfd=0; + Int_t nent=0; + Bool_t ok=false; + gFile = TFile::Open(filePhys); + if(!gFile) { + AliError("No input PHYS data found "); + } + else + { + // gFile->ls(); + ok=true; + TObjArray * TzeroObj = (TObjArray*) gFile->Get("fTzeroObject"); + for (Int_t i=0; i<24; i++) + { + TH1F *cfddiff = (TH1F*)TzeroObj->At(i); + TH1F *cfdtime = (TH1F*)TzeroObj->At(i+24); + if(!cfddiff) AliWarning(Form("no histograms collected by PHYS DA for channel %i", i)); + // printf(" i = %d buf1 = %s\n", i, buf1); + if(cfddiff) { + GetMeanAndSigma(cfddiff,meandiff, sigmadiff); + nent=cfddiff->GetEntries(); + if(nent<500 || cfddiff->GetRMS()>20. ) { + ok=false; + AliWarning(Form("Data is not good enouph in PMT %i - mean %f rsm %f nentries %i", i,meandiff,sigmadiff, nent)); + } else { ok=false; AliWarning(Form("Data is not good enouph in PMT %i , no clean peak", i)); } + if(!cfdtime) AliWarning(Form("no histograms collected by PHYS DA for channel %i", i)); + // printf(" i = %d buf1 = %s\n", i, buf1); + if(cfdtime) { + GetMeanAndSigma(cfdtime,meancfdtime, sigmacfdtime); + nent=cfdtime->GetEntries(); + if(nent<500 || cfdtime->GetRMS()>30. ) { + ok=false; + AliWarning(Form("Data is not good enouph in PMT %i CFD data - mean %f rsm %f nentries %i", i,meancfdtime, sigmacfdtime, nent)); + + } + } + else + { + ok=false; + AliWarning(Form("Data is not good enouph in PMT %i , no clean peak", i)); + } } - SetTimeEq(i,mean); - SetTimeEqRms(i,rms); - if (cfd) delete cfd; + SetTimeEq(i,meandiff); + SetTimeEqRms(i,sigmadiff); + SetCFDvalue(i,0,meancfdtime); + SetCFDvalue(i,0,sigmacfdtime); + if (cfddiff) delete cfddiff; + if (cfdtime) delete cfdtime; + } TH1F *ver = (TH1F*) gFile->Get("hVertex"); if(!ver) AliWarning("no T0 histogram collected by PHYS DA "); @@ -180,4 +263,27 @@ Bool_t AliT0CalibTimeEq::ComputeOnlineParams(const char* filePhys) return ok; } +//________________________________________________________________________ +void AliT0CalibTimeEq::GetMeanAndSigma(TH1F* hist, Float_t &mean, Float_t &sigma) { + + const double threshold = 0.5; // threshold level with respect to maximum bin + const double window = 5.; //fit window + double norm = hist->Integral(); // normalize to one count + hist->Scale(1./norm); + + double meanEstimate, sigmaEstimate; + int maxBin; + maxBin = hist->GetMaximumBin(); //position of maximum + meanEstimate = hist->GetBinCenter( maxBin); // mean of gaussian sitting in maximum + sigmaEstimate = hist->GetRMS(); + TF1* fit= new TF1("fit","gaus", meanEstimate - window*sigmaEstimate, meanEstimate + window*sigmaEstimate); + fit->SetParameters(hist->GetBinContent(maxBin), meanEstimate, sigmaEstimate); + hist->Fit("fit","R"); + + mean = (Float_t) fit->GetParameter(1); + sigma = (Float_t) fit->GetParameter(2); + + delete fit; +} + diff --git a/T0/AliT0CalibTimeEq.h b/T0/AliT0CalibTimeEq.h index 845abbdf5be..12d4bc1fd0b 100644 --- a/T0/AliT0CalibTimeEq.h +++ b/T0/AliT0CalibTimeEq.h @@ -9,6 +9,7 @@ //////////////////////////////////////////////// #include "TNamed.h" +#include "TH1F.h" class AliT0CalibTimeEq: public TNamed { @@ -23,6 +24,8 @@ class AliT0CalibTimeEq: public TNamed { virtual void Print(Option_t* option= "") const; Bool_t ComputeOnlineParams(const char* filePhys); + Bool_t ComputeOfflineParams(const char* filePhys); + Float_t GetCFDvalue(Int_t channel,Int_t number) const {return fCFDvalue[channel][number];} Float_t* GetCFDvalue() const {return (float*)fCFDvalue;} Float_t GetTimeEq(Int_t channel) const {return fTimeEq[channel];} @@ -41,8 +44,10 @@ class AliT0CalibTimeEq: public TNamed { void SetRmsVertex(Float_t rms=0) { fRmsVertex = rms; }; Float_t GetRmsVertex () {return fRmsVertex;}; + protected: + void GetMeanAndSigma(TH1F* hist, Float_t &mean, Float_t &sigma); Float_t fCFDvalue[24][5]; // CFD values Float_t fTimeEq[24]; // Time Equalized for OCDB Float_t fTimeEqRms[24]; // RMS of Time Equalized for OCDB diff --git a/T0/T0Physda.cxx b/T0/T0Physda.cxx index 4e6921fa71f..b643471601d 100644 --- a/T0/T0Physda.cxx +++ b/T0/T0Physda.cxx @@ -125,6 +125,7 @@ int main(int argc, char **argv) { for(Int_t ic=0; ic<24; ic++) { hCFD1minCFD[ic] = new TH1F(Form("CFD1minCFD%d",ic+1),"CFD-CFD",kcbx,kclx,kcmx); + hCFD[ic] = new TH1F(Form("CFD1minCFD%d",ic+1),"CFD",kt0bx,kt0lx,kt0hx); } TH1F *hVertex = new TH1F("hVertex","T0 time",kt0bx,kt0lx,kt0hx); @@ -183,13 +184,13 @@ int main(int argc, char **argv) { AliT0RawReader *start = new AliT0RawReader(reader, kTRUE); // Read raw data - Int_t allData[105][5]; - for(Int_t i0=0;i0<105;i0++) + Int_t allData[110][5]; + for(Int_t i0=0;i0<107;i0++) for(Int_t j0=0;j0<5;j0++) allData[i0][j0] = 0; if(start->Next()){ - for (Int_t i=0; i<105; i++) { + for (Int_t i=0; i<107; i++) { for(Int_t iHit=0;iHit<5;iHit++){ allData[i][iHit]= start->GetData(i,iHit); } @@ -203,6 +204,7 @@ int main(int argc, char **argv) { Float_t meanShift[24]; for (Int_t ik = 0; ik<24; ik++) { + if(allData[ik+1][0]>0 ) hCFD[ik]->Fill(allData[ik+1][0]); if(ik<12 && allData[ik+1][0]>0 && allData[knpmtC][0]>0 ){ hCFD1minCFD[ik]->Fill(allData[ik+1][0]-allData[knpmtC][0]); } @@ -268,6 +270,8 @@ int main(int argc, char **argv) { for(Int_t j=0;j<24;j++){ hCFD1minCFD[j]->SetDirectory(&hist); hCFD1minCFD[j]->Write(); + hCFD[j]->Write(); + } hVertex->Write(); hist.Close(); -- 2.43.0