]> git.uio.no Git - u/mrichter/AliRoot.git/blame - T0/AliT0CalibSeasonTimeShift.cxx
Correction (Raphael)
[u/mrichter/AliRoot.git] / T0 / AliT0CalibSeasonTimeShift.cxx
CommitLineData
455957bc 1/**************************************************************************
2 * Copyright(c) 1998-1999, ALICE Experiment at CERN, All rights reserved. *
3 * *
4 * Author: The ALICE Off-line Project. *
5 * Contributors are mentioned in the code where appropriate. *
6 * *
7 * Permission to use, copy, modify and distribute this software and its *
8 * documentation strictly for non-commercial purposes is hereby granted *
9 * without fee, provided that the above copyright notice appears in all *
10 * copies and that both the copyright notice and this permission notice *
11 * appear in the supporting documentation. The authors make no claims *
12 * about the suitability of this software for any purpose. It is *
13 * provided "as is" without express or implied warranty. *
14 **************************************************************************/
15
16/* $Id: AliT0CalibSeasonTimeShift.cxx 42881 2010-08-16 10:59:14Z alla $ */
17
18///////////////////////////////////////////////////////////////////////////////
19// //
20// class for T0 calibration TM-AC-AM_6-02-2006
21// equalize time shift for each time CFD channel
22// //
23///////////////////////////////////////////////////////////////////////////////
24
25#include "AliT0CalibSeasonTimeShift.h"
26#include "AliLog.h"
27#include <TFile.h>
28#include <TMath.h>
29#include <TF1.h>
30#include <TProfile.h>
31#include <iostream>
32
33ClassImp(AliT0CalibSeasonTimeShift)
34
35//________________________________________________________________
36 AliT0CalibSeasonTimeShift::AliT0CalibSeasonTimeShift():TNamed()
37{
38 //
57b90263 39 for (Int_t i=0; i<4; i++)
40 fMeanPar[i] = fSigmaPar[i] = 0;
455957bc 41}
42
43//________________________________________________________________
44AliT0CalibSeasonTimeShift::AliT0CalibSeasonTimeShift(const char* name):TNamed()
45{
46 //constructor
47
48 TString namst = "Calib_";
49 namst += name;
50 SetName(namst.Data());
57b90263 51 SetTitle(namst.Data());
52
53 for (Int_t i=0; i<4; i++)
54 fMeanPar[i] = fSigmaPar[i] = 0;
55
455957bc 56}
57
58//________________________________________________________________
59AliT0CalibSeasonTimeShift::AliT0CalibSeasonTimeShift(const AliT0CalibSeasonTimeShift& calibda):TNamed(calibda)
60
61{
62// copy constructor
63 SetName(calibda.GetName());
64 SetTitle(calibda.GetName());
57b90263 65 ((AliT0CalibSeasonTimeShift &) calibda).Copy(*this);
455957bc 66
67
68}
69
70//________________________________________________________________
71AliT0CalibSeasonTimeShift &AliT0CalibSeasonTimeShift::operator =(const AliT0CalibSeasonTimeShift& calibda)
72{
73// assignment operator
74 SetName(calibda.GetName());
75 SetTitle(calibda.GetName());
57b90263 76 if (this != &calibda) ((AliT0CalibSeasonTimeShift &) calibda).Copy(*this);
455957bc 77
78 return *this;
79}
80
81//________________________________________________________________
82AliT0CalibSeasonTimeShift::~AliT0CalibSeasonTimeShift()
83{
84 //
85 // destrictor
86}
87
88
89//________________________________________________________________
90void AliT0CalibSeasonTimeShift::Print(Option_t*) const
91{
92 // print time values
93
94 printf("\n ---- T0 results ----\n\n");
95 printf(" (T0A+T0C)/2 = %f; T0A = %f; T0C = %f; resolution = %f \n", fMeanPar[0], fMeanPar[1],fMeanPar[2],fMeanPar[3]);
96 printf(" sigma(T0A+T0C)/2 = %f; sigma(T0 = %f; sigma(T0C) = %f; sigma(resolution) = %f \n" , fSigmaPar[0], fSigmaPar[1], fSigmaPar[2],fSigmaPar[3]);
97
98}
99
100//________________________________________________________________
cd8ca88e 101Bool_t AliT0CalibSeasonTimeShift::SetT0Par(Float_t par[4],Float_t spar[4])
455957bc 102{
cd8ca88e 103 Bool_t ok=false;
104 for (Int_t i=0; i<4; i++)
455957bc 105 {
106 fMeanPar[i] = par[i];
107 fSigmaPar[i] = spar[i];
cd8ca88e 108 if ( fSigmaPar[i] == 0 || fSigmaPar[i] > 500) ok = false;
455957bc 109 }
cd8ca88e 110 return ok;
455957bc 111}
112
113//________________________________________________________________
32d3e409 114Int_t AliT0CalibSeasonTimeShift::SetT0Par(const char* filePhys, Float_t *cdbtime)
455957bc 115{
cd8ca88e 116 // compute online equalized time
117 Float_t mean, sigma;
32d3e409 118 Int_t ok = 0;
bb54f817 119 TH1F *cfd = NULL;
61bd03f6 120 TObjArray * tzeroObj = NULL;
dd92f5b8 121
455957bc 122 gFile = TFile::Open(filePhys);
123 if(!gFile) {
124 AliError("No input PHYS data found ");
32d3e409 125 ok = 2000;
455957bc 126 }
61bd03f6 127 else {
128 // gFile->ls();
129 // TDirectory *dr = (TDirectory*) gFile->Get("T0Calib");
130 tzeroObj = dynamic_cast<TObjArray*>(gFile->Get("T0Calib"));
cd8ca88e 131 TString histname[4]={"fTzeroORAplusORC", "fTzeroORA", "fTzeroORC", "fResolution"};
132 for (Int_t i=0; i<4; i++)
b0ab3f59 133 {
bb54f817 134 if(cfd) cfd->Reset();
61bd03f6 135 if(tzeroObj)
136 cfd = (TH1F*)tzeroObj->FindObject( histname[i].Data());
bb54f817 137 else
dd92f5b8 138 cfd = (TH1F*)gFile ->Get(histname[i].Data());
139
140 if(!cfd) {
32d3e409 141 ok=300;
142 AliError(Form("no histograms collected for %s", histname[i].Data()));
b9dd81a7 143 return ok;
dd92f5b8 144 }
bb54f817 145 if(cfd) {
32d3e409 146 if( cfd->GetEntries() == 0) {
147 ok=300;
148 AliError(Form("%s histogram is empty", histname[i].Data()));
149 return ok;
150 }
cd8ca88e 151 GetMeanAndSigma(cfd, mean, sigma);
84fe4d28 152 if (sigma == 0 || sigma > 500 || cfd->GetEntries()<200 ){
32d3e409 153 ok = -100;
18086258 154 fMeanPar[i] = cdbtime[i];
84fe4d28 155 fSigmaPar[i] = -1;
18086258 156 }
84fe4d28 157 if ( sigma > 0 && sigma < 500 && cfd->GetEntries()>200)
cd8ca88e 158 {
b9dd81a7 159 fMeanPar[i] = mean;
cd8ca88e 160 fSigmaPar[i] = sigma;
cd8ca88e 161 }
455957bc 162 }
b0ab3f59 163 }
61bd03f6 164 }
165 gFile->Close();
166 delete gFile;
cd8ca88e 167 return ok;
168}
169//________________________________________________________________________
170void AliT0CalibSeasonTimeShift::GetMeanAndSigma(TH1F* hist, Float_t &mean, Float_t &sigma) {
171
18086258 172 const double window =2.; //fit window
cd8ca88e 173
174 double meanEstimate, sigmaEstimate;
175 int maxBin;
176 maxBin = hist->GetMaximumBin(); //position of maximum
177 meanEstimate = hist->GetBinCenter( maxBin); // mean of gaussian sitting in maximum
178 sigmaEstimate = hist->GetRMS();
179 TF1* fit= new TF1("fit","gaus", meanEstimate - window*sigmaEstimate, meanEstimate + window*sigmaEstimate);
180 fit->SetParameters(hist->GetBinContent(maxBin), meanEstimate, sigmaEstimate);
84fe4d28 181 hist->Fit("fit","RQ","Q");
cd8ca88e 182
183 mean = (Float_t) fit->GetParameter(1);
184 sigma = (Float_t) fit->GetParameter(2);
185
186 delete fit;
455957bc 187}