Context Navigation

MJCalibrateSignal.cc@ 18752

Visit:

Last change on this file since 18752 was 9519, checked in by tbretz, 15 years ago
* empty log message *
File size: 35.6 KB

Line
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): Thomas Bretz, 1/2004 <mailto:tbretz@astro.uni-wuerzburg.de>
19	!
20	! Copyright: MAGIC Software Development, 2000-2008
21	!
22	!
23	\* ======================================================================== */
24
25	/////////////////////////////////////////////////////////////////////////////
26	//
27	// MJCalibrateSignal
28	//
29	// This class is reading the output written by callisto. It calibrates
30	// signal and time.
31	//
32	// The signal and time extractors are read from the callisto-output. In
33	// pricipal you could overwrite these default using the resource file,
34	// but this is NOT recommended!
35	//
36	/////////////////////////////////////////////////////////////////////////////
37	#include "MJCalibrateSignal.h"
38
39	#include <TEnv.h>
40	#include <TFile.h>
41
42	// Core
43	#include "MLog.h"
44	#include "MLogManip.h"
45
46	#include "MDirIter.h"
47	#include "MTaskList.h"
48	#include "MParList.h"
49	#include "MEvtLoop.h"
50
51	#include "MStatusDisplay.h"
52
53	// General containers
54	#include "MGeomCam.h"
55	#include "MBadPixelsCam.h"
56	#include "MCalibConstCam.h"
57	#include "MPedestalCam.h"
58	#include "MArrivalTimeCam.h"
59
60	// General histograms
61	#include "MH3.h"
62	#include "MHCamEvent.h"
63	#include "MHVsTime.h"
64
65	// Calibration containers
66	#include "MCalibrationQECam.h"
67	#include "MCalibrationBlindCam.h"
68	#include "MCalibrationChargeCam.h"
69	#include "MCalibrationRelTimeCam.h"
70	#include "MCalibrationChargePINDiode.h"
71
72	// Calibration histograms
73	#include "MHCalibrationChargeCam.h"
74	#include "MHCalibrationChargeBlindCam.h"
75	#include "MHCalibrationChargePINDiode.h"
76	#include "MHCalibrationRelTimeCam.h"
77
78	// Tasks
79	#include "MReadMarsFile.h"
80	#include "MRawFileRead.h"
81	#include "MTaskEnv.h"
82	#include "MContinue.h"
83	#include "MFillH.h"
84	#include "MGeomApply.h"
85	#include "MExtractTimeAndCharge.h"
86	#include "MTriggerPatternDecode.h"
87	#include "MCalibrationPatternDecode.h"
88	#include "MCalibrationChargeCalc.h"
89	#include "MCalibrationRelTimeCalc.h"
90	#include "MCalibCalcFromPast.h"
91	#include "MPedestalSubtract.h"
92	#include "MPedCalcFromLoGain.h"
93	#include "MCalibrateData.h"
94	#include "MExtractPINDiode.h"
95	#include "MExtractBlindPixel.h"
96	#include "MCalibrateRelTimes.h"
97	#include "MBadPixelsCalc.h"
98	#include "MBadPixelsTreat.h"
99	#include "MWriteRootFile.h"
100
101	// Filter
102	#include "MFTriggerPattern.h"
103	#include "MFCosmics.h"
104	#include "MFilterList.h"
105	#include "MFDataPhrase.h"
106
107	// Classes for writing movies
108	#include "MFEvtNumber.h"
109	#include "MMoviePrepare.h"
110	#include "MMovieWrite.h"
111	#include "MImgCleanStd.h"
112
113
114	ClassImp(MJCalibrateSignal);
115
116	using namespace std;
117
118	// --------------------------------------------------------------------------
119	//
120	// Default constructor.
121	//
122	// Sets:
123	// - fIsInterlaced to kTRUE
124	// - fIsRelTimesUpdate to kFALSE
125	// - fIsHiLoCalibration to kFALSE
126	//
127	MJCalibrateSignal::MJCalibrateSignal(const char name, const char title)
128	: fExtractor(0), fIsInterlaced(kTRUE), fIsRelTimesUpdate(kTRUE), fIsMovieMode(kFALSE), fIsTestMode(kFALSE)
129	{
130	fName = name ? name : "MJCalibrateSignal";
131	fTitle = title ? title : "Tool to calibrate data";
132	}
133
134	MJCalibrateSignal::~MJCalibrateSignal()
135	{
136	if (fExtractor)
137	delete fExtractor;
138	}
139
140	void MJCalibrateSignal::SetExtractor(const MExtractor *ext)
141	{
142	if (fExtractor)
143	delete fExtractor;
144
145	fExtractor = ext ? (MExtractor*)ext->Clone() : NULL;
146	}
147
148	Bool_t MJCalibrateSignal::WriteResult() const
149	{
150	if (IsNoStorage() \|\| fIsMovieMode)
151	return kTRUE;
152
153	// FIXME: This is not nice because it will update the signal
154	// file always. Which might make the usage outside of
155	// callisto difficult.
156	TObjArray cont;
157	cont.Add(const_cast<TEnv*>(GetEnv()));
158	cont.Add(const_cast<MSequence*>(&fSequence));
159
160	TNamed cmdline("CommandLine", fCommandLine.Data());
161	cont.Add(&cmdline);
162
163	if (fDisplay)
164	{
165	TString title = "-- Calibrate Signal: ";
166	title += fSequence.GetSequence();
167	title += " --";
168	fDisplay->SetTitle(title, kFALSE);
169
170	cont.Add(fDisplay);
171	}
172
173	const TString name(Form("signal%08d.root", fSequence.GetSequence()));
174	return WriteContainer(cont, name, "UPDATE");
175	}
176
177	Bool_t MJCalibrateSignal::ReadCalibration(TObjArray &l, MBadPixelsCam &cam, MExtractor* &ext2, MExtractor* &ext3, TString &geom) const
178	{
179	TString fname = Form("%s/calib%08d.root", fPathIn.Data(), fSequence.GetSequence());
180
181	*fLog << inf << "Reading from file: " << fname << endl;
182
183	TFile file(fname, "READ");
184	if (!file.IsOpen())
185	{
186	*fLog << err << dbginf << "ERROR - Could not open file " << fname << endl;
187	return kFALSE;
188	}
189
190	TObject *o = file.Get("ExtractSignal");
191	if (o && !o->InheritsFrom(MExtractor::Class()))
192	{
193	*fLog << err << dbginf << "ERROR - ExtractSignal read from " << fname << " doesn't inherit from MExtractor!" << endl;
194	return kFALSE;
195	}
196	ext3 = o ? (MExtractor*)o->Clone() : NULL;
197
198	o = file.Get("ExtractTime");
199	if (o && !o->InheritsFrom(MExtractor::Class()))
200	{
201	*fLog << err << dbginf << "ERROR - ExtractTime read from " << fname << " doesn't inherit from MExtractor!" << endl;
202	return kFALSE;
203	}
204	ext2 = o ? (MExtractor*)o->Clone() : NULL;
205	if (!ext3 && !ext2)
206	{
207	*fLog << err << dbginf << "ERROR - Neither ExtractSignal nor ExrtractTime found in " << fname << "!" << endl;
208	return kFALSE;
209	}
210
211	o = file.Get("MGeomCam");
212	if (o && !o->InheritsFrom(MGeomCam::Class()))
213	{
214	*fLog << err << dbginf << "ERROR - MGeomCam read from " << fname << " doesn't inherit from MGeomCam!" << endl;
215	return kFALSE;
216	}
217	geom = o ? o->ClassName() : "";
218
219	TObjArray cont(l);
220	cont.Add(&cam);
221	return ReadContainer(cont);
222	}
223
224	// --------------------------------------------------------------------------
225	//
226	// MJCalibration allows to setup several option by a resource file:
227	// MJCalibrateSignal.RawData: yes,no
228	//
229	// For more details see the class description and the corresponding Getters
230	//
231	Bool_t MJCalibrateSignal::CheckEnvLocal()
232	{
233	SetInterlaced(GetEnv("Interlaced", fIsInterlaced));
234	SetRelTimesUpdate(GetEnv("RelTimesUpdate", fIsRelTimesUpdate));
235	SetMovieMode(GetEnv("MovieMode", fIsMovieMode));
236
237	return MJCalib::CheckEnvLocal();
238	}
239
240	Bool_t MJCalibrateSignal::Process(MPedestalCam &pedcamab, MPedestalCam &pedcambias,
241	MPedestalCam &pedcamextr)
242	{
243	if (!fSequence.IsValid())
244	{
245	*fLog << err << "ERROR - Sequence invalid..." << endl;
246	return kFALSE;
247	}
248
249	// --------------------------------------------------------------------------------
250
251	*fLog << inf;
252	fLog->Separator(GetDescriptor());
253	*fLog << "Calculate calibrated data from Sequence #";
254	*fLog << fSequence.GetSequence() << endl << endl;
255
256
257	if (!CheckEnv())
258	return kFALSE;
259
260	// --------------------------------------------------------------------------------
261
262	MDirIter iter;
263	if (fSequence.IsValid())
264	{
265	if (fSequence.GetRuns(iter, fIsTestMode?MSequence::kRawCal:MSequence::kRawDat)<=0)
266	return kFALSE;
267	}
268
269	// Read File
270	MHCalibrationChargeCam hchacam;
271	MHCalibrationChargeBlindCam hbndcam;
272	MHCalibrationChargePINDiode hpndiod;
273	MHCalibrationRelTimeCam hrelcam;
274	//MHCalibrationHiLoCam hilocam;
275	//MHCalibrationPulseTimeCam hpulcam;
276
277	hchacam.SetOscillations(kFALSE);
278	hbndcam.SetOscillations(kFALSE);
279	hrelcam.SetOscillations(kFALSE);
280
281	MCalibrationChargeCam calcam;
282	MCalibrationQECam qecam;
283	MCalibrationBlindCam bndcam;
284	MCalibrationChargePINDiode pind;
285	MCalibrationRelTimeCam tmcam;
286	//MCalibrationHiLoCam hilcam;
287	//MCalibrationPulseTimeCam pulcam;
288
289	MBadPixelsCam badpix;
290
291	TObjArray interlacedcont;
292	if (fIsInterlaced)
293	{
294	interlacedcont.Add(&hchacam);
295	if (IsUseBlindPixel())
296	interlacedcont.Add(&hbndcam);
297	if (IsUsePINDiode())
298	interlacedcont.Add(&hpndiod);
299	if (fIsRelTimesUpdate)
300	interlacedcont.Add(&hrelcam);
301	}
302
303	MExtractor *extractor1=fExtractor;
304	MExtractor *extractor2=0;
305	MExtractor *extractor3=0;
306	TString geom;
307
308	TObjArray calibcont;
309	calibcont.Add(&calcam);
310	calibcont.Add(&qecam);
311	calibcont.Add(&bndcam);
312	calibcont.Add(&tmcam);
313	if (IsUseBlindPixel())
314	calibcont.Add(&bndcam);
315
316	if (!ReadCalibration(calibcont, badpix, extractor2, extractor3, geom))
317	return kFALSE;
318
319	*fLog << all;
320	if (!geom.IsNull())
321	*fLog << inf << "Camera geometry found in file: " << geom << endl;
322	else
323	*fLog << inf << "No Camera geometry found using default <MGeomCamMagic>" << endl;
324
325	if (extractor3)
326	{
327	*fLog << underline << "Signal Extractor found in calibration file" << endl;
328	extractor3->Print();
329	*fLog << endl;
330	}
331	else
332	*fLog << inf << "No Signal Extractor: ExtractSignal in file." << endl;
333
334
335	if (extractor1)
336	{
337	*fLog << underline << "Modified Signal Extractor set by user." << endl;
338	extractor1->Print();
339	*fLog << endl;
340	}
341	else
342	{
343	*fLog << inf << "No modified Signal Extractor set by user... using ExtractSignal." << endl;
344	extractor1 = extractor3 ? (MExtractor*)extractor3->Clone() : 0;
345	}
346
347	if (extractor2)
348	{
349	*fLog << underline << "Time Extractor found in calibration file" << endl;
350	extractor2->Print();
351	*fLog << endl;
352	}
353	else
354	*fLog << inf << "No Time Extractor: ExtractTime in file." << endl;
355
356	// This is necessary for the case in which it is not in the files
357	MCalibConstCam constcam;
358
359	//MBadPixelsCam badcam;
360	constcam.SetBadPixels(&badpix);
361
362	// Setup Parlist
363	MParList plist;
364	plist.AddToList(this); // take care of fDisplay!
365	plist.AddToList(&badpix);
366	plist.AddToList(&constcam);
367	//plist.AddToList(&hilcam);
368	plist.AddToList(&calibcont); // Using AddToList(TObjArray *)
369	plist.AddToList(&interlacedcont); // Using AddToList(TObjArray *)
370	//plist.AddToList(&pulcam);
371
372	// Setup Tasklist
373	MTaskList tlist;
374	plist.AddToList(&tlist);
375
376	MReadMarsFile readmc("Events");
377	readmc.DisableAutoScheme();
378
379	MRawFileRead rawread(NULL);
380	rawread.SetForceMode(); // Ignore broken time-stamps
381
382	MRead read = fSequence.IsMonteCarlo() ? (MRead)&readmc : (MRead*)&rawread;
383	read->AddFiles(iter);
384
385	const TString fname(Form("s/(([0-9]+_)?(M[12]_)?[0-9.]+)_[CD]_(.*[.])(raw\|raw[.]gz\|root)$/%s\\/$1_Y_$4root/",
386	Esc(fPathOut).Data()));
387
388	// Skips MC which have no contents. This are precisely the
389	// events which fullfilled the MC Lvl1 trigger and an
390	// arbitrary cut (typically at 50phe) to speed up simulation
391	MContinue contmc("MRawEvtData.GetNumPixels<0.5", "ContEmptyMC");
392
393	//MPointingPosInterpolate pextr;
394	//pextr.AddFiles(&iter);
395
396	MGeomApply apply; // Only necessary to create geometry
397	if (!geom.IsNull())
398	apply.SetGeometry(geom);
399	//MBadPixelsMerge merge(&badpix);
400
401	// Make sure that pedcamab has the correct name
402	pedcamab.SetName("MPedestalFundamental");
403	pedcamextr.SetName("MPedestalFromExtractorRndm");
404	pedcambias.SetName("MPedestalFromExtractor");
405	plist.AddToList(&pedcamextr);
406	plist.AddToList(&pedcambias);
407	plist.AddToList(&pedcamab);
408
409	MArrivalTimeCam timecam;
410	plist.AddToList(&timecam);
411
412	// Check for interleaved events
413	MCalibrationPatternDecode caldec;
414	MTriggerPatternDecode decode;
415
416	MH3 hpat("MRawRunHeader.GetFileID", "MTriggerPattern.GetUnprescaled");
417	if (fSequence.IsMonteCarlo())
418	hpat.SetWeight("100./MMcRunHeader.fNumSimulatedShowers");
419	else
420	hpat.SetWeight("1./TMath::Max(MRawRunHeader.GetRunLength,1)");
421	hpat.SetName("TrigPat");
422	hpat.SetTitle("Rate of the trigger pattern [Hz];File Id;Trigger Pattern;Rate [Hz]");
423	hpat.InitLabels(MH3::kLabelsXY);
424	//hpat.DefineLabelsY("1=Lvl1;2=Cal;3=Cal;4=Lvl2;5=Cal;7=Cal;8=Ped;9=Ped+Trig;13=Ped+Trig;16=Pin;32=Sum");
425	hpat.DefaultLabelY("UNKNOWN");
426	hpat.DefineLabelY( 0, "0"); // 0: No pattern
427	hpat.DefineLabelY( 1, "Trig"); // 1: Lvl1
428	hpat.DefineLabelY( 2, "Cal"); // 2: Cal
429	hpat.DefineLabelY( 3, "Cal"); // Cal+Lvl1
430	hpat.DefineLabelY( 4, "Trig"); // 4: Lvl2
431	hpat.DefineLabelY( 5, "Trig"); // Lvl2+Lvl1
432	hpat.DefineLabelY( 7, "Cal"); // Lvl2+Cal+Lvl
433	hpat.DefineLabelY( 8, "Ped"); // 8: Ped
434	hpat.DefineLabelY( 9, "Ped+Trig"); // Ped+Lvl1
435	hpat.DefineLabelY(10, "Ped+Cal"); // Ped+Cal
436	hpat.DefineLabelY(12, "Ped+Trig"); // Ped+Lvl2
437	hpat.DefineLabelY(13, "Ped+Trig"); // Ped+Lvl2+Lvl1
438	hpat.DefineLabelY(16, "Pin"); // 16: Pin
439	hpat.DefineLabelY(32, "Sum"); // 32: Sum
440	hpat.DefineLabelY(33, "Trig"); // Sum+Lvl1
441	hpat.DefineLabelY(34, "Cal"); // Sum+Cal
442	hpat.DefineLabelY(35, "Cal"); // Sum+Cal+Lvl1
443	hpat.DefineLabelY(36, "Trig"); // Sum+Lvl2
444	hpat.DefineLabelY(37, "Trig"); // Sum+Lvl1+Lvl2
445	hpat.DefineLabelY(39, "Cal"); // Sum+Lvl2+Cal+Lvl1
446
447	MFillH fillpat(&hpat, "", "FillPattern");
448	fillpat.SetDrawOption("box");
449
450	// This will make that for data with version less than 5, where
451	// trigger patterns were not yet correct, all the events in the real
452	// data file will be processed. In any case there are no interleaved
453	// calibration events in such data, so this is fine.
454	// We allow only cosmics triggered events to pass (before prescaling)
455	MFTriggerPattern fcalped("SelectTrigEvts");
456	fcalped.SetInverted();
457	fcalped.SetDefault(kTRUE);
458	fcalped.DenyAll();
459	fcalped.AllowTriggerLvl1();
460	fcalped.AllowTriggerLvl2();
461	fcalped.AllowSumTrigger();
462	if (fIsTestMode)
463	fcalped.AllowCalibration();
464
465	// This will skip interleaved events with a cal- or ped-trigger
466	MContinue contcalped(&fcalped, "ContNonTrigger");
467
468	// Create the pedestal subtracted raw-data
469	MPedestalSubtract pedsub;
470	pedsub.SetPedestalCam(&pedcamab);
471
472	// Do signal and pedestal calculation
473	MPedCalcFromLoGain pedlo1("MPedCalcFundamental");
474	pedlo1.SetPedestalUpdate(kTRUE);
475	pedlo1.SetNamePedestalCamOut("MPedestalFundamental");
476
477	MPedCalcFromLoGain pedlo2("MPedCalcWithExtractorRndm");
478	pedlo2.SetPedestalUpdate(kTRUE);
479	pedlo2.SetRandomCalculation(kTRUE);
480	pedlo2.SetNamePedestalCamOut("MPedestalFromExtractorRndm");
481
482	MPedCalcFromLoGain pedlo3("MPedCalcWithExtractor");
483	pedlo3.SetPedestalUpdate(kTRUE);
484	pedlo3.SetRandomCalculation(kFALSE);
485	pedlo3.SetNamePedestalCamOut("MPedestalFromExtractor");
486
487	if (!extractor1)
488	{
489	*fLog << err << "ERROR - extractor1 == NULL" << endl;
490	return kFALSE;
491	}
492
493	// Setup to use the hi-gain extraction window in the lo-gain
494	// range (the start of the lo-gain range is added automatically
495	// by MPedCalcFromLoGain)
496	//
497	// The window size of the extractor is not yet initialized,
498	// so we have to stick to the extraction range
499	//
500	// Even if we would like to use a range comparable to the
501	// hi-gain extraction we use the lo-gain range to make
502	// sure that exclusions (eg. due to switching noise)
503	// are correctly handled.
504	//
505	pedlo1.SetRangeFromExtractor(*extractor1);
506
507	if (extractor1->InheritsFrom("MExtractTimeAndCharge"))
508	{
509	pedlo2.SetExtractor((MExtractTimeAndCharge*)extractor1);
510	pedlo3.SetExtractor((MExtractTimeAndCharge*)extractor1);
511	}
512	else
513	{
514	pedlo2.SetRangeFromExtractor(*extractor1);
515	pedlo3.SetRangeFromExtractor(*extractor1);
516	}
517
518	//------------------------------
519	//
520	// this is the filter to find pedestal events. For the Siegen FADC
521	// these are all events which are not calibration events because
522	// the pedestal is extracted from the lo-gain signal. For MUX
523	// data this are artifiially triggered events with the pedestal
524	// trigger flag, and for safty without Lvl1 or Lvl2 flag)
525	//
526	// For the time "before" the trigger pattern all events (Siegen FADC)
527	// can be considered to be pedestal, because this was also the time
528	// without artifially calibration events
529	//
530	// Deny or allow is done before prescaling.
531	//
532	MFTriggerPattern fped("SelectPedestals");
533	fped.SetDefault(kTRUE);
534	fped.DenyCalibration();
535	if (!extractor1->HasLoGain())
536	{
537	fped.DenyAll();
538	fped.RequirePedestal();
539	}
540
541	//------------------------------
542	//
543	// Apply a filter against cosmics (this is to make sure that the
544	// trigger system was working properly and no empty events survive)
545	// For every event 5% of the pixel must not be empty. In PostProcess
546	// an error is raised if more than 50% of the events were skipped.
547	//
548	// Since unsuitable pixels are not counted in the ideal case
549	// there is no empty pixel at all.
550	//
551	MFCosmics fcosmicscal;
552	fcosmicscal.SetNamePedestalCam("MPedestalFundamental"); //CORRECT?
553	fcosmicscal.SetMaxEmptyPixels(0.05);
554	fcosmicscal.SetMaxAcceptedFraction(0.5);
555
556	MContinue contcoscal(&fcosmicscal, "ContCosmicsCal");
557
558	//------------------------------
559	//
560	// Remove events which are too bright. These events could either
561	// be calibration events or events which produce many many
562	// saturating low-gains and thus spoil the dead pixel plot.
563	// These events don't seem to be physical events so we can fairly
564	// remove them (they are no analysable gammas anyway)
565	// So we remove all events which have less than 5% empty pixels.
566	//
567	MFCosmics fcosmicsbright;
568	fcosmicsbright.SetNamePedestalCam("MPedestalFundamental");
569	fcosmicsbright.SetMaxEmptyPixels(0.10);
570	fcosmicsbright.SetMinAcceptedFraction(0.9);
571
572	MContinue contbright(&fcosmicsbright, "ContBrightEvts");
573	contbright.SetInverted();
574
575	//------------------------------
576	//
577	// Thie signal extractors
578	//
579	MTaskEnv taskenv1("ExtractSignal");
580	MTaskEnv taskenv2("ExtractTime");
581	MTaskEnv taskenv3("ExtractInterlaced");
582	taskenv1.SetDefault(extractor1);
583	taskenv2.SetDefault(extractor2);
584	taskenv3.SetDefault(extractor3);
585
586	//
587	// This is new calibration to photo-electrons, hard-coded
588	// as decided at the Wuerzburg software meeting 26.01.05
589	//
590	MCalibrateData calib;
591	calib.SetSignalType(MCalibrateData::kPhe);
592	//calib.AddPedestal("Fundamental");
593	calib.AddPedestal("FromExtractor");
594	calib.AddPedestal("FromExtractorRndm");
595	calib.SetPedestalFlag(MCalibrateData::kEvent);
596
597	//----------------------------------------------------------
598
599	MExtractPINDiode pinext;
600	MExtractBlindPixel bldext;
601
602	// Execute for all events with the calibration trigger. If no
603	// trigger pattern is available do not execute it
604	// The selection is done before prescaling.
605	MFTriggerPattern fcalib("SelectCalEvts");
606	fcalib.SetDefault(kFALSE);
607	fcalib.DenyAll();
608	fcalib.RequireCalibration();
609	fcalib.AllowTriggerLvl1();
610	fcalib.AllowTriggerLvl2();
611	fcalib.AllowSumTrigger();
612
613	MCalibrationChargeCalc chcalc;
614	chcalc.SetContinousCalibration();
615	chcalc.SetExtractor(extractor3);
616
617	MCalibrationRelTimeCalc recalc;
618	MCalibCalcFromPast pacalc;
619
620	chcalc.SetPedestals(&pedcamextr);
621
622	pacalc.SetChargeCalc(&chcalc);
623	if (fIsRelTimesUpdate)
624	pacalc.SetRelTimeCalc(&recalc);
625	pacalc.SetCalibrate(&calib);
626
627	//
628	// Calibration histogramming
629	//
630	MFillH filpin(&hpndiod, "MExtractedSignalPINDiode", "FillPINDiode");
631	MFillH filbnd(&hbndcam, "MExtractedSignalBlindPixel", "FillBlindCam");
632	MFillH filcam(&hchacam, "MExtractedSignalCam", "FillChargeCam");
633	MFillH filtme(&hrelcam, "MArrivalTimeCam", "FillRelTime");
634	//MFillH filhil(&hilocam, "MExtractedSignalCam", "FillHiLoRatio");
635	//MFillH filpul(&hpulcam, "MRawEvtData", "FillPulseTime");
636	filpin.SetBit(MFillH::kDoNotDisplay);
637	filbnd.SetBit(MFillH::kDoNotDisplay);
638	filcam.SetBit(MFillH::kDoNotDisplay);
639	filtme.SetBit(MFillH::kDoNotDisplay);
640	//filhil.SetBit(MFillH::kDoNotDisplay);
641	//filpul.SetBit(MFillH::kDoNotDisplay);
642
643	MCalibrateRelTimes caltm;
644	MBadPixelsCalc bpcal;
645	MBadPixelsTreat treat;
646
647	//bpcal.SetNamePedPhotCam("MPedPhotFromExtractor");
648	bpcal.SetNamePedPhotCam("MPedPhotFromExtractorRndm");
649
650	//treat.AddNamePedPhotCam("MPedPhotFundamental");
651	treat.AddNamePedPhotCam("MPedPhotFromExtractor");
652	treat.AddNamePedPhotCam("MPedPhotFromExtractorRndm");
653	if (!extractor2 && !extractor1->InheritsFrom("MExtractTimeAndCharge"))
654	treat.SetProcessTimes(kFALSE);
655
656	MHCamEvent evt0( 0, "PedFLG", "Fundamental Pedestal from Lo Gain;;P [cnts/sl]");
657	MHCamEvent evt1( 2, "PedRmsFLG", "RMS from Extractor applied to ped.;;\\sigma_{p} [cnts/sl]");
658	MHCamEvent evt2( 0, "Extra'd", "Extracted Signal;;S [cnts/sl]");
659	// MHCamEvent evt3(4, "PedPhot", "Calibrated Pedestal;;P [phe]");
660	MHCamEvent evt4( 5, "PedRMS", "Calibrated RMS from Extractor applied to ped.;;\\sigma_{p} [phe]");
661	MHCamEvent evt5( 0, "Interp'd", "Interpolated Signal scaled with A/A_{0};;S [phe]");
662	MHCamEvent evt6(102, "Unsuitable", "Fraction of unsuitable events per Pixel;;[1]");
663	// MHCamEvent evt7( 6, "Times", "Calibrated Arrival Time;;T [fadc sl]");
664	MHCamEvent evt8( 0, "Conv", "Calibration Conv. Factors;;[phe/cnts]");
665	MHCamEvent evt9( 7, "PulsePos", "Pulse Position of cosmics (>50phe);;T [ns]");
666	MHCamEvent evtR( 4, "HiLoCal", "Hi-/Lo-Gain ratio;;Ratio");
667	MHCamEvent evtO( 7, "HiLoOff", "Lo-/Hi-Gain Offset;;Offset");
668	MHCamEvent evtC( 4, "CalPos", "Extracted pulse position of calibration pulses;;T [sl]");
669	//MHCamEvent evt2(0, "Extra'd", "Extracted Calibration Signal;;S [cnts/sl]");
670	evt2.SetErrorSpread(kFALSE);
671	evt5.SetErrorSpread(kFALSE);
672	evt6.SetErrorSpread(kFALSE);
673	evt6.SetThreshold();
674
675	MFillH fill0(&evt0, "MPedestalFundamental", "FillPedFLG");
676	MFillH fill1(&evt1, "MPedestalFromExtractorRndm", "FillPedRmsFLG");
677	MFillH fill2(&evt2, "MExtractedSignalCam", "FillExtracted");
678	// MFillH fill3(&evt3, "MPedPhotFundamental", "FillPedPhot");
679	MFillH fill4(&evt4, "MPedPhotFromExtractorRndm", "FillPedRMS");
680	MFillH fill5(&evt5, "MSignalCam", "FillInterpolated");
681	MFillH fill6(&evt6, "MBadPixelsCam", "FillUnsuitable");
682	// MFillH fill7(&evt7, "MSignalCam", "FillTimes");
683	MFillH fill8(&evt8, "MCalibConstCam", "FillConv");
684	MFillH fill9(&evt9, "MSignalCam", "FillPulsePos");
685	MFillH fillR(&evtR, "MExtractedSignalCam", "FillHiLoCal");
686	MFillH fillO(&evtO, "MArrivalTimeCam", "FillHiLoOff");
687	MFillH fillC(&evtC, "MArrivalTimeCam", "FillCalPos");
688	//MFillH fill2(&evt2, "MExtractedSignalCam", "FillExtractedSignal");
689
690	MHVsTime histbp("MBadPixelsCam.GetNumUnsuitable");
691	histbp.SetName("BadPixTm");
692	histbp.SetTitle("Number of unsuitable pixels;;N");
693	histbp.SetMinimum(0);
694
695	MHVsTime histdp("MSignalCam.GetNumPixelsUnmapped");
696	histdp.SetName("DeadPixTm");
697	histdp.SetTitle("Number of dead/unmapped pixels;;N");
698	histdp.SetMinimum(0);
699
700	// Task to fill the histogram
701	MFillH fillB(&histbp, "MTime", "FillBadPixTm");
702	MFillH fillD(&histdp, "MTime", "FillDeadPixTm");
703	fillB.SetNameTab("BadPixTm");
704	fillD.SetNameTab("DeadPixTm");
705
706	/*
707	MFillH fillP("MHPulseShape", "", "FillPulseShape");
708	fillP.SetNameTab("Pulse");
709	*/
710
711	/*
712	MHVsTime hbadpix("MBadPixelsCam.GetNumUnsuitable");
713	hbadpix.SetNumEvents(50);
714	MFillH fillB(&hbadpix, "MTime");
715	*/
716
717	MTaskEnv fillflorian("FinalFantasy");
718	fillflorian.SetDefault();
719
720	// The second rule is for the case reading raw-files!
721	MWriteRootFile write(2, fname, fOverwrite?"RECREATE":"NEW", "Calibrated Data");
722	// Run Header
723	write.AddContainer("MRawRunHeader", "RunHeaders");
724	// write.AddContainer("MBadPixelsCam", "RunHeaders");
725	write.AddContainer("MGeomCam", "RunHeaders");
726	// Monte Carlo Headers
727	write.AddContainer("MMcRunHeader", "RunHeaders", kFALSE);
728	write.AddContainer("MMcFadcHeader", "RunHeaders", kFALSE);
729	write.AddContainer("MMcTrigHeader", "RunHeaders", kFALSE);
730	write.AddContainer("MMcConfigRunHeader", "RunHeaders", kFALSE);
731	write.AddContainer("MMcCorsikaRunHeader", "RunHeaders", kFALSE);
732	write.AddContainer("MCorsikaRunHeader", "RunHeaders", kFALSE);
733	// Monte Carlo
734	write.AddContainer("MMcEvt", "Events", kFALSE);
735	write.AddContainer("MMcTrig", "Events", kFALSE);
736	write.AddContainer("MCorsikaEvtHeader", "Events", kFALSE);
737	// Data tree
738	write.AddContainer("MSignalCam", "Events");
739	// write.AddContainer("MPedPhotFundamental", "Events");
740	write.AddContainer("MPedPhotFromExtractor", "Events");
741	write.AddContainer("MPedPhotFromExtractorRndm", "Events");
742	write.AddContainer("MTime", "Events", kFALSE);
743	write.AddContainer("MRawEvtHeader", "Events");
744	write.AddContainer("MTriggerPattern", "Events");
745
746	// Trees with slow-control information (obsolete
747	// if we don't read merpped root-files)
748	write.AddTree("Trigger", kFALSE);
749	write.AddTree("Drive", kFALSE);
750	write.AddTree("CC", kFALSE);
751	write.AddTree("Pyrometer", kFALSE);
752	write.AddTree("Currents", kFALSE);
753	write.AddTree("Camera", kFALSE);
754	// Slow-Control: Current-tree
755	write.AddContainer("MCameraDC", "Currents", kFALSE);
756	// Slow-Control: Camera-tree
757	write.AddContainer("MCameraAUX", "Camera", kFALSE);
758	write.AddContainer("MCameraCalibration", "Camera", kFALSE);
759	write.AddContainer("MCameraCooling", "Camera", kFALSE);
760	write.AddContainer("MCameraHV", "Camera", kFALSE);
761	write.AddContainer("MCameraLV", "Camera", kFALSE);
762	write.AddContainer("MCameraLids", "Camera", kFALSE);
763
764	// Write the special MC tree
765	MWriteRootFile writemc(2, fname, fOverwrite?"RECREATE":"NEW", "Calibrated Data");
766	writemc.SetName("WriteMC");
767	writemc.AddContainer("MMcEvtBasic", "OriginalMC", kFALSE);
768	if (fSequence.IsMonteCarlo())
769	writemc.AddCopySource("OriginalMC", kFALSE);
770
771	// Write the special calib tree
772	/*
773	MWriteRootFile writecal(2, fname, fOverwrite?"RECREATE":"NEW");
774	writecal.SetName("WriteCalib");
775	writecal.AddContainer("MBadPixelsCam", "Calib");
776	writecal.AddContainer("MCalibrationChargeCam", "Calib");
777	writecal.AddContainer("MCalibrationRelTimeCam", "Calib");
778	*/
779
780	//-----------------------------------------------------------
781	// Build tasklist
782
783	MTaskList tlist2("AllEvents");
784
785	tlist2.AddToList(&caldec);
786	tlist2.AddToList(&decode);
787	tlist2.AddToList(&fillpat);
788	tlist2.AddToList(&apply);
789	//tlist2.AddToList(&merge);
790	tlist2.AddToList(&pedsub);
791
792	//-----------------------------------------------------------
793	// Pedestal extraction
794
795	MTaskList tlist3("PedEvents");
796	tlist3.SetFilter(&fped); // Deny events with cal-trigger
797
798	tlist2.AddToList(&fped); // If no lo-gain require ped-trigger
799	tlist2.AddToList(&tlist3); // and deny cosmics (lvl1/2) trigger
800
801	tlist3.AddToList(&pedlo1); // extract pedestal events
802	tlist3.AddToList(&pedlo2); // extract pedestal events
803	tlist3.AddToList(&pedlo3); // extract pedestal events
804	tlist3.AddToList(&fill0); // fill pedestal events
805	tlist3.AddToList(&fill1); // fill pedestal events
806
807	//-----------------------------------------------------------
808	// Calibration
809
810	MTaskList tlist4("CalEvents");
811	tlist4.SetFilter(&fcalib); // process only events with cal-trigger
812
813	//MFDataPhrase filcalco("MCalibrationConstCam.IsReadyToSave>0.5", "CalibConstFilter");
814	if (fIsInterlaced && !fIsTestMode)
815	{
816	// The task list is executed for all events with the calibration
817	// trigger
818	tlist2.AddToList(&fcalib); // MFTriggerPattern
819	tlist2.AddToList(&tlist4);
820
821	tlist4.AddToList(&taskenv3);
822	tlist4.AddToList(&contcoscal); // MContinue/ContCosmicsCal
823	if (IsUsePINDiode())
824	tlist4.AddToList(&pinext); // MExtractPINDiode
825	if (IsUseBlindPixel())
826	tlist4.AddToList(&bldext); // MExtractBlindPixel
827	tlist4.AddToList(&pacalc); // MCalibCalcFromPast
828	/*
829	tlist3.AddToList(&filcalco); // CalibConstFilter (IsReadyToSave)
830	fill8.SetFilter(&filcalco);
831	tlist3.AddToList(&fill8); // FillConvUpd
832	*/
833
834	tlist4.AddToList(&filcam); // FillChargeCam
835	if (fIsRelTimesUpdate)
836	tlist4.AddToList(&filtme); // FillRelTime
837	if (IsUseBlindPixel())
838	tlist4.AddToList(&filbnd); // FillBlindCam
839	if (IsUsePINDiode())
840	tlist4.AddToList(&filpin); // FillPINDiode
841	tlist4.AddToList(&chcalc); // MCalibrationChargeCalc
842	if (fIsRelTimesUpdate)
843	tlist4.AddToList(&recalc); // MCalibrationRelTimeCam
844
845	tlist4.AddToList(&fillC); // FillCalPos
846
847	//tlist3.AddToList(&writecal); // MWriteRootFile
848	}
849
850	//-----------------------------------------------------------
851	// Cosmics extraction
852
853	// remove all events with a cal- or ped-trigger (no matter
854	// whether they have lvl1 or lvl2 or any other flag
855	tlist2.AddToList(&contcalped);
856
857	// Extract the signal
858	if (extractor1)
859	tlist2.AddToList(&taskenv1);
860
861	// remove all events which definitly don't have a signal
862	// using MFCosmics (ContCosmicsPed)
863	// tlist2.AddToList(&contcosped);
864
865	// Extract arrival time (if a dedicated extrator given)
866	if (extractor2)
867	tlist2.AddToList(&taskenv2);
868
869	// Check if we have an extremely bright event (remove them,
870	// they are presumably errornous events or calibration events)
871	if (!fIsTestMode)
872	tlist2.AddToList(&contbright);
873
874	/*
875	if (fIsHiLoCalibration)
876	{
877	plist.AddToList(&hilocam);
878	tlist2.AddToList(&filhil);
879	}
880
881	if (fIsPulsePosCheck)
882	{
883	plist.AddToList(&hpulcam);
884	tlist2.AddToList(&filpul);
885	}
886	*/
887
888	tlist2.AddToList(&fill2);
889	tlist2.AddToList(&fill8); // FillConv
890	tlist2.AddToList(&calib); // MCalibrateData
891	if (extractor2 \|\| extractor1->InheritsFrom("MExtractTimeAndCharge"))
892	tlist2.AddToList(&caltm);
893
894	tlist2.AddToList(&bpcal); // MBadPixelsCalc
895	tlist2.AddToList(&treat); // MBadPixelsTreat
896	tlist2.AddToList(&fill6);
897	// tlist2.AddToList(&fill3);
898	tlist2.AddToList(&fill4);
899	tlist2.AddToList(&fill5);
900	//if (extractor2 \|\| extractor1->InheritsFrom("MExtractTimeAndCharge"))
901	// tlist2.AddToList(&fill7);
902	tlist2.AddToList(&fill9);
903	if (!fSequence.IsMonteCarlo())
904	{
905	tlist2.AddToList(&fillB);
906	tlist2.AddToList(&fillD);
907	}
908	if (extractor1->HasLoGain())
909	{
910	tlist2.AddToList(&fillR);
911	tlist2.AddToList(&fillO);
912	}
913
914	/*
915	MFillH fillC("MHCleaning", "", "FillClean");
916	tlist2.AddToList(&fillC);
917
918	//tlist2.AddToList(&fillP);
919	*/
920
921	// ----- Start: Code for encoding movies -----
922
923	MMoviePrepare movprep;
924	MMovieWrite movwrite;
925	movprep.SetRangeFromExtractor(*extractor1);
926
927	//MFDataPhrase movfilt("MMovieData.fMax>150");
928	MFDataPhrase movfilt("MMovieData.fMax>5*MMovieData.fMedianPedestalRms", "MovieFilter");
929
930	MImgCleanStd movclean(8.5, 4.0);
931	movclean.SetMethod(MImgCleanStd::kAbsolute);
932
933	//movprep.SetFilter(&evtnum);
934	movclean.SetFilter(&movfilt);
935	movwrite.SetFilter(&movfilt);
936
937	MTaskList tlistmov("MovieEncoder");
938	tlistmov.AddToList(&movprep);
939	tlistmov.AddToList(&movfilt);
940	tlistmov.AddToList(&movclean);
941	tlistmov.AddToList(&movwrite);
942
943	MFEvtNumber evtnum;
944	//evtnum.SetSelector("ThetaSquared.fVal<0.04");
945	//evtnum.SetFileName("ganymed00000001.root");
946	tlistmov.SetFilter(&evtnum);
947
948	if (fIsMovieMode)
949	{
950	tlist2.AddToList(&evtnum);
951	tlist2.AddToList(&tlistmov);
952	}
953
954	// ----- End: Code for encoding movies -----
955
956	tlist2.AddToList(&fillflorian);
957
958	// Setup List for Drive-tree
959	//MPointingPosCalc pcalc;
960
961	// Now setup main tasklist
962	tlist.AddToList(read);
963
964	if (fSequence.IsMonteCarlo())
965	{
966	if (!fIsMovieMode && !HasNullOut())
967	tlist.AddToList(&writemc);
968	tlist.AddToList(&contmc);
969	}
970
971	//if (IsUseRootData())
972	// tlist2.AddToList(&pextr);
973	tlist.AddToList(&tlist2, fSequence.IsMonteCarlo() ? "Events" : "All");
974
975	//if (fSequence.IsMonteCarlo())
976	// tlist.AddToList(&pcalc, "Drive");
977
978	if (!fIsMovieMode && !HasNullOut())
979	tlist.AddToList(&write);
980
981	// Create and setup the eventloop
982	MEvtLoop evtloop(fName);
983	evtloop.SetParList(&plist);
984	evtloop.SetDisplay(fDisplay);
985	evtloop.SetLogStream(fLog);
986	if (!SetupEnv(evtloop))
987	return kFALSE;
988
989	// Execute first analysis
990	const Bool_t rc = evtloop.Eventloop(fMaxEvents);
991
992	// make sure owned object are deleted
993	if (extractor1 && extractor1!=fExtractor)
994	delete extractor1;
995	if (extractor2)
996	delete extractor2;
997	if (extractor3)
998	delete extractor3;
999
1000	// return if job failed
1001	if (!rc)
1002	{
1003	*fLog << err << GetDescriptor() << ": Failed." << endl;
1004	return kFALSE;
1005	}
1006
1007	// if everything went ok write and display result
1008	DisplayResult(plist);
1009
1010	/*
1011	if (fIsPixelCheck)
1012	{
1013	if (fIsPulsePosCheck)
1014	hpulcam[fCheckedPixId].DrawClone("");
1015
1016	//if (fIsHiLoCalibration)
1017	// hilocam[fCheckedPixId].DrawClone("");
1018	}
1019	interlacedcont.Add(&pulcam);
1020	*/
1021
1022	//if (fIsHiLoCalibration)
1023	// interlacedcont.Add(&hilcam);
1024
1025	//if (fIsPulsePosCheck)
1026	// interlacedcont.Add(plist.FindObject("MHCalibrationPulseTimeCam"));
1027
1028	//if (fIsHiLoCalibration)
1029	// interlacedcont.Add(plist.FindObject("MHCalibrationHiLoCam"));
1030
1031	if (!WriteResult())
1032	return kFALSE;
1033
1034	// return if job went ok
1035	*fLog << all << GetDescriptor() << ": Done." << endl;
1036	*fLog << endl << endl;
1037
1038	return kTRUE;
1039	}
1040
1041
1042	void MJCalibrateSignal::DisplayResult(MParList &plist)
1043	{
1044	/*
1045	if (!fDisplay \|\| !fIsHiLoCalibration)
1046	return;
1047
1048	MCalibrationHiLoCam hcam = (MCalibrationHiLoCam)plist.FindObject("MCalibrationHiLoCam");
1049	MGeomCam geom = (MGeomCam)plist.FindObject("MGeomCam");
1050	if (!hcam \|\| !geom)
1051	return;
1052
1053	// Create histograms to display
1054	MHCamera disp1(*geom, "HiLoConv", "Ratio Amplification HiGain vs. LoGain (Charges)");
1055	MHCamera disp2(*geom, "HiLoDiff", "Arrival Time Diff. HiGain vs. LoGain (Times)");
1056
1057	disp1.SetCamContent(*hcam, 0);
1058	disp1.SetCamError( *hcam, 1);
1059	disp2.SetCamContent(*hcam, 5);
1060	disp2.SetCamError( *hcam, 6);
1061
1062	disp1.SetYTitle("R [1]");
1063	disp2.SetYTitle("\\Delta T [FADC sl.]");
1064
1065	TCanvas &c1 = fDisplay->AddTab("HiLoConv");
1066	c1.Divide(2,3);
1067
1068	disp1.CamDraw(c1, 1, 2, 1);
1069	disp2.CamDraw(c1, 2, 2, 1);
1070	*/
1071	}
1072

Note: See TracBrowser for help on using the repository browser.

Context Navigation

source: branches/Corsika7405Compatibility/mjobs/MJCalibrateSignal.cc@ 18752

Download in other formats: