1 | /* ======================================================================== *\
|
---|
2 | !
|
---|
3 | ! *
|
---|
4 | ! * This file is part of MARS, the MAGIC Analysis and Reconstruction
|
---|
5 | ! * Software. It is distributed to you in the hope that it can be a useful
|
---|
6 | ! * and timesaving tool in analysing Data of imaging Cerenkov telescopes.
|
---|
7 | ! * It is distributed WITHOUT ANY WARRANTY.
|
---|
8 | ! *
|
---|
9 | ! * Permission to use, copy, modify and distribute this software and its
|
---|
10 | ! * documentation for any purpose is hereby granted without fee,
|
---|
11 | ! * provided that the above copyright notice appear in all copies and
|
---|
12 | ! * that both that copyright notice and this permission notice appear
|
---|
13 | ! * in supporting documentation. It is provided "as is" without express
|
---|
14 | ! * or implied warranty.
|
---|
15 | ! *
|
---|
16 | !
|
---|
17 | !
|
---|
18 | ! Author(s): Markus Gaug 02/2004 <mailto:markus@ifae.es>
|
---|
19 | ! Author(s): Thomas Bretz <mailto:tbretz@astro.uni-wuerzburg.de>
|
---|
20 | !
|
---|
21 | ! Copyright: MAGIC Software Development, 2000-2007
|
---|
22 | !
|
---|
23 | !
|
---|
24 | \* ======================================================================== */
|
---|
25 |
|
---|
26 | /////////////////////////////////////////////////////////////////////////////
|
---|
27 | //
|
---|
28 | // MHCalibrationRelTimeCam
|
---|
29 | //
|
---|
30 | // Fills the extracted relative arrival times of MArrivalTimeCam into
|
---|
31 | // the MHCalibrationPix-classes MHCalibrationPix for every:
|
---|
32 | //
|
---|
33 | // - Pixel, stored in the TObjArray's MHCalibrationCam::fHiGainArray
|
---|
34 | // or MHCalibrationCam::fHiGainArray, respectively, depending if
|
---|
35 | // MArrivalTimePix::IsLoGainUsed() is set.
|
---|
36 | //
|
---|
37 | // - Average pixel per AREA index (e.g. inner and outer for the MAGIC camera),
|
---|
38 | // stored in the TObjArray's MHCalibrationCam::fAverageHiGainAreas and
|
---|
39 | // MHCalibrationCam::fAverageHiGainAreas
|
---|
40 | //
|
---|
41 | // - Average pixel per camera SECTOR (e.g. sectors 1-6 for the MAGIC camera),
|
---|
42 | // stored in the TObjArray's MHCalibrationCam::fAverageHiGainSectors
|
---|
43 | // and MHCalibrationCam::fAverageHiGainSectors
|
---|
44 | //
|
---|
45 | // Every relative time is calculated as the difference between the individual
|
---|
46 | // pixel arrival time and the one of pixel 1 (hardware number: 2).
|
---|
47 | // The relative times are filled into a histogram and an array, in order to perform
|
---|
48 | // a Fourier analysis (see MHGausEvents). The signals are moreover averaged on an
|
---|
49 | // event-by-event basis and written into the corresponding average pixels.
|
---|
50 | //
|
---|
51 | // The histograms are fitted to a Gaussian, mean and sigma with its errors
|
---|
52 | // and the fit probability are extracted. If none of these values are NaN's and
|
---|
53 | // if the probability is bigger than MHGausEvents::fProbLimit (default: 0.5%),
|
---|
54 | // the fit is declared valid.
|
---|
55 | // Otherwise, the fit is repeated within ranges of the previous mean
|
---|
56 | // - MHCalibrationPix::fPickupLimit (default: 5) sigma (see MHCalibrationPix::RepeatFit())
|
---|
57 | // In case this does not make the fit valid, the histogram means and RMS's are
|
---|
58 | // taken directly (see MHCalibrationPix::BypassFit()) and the following flags are set:
|
---|
59 | // - MBadPixelsPix::SetUncalibrated( MBadPixelsPix::kRelTimeNotFitted ) and
|
---|
60 | // - MBadPixelsPix::SetUnsuitable( MBadPixelsPix::kUnreliableRun )
|
---|
61 | //
|
---|
62 | // Outliers of more than MHCalibrationPix::fPickupLimit (default: 5) sigmas
|
---|
63 | // from the mean are counted as Pickup events (stored in MHCalibrationPix::fPickup)
|
---|
64 | //
|
---|
65 | // The class also fills arrays with the signal vs. event number, creates a fourier
|
---|
66 | // spectrum (see MHGausEvents::CreateFourierSpectrum()) and investigates if the
|
---|
67 | // projected fourier components follow an exponential distribution.
|
---|
68 | // In case that the probability of the exponential fit is less than
|
---|
69 | // MHGausEvents::fProbLimit (default: 0.5%), the following flags are set:
|
---|
70 | // - MBadPixelsPix::SetUncalibrated( MBadPixelsPix::kRelTimeOscillating ) and
|
---|
71 | // - MBadPixelsPix::SetUnsuitable( MBadPixelsPix::kUnreliableRun )
|
---|
72 | //
|
---|
73 | // This same procedure is performed for the average pixels.
|
---|
74 | //
|
---|
75 | // The following results are written into MCalibrationRelTimeCam:
|
---|
76 | //
|
---|
77 | // - MCalibrationPix::SetMean()
|
---|
78 | // - MCalibrationPix::SetMeanErr()
|
---|
79 | // - MCalibrationPix::SetSigma()
|
---|
80 | // - MCalibrationPix::SetSigmaErr()
|
---|
81 | // - MCalibrationPix::SetProb()
|
---|
82 | // - MCalibrationPix::SetNumPickup()
|
---|
83 | //
|
---|
84 | // For all averaged areas, the fitted sigma is multiplied with the square root of
|
---|
85 | // the number involved pixels in order to be able to compare it to the average of
|
---|
86 | // sigmas in the camera.
|
---|
87 | //
|
---|
88 | /////////////////////////////////////////////////////////////////////////////
|
---|
89 | #include "MHCalibrationRelTimeCam.h"
|
---|
90 | #include "MHCalibrationPix.h"
|
---|
91 |
|
---|
92 | #include "MLog.h"
|
---|
93 | #include "MLogManip.h"
|
---|
94 |
|
---|
95 | #include "MParList.h"
|
---|
96 |
|
---|
97 | #include "MCalibrationIntensityRelTimeCam.h"
|
---|
98 |
|
---|
99 | #include "MCalibrationRelTimeCam.h"
|
---|
100 | #include "MCalibrationRelTimePix.h"
|
---|
101 | #include "MCalibrationPix.h"
|
---|
102 |
|
---|
103 | #include "MArrivalTimeCam.h"
|
---|
104 | #include "MArrivalTimePix.h"
|
---|
105 |
|
---|
106 | #include "MGeomCam.h"
|
---|
107 | #include "MGeomPix.h"
|
---|
108 |
|
---|
109 | #include "MBadPixelsIntensityCam.h"
|
---|
110 | #include "MBadPixelsCam.h"
|
---|
111 | #include "MBadPixelsPix.h"
|
---|
112 |
|
---|
113 | #include <TOrdCollection.h>
|
---|
114 | #include <TPad.h>
|
---|
115 | #include <TVirtualPad.h>
|
---|
116 | #include <TCanvas.h>
|
---|
117 | #include <TStyle.h>
|
---|
118 | #include <TF1.h>
|
---|
119 | #include <TLine.h>
|
---|
120 | #include <TLatex.h>
|
---|
121 | #include <TLegend.h>
|
---|
122 | #include <TGraph.h>
|
---|
123 | #include <TEnv.h>
|
---|
124 |
|
---|
125 | ClassImp(MHCalibrationRelTimeCam);
|
---|
126 |
|
---|
127 | using namespace std;
|
---|
128 |
|
---|
129 | const Float_t MHCalibrationRelTimeCam::fgNumHiGainSaturationLimit = 0.25;
|
---|
130 | const UInt_t MHCalibrationRelTimeCam::fgReferencePixel = 1;
|
---|
131 | const Int_t MHCalibrationRelTimeCam::fgNbins = 400;
|
---|
132 | const Axis_t MHCalibrationRelTimeCam::fgFirst = -9.975;
|
---|
133 | const Axis_t MHCalibrationRelTimeCam::fgLast = 10.025;
|
---|
134 | const Float_t MHCalibrationRelTimeCam::fgProbLimit = 0.0;
|
---|
135 | const TString MHCalibrationRelTimeCam::gsHistName = "RelTime";
|
---|
136 | const TString MHCalibrationRelTimeCam::gsHistTitle = "Arr. Times";
|
---|
137 | const TString MHCalibrationRelTimeCam::gsHistXTitle = "Arr. Time [FADC slices]";
|
---|
138 | const TString MHCalibrationRelTimeCam::gsHistYTitle = "Nr. events";
|
---|
139 | const TString MHCalibrationRelTimeCam::fgReferenceFile = "mjobs/calibrationref.rc";
|
---|
140 |
|
---|
141 | // --------------------------------------------------------------------------
|
---|
142 | //
|
---|
143 | // Default Constructor.
|
---|
144 | //
|
---|
145 | // Sets:
|
---|
146 | // - fReferencePixel to fgReferencePixel
|
---|
147 | // - fNbins to fgNbins
|
---|
148 | // - fFirst to fgFirst
|
---|
149 | // - fLast to fgLast
|
---|
150 | //
|
---|
151 | // - fHistName to gsHistName
|
---|
152 | // - fHistTitle to gsHistTitle
|
---|
153 | // - fHistXTitle to gsHistXTitle
|
---|
154 | // - fHistYTitle to gsHistYTitle
|
---|
155 | //
|
---|
156 | MHCalibrationRelTimeCam::MHCalibrationRelTimeCam(const char *name, const char *title)
|
---|
157 | {
|
---|
158 |
|
---|
159 | fName = name ? name : "MHCalibrationRelTimeCam";
|
---|
160 | fTitle = title ? title : "Histogram class for the relative time calibration of the camera";
|
---|
161 |
|
---|
162 | SetNumHiGainSaturationLimit(fgNumHiGainSaturationLimit);
|
---|
163 |
|
---|
164 | SetReferencePixel();
|
---|
165 |
|
---|
166 | SetBinning(fgNbins, fgFirst, fgLast);
|
---|
167 |
|
---|
168 | SetProbLimit(fgProbLimit);
|
---|
169 |
|
---|
170 | SetHistName (gsHistName .Data());
|
---|
171 | SetHistTitle (gsHistTitle .Data());
|
---|
172 | SetHistXTitle(gsHistXTitle.Data());
|
---|
173 | SetHistYTitle(gsHistYTitle.Data());
|
---|
174 |
|
---|
175 | SetReferenceFile();
|
---|
176 |
|
---|
177 | fInnerRefTime = 2.95;
|
---|
178 | fOuterRefTime = 3.6;
|
---|
179 | }
|
---|
180 |
|
---|
181 | // --------------------------------------------------------------------------
|
---|
182 | //
|
---|
183 | // Creates new MHCalibrationRelTimeCam only with the averaged areas:
|
---|
184 | // the rest has to be retrieved directly, e.g. via:
|
---|
185 | // MHCalibrationRelTimeCam *cam = MParList::FindObject("MHCalibrationRelTimeCam");
|
---|
186 | // - cam->GetAverageSector(5).DrawClone();
|
---|
187 | // - (*cam)[100].DrawClone()
|
---|
188 | //
|
---|
189 | TObject *MHCalibrationRelTimeCam::Clone(const char *) const
|
---|
190 | {
|
---|
191 |
|
---|
192 | MHCalibrationRelTimeCam *cam = new MHCalibrationRelTimeCam();
|
---|
193 |
|
---|
194 | //
|
---|
195 | // Copy the data members
|
---|
196 | //
|
---|
197 | cam->fColor = fColor;
|
---|
198 | cam->fRunNumbers = fRunNumbers;
|
---|
199 | cam->fPulserFrequency = fPulserFrequency;
|
---|
200 | cam->fFlags = fFlags;
|
---|
201 | cam->fNbins = fNbins;
|
---|
202 | cam->fFirst = fFirst;
|
---|
203 | cam->fLast = fLast;
|
---|
204 |
|
---|
205 | cam->fReferenceFile = fReferenceFile;
|
---|
206 | cam->fInnerRefTime = fInnerRefTime;
|
---|
207 | cam->fOuterRefTime = fOuterRefTime;
|
---|
208 |
|
---|
209 | //
|
---|
210 | // Copy the MArrays
|
---|
211 | //
|
---|
212 | cam->fAverageAreaRelSigma = fAverageAreaRelSigma;
|
---|
213 | cam->fAverageAreaRelSigmaVar = fAverageAreaRelSigmaVar;
|
---|
214 | cam->fAverageAreaSat = fAverageAreaSat;
|
---|
215 | cam->fAverageAreaSigma = fAverageAreaSigma;
|
---|
216 | cam->fAverageAreaSigmaVar = fAverageAreaSigmaVar;
|
---|
217 | cam->fAverageAreaNum = fAverageAreaNum;
|
---|
218 | cam->fAverageSectorNum = fAverageSectorNum;
|
---|
219 |
|
---|
220 | if (!IsAverageing())
|
---|
221 | return cam;
|
---|
222 |
|
---|
223 | const Int_t navhi = fAverageHiGainAreas->GetSize();
|
---|
224 |
|
---|
225 | for (int i=0; i<navhi; i++)
|
---|
226 | cam->fAverageHiGainAreas->AddAt(GetAverageHiGainArea(i).Clone(),i);
|
---|
227 |
|
---|
228 | if (IsLoGain())
|
---|
229 | {
|
---|
230 |
|
---|
231 | const Int_t navlo = fAverageLoGainAreas->GetSize();
|
---|
232 | for (int i=0; i<navlo; i++)
|
---|
233 | cam->fAverageLoGainAreas->AddAt(GetAverageLoGainArea(i).Clone(),i);
|
---|
234 |
|
---|
235 | }
|
---|
236 |
|
---|
237 | return cam;
|
---|
238 | }
|
---|
239 |
|
---|
240 | // --------------------------------------------------------------------------
|
---|
241 | //
|
---|
242 | // Gets or creates the pointers to:
|
---|
243 | // - MCalibrationRelTimeCam
|
---|
244 | //
|
---|
245 | // Searches pointer to:
|
---|
246 | // - MArrivalTimeCam
|
---|
247 | //
|
---|
248 | // Calls:
|
---|
249 | // - MHCalibrationCam::InitHiGainArrays()
|
---|
250 | // - MHCalibrationCam::InitLoGainArrays()
|
---|
251 | //
|
---|
252 | // Sets:
|
---|
253 | // - fSumareahi to nareas
|
---|
254 | // - fSumarealo to nareas
|
---|
255 | // - fSumsectorhi to nareas
|
---|
256 | // - fSumsectorlo to nareas
|
---|
257 | // - fNumareahi to nareas
|
---|
258 | // - fNumarealo to nareas
|
---|
259 | // - fNumsectorhi to nareas
|
---|
260 | // - fNumsectorlo to nareas
|
---|
261 | //
|
---|
262 | Bool_t MHCalibrationRelTimeCam::ReInitHists(MParList *pList)
|
---|
263 | {
|
---|
264 |
|
---|
265 | if (!InitCams(pList,"RelTime"))
|
---|
266 | return kFALSE;
|
---|
267 |
|
---|
268 | MArrivalTimeCam *signal = (MArrivalTimeCam*)pList->FindObject("MArrivalTimeCam");
|
---|
269 | if (!signal)
|
---|
270 | {
|
---|
271 | *fLog << err << "MArrivalTimeCam not found... abort." << endl;
|
---|
272 | return kFALSE;
|
---|
273 | }
|
---|
274 |
|
---|
275 | const Int_t npixels = fGeom->GetNumPixels();
|
---|
276 | const Int_t nsectors = fGeom->GetNumSectors();
|
---|
277 | const Int_t nareas = fGeom->GetNumAreas();
|
---|
278 |
|
---|
279 | InitHiGainArrays(npixels,nareas,nsectors);
|
---|
280 | InitLoGainArrays(npixels,nareas,nsectors);
|
---|
281 |
|
---|
282 | fSumareahi .Set(nareas);
|
---|
283 | fSumarealo .Set(nareas);
|
---|
284 | fSumsectorhi.Set(nsectors);
|
---|
285 | fSumsectorlo.Set(nsectors);
|
---|
286 | fNumareahi .Set(nareas);
|
---|
287 | fNumarealo .Set(nareas);
|
---|
288 | fNumsectorhi.Set(nsectors);
|
---|
289 | fNumsectorlo.Set(nsectors);
|
---|
290 |
|
---|
291 | return kTRUE;
|
---|
292 | }
|
---|
293 |
|
---|
294 |
|
---|
295 | // -------------------------------------------------------------------------------
|
---|
296 | //
|
---|
297 | // Retrieves pointer to MArrivalTimeCam:
|
---|
298 | //
|
---|
299 | // Retrieves from MGeomCam:
|
---|
300 | // - number of pixels
|
---|
301 | // - number of pixel areas
|
---|
302 | // - number of sectors
|
---|
303 | //
|
---|
304 | // Fills HiGain or LoGain histograms (MHGausEvents::FillHistAndArray()), respectively
|
---|
305 | // depending on MArrivalTimePix::IsLoGainUsed(), with:
|
---|
306 | // - MArrivalTimePix::GetArrivalTime(pixid) - MArrivalTimePix::GetArrivalTime(1);
|
---|
307 | // (i.e. the time difference between pixel i and pixel 1 (hardware number: 2) )
|
---|
308 | //
|
---|
309 | Bool_t MHCalibrationRelTimeCam::FillHists(const MParContainer *par, const Stat_t w)
|
---|
310 | {
|
---|
311 |
|
---|
312 | MArrivalTimeCam *arrtime = (MArrivalTimeCam*)par;
|
---|
313 | if (!arrtime)
|
---|
314 | {
|
---|
315 | gLog << err << "No argument in MArrivalTime::Fill... abort." << endl;
|
---|
316 | return kFALSE;
|
---|
317 | }
|
---|
318 |
|
---|
319 | const Int_t npixels = fGeom->GetNumPixels();
|
---|
320 | const Int_t nareas = fGeom->GetNumAreas();
|
---|
321 | const Int_t nsectors = fGeom->GetNumSectors();
|
---|
322 |
|
---|
323 | fSumareahi .Reset();
|
---|
324 | fSumarealo .Reset();
|
---|
325 | fSumsectorhi.Reset();
|
---|
326 | fSumsectorlo.Reset();
|
---|
327 | fNumareahi .Reset();
|
---|
328 | fNumarealo .Reset();
|
---|
329 | fNumsectorhi.Reset();
|
---|
330 | fNumsectorlo.Reset();
|
---|
331 |
|
---|
332 | const MArrivalTimePix &refpix = (*arrtime)[fReferencePixel];
|
---|
333 | // FIXME: What do we do if pixel is invalid?
|
---|
334 | if (!refpix.IsArrivalTimeValid())
|
---|
335 | {
|
---|
336 | gLog << warn << "WARNING - Arrival time in refrence pixel " << fReferencePixel << " invalid." << endl;
|
---|
337 | return kTRUE;
|
---|
338 | }
|
---|
339 |
|
---|
340 | const Float_t reftime = refpix.GetArrivalTime();
|
---|
341 | // refpix.IsHiGainSaturated()
|
---|
342 | // ? refpix.GetArrivalTimeLoGain() : refpix.GetArrivalTimeHiGain();
|
---|
343 |
|
---|
344 | for (Int_t i=0; i<npixels; i++)
|
---|
345 | {
|
---|
346 | MHCalibrationPix &histhi = (*this)[i];
|
---|
347 | if (histhi.IsExcluded())
|
---|
348 | continue;
|
---|
349 |
|
---|
350 | const Int_t aidx = (*fGeom)[i].GetAidx();
|
---|
351 | const Int_t sector = (*fGeom)[i].GetSector();
|
---|
352 |
|
---|
353 | const MArrivalTimePix &pix = (*arrtime)[i];
|
---|
354 |
|
---|
355 | // If hi-gain arrival time has been extracted successfully
|
---|
356 | // fill hi-gain histograms and arrays
|
---|
357 | if (pix.IsHiGainValid() && !pix.IsHiGainSaturated())
|
---|
358 | {
|
---|
359 | const Float_t time = pix.GetArrivalTimeHiGain();
|
---|
360 |
|
---|
361 | if (IsOscillations())
|
---|
362 | histhi.FillHistAndArray(time-reftime);
|
---|
363 | else
|
---|
364 | histhi.FillHist(time-reftime);
|
---|
365 |
|
---|
366 | fSumareahi [aidx] += time;
|
---|
367 | fNumareahi [aidx] ++;
|
---|
368 | fSumsectorhi[sector] += time;
|
---|
369 | fNumsectorhi[sector] ++;
|
---|
370 | }
|
---|
371 |
|
---|
372 | if (!pix.IsHiGainSaturated())
|
---|
373 | continue;
|
---|
374 |
|
---|
375 | histhi.AddSaturated(1);
|
---|
376 | /*
|
---|
377 | // If lo-gain arrival time has been extracted successfully,
|
---|
378 | // the hi-gain has saturateed and the lo-gain is switched on
|
---|
379 | // fill hi-gain histograms and arrays
|
---|
380 | if (pix.IsLoGainValid() && IsLoGain())
|
---|
381 | {
|
---|
382 | const Float_t time = pix.GetArrivalTimeLoGain();
|
---|
383 |
|
---|
384 | MHCalibrationPix &histlo = (*this)(i);
|
---|
385 | if (IsOscillations())
|
---|
386 | histlo.FillHistAndArray(time-reftime);
|
---|
387 | else
|
---|
388 | histlo.FillHist(time-reftime);
|
---|
389 |
|
---|
390 | fSumarealo [aidx] += time;
|
---|
391 | fNumarealo [aidx] ++;
|
---|
392 | fSumsectorlo[sector] += time;
|
---|
393 | fNumsectorlo[sector] ++;
|
---|
394 | }*/
|
---|
395 | }
|
---|
396 |
|
---|
397 | for (Int_t j=0; j<nareas; j++)
|
---|
398 | {
|
---|
399 | MHCalibrationPix &histhi = GetAverageHiGainArea(j);
|
---|
400 | if (IsOscillations())
|
---|
401 | histhi.FillHistAndArray(fNumareahi[j] == 0 ? 0. : fSumareahi[j]/fNumareahi[j]);
|
---|
402 | else
|
---|
403 | histhi.FillHist(fNumareahi[j] == 0 ? 0. : fSumareahi[j]/fNumareahi[j]);
|
---|
404 |
|
---|
405 | if (IsLoGain())
|
---|
406 | {
|
---|
407 | MHCalibrationPix &histlo = GetAverageLoGainArea(j);
|
---|
408 | if (IsOscillations())
|
---|
409 | histlo.FillHistAndArray(fNumarealo[j] == 0 ? 0. : fSumarealo[j]/fNumarealo[j]);
|
---|
410 | else
|
---|
411 | histlo.FillHist(fNumarealo[j] == 0 ? 0. : fSumarealo[j]/fNumarealo[j]);
|
---|
412 | }
|
---|
413 | }
|
---|
414 |
|
---|
415 | for (Int_t j=0; j<nsectors; j++)
|
---|
416 | {
|
---|
417 | MHCalibrationPix &histhi = GetAverageHiGainSector(j);
|
---|
418 | if (IsOscillations())
|
---|
419 | histhi.FillHistAndArray(fNumsectorhi[j] == 0 ? 0. : fSumsectorhi[j]/fNumsectorhi[j]);
|
---|
420 | else
|
---|
421 | histhi.FillHist(fNumsectorhi[j] == 0 ? 0. : fSumsectorhi[j]/fNumsectorhi[j]);
|
---|
422 |
|
---|
423 | if (IsLoGain())
|
---|
424 | {
|
---|
425 | MHCalibrationPix &histlo = GetAverageLoGainSector(j);
|
---|
426 | if (IsOscillations())
|
---|
427 | histlo.FillHistAndArray(fNumsectorlo[j] == 0 ? 0. : fSumsectorlo[j]/fNumsectorlo[j]);
|
---|
428 | else
|
---|
429 | histlo.FillHist(fNumsectorlo[j] == 0 ? 0. : fSumsectorlo[j]/fNumsectorlo[j]);
|
---|
430 | }
|
---|
431 | }
|
---|
432 |
|
---|
433 | return kTRUE;
|
---|
434 | }
|
---|
435 |
|
---|
436 | // --------------------------------------------------------------------------
|
---|
437 | //
|
---|
438 | // Calls:
|
---|
439 | // - MHCalibrationCam::FitHiGainArrays() with flags:
|
---|
440 | // MBadPixelsPix::kRelTimeNotFitted and MBadPixelsPix::kRelTimeOscillating
|
---|
441 | // - MHCalibrationCam::FitLoGainArrays() with flags:
|
---|
442 | // MBadPixelsPix::kRelTimeNotFitted and MBadPixelsPix::kRelTimeOscillating
|
---|
443 | //
|
---|
444 | Bool_t MHCalibrationRelTimeCam::FinalizeHists()
|
---|
445 | {
|
---|
446 |
|
---|
447 | *fLog << endl;
|
---|
448 |
|
---|
449 | MCalibrationCam *relcam = fIntensCam ? fIntensCam->GetCam() : fCam;
|
---|
450 | MBadPixelsCam *badcam = fIntensBad ? fIntensBad->GetCam() : fBadPixels;
|
---|
451 |
|
---|
452 | const Int_t nareas = fAverageHiGainAreas->GetSize();
|
---|
453 | const Int_t nsectors = fAverageHiGainSectors->GetSize();
|
---|
454 |
|
---|
455 | TArrayI satarea(nareas);
|
---|
456 | TArrayI satsect(nsectors);
|
---|
457 | fNumareahi .Reset();
|
---|
458 | fNumsectorhi.Reset();
|
---|
459 |
|
---|
460 | for (Int_t i=0; i<fHiGainArray->GetSize(); i++)
|
---|
461 | {
|
---|
462 |
|
---|
463 | MHCalibrationPix &histhi = (*this)[i];
|
---|
464 |
|
---|
465 | if (histhi.IsExcluded())
|
---|
466 | continue;
|
---|
467 |
|
---|
468 | const Int_t aidx = (*fGeom)[i].GetAidx();
|
---|
469 | const Int_t sector = (*fGeom)[i].GetSector();
|
---|
470 |
|
---|
471 | MCalibrationRelTimePix &pix = (MCalibrationRelTimePix&)(*relcam)[i] ;
|
---|
472 |
|
---|
473 | fNumareahi[aidx]++;
|
---|
474 | fNumsectorhi[sector]++;
|
---|
475 | //
|
---|
476 | // Check saturation
|
---|
477 | //
|
---|
478 | if (histhi.GetSaturated() > fNumHiGainSaturationLimit*histhi.GetHGausHist()->GetEntries())
|
---|
479 | {
|
---|
480 | pix.SetHiGainSaturation();
|
---|
481 | histhi.SetExcluded();
|
---|
482 | satarea[aidx]++;
|
---|
483 | satsect[sector]++;
|
---|
484 | }
|
---|
485 | else
|
---|
486 | if (IsLoGain())
|
---|
487 | (*this)(i).SetExcluded();
|
---|
488 |
|
---|
489 | //
|
---|
490 | // Check histogram overflow
|
---|
491 | //
|
---|
492 | CheckOverflow(histhi);
|
---|
493 | if (IsLoGain())
|
---|
494 | CheckOverflow((*this)(i));
|
---|
495 |
|
---|
496 | }
|
---|
497 |
|
---|
498 | for (Int_t j=0; j<nareas; j++)
|
---|
499 | {
|
---|
500 |
|
---|
501 | MHCalibrationPix &histhi = GetAverageHiGainArea(j);
|
---|
502 | MCalibrationRelTimePix &pix = (MCalibrationRelTimePix&)relcam->GetAverageArea(j);
|
---|
503 |
|
---|
504 | if (satarea[j] > 0.5*fNumareahi[j])
|
---|
505 | {
|
---|
506 | pix.SetHiGainSaturation();
|
---|
507 | histhi.SetExcluded();
|
---|
508 | }
|
---|
509 | else
|
---|
510 | if (IsLoGain())
|
---|
511 | GetAverageLoGainArea(j).SetExcluded();
|
---|
512 |
|
---|
513 | //
|
---|
514 | // Check histogram overflow
|
---|
515 | //
|
---|
516 | CheckOverflow(histhi);
|
---|
517 | if (IsLoGain())
|
---|
518 | CheckOverflow(GetAverageLoGainArea(j));
|
---|
519 | }
|
---|
520 |
|
---|
521 | for (Int_t j=0; j<fAverageHiGainSectors->GetSize(); j++)
|
---|
522 | {
|
---|
523 |
|
---|
524 | MHCalibrationPix &histhi = GetAverageHiGainSector(j);
|
---|
525 | MCalibrationRelTimePix &pix = (MCalibrationRelTimePix&)relcam->GetAverageSector(j) ;
|
---|
526 |
|
---|
527 | if (satsect[j] > 0.5*fNumsectorhi[j])
|
---|
528 | {
|
---|
529 | pix.SetHiGainSaturation();
|
---|
530 | histhi.SetExcluded();
|
---|
531 | }
|
---|
532 | else
|
---|
533 | if (IsLoGain())
|
---|
534 | GetAverageLoGainSector(j).SetExcluded();
|
---|
535 |
|
---|
536 | //
|
---|
537 | // Check histogram overflow
|
---|
538 | //
|
---|
539 | CheckOverflow(histhi);
|
---|
540 | if (IsLoGain())
|
---|
541 | CheckOverflow(GetAverageLoGainSector(j));
|
---|
542 | }
|
---|
543 |
|
---|
544 | FitHiGainArrays(*relcam,*badcam,
|
---|
545 | MBadPixelsPix::kRelTimeNotFitted,
|
---|
546 | MBadPixelsPix::kRelTimeOscillating);
|
---|
547 |
|
---|
548 | if (IsLoGain())
|
---|
549 | FitLoGainArrays(*relcam,*badcam,
|
---|
550 | MBadPixelsPix::kRelTimeNotFitted,
|
---|
551 | MBadPixelsPix::kRelTimeOscillating);
|
---|
552 |
|
---|
553 | return kTRUE;
|
---|
554 | }
|
---|
555 |
|
---|
556 | // --------------------------------------------------------------------------
|
---|
557 | //
|
---|
558 | // Sets all pixels to MBadPixelsPix::kUnreliableRun, if following flags are set:
|
---|
559 | // - MBadPixelsPix::kRelTimeNotFitted
|
---|
560 | // - MBadPixelsPix::kRelTimeOscillating
|
---|
561 | //
|
---|
562 | void MHCalibrationRelTimeCam::FinalizeBadPixels()
|
---|
563 | {
|
---|
564 |
|
---|
565 | MBadPixelsCam *badcam = fIntensBad ? fIntensBad->GetCam() : fBadPixels;
|
---|
566 |
|
---|
567 | for (Int_t i=0; i<badcam->GetSize(); i++)
|
---|
568 | {
|
---|
569 | MBadPixelsPix &bad = (*badcam)[i];
|
---|
570 |
|
---|
571 | if (bad.IsUncalibrated( MBadPixelsPix::kRelTimeNotFitted ))
|
---|
572 | bad.SetUnsuitable( MBadPixelsPix::kUnreliableRun );
|
---|
573 |
|
---|
574 | if (bad.IsUncalibrated( MBadPixelsPix::kRelTimeOscillating))
|
---|
575 | bad.SetUnsuitable( MBadPixelsPix::kUnreliableRun );
|
---|
576 |
|
---|
577 | }
|
---|
578 | }
|
---|
579 |
|
---|
580 | // --------------------------------------------------------------------------
|
---|
581 | //
|
---|
582 | // The types are as follows:
|
---|
583 | //
|
---|
584 | // Fitted values:
|
---|
585 | // ==============
|
---|
586 | //
|
---|
587 | // 0: Fitted Mean Relative Arrival Time in FADC slices (MHGausEvents::GetMean()
|
---|
588 | // 1: Error Mean Relative Arrival Time in FADC slices (MHGausEvents::GetMeanErr()
|
---|
589 | // 2: Sigma fitted Relative Arrival Time in FADC slices (MHGausEvents::GetSigma()
|
---|
590 | // 3: Error Sigma Relative Arrival Time in FADC slices (MHGausEvents::GetSigmaErr()
|
---|
591 | //
|
---|
592 | // Useful variables derived from the fit results:
|
---|
593 | // =============================================
|
---|
594 | //
|
---|
595 | // 4: Returned probability of Gauss fit (calls: MHGausEvents::GetProb())
|
---|
596 | //
|
---|
597 | // Localized defects:
|
---|
598 | // ==================
|
---|
599 | //
|
---|
600 | // 5: Gaus fit not OK (calls: MHGausEvents::IsGausFitOK())
|
---|
601 | // 6: Fourier spectrum not OK (calls: MHGausEvents::IsFourierSpectrumOK())
|
---|
602 | //
|
---|
603 | Bool_t MHCalibrationRelTimeCam::GetPixelContent(Double_t &val, Int_t idx, const MGeomCam &cam, Int_t type) const
|
---|
604 | {
|
---|
605 |
|
---|
606 | if (fHiGainArray->GetSize() <= idx)
|
---|
607 | return kFALSE;
|
---|
608 |
|
---|
609 | const MHCalibrationPix &pix = (*this)[idx];
|
---|
610 |
|
---|
611 | switch (type)
|
---|
612 | {
|
---|
613 | case 0:
|
---|
614 | val = pix.GetMean();
|
---|
615 | break;
|
---|
616 | case 1:
|
---|
617 | val = pix.GetMeanErr();
|
---|
618 | break;
|
---|
619 | case 2:
|
---|
620 | val = pix.GetSigma();
|
---|
621 | break;
|
---|
622 | case 3:
|
---|
623 | val = pix.GetSigmaErr();
|
---|
624 | break;
|
---|
625 | case 4:
|
---|
626 | val = pix.GetProb();
|
---|
627 | break;
|
---|
628 | case 5:
|
---|
629 | if (!pix.IsGausFitOK())
|
---|
630 | val = 1.;
|
---|
631 | break;
|
---|
632 | case 6:
|
---|
633 | if (!pix.IsFourierSpectrumOK())
|
---|
634 | val = 1.;
|
---|
635 | break;
|
---|
636 | default:
|
---|
637 | return kFALSE;
|
---|
638 | }
|
---|
639 | return kTRUE;
|
---|
640 | }
|
---|
641 |
|
---|
642 | // --------------------------------------------------------------------------
|
---|
643 | //
|
---|
644 | // Calls MHCalibrationPix::DrawClone() for pixel idx
|
---|
645 | //
|
---|
646 | void MHCalibrationRelTimeCam::DrawPixelContent(Int_t idx) const
|
---|
647 | {
|
---|
648 | (*this)[idx].DrawClone();
|
---|
649 | }
|
---|
650 |
|
---|
651 | // -----------------------------------------------------------------------------
|
---|
652 | //
|
---|
653 | // Default draw:
|
---|
654 | //
|
---|
655 | // Displays the averaged areas, both High Gain and Low Gain
|
---|
656 | //
|
---|
657 | // Calls the Draw of the fAverageHiGainAreas and fAverageLoGainAreas objects with options
|
---|
658 | //
|
---|
659 | void MHCalibrationRelTimeCam::Draw(const Option_t *opt)
|
---|
660 | {
|
---|
661 |
|
---|
662 | const Int_t nareas = fAverageHiGainAreas->GetSize();
|
---|
663 | if (nareas == 0)
|
---|
664 | return;
|
---|
665 |
|
---|
666 | TString option(opt);
|
---|
667 | option.ToLower();
|
---|
668 |
|
---|
669 | if (!option.Contains("datacheck"))
|
---|
670 | {
|
---|
671 | MHCalibrationCam::Draw(opt);
|
---|
672 | return;
|
---|
673 | }
|
---|
674 |
|
---|
675 | //
|
---|
676 | // From here on , the datacheck - Draw
|
---|
677 | //
|
---|
678 | TVirtualPad *pad = gPad ? gPad : MH::MakeDefCanvas(this);
|
---|
679 | pad->SetBorderMode(0);
|
---|
680 | pad->Divide(1,nareas);
|
---|
681 |
|
---|
682 | //
|
---|
683 | // Loop over inner and outer pixels
|
---|
684 | //
|
---|
685 | for (Int_t i=0; i<nareas;i++)
|
---|
686 | {
|
---|
687 |
|
---|
688 | pad->cd(i+1);
|
---|
689 |
|
---|
690 | MHCalibrationPix &hipix = GetAverageHiGainArea(i);
|
---|
691 | //
|
---|
692 | // Ask for Hi-Gain saturation
|
---|
693 | //
|
---|
694 | if (hipix.IsExcluded() && IsLoGain())
|
---|
695 | {
|
---|
696 | MHCalibrationPix &lopix = GetAverageLoGainArea(i);
|
---|
697 | DrawDataCheckPixel(lopix,i ? fOuterRefTime+1.5 : fInnerRefTime+1.5);
|
---|
698 | }
|
---|
699 | else
|
---|
700 | DrawDataCheckPixel(hipix,i ? fOuterRefTime : fInnerRefTime);
|
---|
701 | }
|
---|
702 | }
|
---|
703 |
|
---|
704 | void MHCalibrationRelTimeCam::CheckOverflow( MHCalibrationPix &pix ) const
|
---|
705 | {
|
---|
706 | if (pix.IsExcluded())
|
---|
707 | return;
|
---|
708 |
|
---|
709 | const TH1F &hist = *pix.GetHGausHist();
|
---|
710 |
|
---|
711 | const Int_t n = hist.GetNbinsX();
|
---|
712 | const Float_t max = fOverflowLimit*hist.GetEntries();
|
---|
713 |
|
---|
714 | const Stat_t overflow = hist.GetBinContent(n+1);
|
---|
715 | if (overflow > max)
|
---|
716 | {
|
---|
717 | *fLog << warn << overflow << " overflows above " << hist.GetBinLowEdge(n);
|
---|
718 | *fLog << " in " << pix.GetName() << " (w/o saturation!) " << endl;
|
---|
719 | }
|
---|
720 |
|
---|
721 | const Stat_t underflow = hist.GetBinContent(0);
|
---|
722 | if (underflow > max)
|
---|
723 | {
|
---|
724 | *fLog << warn << underflow << " underflows below " << hist.GetBinLowEdge(1);
|
---|
725 | *fLog << " in " << pix.GetName() << " (w/o saturation!) " << endl;
|
---|
726 | }
|
---|
727 | }
|
---|
728 |
|
---|
729 |
|
---|
730 | // -----------------------------------------------------------------------------
|
---|
731 | //
|
---|
732 | // Draw the average pixel for the datacheck:
|
---|
733 | //
|
---|
734 | // Displays the averaged areas, both High Gain and Low Gain
|
---|
735 | //
|
---|
736 | // Calls the Draw of the fAverageHiGainAreas and fAverageLoGainAreas objects with options
|
---|
737 | //
|
---|
738 | void MHCalibrationRelTimeCam::DrawDataCheckPixel(MHCalibrationPix &pix, const Float_t refline)
|
---|
739 | {
|
---|
740 | if (pix.IsEmpty())
|
---|
741 | return;
|
---|
742 |
|
---|
743 | TVirtualPad *pad = gPad;
|
---|
744 | pad->Divide(1,2, 1e-10, 1e-10);
|
---|
745 | pad->cd(1);
|
---|
746 |
|
---|
747 | gPad->SetBorderMode(0);
|
---|
748 |
|
---|
749 | gPad->SetTicks();
|
---|
750 | if (!pix.IsEmpty() && !pix.IsOnlyOverflow() && !pix.IsOnlyUnderflow())
|
---|
751 | gPad->SetLogy();
|
---|
752 |
|
---|
753 | TH1F *hist = pix.GetHGausHist();
|
---|
754 |
|
---|
755 | //
|
---|
756 | // set the labels bigger
|
---|
757 | //
|
---|
758 | TAxis *xaxe = hist->GetXaxis();
|
---|
759 | TAxis *yaxe = hist->GetYaxis();
|
---|
760 | xaxe->CenterTitle();
|
---|
761 | yaxe->CenterTitle();
|
---|
762 | xaxe->SetTitleSize(0.07);
|
---|
763 | yaxe->SetTitleSize(0.07);
|
---|
764 | xaxe->SetTitleOffset(0.65);
|
---|
765 | yaxe->SetTitleOffset(0.55);
|
---|
766 | xaxe->SetLabelSize(0.06);
|
---|
767 | yaxe->SetLabelSize(0.06);
|
---|
768 | xaxe->SetTitle(hist->GetXaxis()->GetTitle());
|
---|
769 | yaxe->SetTitle(hist->GetYaxis()->GetTitle());
|
---|
770 | xaxe->SetRange(hist->GetMaximumBin()-1500, hist->GetMaximumBin()+1500);
|
---|
771 |
|
---|
772 | hist->Draw();
|
---|
773 |
|
---|
774 | gStyle->SetOptFit();
|
---|
775 |
|
---|
776 | TF1 *fit = pix.GetFGausFit();
|
---|
777 |
|
---|
778 | if (fit)
|
---|
779 | {
|
---|
780 | switch (fColor)
|
---|
781 | {
|
---|
782 | case MCalibrationCam::kGREEN:
|
---|
783 | fit->SetLineColor(kGreen);
|
---|
784 | break;
|
---|
785 | case MCalibrationCam::kBLUE:
|
---|
786 | fit->SetLineColor(kBlue);
|
---|
787 | break;
|
---|
788 | case MCalibrationCam::kUV:
|
---|
789 | fit->SetLineColor(106);
|
---|
790 | break;
|
---|
791 | case MCalibrationCam::kCT1:
|
---|
792 | fit->SetLineColor(006);
|
---|
793 | break;
|
---|
794 | default:
|
---|
795 | fit->SetLineColor(kRed);
|
---|
796 | }
|
---|
797 | fit->Draw("same");
|
---|
798 | }
|
---|
799 |
|
---|
800 | DisplayRefLines(hist, refline);
|
---|
801 |
|
---|
802 | pad->cd(2);
|
---|
803 | gPad->SetBorderMode(0);
|
---|
804 | gPad->SetTicks();
|
---|
805 |
|
---|
806 | pix.CreateGraphEvents();
|
---|
807 | TGraph *gr = pix.GetGraphEvents();
|
---|
808 | if (gr)
|
---|
809 | {
|
---|
810 | TH1F *null2 = gr->GetHistogram();
|
---|
811 |
|
---|
812 | null2->SetMinimum(pix.GetMean()-10.*pix.GetSigma());
|
---|
813 | null2->SetMaximum(pix.GetMean()+10.*pix.GetSigma());
|
---|
814 |
|
---|
815 | //
|
---|
816 | // set the labels bigger
|
---|
817 | //
|
---|
818 | TAxis *xaxe2 = null2->GetXaxis();
|
---|
819 | TAxis *yaxe2 = null2->GetYaxis();
|
---|
820 | xaxe2->CenterTitle();
|
---|
821 | yaxe2->CenterTitle();
|
---|
822 | xaxe2->SetTitleSize(0.07);
|
---|
823 | yaxe2->SetTitleSize(0.07);
|
---|
824 | xaxe2->SetTitleOffset(0.65);
|
---|
825 | yaxe2->SetTitleOffset(0.55);
|
---|
826 | xaxe2->SetLabelSize(0.06);
|
---|
827 | yaxe2->SetLabelSize(0.06);
|
---|
828 | xaxe2->SetRangeUser(0.,pix.GetEvents()->GetSize()/pix.GetEventFrequency());
|
---|
829 | }
|
---|
830 |
|
---|
831 | pix.DrawEvents();
|
---|
832 | }
|
---|
833 |
|
---|
834 | void MHCalibrationRelTimeCam::DisplayRefLines(const TH1F *hist, const Float_t refline) const
|
---|
835 | {
|
---|
836 | TLine *line = new TLine(refline, 0, refline, hist->GetMaximum());
|
---|
837 | line->SetLineColor(kGreen);
|
---|
838 | line->SetLineStyle(2);
|
---|
839 | line->SetLineWidth(3);
|
---|
840 | line->SetBit(kCanDelete);
|
---|
841 | line->Draw();
|
---|
842 |
|
---|
843 | TLegend *leg = new TLegend(0.75,0.01,0.99,0.3);
|
---|
844 | leg->SetBit(kCanDelete);
|
---|
845 | leg->AddEntry(line, "Trigger Calibration", "l");
|
---|
846 | leg->Draw();
|
---|
847 | }
|
---|
848 |
|
---|
849 | Int_t MHCalibrationRelTimeCam::ReadEnv(const TEnv &env, TString prefix, Bool_t print)
|
---|
850 | {
|
---|
851 | Bool_t rc = kFALSE;
|
---|
852 |
|
---|
853 | if (IsEnvDefined(env, prefix, "ReferenceFile", print))
|
---|
854 | {
|
---|
855 | SetReferenceFile(GetEnvValue(env, prefix, "ReferenceFile", fReferenceFile.Data()));
|
---|
856 | rc = kTRUE;
|
---|
857 | }
|
---|
858 |
|
---|
859 | if (IsEnvDefined(env, prefix, "ReferencePixel", print))
|
---|
860 | {
|
---|
861 | SetReferencePixel(GetEnvValue(env, prefix, "ReferencePixel", (Int_t)fReferencePixel));
|
---|
862 | rc = kTRUE;
|
---|
863 | }
|
---|
864 |
|
---|
865 | TEnv refenv(fReferenceFile);
|
---|
866 |
|
---|
867 | fInnerRefTime = refenv.GetValue("InnerRefTime", fInnerRefTime);
|
---|
868 | fOuterRefTime = refenv.GetValue("OuterRefTime", fOuterRefTime);
|
---|
869 |
|
---|
870 | return MHCalibrationCam::ReadEnv(env,prefix,print) ? kTRUE : rc;
|
---|
871 | }
|
---|