NuSLoadMacroDefs.h

Go to the documentation of this file.

#pragma once

#include <stdio.h>
#include <string>

#include "CAFAna/Core/Binning.h"
#include "CAFAna/Core/Cut.h"
#include "CAFAna/Core/Spectrum.h"
#include "CAFAna/Core/Utilities.h"
#include "CAFAna/Cuts/Cuts.h"
#include "NuXAna/Cuts/NusCuts.h"
#include "CAFAna/Cuts/TruthCuts.h"
#include "CAFAna/Vars/Vars.h"
#include "OscLib/IOscCalc.h"


#include "TVector3.h"

#define MAXSPEC 5

class TDirectory;

void Preliminary();
void Simulation();
void CenterTitles(TH1* histo);

namespace ana
{
  // Strings to file paths

  // PNFS concat directory
  const std::string diskdir   = "/pnfs/nova/persistent/users/gsdavies/nus/concat/S16-04-08/";
  const std::string diskdirFD = "/pnfs/nova/persistent/users/gsdavies/nus/concat/S16-05-20/";
  // feature_caf_size (fcs)
  const std::string diskdir_fcs = "/pnfs/nova/persistent/production/concat/R16-11-12-feature_caf_size.a/";
  // I put the bluearc location here for reference but one should just 
  // use the /pnfs directory, it's just as quick (if not quicker)
  // The bluearc location is good to have in case of SAM issues.
  //const std::string diskdir = "/nova/ana/steriles/concat/S16-04-08/";

  // MC decafs
  const std::string fnamenear       = "defname: prod_decaf_S16-04-08_nd_genie_nonswap_genierw_fhc_nova_v08_full_nus_or_numu_contain_v1 with stride 2 with limit 10000";
  const std::string fnamefar        = "prod_decaf_S16-04-08_fd_genie_nonswap_fhc_nova_v08_full_nus_contain_v1";
  const std::string fnameswap       = "prod_decaf_S16-04-08_fd_genie_fluxswap_fhc_nova_v08_full_nus_contain_v1";
  const std::string fnametau        = "prod_decaf_S16-04-08_fd_genie_tau_fhc_nova_v08_full_nus_contain_v1";

  // MC cafs
  const std::string fnamenearcaf    = "defname: prod_caf_R16-03-03-prod2reco.a_nd_genie_nonswap_genierw_fhc_nova_v08_full_v1 with limit 10000";
  const std::string fnamefarcaf     = "prod_caf_R16-03-03-prod2reco.a_fd_genie_nonswap_fhc_nova_v08_full_v1";
  const std::string fnameswapcaf    = "prod_caf_R16-03-03-prod2reco.a_fd_genie_fluxswap_fhc_nova_v08_full_v1";
  const std::string fnametaucaf     = "prod_caf_R16-03-03-prod2reco.a_fd_genie_tau_fhc_nova_v08_full_v1";

  // MC concats
  const std::string fnamenear_concat    = diskdir + "prod_decaf_S16-04-08_nd_genie_nonswap_genierw_fhc_nova_v08_full_nus_or_numu_contain_v1/*.root";

  const std::string fnamenear_concat_fcs = diskdir_fcs + "prod_decaf_R16-11-12-feature_decaf_size.a_nd_genie_nonswap_genierw_fhc_nova_v08_full_v1_nue_or_numu/*.root";
  //const std::string fnamenear_concat_fcs = "bckhouse_newcaf_overlap_prod_caf_R16-11-12-feature_caf_size.a_nd_genie_nonswap_genierw_fhc_nova_v08_full_v1";
  // FD nominal MC by period/epoch incl. veto information
  // First a vector of strings of period1+period2+epoch3b
  const std::vector<std::string> fFDMC_non = {diskdirFD + "by_period/prod_decaf_S16-05-20_fd_genie_nonswap_fhc_nova_v08_period1_nus_contain_v1_prod2-snapshot.root",
                                              diskdirFD + "by_period/prod_decaf_S16-05-20_fd_genie_nonswap_fhc_nova_v08_period2_nus_contain_v1_prod2-snapshot.root",
                                              diskdirFD + "by_period/prod_decaf_S16-05-20_fd_genie_nonswap_fhc_nova_v08_epoch3b_nus_contain_v1_prod2-snapshot.root"};

  const std::vector<std::string> fFDMC_swp = {diskdirFD + "by_period/prod_decaf_S16-05-20_fd_genie_fluxswap_fhc_nova_v08_period1_nus_contain_v1_prod2-snapshot.root",
                                              diskdirFD + "by_period/prod_decaf_S16-05-20_fd_genie_fluxswap_fhc_nova_v08_period2_nus_contain_v1_prod2-snapshot.root",
                                              diskdirFD + "by_period/prod_decaf_S16-05-20_fd_genie_fluxswap_fhc_nova_v08_epoch3b_nus_contain_v1_prod2-snapshot.root"};

  const std::vector<std::string> fFDMC_tau = {diskdirFD + "by_period/prod_decaf_S16-05-20_fd_genie_tau_fhc_nova_v08_period1_nus_contain_v1_prod2-snapshot.root",
                                              diskdirFD + "by_period/prod_decaf_S16-05-20_fd_genie_tau_fhc_nova_v08_period2_nus_contain_v1_prod2-snapshot.root",
                                              diskdirFD + "by_period/prod_decaf_S16-05-20_fd_genie_tau_fhc_nova_v08_epoch3b_nus_contain_v1_prod2-snapshot.root"};

  // Then a string for epoch3c alone
  const std::string fFDMC_non_3c = diskdirFD + "by_period/prod_decaf_S16-05-20_fd_genie_nonswap_fhc_nova_v08_epoch3c_nus_contain_v1_prod2-snapshot.root";
  const std::string fFDMC_swp_3c = diskdirFD + "by_period/prod_decaf_S16-05-20_fd_genie_fluxswap_fhc_nova_v08_epoch3c_nus_contain_v1_prod2-snapshot.root";
  const std::string fFDMC_tau_3c = diskdirFD + "by_period/prod_decaf_S16-05-20_fd_genie_tau_fhc_nova_v08_epoch3c_nus_contain_v1_prod2-snapshot.root";

  const std::string fFDMC_non_fcs = "defname: bckhouse_newcaf_overlap_prod_caf_R16-11-12-feature_caf_size.a_fd_genie_nonswap_fhc_nova_v08_period1_v1 OR defname: bckhouse_newcaf_overlap_prod_caf_R16-11-12-feature_caf_size.a_fd_genie_nonswap_fhc_nova_v08_period2_v1 OR defname: bckhouse_newcaf_overlap_prod_caf_R16-11-12-feature_caf_size.a_fd_genie_nonswap_fhc_nova_v08_epoch3b_v1";
  const std::string fFDMC_swp_fcs = "defname: bckhouse_newcaf_overlap_prod_caf_R16-11-12-feature_caf_size.a_fd_genie_fluxswap_fhc_nova_v08_period1_v1 OR defname: bckhouse_newcaf_overlap_prod_caf_R16-11-12-feature_caf_size.a_fd_genie_fluxswap_fhc_nova_v08_period2_v1 OR defname: bckhouse_newcaf_overlap_prod_caf_R16-11-12-feature_caf_size.a_fd_genie_fluxswap_fhc_nova_v08_epoch3b_v1";
  
  const std::string fFDMC_non_3c_fcs = "prod_caf_R16-11-12-feature_caf_size.a_fd_genie_nonswap_fhc_nova_v08_epoch3c_v1";
  const std::string fFDMC_swp_3c_fcs = "prod_caf_R16-11-12-feature_caf_size.a_fd_genie_fluxswap_fhc_nova_v08_epoch3c_v1";

  
  // These are the old versions. In theory, only use these for the 14db
  const std::string fnamefar_concat     = diskdir + "prod_decaf_S16-04-08_fd_genie_nonswap_fhc_nova_v08_full_nus_contain_v1.root";
  const std::string fnameswap_concat    = diskdir + "prod_decaf_S16-04-08_fd_genie_fluxswap_fhc_nova_v08_full_nus_contain_v1.root";
  const std::string fnametau_concat     = diskdir + "prod_decaf_S16-04-08_fd_genie_tau_fhc_nova_v08_full_nus_contain_v1.root";

  // FD data cosmics
  std::vector<std::string> MakeCosOutOfTimeList()
  {
    std::vector<std::string> fnamecos_concat_1 = Wildcard(
      diskdirFD + "by_period/prod_decaf_S16-05-20_fd_numi_fhc_period1_nus_contain_v1_goodruns_prod2-snapshot/"
      + "prod_decaf_S16-05-20_fd_numi_fhc_period1_nus_contain_v1_goodruns_prod2-snapshot*.root"
    );
    std::vector<std::string> fnamecos_concat_2 = Wildcard(
      diskdirFD + "by_period/prod_decaf_S16-05-20_fd_numi_fhc_period2_nus_contain_v1_goodruns_prod2-snapshot/"
      + "prod_decaf_S16-05-20_fd_numi_fhc_period2_nus_contain_v1_goodruns_prod2-snapshot*.root"
    );
    std::vector<std::string> fnamecos_concat_3bc = Wildcard(
      diskdirFD + "by_period/prod_decaf_S16-05-20_fd_numi_fhc_epoch3*_nus_contain_v1_goodruns_prod2-snapshot.root"
    );
    std::vector<std::string> fnamecos_concat_3d = Wildcard(
      diskdirFD + "by_period/prod_decaf_S16-05-20_fd_numi_fhc_epoch3d_nus_contain_v1_goodruns_prod2-snapshot/"
      + "prod_decaf_S16-05-20_fd_numi_fhc_epoch3d_nus_contain_v1_goodruns_prod2-snapshot*.root"
    );
    std::vector<std::string> ret;
    for(const auto& file : fnamecos_concat_1)   { ret.push_back(file); }
    for(const auto& file : fnamecos_concat_2)   { ret.push_back(file); }
    for(const auto& file : fnamecos_concat_3bc) { ret.push_back(file); }
    for(const auto& file : fnamecos_concat_3d)  { ret.push_back(file); }
    return ret;
  }
  const std::vector<std::string> fnameblind_concat(MakeCosOutOfTimeList());
  const std::string fnamecos_concat = diskdir + "nus_decaf_S16-04-08_fd_data_cosmic/*.root"; // SA cosmic trigger concats, includes nus sandbox and nueveto
  const std::string fnamecos_fcs = "prod_caf_R16-11-12-feature_caf_size.a_fd_cosmic_full_nueveto_v1_goodruns";
  const std::string fnamecosnueveto = "defname:nus_decaf_S16-04-08_fd_data_cosmic with limit 10000";   // SA cosmics decafs, nueveto was ran as a preselection.

  // concats are actual ND NuMI data
  const std::string fnameneardata_concat = diskdir + "prod_decaf_S16-04-08_nd_numi_fhc_full_nus_or_numu_contain_v1_goodruns/*.root";
  // Note: nue or numu only - testing only!
  const std::string fnameneardata_concat_fcs = diskdir_fcs + "prod_decaf_R16-11-12-feature_decaf_size.a_nd_numi_fhc_full_v1_goodruns_nue_or_numu/*.root";
  //const std::string fnameneardata_concat_fcs = "bckhouse_newcaf_overlap_prod_caf_R16-11-12-feature_caf_size.a_nd_numi_fhc_full_v1_goodruns";
  // FD NuMI data
  const std::string fnamefardata  = "";
  
  // MC decaf small sample
  const std::string fnamenearsmall = "defname: prod_decaf_S16-04-08_nd_genie_nonswap_genierw_fhc_nova_v08_full_nus_or_numu_contain_v1 with stride 2 with limit 100";
  const std::string fnamefarsmall  = "defname: prod_decaf_S16-04-08_fd_genie_nonswap_fhc_nova_v08_full_nus_contain_v1 with limit 100";
  const std::string fnameswapsmall = "defname: prod_decaf_S16-04-08_fd_genie_fluxswap_fhc_nova_v08_full_nus_contain_v1 with limit 100";
  const std::string fnametausmall  = "defname: prod_decaf_S16-04-08_fd_genie_tau_fhc_nova_v08_full_nus_contain_v1 with limit 100";

  // Data small samples
  const std::string fnamecossmall  = "defname: nus_decaf_S16-04-08_fd_data_cosmic with limit 100";
  const std::string fnameneardatasmall = "defname: prod_decaf_S16-04-08_nd_genie_nonswap_genierw_fhc_nova_v08_full_nus_contain_v1 with stride 2 with offset 1 with limit 100";

  //......................................................................
  const Binning kAngleCosBins = Binning::Simple(200, -1., 1.); // Binning for -1 to 1 in CosTheta
  const Binning kAngleDegBins = Binning::Simple(360, 0.,  360.); // Binning for angles in degrees

  const Binning kBooleanBins = Binning::Simple(3, -1.5, 1.5);

  const Binning kExtendedEBins = Binning::Simple(300, 0., 30.); // Binning for going to higher energies

  const Binning kFinePIDBins = Binning::Simple(100, 0., 1.); // Binning for PID plots (BUT NOT REMID!!!)

  const Binning kRatioBins = Binning::Simple(100, 0., 1.);

  // Short hand for following macros
  const Cut FlavSel[4] = {kNoCut, kIsNC, (kIsNumuCC || kIsNumuApp), (kIsBeamNue || kIsSig)};

  const Cut kNDTruthFull(
    [](const caf::SRProxy* sr)
    {
      if(sr->hdr.ismc == false) { return true; }
      if(sr->mc.nnu < 1) { return false; }
      if(sr->mc.nu[0].vtx.X() < NDL) { return false; }
      if(sr->mc.nu[0].vtx.X() > NDR) { return false; }
      if(sr->mc.nu[0].vtx.Y() < NDB) { return false; }
      if(sr->mc.nu[0].vtx.Y() > NDT) { return false; }
      if(sr->mc.nu[0].vtx.Z() < NDF) { return false; }

      double NDM = 1587.;
      if(sr->mc.nu[0].vtx.Z() > NDM) { return false; }
      if(sr->mc.nu[0].vtx.Z() > NDE) {
        double MCT = (2./3.)*(NDT - NDB) + NDB;
        if(sr->mc.nu[0].vtx.Y() > MCT) { return false; }
      }
      return true;
    });

  // Macros to quickly generate arrays of Spectrum objects
  // The macros make sure the spectra are in a consistent order for plotting functions
#define SPECARR(label, bins, loaderbeam, loadercos, var, cut) {Spectrum(label, bins, loaderbeam, var, cut && FlavSel[0]), Spectrum(label, bins, loaderbeam, var, cut && FlavSel[1]), Spectrum(label, bins, loaderbeam, var, cut && FlavSel[2]), Spectrum(label, bins, loaderbeam, var, cut && FlavSel[3]), Spectrum(label, bins, loadercos, var, cut)}

#define SPECARR2D(label, loaderbeam, loadercos, binsX, varX, binsY, varY, cut) {Spectrum(label, loaderbeam, binsX, varX, binsY, varY, cut && FlavSel[0]), Spectrum(label, loaderbeam, binsX, varX, binsY, varY, cut && FlavSel[1]), Spectrum(label, loaderbeam, binsX, varX, binsY, varY, cut && FlavSel[2]), Spectrum(label, loaderbeam, binsX, varX, binsY, varY, cut && FlavSel[3]), Spectrum(label, loadercos, binsX, varX, binsY, varY, cut)}
}