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