From 35141e85c8a4a377af35b5ae6d5e13c96e872d92 Mon Sep 17 00:00:00 2001 From: Jacob Larkin Date: Tue, 2 Aug 2022 08:40:28 -0500 Subject: [PATCH 1/4] Add new Vars for 202208 Icarus selection --- sbnana/SBNAna/Vars/NumuVarsIcarus202208.cxx | 65 +++++++++++++++++++++ sbnana/SBNAna/Vars/NumuVarsIcarus202208.h | 15 +++++ 2 files changed, 80 insertions(+) create mode 100644 sbnana/SBNAna/Vars/NumuVarsIcarus202208.cxx create mode 100644 sbnana/SBNAna/Vars/NumuVarsIcarus202208.h diff --git a/sbnana/SBNAna/Vars/NumuVarsIcarus202208.cxx b/sbnana/SBNAna/Vars/NumuVarsIcarus202208.cxx new file mode 100644 index 00000000..49fa2ab7 --- /dev/null +++ b/sbnana/SBNAna/Vars/NumuVarsIcarus202208.cxx @@ -0,0 +1,65 @@ +#include "sbnana/SBNAna/Vars/NumuVarsIcarus202208.h" +#include "sbnana/CAFAna/Core/Var.h" + +#include "sbnanaobj/StandardRecord/Proxy/FwdDeclare.h" +#include "sbnanaobj/StandardRecord/Proxy/SRProxy.h" + +namespace ana { + +// More general vars, i.e., not just numu vars. +// Probably should be moved to another file at some point. +const Var kCRLongestTrackDirY = SIMPLEVAR(nuid.crlongtrkdiry); +const Var kFMTimeVar = SIMPLEVAR(fmatch.time); +const Var kFMScoreVar = SIMPLEVAR(fmatch.score); + +// Mostly just copied from NumuVarsIcarus202106.cxx +const Var kIcarus202208MuonIdx([](const caf::SRSliceProxy* slc) -> int { + // The (dis)qualification of a slice is based upon the track level information. + float Longest(0); + int PTrackInd(-1); + for (std::size_t i(0); i < slc->reco.trk.size(); ++i) { + auto const& trk = slc->reco.trk.at(i); + if(trk.bestplane == -1) continue; + const float Atslc = std::hypot(slc->vertex.x - trk.start.x, + slc->vertex.y - trk.start.y, + slc->vertex.z - trk.start.z); + const bool AtSlice = ( Atslc < 10.0 && trk.pfp.parent_is_primary); + + const float Chi2Proton = trk.chi2pid[trk.bestplane].chi2_proton; + const float Chi2Muon = trk.chi2pid[trk.bestplane].chi2_muon; + + const bool Contained = ( !isnan(trk.end.x) && + (( trk.end.x < -71.1 - 5 && trk.end.x > -369.33 + 5 ) || + ( trk.end.x < 71.1 + 5 && trk.end.x > +369.33 - 5 )) && + !isnan(trk.end.y) && + ( trk.end.y > -181.7 + 5 && trk.end.y < 134.8 - 5 ) && + !isnan(trk.end.z) && + ( trk.end.z > -895.95 + 5 && trk.end.z < 895.95 - 5 ) ); + const bool MaybeMuonExiting = ( !Contained && trk.len > 100); + const bool MaybeMuonContained = ( Contained && Chi2Proton > 60 && Chi2Muon < 30 && trk.len > 50 ); + if ( AtSlice && ( MaybeMuonExiting || MaybeMuonContained ) && trk.len > Longest ) { + Longest = trk.len; + PTrackInd = i; + } + } + return PTrackInd; +}); + +static bool proton_cut(const caf::SRTrackProxy& trk) +{ + return trk.chi2pid[trk.bestplane].chi2_proton < 100; +} + +const Var kIcarus202208NumPions([](const caf::SRSliceProxy* slc) +{ + int count = 0; + auto idx = kIcarus202208MuonIdx(slc); + int muID = -1; + if (idx >= 0) muID = slc->reco.trk.at(idx).pfp.id; + for(auto& trk: slc->reco.trk) { + if(trk.pfp.id != muID && !proton_cut(trk) && trk.pfp.parent_is_primary) + ++count; + } + return count; +}); +} diff --git a/sbnana/SBNAna/Vars/NumuVarsIcarus202208.h b/sbnana/SBNAna/Vars/NumuVarsIcarus202208.h new file mode 100644 index 00000000..262fd054 --- /dev/null +++ b/sbnana/SBNAna/Vars/NumuVarsIcarus202208.h @@ -0,0 +1,15 @@ +#pragma once + +#include "sbnana/CAFAna/Core/Var.h" + +namespace ana { + +extern const Var kCRLongestTrackDirY; +extern const Var kFMTimeVar; +extern const Var kFMScoreVar; + +extern const Var kIcarus202208MuonIdx; +extern const Var kIcarus202208NumPions; +} + + From 584d63b56a6050bed1b89071129f380661a170cd Mon Sep 17 00:00:00 2001 From: Jacob Larkin Date: Tue, 2 Aug 2022 08:41:03 -0500 Subject: [PATCH 2/4] Add new Cuts for Icarus 202208 selection --- sbnana/SBNAna/Cuts/NumuCutsIcarus202208.cxx | 51 +++++++++++++++++++++ sbnana/SBNAna/Cuts/NumuCutsIcarus202208.h | 21 +++++++++ 2 files changed, 72 insertions(+) create mode 100644 sbnana/SBNAna/Cuts/NumuCutsIcarus202208.cxx create mode 100644 sbnana/SBNAna/Cuts/NumuCutsIcarus202208.h diff --git a/sbnana/SBNAna/Cuts/NumuCutsIcarus202208.cxx b/sbnana/SBNAna/Cuts/NumuCutsIcarus202208.cxx new file mode 100644 index 00000000..0cca7bfd --- /dev/null +++ b/sbnana/SBNAna/Cuts/NumuCutsIcarus202208.cxx @@ -0,0 +1,51 @@ +#include "sbnana/CAFAna/Core/Cut.h" +#include "sbnana/SBNAna/Cuts/NumuCutsIcarus202208.h" +#include "sbnana/SBNAna/Vars/NumuVarsIcarus202208.h" + +#include "sbnanaobj/StandardRecord/Proxy/FwdDeclare.h" +#include "sbnanaobj/StandardRecord/Proxy/SRProxy.h" + +namespace ana { + +static bool contained(const caf::SRTrackProxy& trk) +{ + return ((trk.end.x < -71.1 - 5 && trk.end.x > -369.33 + 5) + || (trk.end.x > 71.1 + 5 && trk.end.x < 369.33 - 5)) + && trk.end.y > -181.7 + 5 && trk.end.y < 134.8 - 5 + && trk.end.z > -895.95 + 5 && trk.end.z < 895.95 - 5; +} + +const Cut kIcarus202208FMTimeCut = kFMTimeVar > 0 && kFMTimeVar < 1.8; +const Cut kIcarus202208FMScoreCut = kFMScoreVar < 9; +const Cut kIcarus202208LongTrackDirCut = kCRLongestTrackDirY > -0.91; +const Cut kIcarus202208FoundMuon = kIcarus202208MuonIdx >= 0; +const Cut kIcarus202208RecoFiducial([](const caf::SRSliceProxy* slc) { + return ( !isnan(slc->vertex.x) && + ( ( slc->vertex.x < -71.1 - 25 && slc->vertex.x > -369.33 + 25 ) || + ( slc->vertex.x > 71.1 + 25 && slc->vertex.x < 369.33 - 25 ) ) && + !isnan(slc->vertex.y) && + ( slc->vertex.y > -181.7 + 25 && slc->vertex.y < 134.8 - 25 ) && + !isnan(slc->vertex.z) && + ( slc->vertex.z > -895.95 + 30 && slc->vertex.z < 895.95 - 50 ) ); + }); + +const Cut kIcarus202208NumuSelection = kIcarus202208RecoFiducial && kIcarus202208FMScoreCut && kIcarus202208FMTimeCut && kIcarus202208LongTrackDirCut && kIcarus202208FoundMuon; + +const Cut kIcarus202208NoPion = kIcarus202208NumPions == 0; +const Cut kIcarus202208ContainedHadrons([](const caf::SRSliceProxy* slc){ + auto idx = kIcarus202208MuonIdx(slc); + int muID = -1; + if (idx >= 0) muID = slc->reco.trk.at(idx).pfp.id; + for(auto& trk: slc->reco.trk) { + if(trk.pfp.id != muID && trk.pfp.parent_is_primary) + if(!contained(trk)) return false; + } + return true; +}); +const Cut kIcarus202208ContainedMuon([](const caf::SRSliceProxy* slc){ + return kIcarus202208FoundMuon(slc) && contained(slc->reco.trk.at(kIcarus202208MuonIdx(slc))); +}); +const Cut kIcarus202208ContainedMuonAndHadrons = kIcarus202208ContainedMuon && kIcarus202208ContainedHadrons; +const Cut kIcarus202208QELike = kIcarus202208NumuSelection && kIcarus202208NoPion && kIcarus202208ContainedHadrons; +const Cut kIcarus202208QELikeContainedMuon = kIcarus202208QELike && kIcarus202208ContainedMuon; +} diff --git a/sbnana/SBNAna/Cuts/NumuCutsIcarus202208.h b/sbnana/SBNAna/Cuts/NumuCutsIcarus202208.h new file mode 100644 index 00000000..8a107abc --- /dev/null +++ b/sbnana/SBNAna/Cuts/NumuCutsIcarus202208.h @@ -0,0 +1,21 @@ +#pragma once + +#include "sbnana/CAFAna/Core/Cut.h" + +namespace ana { + +extern const Cut kIcarus202208FMTimeCut; +extern const Cut kIcarus202208FMScoreCut; +extern const Cut kIcarus202208LongTrackDirCut; +extern const Cut kIcarus202208FoundMuon; +extern const Cut kIcarus202208RecoFiducial; + +extern const Cut kIcarus202208NumuSelection; + +extern const Cut kIcarus202208NoPion; +extern const Cut kIcarus202208ContainedHadrons; +extern const Cut kIcarus202208ContainedMuon; +extern const Cut kIcarus202208ContainedMuonAndHadrons; +extern const Cut kIcarus202208QELike; +extern const Cut kIcarus202208QELikeContainedMuon; +} From c0dcffc230c600c79d6fce44d615ffd6acf2f00d Mon Sep 17 00:00:00 2001 From: Jacob Larkin Date: Tue, 2 Aug 2022 08:41:55 -0500 Subject: [PATCH 3/4] Add example file using Icarus 202208 selection --- .../Analysis/202208Scanning/selected_events.C | 38 +++++++++++++++++++ 1 file changed, 38 insertions(+) create mode 100644 sbnana/SBNAna/Analysis/202208Scanning/selected_events.C diff --git a/sbnana/SBNAna/Analysis/202208Scanning/selected_events.C b/sbnana/SBNAna/Analysis/202208Scanning/selected_events.C new file mode 100644 index 00000000..75f0095a --- /dev/null +++ b/sbnana/SBNAna/Analysis/202208Scanning/selected_events.C @@ -0,0 +1,38 @@ +#include "sbnana/CAFAna/Core/Binning.h" +#include "sbnana/CAFAna/Core/Cut.h" +#include "sbnana/CAFAna/Core/Spectrum.h" +#include "sbnana/CAFAna/Core/SpectrumLoader.h" +#include "sbnana/CAFAna/Core/Var.h" + +#include "sbnana/SBNAna/Cuts/NumuCutsIcarus202208.h" +#include "sbnana/SBNAna/Vars/NumuVarsIcarus202208.h" + +#include "sbnanaobj/StandardRecord/Proxy/SRProxy.h" + +#include + +using namespace ana; + +void selected_events() +{ + SpectrumLoader loader("icarus_BNB_Nu_Cosmics_v09_37_02_04_flatcaf"); + + std::ofstream out("selected_events.txt"); + + const SpillVar dummy_var([&out](const caf::SRSpillProxy* spill){ + for(const auto& slc: spill->slc) { + if(kIcarus202208QELikeContainedMuon(&slc)) { + out << "Run: " << spill->hdr.run << "\tSubrun: " << spill->hdr.subrun + << "\tEvent: " << spill->hdr.evt << "\tVertex: (" + << slc.vertex.x << ", " << slc.vertex.y << ", " << slc.vertex.z << ")\n"; + return 1; + } + } + return 0; + }); + + const Binning dummy_bins = Binning::Simple(2, 0, 2); + Spectrum dummy_spec("Dummy Label", dummy_bins, loader, dummy_var, kNoSpillCut); + + loader.Go(); +} From b8e0e466b411a1beb5f2864fdaffdca096bd4063 Mon Sep 17 00:00:00 2001 From: Jacob Larkin Date: Thu, 4 Aug 2022 09:54:43 -0500 Subject: [PATCH 4/4] Add Icarus202208 to two places where it was missing --- sbnana/SBNAna/Cuts/NumuCutsIcarus202208.cxx | 7 +++---- sbnana/SBNAna/Vars/NumuVarsIcarus202208.cxx | 4 ++-- 2 files changed, 5 insertions(+), 6 deletions(-) diff --git a/sbnana/SBNAna/Cuts/NumuCutsIcarus202208.cxx b/sbnana/SBNAna/Cuts/NumuCutsIcarus202208.cxx index 0cca7bfd..a69686ad 100644 --- a/sbnana/SBNAna/Cuts/NumuCutsIcarus202208.cxx +++ b/sbnana/SBNAna/Cuts/NumuCutsIcarus202208.cxx @@ -2,12 +2,11 @@ #include "sbnana/SBNAna/Cuts/NumuCutsIcarus202208.h" #include "sbnana/SBNAna/Vars/NumuVarsIcarus202208.h" -#include "sbnanaobj/StandardRecord/Proxy/FwdDeclare.h" #include "sbnanaobj/StandardRecord/Proxy/SRProxy.h" namespace ana { -static bool contained(const caf::SRTrackProxy& trk) +static bool Icarus202208contained(const caf::SRTrackProxy& trk) { return ((trk.end.x < -71.1 - 5 && trk.end.x > -369.33 + 5) || (trk.end.x > 71.1 + 5 && trk.end.x < 369.33 - 5)) @@ -38,12 +37,12 @@ const Cut kIcarus202208ContainedHadrons([](const caf::SRSliceProxy* slc){ if (idx >= 0) muID = slc->reco.trk.at(idx).pfp.id; for(auto& trk: slc->reco.trk) { if(trk.pfp.id != muID && trk.pfp.parent_is_primary) - if(!contained(trk)) return false; + if(!Icarus202208contained(trk)) return false; } return true; }); const Cut kIcarus202208ContainedMuon([](const caf::SRSliceProxy* slc){ - return kIcarus202208FoundMuon(slc) && contained(slc->reco.trk.at(kIcarus202208MuonIdx(slc))); + return kIcarus202208FoundMuon(slc) && Icarus202208contained(slc->reco.trk.at(kIcarus202208MuonIdx(slc))); }); const Cut kIcarus202208ContainedMuonAndHadrons = kIcarus202208ContainedMuon && kIcarus202208ContainedHadrons; const Cut kIcarus202208QELike = kIcarus202208NumuSelection && kIcarus202208NoPion && kIcarus202208ContainedHadrons; diff --git a/sbnana/SBNAna/Vars/NumuVarsIcarus202208.cxx b/sbnana/SBNAna/Vars/NumuVarsIcarus202208.cxx index 49fa2ab7..5e76addf 100644 --- a/sbnana/SBNAna/Vars/NumuVarsIcarus202208.cxx +++ b/sbnana/SBNAna/Vars/NumuVarsIcarus202208.cxx @@ -45,7 +45,7 @@ const Var kIcarus202208MuonIdx([](const caf::SRSliceProxy* slc) -> int { return PTrackInd; }); -static bool proton_cut(const caf::SRTrackProxy& trk) +static bool Icarus202208_proton_cut(const caf::SRTrackProxy& trk) { return trk.chi2pid[trk.bestplane].chi2_proton < 100; } @@ -57,7 +57,7 @@ const Var kIcarus202208NumPions([](const caf::SRSliceProxy* slc) int muID = -1; if (idx >= 0) muID = slc->reco.trk.at(idx).pfp.id; for(auto& trk: slc->reco.trk) { - if(trk.pfp.id != muID && !proton_cut(trk) && trk.pfp.parent_is_primary) + if(trk.pfp.id != muID && !Icarus202208_proton_cut(trk) && trk.pfp.parent_is_primary) ++count; } return count;