Context Navigation

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

Visit:

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

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

Download in other formats: