nuCrystalBall.C

Go to the documentation of this file.

#ifdef __CINT__
void nuCrystalBall(const TString filename, const TString opt)
{
  std::cout << "Sorry, you must run in compiled mode" << std::endl;
}
#else

#include "CAFAna/Core/SpectrumLoader.h"
#include "CAFAna/Vars/GenieWeights.h"
#include "CAFAna/Vars/PPFXWeights.h"
#include "CAFAna/Cuts/Cuts.h"
#include "CAFAna/Cuts/SpillCuts.h"
#include "CAFAna/Prediction/PredictionNoExtrap.h"

#include "TString.h"
#include "TFile.h"

#include <iostream>
#include <vector>

#include "Utilities/func/MathUtil.h"

#include "CVN/macros/nuCrystalBallDefs.h"

using namespace ana;

void nuCrystalBall(const TString filename, const TString opt)
{
  // Create loaders
  // Set miniprod file definitions
  std::vector<SpectrumLoader> loaders;
  loaders.push_back(SpectrumLoader("prod_caf_R17-09-05-prod4recopreview.f_fd_genie_nonswap_fhc_nova_v08_full_v1_addShortSimpleCVN"));
  loaders.push_back(SpectrumLoader("prod_caf_R17-09-05-prod4recopreview.f_fd_genie_nonswap_rhc_nova_v08_full_v1_addShortSimpleCVN"));
  std::vector<SpectrumLoader> loadersSwap;
  loadersSwap.push_back(SpectrumLoader("prod_caf_R17-09-05-prod4recopreview.f_fd_genie_fluxswap_fhc_nova_v08_full_v1_addShortSimpleCVN"));
  loadersSwap.push_back(SpectrumLoader("prod_caf_R17-09-05-prod4recopreview.f_fd_genie_fluxswap_rhc_nova_v08_full_v1_addShortSimpleCVN"));

  // Define weight
  const Var weight = kPPFXFluxCVWgt * kXSecCVWgt2017;

  // Get vectors of vars, cuts, cut names, pid cuts and pid cut names
  std::vector<VarDef> vars = GetVars(opt);
  std::vector<CutDef> cuts = GetCuts(opt);
  std::vector<std::vector<PIDCutDef> > pidCuts = GetPIDCuts(opt);

  // Create vector to fill with predictions
  std::vector<PredDef> preds;

  // Loop over FHC and RHC
  for(auto rhc = 0; rhc < 2; ++rhc)
  {
    loaders[rhc].SetSpillCut(kStandardSpillCuts);
    loadersSwap[rhc].SetSpillCut(kStandardSpillCuts);

    TString fhcName = "fhc";
    if(rhc)
      fhcName = "rhc";

    // Make predictions for var distributions to do the fom optimisation
    // PIDS:
    // 1. short Simple
    // 2. ELU
    // 3. 2017
    // Loop over vars
    for(const auto var:vars)
    {
      // Loop over cuts
      for(auto cut:cuts)
      {
        // Add No PID cut prediction
        preds.push_back({
          new PredictionNoExtrap(loaders[rhc], loadersSwap[rhc],
                                 var.name.Data(), var.binning, var.var,
                                 cut.cut, kNoShift, weight),
          cut.name.Data(), "noPIDCut", var.name.Data(), fhcName});
        std::cout << "Making prediction: "
          << var.name.Data() << " "
          << "noPIDCut" << " "
          << cut.name.Data() << " "
          << fhcName << std::endl;

        // Now lets make a PID cut
        // Loop over PID cuts
        for(auto pidCut:pidCuts[rhc])
        {
          // Select if we want to create a prediction for combo
          if(var.type.Contains("pid", TString::kIgnoreCase))
            // only apply cut if relevant for PID
          {
            if(!((pidCut.name.Contains(var.name.Data())) ||
                 (pidCut.name.Contains("allPIDs"))))
              continue;  // Don't apply PID cut
          }
          // Add prediciton
          preds.push_back({
            new PredictionNoExtrap(loaders[rhc], loadersSwap[rhc],
                                   var.name.Data(), var.binning, var.var,
                                   cut.cut && pidCut.cut, kNoShift, weight),
            cut.name.Data(), pidCut.name.Data(), var.name.Data(), fhcName});
          std::cout << "Making prediction: "
            << var.name.Data() << " "
            << pidCut.name.Data() << " "
            << cut.name.Data() << " "
            << fhcName << std::endl;
        }// loop over PID cuts
      }// loop over cuts
    }// loop over vars
  }// loop over FHC/RHC

  // Loop over FHC and RHC
  for(auto rhc = 0; rhc < 1; ++rhc)
  {
    std::cout << "Loaders Go!" << std::endl;
    loaders[rhc].Go();
    loadersSwap[rhc].Go();
  }

  // Write predictions to file
  auto outFile = new TFile(filename, "recreate");
  std::cout << "writing predictions..." << std::endl;
  for(auto pred:preds)
  {
    std::cout << "pred_nxp_"
      << pred.cutName << "_"
      << pred.pidCutName << "_"
      << pred.varName << "_"
      << pred.fhcName << std::endl;
    pred.prediction->SaveTo(outFile->mkdir(
      "pred_nxp_" + pred.cutName + "_" +
      pred.pidCutName + "_" + pred.varName + "_" + pred.fhcName));
  }
std::cout << "Ending script..." << std::endl;
}// end of nuCrystalBall
#endif