adding some cuts for systematic study

[u/mrichter/AliRoot.git] / T0 / AliT0CalibSeasonTimeShift.cxx

diff --git a/T0/AliT0CalibSeasonTimeShift.cxx b/T0/AliT0CalibSeasonTimeShift.cxx
index e67bd7dea24716593ceaf948eb60e445ae294601..2d45e5bdf40bc918a728010e2ada80485b817c9c 100644 (file)
--- a/T0/AliT0CalibSeasonTimeShift.cxx
+++ b/T0/AliT0CalibSeasonTimeShift.cxx
@@ -36,7 +36,8 @@ ClassImp(AliT0CalibSeasonTimeShift)
   AliT0CalibSeasonTimeShift::AliT0CalibSeasonTimeShift():TNamed()
 {
   //
-
+  for (Int_t i=0; i<4; i++)
+    fMeanPar[i] = fSigmaPar[i] = 0; 
 }
 
 //________________________________________________________________
@@ -47,7 +48,11 @@ AliT0CalibSeasonTimeShift::AliT0CalibSeasonTimeShift(const char* name):TNamed()
   TString namst = "Calib_";
   namst += name;
   SetName(namst.Data());
-  SetTitle(namst.Data());
+  SetTitle(namst.Data()); 
+  
+  for (Int_t i=0; i<4; i++)
+    fMeanPar[i] = fSigmaPar[i] = 0; 
+ 
 }
 
 //________________________________________________________________
@@ -57,6 +62,7 @@ AliT0CalibSeasonTimeShift::AliT0CalibSeasonTimeShift(const AliT0CalibSeasonTimeS
 // copy constructor
   SetName(calibda.GetName());
   SetTitle(calibda.GetName());
+  ((AliT0CalibSeasonTimeShift &) calibda).Copy(*this);
 
 
 }
@@ -67,6 +73,7 @@ AliT0CalibSeasonTimeShift &AliT0CalibSeasonTimeShift::operator =(const AliT0Cali
 // assignment operator
   SetName(calibda.GetName());
   SetTitle(calibda.GetName());
+  if (this != &calibda) ((AliT0CalibSeasonTimeShift &) calibda).Copy(*this);
  
   return *this;
 }
@@ -91,47 +98,93 @@ void  AliT0CalibSeasonTimeShift::Print(Option_t*) const
 } 
 
 //________________________________________________________________
-void  AliT0CalibSeasonTimeShift::SetT0Par(Float_t par[4],Float_t spar[4])
+Bool_t  AliT0CalibSeasonTimeShift::SetT0Par(Float_t par[4],Float_t spar[4])
 {
-  for (Int_t i=0; i<4; i++)
+  Bool_t ok=false;
+ for (Int_t i=0; i<4; i++)
     {
       fMeanPar[i] = par[i];
       fSigmaPar[i] = spar[i];
+      if ( fSigmaPar[i] == 0 ||  fSigmaPar[i] > 500) ok = false;
     }
-
+ return ok;
 }
 
 //________________________________________________________________
-void AliT0CalibSeasonTimeShift::SetT0Par(const char* filePhys)
+Int_t AliT0CalibSeasonTimeShift::SetT0Par(const char* filePhys, Float_t *cdbtime)
 {
+  // compute shifts fo T0A, T0C, T0AC and resolution
+  //Status : writeok:
+  // 0 OK
+  // 2000 - no data file
+  // 300 no one histogram or it is empty
+  //-100 peak is very narrow
+
+  Float_t mean, sigma;
+  Int_t ok = 0;
+  TH1F *cfd = NULL;
+  TObjArray * tzeroObj = NULL;
 
   gFile = TFile::Open(filePhys);
   if(!gFile) {
     AliError("No input PHYS data found ");
+    return 2000;
   }
-  else
-    {
-      gFile->ls();
- 
-      TString histname[4]={"meanAC", "meanA", "meanC", "resolution"};
-       for (Int_t i=0; i<4; i++)
-       {
-         TH1F *cfd = (TH1F*) gFile->Get(histname[i].Data());
-         if(!cfd) AliWarning(Form("no histograms collected for %s", histname[i].Data()));
-         if(cfd) {
-           TF1 *g = new TF1("g", "gaus",-2,2);
-           cfd->Fit("g"," ","Q",-2,2);
-           Double_t par[3];
-           g->GetParameters(&par[0]);
-           fMeanPar[i] = par[1];
-           fSigmaPar[i]=par[2];
-
+  else {
+    //    gFile->ls();
+    //    TDirectory *dr = (TDirectory*) gFile->Get("T0Calib");
+    tzeroObj = dynamic_cast<TObjArray*>(gFile->Get("T0Calib"));
+    TString histname[4]={"fTzeroORAplusORC", "fTzeroORA", "fTzeroORC",  "fResolution"};
+    for (Int_t i=0; i<4; i++)
+      {
+       if(cfd) cfd->Reset();
+       if(tzeroObj) 
+         cfd = (TH1F*)tzeroObj->FindObject( histname[i].Data());
+       else
+         cfd =  (TH1F*)gFile ->Get(histname[i].Data());
+
+       if(!cfd) {
+         AliError(Form("no histograms collected for %s", histname[i].Data()));
+         return 300;
+       }
+       if(cfd) {
+         if( cfd->GetEntries() == 0) {
+         AliError(Form("%s histogram is empty", histname[i].Data()));
+         return 300;
          }
-       } 
+         GetMeanAndSigma(cfd, mean, sigma);
+         if (sigma == 0 || sigma > 600 || cfd->GetEntries()<50 ){ //!!!
+           AliError(Form("%s low statsitics or bad histogram, OCDB value is = %f", histname[i].Data(), cdbtime[i]) );
+           return 400;
+         }
+         if ( sigma > 0 && sigma < 600 && cfd->GetEntries()>=50) //!!!
+           { 
+             fMeanPar[i] =   mean;
+             fSigmaPar[i] = sigma;
+           }
+       }
+      } 
+  }
+  gFile->Close();
+  delete gFile;
+  return ok;
+}
+//________________________________________________________________________
+void AliT0CalibSeasonTimeShift::GetMeanAndSigma(TH1F* hist,  Float_t &mean, Float_t &sigma) {
 
-         gFile->Close();
-         delete gFile;
-         
-    }
-  
+  const double window =3.;  //fit window 
+ 
+  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;
 }