Fixed compiler warnings
[u/mrichter/AliRoot.git] / FMD / AliFMDGainDA.cxx
CommitLineData
f631f26d 1/**************************************************************************
2 * Copyright(c) 2004, 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/** @file AliFMDGainDA.cxx
17 @author Hans Hjersing Dalsgaard <canute@nbi.dk>
18 @date Mon Mar 13 13:46:05 2008
19 @brief Derived class for the pulse gain detector algorithm.
20*/
cf291b42 21//____________________________________________________________________
22//
23// This class contains the implementation of the gain detector
24// algorithms (DA) for the FMD. The data is collected in histograms
25// that are reset for each pulse length after the mean and standard
26// deviation are put into a TGraphErrors object. After a certain
27// number of pulses (usually 8) the graph is fitted to a straight
28// line. The gain is then slope of this line as it combines the known
29// pulse and the response of the detector.
30//
f631f26d 31#include "AliFMDGainDA.h"
32#include "iostream"
33#include "fstream"
34#include "AliLog.h"
35#include "TF1.h"
36#include "TH1.h"
37#include "TMath.h"
38#include "TGraphErrors.h"
39#include "AliFMDParameters.h"
40
41//_____________________________________________________________________
42ClassImp(AliFMDGainDA)
cf291b42 43#if 0 // Do not delete - here to let Emacs indent properly
44;
45#endif
f631f26d 46
47//_____________________________________________________________________
cf291b42 48AliFMDGainDA::AliFMDGainDA()
49 : AliFMDBaseDA(),
50 fGainArray(),
51 fPulseSize(32),
52 fHighPulse(256),
53 fPulseLength(100),
54 fEventsPerChannel(0),
55 fCurrentPulse(0),
56 fCurrentChannel(0),
57 fNumberOfStripsPerChip(128)
f631f26d 58{
59 fOutputFile.open("gains.csv");
cf291b42 60 fGainArray.SetOwner();
f631f26d 61}
62
63//_____________________________________________________________________
cf291b42 64AliFMDGainDA::AliFMDGainDA(const AliFMDGainDA & gainDA)
65 : AliFMDBaseDA(gainDA),
66 fGainArray(gainDA.fGainArray),
67 fPulseSize(gainDA.fPulseSize),
68 fHighPulse(gainDA.fHighPulse),
69 fPulseLength(gainDA.fPulseLength),
70 fEventsPerChannel(gainDA.fEventsPerChannel),
71 fCurrentPulse(gainDA.fCurrentPulse),
72 fCurrentChannel(gainDA.fCurrentChannel),
73 fNumberOfStripsPerChip(gainDA.fNumberOfStripsPerChip)
74{
f631f26d 75}
76
77//_____________________________________________________________________
cf291b42 78AliFMDGainDA::~AliFMDGainDA()
79{
f631f26d 80}
81
82//_____________________________________________________________________
cf291b42 83void AliFMDGainDA::Init()
84{
f631f26d 85 fEventsPerChannel = (fPulseLength*fHighPulse) / fPulseSize ;
86
cf291b42 87 //8 pulser values * 128 strips * 100 samples
88 SetRequiredEvents(fEventsPerChannel*fNumberOfStripsPerChip);
f631f26d 89 TObjArray* detArray;
90 TObjArray* ringArray;
91 TObjArray* sectorArray;
92
93 for(UShort_t det=1;det<=3;det++) {
94 detArray = new TObjArray();
95 detArray->SetOwner();
96 fGainArray.AddAtAndExpand(detArray,det);
97 for (UShort_t ir = 0; ir < 2; ir++) {
98 Char_t ring = (ir == 0 ? 'O' : 'I');
99 UShort_t nsec = (ir == 0 ? 40 : 20);
100 UShort_t nstr = (ir == 0 ? 2 : 4);
101 ringArray = new TObjArray();
102 ringArray->SetOwner();
103 detArray->AddAtAndExpand(ringArray,ir);
104 for(UShort_t sec =0; sec < nsec; sec++) {
105 sectorArray = new TObjArray();
106 sectorArray->SetOwner();
107 ringArray->AddAtAndExpand(sectorArray,sec);
108 for(UShort_t strip = 0; strip < nstr; strip++) {
cf291b42 109 TH1S* hChannel = new TH1S(Form("hFMD%d%c_%d_%d",det,ring,sec,strip),
110 Form("hFMD%d%c_%d_%d",det,ring,sec,strip),
111 1024,0,1023);
f631f26d 112 hChannel->SetDirectory(0);
113 sectorArray->AddAtAndExpand(hChannel,strip);
114 }
115 }
116 }
117 }
118}
119
120//_____________________________________________________________________
121void AliFMDGainDA::AddChannelContainer(TObjArray* sectorArray,
cf291b42 122 UShort_t det ,
123 Char_t ring,
f631f26d 124 UShort_t sec,
cf291b42 125 UShort_t strip)
126{
f631f26d 127 TGraphErrors* hChannel = new TGraphErrors();
cf291b42 128 hChannel->SetName(Form("FMD%d%c[%02d,%03d]", det, ring, sec, strip));
129 hChannel->SetTitle(Form("FMD%d%c[%02d,%03d] ADC vs DAC",
130 det, ring, sec, strip));
f631f26d 131 sectorArray->AddAtAndExpand(hChannel,strip);
132}
133
134//_____________________________________________________________________
135void AliFMDGainDA::FillChannels(AliFMDDigit* digit) {
136
137 UShort_t det = digit->Detector();
138 Char_t ring = digit->Ring();
139 UShort_t sec = digit->Sector();
140 UShort_t strip = digit->Strip();
141
142
cf291b42 143 if(strip % fNumberOfStripsPerChip) return;
144 Int_t vaChip = strip / fNumberOfStripsPerChip;
145 TH1S* hChannel = GetChannelHistogram(det, ring, sec, vaChip);
f631f26d 146 hChannel->Fill(digit->Counts());
147 UpdatePulseAndADC(det,ring,sec,strip);
148}
149
150//_____________________________________________________________________
151void AliFMDGainDA::Analyse(UShort_t det,
152 Char_t ring,
153 UShort_t sec,
154 UShort_t strip) {
155 TGraphErrors* grChannel = GetChannel(det,ring,sec,strip);
156 if(!grChannel->GetN()) {
cf291b42 157 // AliWarning(Form("No entries for FMD%d%c, sector %d, strip %d",
158 // det, ring , sec, strip));
f631f26d 159 return;
160 }
9c958f2b 161 TF1 fitFunc("fitFunc","pol1",-10,280);
162 fitFunc.SetParameters(100,3);
f631f26d 163 grChannel->Fit("fitFunc","Q","Q",0,fHighPulse);
164 AliFMDParameters* pars = AliFMDParameters::Instance();
165 UInt_t ddl, board,chip,channel;
166 pars->Detector2Hardware(det,ring,sec,strip,ddl,board,chip,channel);
167
168 Float_t chi2ndf = 0;
9c958f2b 169 if(fitFunc.GetNDF())
170 chi2ndf = fitFunc.GetChisquare() / fitFunc.GetNDF();
f631f26d 171 ddl = ddl + kBaseDDL;
172
173 Int_t relStrip = strip%fNumberOfStripsPerChip;
174
175 fOutputFile << ddl << ','
176 << board << ','
177 << chip << ','
178 << channel << ','
179 << relStrip << ','
9c958f2b 180 << fitFunc.GetParameter(1) << ','
181 << fitFunc.GetParError(1) << ','
f631f26d 182 << chi2ndf <<"\n";
183
184
185 if(fSaveHistograms) {
cf291b42 186 gDirectory->cd(Form("%s:FMD%d%c/sector_%d/strip_%d",
187 fDiagnosticsFilename,det,ring,sec,strip));
f631f26d 188 grChannel->Write(Form("grFMD%d%c_%d_%d",det,ring,sec,strip));
189 }
190
9c958f2b 191
192
f631f26d 193}
194
195//_____________________________________________________________________
196void AliFMDGainDA::WriteHeaderToFile() {
9c958f2b 197 AliFMDParameters* pars = AliFMDParameters::Instance();
198 fOutputFile.write(Form("# %s \n",pars->GetGainShuttleID()),9);
f631f26d 199 fOutputFile.write("# Rcu, Board, Chip, Channel, Strip, Gain, Error, Chi2/NDF \n",59);
200
201}
202
203//_____________________________________________________________________
204TH1S* AliFMDGainDA::GetChannelHistogram(UShort_t det,
205 Char_t ring,
206 UShort_t sec,
207 UShort_t strip) {
208
209 UShort_t Ring = 1;
210 if(ring == 'O')
211 Ring = 0;
212
213
214 TObjArray* detArray = static_cast<TObjArray*>(fGainArray.At(det));
215 TObjArray* ringArray = static_cast<TObjArray*>(detArray->At(Ring));
216 TObjArray* secArray = static_cast<TObjArray*>(ringArray->At(sec));
217 TH1S* hChannel = static_cast<TH1S*>(secArray->At(strip));
218
219 return hChannel;
220}
221
222//_____________________________________________________________________
cf291b42 223TGraphErrors* AliFMDGainDA::GetChannel(UShort_t det,
224 Char_t ring,
225 UShort_t sec,
226 UShort_t strip) {
f631f26d 227
cf291b42 228 UShort_t iring = (ring == 'O' ? 0 : 1);
229 TObjArray* detArray = static_cast<TObjArray*>(fDetectorArray.At(det));
230 TObjArray* ringArray = static_cast<TObjArray*>(detArray->At(iring));
231 TObjArray* secArray = static_cast<TObjArray*>(ringArray->At(sec));
232 TGraphErrors* hChannel = static_cast<TGraphErrors*>(secArray->At(strip));
f631f26d 233
234 return hChannel;
235}
236
237//_____________________________________________________________________
cf291b42 238void AliFMDGainDA::UpdatePulseAndADC(UShort_t det,
239 Char_t ring,
240 UShort_t sec,
241 UShort_t strip)
242{
243 if(strip % fNumberOfStripsPerChip) return;
244 if(((GetCurrentEvent()) % fPulseLength) && GetCurrentEvent() > 0) return;
f631f26d 245
cf291b42 246 Int_t vaChip = strip/fNumberOfStripsPerChip;
247 TH1S* hChannel = GetChannelHistogram(det,ring,sec,vaChip);
f631f26d 248
249 if(!hChannel->GetEntries()) {
cf291b42 250 AliWarning(Form("No entries for FMD%d%c, sector %d, strip %d",
251 det, ring , sec, strip));
f631f26d 252 return;
253 }
254 Double_t mean = hChannel->GetMean();
255 Double_t rms = hChannel->GetRMS();
cf291b42 256 Double_t pulse = Double_t(fCurrentPulse) * fPulseSize;
257 Int_t firstBin = hChannel->GetXaxis()->GetFirst();
258 Int_t lastBin = hChannel->GetXaxis()->GetLast();
9c958f2b 259 hChannel->GetXaxis()->SetRangeUser(mean-4*rms,mean+4*rms);
f631f26d 260
cf291b42 261 mean = hChannel->GetMean();
262 rms = hChannel->GetRMS();
efd39294 263
264 hChannel->GetXaxis()->SetRange(firstBin,lastBin);
265
f631f26d 266 Int_t channelNumber = strip + (GetCurrentEvent()-1)/800;
267
268 TGraphErrors* channel = GetChannel(det,ring,sec,channelNumber);
269
270 channel->SetPoint(fCurrentPulse,pulse,mean);
271 channel->SetPointError(fCurrentPulse,0,rms);
272
273 if(fSaveHistograms) {
cf291b42 274 gDirectory->cd(Form("%s:FMD%d%c/sector_%d/strip_%d",
275 fDiagnosticsFilename,det,ring,sec,channelNumber));
276 hChannel->Write(Form("hFMD%d%c_%d_%d_pulse_%d",
277 det,ring,sec,channelNumber,fCurrentPulse));
f631f26d 278 }
279
280 hChannel->Reset();
281
282}
283
284//_____________________________________________________________________
cf291b42 285void AliFMDGainDA::ResetPulseAndUpdateChannel()
286{
287 fCurrentPulse = 0;
f631f26d 288}
289
290//_____________________________________________________________________
cf291b42 291void AliFMDGainDA::FinishEvent()
292{
293 if(GetCurrentEvent()>0 && (GetCurrentEvent() % fPulseLength == 0))
f631f26d 294 fCurrentPulse++;
295
cf291b42 296 if(GetCurrentEvent()>0 && (GetCurrentEvent()) % fEventsPerChannel == 0)
f631f26d 297 fCurrentPulse = 0;
f631f26d 298}
299//_____________________________________________________________________
300//
301//EOF
302//