Index: trunk/MagicSoft/Mars/macros/ONOFFCT1Analysis.C
===================================================================
--- trunk/MagicSoft/Mars/macros/ONOFFCT1Analysis.C	(revision 2167)
+++ trunk/MagicSoft/Mars/macros/ONOFFCT1Analysis.C	(revision 2169)
@@ -153,25 +153,15 @@
     //  - write root file for ON (or OFF or MC) data (ON1.root, ...);
 
-    Bool_t JobA    = kTRUE;  
+    Bool_t JobA    = kFALSE;  
     Bool_t WPad    = kFALSE;   // write out padding histograms ?
-    Bool_t RPad    = kTRUE;   // read in padding histograms ?
-    Bool_t Wout    = kTRUE;   // write out root file ON1.root 
+    Bool_t RPad    = kFALSE;   // read in padding histograms ?
+    Bool_t Wout    = kFALSE;   // write out root file ON1.root 
                                // (or OFF1.root or MC1.root)?
 
-
-    // Job A_MC : read MC gamma data, 
-    //  - read sigmabar vs. Theta plot from ON data  
-    //  - do padding; 
-    //  - write root file for MC gammas (MC1.root);
-
-    Bool_t JobA_MC  = kFALSE;  
-    Bool_t WMC1     = kFALSE;  // write out root file MC1.root ?
-
-
-    // Job B_NN_UP : read ON1.root (or MC1.root) file 
+    // Job B_NN_UP : read ON1.root (OFF1.root or MC1.root) file 
     //  - depending on RlookNN : create (or read in) hadron and gamma matrix
     //  - calculate hadroness for method of NEAREST NEIGHBORS
-    //    and for the SUPERCUTS
-    //  - update the input files with the hadronesses (ON1.root or MC1.root)
+    //  - update the input files with the hadroness (ON1.root, OFF1.root
+    //     or MC1.root)
 
     Bool_t JobB_NN_UP  = kFALSE;  
@@ -180,9 +170,11 @@
 
 
-    // Job B_RF_UP : read ON1.root (or MC1.root) file 
+    // Job B_RF_UP : read ON1.root (OFF1.root or MC1.root) file 
     //  - depending on RLook : create (or read in) hadron and gamma matrix
     //  - depending on RTree : create (or read in) trees
     //  - calculate hadroness for method of RANDOM FOREST
-    //  - update the input files with the hadroness (ON1.root or MC1.root)
+    //    and for the SUPERCUTS
+    //  - update the input files with the hadronesses (ON1.root, OFF1.root
+    //     or MC1.root)
 
     Bool_t JobB_RF_UP  = kFALSE;  
@@ -211,5 +203,5 @@
     //  - make plots
 
-    Bool_t JobD  = kFALSE;  
+    Bool_t JobD  = kTRUE;  
 
 
@@ -220,6 +212,6 @@
     //  - optimize energy estimation for events passing the final cuts
     //  - write parameters of energy estimator onto file
-    //  - update ON1.root and MC1.root files with estimated energy
-    //    (ON_XX1.root and MC_XX1.root)
+    //  - update ON1.root, OFF1.root and MC1.root files with estimated energy
+    //    (ON_XX1.root, OFF_XX1.root and MC_XX1.root)
 
     Bool_t JobE_EST_UP  = kFALSE;  
@@ -693,295 +685,4 @@
 
 
-  //---------------------------------------------------------------------
-   // Job A_MC
-   //=========
-
-    // read MC gamma data  
-
-    //    - to pad them
-    //      (using the 2D-histogram "sigmabar versus Theta" 
-    //       (SigmaTheta_ON.root)  of the ON data)
-
-    //    - to write a file of padded MC gamma events (MC1.root)
-    //      (after the standard cuts, before the g/h separation)
-    //      (to be used together with the corresponding hadron file
-    //       for the optimization of the g/h separation)
-
-
- if (JobA_MC)
- {
-    gLog << "=====================================================" << endl;
-    gLog << "Macro CT1Analysis : Start of Job A_MC" << endl;
-
-    gLog << "" << endl;
-    gLog << "Macro CT1Analysis : JobA_MC, WMC1 = " 
-         << JobA_MC  << ",  " << WMC1 << endl;
-
-
-    // name of input root file
-    TString filenamein(mcfile);
-
-    // name of output root file
-    TString outNameImage = outPath;
-    outNameImage += "MC";
-    outNameImage += "1.root";
-
-    //------------------------------------------------
-    // use for padding sigmabar vs. Theta from ON data
-    TString typeHist = "ON";
-    gLog << "typeHist = " << typeHist << endl;
-
-    // name of file containing the histograms for the padding
-    TString outNameSigTh = outPath;
-    outNameSigTh += "SigmaTheta_";
-    outNameSigTh += typeHist;
-    outNameSigTh += ".root";
-
-
-    //------------------------------------
-    // Get the histograms "2D-ThetaSigmabar"
-    // and                "3D-ThetaPixSigma"
-    // and                "3D-ThetaPixDiff"
-    // and                "2D-IdBlindPixels"
-    // and                "2D-NBlindPixels"
-
-
-      gLog << "Reading in file " << outNameSigTh << endl;
-
-      TFile *infile = new TFile(outNameSigTh);
-      infile->ls();
-
-      TH2D *fHSigmaTheta = 
-      (TH2D*) gROOT->FindObject("2D-ThetaSigmabar");
-      if (!fHSigmaTheta)
-	{
-          gLog << "Object '2D-ThetaSigmabar' not found on root file" << endl;
-          return;
-	}
-      gLog << "Object '2D-ThetaSigmabar' was read in" << endl;
-
-      TH3D *fHSigmaPixTheta = 
-      (TH3D*) gROOT->FindObject("3D-ThetaPixSigma");
-      if (!fHSigmaPixTheta)
-	{
-          gLog << "Object '3D-ThetaPixSigma' not found on root file" << endl;
-          return;
-	}
-      gLog << "Object '3D-ThetaPixSigma' was read in" << endl;
-
-      TH3D *fHDiffPixTheta = 
-      (TH3D*) gROOT->FindObject("3D-ThetaPixDiff");
-      if (!fHDiffPixTheta)
-	{
-          gLog << "Object '3D-ThetaPixDiff' not found on root file" << endl;
-          return;
-	}
-      gLog << "Object '3D-ThetaPixDiff' was read in" << endl;
-
-
-      TH2D *fHIdBlindPixels = 
-      (TH2D*) gROOT->FindObject("2D-IdBlindPixels");
-      if (!fHIdBlindPixels)
-	{
-          gLog << "Object '2D-IdBlindPixels' not found on root file" << endl;
-          return;
-	}
-      gLog << "Object '2D-IdBlindPixels' was read in" << endl;
-
-      TH2D *fHNBlindPixels = 
-      (TH2D*) gROOT->FindObject("2D-NBlindPixels");
-      if (!fHNBlindPixels)
-	{
-          gLog << "Object '2D-NBlindPixels' not found on root file" << endl;
-          return;
-	}
-      gLog << "Object '2D-NBlindPixels' was read in" << endl;
-
-    //------------------------------------
-
-    MTaskList tlist;
-    MParList plist;
-
-    char *sourceName = "MSrcPosCam";
-    MSrcPosCam source(sourceName);
-
-
-    // geometry is needed in  MHHillas... classes 
-    MGeomCam *fGeom = 
-             (MGeomCam*)plist->FindCreateObj("MGeomCamCT1", "MGeomCam");
-
-    //-------------------------------------------
-    // create the tasks which should be executed 
-    //
-
-    MCT1ReadPreProc read(filenamein);
-
-    MBlindPixelCalc blind;
-    blind.SetUseBlindPixels();
-
-    MFCT1SelBasic selbasic;
-    MContinue contbasic(&selbasic);
-    contbasic.SetName("SelBasic");
-
-
-    // There are 2 options for Thomas Schweizer's padding
-    //     fPadFlag = 1   get Sigmabar from fHSigmaTheta
-    //                    and Sigma    from fHDiffPixTheta
-    //     fPadFlag = 2   get Sigma    from fHSigmaPixTheta
-    
-    MPadSchweizer padthomas("MPadSchweizer","Task for the padding (Schweizer)");
-    padthomas.SetHistograms(fHSigmaTheta, fHSigmaPixTheta, fHDiffPixTheta,
-                            fHIdBlindPixels, fHNBlindPixels);
-    padthomas.SetPadFlag(1);
-
-    MFillH fillblind("MCBlindPixels[MHBlindPixels]", "MBlindPixels");
-    fillblind.SetName("HBlind");
-
-
-    //...........................................
-
-    MSigmabarCalc sigbarcalc;
-
-    MFillH fillsigtheta ("MCSigmaTheta[MHSigmaTheta]", "MMcEvt");
-    fillsigtheta.SetName("HSigmaTheta");
-
-    MImgCleanStd    clean; 
-
-    // calculation of  image parameters ---------------------
-    TString fHilName    = "MHillas";
-    TString fHilNameExt = "MHillasExt";
-    TString fHilNameSrc = "MHillasSrc";
-    TString fImgParName = "MNewImagePar";
-
-    MHillasCalc    hcalc;
-    hcalc.SetNameHillas(fHilName);
-    hcalc.SetNameHillasExt(fHilNameExt);
-    hcalc.SetNameNewImgPar(fImgParName);
-
-    MHillasSrcCalc hsrccalc(sourceName, fHilNameSrc);
-    hsrccalc.SetInput(fHilName);
-
-
-    MFillH hfill1("MHHillas",    fHilName);
-    hfill1.SetName("HHillas");
-
-    MFillH hfill2("MHStarMap",   fHilName);
-    hfill2.SetName("HStarMap");
-
-    MFillH hfill3("MHHillasExt",   fHilNameSrc);
-    hfill3.SetName("HHillasExt");
-
-    MFillH hfill4("MHHillasSrc",   fHilNameSrc);
-    hfill4.SetName("HHillasSrc");
-
-    MFillH hfill5("MHNewImagePar", fImgParName);
-    hfill5.SetName("HNewImagePar");
-    // --------------------------------------------------
-
-    MFCT1SelStandard selstandard(fHilNameSrc);
-    selstandard.SetHillasName(fHilName);
-    selstandard.SetImgParName(fImgParName);
-    selstandard.SetCuts(92, 4, 60, 0.4, 1.05, 0.0, 0.0);
-    MContinue contstandard(&selstandard);
-    contstandard.SetName("SelStandard");
-
-
-
-      MWriteRootFile write(outNameImage);
-
-      write.AddContainer("MRawRunHeader", "RunHeaders");
-      write.AddContainer("MTime",         "Events");
-      write.AddContainer("MMcEvt",        "Events");
-      write.AddContainer("ThetaOrig",     "Events");
-      write.AddContainer("MSrcPosCam",    "Events");
-      write.AddContainer("MSigmabar",     "Events");
-      write.AddContainer("MHillas",       "Events");
-      write.AddContainer("MHillasExt",    "Events");
-      write.AddContainer("MHillasSrc",    "Events");
-      write.AddContainer("MNewImagePar",  "Events");
-
-
-
-    //*****************************
-    // entries in MParList
-
-    plist.AddToList(&tlist);
-    InitBinnings(&plist);
-
-    plist.AddToList(&source);
-
-
-    //*****************************
-    // entries in MTaskList
-
-    tlist.AddToList(&read);
-    tlist.AddToList(&padthomas);
-    tlist.AddToList(&blind);
-
-    tlist.AddToList(&contbasic);
-    tlist.AddToList(&fillblind);
-    tlist.AddToList(&sigbarcalc);
-    tlist.AddToList(&fillsigtheta);
-    tlist.AddToList(&clean);
-
-    tlist.AddToList(&hcalc);
-    tlist.AddToList(&hsrccalc);
-
-    tlist.AddToList(&hfill1);
-    tlist.AddToList(&hfill2);
-    tlist.AddToList(&hfill3);
-    tlist.AddToList(&hfill4);
-    tlist.AddToList(&hfill5);
-
-    tlist.AddToList(&contstandard);
-    if (WMC1)
-      tlist.AddToList(&write);
-
-    //*****************************
-
-
-    //-------------------------------------------
-    // Execute event loop
-    //
-    MProgressBar bar;
-    MEvtLoop evtloop;
-    evtloop.SetParList(&plist);
-    evtloop.ReadEnv(env, "", printEnv);
-    evtloop.SetProgressBar(&bar);
-    //if (WMC1)    
-    //  evtloop.Write();
-
-    Int_t maxevents = -1;
-    //Int_t maxevents = 1000;
-    if ( !evtloop.Eventloop(maxevents) )
-        return;
-
-    tlist.PrintStatistics(0, kTRUE);
-
-
-    //-------------------------------------------
-    // Display the histograms
-    //
-
-    plist.FindObject("MCSigmaTheta",  "MHSigmaTheta")->DrawClone();
-    plist.FindObject("MCBlindPixels", "MHBlindPixels")->DrawClone();
-
-    plist.FindObject("MHHillas")->DrawClone();
-    plist.FindObject("MHHillasExt")->DrawClone();
-    plist.FindObject("MHHillasSrc")->DrawClone();
-    plist.FindObject("MHNewImagePar")->DrawClone();
-    plist.FindObject("MHStarMap")->DrawClone();
-
-
-
-    DeleteBinnings(&plist);
-
-    gLog << "Macro CT1Analysis : End of Job A_MC" 
-         << endl;
-    gLog << "=========================================================" 
-         << endl;
- }
-
-
 
   //---------------------------------------------------------------------
@@ -1010,8 +711,9 @@
 
     //--------------------------------------------
-    // file to be updated (either ON or MC)
-
-    TString typeInput = "ON";
-    //TString typeInput = "MC";
+    // file to be updated (ON, OFF or MC)
+
+    //TString typeInput = "ON";
+    //TString typeInput = "OFF";
+    TString typeInput = "MC";
     gLog << "typeInput = " << typeInput << endl;
 
@@ -1034,9 +736,9 @@
     // files to be read for generating the look-alike events
     // "hadrons" :
-    TString filenameON = outPath;
-    filenameON += "ON";
-    filenameON += "1.root";
+    TString filenameHad = outPath;
+    filenameHad += "OFF";
+    filenameHad += "1.root";
     Int_t howManyHadrons = 500000;
-    gLog << "filenameON = " << filenameON << ",   howManyHadrons = "
+    gLog << "filenameHad = " << filenameHad << ",   howManyHadrons = "
          << howManyHadrons  << endl; 
     
@@ -1058,5 +760,5 @@
     outNameGammas += ".root";
 
-    TString typeMatrixHadrons = "ON";
+    TString typeMatrixHadrons = "OFF";
     gLog << "typeMatrixHadrons = " << typeMatrixHadrons << endl;
 
@@ -1134,5 +836,5 @@
     readg.DisableAutoScheme();
 
-    MReadMarsFile  readh("Events", filenameON);
+    MReadMarsFile  readh("Events", filenameHad);
     readh.DisableAutoScheme();
 
@@ -1383,12 +1085,4 @@
 
     //.......................................................................
-    // calculation of hadroness for the supercuts
-    // (=0.25 if fullfilled, =0.75 otherwise)
-
-    TString hadSCName = "HadSC";
-    MCT1SupercutsCalc sccalc(fHilName, fHilNameSrc);
-    sccalc.SetHadronnessName(hadSCName);
-
-    //.......................................................................
 
       //MWriteRootFile write(outNameImage, "UPDATE");
@@ -1406,6 +1100,7 @@
       write.AddContainer("MNewImagePar",  "Events");
 
+      write.AddContainer("HadRF",         "Events");
+      write.AddContainer("HadSC",         "Events");
       write.AddContainer("HadNN",         "Events");
-      write.AddContainer("HadSC",         "Events");
 
 
@@ -1428,6 +1123,4 @@
     MFillH fillhadnn("hadNN[MHHadronness]", hadNNName);
     fillhadnn.SetName("HhadNN");
-    MFillH fillhadsc("hadSC[MHHadronness]", hadSCName);
-    fillhadsc.SetName("HhadSC");
 
     MFCT1SelFinal selfinal(fHilNameSrc);
@@ -1470,7 +1163,5 @@
 
     tliston.AddToList(&nncalc);
-    tliston.AddToList(&sccalc);
     tliston.AddToList(&fillhadnn);
-    tliston.AddToList(&fillhadsc);
 
     tliston.AddToList(&hfill1);
@@ -1505,5 +1196,4 @@
     //
     pliston.FindObject("hadNN", "MHHadronness")->DrawClone();
-    pliston.FindObject("hadSC", "MHHadronness")->DrawClone();
 
     pliston.FindObject("MHHillas")->DrawClone();
@@ -1557,7 +1247,8 @@
 
     //--------------------------------------------
-    // file to be updated (either ON or MC)
+    // file to be updated (ON, OFF or MC)
 
     TString typeInput = "ON";
+    //TString typeInput = "OFF";
     //TString typeInput = "MC";
     gLog << "typeInput = " << typeInput << endl;
@@ -1566,5 +1257,5 @@
     TString filenameData = outPath;
     filenameData += typeInput;
-    filenameData += "2.root";
+    filenameData += "1.root";
     gLog << "filenameData = " << filenameData << endl; 
 
@@ -1572,5 +1263,5 @@
     TString outNameImage = outPath;
     outNameImage += typeInput;
-    outNameImage += "3.root";
+    outNameImage += "2.root";
     //TString outNameImage = filenameData;
 
@@ -1580,9 +1271,9 @@
     // files to be read for generating the look-alike events
     // "hadrons" :
-    TString filenameON = outPath;
-    filenameON += "ON";
-    filenameON += "1.root";
+    TString filenameHad = outPath;
+    filenameHad += "OFF";
+    filenameHad += "1.root";
     Int_t howManyHadrons = 1000000;
-    gLog << "filenameON = " << filenameON << ",   howManyHadrons = "
+    gLog << "filenameHad = " << filenameHad << ",   howManyHadrons = "
          << howManyHadrons  << endl; 
     
@@ -1604,5 +1295,5 @@
     outNameGammas += ".root";
 
-    TString typeMatrixHadrons = "ON";
+    TString typeMatrixHadrons = "OFF";
     gLog << "typeMatrixHadrons = " << typeMatrixHadrons << endl;
 
@@ -1687,5 +1378,5 @@
     readg.DisableAutoScheme();
 
-    MReadMarsFile  readh("Events", filenameON);
+    MReadMarsFile  readh("Events", filenameHad);
     readh.DisableAutoScheme();
 
@@ -2046,4 +1737,13 @@
 
     //.......................................................................
+    // calculation of hadroness for the supercuts
+    // (=0.25 if fullfilled, =0.75 otherwise)
+
+    TString hadSCName = "HadSC";
+    MCT1SupercutsCalc sccalc(fHilName, fHilNameSrc);
+    sccalc.SetHadronnessName(hadSCName);
+
+
+    //.......................................................................
 
       //MWriteRootFile write(outNameImage, "UPDATE");
@@ -2061,9 +1761,6 @@
       write.AddContainer("MNewImagePar",  "Events");
 
-      write.AddContainer("HadNN",         "Events");
+      write.AddContainer("HadRF",         "Events");
       write.AddContainer("HadSC",         "Events");
-
-      write.AddContainer(hadRFName,       "Events");
-
 
     //-----------------------------------------------------------------
@@ -2089,4 +1786,7 @@
     MFillH fillhadrf("hadRF[MHHadronness]", hadRFName);
     fillhadrf.SetName("HhadRF");
+    MFillH fillhadsc("hadSC[MHHadronness]", hadSCName);
+    fillhadsc.SetName("HhadSC");
+
 
     MFCT1SelFinal selfinal(fHilNameSrc);
@@ -2128,6 +1828,8 @@
 
     tliston.AddToList(&rfcalc);
+    tliston.AddToList(&sccalc);
     tliston.AddToList(&fillranfor);
     tliston.AddToList(&fillhadrf);
+    tliston.AddToList(&fillhadsc);
 
     tliston.AddToList(&hfill1);
@@ -2164,4 +1866,5 @@
     pliston.FindObject("MHRanForest")->DrawClone();
     pliston.FindObject("hadRF", "MHHadronness")->DrawClone();
+    pliston.FindObject("hadSC", "MHHadronness")->DrawClone();
 
     pliston.FindObject("MHHillas")->DrawClone();
@@ -2201,6 +1904,6 @@
     // name of input data file
     TString filenameData = outPath;
-    filenameData += "ON";
-    filenameData += "3.root";
+    filenameData += "OFF";
+    filenameData += "2.root";
     gLog << "filenameData = " << filenameData << endl;
 
@@ -2208,5 +1911,5 @@
     TString filenameMC = outPath;
     filenameMC += "MC";
-    filenameMC += "3.root";
+    filenameMC += "2.root";
     gLog << "filenameMC   = " << filenameMC   << endl;
 
@@ -2252,11 +1955,11 @@
     MFCT1SelFinal selfinalgh(fHilNameSrc);
     selfinalgh.SetCuts(maxhadronness, 100.0, maxdist);
-    selfinalgh.SetHadronnessName(hadNNName);
+    selfinalgh.SetHadronnessName(hadRFName);
     selfinalgh.SetName("SelFinalgh");
     MContinue contfinalgh(&selfinalgh);
     contfinalgh.SetName("ContSelFinalgh");
 
-    MFillH fillhadnn("hadNN[MHHadronness]", hadNNName);
-    fillhadnn.SetName("HhadNN");
+    //MFillH fillhadnn("hadNN[MHHadronness]", hadNNName);
+    //fillhadnn.SetName("HhadNN");
     MFillH fillhadsc("hadSC[MHHadronness]", hadSCName);
     fillhadsc.SetName("HhadSC");
@@ -2266,5 +1969,5 @@
     MFCT1SelFinal selfinal(fHilNameSrc);
     selfinal.SetCuts(maxhadronness, maxalpha, maxdist);
-    selfinal.SetHadronnessName(hadNNName);
+    selfinal.SetHadronnessName(hadRFName);
     selfinal.SetName("SelFinal");
     MContinue contfinal(&selfinal);
@@ -2300,5 +2003,5 @@
     tliston.AddToList(&read);
 
-    tliston.AddToList(&fillhadnn);
+    //tliston.AddToList(&fillhadnn);
     tliston.AddToList(&fillhadsc);
     tliston.AddToList(&fillhadrf);
@@ -2335,5 +2038,5 @@
     //
 
-    pliston.FindObject("hadNN", "MHHadronness")->DrawClone();
+    //pliston.FindObject("hadNN", "MHHadronness")->DrawClone();
     pliston.FindObject("hadSC", "MHHadronness")->DrawClone();
     pliston.FindObject("hadRF", "MHHadronness")->DrawClone();
@@ -2372,10 +2075,10 @@
 
     // type of data to be analysed
-    //TString typeData = "ON";
+    TString typeData = "ON";
     //TString typeData = "OFF";
-    TString typeData = "MC";
+    //TString typeData = "MC";
     gLog << "typeData = " << typeData << endl;
 
-    TString ext      = "3.root";
+    TString ext      = "2.root";
 
 
@@ -2392,5 +2095,5 @@
 
     // maximum values of the hadronness, |ALPHA| and DIST
-    Float_t maxhadronness   = 0.20;
+    Float_t maxhadronness   = 0.30;
     Float_t maxalpha        = 20.0;
     Float_t maxdist         = 10.0;
@@ -2447,6 +2150,6 @@
     contfinalgh.SetName("ContSelFinalgh");
 
-    MFillH fillhadnn("hadNN[MHHadronness]", "HadNN");
-    fillhadnn.SetName("HhadNN");
+    //MFillH fillhadnn("hadNN[MHHadronness]", "HadNN");
+    //fillhadnn.SetName("HhadNN");
     MFillH fillhadsc("hadSC[MHHadronness]", "HadSC");
     fillhadsc.SetName("HhadSC");
@@ -2492,7 +2195,5 @@
     tliston.AddToList(&contfinalgh);
 
-    tliston.AddToList(&contfinal);
-
-    tliston.AddToList(&fillhadnn);
+    //tliston.AddToList(&fillhadnn);
     tliston.AddToList(&fillhadsc);
     tliston.AddToList(&fillhadrf);
@@ -2504,5 +2205,5 @@
     tliston.AddToList(&hfill5);
 
-
+    tliston.AddToList(&contfinal);
 
     //*****************************
@@ -2524,9 +2225,10 @@
 
 
+
     //-------------------------------------------
     // Display the histograms
     //
 
-    pliston.FindObject("hadNN", "MHHadronness")->DrawClone();
+    //pliston.FindObject("hadNN", "MHHadronness")->DrawClone();
     pliston.FindObject("hadRF", "MHHadronness")->DrawClone();
     pliston.FindObject("hadSC", "MHHadronness")->DrawClone();
@@ -2538,4 +2240,16 @@
     pliston.FindObject("MHStarMap")->DrawClone();
 
+    //-------------------------------------------
+    // fit alpha distribution to get the number of excess events
+    //
+
+    MHHillasSrc* hillasSrc = 
+      (MHHillasSrc*)(pliston->FindObject("MHHillasSrc"));
+    TH1F* alphaHist = (TH1F*)(hillasSrc->GetHistAlpha());
+  
+    MHOnSubtraction onsub;
+    onsub.Calc(alphaHist, &pliston, kTRUE, 25);
+    //-------------------------------------------
+
     DeleteBinnings(&pliston);
 
@@ -2570,8 +2284,9 @@
 
 
-    TString typeON = "ON";
-    TString typeMC = "MC";
-    TString ext    = "3.root";
-    TString extout = "4.root";
+    TString typeON  = "ON";
+    TString typeOFF = "OFF";
+    TString typeMC  = "MC";
+    TString ext    = "2.root";
+    TString extout = "3.root";
 
     //------------------------------
@@ -2586,7 +2301,7 @@
     // and definition of final selections
 
-    TString XX("NN");
+    //TString XX("NN");
     //TString XX("SC");
-    //TString XX("RF");
+    TString XX("RF");
     TString fhadronnessName("Had");
     fhadronnessName += XX;
@@ -2615,4 +2330,10 @@
     filenameON += ext;
     gLog << "filenameON = " << filenameON << endl;
+
+    // name of OFF file to be updated
+    TString filenameOFF(outPath);
+    filenameOFF += typeOFF;
+    filenameOFF += ext;
+    gLog << "filenameOFF = " << filenameOFF << endl;
 
     // name of MC file to be updated
@@ -2631,4 +2352,12 @@
     gLog << "filenameONup = " << filenameONup << endl;
 
+    // name of updated OFF file 
+    TString filenameOFFup(outPath);
+    filenameOFFup += typeOFF;
+    filenameOFFup += "_";
+    filenameOFFup += XX;
+    filenameOFFup += extout;
+    gLog << "filenameOFFup = " << filenameOFFup << endl;
+
     // name of updated MC file 
     TString filenameMCup(outPath);
@@ -2683,21 +2412,21 @@
     //==========   start update   ============================================
     //
-    // Update ON and MC root files with the estimated energy
+    // Update ON, OFF and MC root files with the estimated energy
 
     //---------------------------------------------------
-    // Update ON data
+    // Update OFF data 
     //
     gLog << "============================================================"
          << endl;
-    gLog << "Macro CT1Analysis.C : update file '" << filenameON
+    gLog << "Macro CT1Analysis.C : update file '" << filenameOFF
          << "'" << endl;
 
-    MTaskList tliston;
-    MParList pliston;
+    MTaskList tlistoff;
+    MParList plistoff;
 
 
     // geometry is needed in  MHHillas... classes 
     MGeomCam *fGeom = 
-             (MGeomCam*)pliston->FindCreateObj("MGeomCamCT1", "MGeomCam");
+             (MGeomCam*)plistoff->FindCreateObj("MGeomCamCT1", "MGeomCam");
 
     //-------------------------------------------
@@ -2705,5 +2434,5 @@
     //
 
-    MReadMarsFile read("Events", filenameON);
+    MReadMarsFile read("Events", filenameOFF);
     read.DisableAutoScheme();
 
@@ -2719,5 +2448,5 @@
     //.......................................................................
 
-      MWriteRootFile write(filenameONup);
+      MWriteRootFile write(filenameOFFup);
 
       write.AddContainer("MRawRunHeader", "RunHeaders");
@@ -2732,5 +2461,5 @@
       write.AddContainer("MNewImagePar",  "Events");
 
-      write.AddContainer("HadNN",         "Events");
+      //write.AddContainer("HadNN",         "Events");
       write.AddContainer("HadSC",         "Events");
       write.AddContainer("HadRF",         "Events");
@@ -2749,6 +2478,6 @@
     // entries in MParList
 
-    pliston.AddToList(&tliston);
-    InitBinnings(&pliston);
+    plistoff.AddToList(&tlistoff);
+    InitBinnings(&plistoff);
 
 
@@ -2756,8 +2485,8 @@
     // entries in MTaskList
     
-    tliston.AddToList(&read);
-    tliston.AddToList(&eest2);
-    tliston.AddToList(&write);
-    tliston.AddToList(&contfinal);
+    tlistoff.AddToList(&read);
+    tlistoff.AddToList(&eest2);
+    tlistoff.AddToList(&write);
+    tlistoff.AddToList(&contfinal);
 
     //*****************************
@@ -2768,5 +2497,5 @@
     MProgressBar bar;
     MEvtLoop evtloop;
-    evtloop.SetParList(&pliston);
+    evtloop.SetParList(&plistoff);
     evtloop.SetProgressBar(&bar);
 
@@ -2776,18 +2505,24 @@
         return;
 
-    tliston.PrintStatistics(0, kTRUE);
-    DeleteBinnings(&pliston);
+    tlistoff.PrintStatistics(0, kTRUE);
+    DeleteBinnings(&plistoff);
 
     //---------------------------------------------------
+
     //---------------------------------------------------
-    // Update MC data
+    // Update ON data
     //
     gLog << "============================================================"
          << endl;
-    gLog << "Macro CT1Analysis.C : update file '" << filenameMC
+    gLog << "Macro CT1Analysis.C : update file '" << filenameON
          << "'" << endl;
 
-    MTaskList tlistmc;
-    MParList plistmc;
+    MTaskList tliston;
+    MParList pliston;
+
+
+    // geometry is needed in  MHHillas... classes 
+    MGeomCam *fGeom = 
+             (MGeomCam*)pliston->FindCreateObj("MGeomCamCT1", "MGeomCam");
 
     //-------------------------------------------
@@ -2795,5 +2530,5 @@
     //
 
-    MReadMarsFile read("Events", filenameMC);
+    MReadMarsFile read("Events", filenameON);
     read.DisableAutoScheme();
 
@@ -2809,5 +2544,5 @@
     //.......................................................................
 
-      MWriteRootFile write(filenameMCup);
+      MWriteRootFile write(filenameONup);
 
       write.AddContainer("MRawRunHeader", "RunHeaders");
@@ -2822,5 +2557,5 @@
       write.AddContainer("MNewImagePar",  "Events");
 
-      write.AddContainer("HadNN",         "Events");
+      //write.AddContainer("HadNN",         "Events");
       write.AddContainer("HadSC",         "Events");
       write.AddContainer("HadRF",         "Events");
@@ -2839,6 +2574,6 @@
     // entries in MParList
 
-    plistmc.AddToList(&tlistmc);
-    InitBinnings(&plistmc);
+    pliston.AddToList(&tliston);
+    InitBinnings(&pliston);
 
 
@@ -2846,8 +2581,8 @@
     // entries in MTaskList
     
-    tlistmc.AddToList(&read);
-    tlistmc.AddToList(&eest2);
-    tlistmc.AddToList(&write);
-    tlistmc.AddToList(&contfinal);
+    tliston.AddToList(&read);
+    tliston.AddToList(&eest2);
+    tliston.AddToList(&write);
+    tliston.AddToList(&contfinal);
 
     //*****************************
@@ -2858,5 +2593,5 @@
     MProgressBar bar;
     MEvtLoop evtloop;
-    evtloop.SetParList(&plistmc);
+    evtloop.SetParList(&pliston);
     evtloop.SetProgressBar(&bar);
 
@@ -2866,199 +2601,19 @@
         return;
 
-    tlistmc.PrintStatistics(0, kTRUE);
-    DeleteBinnings(&plistmc);
-
-
-    //==========   end update   ============================================
-   }
-    
-    enparam.Close();
-
-    gLog << "Macro CT1Analysis : End of Job E_EST_UP" << endl;
-    gLog << "=======================================================" << endl;
- }
-  //---------------------------------------------------------------------
-
-
-  //---------------------------------------------------------------------
-  // Job F_XX
-  //=========
-
-    //  - select g/h separation method XX
-    //  - read MC_XX2.root file 
-    //  - calculate eff. collection area
-    //  - read ON_XX2.root file 
-    //  - apply final cuts
-    //  - calculate flux
-    //  - write root file for ON data after final cuts (ON_XX3.root))
-
-
- if (JobF_XX)
- {
-    gLog << "=====================================================" << endl;
-    gLog << "Macro CT1Analysis : Start of Job F_XX" << endl;
-
-    gLog << "" << endl;
-    gLog << "Macro CT1Analysis : JobF_XX, WXX = " 
-         << JobF_XX  << ",  " << WXX << endl;
-
-    // type of data to be analysed
-    TString typeData = "ON";
-    //TString typeData = "OFF";
-    //TString typeData = "MC";
-    gLog << "typeData = " << typeData << endl;
-
-    TString typeMC   = "MC";
-    TString ext      = "2.root";
-    TString extout   = "3.root";
-
-    //------------------------------
-    // selection of g/h separation method
-    // and definition of final selections
-
-    //TString XX("NN");
-    TString XX("SC");
-    //TString XX("RF");
-    TString fhadronnessName("Had");
-    fhadronnessName += XX;
-    gLog << "fhadronnessName = " << fhadronnessName << endl;
-
-    // maximum values of the hadronness, |ALPHA| and DIST
-    Float_t maxhadronness   = 0.40;
-    Float_t maxalpha        = 20.0;
-    Float_t maxdist         = 10.0;
-    gLog << "Maximum values of hadronness, |ALPHA| and DIST = "
-         << maxhadronness << ",  " << maxalpha << ",  " 
-         << maxdist << endl;
-
-
-    //------------------------------
-    // name of MC file to be used for calculating the eff. collection areas
-    TString filenameArea(outPath);
-    filenameArea += typeMC;
-    filenameArea += "_";
-    filenameArea += XX;
-    filenameArea += ext; 
-    gLog << "filenameArea = " << filenameArea << endl;
-
-    //------------------------------
-    // name of file containing the eff. collection areas
-    TString collareaName(outPath);
-    collareaName += "area_";
-    collareaName += XX;
-    collareaName += ".root";
-    gLog << "collareaName = " << collareaName << endl;
-
-    //------------------------------
-    // name of data file to be analysed
-    TString filenameData(outPath);
-    filenameData += typeData;
-    filenameData += "_";
-    filenameData += XX;
-    filenameData += ext;
-    gLog << "filenameData = " << filenameData << endl;
-
-    //------------------------------
-    // name of output data file (after the final cuts)
-    TString filenameDataout(outPath);
-    filenameDataout += typeData;
-    filenameDataout += "_";
-    filenameDataout += XX;
-    filenameDataout += extout;
-    gLog << "filenameDataout = " << filenameDataout << endl;
-
-
-    //====================================================================
-    gLog << "-----------------------------------------------" << endl;
-    gLog << "Start calculation of effective collection areas" << endl;
-    MParList  parlist;
-    MTaskList tasklist;
-
-    //---------------------------------------
-    // Setup the tasks to be executed
-    //
-    MReadMarsFile reader("Events", filenameArea);
-    reader.DisableAutoScheme();
-
-    MFCT1SelFinal cuthadrons;
-    cuthadrons.SetHadronnessName(fhadronnessName);
-    cuthadrons.SetCuts(maxhadronness, maxalpha, maxdist);
-
-    MContinue conthadrons(&cuthadrons);
-
-    MHMcCT1CollectionArea* collarea = new MHMcCT1CollectionArea();
-
-    MFillH filler("MHMcCT1CollectionArea", "MMcEvt");
-    filler.SetName("CollectionArea");
-
-    //********************************
-    // entries in MParList
-
-    parlist.AddToList(&tasklist);
-    InitBinnings(&parlist);
-    parlist.AddToList(collarea);
-
-    //********************************
-    // entries in MTaskList
-
-    tasklist.AddToList(&reader);   
-    tasklist.AddToList(&conthadrons);
-    tasklist.AddToList(&filler);
-
-    //********************************
-
-    //-----------------------------------------
-    // Execute event loop
-    //
-    MEvtLoop magic;
-    magic.SetParList(&parlist);
-
-    MProgressBar bar;
-    magic.SetProgressBar(&bar);
-    if (!magic.Eventloop())
-        return;
-
-    tasklist.PrintStatistics(0, kTRUE);
-
-    // Calculate effective collection areas 
-    // and display the histograms
-    //
-    //MHMcCT1CollectionArea *collarea = parlist.FindObject("MHMcCT1CollectionArea");
-    collarea->CalcEfficiency();
-    collarea->DrawClone("lego");
-
-    //---------------------------------------------
-    // Write histograms to a file 
-    //
-
-    TFile f(collareaName, "RECREATE");
-    collarea->GetHist()->Write();
-    collarea->GetHAll()->Write();
-    collarea->GetHSel()->Write();
-    f.Close();
-
-
-    gLog << "Calculation of effective collection areas done" << endl;
-    gLog << "-----------------------------------------------" << endl;    
-    //------------------------------------------------------------------
-
-
-    //*************************************************************************
-    //
-    // Analyse the data
-    //
-
-    MTaskList tliston;
-    MParList pliston;
-
-    // geometry is needed in  MHHillas... classes 
-    MGeomCam *fGeom = 
-             (MGeomCam*)pliston->FindCreateObj("MGeomCamCT1", "MGeomCam");
-
-
-    TString fHilName    = "MHillas"; 
-    TString fHilNameExt = "MHillasExt"; 
-    TString fHilNameSrc = "MHillasSrc"; 
-    TString fImgParName = "MNewImagePar"; 
+    tliston.PrintStatistics(0, kTRUE);
+    DeleteBinnings(&pliston);
+
+    //---------------------------------------------------
+
+    //---------------------------------------------------
+    // Update MC data
+    //
+    gLog << "============================================================"
+         << endl;
+    gLog << "Macro CT1Analysis.C : update file '" << filenameMC
+         << "'" << endl;
+
+    MTaskList tlistmc;
+    MParList plistmc;
 
     //-------------------------------------------
@@ -3066,11 +2621,19 @@
     //
 
-    MReadMarsFile read("Events", filenameData);
+    MReadMarsFile read("Events", filenameMC);
     read.DisableAutoScheme();
 
+    //---------------------------
+    // calculate estimated energy
+
+    MEnergyEstParam eest2(fHilName);
+    eest2.Add(fHilNameSrc);
+
+    eest2.SetCoeffA(parA);
+    eest2.SetCoeffB(parB);
+
     //.......................................................................
 
-
-      MWriteRootFile write(filenameDataout);
+      MWriteRootFile write(filenameMCup);
 
       write.AddContainer("MRawRunHeader", "RunHeaders");
@@ -3085,5 +2648,272 @@
       write.AddContainer("MNewImagePar",  "Events");
 
-      write.AddContainer("HadNN",         "Events");
+      //write.AddContainer("HadNN",         "Events");
+      write.AddContainer("HadSC",         "Events");
+      write.AddContainer("HadRF",         "Events");
+
+      write.AddContainer("MEnergyEst",    "Events");
+
+    //-----------------------------------------------------------------
+
+    MFCT1SelFinal selfinal(fHilNameSrc);
+    selfinal.SetCuts(maxhadronness, maxalpha, maxdist);
+    selfinal.SetHadronnessName(fhadronnessName);
+    MContinue contfinal(&selfinal);
+
+
+    //*****************************
+    // entries in MParList
+
+    plistmc.AddToList(&tlistmc);
+    InitBinnings(&plistmc);
+
+
+    //*****************************
+    // entries in MTaskList
+    
+    tlistmc.AddToList(&read);
+    tlistmc.AddToList(&eest2);
+    tlistmc.AddToList(&write);
+    tlistmc.AddToList(&contfinal);
+
+    //*****************************
+
+    //-------------------------------------------
+    // Execute event loop
+    //
+    MProgressBar bar;
+    MEvtLoop evtloop;
+    evtloop.SetParList(&plistmc);
+    evtloop.SetProgressBar(&bar);
+
+    Int_t maxevents = -1;
+    //Int_t maxevents = 1000;
+    if ( !evtloop.Eventloop(maxevents) )
+        return;
+
+    tlistmc.PrintStatistics(0, kTRUE);
+    DeleteBinnings(&plistmc);
+
+
+    //==========   end update   ============================================
+   }
+    
+    enparam.Close();
+
+    gLog << "Macro CT1Analysis : End of Job E_EST_UP" << endl;
+    gLog << "=======================================================" << endl;
+ }
+  //---------------------------------------------------------------------
+
+
+  //---------------------------------------------------------------------
+  // Job F_XX
+  //=========
+
+    //  - select g/h separation method XX
+    //  - read MC_XX2.root file 
+    //  - calculate eff. collection area
+    //  - read ON_XX2.root file 
+    //  - apply final cuts
+    //  - calculate flux
+    //  - write root file for ON data after final cuts (ON_XX3.root))
+
+
+ if (JobF_XX)
+ {
+    gLog << "=====================================================" << endl;
+    gLog << "Macro CT1Analysis : Start of Job F_XX" << endl;
+
+    gLog << "" << endl;
+    gLog << "Macro CT1Analysis : JobF_XX, WXX = " 
+         << JobF_XX  << ",  " << WXX << endl;
+
+    // type of data to be analysed
+    //TString typeData = "ON";
+    //TString typeData = "OFF";
+    TString typeData = "MC";
+    gLog << "typeData = " << typeData << endl;
+
+    TString typeMC   = "MC";
+    TString ext      = "3.root";
+    TString extout   = "4.root";
+
+    //------------------------------
+    // selection of g/h separation method
+    // and definition of final selections
+
+    //TString XX("NN");
+    //TString XX("SC");
+    TString XX("RF");
+    TString fhadronnessName("Had");
+    fhadronnessName += XX;
+    gLog << "fhadronnessName = " << fhadronnessName << endl;
+
+    // maximum values of the hadronness, |ALPHA| and DIST
+    Float_t maxhadronness   = 0.40;
+    Float_t maxalpha        = 20.0;
+    Float_t maxdist         = 10.0;
+    gLog << "Maximum values of hadronness, |ALPHA| and DIST = "
+         << maxhadronness << ",  " << maxalpha << ",  " 
+         << maxdist << endl;
+
+
+    //------------------------------
+    // name of MC file to be used for calculating the eff. collection areas
+    TString filenameArea(outPath);
+    filenameArea += typeMC;
+    filenameArea += "_";
+    filenameArea += XX;
+    filenameArea += ext; 
+    gLog << "filenameArea = " << filenameArea << endl;
+
+    //------------------------------
+    // name of file containing the eff. collection areas
+    TString collareaName(outPath);
+    collareaName += "area_";
+    collareaName += XX;
+    collareaName += ".root";
+    gLog << "collareaName = " << collareaName << endl;
+
+    //------------------------------
+    // name of data file to be analysed
+    TString filenameData(outPath);
+    filenameData += typeData;
+    filenameData += "_";
+    filenameData += XX;
+    filenameData += ext;
+    gLog << "filenameData = " << filenameData << endl;
+
+    //------------------------------
+    // name of output data file (after the final cuts)
+    TString filenameDataout(outPath);
+    filenameDataout += typeData;
+    filenameDataout += "_";
+    filenameDataout += XX;
+    filenameDataout += extout;
+    gLog << "filenameDataout = " << filenameDataout << endl;
+
+
+    //====================================================================
+    gLog << "-----------------------------------------------" << endl;
+    gLog << "Start calculation of effective collection areas" << endl;
+    MParList  parlist;
+    MTaskList tasklist;
+
+    //---------------------------------------
+    // Setup the tasks to be executed
+    //
+    MReadMarsFile reader("Events", filenameArea);
+    reader.DisableAutoScheme();
+
+    MFCT1SelFinal cuthadrons;
+    cuthadrons.SetHadronnessName(fhadronnessName);
+    cuthadrons.SetCuts(maxhadronness, maxalpha, maxdist);
+
+    MContinue conthadrons(&cuthadrons);
+
+    //MHMcCT1CollectionArea* collarea = new MHMcCT1CollectionArea();
+    //MHMcCT1CollectionArea* collarea;
+
+    MFillH filler("MHMcCT1CollectionArea", "MMcEvt");
+    filler.SetName("CollectionArea");
+
+    //********************************
+    // entries in MParList
+
+    parlist.AddToList(&tasklist);
+    InitBinnings(&parlist);
+    //parlist.AddToList(collarea);
+
+    //********************************
+    // entries in MTaskList
+
+    tasklist.AddToList(&reader);   
+    tasklist.AddToList(&conthadrons);
+    tasklist.AddToList(&filler);
+
+    //********************************
+
+    //-----------------------------------------
+    // Execute event loop
+    //
+    MEvtLoop magic;
+    magic.SetParList(&parlist);
+
+    MProgressBar bar;
+    magic.SetProgressBar(&bar);
+    if (!magic.Eventloop())
+        return;
+
+    tasklist.PrintStatistics(0, kTRUE);
+
+    // Calculate effective collection areas 
+    // and display the histograms
+    //
+
+    MHMcCT1CollectionArea *collarea = 
+        (MHMcCT1CollectionArea*)parlist.FindObject("MHMcCT1CollectionArea");
+
+    collarea->CalcEfficiency();
+    collarea->DrawClone("lego");
+
+    //---------------------------------------------
+    // Write histograms to a file 
+    //
+
+    TFile f(collareaName, "RECREATE");
+    collarea->GetHist()->Write();
+    collarea->GetHAll()->Write();
+    collarea->GetHSel()->Write();
+    f.Close();
+
+    //delete collarea;
+
+    gLog << "Calculation of effective collection areas done" << endl;
+    gLog << "-----------------------------------------------" << endl;    
+    //------------------------------------------------------------------
+
+
+    //*************************************************************************
+    //
+    // Analyse the data
+    //
+
+    MTaskList tliston;
+    MParList pliston;
+
+    // geometry is needed in  MHHillas... classes 
+    MGeomCam *fGeom = 
+             (MGeomCam*)pliston->FindCreateObj("MGeomCamCT1", "MGeomCam");
+
+    TString fHilName    = "MHillas"; 
+    TString fHilNameExt = "MHillasExt"; 
+    TString fHilNameSrc = "MHillasSrc"; 
+    TString fImgParName = "MNewImagePar"; 
+
+    //-------------------------------------------
+    // create the tasks which should be executed 
+    //
+
+    MReadMarsFile read("Events", filenameData);
+    read.DisableAutoScheme();
+
+    //.......................................................................
+
+
+      MWriteRootFile write(filenameDataout);
+
+      write.AddContainer("MRawRunHeader", "RunHeaders");
+      write.AddContainer("MTime",         "Events");
+      write.AddContainer("MMcEvt",        "Events");
+      write.AddContainer("ThetaOrig",     "Events");
+      write.AddContainer("MSrcPosCam",    "Events");
+      write.AddContainer("MSigmabar",     "Events");
+      write.AddContainer("MHillas",       "Events");
+      write.AddContainer("MHillasExt",    "Events");
+      write.AddContainer("MHillasSrc",    "Events");
+      write.AddContainer("MNewImagePar",  "Events");
+
+      //write.AddContainer("HadNN",         "Events");
       write.AddContainer("HadSC",         "Events");
       write.AddContainer("HadRF",         "Events");
@@ -3104,6 +2934,6 @@
     contfinalgh.SetName("ContSelFinalgh");
 
-    MFillH fillhadnn("hadNN[MHHadronness]", "HadNN");
-    fillhadnn.SetName("HhadNN");
+    //MFillH fillhadnn("hadNN[MHHadronness]", "HadNN");
+    //fillhadnn.SetName("HhadNN");
     MFillH fillhadsc("hadSC[MHHadronness]", "HadSC");
     fillhadsc.SetName("HhadSC");
@@ -3152,5 +2982,5 @@
       tliston.AddToList(&write);
 
-    tliston.AddToList(&fillhadnn);
+    //tliston.AddToList(&fillhadnn);
     tliston.AddToList(&fillhadsc);
     tliston.AddToList(&fillhadrf);
@@ -3186,5 +3016,5 @@
     //
 
-    pliston.FindObject("hadNN", "MHHadronness")->DrawClone();
+    //pliston.FindObject("hadNN", "MHHadronness")->DrawClone();
     pliston.FindObject("hadRF", "MHHadronness")->DrawClone();
     pliston.FindObject("hadSC", "MHHadronness")->DrawClone();