source: trunk/MagicSoft/Mars/mjobs/MJOptimizeEnergy.cc@ 7126

/* ======================================================================== *\
!
! *
! * 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 <mailto:tbretz@astro.uni-wuerzburg.de>
!
!   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);
}