Index: trunk/MagicSoft/Mars/macros/train/traindisp.C
===================================================================
--- trunk/MagicSoft/Mars/macros/train/traindisp.C	(revision 8634)
+++ trunk/MagicSoft/Mars/macros/train/traindisp.C	(revision 8706)
@@ -1,5 +1,54 @@
+/* ======================================================================== *\
+!
+! *
+! * 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, 11/2005 <mailto:tbretz@astro.uni-wuerzburg.de>
+!
+!   Copyright: MAGIC Software Development, 2000-2007
+!
+!
+\* ======================================================================== *//////////////////////////////////////////////////////////////////////////////
+
+/////////////////////////////////////////////////////////////////////////////
+//
+// traindisp.C
+// ===========
+//
+// Trains a random forest for disp estimation.
+//
+// The additional code below shows how the MagicCuts can be used in
+// training and testing. There is also code which allows to change the
+// slope of the input spectrum. If a min- or max-break energy is given
+// the new slope is only applied in this region. Further possibilities
+// for additional cuts are shown.
+//
+// SequencesOn:   Monte Carlo Sequences used for training
+// SequencesOff:  Monte Carlo Sequences used for testing
+//
+// Use:
+//  opt.AddPreCut    to use cut for test and training
+//  opt.AddTestCut   to use cut for test only
+//  opt.AddTrainCut  to use cut for train only
+//  opt.AddPreTask   a task executed before filling the matrix/evaluating the RF
+//  opt.AddPostTask  a task executed before training, after evaluating the RF
+//
+/////////////////////////////////////////////////////////////////////////////
+
 void traindisp()
 {
-    MDataSet set("mctesttrain.txt");
+    MDataSet set("mcsponde/mcdataset-for-training.txt");
     set.SetNumAnalysis(1);            // Necessary
 
@@ -10,4 +59,9 @@
     //opt.SetDebug();
 
+    // These are control parameters how to train the random forest (use with care!)
+    //opt.SetNumTrees(100);
+    //opt.SetNdSize(5);
+    //opt.SetNumTry(0);
+
     // ------- Parameters to train Random Forest --------
     opt.AddParameter("MHillas.fLength");
@@ -15,15 +69,16 @@
     opt.AddParameter("MHillas.fSize");
     opt.AddParameter("MNewImagePar.fLeakage1");
-    opt.AddParameter("MNewImagePar.fLeakage2");
-    opt.AddParameter("MNewImagePar.fConc");
-    opt.AddParameter("MNewImagePar.fConc1");
     opt.AddParameter("MPointingPos.fZd");
-    opt.AddParameter("MNewImagePar.fUsedArea");
-    opt.AddParameter("MNewImagePar.fCoreArea");
-    opt.AddParameter("MNewImagePar.fNumUsedPixels");
-    opt.AddParameter("MNewImagePar.fNumCorePixels");
-    opt.AddParameter("MImagePar.fSizeSinglePixels");
-    opt.AddParameter("MHillasExt.fM3Long*sign(MHillasSrc.fCosDeltaAlpha)");
-    opt.AddParameter("MHillasExt.fAsym*sign(MHillasSrc.fCosDeltaAlpha)");
+    opt.AddParameter("MHillasExt.fSlopeLong*sign(MHillasSrc.fCosDeltaAlpha)");
+    //opt.AddParameter("MNewImagePar.fLeakage2");
+    //opt.AddParameter("MNewImagePar.fConc");
+    //opt.AddParameter("MNewImagePar.fConc1");
+    //opt.AddParameter("MNewImagePar.fUsedArea");
+    //opt.AddParameter("MNewImagePar.fCoreArea");
+    //opt.AddParameter("MNewImagePar.fNumUsedPixels");
+    //opt.AddParameter("MNewImagePar.fNumCorePixels");
+    //opt.AddParameter("MImagePar.fSizeSinglePixels");
+    //opt.AddParameter("MHillasExt.fM3Long*sign(MHillasSrc.fCosDeltaAlpha)");
+    //opt.AddParameter("MHillasExt.fAsym*sign(MHillasSrc.fCosDeltaAlpha)");
 
     // -------------------- Run ----------------------------
@@ -38,13 +93,17 @@
      cuts.SetThetaCut(MFMagicCuts::kOn);
 
-     TArrayD arr(10);
-     arr[0]=  1.3245;
-     arr[1]=  0.208700;
-     arr[2]=  0.229200;
-     arr[3]=  5.305200;
-     arr[4]=  0.098930;
-     arr[5]= -0.082950;
-     arr[6]=  8.2957;
-     arr[7]=  0.8677;
+     TArrayD arr(13);
+     arr[0]  =  1.15136;
+     arr[1]  =  0.215;
+     arr[2]  =  0.215468;
+     arr[3]  =  5.63973;
+     arr[4]  =  0.0836169;
+     arr[5]  = -0.07;
+     arr[6]  =  7.2;
+     arr[7]  =  0.5;
+     arr[8]  =  0.0681437;
+     arr[9]  =  2.62932;
+     arr[10] =  1.51279;
+     arr[11] =  0.0507821;
 
      cuts.SetVariables(arr);
@@ -58,8 +117,46 @@
      -------------------- Other cuts ----------------------
      opt.AddPreCut("MNewImagePar.fLeakage1<0.0001");
-    */
-
      opt.AddPreCut("MHillas.fSize>200");
      opt.AddPreCut("MHillasSrc.fDCA*MGeomCam.fConvMm2Deg<0.3");
+     opt.AddPreCut("MPointingPos.fZd<25");
+
+     ------------------ Zd distribution -------------------
+     TFile file("ganymed00001111.root");
+
+     MStatusArray arr;
+     if (arr.Read()<=0)
+        return;
+     TH1D *vstime = (TH1D*)arr.FindObjectInCanvas("Theta",  "TH1D", "OnTime");
+     if (!vstime)
+         return;
+
+     MMcSpectrumWeight weight;
+     weight.SetWeightsZd(vstime);
+     opt.AddPreTask(&weight);
+
+     ---------------------- Histogram  --------------------
+
+     MHn hist("MyHist", "Energy Residual (lg E_{est} - lg E_{mc})");
+
+     MBinning bins(50, -0.5, 0.5);
+
+     hist.AddHist("MHillasExt.fM3Long*sign(MHillasSrc.fCosDeltaAlpha)*3.37e-3", "Disp.fVal-MHillasSrc.fDist*3.37e-3");
+     hist.InitName("ResM3l;M3Long;ResidualDist");
+     hist.InitTitle(";M3l [\\circ];Disp-Dist [\\circ];");
+     hist.SetDrawOption("colz profx");
+     hist.SetBinnings(&bins);
+
+     hist.AddHist("MHillasExt.fAsym*sign(MHillasSrc.fCosDeltaAlpha)*3.37e-3", "Disp.fVal-MHillasSrc.fDist*3.37e-3");
+     hist.InitName("ResAsym;M3Long;ResidualDist");
+     hist.InitTitle(";Asym [\\circ];Disp-Dist [\\circ];");
+     hist.SetDrawOption("colz profx");
+     hist.SetBinnings(&bins);
+
+     MFillH fill(&hist, "", "FillMyHist");
+     //fill.SetWeight(); // Enable weights to be used to fill this histogram
+     opt.AddTestTask(&fill);
+
+     ------------------------------------------------------
+    */
 
      // The number is the number of events read from the file
@@ -67,11 +164,2 @@
      opt.Train("rf-disp.root", set, 30000);
 }
-/*
- // SequencesOn:   Monte Carlo Sequences used for training
- // SequencesOff:  Monte Carlo Sequences used for testing
-
- // Use:
- // opt.AddPreCut    to use cut for test and training
- // opt.AddTestCut   to use cut for test only
- // opt.AddTrainCut  to use cut for train only
- */
Index: trunk/MagicSoft/Mars/macros/train/trainenergy.C
===================================================================
--- trunk/MagicSoft/Mars/macros/train/trainenergy.C	(revision 8634)
+++ trunk/MagicSoft/Mars/macros/train/trainenergy.C	(revision 8706)
@@ -18,5 +18,5 @@
 !   Author(s): Thomas Bretz, 11/2005 <mailto:tbretz@astro.uni-wuerzburg.de>
 !
-!   Copyright: MAGIC Software Development, 2000-2006
+!   Copyright: MAGIC Software Development, 2000-2007
 !
 !
@@ -36,8 +36,18 @@
 // for additional cuts are shown.
 //
+// SequencesOn:   Monte Carlo Sequences used for training
+// SequencesOff:  Monte Carlo Sequences used for testing
+//
+// Use:
+//  opt.AddPreCut    to use cut for test and training
+//  opt.AddTestCut   to use cut for test only
+//  opt.AddTrainCut  to use cut for train only
+//  opt.AddPreTask   a task executed before filling the matrix/evaluating the RF
+//  opt.AddPostTask  a task executed before training, after evaluating the RF
+//
 /////////////////////////////////////////////////////////////////////////////
 void trainenergy()
 {
-    MDataSet set("mcdataset.txt");
+    MDataSet set("mcsponde/mcdataset-for-training.txt");
     set.SetNumAnalysis(1);            // Necessary
 
@@ -48,15 +58,30 @@
     //opt.SetDebug();
 
+    // These are control parameters how to train the random forest (use with care!)
+    //opt.SetNumTrees(100);
+    //opt.SetNdSize(5);
+    //opt.SetNumTry(0);
+
     // ------- Parameters to train Random Forest --------
+    // The following parameters are the best parameters
     opt.AddParameter("MHillas.fSize");
     opt.AddParameter("MHillasSrc.fDist");
     opt.AddParameter("MPointingPos.fZd");
-    opt.AddParameter("MHillas.GetArea");
-    opt.AddParameter("MNewImagePar.fUsedArea");
-    opt.AddParameter("MNewImagePar.fCoreArea");
     opt.AddParameter("MNewImagePar.fLeakage1");
-    opt.AddParameter("MNewImagePar.fLeakage2");
-    opt.AddParameter("MNewImagePar.fConc");
-    opt.AddParameter("MNewImagePar.fConc1");
+    opt.AddParameter("MHillasExt.fSlopeLong*sign(MHillasSrc.fCosDeltaAlpha)");
+    // The influence of these parameters is unclear
+    //opt.AddParameter("MHillasExt.fM3Long*sign(MHillasSrc.fCosDeltaAlpha)");
+    //opt.AddParameter("MNewImagePar.fConcCOG");
+    //opt.AddParameter("MHillas.GetArea");
+    //opt.AddParameter("MNewImagePar.fConc1");
+    //opt.AddParameter("MNewImagePar.fConc");
+    //opt.AddParameter("MNewImagePar.fUsedArea");
+    //opt.AddParameter("MNewImagePar.fCoreArea");
+    //opt.AddParameter("MNewImagePar.fLeakage2");
+
+    // Setup how to train the RF. This one gives best results so far
+    opt.SetTrainFunc("log(MMcEvt.fEnergy)/log(MHillas.fSize)", "MHillas.fSize^x");
+    // opt.SetTrainLog(); // Train in log-energy
+    // opt.SetTrainLin(); // Train just in energy
 
     // -------------------- Run ----------------------------
@@ -65,50 +90,78 @@
     opt.SetDisplay(d);
 
+    // -------------------- Magic-Cuts ----------------------
+    // It is recommended to test with your cuts, but train with all events
+
+    MFMagicCuts cuts;
+    cuts.SetHadronnessCut(MFMagicCuts::kArea);
+    cuts.SetThetaCut(MFMagicCuts::kOn);
+
+    TArrayD arr(13);
+    arr[0]  =  1.15136;
+    arr[1]  =  0.215;
+    arr[2]  =  0.215468;
+    arr[3]  =  5.63973;
+    arr[4]  =  0.0836169;
+    arr[5]  = -0.07;
+    arr[6]  =  7.2;
+    arr[7]  =  0.5;
+    arr[8]  =  0.0681437;
+    arr[9]  =  2.62932;
+    arr[10] =  1.51279;
+    arr[11] =  0.0507821;
+
+    cuts.SetVariables(arr);
+
+    // opt.AddPreCut(&cuts);    // add cut for training and testing
+    // opt.AddTrainCut(&cuts);  // add cut for training only
+    opt.AddTestCut(&cuts);   // add cut for testing only
+
     /*
-     -------------------- Magic-Cuts ----------------------
-     MFMagicCuts cuts;
-     cuts.SetHadronnessCut(MFMagicCuts::kArea);
-     cuts.SetThetaCut(MFMagicCuts::kOn);
+     -------------------- Energy Slope --------------------
+     // This is a way to throw away events to a different slope
+     MFEnergySlope slope(-4.0); // New slope for mc spectrum
+     opt.AddPreCut(&slope);     // throw away events to change slope
 
-     TArrayD arr(10);
-     arr[0]=  1.3245;
-     arr[1]=  0.208700;
-     arr[2]=  0.229200;
-     arr[3]=  5.305200;
-     arr[4]=  0.098930;
-     arr[5]= -0.082950;
-     arr[6]=  8.2957;
-     arr[7]=  0.8677;
-
-     cuts.SetVariables(arr);
-
-     opt.AddPreCut(&cuts);
-
-     -------------------- Energy Slope --------------------
-     MFEnergySlope slope(-2.8); // New slope for mc spectrum
-     slope.SetMcMinEnergy(80);  // Set break energy from -2.6 to -2.8
-     opt.AddPreCut(&slope);     // throw away events to change slope
+     // This is a way to weight the events to a different spectrum
+     MMcSpectrumWeight weight;
+     weight.SetFormula("pow(X/300, -2.31-0.26*log10(X/300))");
+     opt.SetWeights(&weight);
 
      -------------------- Other cuts ----------------------
      opt.AddPreCut("MHillasSrc.fDist*MGeomCam.fConvMm2Deg<1.0");
      opt.AddPreCut("MHillas.fSize>200");
+
+     ------------------ Zd distribution -------------------
+     TFile file("ganymed00001111.root");
+
+     MStatusArray arr;
+     if (arr.Read()<=0)
+        return;
+     TH1D *vstime = (TH1D*)arr.FindObjectInCanvas("Theta",  "TH1D", "OnTime");
+     if (!vstime)
+         return -1;
+
+     MMcSpectrumWeight weight;
+     weight.SetWeightsZd(vstime);
+     opt.AddPreTask(&weight);
+
+     ---------------------- Histogram  --------------------
+     MHn hist("MyHist", "Energy Residual (lg E_{est} - lg E_{mc})");
+
+     hist.AddHist("MNewImagePar.fConcCOG", "log10(MEnergyEst.fVal)-log10(MMcEvt.fEnergy)");
+     hist.InitName("ResConc;Conc;EnergyResidual");
+     hist.InitTitle(";C;\\Delta lg E;");
+     hist.SetDrawOption("colz profx");
+
+     MBinning bins(50, 0, 1);
+     hist.SetBinnings(&bins);
+
+     MFillH fill(&hist, "", "FillMyHist");
+     //fill.SetWeight(); // Enable weights to be used to fill this histogram
+     opt.AddTestTask(&fill);
+
+     ------------------------------------------------------
     */
-
-    // Things which could be implemented in a future version
-    // PreCut in Energy oder Size?
-    // PreCut in Max Dist
-    // PreCut in Theta^2 (0.25 oder 0.3)
-    // New spectral slope
-    // Zenith angle distribution in OnTime
 
      opt.Train("rf-energy.root", set, 30000);
 }
-/*
- // SequencesOn:   Monte Carlo Sequences used for training
- // SequencesOff:  Monte Carlo Sequences used for testing
-
- // Use:
- // opt.AddPreCut    to use cut for test and training
- // opt.AddTestCut   to use cut for test only
- // opt.AddTrainCut  to use cut for train only
- */
