resolutionplotter.C

Go to the documentation of this file.

/* Gets ND Data/MC Numu varaibles spectra and converts them into POT normalized histograms.
 */
#ifdef __CINT__
void resolutionplotter()
{
  std::cout<<" Sorry, you must run in compiled mode" <<std::endl;
}
#else
#include "TLine.h"
#include "Rtypes.h"

#include "TCanvas.h"
#include "TFile.h"
#include "TGraph.h"
#include "TH1.h"
#include "TH2.h"
#include "TMath.h"
#include "TGaxis.h"
#include "TMultiGraph.h"
#include "TLegend.h"
#include "TLegendEntry.h"
#include "TLatex.h"
#include "TStyle.h"
#include "THStack.h"
#include "TPaveText.h"
#include "TList.h"
#include "TGaxis.h"
#include "TAttLine.h"
#include "TAttMarker.h"

#include <cmath>
#include <iostream>
#include <vector>
#include <list>
#include <sstream>
#include <string>
#include <sstream>
#include <fstream>
#include <iomanip>

#include "CAFAna/Core/Binning.h"
#include "CAFAna/Cuts/Cuts.h"
#include "CAFAna/Core/Spectrum.h"
#include "CAFAna/Core/SpectrumLoader.h"
#include "CAFAna/Vars/Vars.h"
#include "CAFAna/Core/Utilities.h"
#include "CAFAna/Analysis/Plots.h"
#include "CAFAna/Cuts/SpillCuts.h"
#include "3FlavorAna/Cuts/NumuCuts.h"
#include "CAFAna/Cuts/TruthCuts.h"
#include "3FlavorAna/Vars/NumuVars.h"
#include "CAFAna/Core/SpectrumLoaderBase.h"
#include "CAFAna/Analysis/Style.h"
#include "TColor.h"
#include "Utilities/rootlogon.C"



using namespace ana;
void NDSimulation()
{
  TLatex* prelim = new TLatex(.9, .95, "NOvA Simulation");
  prelim->SetTextColor(kGray+1);
  prelim->SetNDC();
  prelim->SetTextSize(2/30.);
  prelim->SetTextAlign(32);
   TGaxis::SetMaxDigits(3);
  prelim->Draw();
}

void resolutionplotter()
{



/*********************************************************************************************************************/

///Read the Spectrum Root file /////

  TFile *f1= new TFile ("/nova/app/users/acedeno/tag_releaseS17-09-24/NDAna/ncpi0_Systematics/ResolutionstudyOptimizedProng1_2_cuts_kPi0updated.root","READ");
  f1->cd(); 
 

   ///load delta energy 2 prong CVN cuts spectrums (stacked)/////
  auto sndMC=*Spectrum::LoadFrom(f1->GetDirectory("ndMC"));
  auto sndMC3P=*Spectrum::LoadFrom(f1->GetDirectory("ndMC3P"));
  auto sndNCNonPi0=*Spectrum::LoadFrom(f1->GetDirectory("ndNCNonPi0"));
  //auto sndNCMCSig=*Spectrum::LoadFrom(f1->GetDirectory("ndNCMCSig"));
  auto sndNCPi0=*Spectrum::LoadFrom(f1->GetDirectory("ndNCPi0"));
  auto sndCCPi0=*Spectrum::LoadFrom(f1->GetDirectory("ndCCPi0"));
  auto sndNCPi0L=*Spectrum::LoadFrom(f1->GetDirectory("ndNCPi0L"));
  auto sndCCNonPi0=*Spectrum::LoadFrom(f1->GetDirectory("ndCCNonPi0"));
  
  ///load delta energy 2 prong CVN cuts spectrums (unstacked)/////
  auto sndMC1=*Spectrum::LoadFrom(f1->GetDirectory("ndMC"));
  auto sndMC3P1=*Spectrum::LoadFrom(f1->GetDirectory("ndMC3P"));
  auto sndNCNonPi01=*Spectrum::LoadFrom(f1->GetDirectory("ndNCNonPi0"));
  //auto sndNCMCSig=*Spectrum::LoadFrom(f1->GetDirectory("ndNCMCSig"));
  auto sndNCPi01=*Spectrum::LoadFrom(f1->GetDirectory("ndNCPi0"));
  auto sndCCPi01=*Spectrum::LoadFrom(f1->GetDirectory("ndCCPi0"));
  auto sndNCPi0L1=*Spectrum::LoadFrom(f1->GetDirectory("ndNCPi0L"));
  auto sndCCNonPi01=*Spectrum::LoadFrom(f1->GetDirectory("ndCCNonPi0"));
  
   ///load delta energy 1 prong CVN cut spectrums (stacked)/////
  auto sndMC1Cut=*Spectrum::LoadFrom(f1->GetDirectory("ndMC1Cut"));
  auto sndMC3P1Cut=*Spectrum::LoadFrom(f1->GetDirectory("ndMC3P1Cut"));
  auto sndNCNonPi01Cut=*Spectrum::LoadFrom(f1->GetDirectory("ndNCNonPi01Cut"));
  //auto sndNCMCSig=*Spectrum::LoadFrom(f1->GetDirectory("ndNCMCSig"));
  auto sndNCPi01Cut=*Spectrum::LoadFrom(f1->GetDirectory("ndNCPi01Cut"));
  auto sndCCPi01Cut=*Spectrum::LoadFrom(f1->GetDirectory("ndCCPi01Cut"));
  auto sndNCPi0L1Cut=*Spectrum::LoadFrom(f1->GetDirectory("ndNCPi0L1Cut"));
  auto sndCCNonPi01Cut=*Spectrum::LoadFrom(f1->GetDirectory("ndCCNonPi01Cut"));

  ///load reco mass 1 prong CVN cut spectrums (unstacked)/////
  auto sndMC1Cut2=*Spectrum::LoadFrom(f1->GetDirectory("ndMC1Cut"));
  auto sndMC3P1Cut2=*Spectrum::LoadFrom(f1->GetDirectory("ndMC3P1Cut"));
  auto sndNCNonPi01Cut2=*Spectrum::LoadFrom(f1->GetDirectory("ndNCNonPi01Cut"));
  //auto sndNCMCSig=*Spectrum::LoadFrom(f1->GetDirectory("ndNCMCSig"));
  auto sndNCPi01Cut2=*Spectrum::LoadFrom(f1->GetDirectory("ndNCPi01Cut"));
  auto sndCCPi01Cut2=*Spectrum::LoadFrom(f1->GetDirectory("ndCCPi01Cut"));
  auto sndNCPi0L1Cut2=*Spectrum::LoadFrom(f1->GetDirectory("ndNCPi0L1Cut"));
  auto sndCCNonPi01Cut2=*Spectrum::LoadFrom(f1->GetDirectory("ndCCNonPi01Cut"));
   
  // ////2D Plots 2 prong & 3 prong///////
  auto snd2dplot2prong1cut=*Spectrum::LoadFrom(f1->GetDirectory("nd2dplot2prong1cut"));
  auto snd2dplot2prong2cuts=*Spectrum::LoadFrom(f1->GetDirectory("nd2dplot2prong2cuts"));

  /*********************************************************************************************************************/

  ///How to Scale Histogram////
 
  // double myPOT=sndMC.POT();
  const double myPOT = 8e20;
 
/*********************************************************************************************************************/

   //------histograms--------//

  ///load delta energy  spectrums 2 prong CVN cuts (stacked) /////
  auto hndMC = sndMC.ToTH1(myPOT);
  auto hndMC3P = sndMC3P.ToTH1(myPOT);
  auto hndNCNonPi0 = sndNCNonPi0.ToTH1(myPOT);
  auto hndCCPi0 = sndCCPi0.ToTH1(myPOT);
  auto hndNCPi0 = sndNCPi0.ToTH1(myPOT);
  auto hndNCPi0L = sndNCPi0L.ToTH1(myPOT);
  auto hndCCNonPi0 = sndCCNonPi0.ToTH1(myPOT);

  ///load delta energy spectrums 2 prong CVN cuts (unstacked) /////
  auto hndMC1 = sndMC.ToTH1(myPOT);
  auto hndMC3P1 = sndMC3P.ToTH1(myPOT);
  auto hndNCNonPi01 = sndNCNonPi0.ToTH1(myPOT);
  auto hndCCPi01 = sndCCPi0.ToTH1(myPOT);
  auto hndNCPi01 = sndNCPi0.ToTH1(myPOT);
  auto hndNCPi0L1 = sndNCPi0L.ToTH1(myPOT);
  auto hndCCNonPi01 = sndCCNonPi0.ToTH1(myPOT);


  ///load delta energy spectrums 1 prong CVN cut (Stacked) /////
  auto hndMC1Cut = sndMC1Cut.ToTH1(myPOT);
  auto hndMC3P1Cut = sndMC3P1Cut.ToTH1(myPOT);
  auto hndNCNonPi01Cut = sndNCNonPi01Cut.ToTH1(myPOT);
  auto hndCCPi01Cut = sndCCPi01Cut.ToTH1(myPOT);
  auto hndNCPi01Cut = sndNCPi01Cut.ToTH1(myPOT);
  auto hndNCPi0L1Cut = sndNCPi0L1Cut.ToTH1(myPOT);
  auto hndCCNonPi01Cut = sndCCNonPi01Cut.ToTH1(myPOT);

  ///load delta energy spectrums 1 prong CVN cut (Unstacked)/////
  auto hndMC1Cut2 = sndMC1Cut2.ToTH1(myPOT);
  auto hndMC3P1Cut2 = sndMC3P1Cut2.ToTH1(myPOT);
  auto hndNCNonPi01Cut2 = sndNCNonPi01Cut2.ToTH1(myPOT);
  auto hndCCPi01Cut2 = sndCCPi01Cut2.ToTH1(myPOT);
  auto hndNCPi01Cut2 = sndNCPi01Cut2.ToTH1(myPOT);
  auto hndNCPi0L1Cut2 = sndNCPi0L1Cut2.ToTH1(myPOT);
  auto hndCCNonPi01Cut2 = sndCCNonPi01Cut2.ToTH1(myPOT);
 
  ////2D Plots 2 prong & 3 prong///////

  auto hnd2dplot2prong1cut = snd2dplot2prong1cut.ToTH2(myPOT);
  auto hnd2dplot2prong2cuts = snd2dplot2prong2cuts.ToTH2(myPOT);

//****************************Updated Date and string path for saving files **********************************************************///
  TString outDir2 = "/nova/app/users/acedeno/tag_releaseS17-09-24/NDAna/ncpi0_Systematics/Resolution";
  struct FileUp{
    std::string name;
  };
  struct DateUp{
    std::string date;
  };
 
  const int File =6;
  const FileUp CurrentFile[File] = {
    {"/2cuts_resolution_deltaenergy_Stacked"},
    {"/2cuts_resolution_deltaenergy_unstacked"},
    {"/1cut_resolution_deltaenergy_stacked"},
    {"/1cut_resolution_deltaenergy_unstacked"},
    {"/1cut_resolution_2DPlot"},
    {"/2cut_resolution_2DPlot"},
    
  };

   const int Date = 1;
   const DateUp CurrentDate[Date] = {
     {"_03-16-18"},
  };

//------------------------------------------------------------


/*********************************************************************************************************************/
  /// Plot 2 /// Stacked Delta Energy 2 P-CVN
 
   //gStyle->SetTitleY(0.9);
   TCanvas *cs2 = new TCanvas("cs2","cs2",10,10,855,500);
   THStack *hs2 = new THStack("hs2"," ");
 //TCanvas *c2_ndMC =new TCanvas("c2_ndMC","ndMC",800,500);
   //TH1D *c2_ndNCMCSig = new TH1D("c2_ndNCMCSig","ndNCMCSig",100,-4,4);
   // hndNCMCSig->SetLineColor(kBlack);
   // hndNCMCSig->SetFillColor(kGray);
   // hndNCMCSig->SetLineWidth(3);
   // hndNCMCSig->SetMarkerStyle(21);
   // hndNCMCSig->SetMarkerColor(kGray);
   
   // hs2->Add(hndNCMCSig);
  
   //TH1D *c2_ndCCMCSig = new TH1D("c2_ndCCMCSig","ndCCMCSig",100,-4,4);
   hndNCPi0L->SetLineColor(kBlack);
   hndNCPi0L->SetFillColor(kOrange);
   hndNCPi0L->SetLineWidth(3);
   hndNCPi0L->SetMarkerStyle(21);
   hndNCPi0L->SetMarkerColor(kOrange);
   hs2->Add(hndNCPi0L);

   //TH1D *c2_hndCCPi0 = new TH1D("c2_ndCCPi0","ndCCPi0",100,-4,4);
   hndCCPi0->SetLineColor(kBlack);
   hndCCPi0->SetFillColor(kMagenta);
   hndCCPi0->SetLineWidth(3);
   hndCCPi0->SetMarkerStyle(21);
   hndCCPi0->SetMarkerColor(kMagenta);

   hs2->Add(hndCCPi0);
   //TH1D *c2_hndCCNonPi0 = new TH1D("c2_hndCCNonPi0","hndCCNonPi0",100,-4,4);
   hndCCNonPi0->SetLineColor(kBlack);
   hndCCNonPi0->SetFillColor(kAzure);
   hndCCNonPi0->SetLineWidth(3);
   hndCCNonPi0->SetMarkerStyle(21);
   hndCCNonPi0->SetMarkerColor(kAzure);
   hs2->Add(hndCCNonPi0);
   
   //TH1D *c2_hndNCPi0 = new TH1D("c2_ndNCPi0","ndNCPi0",100,-4,4);
   hndNCNonPi0->SetLineColor(kBlack);
   hndNCNonPi0->SetFillColor(kCyan);
   hndNCNonPi0->SetLineWidth(3);
   hndNCNonPi0->SetMarkerStyle(21);
   hndNCNonPi0->SetMarkerColor(kCyan);
   hs2->Add(hndNCNonPi0);
   
   //TH1D *c2_hndNCPi0L = new TH1D("c2_ndNCPi0L","ndNCPi0L",100,-4,4);
   hndNCPi0->SetLineColor(kBlack);
   hndNCPi0->SetFillColor(kGreen);
   hndNCPi0->SetLineWidth(3);
   hndNCPi0->SetMarkerStyle(21);
   hndNCPi0->SetMarkerColor(kGreen);
   hs2->Add(hndNCPi0);

   //TCanvas *cs2 = new TCanvas("cs2","cs2",10,10,855,500);
   TText T; T.SetTextFont(42); T.SetTextAlign(21);
   hs2->Draw("hist"); T.DrawTextNDC(.5,.95,"");
   hs2->GetYaxis()->SetTitle("Events/ 8e20 POT");
   hs2->GetXaxis()->SetTitle("(Reco-True) #Delta E [GeV]");
   hs2->SetMaximum(450);
   hs2->GetXaxis()->SetRangeUser(-1,1);
   gPad->Modified();
   cs2->cd();
   //cs2->SetTopMargin(0.18);
   hndMC3P->SetLineColor(kRed);
   hndMC3P->SetLineWidth(3);
   hndMC3P->Draw("same");
   TLegend * leg2 = new TLegend(0.62,0.5,0.73,0.82);
   //leg2->SetTextAlign(40);
   leg2->SetTextFont(42);
   leg2->SetTextSize(0.05);
   leg2->SetFillStyle(0);
   leg2->AddEntry(hndMC3P,"Total PreSel","kFullDotSmall");
   // leg2->AddEntry(hndNCMCSig,"Signal","kFullDotSmall");
   leg2->AddEntry(hndCCPi0,"CCPi0","kFullDotSmall");
   leg2->AddEntry(hndCCNonPi0,"CCNonPi0","kFullDotSmall");
   leg2->AddEntry(hndNCPi0L,"NC w/Pi0 < 0.3 E","kFullDotSmall");
   leg2->AddEntry(hndNCPi0,"NC w/Pi0 > 0.3 E","kFullDotSmall");
   leg2->AddEntry(hndNCNonPi0,"NCNonPi0","kFullDotSmall");
   
   leg2->Draw("same"); 
    NDSimulation();
    //cs2->SaveAs("/nova/app/users/acedeno/tag_releas9-24/NDAna/ncpi0_Systematics/NcPi0Plots/plot3_2Pr_0_3_CVN07__11-02-17.pdf");
    cs2->SaveAs(outDir2+CurrentFile[0].name.c_str()+CurrentDate[0].date.c_str()+".pdf");

// /*********************************************************************************************************************/
/// Plot 2 /// Unstacked Delta energy 2 P-CVN cuts
 
   //gStyle->SetTitleY(0.9);
   TCanvas *cs13 = new TCanvas("cs13","cs13",10,10,855,500);
   cs13->cd();

   hndMC3P1->SetLineColor(kRed);
   hndMC3P1->SetLineWidth(3);
   // hndMC3P1->SetMarkerStyle(21);
   //hndMC3P1->SetMarkerColor(kRed);
  
   hndCCPi01->SetLineColor(kMagenta);
   hndCCPi01->SetLineWidth(3);
   // hndCCPi01->SetMarkerStyle(21);
   // hndCCPi01->SetMarkerColor(kMagenta);

   hndCCNonPi01->SetLineColor(kAzure);
   hndCCNonPi01->SetLineWidth(3);
   // hndCCNonPi01->SetMarkerStyle(21);
   //hndCCNonPi01->SetMarkerColor(kAzure);

   hndNCPi0L1->SetLineColor(kOrange);
   hndNCPi0L1->SetLineWidth(3);
   // hndNCPi0L1->SetMarkerStyle(21);
   //hndNCPi0L1->SetMarkerColor(kOrange);
  
   hndNCPi01->SetLineColor(kGreen); 
   hndNCPi01->SetLineWidth(3);
   // hndNCPi01->SetMarkerStyle(21);
   //hndNCPi01->SetMarkerColor(kGreen);
  

   hndNCNonPi01->SetLineColor(kCyan);
   hndNCNonPi01->SetLineWidth(3);
   // hndNCNonPi01->SetMarkerStyle(21);
   //hndNCNonPi01->SetMarkerColor(kCyan);

   hndMC3P1->Draw("hist E1");
   hndCCPi01->Draw("hist E1 same");
   hndCCNonPi01->Draw("hist E1 same");
   hndNCPi0L1->Draw("hist E1 same");
   hndNCPi01->Draw("hist E1 same");
   hndNCNonPi01->Draw("hist E1 same");
   hndMC3P1->GetYaxis()->SetTitle("Events/8e20 POT");
   hndMC3P1->GetXaxis()->SetRangeUser(-.5,.5);
   hndMC3P1->SetTitleOffset(0.95,"Y");
   // hndMC3P1->SetTitle("2 Prong Event Preselection + RemID < 0.36 + Prong 1 & 2 CVN Gamma ID > 0.7 ");
   
   TLegend * leg13 = new TLegend(0.62,0.5,0.73,0.82);
   leg13->SetTextFont(42);
   leg13->SetTextSize(0.05);
   leg13->AddEntry(hndMC3P1,"Total PreSel","l");
   leg13->AddEntry(hndCCPi01,"CCPi0","l");
   leg13->AddEntry(hndCCNonPi01,"CCNonPi0","l");
   leg13->AddEntry(hndNCPi0L1,"NC w/Pi0 < 0.3 E","l");
   leg13->AddEntry(hndNCPi01,"NC w/Pi0 > 0.3 E","l");
   leg13->AddEntry(hndNCNonPi01,"NCNonPi0","l");
   
   leg13->Draw("same"); 
    NDSimulation();
    //cs13->SaveAs("/nova/app/users/acedeno/tag_releaseS17-09-24/NDAna/ncpi0_Systematics/NcPi0Plots/plot13_2Pr_0_3_CVN07__11-02-17.pdf");
   cs13->SaveAs(outDir2+CurrentFile[1].name.c_str()+CurrentDate[0].date.c_str()+".pdf");

/*********************************************************************************************************************/

/// Plot 2b /// delta energy stacked 1 P-CVN cut 
 
   //gStyle->SetTitleY(0.9);
   TCanvas *cs2b = new TCanvas("cs2b","cs2b",10,10,855,500);
   THStack *hs2b = new THStack("hs2b"," ");
 //TCanvas *c2_ndMC =new TCanvas("c2_ndMC","ndMC",800,500);
   //TH1D *c2_ndNCMCSig = new TH1D("c2_ndNCMCSig","ndNCMCSig",100,-4,4);
   // hndNCMCSig->SetLineColor(kBlack);
   // hndNCMCSig->SetFillColor(kGray);
   // hndNCMCSig->SetLineWidth(3);
   // hndNCMCSig->SetMarkerStyle(21);
   // hndNCMCSig->SetMarkerColor(kGray);
   
   // hs2->Add(hndNCMCSig);
  
   //TH1D *c2_ndCCMCSig = new TH1D("c2_ndCCMCSig","ndCCMCSig",100,-4,4);
   hndNCPi0L1Cut->SetLineColor(kBlack);
   hndNCPi0L1Cut->SetFillColor(kOrange);
   hndNCPi0L1Cut->SetLineWidth(3);
   hndNCPi0L1Cut->SetMarkerStyle(21);
   hndNCPi0L1Cut->SetMarkerColor(kOrange);
   hs2b->Add(hndNCPi0L1Cut);

   //TH1D *c2_hndCCPi0 = new TH1D("c2_ndCCPi0","ndCCPi0",100,-4,4);
   hndCCPi01Cut->SetLineColor(kBlack);
   hndCCPi01Cut->SetFillColor(kMagenta);
   hndCCPi01Cut->SetLineWidth(3);
   hndCCPi01Cut->SetMarkerStyle(21);
   hndCCPi01Cut->SetMarkerColor(kMagenta);

   hs2b->Add(hndCCPi01Cut);
   //TH1D *c2_hndCCNonPi0 = new TH1D("c2_hndCCNonPi0","hndCCNonPi0",100,-4,4);
   hndCCNonPi01Cut->SetLineColor(kBlack);
   hndCCNonPi01Cut->SetFillColor(kAzure);
   hndCCNonPi01Cut->SetLineWidth(3);
   hndCCNonPi01Cut->SetMarkerStyle(21);
   hndCCNonPi01Cut->SetMarkerColor(kAzure);
   hs2b->Add(hndCCNonPi01Cut);
   
   //TH1D *c2_hndNCPi0 = new TH1D("c2_ndNCPi0","ndNCPi0",100,-4,4);
   hndNCNonPi01Cut->SetLineColor(kBlack);
   hndNCNonPi01Cut->SetFillColor(kCyan);
   hndNCNonPi01Cut->SetLineWidth(3);
   hndNCNonPi01Cut->SetMarkerStyle(21);
   hndNCNonPi01Cut->SetMarkerColor(kCyan);
   hs2b->Add(hndNCNonPi01Cut);
   
   //TH1D *c2_hndNCPi0L = new TH1D("c2_ndNCPi0L","ndNCPi0L",100,-4,4);
   hndNCPi01Cut->SetLineColor(kBlack);
   hndNCPi01Cut->SetFillColor(kGreen);
   hndNCPi01Cut->SetLineWidth(3);
   hndNCPi01Cut->SetMarkerStyle(21);
   hndNCPi01Cut->SetMarkerColor(kGreen);
   hs2b->Add(hndNCPi01Cut);

   //TCanvas *cs2 = new TCanvas("cs2","cs2",10,10,855,500);
   TText T99; T99.SetTextFont(42); T99.SetTextAlign(21);
   hs2b->Draw("hist"); T99.DrawTextNDC(.5,.95,"");
   hs2b->GetYaxis()->SetTitle("Events/ 8e20 POT");
   hs2b->GetXaxis()->SetTitle("(Reco-True) #Delta E [GeV]");
   hs2b->SetMaximum(1200);
   hs2b->GetXaxis()->SetRangeUser(-1,1);
   gPad->Modified();
   cs2b->cd();
   //cs2->SetTopMargin(0.18);
   hndMC3P1Cut->SetLineColor(kRed);
   hndMC3P1Cut->SetLineWidth(3);
   hndMC3P1Cut->Draw("same");
   TLegend * leg2b = new TLegend(0.62,0.5,0.73,0.82);
   //leg2->SetTextAlign(40);
   leg2b->SetTextFont(42);
   leg2b->SetTextSize(0.05);
   leg2b->SetFillStyle(0);
   leg2b->AddEntry(hndMC3P1Cut,"Total PreSel","kFullDotSmall");
   // leg2->AddEntry(hndNCMCSig,"Signal","kFullDotSmall");
   leg2b->AddEntry(hndCCPi01Cut,"CCPi0","kFullDotSmall");
   leg2b->AddEntry(hndCCNonPi01Cut,"CCNonPi0","kFullDotSmall");
   leg2b->AddEntry(hndNCPi0L1Cut,"NC w/Pi0 < 0.3 E","kFullDotSmall");
   leg2b->AddEntry(hndNCPi01Cut,"NC w/Pi0 > 0.3 E","kFullDotSmall");
   leg2b->AddEntry(hndNCNonPi01Cut,"NCNonPi0","kFullDotSmall");
   
   leg2b->Draw("same"); 
    NDSimulation();
    //cs2->SaveAs("/nova/app/users/acedeno/tag_releaseS17-09-24/NDAna/ncpi0_Systematics/NcPi0Plots/plot3_2Pr_0_3_CVN07__11-02-17.pdf");
    cs2b->SaveAs(outDir2+CurrentFile[2].name.c_str()+CurrentDate[0].date.c_str()+".pdf");

/*********************************************************************************************************************/

/// Plot 2 /// delta energy unstacked 1 P-CVN cut
 
   //gStyle->SetTitleY(0.9);
   TCanvas *cs13b = new TCanvas("cs13b","cs13b",10,10,855,500);
   cs13b->cd();

   hndMC3P1Cut2->SetLineColor(kRed);
   hndMC3P1Cut2->SetLineWidth(3);
   // hndMC3P1->SetMarkerStyle(21);
   //hndMC3P1->SetMarkerColor(kRed);
  
   hndCCPi01Cut2->SetLineColor(kMagenta);
   hndCCPi01Cut2->SetLineWidth(3);
   // hndCCPi01->SetMarkerStyle(21);
   // hndCCPi01->SetMarkerColor(kMagenta);

   hndCCNonPi01Cut2->SetLineColor(kAzure);
   hndCCNonPi01Cut2->SetLineWidth(3);
   // hndCCNonPi01->SetMarkerStyle(21);
   //hndCCNonPi01->SetMarkerColor(kAzure);

   hndNCPi0L1Cut2->SetLineColor(kOrange);
   hndNCPi0L1Cut2->SetLineWidth(3);
   // hndNCPi0L1->SetMarkerStyle(21);
   //hndNCPi0L1->SetMarkerColor(kOrange);
  
   hndNCPi01Cut2->SetLineColor(kGreen); 
   hndNCPi01Cut2->SetLineWidth(3);
   // hndNCPi01->SetMarkerStyle(21);
   //hndNCPi01->SetMarkerColor(kGreen);
  

   hndNCNonPi01Cut2->SetLineColor(kCyan);
   hndNCNonPi01Cut2->SetLineWidth(3);
   // hndNCNonPi01->SetMarkerStyle(21);
   //hndNCNonPi01->SetMarkerColor(kCyan);

   hndMC3P1Cut2->Draw("hist E1");
   hndCCPi01Cut2->Draw("hist E1 same");
   hndCCNonPi01Cut2->Draw("hist E1 same");
   hndNCPi0L1Cut2->Draw("hist E1 same");
   hndNCPi01Cut2->Draw("hist E1 same");
   hndNCNonPi01Cut2->Draw("hist E1 same");
   hndMC3P1Cut2->GetYaxis()->SetTitle("Events/8e20 POT");
   hndMC3P1Cut2->GetXaxis()->SetRangeUser(-.5,.5);
   hndMC3P1Cut2->SetTitleOffset(0.95,"Y");
   // hndMC3P1->SetTitle("2 Prong Event Preselection + RemID < 0.36 + Prong 1 & 2 CVN Gamma ID > 0.7 ");
   
   TLegend * leg13b = new TLegend(0.62,0.5,0.73,0.82);
   leg13b->SetTextFont(42);
   leg13b->SetTextSize(0.05);
   leg13b->AddEntry(hndMC3P1Cut2,"Total PreSel","l");
   leg13b->AddEntry(hndCCPi01Cut2,"CCPi0","l");
   leg13b->AddEntry(hndCCNonPi01Cut2,"CCNonPi0","l");
   leg13b->AddEntry(hndNCPi0L1Cut2,"NC w/Pi0 < 0.3 E","l");
   leg13b->AddEntry(hndNCPi01Cut2,"NC w/Pi0 > 0.3 E","l");
   leg13b->AddEntry(hndNCNonPi01Cut2,"NCNonPi0","l");
   
   leg13b->Draw("same"); 
    NDSimulation();
    //cs13->SaveAs("/nova/app/users/acedeno/tag_releaseS17-09-24/NDAna/ncpi0_Systematics/NcPi0Plots/plot13_2Pr_0_3_CVN07__11-02-17.pdf");
   cs13b->SaveAs(outDir2+CurrentFile[3].name.c_str()+CurrentDate[0].date.c_str()+".pdf");

/*********************************************************************************************************************/
   /// Plot 16 ///
   TCanvas *cs10 = new TCanvas("cs10","cs10",10,10,855,500);
   cs10->cd();
   hnd2dplot2prong1cut->Draw("colz");
   hnd2dplot2prong1cut->GetXaxis()->SetTitle("#pi^{0} Reco Energy [GeV]");
   hnd2dplot2prong1cut->GetYaxis()->SetTitle("(Reco-True) #Delta E [GeV]");
   hnd2dplot2prong1cut->SetTitleOffset(0.9,"Y");
   // hnd2dplot2prong->SetTitle("2 Prong Event Preselection + RemID < 0.36 + Prong 1 & 2 CVN Gamma ID > 0.7 Cut");
 NDSimulation();
 // cs10->SaveAs("/nova/app/users/acedeno/tag_releaseS17-09-24/NDAna/ncpi0_Systematics/NcPi0Plots/plot11_2Pr_2dplot_0_3_CVN07__11-02-17.pdf");
   cs10->SaveAs(outDir2+CurrentFile[4].name.c_str()+CurrentDate[0].date.c_str()+".pdf");
  
// /*********************************************************************************************************************/

   /// Plot 16 ///
   TCanvas *cs14 = new TCanvas("cs14","cs14",10,10,855,500);
   cs14->cd();
   hnd2dplot2prong2cuts->Draw("colz");
   hnd2dplot2prong2cuts->GetXaxis()->SetTitle("#pi^{0} Reco Energy [GeV]");
   hnd2dplot2prong2cuts->GetYaxis()->SetTitle("(Reco-True) #Delta E [GeV]");
   hnd2dplot2prong2cuts->SetTitleOffset(0.9,"Y");
   // hnd2dplot2prong->SetTitle("2 Prong Event Preselection + RemID < 0.36 + Prong 1 & 2 CVN Gamma ID > 0.7 Cut");
 NDSimulation();
 // cs10->SaveAs("/nova/app/users/acedeno/tag_releaseS17-09-24/NDAna/ncpi0_Systematics/NcPi0Plots/plot11_2Pr_2dplot_0_3_CVN07__11-02-17.pdf");
   cs14->SaveAs(outDir2+CurrentFile[5].name.c_str()+CurrentDate[0].date.c_str()+".pdf");
  
// /*********************************************************************************************************************/
   //////// Bias vs pion total reco/truth energy and Absolute Resolution vs pion total reco/truth energy/////


   //h=(TH2D*) f->Get("nd2dplot2prong2cuts");
   auto histo=*Spectrum::LoadFrom(f1->GetDirectory("nd2dplot2prong2cuts"));
   auto h = histo.ToTH2(myPOT);
   float datapot = 8.09e20/3.54e21;
   h->Scale(datapot);
   h->Sumw2();
   int nbins = h->GetNbinsX();
   double xmin = h->GetXaxis()->GetXmin();
   
   double xmax = h->GetXaxis()->GetXmax();
  
   TH1D *hrms = new TH1D("hrms","hrms",nbins,xmin,xmax);
   for( int i = 1; i <= nbins; ++i){
     TH1D *projection = (TH1D*) h->ProjectionY("",i,i);
     projection->Sumw2();
     double rms = projection->GetMean();
     cout<<"rms:"<<rms<<endl;
     double rmserr = projection->GetMeanError();
     
     hrms-> SetBinContent(i,rms);
     cout<<"hrms:"<<hrms<<endl;
     hrms-> SetBinError(i,rmserr);
   }
   TCanvas *c = new TCanvas();
   c->SetLeftMargin(0.15);
   //c->SetLogz();
   gROOT->SetStyle("novaStyle");
   c->SetTickx(1);
   c->SetTicky(1);
   gStyle->SetTitleX(.39);
   gStyle->SetTitleY(0.95);
   hrms->SetLineColor(kBlack);
   hrms->Draw("e1");
   //h->Draw();
   //c->Close();
 f1->Close();

}
#endif