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