reduce_numu_ana2017_prod.C
Go to the documentation of this file.
1 #include "CAFAna/Cuts/Cuts.h"
4 #include "CAFAna/Vars/Vars.h"
7 #include "CAFAna/Core/Cut.h"
8 
9 using namespace ana;
10 
12 {
13  FileReducer reduce(fname, outname);
14 
15  reduce.SetMetadata("nova.decaf_skim","numu_contain");
16 
17  const Cut kIsND({"hdr.det"},
18  [](const caf::StandardRecord* sr)
19  {
20  return sr->hdr.det==caf::kNEARDET;
21  });
22 
23  const Cut kIsFD({"hdr.det"},
24  [](const caf::StandardRecord* sr)
25  {
26  return sr->hdr.det==caf::kFARDET;
27  });
28 
29  const Cut kNumuContainND2017({"vtx.elastic.IsValid","vtx.elastic.*","trk.kalman.*","slc.firstplane","slc.lastplane","sel.contain.*"},
30  [](const caf::StandardRecord* sr)
31  { if( !sr->vtx.elastic.IsValid ) return false;
32  // reconstructed showers all contained
33  for( unsigned int i = 0; i < sr->vtx.elastic.fuzzyk.nshwlid; ++i ) {
34  TVector3 start = sr->vtx.elastic.fuzzyk.png[i].shwlid.start;
35  TVector3 stop = sr->vtx.elastic.fuzzyk.png[i].shwlid.stop;
36  if( std::min( start.X(), stop.X() ) < -180.0 ) return false;
37  if( std::max( start.X(), stop.X() ) > 180.0 ) return false;
38  if( std::min( start.Y(), stop.Y() ) < -180.0 ) return false;
39  if( std::max( start.Y(), stop.Y() ) > 180.0 ) return false;
40  if( std::min( start.Z(), stop.Z() ) < 20.0 ) return false;
41  if( std::max( start.Z(), stop.Z() ) > 1525.0 ) return false;
42  }
43 
44  // only primary muon track present in muon catcher
45  if( sr->trk.kalman.ntracks < 1 ) return false;
46  for( unsigned int i = 0; i < sr->trk.kalman.ntracks; ++i ) {
47  if( i == sr->trk.kalman.idxremid ) continue;
48  else if( sr->trk.kalman.tracks[i].start.Z() > 1275 ||
49  sr->trk.kalman.tracks[i].stop.Z() > 1275 )
50  return false;
51  }
52 
53  return ( sr->trk.kalman.ntracks > sr->trk.kalman.idxremid
54  && sr->slc.firstplane > 1 // skip 0 and 1
55  && sr->slc.lastplane < 212 // skip 212 and 213
56  && sr->trk.kalman.tracks[0].start.Z() < 1100
57  // vertex definitely outside mC
58  && ( sr->trk.kalman.tracks[0].stop.Z() < 1275
59  || sr->sel.contain.kalyposattrans < 55 ) // air gap
60  && sr->sel.contain.kalfwdcellnd > 5
61  && sr->sel.contain.kalbakcellnd > 10 );
62  }
63  );
64 
65  const Cut kNumuDecafPresel2017({"slc.nhit","sel.remid.pid","sel.cvn.numuid"},
66  [](const caf::StandardRecord* sr)
67  { return (sr->slc.nhit > 20
68  && sr->sel.remid.pid > 0
69  && sr->sel.cvn.numuid > 0);
70  }
71  );
72 
73  const Cut kNumuCosmicRej2017({"sel.cosrej.anglekal","sel.cosrej.numucontpid","slc.nhit"},
74  [](const caf::StandardRecord* sr)
75  { return ( sr->sel.cosrej.anglekal > 0.5 &&
76  sr->sel.cosrej.numucontpid > 0.5 &&
77  sr->slc.nhit < 400);
78  }
79  );
80 
81  const Cut kNumuOptimizedContainFD2017({"hdr.dibmask","slc.firstplane","slc.lastplane","sel.contain.*"},
82  [](const caf::StandardRecord* sr)
83  {
84  std::pair<int, int> planes = calcFirstLastLivePlane(sr->slc.firstplane, std::bitset<14>(sr->hdr.dibmask));
85  int planestofront = sr->slc.firstplane - planes.first;
86  int planestoback = planes.second - sr->slc.lastplane;
87  return (
88  sr->sel.contain.kalfwdcell > 6
89  && sr->sel.contain.kalbakcell > 6
90  && sr->sel.contain.cosfwdcell > 0
91  && sr->sel.contain.cosbakcell > 7
92  && planestofront > 1
93  && planestoback > 1
94  );
95  }
96  );
97  //Prongs cut for FD
99  //Numu Containment for FD
100  const Cut kNumuContainFD2017 = kNumuProngsContainFD2017 && kNumuOptimizedContainFD2017;
101 
102  const Cut kNumuSACosmicRejLoose({"sel.cosrej.numuSAcontpid"},
103  [](const caf::StandardRecord* sr)
104  { return ( sr->sel.cosrej.numuSAcontpid > 0.5 );
105  }
106  );
107 
108  const Cut kNumuDecafCutFD2017 =
109  kNumuDecafPresel2017 && kNumuContainFD2017 && ( kNumuCosmicRejLoose || kNumuSACosmicRejLoose );
110  const Cut kNumuDecafCutND2017 =
112 
113  const Cut kNumuDecafCut = (kIsND && kNumuDecafCutND2017) || (kIsFD && kNumuDecafCutFD2017);
114 
115  reduce.AddEventCut(kNumuDecafCut);
116 
118 
119  reduce.Go();
120 }
int cosfwdcell
cosmictrack projected # cells from end point forwards to det edge
Definition: SRContain.h:30
Det_t det
Detector, ND = 1, FD = 2, NDOS = 3.
Definition: SRHeader.h:28
Near Detector underground.
Definition: SREnums.h:10
T max(const caf::Proxy< T > &a, T b)
Far Detector at Ash River.
Definition: SREnums.h:11
Cuts and Vars for the 2020 FD DiF Study.
Definition: vars.h:6
SRHeader hdr
Header branch: run, subrun, etc.
const Var kDistAllBottom([](const caf::SRProxy *sr){if(std::isnan(1.*sr->sel.nuecosrej.distallpngbottom)) return-1000.0f;return float(sr->sel.nuecosrej.distallpngbottom);})
Distance of all showers in slice from the bottom edge of detector.
Definition: NueVars.h:33
const Var kDistAllWest([](const caf::SRProxy *sr){if(std::isnan(1.*sr->sel.nuecosrej.distallpngwest)) return-1000.0f;return float(sr->sel.nuecosrej.distallpngwest);})
Distance of all showers in slice from the west edge of detector.
Definition: NueVars.h:36
unsigned int lastplane
last plane
Definition: SRSlice.h:27
SRCosRej cosrej
Output from CosRej (Cosmic Rejection)
Definition: SRIDBranch.h:47
unsigned int firstplane
first plane
Definition: SRSlice.h:26
size_t ntracks
Definition: SRKalman.h:23
const Var kDistAllTop([](const caf::SRProxy *sr){if(std::isnan(1.*sr->sel.nuecosrej.distallpngtop)) return-1000.0f;return float(sr->sel.nuecosrej.distallpngtop);})
Distance of all showers in slice from the top edge of detector.
Definition: NueVars.h:30
unsigned int idxremid
index number of the best ReMId track
Definition: SRKalman.h:32
const Cut kNumuContainND2017([](const caf::SRProxy *sr){if(!sr->vtx.elastic.IsValid) return false;for(unsigned int i=0;i< sr->vtx.elastic.fuzzyk.nshwlid;++i){TVector3 start=sr->vtx.elastic.fuzzyk.png[i].shwlid.start;TVector3 stop=sr->vtx.elastic.fuzzyk.png[i].shwlid.stop;if(std::min(start.X(), stop.X())< -180.0) return false;if(std::max(start.X(), stop.X()) > 180.0) return false;if(std::min(start.Y(), stop.Y())< -180.0) return false;if(std::max(start.Y(), stop.Y()) > 180.0) return false;if(std::min(start.Z(), stop.Z())< 20.0) return false;if(std::max(start.Z(), stop.Z()) > 1525.0) return false;}if(sr->trk.kalman.ntracks< 1) return false;for(unsigned int i=0;i< sr->trk.kalman.ntracks;++i){if(i==sr->trk.kalman.idxremid) continue;else if(sr->trk.kalman.tracks[i].start.Z() > 1275||sr->trk.kalman.tracks[i].stop.Z() > 1275) return false;}return(sr->trk.kalman.ntracks > sr->trk.kalman.idxremid &&sr->slc.firstplane > 1 &&sr->slc.lastplane< 212 &&sr->trk.kalman.tracks[0].start.Z()< 1100 &&(sr->trk.kalman.tracks[0].stop.Z()< 1275 ||sr->sel.contain.kalyposattrans< 55) &&sr->sel.contain.kalfwdcellnd > 5 &&sr->sel.contain.kalbakcellnd > 10);})
Definition: NumuCuts2017.h:11
std::vector< SRFuzzyKProng > png
Vector of 3D prong objects.
Definition: SRFuzzyK.h:19
const Var kDistAllBack([](const caf::SRProxy *sr){if(std::isnan(1.*sr->sel.nuecosrej.distallpngback)) return-1000.0f;return float(sr->sel.nuecosrej.distallpngback);})
Distance of all showers in slice from the back edge of detector.
Definition: NueVars.h:42
float kalyposattrans
Y position of Kalman track and transition (ND only, use to check if went through air gap) ...
Definition: SRContain.h:50
const Cut kNumuSACosmicRejLoose([](const caf::SRProxy *sr){return( sr->sel.cosrej.numucontpid2019 > 0.5);})
Definition: NumuCuts.h:36
void SetMetadata(const std::string &key, const std::string &val)
Cause the output file to be systematically shifted.
Definition: FileReducer.h:102
int cosbakcell
cosmictrack projected # cells from start point backwards to det edge
Definition: SRContain.h:35
int kalbakcellnd
Kalmantrack projected # cells from start point backwards to det edge, including muon catcher...
Definition: SRContain.h:46
const Var kDistAllEast([](const caf::SRProxy *sr){if(std::isnan(1.*sr->sel.nuecosrej.distallpngeast)) return-1000.0f;return float(sr->sel.nuecosrej.distallpngeast);})
Distance of all showers in slice from the east edge of detector.
Definition: NueVars.h:39
void AddReductionStep(const std::function< ReductionFunc > &f)
Run the specified reduction function over each event.
Definition: FileReducer.h:85
SRKalman kalman
Tracks produced by KalmanTrack.
Definition: SRTrackBranch.h:24
float pid
PID value output by kNN.
Definition: SRRemid.h:25
int kalfwdcellnd
Kalmantrack projected # cells from end point forwards to det edge, including muon catcher...
Definition: SRContain.h:41
unsigned int nshwlid
number of shwlid showers - either 0 or number of 3d prongs
Definition: SRFuzzyK.h:24
const Cut kNumuContainFD2017
Definition: NumuCuts2017.h:21
float anglekal
cos of angle of best ReMId Kalman track
Definition: SRCosRej.h:23
SRCVNResult cvn
Horrible hack to appease CAFAna.
Definition: SRIDBranch.h:53
const Cut kNumuOptimizedContainFD2017([](const caf::SRProxy *sr){std::pair< int, int > planes=calcFirstLastLivePlane(sr->slc.firstplane, std::bitset< 14 >(sr->hdr.dibmask));int planestofront=sr->slc.firstplane-planes.first;int planestoback=planes.second-sr->slc.lastplane;return(sr->sel.contain.kalfwdcell > 6 &&sr->sel.contain.kalbakcell > 6 &&sr->sel.contain.cosfwdcell > 0 &&sr->sel.contain.cosbakcell > 7 &&planestofront > 1 &&planestoback > 1);})
Definition: NumuCuts2017.h:16
SRRemid remid
Output from RecoMuonID (ReMId) package.
Definition: SRIDBranch.h:39
void AddEventCut(const Cut &cut)
Only copy records to the output file if they pass this cut.
const Cut kNumuProngsContainFD2017
Definition: NumuCuts2017.h:19
caf::StandardRecord * sr
const Cut kNumuDecafCutFD2017
Definition: NumuCuts2017.h:33
const Cut kNumuDecafPresel2017([](const caf::SRProxy *sr){return(sr->slc.nhit > 20 &&sr->sel.remid.pid > 0 &&sr->sel.cvn.numuid > 0);})
Definition: NumuCuts2017.h:13
unsigned int nhit
number of hits
Definition: SRSlice.h:22
virtual void Go() override
Load all the registered spectra.
void reduce_numu_ana2017_prod(std::string fname, std::string outname)
Create smaller CAFs.
Definition: FileReducer.h:56
const Cut kNumuCosmicRej2017([](const caf::SRProxy *sr){return(sr->sel.cosrej.anglekal > 0.5 && sr->sel.cosrej.numucontpid2020 > 0.5 && sr->slc.nhit< 400);})
Definition: NumuCuts2017.h:24
The StandardRecord is the primary top-level object in the Common Analysis File trees.
const Var kDistAllFront([](const caf::SRProxy *sr){if(std::isnan(1.*sr->sel.nuecosrej.distallpngfront)) return-1000.0f;return float(sr->sel.nuecosrej.distallpngfront);})
Distance of all showers in slice from the front edge of detector.
Definition: NueVars.h:45
int kalfwdcell
Kalmantrack projected # cells from end point forwards to det edge.
Definition: SRContain.h:40
SRIDBranch sel
Selector (PID) branch.
SRElastic elastic
Single vertex found by Elastic Arms.
float numuid
Likelihood Charge Current NuMu.
Definition: SRCVNResult.h:20
std::pair< int, int > calcFirstLastLivePlane(int plane, std::bitset< 14 > binary, caf::Det_t det)
int kalbakcell
Kalmantrack projected # cells from start point backwards to det edge.
Definition: SRContain.h:45
SRSlice slc
Slice branch: nhit, extents, time, etc.
SRFuzzyK fuzzyk
Primary 3D prong object.
Definition: SRElastic.h:44
T min(const caf::Proxy< T > &a, T b)
SRTrackBranch trk
Track branch: nhit, len, etc.
SRContain contain
Output from SRContain (containment related variables)
Definition: SRIDBranch.h:51
const Cut kNumuDecafCutND2017
Definition: NumuCuts2017.h:35
void ClearMultiNuInfo(caf::StandardRecord *sr)
Definition: FileReducer.cxx:25
std::vector< SRKalmanTrack > tracks
3D Tracks produced by KalmanTrack
Definition: SRKalman.h:16
const Cut kNumuCosmicRejLoose([](const caf::SRProxy *sr){return(sr->sel.cosrej.numucontpid2020 > 0.48);})
Definition: NumuCuts.h:34
SRVertexBranch vtx
Vertex branch: location, time, etc.
unsigned short dibmask
diblock mask (bitfield, lowest bit = diblock 1)
Definition: SRHeader.h:34
enum BeamMode string