source: trunk/MagicSoft/Mars/mtemp/mmpi/macros/SuperCutsONOFFMacro.C@ 6985

// Silly macro to run the classes that optimize supercuts 
// using ON and OFF data.

// The user only needs to fill/change the variables that control 
// the optimization procedure.


gROOT -> Reset();

void SuperCutsONOFFMacroNew()
{
  gLog.SetNoColors();

    // File containing the data (ON/OFF DATA and path for files (root/ps))

    // From magicserv01
//    TString  ONDataFilename("/.magic/data16a/mazin/data/mcdata/2004_06_30/*HillasON.root");
                                       
//    TString  OFFDataFilename("/.magic/data16a/mazin/data/mcdata/2004_06_30/*HillasOFF.root");

//   TString  ONDataFilename("/.magic/magicserv01/scratch/Period21/HillasFiles/CrabNebula_20040914.root");
//   TString  ONDataFilename("/.magic/magicserv01/scratch/Period21/HillasFiles/CrabNebula_20040921.root");
//   TString  ONDataFilename("/.magic/magicserv01/scratch/Period21/HillasFiles/CrabNebula_20040922.root");

//   TString  ONDataFilename("/.magic/magicserv01/scratch/Period21/HillasFiles/CrabNebula_2004*.root");
  
   //TString  ONDataFilename("/.magic/data22a/mazin/HillasParam/Period24/2004_12_18/Crab_40_25_1_AbsHillas.root");
   TString  ONDataFilename("/.magic/data04a/nadia/Analysis/starfiles/DF/40_25_1/CrabNebula_*_Hillas_Abs.root");
   TString  OFFDataFilename("/.magic/data04a/nadia/Analysis/starfiles/DF/40_25_1/OffCrab1*Hillas_Abs.root");
   //TString ONDataFilename("/.magic/data03a/mazin/results/2004_09_21/CrabNebulaNadiaHillas.root");  
   //TString OFFDataFilename("/.magic/data03a/mazin/results/2004_09_21/CrabNebulaNadiaHillas.root");  


    TString PathForFiles ("/.magic/data03a/mazin/results/Crab/DF/CrabNadia/SuperCuts/Size2000/");  

//    TString PathForFiles ("/.magic/magicserv01/scratch/hbartko/SC/"); // was: Daniel/SuperCuts/Mrk421/2004_04_22/4slices_3520_nc/E800_1200_Opt_MC_Test/

    

    // Boolean variables defining the job of the 
    // macro
 
    // Boolean variable that decides wether data is read from files specified above
    // (ON/OFF) or read from already existing Matrices (which are obviously stored
    // in a root file). The names of the files storing those matrices are produced 
    // automatically using information provided by some of the next variables whose 
    // values must be specified by user.

    Bool_t ReadMatrixFromRootFiles = kTRUE;  
  

    // Boolean variable that controls the supercuts 
    // optimization using the training sample
    // The optimized cuts will be  written in root file
    // located in directory specified before. Name of 
    // the root files is created automatically.

    Bool_t TrainParams   = kTRUE;  

    // Variable that allows the user to skip the optimization on the 
    // train sample. If optimization is skipped (value kTRUE), the
    // previously optimized supercuts (stored in root file) are used 
    // on the train sample.
    
    Bool_t SkipOptimization = kFALSE;

    // Boolean variable that allows the user to write the initial parameters 
    // into the root file that will be used to store the optimum cuts.
    // If fUseInitialSCParams = kTRUE , parameters are written.
    // In this way, the initial SC parameters can be applied on the data (train/test) 
    
    // The initial parameters are ONLY written to the root file if 
    // there is NO SC params optimization, i.e., if variable 
    // fSkipOptimization = kTRUE;
  
    // NOTE: be aware that, if there was a file with the optimized cuts 
    // (previously computed), it will be overwritten with the initial 
    // SC parameters. This is something that I WILL HAVE TO CHANGE IN 
    // future. Yet for the time being...

    Bool_t UseInitialSCParams = kTRUE;  // was kFALSE 
    


    // Variable that decides whether the optimized cuts are used 
    // in the test sample.

    Bool_t TestParams = kTRUE;  


    // Boolean variable that controls wether to combine, OR not, the 
    // alpha distributions computed (after cuts) for the several theta bins
    // in which the TRAIN sample was divided.

    Bool_t CombineCosThetaBinsForTrainSample = kFALSE; 
    
    // Boolean variable that controls wether to combine, OR not, the 
    // alpha distribution computed (after cuts) for the several theta bins
    // in which the TEST sample was divided.

    Bool_t CombineCosThetaBinsForTestSample = kFALSE; 
    

    // Fraction of ON events used for the training/testing
    Double_t whichfractiontrain = 0.5;
    Double_t whichfractiontest = 0.5;

    // Fraction of OFF events used for the training/testing 
    Double_t whichfractiontrainOFF = 0.5;
    Double_t whichfractiontestOFF = 0.5;


    // Efficiency for gammas when using this set of dynamical cuts
    // (i.e., fraction of initial gammas that remain after cuts)

    // Current value is the first estimation of the efficiency of cuts
    // on Mkn421 at a SIZE > 2000 photons

    Double_t gammaeff = 0.6;  
                         

    // Alpha value (degrees) below which signal is expected
    
    Double_t alphasig = 9; 

    // Definition of alpha bkg region (where no signal is expected)

    Double_t alphabkgmin = 30;
    Double_t alphabkgmax = 90;
    
    // Definition of the Size range

    Double_t SizeLow = 2000;
//    Double_t SizeUp =  2000;
    Double_t SizeUp = 1000000;
    
    Double_t LeakageMax = 0.04;
    Double_t DistMax    = 1.5;
    Double_t DistMin    = 0.1;

     // Definition of binning of alpha plots
    Int_t NAlphaBins = 35;
    Double_t AlphaBinLow = -9;
    Double_t AlphaBinUp = 96;

 
    // Boolean variable used to determine wether the normalization factor is 
    // computed from method 1) or 2)
    // 1) Using total number of ON and OFF events before cuts, and tuning the factor 
    //    correcting for "contamination" of gamma events in ON sample 
    // 2) Using number of ON and OFF events after cuts in the background 
    //    region determined by variables fAlphaBkgMin-fAlphaBkgMax
    
    Bool_t NormFactorFromAlphaBkg = kTRUE; // if kTRUE, method 2) is used


    // Boolean variable used to disable the usage ("serious" usage) of the 
    // quantities computed from fits. This will be useful in those cases 
    // where there is too few events to perform a decent fit to the 
    // alpha histograms.
    
    Bool_t UseFittedQuantities = kTRUE;


    // Boolean variable used to control wether to use theta information 
    // in the computation of teh dynamical cuts that take place within
    // class MCT1SupercutsCalc
    //Bool_t NotUseTheta = kFALSE; // kTRUE removes theta from the parameterization of cuts
    Bool_t NotUseTheta = kTRUE; // kTRUE removes theta from the parameterization of cuts

    // Boolean variable used to decide wether to use dynamical cuts or static cuts
    // kTRUE means that static cuts are used.
    Bool_t UseStaticCuts = kFALSE;





    // Name of the Postscript document where all plots
    // will be saved.
    // STORAGE OF PSFILE IS NOT WORKING PROPERLY
    // For the time being, several ps files are produced 
    // and saved in the directory specified by PathForFiles

    /*
    TString PsFileName = ("PsTest23.ps");
    TString CompletePsFileName = (PathForFiles);
    CompletePsFileName += PsFileName;
    TPostScript* PsFile = new TPostScript(CompletePsFileName, 111);
    */

    // Boolean variable used to decide wether initial parameters are 
    // read from ascii file or not. If kTRUE, parameters are retrieved 
    // from ascii file. 

    Bool_t ReadInitParamsFromAsciiFile = kTRUE;

    // Number of SC parameters. The aim of this variable is to cross check 
    // that the number of parameters read from an ascii file 
    // is teh one the user wants. 

    Int_t NInitSCPar = 104;

    // Name of  the ascii file containing the 2 columns, the first one 
    // for initial parameters and the second one for the steps
    // Name must contain also the path.
 
    const char* InitSCParamAsciiFile = 
      // {"../InitialSCParametersSteps/InitSCParamsAndStepsDanielModified1.txt"};
      // {"../InitialSCParametersSteps/FixedStaticCutsInLengthWidthDist.txt"};
      // {"../InitialSCParametersSteps/FixedStaticCutsInLengthWidthDist11.txt"};
      // {"../InitialSCParametersSteps/InitSCParamsAndStepsDanielModified1.txt"};
      // {"../InitialSCParametersSteps/StartingValuesForOptimizationMkn421.txt"};
      // {"../InitialSCParametersSteps/StartingValuesForOptimizationMkn421DynStaticCutsFixedPol2SizeCut3000.txt"};
      // {"../InitialSCParametersSteps/StartingValuesForOptimizationMkn421DynWithDistParametersFixed.txt"};
      // {"../InitialSCParametersSteps/StartingValuesForOptimizationMkn421DynStaticCutsVariablePol2.txt"};
      // {"../InitialSCParametersSteps/StartingValuesForOptimizationMkn421DynStaticCutsVariablePol2WidthCutLowFixed.txt"};
      // {"../InitialSCParametersSteps/StartingValuesForOptimizationMkn421DynCutsOnSize.txt"};
      // {"../InitialSCParametersSteps/StartingValuesForOptimizationMkn421DynCutsOnSizeAndDist.txt"};
      //{"mtemp/mmpi/asciifiles/OptimizedMkn421DynCutsGridWithSelected22pointsMay19.txt"};
      //{"/home/pcmagic16/mazin/mars/Mars260804/mtemp/mmpi/asciifiles/StartingValuesForSmallSizes1.txt"};
//      {"/home/pcmagic16/mazin/mars/Mars260804/mtemp/mmpi/asciifiles/SmallSizeRFStartValues.txt"};
      {"/home/pcmagic16/mazin/mars/Mars260804/mtemp/mmpi/asciifiles/OptimizedCrabFrom2000_a.txt"};
      //{"mtemp/mmpi/asciifiles/VeryGoodDummyValuesWithConcCut.txt"};


    // Name of the root file where alpha distributions, TTree objects
    // with info about the events and cuts applied and  info support histograms 
    // will be stored. 
    // Write only the name of the file. The Path 
    // is the one defined previously

    TString RootFilename = ("RootFileDynCuts.root");

    

    
    

    // Vector containing the theta bins in which data (ON/OFF train/test)
    // will be divided. Actually this vector contains the cosinus of 
    // these theta bins. The dimension of the vector is N+1, where 
    // N is the number of theta bins intended. The first component of the 
    // vector is the low bin edge of the first bin, and the last 
    // vector component the upper bin edge of the last bin.


    TArrayD CosThetaRangeVector(2);
    CosThetaRangeVector[1] = 0.866;  //30
    //CosThetaRangeVector[1] = 0.74;  //42
    //CosThetaRangeVector[1] = 0.82;  //35
    //CosThetaRangeVector[1] = 0;
    //CosThetaRangeVector[1] = 0.825;
    //CosThetaRangeVector[2] = 0.921; 
    //CosThetaRangeVector[3] = 0.961; 
    CosThetaRangeVector[0] = 0.5;  // 60
    //CosThetaRangeVector[0] = 1;


    /*
    TArrayD CosThetaRangeVector(2);
    CosThetaRangeVector[0] = 0.622;
    // CosThetaRangeVector[1] = 0.825;
    //CosThetaRangeVector[2] = 0.921; 
    //CosThetaRangeVector[0] = 0.961; 
    CosThetaRangeVector[1] = 0.984;
    
    */
    // Object of MCT1FindSupercutsONOFFThetaLoop created, data that was specified 
    // above is introduced and ... and the party starts.
     
    MFindSupercutsONOFFThetaLoop FindSupercuts("MFindSupercutsONOFFThetaLoop", 
                                               "Optimizer for the supercuts");

    
    FindSupercuts.SetPathForFiles(PathForFiles);

    FindSupercuts.SetDataONOFFRootFilenames(ONDataFilename, OFFDataFilename);

    FindSupercuts.SetFractionTrainTestOnOffEvents(whichfractiontrain, 
                                                  whichfractiontest, 
                                                  whichfractiontrainOFF, 
                                                  whichfractiontestOFF);


    FindSupercuts.SetGammaEfficiency(gammaeff);


    FindSupercuts.SetAlphaSig(alphasig);

    // Bkg alpha region is set 
    FindSupercuts.SetAlphaBkgMin(alphabkgmin);
    FindSupercuts.SetAlphaBkgMax(alphabkgmax);

    // alpha bkg and signal region set in object FindSupercuts
    // are re-checked in order to be sure that make sense
        
    FindSupercuts.CheckAlphaSigBkg();


    // binning for alpha plots is defined

    FindSupercuts.SetAlphaPlotBinining(NAlphaBins, AlphaBinLow, 
                                       AlphaBinUp);




    // Size range is defined

    FindSupercuts.SetSizeRange(SizeLow, SizeUp);
    FindSupercuts.SetFilters(LeakageMax, DistMax, DistMin);
//    FindSupercuts.SetFilters(0.001, 1.5, 0.2);



    FindSupercuts.SetNormFactorFromAlphaBkg(NormFactorFromAlphaBkg);

    FindSupercuts.SetUseFittedQuantities(UseFittedQuantities);

    FindSupercuts.SetVariableUseStaticCuts(UseStaticCuts);

    FindSupercuts.SetVariableNotUseTheta(NotUseTheta);

    FindSupercuts.SetReadMatricesFromFile(ReadMatrixFromRootFiles);

    FindSupercuts.SetTrainParameters(TrainParams);
    FindSupercuts.SetSkipOptimization(SkipOptimization);
    FindSupercuts.SetUseInitialSCParams(UseInitialSCParams);

    FindSupercuts.SetTestParameters(TestParams);

   
    
    FindSupercuts.SetHadronnessName("MHadSC");
    FindSupercuts.SetHadronnessNameOFF("MHadOFFSC");
    
    FindSupercuts.SetAlphaDistributionsRootFilename (RootFilename);
    
    // FindSupercuts.SetPostScriptFile (PsFile);

    FindSupercuts.SetCosThetaRangeVector (CosThetaRangeVector);
    

    // Names for all root files (matrices, alpha distributions...)
    // are created 
    FindSupercuts.SetALLNames();

    if(ReadInitParamsFromAsciiFile)
    {
        // Initial SC Parameters and steps are retrieved from 
        // Ascii file
        if(!FindSupercuts.ReadSCParamsFromAsciiFile(InitSCParamAsciiFile,
                                                    NInitSCPar))
        {
            cout << "Initial SC Parameters could not be read from Ascii file "
                 << InitSCParamAsciiFile << endl
                 << "Aborting execution of macro... " << endl;
            return;

        }
    }
 
   


    // Finally loop over all theta bins defined is executed
    
    if (!FindSupercuts.LoopOverThetaRanges())
      {
        cout << "Function MFindSupercutsONOFFThetaLoop::LoopOverThetaRanges()" << endl
             << "could not be performed" << endl;

      }

    

    // Nex and Significance are computed vs alphasig
    
    if (!FindSupercuts.ComputeNexSignificanceVSAlphaSig())
      {
        cout << "Function MFindSupercutsONOFFThetaLoop::ComputeNexSignificanceVSAlphaSig()" << endl
             << "could not be performed" << endl;
        

      }
    



    // Option to store ps files in a single ps document is still not working
    /*
    PsFile -> Close();
    PsFile = NULL;
    */

    // Several theta bins are combined to produced a single alpha plot (for train and test)
    // with single Nex and significances

    if (CombineCosThetaBinsForTrainSample || CombineCosThetaBinsForTestSample)
    {
        if(!FindSupercuts.ComputeOverallSignificance(CombineCosThetaBinsForTrainSample, 
                                                     CombineCosThetaBinsForTestSample))
          {
            cout << "Function MFindSupercutsONOFFThetaLoop::ComputeOverallSignificance" << endl
                 << "could not be performed" << endl;
          }
       

    }


   
}