getTimePeakTotal.C

Go to the documentation of this file.

#include "TCanvas.h"
#include "TF2.h"
#include "TFile.h"
#include "TGraph.h"
#include "TGraphErrors.h"
#include "TGraphAsymmErrors.h"
#include "TH1D.h"
#include "TH2D.h"
#include "THStack.h"
#include "TLegend.h"
#include "TLine.h"
#include "TMath.h"
#include "TPad.h"
#include "TPaveLabel.h"
#include "TProfile.h"
#include "TStyle.h"
#include "TTimeStamp.h"

#include "OscLib/OscCalcPMNSOpt.h"
#include "StandardRecord/Proxy/SRProxy.h"

#include "CAFAna/Analysis/Calcs.h"
#include "CAFAna/Analysis/Plots.h"
#include "CAFAna/Core/Binning.h"
#include "CAFAna/Core/EventList.h"
#include "CAFAna/Core/LoadFromFile.h"
#include "CAFAna/Core/ReweightableSpectrum.h"
#include "CAFAna/Core/Spectrum.h"
#include "CAFAna/Core/SpectrumLoader.h"
#include "CAFAna/Cuts/Cuts.h"
#include "3FlavorAna/Cuts/NumuCuts.h"
#include "3FlavorAna/Cuts/NumuCuts2018.h"
#include "3FlavorAna/Vars/NumuVars.h"




#include "CAFAna/Cuts/TimingCuts.h"
#include "CAFAna/Prediction/PredictionNoExtrap.h"
#include "CAFAna/Vars/GenieWeights.h"
#include "CAFAna/Vars/PPFXWeights.h"
#include "CAFAna/Vars/Vars.h"
#include "CAFAna/Vars/TruthVars.h"
#include "CAFAna/Vars/XsecTunes.h"
#include "CAFAna/Vars/GenieWeights.h"
#include "CAFAna/Vars/PPFXWeights.h"

#include "DataQuality/TimePeak/Plots/VarsAndCuts.h"
#include "DataQuality/TimePeak/Plots/Includes.h"

using namespace ana;


struct Plot
{
  std::string label;
  std::string fname;
  Binning bins;
  Var var;
};

void getTimePeakTotal(std::string samdef = "" , std::string s_Period = "" )
{
  
  SpectrumLoader loader(samdef);
  int timestamp = 1610474064;
  int nDays = 41;  
  const int kNumPlots = 11;

  
  //ADDED THIS SO THAT CAN COUNT NUMBER OF SPILLS.
  const SpillVar kCountMySpills([](const caf::SRSpillProxy* sp){
      return sp->isgoodspill;
    });
  TH1D *hSpill = new TH1D("hSpill", "hSpill", 2, -0.5, 1.5);
  loader.AddSpillHistogram(hSpill, kCountMySpills, kNoSpillCut);

 
  Spectrum* hData[kNumPlots];
  Spectrum* hBkg[kNumPlots];
  
  Plot plots[kNumPlots] =
    {
      {"Slice time [#mus]", "tus", Binning::Simple(45,1,541), ana::kSliceTime},
      {"Run", "run", Binning::Simple(15000, 15000, 30000), SIMPLEVAR(hdr.run)},
      {"Time", "time", Binning::Simple(nDays<3 ? 24*nDays : nDays,timestamp - 86400*nDays,timestamp), SIMPLEVAR(hdr.unixtime)},
      {"Max y [m]", "maxy", kXYBins, kBoxMaxY},
      {"Slice time [#mus]", "tighttus", Binning::Simple(360,205,241), ana::kSliceTime},
      {"No. Kalman tracks", "nkal", Binning::Simple(10, 0, 10), SIMPLEVAR(trk.kalman.ntracks)},
      {"No. Hits", "nhit", Binning::Simple(10, 0, 500), ana::kNHit},
      {"Kalman angle", "angkal", Binning::Simple(5, 0, 1), SIMPLEVAR(sel.cosrej.anglekal)},
      {"Z direction", "dirz", Binning::Simple(5, 0.5, 1), ana::kDirZ},
      {"BDT", "bdt", Binning::Simple(20, 0, 1), ana::kNumuContPID},
      {"p_{T}/p", "ptp", Binning::Simple(20, 0, 1), ana::kPtP}
    };
  
  for(int i = 0; i < kNumPlots; ++i){
    Plot p = plots[i];
    hData[i] = new Spectrum(p.label, p.bins, loader, p.var, i <= 1 ? kCut : kCut && kInTimeTight, kNoShift, kUnweighted);
    hBkg[i]  = new Spectrum(p.label, p.bins, loader, p.var, kCut && !kInTimeTight, kNoShift, kUnweighted);
  }
  
  Spectrum potTime("Date","POT #times Mass", loader, Binning::Simple(nDays<3 ? 24*nDays : nDays,timestamp - 86400*nDays,timestamp), SIMPLEVAR(hdr.unixtime), 
                   Binning::Simple(500,0,50e12), kMassFracPOT, kNoCut, kNoShift, kUnweighted);
  Spectrum massDenom("", Binning::Simple(1, 0, 1), loader, kUnweighted, kNoCut, kNoShift, kSpillPOT);
  Spectrum massNum("", Binning::Simple(1, 0, 1), loader, kUnweighted, kNoCut, kNoShift, kMassFracPOT);
  Spectrum *hSig_time = new Spectrum("time", Binning::Simple(nDays<3 ? 24*nDays : nDays,timestamp - 86400*nDays,timestamp), loader,
                                     SIMPLEVAR(hdr.unixtime),kCut && kInTimeTight,kNoShift,kUnweighted);
  
  Spectrum *hBkg_time = new Spectrum("time", Binning::Simple(nDays<3 ? 24*nDays : nDays,timestamp - 86400*nDays,timestamp), loader,
                                     SIMPLEVAR(hdr.unixtime),kCut && !kInTimeTight, kNoShift,kUnweighted);
  
  loader.Go();
  
  
  
  std::string fileName = "dataspectra_"+s_Period+".root";
  TFile *file = TFile::Open(fileName.c_str(), "RECREATE");
  
  for(int i = 0; i < kNumPlots; ++i){
    Plot p = plots[i];
    hData[i]->SaveTo(file, ("data_"+p.fname).c_str() );
    hBkg[i]->SaveTo(file, ("bkg_"+p.fname).c_str() );
  }
  
  potTime.SaveTo(file, "potTime");
  massDenom.SaveTo(file, "massDenom");
  massNum.SaveTo(file, "massNum");
  hSig_time->SaveTo(file, "sig_time");
  hBkg_time->SaveTo(file,"BKG_time");
  hSpill->Write();
  
  
  file->Close();
  
  std::cout << "Spectra made! \n" << std::endl;
  
}