#include <NuCutImps.h>
Public Member Functions | |
CC2014Resolution (const NuPlots *plots, TString ana, int quantile) | |
virtual | ~CC2014Resolution () |
virtual Bool_t | InFidVol (const NuEvent &nu) const |
void | Preselection (const NuEvent &nu) |
void | Selection (const NuEvent &nu) |
Protected Member Functions | |
void | Resolution (const NuEvent &nu) |
void | SetCutFile (const char *filename) |
Protected Attributes | |
TTree * | fResCuts_r1 |
TTree * | fResCuts_r2 |
TTree * | fResCuts_r3 |
TTree * | fResCuts_phe |
int | fQuantile |
Private Attributes | |
CC2014 | fCC2014 |
Definition at line 288 of file NuCutImps.h.
NuCutImps::CC2014Resolution::CC2014Resolution | ( | const NuPlots * | plots, | |
TString | ana, | |||
int | quantile | |||
) |
Definition at line 1939 of file NuCutImps.cxx.
References fQuantile, NuLibrary::general, NuGeneral::GetReleaseDirToUse(), NuLibrary::Instance(), NuCuts::kCC2014, NuCut::SetAnaVersion(), SetCutFile(), and NuCut::SetFidVol().
01942 : NuCut(ana, plots), fResCuts_r1(0), fResCuts_r2(0), 01943 fResCuts_r3(0), fResCuts_phe(0), fQuantile(quantile) 01944 { 01945 assert(fQuantile >= 0 && fQuantile <= 4); 01946 01947 SetFidVol("cc2008"); 01948 SetAnaVersion(NuCuts::kCC2014); 01949 01950 const NuGeneral& gen = NuLibrary::Instance().general; 01951 const TString releaseDir = gen.GetReleaseDirToUse("NtupleUtils"); 01952 const TString fullPath = releaseDir+"/NtupleUtils/data/cuts_numu.root"; 01953 SetCutFile(fullPath); 01954 }
NuCutImps::CC2014Resolution::~CC2014Resolution | ( | ) | [virtual] |
Definition at line 1956 of file NuCutImps.cxx.
Bool_t NuCutImps::CC2014Resolution::InFidVol | ( | const NuEvent & | nu | ) | const [virtual] |
Fiducial volume calculation function. This is the basic 'infid' derived version, and should be overridden/passed to for any more complicated evaluations.
Reimplemented from NuCut.
Definition at line 1960 of file NuCutImps.cxx.
References fCC2014, and NuCutImps::CC2014::InFidVol().
01961 { 01962 return fCC2014.InFidVol(nu); 01963 }
void NuCutImps::CC2014Resolution::Preselection | ( | const NuEvent & | nu | ) | [virtual] |
Implements NuCut.
Definition at line 1965 of file NuCutImps.cxx.
References NuCut::Defer_Preselection(), and fCC2014.
01966 { 01967 Defer_Preselection(fCC2014, nu); 01968 }
void NuCutImps::CC2014Resolution::Resolution | ( | const NuEvent & | nu | ) | [protected] |
Definition at line 1998 of file NuCutImps.cxx.
References NuEvent::beamTypeDB, NuEvent::energy, fQuantile, fResCuts_phe, fResCuts_r1, fResCuts_r2, fResCuts_r3, NuCut::Keep_If(), Msg::kInfo, BeamType::kL250z200i, MAXMSG, NuEvent::resolution, and NuEvent::runPeriod.
Referenced by Selection().
01999 { 02000 Double_t low_en_edge = 0; 02001 Double_t high_en_edge = 0; 02002 Double_t low_res_edge = 0; 02003 Double_t high_res_edge = 0; 02004 Double_t quantiles[6]; 02005 TTree *resCuts; 02006 resCuts=fResCuts_r3; 02007 //Different cuts per run period! 02008 switch(nu.runPeriod){ 02009 case 1: 02010 resCuts=fResCuts_r1; 02011 break; 02012 case 2: 02013 resCuts=fResCuts_r2; 02014 break; 02015 case 3: 02016 02017 resCuts=fResCuts_r3; 02018 break; 02019 //default: 02020 //assert(0 && "Badly configured runPeriod"); 02021 } 02022 //Whatever run period this is nowadays... 02023 if(nu.beamTypeDB == BeamType::kL250z200i) { 02024 resCuts=fResCuts_phe; 02025 MAXMSG("NuCutImps",Msg::kInfo,10)<<"Using pHE resolution cuts"<<endl; 02026 } 02027 02028 resCuts->SetBranchAddress("quantiles", quantiles); 02029 resCuts->SetBranchAddress("lowedge", &low_en_edge); 02030 resCuts->SetBranchAddress("highedge", &high_en_edge); 02031 02032 for(Int_t i = 0; i < resCuts->GetEntries(); i++){ 02033 resCuts->GetEntry(i); 02034 02035 if(fQuantile == 0) low_res_edge = 0.0; 02036 else low_res_edge = quantiles[fQuantile]; 02037 02038 if(fQuantile == 4) high_res_edge = 9999999.; 02039 else high_res_edge = quantiles[fQuantile+1]; 02040 02041 if(nu.energy > low_en_edge && nu.energy <= high_en_edge) 02042 { 02043 MAXMSG("NuCutImps", Msg::kInfo, 20) 02044 << "Entry:" << i << " Energy:" << nu.energy << " sigma/E:" 02045 << (nu.resolution/nu.energy) << endl 02046 << "Keep if between " << low_res_edge << " and " 02047 << high_res_edge << endl; 02048 02049 Keep_If((nu.resolution/nu.energy) > low_res_edge && 02050 (nu.resolution/nu.energy) <= high_res_edge, "resolution"); 02051 return; 02052 } 02053 } 02054 // God knows how we got here. We'd better keep the event, but not 02055 // double-count it, so bung it in quantile 4. 02056 Keep_If(fQuantile == 4, "Overflow"); 02057 }
void NuCutImps::CC2014Resolution::Selection | ( | const NuEvent & | nu | ) | [virtual] |
Implements NuCut.
Definition at line 1970 of file NuCutImps.cxx.
References NuCut::Defer_Selection(), NuEvent::detector, fCC2014, Detector::kFar, and Resolution().
01971 { 01972 Defer_Selection(fCC2014, nu); 01973 01974 if(nu.detector == Detector::kFar/* && nu.charge<0*/) Resolution(nu); 01975 }
void NuCutImps::CC2014Resolution::SetCutFile | ( | const char * | filename | ) | [protected] |
Definition at line 1977 of file NuCutImps.cxx.
References fResCuts_phe, fResCuts_r1, fResCuts_r2, fResCuts_r3, Msg::kInfo, Msg::kWarning, and MSG.
Referenced by CC2014Resolution().
01978 { 01979 TDirectory* tmpd = gDirectory; //get the right path 01980 TFile* cutfile = new TFile(filename, "READ"); 01981 gDirectory = tmpd; 01982 if(cutfile->IsZombie()){ 01983 MSG("NuCutImps", Msg::kWarning) << "cutfile does not exist" << endl; 01984 return; 01985 } 01986 //I know, I know, this isn't a great way to do it... 01987 01988 fResCuts_r1 = (TTree*)cutfile->Get("cuts_r1"); 01989 fResCuts_r2 = (TTree*)cutfile->Get("cuts_r2"); 01990 fResCuts_r3 = (TTree*)cutfile->Get("cuts_r3"); 01991 fResCuts_phe = (TTree*)cutfile->Get("cuts_phe"); 01992 01993 MSG("NuCutImps", Msg::kInfo) << "cutfile exists: " 01994 << fResCuts_r1->GetEntries() 01995 << " R1 entries" << endl; 01996 }
CC2014 NuCutImps::CC2014Resolution::fCC2014 [private] |
Definition at line 290 of file NuCutImps.h.
Referenced by InFidVol(), Preselection(), and Selection().
int NuCutImps::CC2014Resolution::fQuantile [protected] |
Definition at line 309 of file NuCutImps.h.
Referenced by CC2014Resolution(), and Resolution().
TTree* NuCutImps::CC2014Resolution::fResCuts_phe [protected] |
Definition at line 308 of file NuCutImps.h.
Referenced by Resolution(), and SetCutFile().
TTree* NuCutImps::CC2014Resolution::fResCuts_r1 [protected] |
Definition at line 305 of file NuCutImps.h.
Referenced by Resolution(), and SetCutFile().
TTree* NuCutImps::CC2014Resolution::fResCuts_r2 [protected] |
Definition at line 306 of file NuCutImps.h.
Referenced by Resolution(), and SetCutFile().
TTree* NuCutImps::CC2014Resolution::fResCuts_r3 [protected] |
Definition at line 307 of file NuCutImps.h.
Referenced by Resolution(), and SetCutFile().