Context Navigation

source: trunk/MagicSoft/Mars/mjobs/MJCut.cc@ 7333

Visit:

Last change on this file since 7333 was 7287, checked in by tbretz, 19 years ago
* empty log message *
File size: 24.9 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/2005 <mailto:tbretz@astro.uni-wuerzburg.de>
19	!
20	! Copyright: MAGIC Software Development, 2000-2005
21	!
22	!
23	\* ======================================================================== */
24
25	/////////////////////////////////////////////////////////////////////////////
26	//
27	// MJCut
28	//
29	// FIXME: Preparation for wobble mode missing
30	//
31	/////////////////////////////////////////////////////////////////////////////
32	#include "MJCut.h"
33
34	// Root
35	#include <TEnv.h>
36	#include <TFile.h>
37
38	// Environment
39	#include "MLog.h"
40	#include "MLogManip.h"
41
42	// Eventloop
43	#include "MParList.h"
44	#include "MTaskList.h"
45	#include "MEvtLoop.h"
46
47	// Display
48	#include "MStatusDisplay.h"
49
50	// Tasks
51	#include "MReadReports.h"
52	#include "MReadMarsFile.h"
53	#include "MPrint.h"
54	#include "MContinue.h"
55	#include "MTaskEnv.h"
56	#include "MPointingDevCalc.h"
57	#include "MSrcPosCalc.h"
58	#include "MSrcPosCorrect.h"
59	#include "MHillasCalc.h"
60	#include "MFillH.h"
61	#include "MWriteRootFile.h"
62
63	// Fit signal environment
64	#include "../mhflux/MAlphaFitter.h"
65	#include "../mhflux/MHAlpha.h"
66
67	// Containers
68	#include "MH3.h"
69	#include "MBinning.h"
70	#include "MDataSet.h"
71	#include "MParameters.h"
72	#include "MPointingPos.h"
73	#include "MObservatory.h"
74
75	ClassImp(MJCut);
76
77	using namespace std;
78
79	// --------------------------------------------------------------------------
80	//
81	// Default constructor. Set defaults for fStoreSummary, fStoreresult,
82	// fWriteOnly, fIsWobble and fFullDisplay to kFALSE and initialize
83	// /fEstimateEnergy and/ fCalcHadronness with NULL.
84	//
85	MJCut::MJCut(const char name, const char title)
86	: fStoreSummary(kFALSE), fStoreResult(kTRUE), fWriteOnly(kFALSE),
87	fIsWobble(kFALSE), fIsMonteCarlo(kFALSE), fFullDisplay(kTRUE),
88	fNameHist("MHThetaSq"), fCalcHadronness(0)
89	{
90	fName = name ? name : "MJCut";
91	fTitle = title ? title : "Standard program to perform g/h-separation cuts";
92	}
93
94	// --------------------------------------------------------------------------
95	//
96	// Destructor. Delete fEstimateEnergy and fCalcHadronness if != NULL
97	//
98	MJCut::~MJCut()
99	{
100	//if (fEstimateEnergy)
101	// delete fEstimateEnergy;
102	if (fCalcHadronness)
103	delete fCalcHadronness;
104	}
105
106	// --------------------------------------------------------------------------
107	//
108	// Set the name of the summary file (events after cut0)
109	// If you give a name the storage of this file is enabled implicitly.
110	// If you give no filename the storage is neither enabled nor disabled,
111	// but the storage file name is reset.
112	// If no filename is set the default filename is used.
113	// You can explicitly enable or disable the storage using EnableStoreOf*()
114	// The default argument is no filename.
115	//
116	void MJCut::SetNameSummaryFile(const char *name)
117	{
118	fNameSummary=name;
119	if (!fNameSummary.IsNull())
120	EnableStorageOfSummary();
121	}
122
123	// --------------------------------------------------------------------------
124	//
125	// Set the name of the summary file (events after cut3)
126	// If you give a name the storage of this file is enabled implicitly.
127	// If you give no filename the storage is neither enabled nor disabled,
128	// but the storage file name is reset.
129	// If no filename is set the default filename is used.
130	// You can explicitly enable or disable the storage using EnableStoreOf*()
131	// The default argument is no filename.
132	//
133	void MJCut::SetNameResultFile(const char *name)
134	{
135	fNameResult=name;
136	if (!fNameResult.IsNull())
137	EnableStorageOfResult();
138	}
139
140	// --------------------------------------------------------------------------
141	//
142	// Setup a task estimating the energy. The given task is cloned.
143	//
144	/*
145	void MJCut::SetEnergyEstimator(const MTask *task)
146	{
147	if (fEstimateEnergy)
148	delete fEstimateEnergy;
149	fEstimateEnergy = task ? (MTask*)task->Clone() : 0;
150	}
151	*/
152
153	// --------------------------------------------------------------------------
154	//
155	// Setup a task calculating the hadronness. The given task is cloned.
156	//
157	void MJCut::SetHadronnessCalculator(const MTask *task)
158	{
159	if (fCalcHadronness)
160	delete fCalcHadronness;
161	fCalcHadronness = task ? (MTask*)task->Clone() : 0;
162	}
163
164	// --------------------------------------------------------------------------
165	//
166	// return fOutputPath+"/ganymed%08d.root", num
167	//
168	TString MJCut::GetOutputFile(UInt_t num) const
169	{
170	TString p(fPathOut);
171	p += "/";
172	p += fNameOutput.IsNull() ? Form("ganymed%08d.root", num) : fNameOutput.Data();
173	return p;
174	}
175
176	/*
177	Bool_t MJCut::ReadTasks(const char fname, MTask &env1, MTask* &env2) const
178	{
179	// const 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("EstimateEnergy");
191	if (o && !o->InheritsFrom(MTask::Class()))
192	{
193	*fLog << err << dbginf << "ERROR - EstimateEnergy read from " << fname << " doesn't inherit from MTask!" << endl;
194	return kFALSE;
195	}
196	env1 = o ? (MTask*)o->Clone() : FNULL;
197
198	o = file.Get("CalcHadronness");
199	if (o && !o->InheritsFrom(MTask::Class()))
200	{
201	*fLog << err << dbginf << "ERROR - CalcHadronness read from " << fname << " doesn't inherit from MTask!" << endl;
202	return kFALSE;
203	}
204	env2 = o ? (MTask*)o->Clone() : NULL;
205
206	return kTRUE;
207	}
208	*/
209
210	// --------------------------------------------------------------------------
211	//
212	// Write the tasks in cont to the file corresponding to analysis number num,
213	// see GetOutputFile()
214	//
215	Bool_t MJCut::WriteTasks(UInt_t num, TObjArray &cont) const
216	{
217	if (fPathOut.IsNull())
218	{
219	*fLog << inf << "No output path specified via SetPathOut - no output written." << endl;
220	return kTRUE;
221	}
222
223	const TString oname(GetOutputFile(num));
224
225	*fLog << inf << "Writing to file: " << oname << endl;
226
227	TFile *file = 0;
228	if (fNameResult.IsNull() && fStoreResult)
229	{
230	file = (TFile*)gROOT->GetListOfFiles()->FindObject(oname);
231	if (file)
232	file->cd();
233	}
234	else
235	file = TFile::Open(oname, fOverwrite?"RECREATE":"NEW", "File created by MJCut", 9);
236
237	if (!file)
238	{
239	*fLog << err << "ERROR - Couldn't open file " << oname << " for writing..." << endl;
240	return kFALSE;
241	}
242
243	const Bool_t rc = WriteContainer(cont);
244
245	if (!(fNameResult.IsNull() && fStoreResult))
246	delete file;
247
248	return rc;
249	}
250
251	// --------------------------------------------------------------------------
252	//
253	// Write the result plots and other results to the file corresponding to
254	// analysis number num, see GetOutputFile()
255	//
256	Bool_t MJCut::WriteResult(const MParList &plist, UInt_t num) const
257	{
258	TObjArray arr;
259
260	// Save all MBinnings
261	TIter Next(plist);
262	TObject *o=0;
263	while ((o=Next()))
264	if (o->InheritsFrom(MBinning::Class()))
265	arr.Add(o);
266
267	// Save also the result, not only the setup
268	const MHAlpha halpha = (MHAlpha)plist.FindObject("Hist", "MHAlpha");
269	if (halpha)
270	arr.Add((TObject*)(&halpha->GetAlphaFitter()));
271
272	const TString fname(fNameOutput.IsNull() ? Form("ganymed%08d.root", num) : fNameOutput.Data());
273
274	// If requested, write to already open output file
275	if (fNameResult.IsNull() && fStoreResult)
276	{
277	TFile file = (TFile)gROOT->GetListOfFiles()->FindObject(fname);
278	if (file)
279	{
280	file->cd();
281	return WriteContainer(arr);
282	}
283	}
284
285	return WriteContainer(arr, fname, "UPDATE");
286	}
287
288	// --------------------------------------------------------------------------
289	//
290	// MJCut allows to setup several option by a resource file:
291	// MJCut.WriteSummary: yes, no
292	// MJCut.SummaryFile: filename
293	// MJCut.WriteResult: yes, no
294	// MJCut.ResultFile: filename
295	// MJCut.HistName: MHAlpha
296	//
297	Bool_t MJCut::CheckEnvLocal()
298	{
299	const TString f0(GetEnv("SummaryFile", ""));
300	const TString f1(GetEnv("ResultFile", ""));
301	if (!f0.IsNull())
302	SetNameSummaryFile(f0);
303	if (!f1.IsNull())
304	SetNameResultFile(f1);
305
306	EnableStorageOfSummary(GetEnv("SummaryFile", fStoreSummary));
307	EnableStorageOfResult(GetEnv("ResultFile", fStoreResult));
308	EnableWobbleMode(GetEnv("WobbleMode", fIsWobble));
309	EnableMonteCarloMode(GetEnv("MonteCarlo", fIsMonteCarlo));
310	EnableFullDisplay(GetEnv("FullDisplay", fFullDisplay));
311	//EnableSubstraction(GetEnv("HistogramSubstraction", fSubstraction));
312
313	SetNameHist(GetEnv("NameHist", fNameHist));
314	SetNameHistFS(GetEnv("NameHistFS", fNameHistFS));
315
316	return kTRUE;
317	}
318
319	// --------------------------------------------------------------------------
320	//
321	// Setup write to write:
322	// container tree optional?
323	// -------------- ---------- -----------
324	// "MHillas" to "Events"
325	// "MHillasSrc" to "Events"
326	// "MHadronness" to "Events" yes
327	// "MEnergyEst" to "Events" yes
328	// "DataType" to "Events"
329	//
330	void MJCut::SetupWriter(MWriteRootFile write, const char name) const
331	{
332	if (!write)
333	return;
334
335	write->SetName(name);
336	write->AddContainer("MHillas", "Events");
337	write->AddContainer("MHillasSrc", "Events");
338	write->AddContainer("MHillasExt", "Events");
339	write->AddContainer("MPointingPos", "Events");
340	write->AddContainer("MHillasSrcAnti", "Events", kFALSE);
341	write->AddContainer("MImagePar", "Events", kFALSE);
342	write->AddContainer("MNewImagePar", "Events", kFALSE);
343	write->AddContainer("MNewImagePar2", "Events", kFALSE);
344	write->AddContainer("MHadronness", "Events", kFALSE);
345	write->AddContainer("MSrcPosCam", "Events", kFALSE);
346	write->AddContainer("MSrcPosAnti", "Events", kFALSE);
347	write->AddContainer("ThetaSquared", "Events", kFALSE);
348	write->AddContainer("OpticalAxis", "Events", kFALSE);
349	write->AddContainer("Disp", "Events", kFALSE);
350	write->AddContainer("MTime", "Events", kFALSE);
351	write->AddContainer("MMcEvt", "Events", kFALSE);
352	write->AddContainer("DataType", "Events");
353	}
354
355	// --------------------------------------------------------------------------
356	//
357	// Create a new instance of an object with name name of class
358	// type fNameHist in parlist. It must derive from MHAlpha.
359	// Call ForceUsingSize for it and return its pointer.
360	// If something fails NULL is returned.
361	//
362	MHAlpha MJCut::CreateNewHist(MParList &plist, const char name) const
363	{
364	TClass *cls = gROOT->GetClass(fNameHist);
365	if (!cls)
366	{
367	*fLog << err << "Class " << fNameHist << " not found in dictionary... abort." << endl;
368	return NULL;
369	}
370	if (!cls->InheritsFrom(MHAlpha::Class()))
371	{
372	*fLog << err << "Class " << fNameHist << " doesn't inherit from MHAlpha... abort." << endl;
373	return NULL;
374	}
375
376	const TString objname(Form("Hist%s", name));
377	MHAlpha h = (MHAlpha)plist.FindCreateObj(fNameHist, objname);
378	if (!h)
379	return NULL;
380
381	h->ForceUsingSize();
382
383	return h;
384	}
385
386	// --------------------------------------------------------------------------
387	//
388	// Create a new instance of an object with name name of class
389	// type fNameHistFS in parlist. It must derive from MHFalseSource
390	// If something fails NULL is returned.
391	//
392	MH MJCut::CreateNewHistFS(MParList &plist, const char name) const
393	{
394	const TString cname(fNameHistFS.IsNull()?"MHFalseSource":fNameHistFS);
395
396	TClass *cls = gROOT->GetClass(cname);
397	if (!cls)
398	{
399	*fLog << err << "Class " << cname << " not found in dictionary... abort." << endl;
400	return NULL;
401	}
402	if (!cls->InheritsFrom("MHFalseSource"))
403	{
404	*fLog << err << "Class " << cname << " doesn't inherit from MHFalseSource... abort." << endl;
405	return NULL;
406	}
407
408	const TString objname(Form("FS%s", name));
409	return (MH*)plist.FindCreateObj(cname, objname);
410	}
411
412	Bool_t MJCut::Process(const MDataSet &set)
413	{
414	if (!set.IsValid())
415	{
416	*fLog << err << "ERROR - DataSet invalid!" << endl;
417	return kFALSE;
418	}
419
420	CheckEnv();
421
422	// --------------------------------------------------------------------------------
423
424	*fLog << inf;
425	fLog->Separator(GetDescriptor());
426	*fLog << "Perform cuts for data set " << set.GetName() << endl;
427	*fLog << endl;
428
429	// --------------------------------------------------------------------------------
430
431	// Setup Parlist
432	MParList plist;
433	plist.AddToList(this); // take care of fDisplay!
434
435	MParameterI par("DataType");
436	plist.AddToList(&par);
437
438	// Setup Tasklist
439	MTaskList tlist;
440	plist.AddToList(&tlist);
441
442	// La Palma Magic1
443	MObservatory obs;
444	plist.AddToList(&obs);
445
446	// Possible source position (eg. Wobble Mode)
447	MPointingPos source("MSourcePos");
448	if (set.HasSource())
449	{
450	if (!set.GetSourcePos(source))
451	return kFALSE;
452	plist.AddToList(&source);
453	*fLog << inf << "Using Source Position: " << source.GetTitle() << endl;
454	}
455	else
456	*fLog << inf << "No source position applied..." << endl;
457
458	// Initialize default binnings
459	MBinning bins1(18, 0, 90, "BinningAlpha", "lin");
460	MBinning bins2(15, 10, 1e6 , "BinningSize", "log");
461	MBinning bins3(67, -0.005, 0.665, "BinningTheta", "asin");
462	MBinning bins4("BinningFalseSource");
463	MBinning bins5("BinningWidth");
464	MBinning bins6("BinningLength");
465	MBinning bins7("BinningDist");
466	MBinning bins8("BinningMaxDist");
467	MBinning bins9("BinningM3Long");
468	MBinning bins0("BinningConc1");
469	plist.AddToList(&bins1);
470	plist.AddToList(&bins2);
471	plist.AddToList(&bins3);
472	plist.AddToList(&bins4);
473	plist.AddToList(&bins5);
474	plist.AddToList(&bins6);
475	plist.AddToList(&bins7);
476	plist.AddToList(&bins8);
477	plist.AddToList(&bins9);
478	plist.AddToList(&bins0);
479	//plist.AddToList(&binsa);
480
481	// --------------------------------------------------------------------------------
482
483	// Setup fitter and histograms
484	MAlphaFitter fit;
485	plist.AddToList(&fit);
486	if (fIsWobble)
487	fit.SetScaleMode(MAlphaFitter::kNone);
488
489	MHAlpha *halphaoff = CreateNewHist(plist, "Off");
490	MFillH falpha(halphaoff, "", "FillHist");
491	MH *hfsoff = CreateNewHistFS(plist, "Off");
492	MFillH ffs(hfsoff, "MHillas", "FillFS");
493
494	// FIXME: If fPathIn read cuts and energy estimator from file!
495	MContinue cont0("", "Cut0");
496	MContinue cont1("", "Cut1");
497	MContinue cont2("", "Cut2");
498	MContinue cont3("", "Cut3");
499	cont0.SetAllowEmpty();
500	cont1.SetAllowEmpty();
501	cont2.SetAllowEmpty();
502	cont3.SetAllowEmpty();
503
504	// ------------- Loop Off Data --------------------
505	MReadReports readoffdata;
506	readoffdata.AddTree("Events", "MTime.", MReadReports::kMaster);
507	readoffdata.AddTree("Drive", MReadReports::kRequired);
508	readoffdata.AddTree("Starguider", MReadReports::kRequired);
509	readoffdata.AddTree("EffectiveOnTime");
510
511	MReadMarsFile readoffmc("Events");
512	readoffmc.DisableAutoScheme();
513
514	MRead &readoff = fIsMonteCarlo ? (MRead&)readoffmc : (MRead&)readoffdata;
515	if (fIsWobble)
516	set.AddFilesOn(readoff);
517	else
518	set.AddFilesOff(readoff);
519
520	const TString path(Form("%s/", fPathOut.Data()));
521	TString fname0(path);
522	TString fname1(path);
523	fname0 += fNameSummary.IsNull() ? (TString) Form("ganymed%08d-summary.root", set.GetNumAnalysis()) : fNameSummary;
524	fname1 += fNameResult.IsNull() ? (TString) Form("ganymed%08d.root", set.GetNumAnalysis()) : fNameResult;
525
526	MWriteRootFile *write0 = CanStoreSummary() ? new MWriteRootFile(fPathOut.IsNull()?0:fname0.Data(), fOverwrite?"RECREATE":"NEW") : 0;
527	MWriteRootFile *write1 = CanStoreResult() ? new MWriteRootFile(fPathOut.IsNull()?0:fname1.Data(), fOverwrite?"RECREATE":"NEW") : 0;
528	SetupWriter(write0, "WriteAfterCut0");
529	SetupWriter(write1, "WriteAfterCut3");
530
531	/*
532	MEnergyEstimate est;
533
534	MTaskEnv taskenv1("EstimateEnergy");
535	taskenv1.SetDefault(fEstimateEnergy ? fEstimateEnergy : &est);
536	*/
537	MTaskEnv taskenv2("CalcHadronness");
538	taskenv2.SetDefault(fCalcHadronness);
539
540	MFillH fill1a("MHHillasOffPre [MHHillas]", "MHillas", "FillHillasPre");
541	MFillH fill2a("MHHillasOffPost [MHHillas]", "MHillas", "FillHillasPost");
542	MFillH fill3a("MHVsSizeOffPost [MHVsSize]", "MHillasSrc", "FillVsSizePost");
543	MFillH fill4a("MHHilExtOffPost [MHHillasExt]", "MHillasSrc", "FillHilExtPost");
544	MFillH fill5a("MHHilSrcOffPost [MHHillasSrc]", "MHillasSrc", "FillHilSrcPost");
545	MFillH fill6a("MHImgParOffPost [MHImagePar]", "MImagePar", "FillImgParPost");
546	MFillH fill7a("MHNewParOffPost [MHNewImagePar]", "MNewImagePar", "FillNewParPost");
547	fill1a.SetNameTab("PreCut");
548	fill2a.SetNameTab("PostCut");
549	fill3a.SetNameTab("VsSize");
550	fill4a.SetNameTab("HilExt");
551	fill5a.SetNameTab("HilSrc");
552	fill6a.SetNameTab("ImgPar");
553	fill7a.SetNameTab("NewPar");
554
555	MPrint print2("MEffectiveOnTime");
556	print2.EnableSkip();
557
558	// How to get source position from off- and on-data?
559	MSrcPosCalc scalc;
560	scalc.SetMode(fIsWobble?MSrcPosCalc::kWobble:MSrcPosCalc::kOffData); /********************/
561
562	MSrcPosCorrect scor;
563
564	MHillasCalc hcalc;
565	MHillasCalc hcalc2("MHillasCalcAnti");
566	hcalc.SetFlags(MHillasCalc::kCalcHillasSrc);
567	hcalc2.SetFlags(MHillasCalc::kCalcHillasSrc);
568	hcalc2.SetNameHillasSrc("MHillasSrcAnti");
569	hcalc2.SetNameSrcPosCam("MSrcPosAnti");
570
571	MTaskList tlist2;
572	tlist2.AddToList(&scalc);
573	tlist2.AddToList(&scor);
574	tlist2.AddToList(&hcalc);
575	if (fIsWobble)
576	tlist2.AddToList(&hcalc2);
577	//tlist2.AddToList(&taskenv1);
578	tlist2.AddToList(&taskenv2);
579	tlist2.AddToList(&cont0);
580	if (write0)
581	tlist2.AddToList(write0);
582	if (!fWriteOnly)
583	tlist2.AddToList(&fill1a);
584	tlist2.AddToList(&cont1);
585	if (!fWriteOnly && (!fIsWobble \|\| !fNameHistFS.IsNull()))
586	tlist2.AddToList(&ffs);
587	tlist2.AddToList(&cont2);
588	if (!fWriteOnly)
589	{
590	tlist2.AddToList(&fill2a);
591	if (fFullDisplay)
592	{
593	tlist2.AddToList(&fill3a);
594	tlist2.AddToList(&fill4a);
595	tlist2.AddToList(&fill5a);
596	tlist2.AddToList(&fill6a);
597	tlist2.AddToList(&fill7a);
598	}
599	}
600	if (!fWriteOnly)
601	tlist2.AddToList(&falpha);
602	tlist2.AddToList(&cont3);
603	if (write1)
604	tlist2.AddToList(write1);
605
606	MPointingDevCalc devcalc;
607
608	tlist.AddToList(&readoff);
609	if (gLog.GetDebugLevel()>4)
610	tlist.AddToList(&print2, "EffectiveOnTime");
611	tlist.AddToList(&devcalc, "Starguider");
612	tlist.AddToList(&tlist2, "Events");
613
614	par.SetVal(0);
615
616	// Create and setup the eventloop
617	MEvtLoop evtloop(fName);
618	evtloop.SetParList(&plist);
619	evtloop.SetDisplay(fDisplay);
620	evtloop.SetLogStream(fLog);
621	if (!SetupEnv(evtloop))
622	return kFALSE;
623
624	TObjArray cont;
625	cont.Add(&cont0);
626	cont.Add(&cont1);
627	cont.Add(&cont2);
628	cont.Add(&cont3);
629	//if (taskenv1.GetTask())
630	// cont.Add(taskenv1.GetTask());
631	if (taskenv2.GetTask())
632	cont.Add(taskenv2.GetTask());
633
634	if (!WriteTasks(set.GetNumAnalysis(), cont))
635	return kFALSE;
636
637	if (set.HasOffSequences() \|\| fIsWobble)
638	{
639	// Execute first analysis
640	if (!evtloop.Eventloop(fMaxEvents))
641	{
642	*fLog << err << GetDescriptor() << ": Processing of off-sequences failed." << endl;
643	return kFALSE;
644	}
645
646	if (!evtloop.GetDisplay())
647	{
648	*fLog << err << GetDescriptor() << ": Execution stopped by user." << endl;
649	return kFALSE;
650	}
651
652	//plist.FindObject("MTimeEffectiveOnTime")->Clear();
653	}
654	else
655	{
656	// This is the simplest way to remove the two object from the parlist
657	delete halphaoff;
658	delete hfsoff;
659	}
660
661	// ------------- Loop On Data --------------------
662	MReadReports readondata;
663	readondata.AddTree("Events", "MTime.", MReadReports::kMaster);
664	readondata.AddTree("Drive", MReadReports::kRequired);
665	readondata.AddTree("Starguider", MReadReports::kRequired);
666	readondata.AddTree("EffectiveOnTime");
667
668	MReadMarsFile readonmc("Events");
669	readonmc.DisableAutoScheme();
670
671	MRead &readon = fIsMonteCarlo ? (MRead&)readonmc : (MRead&)readondata;
672	set.AddFilesOn(readon);
673
674	scalc.SetMode(MSrcPosCalc::kDefault);
675
676	MFillH fill1b("MHHillasOnPre [MHHillas]", "MHillas", "FillHillasPre");
677	MFillH fill2b("MHHillasOnPost [MHHillas]", "MHillas", "FillHillasPost");
678	MFillH fill3b("MHVsSizeOnPost [MHVsSize]", "MHillasSrc", "FillVsSizePost");
679	MFillH fill4b("MHHilExtOnPost [MHHillasExt]", "MHillasSrc", "FillHilExtPost");
680	MFillH fill5b("MHHilSrcOnPost [MHHillasSrc]", "MHillasSrc", "FillHilSrcPost");
681	MFillH fill6b("MHImgParOnPost [MHImagePar]", "MImagePar", "FillImgParPost");
682	MFillH fill7b("MHNewParOnPost [MHNewImagePar]", "MNewImagePar", "FillNewParPost");
683	fill1b.SetNameTab("PreCut");
684	fill2b.SetNameTab("PostCut");
685	fill3b.SetNameTab("VsSize");
686	fill4b.SetNameTab("HilExt");
687	fill5b.SetNameTab("HilSrc");
688	fill6b.SetNameTab("ImgPar");
689	fill7b.SetNameTab("NewPar");
690	fill1b.SetDrawOption(set.HasOffSequences()\|\|fIsWobble?"same":"");
691	fill2b.SetDrawOption(set.HasOffSequences()\|\|fIsWobble?"same":"");
692	fill3b.SetDrawOption(set.HasOffSequences()\|\|fIsWobble?"same":"");
693	fill4b.SetDrawOption(set.HasOffSequences()\|\|fIsWobble?"same":"");
694	fill5b.SetDrawOption(set.HasOffSequences()\|\|fIsWobble?"same":"");
695	fill6b.SetDrawOption(set.HasOffSequences()\|\|fIsWobble?"same":"");
696	fill7b.SetDrawOption(set.HasOffSequences()\|\|fIsWobble?"same":"");
697
698	/*
699	MHVsTime hvs("MEffectiveOnTime.fVal");
700	hvs.SetTitle("Effective On-Time vs. Time;;T_{on}");
701	MFillH fillvs(&hvs, "MTimeEffectiveOnTime", "FillOnTime");
702	fillvs.SetNameTab("OnTime");
703	*/
704	MH3 hvs("MPointingPos.fZd");
705	hvs.SetName("Theta");
706	hvs.SetTitle("Effective On-Time vs. Zenith Angle;\\Theta [\\circ];T_{on} [s]");
707	MFillH fillvs(&hvs, "", "FillOnTime");
708	if (!fIsMonteCarlo)
709	fillvs.SetWeight("MEffectiveOnTime");
710	fillvs.SetNameTab("OnTime");
711
712	/*
713	MParameterD weight;
714	weight.SetVal(-1);
715	fill2a.SetWeight(&weight);
716	fill3a.SetWeight(&weight);
717	fill4a.SetWeight(&weight);
718	fill5a.SetWeight(&weight);
719	fill6a.SetWeight(&weight);
720	fill7a.SetWeight(&weight);
721	if (fSubstraction)
722	{
723	fill2a.SetNameTab("PostCut-");
724	fill3a.SetNameTab("VsSize-");
725	fill4a.SetNameTab("HilExt-");
726	fill5a.SetNameTab("HilSrc-");
727	fill6a.SetNameTab("ImgPar-");
728	fill7a.SetNameTab("NewPar-");
729	}
730	*/
731	MHAlpha *halphaon=CreateNewHist(plist);
732	MFillH falpha2(halphaon, "", "FillHist");
733	MH *hfs=CreateNewHistFS(plist);
734	MFillH ffs2(hfs, "MHillas", "FillFS");
735	MFillH fillphi("MHPhi", "", "FillPhi");
736	fillphi.SetDrawOption("anticut");
737
738	tlist.Replace(&readon);
739	if (!fWriteOnly)
740	{
741	tlist2.Replace(&fill1b);
742	/* if (fIsWobble)
743	{
744	tlist2.AddToListAfter(&fill2b, &fill2a);
745	tlist2.AddToListAfter(&fill3b, &fill3a);
746	}
747	else
748	*/
749	tlist2.Replace(&fill2b);
750	if (fFullDisplay)
751	{
752	tlist2.Replace(&fill3b);
753	tlist2.Replace(&fill4b);
754	tlist2.Replace(&fill5b);
755	tlist2.Replace(&fill6b);
756	tlist2.Replace(&fill7b);
757	}
758	tlist2.Replace(&falpha2);
759	if (!fIsWobble \|\| !fNameHist.IsNull())
760	tlist2.Replace(&ffs2);
761	if (fIsWobble)
762	{
763	tlist2.AddToListAfter(&fillphi, &falpha2);
764	if (!fNameHist.IsNull())
765	tlist2.RemoveFromList(&ffs);
766	}
767
768	if (!fIsMonteCarlo)
769	tlist.AddToList(&fillvs, "EffectiveOnTime");
770	else
771	tlist2.AddToListBefore(&fillvs, &scalc);
772	}
773
774	par.SetVal(1);
775
776	// Execute first analysis
777	if (!evtloop.Eventloop(fMaxEvents))
778	{
779	*fLog << err << GetDescriptor() << ": Processing of on-sequences failed." << endl;
780	return kFALSE;
781	}
782
783	if (write0)
784	delete write0;
785	if (write1)
786	delete write1;
787
788	// FIXME: Perform fit and plot energy dependant alpha plots
789	// and fit result to new tabs!
790	if (!WriteResult(plist, set.GetNumAnalysis()))
791	return kFALSE;
792
793	*fLog << all << GetDescriptor() << ": Done." << endl;
794	*fLog << endl << endl;
795
796	return kTRUE;
797	}

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

Download in other formats: