/* ======================================================================== *\
!
! *
! * 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, 5/2002 <mailto:tbretz@astro.uni-wuerzburg.de>
!   Author(s): Rudy Bock, 5/2002 <mailto:rkb@mppmu.mpg.de>
!
!   Copyright: MAGIC Software Development, 2000-2002
!
!
\* ======================================================================== */

void multidimdist()
{
    //
    // This is a demonstration program which calculates the Hillas
    // parameter out of a Magic root file.

    //
    // Create a empty Parameter List and an empty Task List
    // The tasklist is identified in the eventloop by its name
    //
    MParList  plist;

    MTaskList tlist;
    plist.AddToList(&tlist);

    //
    // Now setup the tasks and tasklist:
    // ---------------------------------
    //
    // The first argument is the tree you want to read.
    //   Events:     Cosmic ray events
    //   PedEvents:  Pedestal Events
    //   CalEvents:  Calibration Events
    //
    MReadMarsFile  read("Events", "~/Mars/star.root");
    read.DisableAutoScheme();
    tlist.AddToList(&read);

    MFParticleId fgamma("MMcEvt", '=', kGAMMA);
    tlist.AddToList(&fgamma);

    MFParticleId fhadrons("MMcEvt", '!', kGAMMA);
    tlist.AddToList(&fhadrons);

    MHMatrix matrix("MatrixGammas");
    matrix.AddColumn("MHillas.fWidth");
    matrix.AddColumn("MHillas.fLength");
    matrix.AddColumn("MHillas.fWidth*MHillas.fLength/MHillas.fSize");
    matrix.AddColumn("abs(MHillas.fAsym)");
    matrix.AddColumn("abs(MHillas.fM3Long)");
    matrix.AddColumn("abs(MHillas.fM3Trans)");
    matrix.AddColumn("abs(HillasSource.fHeadTail)");
    matrix.AddColumn("MHillas.fConc");
    matrix.AddColumn("MHillas.fConc1");
    matrix.AddColumn("HillasSource.fDist");
    matrix.AddColumn("log10(MHillas.fSize)");
    plist.AddToList(&matrix);

    MHMatrix matrix2("MatrixHadrons");
    matrix2.AddColumns(matrix.GetColumns());
    plist.AddToList(&matrix2);

    MFillH fillmat("MatrixGammas");
    fillmat.SetFilter(&fgamma);
    tlist.AddToList(&fillmat);

    MFillH fillmat2("MatrixHadrons");
    fillmat2.SetFilter(&fhadrons);
    tlist.AddToList(&fillmat2);

    //
    // Create and setup the eventloop
    //
    MEvtLoop evtloop;
    evtloop.SetParList(&plist);

    //
    // Execute your analysis
    //
    if (!evtloop.Eventloop())
        return;

    tlist.PrintStatistics();

    matrix.Print("size");
    matrix2.Print("size");

    // ---------------------------------------------------------

    MTaskList tlist2;

    plist.Replace(&tlist2);

    MReadMarsFile read2("Events", "~/Mars/star.root");
//    read2.AddFile("~/Mars/star2.root");
    read2.DisableAutoScheme();
    tlist2.AddToList(&read2);

    MMultiDimDistCalc calc;
    calc.SetUseNumRows(15);
    //calc.SetUseKernelMethod(kTRUE);
    tlist2.AddToList(&calc);

    MFillH fillh("MHHadronness");

    /*
     MF filter("MMcEvt.fEnergy < 100");
     fillh.SetFilter(&filter);
     tlist2.AddToList(&filter);
     */

    tlist2.AddToList(&fillh);

    //
    // Execute your analysis
    //
    if (!evtloop.Eventloop())
        return;

    tlist2.PrintStatistics();

    plist.FindObject("MHHadronness")->DrawClone();
    plist.FindObject("MHHadronness")->Print();
}