1 #include "AliPHOSDA1.h"
6 //----------------------------------------------------------------
7 AliPHOSDA1::AliPHOSDA1(int module) : TNamed(),
8 fHistoFile(0),fMod(module),fWriteToFile(kTRUE),fHistoArray()
11 // Create DA1 ("Calibration DA") object.
12 // module is the PHOS module number (0..4).
13 // Checks existence of histograms which might have been left
14 // from the previous runs to continue their filling.
15 // Histogram names: module_iX_iZ_gain for TH2 and module_iX_iZ for TH1.
16 // Root file name: PHOS_ModuleX_Calib.root, where X - module number.
19 sprintf(name,"PHOS_Module%d_Calib",fMod);
22 SetTitle("Calibration Detector Algorithm");
25 sprintf(rootname,"%s.root",GetName());
27 fHistoFile = new TFile(rootname,"update");
28 fHistoArray.SetName("histo_container");
34 for(Int_t iX=0; iX<64; iX++) {
35 for(Int_t iZ=0; iZ<56; iZ++) {
37 sprintf(hname,"%d_%d_%d",fMod,iX,iZ);
38 hist1 = (TH1F*)fHistoFile->Get(hname);
40 fHgLgRatio[iX][iZ] = hist1;
41 fHistoArray.Add(hist1);
44 fHgLgRatio[iX][iZ] = 0;
46 for(Int_t iGain=0; iGain<2; iGain++) {
47 sprintf(hname,"%d_%d_%d_%d",fMod,iX,iZ,iGain);
48 hist2 = (TH2F*)fHistoFile->Get(hname);
50 fTimeEnergy[iX][iZ][iGain] = hist2;
51 fHistoArray.Add(hist2);
54 fTimeEnergy[iX][iZ][iGain] = 0;
62 //-------------------------------------------------------------------
63 AliPHOSDA1::AliPHOSDA1(Int_t module, TH2F oldTimeEnergy[64][56][2]) :
64 TNamed(),fHistoFile(0),fMod(module),fWriteToFile(kFALSE),fHistoArray()
67 // oldTimeEnergy is an array of histograms kept from the previous run.
68 // By default the final histograms will not be saved to the root file.
71 sprintf(name,"PHOS_Module%d_Calib",fMod);
74 SetTitle("Calibration Detector Algorithm");
75 fHistoArray.SetName("histo_container");
79 for(Int_t iX=0; iX<64; iX++) {
80 for(Int_t iZ=0; iZ<56; iZ++) {
82 fHgLgRatio[iX][iZ] = 0;
83 for(Int_t iGain=0; iGain<2; iGain++) {
85 if((&(oldTimeEnergy[iX][iZ][iGain]))->GetEntries()) {
86 fTimeEnergy[iX][iZ][iGain] = &(oldTimeEnergy[iX][iZ][iGain]);
87 fHistoArray.Add(&(oldTimeEnergy[iX][iZ][iGain]));
90 fTimeEnergy[iX][iZ][iGain] = 0;
98 for(Int_t iX=0; iX<64; iX++) {
99 for(Int_t iZ=0; iZ<56; iZ++) {
100 fHgLgRatio[iX][iZ] = 0;
101 for(Int_t iGain=0; iGain<2; iGain++) {
102 fTimeEnergy[iX][iZ][iGain] = 0;
111 //-------------------------------------------------------------------
112 AliPHOSDA1::AliPHOSDA1(const AliPHOSDA1& da) : TNamed(da),
113 fHistoFile(0),fMod(da.fMod),fWriteToFile(da.fWriteToFile),fHistoArray(da.fHistoArray)
117 fHistoFile = new TFile(da.GetName(),"update");
123 for(Int_t iX=0; iX<64; iX++) {
124 for(Int_t iZ=0; iZ<56; iZ++) {
126 sprintf(hname,"%d_%d_%d",fMod,iX,iZ);
127 hist1 = (TH1F*)da.fHistoFile->Get(hname);
128 if(hist1) fHgLgRatio[iX][iZ] = hist1;
130 fHgLgRatio[iX][iZ] = 0;
132 for(Int_t iGain=0; iGain<2; iGain++) {
133 sprintf(hname,"%d_%d_%d_%d",fMod,iX,iZ,iGain);
134 hist2 = (TH2F*)da.fHistoFile->Get(hname);
135 if(hist2) fTimeEnergy[iX][iZ][iGain] = hist2;
137 fTimeEnergy[iX][iZ][iGain] = 0;
145 //-------------------------------------------------------------------
146 AliPHOSDA1& AliPHOSDA1::operator= (const AliPHOSDA1& da)
148 //Assignment operator.
152 TString oldname(fHistoFile->GetName());
153 TString newname(da.fHistoFile->GetName());
155 if(oldname != newname) {
157 fHistoFile = new TFile(da.fHistoFile->GetName(),"update");
162 SetName(da.GetName());
163 SetTitle(da.GetTitle());
165 for(Int_t iX=0; iX<64; iX++) {
166 for(Int_t iZ=0; iZ<56; iZ++) {
168 if(fHgLgRatio[iX][iZ]) delete fHgLgRatio[iX][iZ];
169 fHgLgRatio[iX][iZ] = da.fHgLgRatio[iX][iZ];
171 for(Int_t iGain=0; iGain<2; iGain++) {
172 if (fTimeEnergy[iX][iZ][iGain]) delete fTimeEnergy[iX][iZ][iGain];
173 fTimeEnergy[iX][iZ][iGain] = da.fTimeEnergy[iX][iZ][iGain];
179 fHistoArray = da.fHistoArray;
187 //-------------------------------------------------------------------
188 AliPHOSDA1::~AliPHOSDA1()
193 if(fHistoFile) delete fHistoFile;
197 //-------------------------------------------------------------------
198 void AliPHOSDA1::FillHistograms(Float_t e[64][56][2], Float_t t[64][56][2])
200 // Fill energy vs time-of-flight and HG/LG ratio histograms in one event.
201 // Energy and TOF data are encoded as e[X][Z][gain] and t[X][Z][gain],
202 // where X(0..63) and Z(0..55) - crystal position in the module,
203 // gain=0 - low gain, gain=1 - high gain.
204 // Energy in ADC counts, time in samples.
205 // If no energy or time read for particular channel,
206 // the correspondent array entry should be filled by zero.
207 // WARNING: this function should be called once per event!
212 for(Int_t iX=0; iX<64; iX++) {
213 for (Int_t iZ=0; iZ<56; iZ++) {
216 if(e[iX][iZ][0] && e[iX][iZ][1]) {
217 if(fHgLgRatio[iX][iZ]) {
218 // printf("iX=%d iZ=%d,e[iX][iZ][1]=%.3f,e[iX][iZ][0]=%.3f, t1=%.3f,t0=%.3f\n",
219 // iX,iZ,e[iX][iZ][1],e[iX][iZ][0], t[iX][iZ][1],t[iX][iZ][0]);
220 fHgLgRatio[iX][iZ]->Fill(e[iX][iZ][1]/e[iX][iZ][0]);
224 sprintf(hname,"%d_%d_%d",fMod,iX,iZ);
225 sprintf(htitl,"HG/LG ratio for crystal %d_%d_%d",fMod,iX,iZ);
226 fHgLgRatio[iX][iZ] = new TH1F(hname,htitl,400,14.,18.);
227 // printf("iX=%d iZ=%d,e[iX][iZ][1]=%.3f,e[iX][iZ][0]=%.3f\n",
228 // iX,iZ,e[iX][iZ][1],e[iX][iZ][0]);
229 fHgLgRatio[iX][iZ]->Fill(e[iX][iZ][1]/e[iX][iZ][0]);
230 fHistoArray.Add(fHgLgRatio[iX][iZ]);
235 for(Int_t iGain=0; iGain<2; iGain++) {
236 if(e[iX][iZ][iGain]<10) continue;
238 if(fTimeEnergy[iX][iZ][iGain])
239 fTimeEnergy[iX][iZ][iGain]->Fill(e[iX][iZ][iGain],t[iX][iZ][iGain]);
241 sprintf(hname,"%d_%d_%d_%d",fMod,iX,iZ,iGain);
242 sprintf(htitl,"Energy vs TOF for crystal %d_%d_%d and gain %d",fMod,iX,iZ,iGain);
243 fTimeEnergy[iX][iZ][iGain] = new TH2F(hname,htitl,1024,0.,1024.,100,0.,100);
244 fTimeEnergy[iX][iZ][iGain]->Fill(e[iX][iZ][iGain],t[iX][iZ][iGain]);
245 fHistoArray.Add(fTimeEnergy[iX][iZ][iGain]);
254 //-------------------------------------------------------------------
255 void AliPHOSDA1::UpdateHistoFile()
257 // Write histograms to file
259 if(!fWriteToFile) return;
260 if(!fHistoFile) return;
261 if(!fHistoFile->IsOpen()) return;
266 for(Int_t iX=0; iX<64; iX++) {
267 for(Int_t iZ=0; iZ<56; iZ++) {
269 hist1 = fHgLgRatio[iX][iZ];
270 if(hist1) hist1->Write(hist1->GetName(),TObject::kWriteDelete);
272 for(Int_t iGain=0; iGain<2; iGain++) {
273 hist2 = fTimeEnergy[iX][iZ][iGain];
274 if(hist2) hist2->Write(hist2->GetName(),TObject::kWriteDelete);
282 //-----------------------------------------------------------------------
283 void AliPHOSDA1::SetWriteToFile(Bool_t write)
286 fWriteToFile = write;
292 sprintf(rootname,"%s.root",GetName());
293 fHistoFile = new TFile(rootname,"update");
296 fWriteToFile = write;