1 #include "AliPHOSRcuDA1.h"
4 ClassImp(AliPHOSRcuDA1)
6 //----------------------------------------------------------------
7 AliPHOSRcuDA1::AliPHOSRcuDA1(Int_t module, Int_t rcu) : TNamed(),
8 fHistoFile(0),fMod(module),fRCU(rcu),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_RCUY_Calib.root, where X - module number,
17 // Y - RCU number. If no RCU specified (rcu<0), file name is simply
18 // PHOS_ModuleX_Calib.root.
21 if(rcu<0) sprintf(name,"PHOS_Module%d_Calib",fMod);
22 else sprintf(name,"PHOS_Module%d_RCU%d_Calib",fMod,fRCU);
25 SetTitle("Calibration Detector Algorithm");
28 sprintf(rootname,"%s.root",GetName());
30 fHistoFile = new TFile(rootname,"update");
31 fHistoArray.SetName("histo_container");
37 for(Int_t iX=0; iX<64; iX++) {
38 for(Int_t iZ=0; iZ<56; iZ++) {
40 sprintf(hname,"%d_%d_%d",fMod,iX,iZ);
41 hist1 = (TH1F*)fHistoFile->Get(hname);
43 fHgLgRatio[iX][iZ] = hist1;
44 fHistoArray.Add(hist1);
47 fHgLgRatio[iX][iZ] = 0;
49 for(Int_t iGain=0; iGain<2; iGain++) {
50 sprintf(hname,"%d_%d_%d_%d",fMod,iX,iZ,iGain);
51 hist2 = (TH2F*)fHistoFile->Get(hname);
53 fTimeEnergy[iX][iZ][iGain] = hist2;
54 fHistoArray.Add(hist2);
57 fTimeEnergy[iX][iZ][iGain] = 0;
65 //-------------------------------------------------------------------
66 AliPHOSRcuDA1::AliPHOSRcuDA1(Int_t module, Int_t rcu, TObjArray* oldTimeEnergy) :
67 TNamed(),fHistoFile(0),fMod(module),fRCU(rcu),fWriteToFile(kFALSE),
71 // oldTimeEnergy is an array of histograms kept from the previous run.
72 // By default the final histograms will not be saved to the root file.
75 if(rcu<0) sprintf(name,"PHOS_Module%d_Calib",fMod);
76 else sprintf(name,"PHOS_Module%d_RCU%d_Calib",fMod,fRCU);
79 SetTitle("Calibration Detector Algorithm");
80 fHistoArray.SetName("histo_container");
86 for(Int_t iX=0; iX<64; iX++) {
87 for(Int_t iZ=0; iZ<56; iZ++) {
89 sprintf(hname,"%d_%d_%d",fMod,iX,iZ);
91 hist1 = (TH1F*)oldTimeEnergy->FindObject(hname);
93 fHgLgRatio[iX][iZ] = hist1;
94 fHistoArray.Add(hist1);
97 fHgLgRatio[iX][iZ] = 0;
99 for(Int_t iGain=0; iGain<2; iGain++) {
100 sprintf(hname,"%d_%d_%d_%d",fMod,iX,iZ,iGain);
102 hist2 = (TH2F*)oldTimeEnergy->FindObject(hname);
104 fTimeEnergy[iX][iZ][iGain] = hist2;
105 fHistoArray.Add(hist2);
108 fTimeEnergy[iX][iZ][iGain] = 0;
116 //-------------------------------------------------------------------
117 AliPHOSRcuDA1::~AliPHOSRcuDA1()
122 if(fHistoFile) delete fHistoFile;
126 //-------------------------------------------------------------------
127 void AliPHOSRcuDA1::FillHistograms(Float_t e[64][56][2], Float_t t[64][56][2])
129 // Fill energy vs time-of-flight and HG/LG ratio histograms in one event.
130 // Energy and TOF data are encoded as e[X][Z][gain] and t[X][Z][gain],
131 // where X(0..63) and Z(0..55) - crystal position in the module,
132 // gain=0 - low gain, gain=1 - high gain.
133 // Energy in ADC counts, time in samples.
134 // If no energy or time read for particular channel,
135 // the correspondent array entry should be filled by zero.
136 // WARNING: this function should be called once per event!
141 for(Int_t iX=0; iX<64; iX++) {
142 for (Int_t iZ=0; iZ<56; iZ++) {
145 if(e[iX][iZ][0]>1 && e[iX][iZ][1]>1 ) {
147 if(fHgLgRatio[iX][iZ]) {
148 //printf("iX=%d iZ=%d,e[iX][iZ][1]=%.3f,e[iX][iZ][0]=%.3f, t1=%.3f,t0=%.3f\n",
149 // iX,iZ,e[iX][iZ][1],e[iX][iZ][0], t[iX][iZ][1],t[iX][iZ][0]);
150 fHgLgRatio[iX][iZ]->Fill(e[iX][iZ][1]/e[iX][iZ][0]);
154 sprintf(hname,"%d_%d_%d",fMod,iX,iZ);
155 sprintf(htitl,"HG/LG ratio for crystal %d_%d_%d",fMod,iX,iZ);
156 fHgLgRatio[iX][iZ] = new TH1F(hname,htitl,400,14.,18.);
157 // printf("iX=%d iZ=%d,e[iX][iZ][1]=%.3f,e[iX][iZ][0]=%.3f\n",
158 // iX,iZ,e[iX][iZ][1],e[iX][iZ][0]);
159 fHgLgRatio[iX][iZ]->Fill(e[iX][iZ][1]/e[iX][iZ][0]);
160 fHistoArray.Add(fHgLgRatio[iX][iZ]);
165 for(Int_t iGain=0; iGain<2; iGain++) {
167 if(e[iX][iZ][iGain]<10) continue;
168 if(!(t[iX][iZ][iGain]>0)) continue;
170 if(fTimeEnergy[iX][iZ][iGain])
171 fTimeEnergy[iX][iZ][iGain]->Fill(e[iX][iZ][iGain],t[iX][iZ][iGain]);
173 sprintf(hname,"%d_%d_%d_%d",fMod,iX,iZ,iGain);
174 sprintf(htitl,"Energy vs TOF for crystal %d_%d_%d and gain %d",fMod,iX,iZ,iGain);
175 fTimeEnergy[iX][iZ][iGain] = new TH2F(hname,htitl,100,0.,1024.,100,0.,100);
176 fTimeEnergy[iX][iZ][iGain]->Fill(e[iX][iZ][iGain],t[iX][iZ][iGain]);
177 fHistoArray.Add(fTimeEnergy[iX][iZ][iGain]);
186 //-------------------------------------------------------------------
187 void AliPHOSRcuDA1::UpdateHistoFile()
189 // Write histograms to file
191 if(!fWriteToFile) return;
192 if(!fHistoFile) return;
193 if(!fHistoFile->IsOpen()) return;
200 //-----------------------------------------------------------------------
201 void AliPHOSRcuDA1::SetWriteToFile(Bool_t write)
204 fWriteToFile = write;
210 sprintf(rootname,"%s.root",GetName());
211 fHistoFile = new TFile(rootname,"update");
214 fWriteToFile = write;