Index: trunk/MagicSoft/Mars/Changelog =================================================================== --- trunk/MagicSoft/Mars/Changelog (revision 7120) +++ trunk/MagicSoft/Mars/Changelog (revision 7121) @@ -21,4 +21,19 @@ -*-*- END OF LINE -*-*- + 2005/06/01 Thomas Bretz + + * mjobs/MJOptimize.[h,cc]: + - moved code for cut optimization to MJOptimizeCuts + - moved code for energy optimization to MJOptimizeEnergy + + * mjobs/MJOptimizeCuts.[h,cc], mjobs/MJOptimizeEnergy.[h,cc]: + - added + + * mjobs/Makefile, mjobs/JobsLinkDef.h: + - added MJOptimizeCuts + - added MJOptimizeEnergy + + + 2005/06/01 Daniela Dorner Index: trunk/MagicSoft/Mars/NEWS =================================================================== --- trunk/MagicSoft/Mars/NEWS (revision 7120) +++ trunk/MagicSoft/Mars/NEWS (revision 7121) @@ -47,4 +47,12 @@ For details about the requirements of these plots look at the corresponding classes + + - MJOptimize: For On/Off and Wobble optimization MJOptimize + has been replaced by MJOptimizeCuts. To use MHThetaSq instead + of MHAlpha call MJOptimizeCuts("MHThetaSq"). Make sure in this + case MAlphaFitter is setup correctly. For On-Only optimization + MJOptimize has been Replaced by MJOptimizeCuts, too. Call + RunOn() instead of Run(). For energy estimation MJOptimize + has been replaced by MJOptimizeEnergy - sponde: The input MC spectrum can now be weighted to fake a Index: trunk/MagicSoft/Mars/mjobs/JobsLinkDef.h =================================================================== --- trunk/MagicSoft/Mars/mjobs/JobsLinkDef.h (revision 7120) +++ trunk/MagicSoft/Mars/mjobs/JobsLinkDef.h (revision 7121) @@ -17,6 +17,9 @@ #pragma link C++ class MJStar+; #pragma link C++ class MJCut+; -#pragma link C++ class MJOptimize+; #pragma link C++ class MJSpectrum+; +#pragma link C++ class MJOptimize+; +#pragma link C++ class MJOptimizeCuts+; +#pragma link C++ class MJOptimizeEnergy+; + #endif Index: trunk/MagicSoft/Mars/mjobs/MJOptimize.cc =================================================================== --- trunk/MagicSoft/Mars/mjobs/MJOptimize.cc (revision 7120) +++ trunk/MagicSoft/Mars/mjobs/MJOptimize.cc (revision 7121) @@ -90,34 +90,29 @@ #include +#include "MHMatrix.h" + +// environment +#include "MLog.h" +#include "MLogManip.h" + +#include "MDirIter.h" +#include "MStatusDisplay.h" + +// eventloop #include "MParList.h" #include "MTaskList.h" -#include "MGeomCamCT1.h" -#include "MFEventSelector.h" -#include "MFDataMember.h" +#include "MEvtLoop.h" + +// parameters +#include "MParameters.h" + +// tasks #include "MReadTree.h" -#include "MHMatrix.h" -#include "MEnergyEstimate.h" #include "MMatrixLoop.h" -#include "MChisqEval.h" -#include "MEvtLoop.h" -#include "MDataElement.h" -#include "MDataMember.h" -#include "MLog.h" -#include "MLogManip.h" -#include "MParameters.h" #include "MFillH.h" + +// filters #include "MF.h" #include "MFilterList.h" -#include "../mfilter/MFMagicCuts.h" -#include "MContinue.h" -#include "MGeomCamMagic.h" -#include "MHillasSrcCalc.h" -#include "MHMatrix.h" -#include "MMatrixLoop.h" -#include "../mhflux/MHAlpha.h" -#include "MStatusDisplay.h" -#include "../mhflux/MHEnergyEst.h" -#include "MDirIter.h" -#include "MSequence.h" using namespace std; @@ -507,4 +502,10 @@ Bool_t MJOptimize::Optimize(MEvtLoop &evtloop) { + if (fParameters.GetSize()==0) + { + *fLog << err << GetDescriptor() << "::Optimize: ERROR - Sorry, no parameters defined." << endl; + return kFALSE; + } + if (fType==kNone) return kTRUE; @@ -806,270 +807,2 @@ *fLog << endl; } - -//------------------------------------------------------------------------ -// -Bool_t MJOptimize::Run(const char *fname, MFilter *filter, MAlphaFitter *fit) -{ - if (fParameters.GetSize()==0) - { - *fLog << err << "Sorry, no parameters defined." << endl; - return kFALSE; - } - - fLog->Separator("Preparing On-only-optimization"); - - MParList parlist; - - MGeomCamMagic geom; // For GetConvMm2Deg - parlist.AddToList(&geom); - - MHMatrix m("M"); - AddRulesToMatrix(m); - parlist.AddToList(&m); - - MHAlpha hist; - hist.SkipHistTime(); - hist.SkipHistTheta(); - hist.SkipHistEnergy(); - hist.InitMapping(&m); - - if (filter && filter->InheritsFrom(MFMagicCuts::Class())) - ((MFMagicCuts*)filter)->InitMapping(&m); - - MReadTree read("Events"); - read.DisableAutoScheme(); // AutoScheme doesn't seem to be faster! - if (fname) - read.AddFile(fname); - else - AddSequences(read, fNamesOn); - if (!FillMatrix(read, parlist)) - return kFALSE; - - MTaskList tasklist; - parlist.Replace(&tasklist); - if (fit) - parlist.AddToList(fit); - - MFilterList list; - SetupFilters(list, filter); - - MContinue contin(&list); - parlist.AddToList(&list); - - MFillH fill(&hist); - - MMatrixLoop loop(&m); - - tasklist.AddToList(&loop); - tasklist.AddToList(&list); - tasklist.AddToList(&contin); - tasklist.AddToList(&fill); - - // Optimize with the tasklist in this parameterlist - if (!Optimize(parlist)) - return kFALSE; - - // Copy the result back to be accessible by the user - if (fit) - hist.GetAlphaFitter().Copy(*fit); - // Print the result - hist.GetAlphaFitter().Print("result"); - - // Store result if requested - TObjArray cont; - cont.Add(&contin); - return WriteContainer(cont, fNameOut); -} - -//------------------------------------------------------------------------ -// -// Make sure, that filter->GetDataMember is correctly implemented!!!! -// -Bool_t MJOptimize::RunOnOff(const char *fname, MFilter *filter, MAlphaFitter *fit, const char *tree) -{ - if (fParameters.GetSize()==0) - { - *fLog << err << "Sorry, no parameters defined." << endl; - return kFALSE; - } - - fLog->Separator("Preparing On/Off-optimization"); - - MParList parlist; - - MGeomCamMagic geom; // For GetConvMm2Deg - parlist.AddToList(&geom); - - MHMatrix m("M"); - AddRulesToMatrix(m); - parlist.AddToList(&m); - - const Int_t idxdatatype = m.AddColumn("DataType.fVal"); - - MHAlpha histon, histof("MHAlphaOff"); - histon.SkipHistTime(); - histon.SkipHistTheta(); - //histon.SkipHistEnergy(); - histof.SkipHistTime(); - histof.SkipHistTheta(); - //histof.SkipHistEnergy(); - histon.ForceUsingSize(); - histof.ForceUsingSize(); - histon.InitMapping(&m, 1); - histof.InitMapping(&m, 1); - - if (filter && filter->InheritsFrom(MFMagicCuts::Class())) - ((MFMagicCuts*)filter)->InitMapping(&m); - - parlist.AddToList(&histon); - parlist.AddToList(&histof); - - if (fname) - { - MReadTree read(tree); - read.DisableAutoScheme(); // AutoScheme doesn't seem to be faster! - read.AddFile(fname); - if (!FillMatrix(read, parlist)) - return kFALSE; - } - else - { - MParameterI par("DataType"); - parlist.AddToList(&par); - - gLog.Separator("Reading On-Data"); - par.SetVal(1); - MReadTree readon(tree); - readon.DisableAutoScheme(); // AutoScheme doesn't seem to be faster! - AddSequences(readon, fNamesOn); - if (!FillMatrix(readon, parlist)) - return kFALSE; - - gLog.Separator("Reading Off-Data"); - par.SetVal(0); - MReadTree readoff(tree); - readoff.DisableAutoScheme(); // AutoScheme doesn't seem to be faster! - AddSequences(readoff, fNamesOff); - if (!FillMatrix(readoff, parlist)) - return kFALSE; - } - - MTaskList tasklist; - parlist.Replace(&tasklist); - if (fit) - parlist.AddToList(fit); - - MFilterList list; - SetupFilters(list, filter); - - MContinue contin(&list); - parlist.AddToList(&list); - - MFillH fillof(&histof, "", "FillHistOff"); - MFillH fillon(&histon, "", "FillHistOn"); - - MF f0(Form("M[%d]<0.5", idxdatatype), "FilterOffData"); - MF f1(Form("M[%d]>0.5", idxdatatype), "FilterOnData"); - - fillof.SetFilter(&f0); - fillon.SetFilter(&f1); - - MMatrixLoop loop(&m); - - tasklist.AddToList(&loop); - tasklist.AddToList(&list); - tasklist.AddToList(&contin); - tasklist.AddToList(&f0); - tasklist.AddToList(&f1); - tasklist.AddToList(&fillof); - tasklist.AddToList(&fillon); - - // Optimize with the tasklist in this parameterlist - if (!Optimize(parlist)) - return kFALSE; - - // Copy the result back to be accessible by the user - if (fit) - histon.GetAlphaFitter().Copy(*fit); - // Print the result - histon.GetAlphaFitter().Print("result"); - - // Store result if requested - TObjArray cont; - cont.Add(&contin); - return WriteContainer(cont, fNameOut); -} - -//------------------------------------------------------------------------ -// -// Read all events from file which do match rules and optimize -// energy estimator. -// -Bool_t MJOptimize::RunEnergy(const char *fname, const char *rule) -{ - if (fParameters.GetSize()==0) - { - *fLog << err << "Sorry, no parameters defined." << endl; - return kFALSE; - } - - fLog->Separator("Preparing Energy optimization"); - - MParList parlist; - - MParameterI par("DataType"); - par.SetVal(1); - parlist.AddToList(&par); - - MFDataMember filter("DataType.fVal", '>', 0.5); - fPreCuts.Add(&filter); - - MGeomCamMagic geom; // For GetConvMm2Deg - parlist.AddToList(&geom); - - MHMatrix m("M"); - AddRulesToMatrix(m); - parlist.AddToList(&m); - - MHEnergyEst hist; - hist.InitMapping(&m); - - MEnergyEstimate est("MParameters"); - est.SetRule(rule); - parlist.AddToList(&est); - - MReadTree read("Events"); - // NECESSARY BECAUSE OF MDataFormula GetRules missing - read.DisableAutoScheme(); - if (fname) - read.AddFile(fname); - else - AddSequences(read, fNamesOn); - if (!FillMatrix(read, parlist, kTRUE)) - return kFALSE; - - fPreCuts.Remove(&filter); - - MTaskList tasklist; - parlist.Replace(&tasklist); - - MFillH fill(&hist); - - MMatrixLoop loop(&m); - - tasklist.AddToList(&loop); - tasklist.AddToList(&est); - tasklist.AddToList(&fill); - - // Optimize with the tasklist in this parameterlist - if (!Optimize(parlist)) - return kFALSE; - - // Print the result - hist.Print(); - - // Store result if requested - TObjArray cont; - cont.Add(&est); - return WriteContainer(cont, fNameOut); -} Index: trunk/MagicSoft/Mars/mjobs/MJOptimize.h =================================================================== --- trunk/MagicSoft/Mars/mjobs/MJOptimize.h (revision 7120) +++ trunk/MagicSoft/Mars/mjobs/MJOptimize.h (revision 7121) @@ -46,21 +46,9 @@ TList fRules; TList fFilter; - TList fPreCuts; - - TList fNamesOn; - TList fNamesOff; - - TString fNameOut; void AddPoint(TList *l, Int_t idx, Float_t val) const; TList *GetPlots() const; - void AddRulesToMatrix(MHMatrix &m) const; - void SetupFilters(MFilterList &list, MFilter *filter=0) const; - Bool_t FillMatrix(MReadTree &read, MParList &l, Bool_t userules=kFALSE); - MEvtLoop *fEvtLoop; //! - - Bool_t AddSequences(MRead &read, TList &list) const; // Minuit Interface @@ -87,4 +75,17 @@ Bool_t Optimize(MEvtLoop &evtloop); + +protected: + TList fPreCuts; + + TList fNamesOn; + TList fNamesOff; + + TString fNameOut; + + void AddRulesToMatrix(MHMatrix &m) const; + void SetupFilters(MFilterList &list, MFilter *filter=0) const; + Bool_t AddSequences(MRead &read, TList &list) const; + Bool_t FillMatrix(MReadTree &read, MParList &l, Bool_t userules=kFALSE); public: @@ -153,40 +154,4 @@ Bool_t Optimize(MParList &list); - // Special optimizing routines - Bool_t Run(const char *fname, MFilter *filter, MAlphaFitter *fit=0); - Bool_t Run(const char *fname, MAlphaFitter *fit=0) - { - return Run(fname, 0, fit); - } - - Bool_t Run(MFilter *filter, MAlphaFitter *fit=0) - { - return Run(0, filter, fit); - } - Bool_t Run(MAlphaFitter *fit=0) - { - return Run(0, 0, fit); - } - - Bool_t RunOnOff(const char *fname, MFilter *filter, MAlphaFitter *fit=0, const char *tree="Events"); - Bool_t RunOnOff(const char *fname, MAlphaFitter *fit=0, const char *tree="Events") - { - return RunOnOff(fname, 0, fit, tree); - } - Bool_t RunEnergy(const char *fname, const char *rule); - - Bool_t RunOnOff(MFilter *filter, MAlphaFitter *fit=0, const char *tree="Events") - { - return RunOnOff(0, filter, fit, tree); - } - Bool_t RunOnOff(MAlphaFitter *fit=0, const char *tree="Events") - { - return RunOnOff(fit, tree); - } - Bool_t RunEnergy(const char *rule) - { - return RunEnergy(0, rule); - } - ClassDef(MJOptimize, 0) // Class for optimization of the Supercuts }; Index: trunk/MagicSoft/Mars/mjobs/MJOptimizeCuts.cc =================================================================== --- trunk/MagicSoft/Mars/mjobs/MJOptimizeCuts.cc (revision 7121) +++ trunk/MagicSoft/Mars/mjobs/MJOptimizeCuts.cc (revision 7121) @@ -0,0 +1,338 @@ +/* ======================================================================== *\ +! +! * +! * This file is part of MARS, the MAGIC Analysis and Reconstruction +! * Software. It is distributed to you in the hope that it can be a useful +! * and timesaving tool in analysing Data of imaging Cerenkov telescopes. +! * It is distributed WITHOUT ANY WARRANTY. +! * +! * Permission to use, copy, modify and distribute this software and its +! * documentation for any purpose is hereby granted without fee, +! * provided that the above copyright notice appear in all copies and +! * that both that copyright notice and this permission notice appear +! * in supporting documentation. It is provided "as is" without express +! * or implied warranty. +! * +! +! +! Author(s): Thomas Bretz, 9/2004 +! +! Copyright: MAGIC Software Development, 2000-2004 +! +! +\* ======================================================================== */ + +///////////////////////////////////////////////////////////////////////////// +// +// MJOptimizeCuts +// +// Class for otimizing the parameters of the supercuts +// +// Minimization Control +// ==================== +// +// To choose the minimization algorithm use: +// void SetOptimizer(Optimizer_t o); +// +// Allowed options are: +// enum Optimizer_t +// { +// kMigrad, // Minimize by the method of Migrad +// kSimplex, // Minimize by the method of Simplex +// kMinimize, // Migrad + Simplex (if Migrad fails) +// kMinos, // Minos error determination +// kImprove, // Local minimum search +// kSeek, // Minimize by the method of Monte Carlo +// kNone // Skip optimization +// }; +// +// For more details on the methods see TMinuit. +// +// +// You can change the behaviour of the minimization using +// +// void SetNumMaxCalls(UInt_t num=0); +// void SetTolerance(Float_t tol=0); +// +// While NumMaxCalls is the first, Tolerance the second arguement. +// For more details start root and type +// +// gMinuit->mnhelp("command") +// +// while command can be +// * MIGRAD +// * SIMPLEX +// * MINIMIZE +// * MINOS +// * IMPROVE +// * SEEK +// +// The default (num==0 and tol==0) should always give you the +// corresponding defaults used in Minuit. +// +// +// FIXME: Implement changing cut in hadronness... +// FIXME: Show MHSignificance on MStatusDisplay during filling... +// FIXME: Choose step-size percentage as static data membewr +// FIXME: Choose minimization method +// +///////////////////////////////////////////////////////////////////////////// +#include "MJOptimizeCuts.h" + +#include "MHMatrix.h" + +// Environment +#include "MLog.h" +#include "MLogManip.h" + +// Eventloop +#include "MReadTree.h" +#include "MParList.h" +#include "MTaskList.h" + +// Parameter container +#include "MGeomCamMagic.h" +#include "MParameters.h" +#include "MFilterList.h" + +// histograms +#include "../mhflux/MHAlpha.h" + +// Tasks +#include "MF.h" +#include "MFillH.h" +#include "MContinue.h" +#include "MMatrixLoop.h" + +#include "MFMagicCuts.h" + +ClassImp(MJOptimizeCuts); + +using namespace std; + +//------------------------------------------------------------------------ +// +MHAlpha *MJOptimizeCuts::CreateNewHist(const char *name) const +{ + TClass *cls = gROOT->GetClass(fNameHist); + if (!cls) + { + *fLog << err << "Class " << fNameHist << " not found in dictionary... abort." << endl; + return NULL; + } + if (!cls->InheritsFrom(MHAlpha::Class())) + { + *fLog << err << "Class " << fNameHist << " doesn't inherit from MHAlpha... abort." << endl; + return NULL; + } + + MHAlpha *h = (MHAlpha*)cls->New(); + if (h && name) + h->SetName(name); + + return h; +} + +//------------------------------------------------------------------------ +// +Bool_t MJOptimizeCuts::RunOnOffCore(MHAlpha &histon, MHAlpha &histof, const char *fname, MFilter *filter, MAlphaFitter *fit, const char *tree) +{ + fLog->Separator("Preparing On/Off-optimization"); + + MParList parlist; + + MGeomCamMagic geom; // For GetConvMm2Deg + parlist.AddToList(&geom); + + MHMatrix m("M"); + AddRulesToMatrix(m); + parlist.AddToList(&m); + + const Int_t idxdatatype = m.AddColumn("DataType.fVal"); + + histon.SkipHistTime(); + histon.SkipHistTheta(); + //histon.SkipHistEnergy(); + histof.SkipHistTime(); + histof.SkipHistTheta(); + //histof.SkipHistEnergy(); + histon.ForceUsingSize(); + histof.ForceUsingSize(); + histon.InitMapping(&m, 1); + histof.InitMapping(&m, 1); + + if (filter && filter->InheritsFrom(MFMagicCuts::Class())) + ((MFMagicCuts*)filter)->InitMapping(&m); + + parlist.AddToList(&histon); + parlist.AddToList(&histof); + + if (fname) + { + MReadTree read(tree); + read.DisableAutoScheme(); // AutoScheme doesn't seem to be faster! + read.AddFile(fname); + if (!FillMatrix(read, parlist)) + return kFALSE; + } + else + { + MParameterI par("DataType"); + parlist.AddToList(&par); + + gLog.Separator("Reading On-Data"); + par.SetVal(1); + MReadTree readon(tree); + readon.DisableAutoScheme(); // AutoScheme doesn't seem to be faster! + AddSequences(readon, fNamesOn); + if (!FillMatrix(readon, parlist)) + return kFALSE; + + gLog.Separator("Reading Off-Data"); + par.SetVal(0); + MReadTree readoff(tree); + readoff.DisableAutoScheme(); // AutoScheme doesn't seem to be faster! + AddSequences(readoff, fNamesOff); + if (!FillMatrix(readoff, parlist)) + return kFALSE; + } + + MTaskList tasklist; + parlist.Replace(&tasklist); + if (fit) + parlist.AddToList(fit); + + MFilterList list; + SetupFilters(list, filter); + + MContinue contin(&list); + parlist.AddToList(&list); + + MFillH fillof(&histof, "", "FillHistOff"); + MFillH fillon(&histon, "", "FillHistOn"); + + MF f0(Form("M[%d]<0.5", idxdatatype), "FilterOffData"); + MF f1(Form("M[%d]>0.5", idxdatatype), "FilterOnData"); + + fillof.SetFilter(&f0); + fillon.SetFilter(&f1); + + MMatrixLoop loop(&m); + + tasklist.AddToList(&loop); + tasklist.AddToList(&list); + tasklist.AddToList(&contin); + tasklist.AddToList(&f0); + tasklist.AddToList(&f1); + tasklist.AddToList(&fillof); + tasklist.AddToList(&fillon); + + // Optimize with the tasklist in this parameterlist + if (!Optimize(parlist)) + return kFALSE; + + // Copy the result back to be accessible by the user + if (fit) + histon.GetAlphaFitter().Copy(*fit); + + // Print the result + histon.GetAlphaFitter().Print("result"); + + // Store result if requested + TObjArray cont; + cont.Add(&contin); + return WriteContainer(cont, fNameOut); +} + +Bool_t MJOptimizeCuts::RunOnCore(MHAlpha &hist, const char *fname, MFilter *filter, MAlphaFitter *fit) +{ + fLog->Separator("Preparing On-only-optimization"); + + MParList parlist; + + MGeomCamMagic geom; // For GetConvMm2Deg + parlist.AddToList(&geom); + + MHMatrix m("M"); + AddRulesToMatrix(m); + parlist.AddToList(&m); + + hist.SkipHistTime(); + hist.SkipHistTheta(); + hist.SkipHistEnergy(); + hist.InitMapping(&m); + + if (filter && filter->InheritsFrom(MFMagicCuts::Class())) + ((MFMagicCuts*)filter)->InitMapping(&m); + + MReadTree read("Events"); + read.DisableAutoScheme(); // AutoScheme doesn't seem to be faster! + if (fname) + read.AddFile(fname); + else + AddSequences(read, fNamesOn); + if (!FillMatrix(read, parlist)) + return kFALSE; + + MTaskList tasklist; + parlist.Replace(&tasklist); + if (fit) + parlist.AddToList(fit); + + MFilterList list; + SetupFilters(list, filter); + + MContinue contin(&list); + parlist.AddToList(&list); + + MFillH fill(&hist); + + MMatrixLoop loop(&m); + + tasklist.AddToList(&loop); + tasklist.AddToList(&list); + tasklist.AddToList(&contin); + tasklist.AddToList(&fill); + + // Optimize with the tasklist in this parameterlist + if (!Optimize(parlist)) + return kFALSE; + + // Copy the result back to be accessible by the user + if (fit) + hist.GetAlphaFitter().Copy(*fit); + // Print the result + hist.GetAlphaFitter().Print("result"); + + // Store result if requested + TObjArray cont; + cont.Add(&contin); + return WriteContainer(cont, fNameOut); +} + +//------------------------------------------------------------------------ +// +Bool_t MJOptimizeCuts::RunOnOff(const char *fname, MFilter *filter, MAlphaFitter *fit, const char *tree) +{ + MHAlpha *histon = CreateNewHist(); + MHAlpha *histof = CreateNewHist("MHAlphaOff"); + + const Bool_t rc = RunOnOffCore(*histon, *histof, fname, filter, fit, tree); + + delete histon; + delete histof; + + return rc; +} + +//------------------------------------------------------------------------ +// +Bool_t MJOptimizeCuts::RunOn(const char *fname, MFilter *filter, MAlphaFitter *fit) +{ + MHAlpha *histon = CreateNewHist(); + + const Bool_t rc = RunOnCore(*histon, fname, filter, fit); + + delete histon; + return rc; +} Index: trunk/MagicSoft/Mars/mjobs/MJOptimizeCuts.h =================================================================== --- trunk/MagicSoft/Mars/mjobs/MJOptimizeCuts.h (revision 7121) +++ trunk/MagicSoft/Mars/mjobs/MJOptimizeCuts.h (revision 7121) @@ -0,0 +1,57 @@ +#ifndef MARS_MJOptimizeCuts +#define MARS_MJOptimizeCuts + +#ifndef MARS_MJOptimize +#include "MJOptimize.h" +#endif + +class MHAlpha; + +class MJOptimizeCuts : public MJOptimize +{ +private: + TString fNameHist; + + MHAlpha *CreateNewHist(const char *name=0) const; + + Bool_t RunOnCore(MHAlpha &hist, const char *fname, MFilter *filter, MAlphaFitter *fit); + Bool_t RunOnOffCore(MHAlpha &histon, MHAlpha &histoff, const char *fname, MFilter *filter, MAlphaFitter *fit, const char *tree); + +public: + MJOptimizeCuts(const char *name="MHAlpha") : MJOptimize(), fNameHist(name) { } + + // Optimization On-Off/Wobble + Bool_t RunOnOff(const char *fname, MFilter *filter, MAlphaFitter *fit=0, const char *tree="Events"); + Bool_t RunOnOff(const char *fname, MAlphaFitter *fit=0, const char *tree="Events") + { + return RunOnOff(fname, 0, fit, tree); + } + Bool_t RunOnOff(MFilter *filter, MAlphaFitter *fit=0, const char *tree="Events") + { + return RunOnOff(0, filter, fit, tree); + } + Bool_t RunOnOff(MAlphaFitter *fit=0, const char *tree="Events") + { + return RunOnOff(fit, tree); + } + + // Optimization On-only + Bool_t RunOn(const char *fname, MFilter *filter, MAlphaFitter *fit=0); + Bool_t RunOn(const char *fname, MAlphaFitter *fit=0) + { + return RunOn(fname, 0, fit); + } + + Bool_t RunOn(MFilter *filter, MAlphaFitter *fit=0) + { + return RunOn(0, filter, fit); + } + Bool_t RunOn(MAlphaFitter *fit=0) + { + return RunOn(0, 0, fit); + } + + ClassDef(MJOptimizeCuts, 0) // Class for optimization of the Supercuts +}; + +#endif Index: trunk/MagicSoft/Mars/mjobs/MJOptimizeEnergy.cc =================================================================== --- trunk/MagicSoft/Mars/mjobs/MJOptimizeEnergy.cc (revision 7121) +++ trunk/MagicSoft/Mars/mjobs/MJOptimizeEnergy.cc (revision 7121) @@ -0,0 +1,179 @@ +/* ======================================================================== *\ +! +! * +! * This file is part of MARS, the MAGIC Analysis and Reconstruction +! * Software. It is distributed to you in the hope that it can be a useful +! * and timesaving tool in analysing Data of imaging Cerenkov telescopes. +! * It is distributed WITHOUT ANY WARRANTY. +! * +! * Permission to use, copy, modify and distribute this software and its +! * documentation for any purpose is hereby granted without fee, +! * provided that the above copyright notice appear in all copies and +! * that both that copyright notice and this permission notice appear +! * in supporting documentation. It is provided "as is" without express +! * or implied warranty. +! * +! +! +! Author(s): Thomas Bretz, 9/2004 +! +! Copyright: MAGIC Software Development, 2000-2004 +! +! +\* ======================================================================== */ + +///////////////////////////////////////////////////////////////////////////// +// +// MJOptimize +// +// Class for otimizing the parameters of the supercuts +// +// Minimization Control +// ==================== +// +// To choose the minimization algorithm use: +// void SetOptimizer(Optimizer_t o); +// +// Allowed options are: +// enum Optimizer_t +// { +// kMigrad, // Minimize by the method of Migrad +// kSimplex, // Minimize by the method of Simplex +// kMinimize, // Migrad + Simplex (if Migrad fails) +// kMinos, // Minos error determination +// kImprove, // Local minimum search +// kSeek, // Minimize by the method of Monte Carlo +// kNone // Skip optimization +// }; +// +// For more details on the methods see TMinuit. +// +// +// You can change the behaviour of the minimization using +// +// void SetNumMaxCalls(UInt_t num=0); +// void SetTolerance(Float_t tol=0); +// +// While NumMaxCalls is the first, Tolerance the second arguement. +// For more details start root and type +// +// gMinuit->mnhelp("command") +// +// while command can be +// * MIGRAD +// * SIMPLEX +// * MINIMIZE +// * MINOS +// * IMPROVE +// * SEEK +// +// The default (num==0 and tol==0) should always give you the +// corresponding defaults used in Minuit. +// +// +// FIXME: Implement changing cut in hadronness... +// FIXME: Show MHSignificance on MStatusDisplay during filling... +// FIXME: Choose step-size percentage as static data membewr +// FIXME: Choose minimization method +// +///////////////////////////////////////////////////////////////////////////// +#include "MJOptimizeEnergy.h" + +#include "MHMatrix.h" + +// environment +#include "MLog.h" +#include "MLogManip.h" + +// eventloop +#include "MParList.h" +#include "MTaskList.h" +#include "MEvtLoop.h" + +// parameters +#include "MParameters.h" +#include "MGeomCamMagic.h" + +// histograms +#include "../mhflux/MHEnergyEst.h" + +// tasks +#include "MReadTree.h" +#include "MMatrixLoop.h" +#include "MEnergyEstimate.h" +#include "MFillH.h" + +// filters +#include "MFDataMember.h" + +using namespace std; + +ClassImp(MJOptimizeEnergy); + +//------------------------------------------------------------------------ +// +// Read all events from file which do match rules and optimize +// energy estimator. +// +Bool_t MJOptimizeEnergy::RunEnergy(const char *fname, const char *rule) +{ + fLog->Separator("Preparing Energy optimization"); + + MParList parlist; + + MParameterI par("DataType"); + par.SetVal(1); + parlist.AddToList(&par); + + MFDataMember filter("DataType.fVal", '>', 0.5); + fPreCuts.Add(&filter); + + MGeomCamMagic geom; // For GetConvMm2Deg + parlist.AddToList(&geom); + + MHMatrix m("M"); + AddRulesToMatrix(m); + parlist.AddToList(&m); + + MHEnergyEst hist; + hist.InitMapping(&m); + + MEnergyEstimate est("MParameters"); + est.SetRule(rule); + parlist.AddToList(&est); + + MReadTree read("Events"); + // NECESSARY BECAUSE OF MDataFormula GetRules missing + read.DisableAutoScheme(); + if (fname) + read.AddFile(fname); + else + AddSequences(read, fNamesOn); + if (!FillMatrix(read, parlist, kTRUE)) + return kFALSE; + + fPreCuts.Remove(&filter); + + MTaskList tasklist; + parlist.Replace(&tasklist); + + MFillH fill(&hist); + + MMatrixLoop loop(&m); + + tasklist.AddToList(&loop); + tasklist.AddToList(&est); + tasklist.AddToList(&fill); + + // Optimize with the tasklist in this parameterlist + if (!Optimize(parlist)) + return kFALSE; + + // Print the result + hist.Print(); + + // Store result if requested + TObjArray cont; + cont.Add(&est); + return WriteContainer(cont, fNameOut); +} Index: trunk/MagicSoft/Mars/mjobs/MJOptimizeEnergy.h =================================================================== --- trunk/MagicSoft/Mars/mjobs/MJOptimizeEnergy.h (revision 7121) +++ trunk/MagicSoft/Mars/mjobs/MJOptimizeEnergy.h (revision 7121) @@ -0,0 +1,21 @@ +#ifndef MARS_MJOptimizeEnergy +#define MARS_MJOptimizeEnergy + +#ifndef MARS_MJOptimize +#include "MJOptimize.h" +#endif + +class MJOptimizeEnergy : public MJOptimize +{ +public: + // Special optimizing routines + Bool_t RunEnergy(const char *fname, const char *rule); + Bool_t RunEnergy(const char *rule) + { + return RunEnergy(0, rule); + } + + ClassDef(MJOptimizeEnergy, 0) // Class for optimization of the Supercuts +}; + +#endif Index: trunk/MagicSoft/Mars/mjobs/Makefile =================================================================== --- trunk/MagicSoft/Mars/mjobs/Makefile (revision 7120) +++ trunk/MagicSoft/Mars/mjobs/Makefile (revision 7121) @@ -36,4 +36,6 @@ MJCut.cc \ MJOptimize.cc \ + MJOptimizeCuts.cc \ + MJOptimizeEnergy.cc \ MJSpectrum.cc