1 | /* ======================================================================== *\
|
---|
2 | ! $Name: not supported by cvs2svn $:$Id: MExtractTimeAndChargeDigitalFilter.cc,v 1.77 2007-02-04 15:30:16 tbretz Exp $
|
---|
3 | ! --------------------------------------------------------------------------
|
---|
4 | !
|
---|
5 | ! *
|
---|
6 | ! * This file is part of MARS, the MAGIC Analysis and Reconstruction
|
---|
7 | ! * Software. It is distributed to you in the hope that it can be a useful
|
---|
8 | ! * and timesaving tool in analysing Data of imaging Cerenkov telescopes.
|
---|
9 | ! * It is distributed WITHOUT ANY WARRANTY.
|
---|
10 | ! *
|
---|
11 | ! * Permission to use, copy, modify and distribute this software and its
|
---|
12 | ! * documentation for any purpose is hereby granted without fee,
|
---|
13 | ! * provided that the above copyright notice appear in all copies and
|
---|
14 | ! * that both that copyright notice and this permission notice appear
|
---|
15 | ! * in supporting documentation. It is provided "as is" without express
|
---|
16 | ! * or implied warranty.
|
---|
17 | ! *
|
---|
18 | !
|
---|
19 | !
|
---|
20 | ! Author(s): Hendrik Bartko, 09/2004 <mailto:hbartko@mppmu.mpg.de>
|
---|
21 | ! Author(s): Markus Gaug, 05/2004 <mailto:markus@ifae.es>
|
---|
22 | ! Author(s): Diego Tescaro, 05/2004 <mailto:tescaro@pd.infn.it>
|
---|
23 | !
|
---|
24 | ! Copyright: MAGIC Software Development, 2000-2004
|
---|
25 | !
|
---|
26 | !
|
---|
27 | \* ======================================================================== */
|
---|
28 |
|
---|
29 | //////////////////////////////////////////////////////////////////////////////
|
---|
30 | //
|
---|
31 | // MExtractTimeAndChargeDigitalFilter
|
---|
32 | //
|
---|
33 | // Hendrik has promised to write more documentation
|
---|
34 | //
|
---|
35 | // The following variables have to be set by the derived class and
|
---|
36 | // do not have defaults:
|
---|
37 | // - fNumHiGainSamples
|
---|
38 | // - fNumLoGainSamples
|
---|
39 | // - fSqrtHiGainSamples
|
---|
40 | // - fSqrtLoGainSamples
|
---|
41 | //
|
---|
42 | // The reading of automatic weights files (color, type) can be switched
|
---|
43 | // off using EnableAutomaticWeights(kFALSE).
|
---|
44 | //
|
---|
45 | // An empty name or "-" as the weights file name is a synonym for
|
---|
46 | // setting all weights to 1
|
---|
47 | //
|
---|
48 | // The digital filter can even do a little extrapolation around the
|
---|
49 | // extraction window. For real Pulses having their maximum around
|
---|
50 | // slice 0 the extracted time is similar to a gaussian around 0 with
|
---|
51 | // sigma 0.4.
|
---|
52 | //
|
---|
53 | //
|
---|
54 | // Input Containers:
|
---|
55 | // MRawEvtData
|
---|
56 | // MRawRunHeader
|
---|
57 | // MPedestalCam
|
---|
58 | // [MCalibrationPattern]
|
---|
59 | //
|
---|
60 | // Output Containers:
|
---|
61 | // MArrivalTimeCam
|
---|
62 | // MExtractedSignalCam
|
---|
63 | //
|
---|
64 | //////////////////////////////////////////////////////////////////////////////
|
---|
65 | #include "MExtractTimeAndChargeDigitalFilter.h"
|
---|
66 |
|
---|
67 | #include <errno.h>
|
---|
68 | #include <fstream>
|
---|
69 |
|
---|
70 | #include <TRandom.h>
|
---|
71 |
|
---|
72 | #include "MLog.h"
|
---|
73 | #include "MLogManip.h"
|
---|
74 |
|
---|
75 | #include "MParList.h"
|
---|
76 |
|
---|
77 | #include "MRawRunHeader.h"
|
---|
78 | #include "MCalibrationPattern.h"
|
---|
79 | #include "MExtractedSignalCam.h"
|
---|
80 | #include "MExtralgoDigitalFilter.h"
|
---|
81 |
|
---|
82 | ClassImp(MExtractTimeAndChargeDigitalFilter);
|
---|
83 |
|
---|
84 | using namespace std;
|
---|
85 |
|
---|
86 | const Byte_t MExtractTimeAndChargeDigitalFilter::fgHiGainFirst = 0;
|
---|
87 | const Byte_t MExtractTimeAndChargeDigitalFilter::fgHiGainLast = 16;
|
---|
88 | const Int_t MExtractTimeAndChargeDigitalFilter::fgLoGainFirst = 1;
|
---|
89 | const Byte_t MExtractTimeAndChargeDigitalFilter::fgLoGainLast = 14;
|
---|
90 | const Int_t MExtractTimeAndChargeDigitalFilter::fgBinningResolutionHiGain = 10;
|
---|
91 | const Int_t MExtractTimeAndChargeDigitalFilter::fgBinningResolutionLoGain = 10;
|
---|
92 | const Float_t MExtractTimeAndChargeDigitalFilter::fgOffsetLoGain = 0.95;
|
---|
93 |
|
---|
94 | // --------------------------------------------------------------------------
|
---|
95 | //
|
---|
96 | // Default constructor.
|
---|
97 | //
|
---|
98 | // Calls:
|
---|
99 | // - SetWindowSize();
|
---|
100 | // - SetRange(fgHiGainFirst, fgHiGainLast, fgLoGainFirst, fgLoGainLast)
|
---|
101 | // - SetBinningResolution();
|
---|
102 | //
|
---|
103 | // Sets all weights to 1.
|
---|
104 | //
|
---|
105 | MExtractTimeAndChargeDigitalFilter::MExtractTimeAndChargeDigitalFilter(const char *name, const char *title)
|
---|
106 | : fBinningResolutionHiGain(fgBinningResolutionHiGain),
|
---|
107 | fBinningResolutionLoGain(fgBinningResolutionLoGain),
|
---|
108 | fAutomaticWeights(kTRUE)
|
---|
109 | {
|
---|
110 | fName = name ? name : "MExtractTimeAndChargeDigitalFilter";
|
---|
111 | fTitle = title ? title : "Digital Filter";
|
---|
112 |
|
---|
113 | SetRange(fgHiGainFirst, fgHiGainLast, fgLoGainFirst, fgLoGainLast);
|
---|
114 | SetWindowSize(3, 5);
|
---|
115 | SetOffsetLoGain(fgOffsetLoGain);
|
---|
116 | }
|
---|
117 |
|
---|
118 | // ---------------------------------------------------------------------------------------
|
---|
119 | //
|
---|
120 | // Checks:
|
---|
121 | // - if a window is bigger than the one defined by the ranges, set it
|
---|
122 | // to the available range
|
---|
123 | //
|
---|
124 | // Sets:
|
---|
125 | // - fNumHiGainSamples to: (Float_t)fWindowSizeHiGain
|
---|
126 | // - fNumLoGainSamples to: (Float_t)fWindowSizeLoGain
|
---|
127 | //
|
---|
128 | // This function might be used to turn the digital filter into a
|
---|
129 | // sliding window extractor by setting the filename to NULL
|
---|
130 | //
|
---|
131 | void MExtractTimeAndChargeDigitalFilter::SetWindowSize(Int_t windowh, Int_t windowl)
|
---|
132 | {
|
---|
133 | if (windowh > fHiGainLast-fHiGainFirst+1)
|
---|
134 | {
|
---|
135 | *fLog << err << "ERROR - The new hi-gain window size exceeds the extraction range." << endl;
|
---|
136 | return;
|
---|
137 | }
|
---|
138 |
|
---|
139 | if (windowl > (int)fLoGainLast-fLoGainFirst+1)
|
---|
140 | {
|
---|
141 | *fLog << err << "ERROR - The new hi-gain window size exceeds the extraction range." << endl;
|
---|
142 | return;
|
---|
143 | }
|
---|
144 |
|
---|
145 | fWindowSizeHiGain = windowh;
|
---|
146 | fWindowSizeLoGain = windowl;
|
---|
147 |
|
---|
148 | /*
|
---|
149 | const Int_t availhirange = (Int_t)();
|
---|
150 |
|
---|
151 | if (fWindowSizeHiGain > availhirange)
|
---|
152 | {
|
---|
153 | *fLog << warn << GetDescriptor() << ": Hi Gain window size: " << Form("%2i",fWindowSizeHiGain);
|
---|
154 | *fLog << " is bigger than available range: [" << Form("%2i", (int)fHiGainFirst);
|
---|
155 | *fLog << "," << Form("%21", (int)fHiGainLast) << "]" << endl;
|
---|
156 |
|
---|
157 | fHiGainLast = fHiGainFirst + fWindowSizeHiGain;
|
---|
158 |
|
---|
159 | *fLog << warn << GetDescriptor() << ": Will set the upper range to: " << (int)fHiGainLast << endl;
|
---|
160 | }
|
---|
161 |
|
---|
162 | if (fWindowSizeHiGain < 2)
|
---|
163 | {
|
---|
164 | fWindowSizeHiGain = 2;
|
---|
165 | *fLog << warn << GetDescriptor() << ": High Gain window size set to two samples" << endl;
|
---|
166 | }
|
---|
167 |
|
---|
168 | if (fLoGainLast != 0 && fWindowSizeLoGain != 0)
|
---|
169 | {
|
---|
170 | const Int_t availlorange = (Int_t)(fLoGainLast-fLoGainFirst+1);
|
---|
171 |
|
---|
172 | if (fWindowSizeLoGain > availlorange)
|
---|
173 | {
|
---|
174 | *fLog << warn << GetDescriptor() << ": Lo Gain window size: " << Form("%2i",fWindowSizeLoGain);
|
---|
175 | *fLog << " is bigger than available range: [" << Form("%2i", (int)fLoGainFirst);
|
---|
176 | *fLog << "," << Form("%21", (int)fLoGainLast) << "]" << endl;
|
---|
177 |
|
---|
178 | fLoGainLast = fLoGainFirst + fWindowSizeLoGain;
|
---|
179 |
|
---|
180 | *fLog << warn << GetDescriptor() << ": Will set the upper range to: " << (int)fLoGainLast << endl;
|
---|
181 | }
|
---|
182 |
|
---|
183 | if (fWindowSizeLoGain<2)
|
---|
184 | {
|
---|
185 | fWindowSizeLoGain = 2;
|
---|
186 | *fLog << warn << GetDescriptor() << ": Low Gain window size set to two samples" << endl;
|
---|
187 | }
|
---|
188 | }*/
|
---|
189 | //
|
---|
190 | // We need here the effective number of samples which is about 2.5 in the case of a window
|
---|
191 | // size of 6. The exact numbers have to be found still.
|
---|
192 | //
|
---|
193 | fNumHiGainSamples = fWindowSizeHiGain;
|
---|
194 | fNumLoGainSamples = fWindowSizeLoGain;
|
---|
195 | fSqrtHiGainSamples = TMath::Sqrt(fNumHiGainSamples);
|
---|
196 | fSqrtLoGainSamples = TMath::Sqrt(fNumLoGainSamples);
|
---|
197 | }
|
---|
198 |
|
---|
199 |
|
---|
200 | // --------------------------------------------------------------------------
|
---|
201 | //
|
---|
202 | // Executing MExtractTimeAndCharge::PreProcess and searching for
|
---|
203 | // MCalibrationPattern
|
---|
204 | //
|
---|
205 | Int_t MExtractTimeAndChargeDigitalFilter::PreProcess(MParList *pList)
|
---|
206 | {
|
---|
207 | if (!MExtractTimeAndCharge::PreProcess(pList))
|
---|
208 | return kFALSE;
|
---|
209 |
|
---|
210 | fCalibPattern = (MCalibrationPattern*)pList->FindObject("MCalibrationPattern");
|
---|
211 | return kTRUE;
|
---|
212 | }
|
---|
213 |
|
---|
214 | // --------------------------------------------------------------------------
|
---|
215 | //
|
---|
216 | // The weights are determined using GetAutimaticWeights().
|
---|
217 | //
|
---|
218 | // kFALSE is returned if it returned an error.
|
---|
219 | // kTRUE is returned if no new weights were set.
|
---|
220 | //
|
---|
221 | // If new weights are set
|
---|
222 | // fNumHiGainSamples
|
---|
223 | // fNumLoGainSamples
|
---|
224 | // fSqrtHiGainSamples
|
---|
225 | // fSqrtLoGainSamples
|
---|
226 | // and
|
---|
227 | // fSignals->SetUsedFADCSlices(...)
|
---|
228 | // is updated accordingly.
|
---|
229 | //
|
---|
230 | Bool_t MExtractTimeAndChargeDigitalFilter::GetWeights()
|
---|
231 | {
|
---|
232 | switch (GetAutomaticWeights())
|
---|
233 | {
|
---|
234 | case kERROR: // An error occured
|
---|
235 | return kFALSE;
|
---|
236 | case kFALSE: // No new weights set
|
---|
237 | return kTRUE;
|
---|
238 | }
|
---|
239 |
|
---|
240 | //
|
---|
241 | // We need here the effective number of samples. In pricipal the number
|
---|
242 | // is different depending on the weights used and must be set
|
---|
243 | // event by event.
|
---|
244 | //
|
---|
245 | fNumHiGainSamples = fAmpWeightsHiGain.GetSum()/fBinningResolutionHiGain;
|
---|
246 | fNumLoGainSamples = fAmpWeightsLoGain.GetSum()/fBinningResolutionLoGain;
|
---|
247 | fSqrtHiGainSamples = TMath::Sqrt(fNumHiGainSamples);
|
---|
248 | fSqrtLoGainSamples = TMath::Sqrt(fNumLoGainSamples);
|
---|
249 |
|
---|
250 | // From MExtractTimeAndCharge::ReInit
|
---|
251 | if (fSignals)
|
---|
252 | fSignals->SetUsedFADCSlices(fHiGainFirst, fHiGainLast/*+fHiLoLast*/, fNumHiGainSamples,
|
---|
253 | fLoGainFirst, fLoGainLast, fNumLoGainSamples);
|
---|
254 | return kTRUE;
|
---|
255 | }
|
---|
256 |
|
---|
257 | // --------------------------------------------------------------------------
|
---|
258 | //
|
---|
259 | // InitArrays
|
---|
260 | //
|
---|
261 | // Gets called in the ReInit() and initialized the arrays
|
---|
262 | //
|
---|
263 | Bool_t MExtractTimeAndChargeDigitalFilter::InitArrays(Int_t n)
|
---|
264 | {
|
---|
265 | if (!fRunHeader)
|
---|
266 | return kFALSE;
|
---|
267 |
|
---|
268 | return GetWeights();
|
---|
269 | }
|
---|
270 |
|
---|
271 | // --------------------------------------------------------------------------
|
---|
272 | //
|
---|
273 | // Check if reading a new weights file is necessary because the calibration
|
---|
274 | // pattern has changed. (Cannot be done in ReInit, because at this time
|
---|
275 | // the calibration pattern is not available.
|
---|
276 | // Then process the event.
|
---|
277 | //
|
---|
278 | Int_t MExtractTimeAndChargeDigitalFilter::Process()
|
---|
279 | {
|
---|
280 | // Change Weights if the calibration patter changes
|
---|
281 | if (!GetWeights())
|
---|
282 | return kERROR;
|
---|
283 |
|
---|
284 | // Process event
|
---|
285 | return MExtractTimeAndCharge::Process();
|
---|
286 | }
|
---|
287 |
|
---|
288 | // --------------------------------------------------------------------------
|
---|
289 | //
|
---|
290 | // Apply the digital filter algorithm to the high-gain slices.
|
---|
291 | //
|
---|
292 | void MExtractTimeAndChargeDigitalFilter::FindTimeAndChargeHiGain2(const Float_t *ptr, Int_t num,
|
---|
293 | Float_t &sum, Float_t &dsum,
|
---|
294 | Float_t &time, Float_t &dtime,
|
---|
295 | Byte_t sat, Int_t maxpos) const
|
---|
296 | {
|
---|
297 | // Do some handling if maxpos is last slice!
|
---|
298 |
|
---|
299 | MExtralgoDigitalFilter df(fBinningResolutionHiGain, fWindowSizeHiGain,
|
---|
300 | fAmpWeightsHiGain.GetArray(),
|
---|
301 | fTimeWeightsHiGain.GetArray(),
|
---|
302 | fPulseHiGain.GetArray());
|
---|
303 | df.SetData(num, ptr);
|
---|
304 |
|
---|
305 | if (IsNoiseCalculation())
|
---|
306 | {
|
---|
307 | sum = df.ExtractNoise(gRandom->Integer(fBinningResolutionHiGain));
|
---|
308 | return;
|
---|
309 | }
|
---|
310 |
|
---|
311 | df.Extract(/*maxpos*/);
|
---|
312 | df.GetSignal(sum, dsum);
|
---|
313 | df.GetTime(time, dtime);
|
---|
314 | }
|
---|
315 |
|
---|
316 | void MExtractTimeAndChargeDigitalFilter::FindTimeAndChargeLoGain2(const Float_t *ptr, Int_t num,
|
---|
317 | Float_t &sum, Float_t &dsum,
|
---|
318 | Float_t &time, Float_t &dtime,
|
---|
319 | Byte_t sat, Int_t maxpos) const
|
---|
320 | {
|
---|
321 | MExtralgoDigitalFilter df(fBinningResolutionLoGain, fWindowSizeLoGain,
|
---|
322 | fAmpWeightsLoGain.GetArray(),
|
---|
323 | fTimeWeightsLoGain.GetArray(),
|
---|
324 | fPulseLoGain.GetArray());
|
---|
325 |
|
---|
326 | df.SetData(num, ptr);
|
---|
327 |
|
---|
328 | if (IsNoiseCalculation())
|
---|
329 | {
|
---|
330 | sum = df.ExtractNoise(gRandom->Integer(fBinningResolutionHiGain));
|
---|
331 | return;
|
---|
332 | }
|
---|
333 |
|
---|
334 | df.Extract(/*maxpos*/);
|
---|
335 | df.GetSignal(sum, dsum);
|
---|
336 | df.GetTime(time, dtime);
|
---|
337 | }
|
---|
338 |
|
---|
339 |
|
---|
340 | // --------------------------------------------------------------------------
|
---|
341 | //
|
---|
342 | // Read the setup from a TEnv, eg:
|
---|
343 | // MJPedestal.MExtractor.WeightsFile: filename
|
---|
344 | // MJPedestal.MExtractor.AutomaticWeights: off
|
---|
345 | //
|
---|
346 | Int_t MExtractTimeAndChargeDigitalFilter::ReadEnv(const TEnv &env, TString prefix, Bool_t print)
|
---|
347 | {
|
---|
348 |
|
---|
349 | Bool_t rc = kFALSE;
|
---|
350 |
|
---|
351 | if (IsEnvDefined(env, prefix, "AutomaticWeights", print))
|
---|
352 | {
|
---|
353 | EnableAutomaticWeights(GetEnvValue(env, prefix, "AutomaticWeights", fAutomaticWeights));
|
---|
354 | rc = kTRUE;
|
---|
355 | }
|
---|
356 |
|
---|
357 | if (IsEnvDefined(env, prefix, "WeightsFile", print))
|
---|
358 | {
|
---|
359 | SetNameWeightsFile(GetEnvValue(env, prefix, "WeightsFile", ""));
|
---|
360 | rc = kTRUE;
|
---|
361 | }
|
---|
362 |
|
---|
363 | return MExtractTimeAndCharge::ReadEnv(env, prefix, print) ? kTRUE : rc;
|
---|
364 | }
|
---|
365 |
|
---|
366 | //----------------------------------------------------------------------------
|
---|
367 | //
|
---|
368 | // If automatic weights are requested, no default weights (name.IsNull())
|
---|
369 | // are requested, fRunHeader is available and fRunHeader->IsMonteCarloRun()
|
---|
370 | // is true prepend "MC_" in from of the name.
|
---|
371 | //
|
---|
372 | // return poth+name;
|
---|
373 | //
|
---|
374 | TString MExtractTimeAndChargeDigitalFilter::CompileWeightFileName(TString path, const TString &name) const
|
---|
375 | {
|
---|
376 | if (fAutomaticWeights && !name.IsNull() && fRunHeader && fRunHeader->IsMonteCarloRun())
|
---|
377 | path += "MC_";
|
---|
378 |
|
---|
379 | path += name;
|
---|
380 |
|
---|
381 | return path;
|
---|
382 | }
|
---|
383 |
|
---|
384 | //----------------------------------------------------------------------------
|
---|
385 | //
|
---|
386 | // Read a pre-defined weights file into the class.
|
---|
387 | // This is mandatory for the extraction
|
---|
388 | //
|
---|
389 | // If filenname is empty, then all weights will be set to 1.
|
---|
390 | //
|
---|
391 | // Returns:
|
---|
392 | // kTRUE: new weights set
|
---|
393 | // kFALSE: no weights set
|
---|
394 | // kERROR: error
|
---|
395 | //
|
---|
396 | Int_t MExtractTimeAndChargeDigitalFilter::ReadWeightsFile(TString filename, TString path)
|
---|
397 | {
|
---|
398 | if (filename.IsNull())
|
---|
399 | {
|
---|
400 | fAmpWeightsHiGain .Set(fBinningResolutionHiGain*fWindowSizeHiGain);
|
---|
401 | fAmpWeightsLoGain .Set(fBinningResolutionLoGain*fWindowSizeLoGain);
|
---|
402 | fTimeWeightsHiGain.Set(fBinningResolutionHiGain*fWindowSizeHiGain);
|
---|
403 | fTimeWeightsLoGain.Set(fBinningResolutionLoGain*fWindowSizeLoGain);
|
---|
404 |
|
---|
405 | fAmpWeightsHiGain.Reset(1);
|
---|
406 | fTimeWeightsHiGain.Reset(1);
|
---|
407 | fAmpWeightsLoGain.Reset(1);
|
---|
408 | fTimeWeightsLoGain.Reset(1);
|
---|
409 | return kTRUE;
|
---|
410 | }
|
---|
411 |
|
---|
412 | // Add "MC_" in front of the filename if necessary
|
---|
413 | filename = CompileWeightFileName(path, filename);
|
---|
414 |
|
---|
415 | //filename = MJob::ExpandPath(filename);
|
---|
416 |
|
---|
417 | if (fNameWeightsFileSet==filename)
|
---|
418 | return kFALSE; // No file read
|
---|
419 |
|
---|
420 | ifstream fin(filename.Data());
|
---|
421 | if (!fin)
|
---|
422 | {
|
---|
423 | *fLog << err << GetDescriptor() << ": ERROR - Cannot open file " << filename << ": ";
|
---|
424 | *fLog << strerror(errno) << endl;
|
---|
425 | return kERROR;
|
---|
426 | }
|
---|
427 |
|
---|
428 | *fLog << all << GetDescriptor() << ": Reading weights in " << filename << "..." << flush;
|
---|
429 |
|
---|
430 | Int_t len = 0;
|
---|
431 | Int_t cnt = 0;
|
---|
432 | Int_t line = 0;
|
---|
433 | Bool_t hi = kFALSE;
|
---|
434 | Bool_t lo = kFALSE;
|
---|
435 |
|
---|
436 | TString str;
|
---|
437 |
|
---|
438 | while (1)
|
---|
439 | {
|
---|
440 | str.ReadLine(fin);
|
---|
441 | if (!fin)
|
---|
442 | break;
|
---|
443 |
|
---|
444 | line++;
|
---|
445 |
|
---|
446 | if (str.Contains("# High Gain Weights:"))
|
---|
447 | {
|
---|
448 | if (hi)
|
---|
449 | {
|
---|
450 | *fLog << err << "ERROR - 'High Gain Weights' found twice in line #" << line << "." << endl;
|
---|
451 | return kERROR;
|
---|
452 | }
|
---|
453 |
|
---|
454 | if (2!=sscanf(str.Data(), "# High Gain Weights: %2i %2i", &fWindowSizeHiGain, &fBinningResolutionHiGain))
|
---|
455 | {
|
---|
456 | *fLog << err << "ERROR - Wrong number of arguments in line #" << line << ":" << endl;
|
---|
457 | *fLog << str << endl;
|
---|
458 | return kERROR;
|
---|
459 | }
|
---|
460 |
|
---|
461 | len = fBinningResolutionHiGain*fWindowSizeHiGain;
|
---|
462 | fAmpWeightsHiGain .Set(len);
|
---|
463 | fTimeWeightsHiGain.Set(len);
|
---|
464 | fPulseHiGain.Set(len);
|
---|
465 | hi = kTRUE;
|
---|
466 | continue;
|
---|
467 | }
|
---|
468 |
|
---|
469 | if (str.Contains("# Low Gain Weights:"))
|
---|
470 | {
|
---|
471 | if (lo)
|
---|
472 | {
|
---|
473 | *fLog << err << "ERROR - 'Lo Gain Weights' found twice in line #" << line << "." << endl;
|
---|
474 | return kERROR;
|
---|
475 | }
|
---|
476 |
|
---|
477 | if (2!=sscanf(str.Data(),"# Low Gain Weights: %2i %2i", &fWindowSizeLoGain, &fBinningResolutionLoGain))
|
---|
478 | {
|
---|
479 | *fLog << err << "ERROR - Wrong number of arguments in line #" << line << ":" << endl;
|
---|
480 | *fLog << str << endl;
|
---|
481 | return kERROR;
|
---|
482 | }
|
---|
483 |
|
---|
484 | len = fBinningResolutionLoGain*fWindowSizeLoGain;
|
---|
485 | fAmpWeightsLoGain .Set(len);
|
---|
486 | fTimeWeightsLoGain.Set(len);
|
---|
487 | fPulseLoGain.Set(len);
|
---|
488 | lo = kTRUE;
|
---|
489 | continue;
|
---|
490 | }
|
---|
491 |
|
---|
492 | // Handle lines with comments
|
---|
493 | if (str.Contains("#"))
|
---|
494 | continue;
|
---|
495 |
|
---|
496 | // Nothing found so far
|
---|
497 | if (len == 0)
|
---|
498 | continue;
|
---|
499 |
|
---|
500 | if (3!=sscanf(str.Data(), "%f %f %f",
|
---|
501 | lo ? &fAmpWeightsLoGain [cnt] : &fAmpWeightsHiGain [cnt],
|
---|
502 | lo ? &fTimeWeightsLoGain[cnt] : &fTimeWeightsHiGain[cnt],
|
---|
503 | lo ? &fPulseLoGain[cnt] : &fPulseHiGain[cnt]))
|
---|
504 | {
|
---|
505 | *fLog << err << "ERROR - Wrong number of arguments in line #" << line << ":" << endl;
|
---|
506 | *fLog << str << endl;
|
---|
507 | return kERROR;
|
---|
508 | }
|
---|
509 |
|
---|
510 | if (++cnt == len)
|
---|
511 | {
|
---|
512 | len = 0;
|
---|
513 | cnt = 0;
|
---|
514 | }
|
---|
515 | }
|
---|
516 |
|
---|
517 | if (cnt != len)
|
---|
518 | {
|
---|
519 | *fLog << err << "ERROR - Size mismatch in weights file " << filename << endl;
|
---|
520 | return kERROR;
|
---|
521 | }
|
---|
522 |
|
---|
523 | if (!hi)
|
---|
524 | {
|
---|
525 | *fLog << err << "ERROR - No correct header found in weights file " << filename << endl;
|
---|
526 | return kERROR;
|
---|
527 | }
|
---|
528 |
|
---|
529 | *fLog << "done." << endl;
|
---|
530 |
|
---|
531 | *fLog << inf << " File contains " << fWindowSizeHiGain << " hi-gain slices ";
|
---|
532 | *fLog << "with a resolution of " << fBinningResolutionHiGain << endl;
|
---|
533 |
|
---|
534 | *fLog << inf << " File contains " << fWindowSizeLoGain << " lo-gain slices ";
|
---|
535 | *fLog << "with a resolution of " << fBinningResolutionLoGain << endl;
|
---|
536 |
|
---|
537 | //CalcBinningResArrays();
|
---|
538 |
|
---|
539 | switch (fWindowSizeHiGain)
|
---|
540 | {
|
---|
541 | case 4:
|
---|
542 | SetResolutionPerPheHiGain(0.036);
|
---|
543 | break;
|
---|
544 | case 6:
|
---|
545 | SetResolutionPerPheHiGain(0.021);
|
---|
546 | break;
|
---|
547 | default:
|
---|
548 | *fLog << warn << "Could not set the high-gain extractor resolution per phe for window size "
|
---|
549 | << fWindowSizeHiGain << endl;
|
---|
550 | }
|
---|
551 |
|
---|
552 | switch (fWindowSizeLoGain)
|
---|
553 | {
|
---|
554 | case 4:
|
---|
555 | SetResolutionPerPheLoGain(0.005);
|
---|
556 | break;
|
---|
557 | case 6:
|
---|
558 | SetResolutionPerPheLoGain(0.004);
|
---|
559 | break;
|
---|
560 | default:
|
---|
561 | *fLog << warn << "Could not set the low-gain extractor resolution per phe for window size "
|
---|
562 | << fWindowSizeLoGain << endl;
|
---|
563 | }
|
---|
564 |
|
---|
565 | fNameWeightsFileSet = filename;
|
---|
566 |
|
---|
567 | return kTRUE;
|
---|
568 | }
|
---|
569 |
|
---|
570 |
|
---|
571 | //----------------------------------------------------------------------------
|
---|
572 | //
|
---|
573 | // The default (+ prepending possible "MC_") is read for:
|
---|
574 | //
|
---|
575 | // - RunType: Pedestal (independant of fAutomaticWeights)
|
---|
576 | // - fAutomaticWeights disabled
|
---|
577 | //
|
---|
578 | // if fAutomaticWeights enabled:
|
---|
579 | // - fNameWeightsFile.IsNull()
|
---|
580 | // - !fCalibPattern
|
---|
581 | // - fCalibPattern->GetPulserColor()==MCalibrationCam::kNONE
|
---|
582 | //
|
---|
583 | // If automatic weights are enabled, the case above didn't take place and
|
---|
584 | // fNameWeightsFile starts with "calibration_weights_"
|
---|
585 | // - the color (blue, UV) is replaced by the appropriate one
|
---|
586 | // taken from the calibration pattern
|
---|
587 | //
|
---|
588 | // In most cases a debug output is printed. Further output about the color
|
---|
589 | // determination can be switched on with debug level > 5;
|
---|
590 | //
|
---|
591 | // Returns:
|
---|
592 | // kFALSE: No new weights set
|
---|
593 | // kTRUE: New weights set
|
---|
594 | // kERROR: Error
|
---|
595 | //
|
---|
596 | Int_t MExtractTimeAndChargeDigitalFilter::GetAutomaticWeights()
|
---|
597 | {
|
---|
598 | const Ssiz_t pos = fNameWeightsFile.Last('/')+1;
|
---|
599 | const Ssiz_t len = fNameWeightsFile.Length();
|
---|
600 |
|
---|
601 | // Split file name in path and name
|
---|
602 | TString path = fNameWeightsFile(0, pos>=0?pos:len);
|
---|
603 | TString name = fNameWeightsFile(pos>=0?pos:0, len);
|
---|
604 |
|
---|
605 | // Remove trailing "MC_" for automatic weights
|
---|
606 | if (fAutomaticWeights && name.BeginsWith("MC_"))
|
---|
607 | name.Remove(0, 3);
|
---|
608 |
|
---|
609 | // In case of a pedetsal run no calibration pattern can be available
|
---|
610 | // the default weights are always used.
|
---|
611 | if (fRunHeader->GetRunType()==MRawRunHeader::kRTPedestal)
|
---|
612 | {
|
---|
613 | *fLog << dbg << "Pedestal file... using default weights: " << fNameWeightsFile << endl;
|
---|
614 | return ReadWeightsFile(name, path);
|
---|
615 | }
|
---|
616 |
|
---|
617 | // If automatic weights are switched off use default weights
|
---|
618 | if (!fAutomaticWeights)
|
---|
619 | {
|
---|
620 | *fLog << dbg << "Automatic weights switched off... using default weights: " << fNameWeightsFile << endl;
|
---|
621 | return ReadWeightsFile(name, path);
|
---|
622 | }
|
---|
623 |
|
---|
624 | // If automatic weights are switched on but no filename is given raise error
|
---|
625 | if (fNameWeightsFile.IsNull())
|
---|
626 | {
|
---|
627 | *fLog << err << "ERROR - Cannot get automatic weights without default filename." << endl;
|
---|
628 | return kERROR;
|
---|
629 | }
|
---|
630 |
|
---|
631 | // If this is no pedestal run, automatic weights are requested and a
|
---|
632 | // filename for the weights file is given pedestal-extraction from
|
---|
633 | // cosmics data is assumed.
|
---|
634 | if (!fCalibPattern)
|
---|
635 | {
|
---|
636 | *fLog << dbg << "No decoded calibration pattern available... using default weights: " << fNameWeightsFile << endl;
|
---|
637 | return ReadWeightsFile(name, path);
|
---|
638 | }
|
---|
639 |
|
---|
640 | const Bool_t debug = gLog.GetDebugLevel()>5;
|
---|
641 |
|
---|
642 | // If no calibration pattern is available do not change the
|
---|
643 | // current weighs or current weights file name.
|
---|
644 | if (fCalibPattern->GetPulserColor()==MCalibrationCam::kNONE)
|
---|
645 | {
|
---|
646 | // If we are extracting data and the calibration pattern is kNONE
|
---|
647 | // we assume that it is a data file without interleaved events
|
---|
648 | // and calibration pattern information available.
|
---|
649 | if ((fRunHeader->GetRunType()!=MRawRunHeader::kRTData && !fRunHeader->IsMonteCarloRun()) || debug)
|
---|
650 | *fLog << dbg << "No calibration color set so far... guessing default: " << fNameWeightsFile << endl;
|
---|
651 |
|
---|
652 | return ReadWeightsFile(name, path);
|
---|
653 | }
|
---|
654 |
|
---|
655 | if (debug)
|
---|
656 | {
|
---|
657 | *fLog << dbg << endl;
|
---|
658 | *fLog << underline << GetDescriptor() << endl;
|
---|
659 | *fLog << " Trying to get automatic weight for " << fNameWeightsFile << endl;
|
---|
660 | *fLog << " Run type: ";
|
---|
661 | }
|
---|
662 |
|
---|
663 | if (name.BeginsWith("calibration_weights_") && fCalibPattern)
|
---|
664 | {
|
---|
665 | if (debug)
|
---|
666 | *fLog << " Calibration with color " << fCalibPattern->GetPulserColorStr() << ", setting ";
|
---|
667 | switch (fCalibPattern->GetPulserColor())
|
---|
668 | {
|
---|
669 | case MCalibrationCam::kBLUE: // 2
|
---|
670 | case MCalibrationCam::kGREEN: // 1
|
---|
671 | if (debug)
|
---|
672 | *fLog << "blue/green, ";
|
---|
673 | name.ReplaceAll("UV", "blue");
|
---|
674 | break;
|
---|
675 |
|
---|
676 | case MCalibrationCam::kUV: // 3
|
---|
677 | case MCalibrationCam::kCT1: // 0
|
---|
678 | if (debug)
|
---|
679 | *fLog << "UV/CT1, ";
|
---|
680 | name.ReplaceAll("blue", "UV");
|
---|
681 | break;
|
---|
682 | case MCalibrationCam::kNONE:
|
---|
683 | break;
|
---|
684 | default: // kNone + etc
|
---|
685 | *fLog << err << "ERROR - Cannot get automatic weights for " << fCalibPattern->GetPulserColorStr() << endl;
|
---|
686 | return kERROR;
|
---|
687 | }
|
---|
688 | }
|
---|
689 |
|
---|
690 | return ReadWeightsFile(name, path);
|
---|
691 | }
|
---|
692 |
|
---|
693 | //----------------------------------------------------------------------------
|
---|
694 | //
|
---|
695 | // Print the setup of the digital filter extraction used. Use
|
---|
696 | // the option "weights" if you want to print also all weights.
|
---|
697 | //
|
---|
698 | void MExtractTimeAndChargeDigitalFilter::Print(Option_t *o) const
|
---|
699 | {
|
---|
700 | if (IsA()==Class())
|
---|
701 | *fLog << GetDescriptor() << ":" << endl;
|
---|
702 |
|
---|
703 | MExtractTimeAndCharge::Print(o);
|
---|
704 | *fLog << " Window Size HiGain: " << setw(2) << fWindowSizeHiGain << " LoGain: " << setw(2) << fWindowSizeLoGain << endl;
|
---|
705 | *fLog << " Binning Res HiGain: " << setw(2) << fBinningResolutionHiGain << " LoGain: " << setw(2) << fBinningResolutionHiGain << endl;
|
---|
706 | *fLog << " Weights File desired: " << (fNameWeightsFile.IsNull()?"-":fNameWeightsFile.Data()) << endl;
|
---|
707 | if (!fNameWeightsFileSet.IsNull())
|
---|
708 | *fLog << " Weights File set: " << fNameWeightsFileSet << endl;
|
---|
709 |
|
---|
710 | TString opt(o);
|
---|
711 | if (!opt.Contains("weights"))
|
---|
712 | return;
|
---|
713 |
|
---|
714 | *fLog << endl;
|
---|
715 | *fLog << inf << "Using the following weights: " << endl;
|
---|
716 | *fLog << "Hi-Gain:" << endl;
|
---|
717 | for (Int_t i=0; i<fBinningResolutionHiGain*fWindowSizeHiGain; i++)
|
---|
718 | *fLog << " " << fAmpWeightsHiGain[i] << " \t " << fTimeWeightsHiGain[i] << endl;
|
---|
719 |
|
---|
720 | *fLog << "Lo-Gain:" << endl;
|
---|
721 | for (Int_t i=0; i<fBinningResolutionLoGain*fWindowSizeLoGain; i++)
|
---|
722 | *fLog << " " << fAmpWeightsLoGain[i] << " \t " << fTimeWeightsLoGain[i] << endl;
|
---|
723 | }
|
---|