67 std::set<int>
ids)
const;
114 mf::LogWarning(
"RecoCheckAna") <<
"attempting to run MC truth check on " 115 <<
"real data, bail";
133 for(
unsigned int islc = 0; islc < hslices->size(); ++islc){
136 if(fmtrack.isValid()){
137 std::vector< art::Ptr<rb::Track> > trkcol = fmtrack.at(islc);
138 std::vector< art::Ptr<rb::Cluster> > tracks;
139 for(
unsigned int itrk = 0; itrk < trkcol.size(); ++itrk){
140 tracks.push_back(trkcol[itrk]);
142 for(
unsigned int n = 0;
n <
fPlots.size(); ++
n){
152 if(fmprong.isValid()){
153 std::vector< art::Ptr<rb::Prong> > shwcol = fmprong.at(islc);
154 std::vector< art::Ptr<rb::Cluster> > shws;
155 for(
unsigned int ishw = 0; ishw < shwcol.size(); ++ishw){
156 shws.push_back(shwcol[ishw]);
187 for(
int n = 0;
n < 4; ++
n){
190 std::set<int> allowedDaughters;
198 pdgs.insert(-13); pdgs.insert(+13);
199 allowedDaughters.insert(-11); allowedDaughters.insert(+11);
203 pdgs.insert(-2212); pdgs.insert(+2212);
206 name =
"piPlusMinus";
207 pdgs.insert(-211); pdgs.insert(+211);
212 p.
purity = tfs->
make<TH1D>((name+
"Purity").c_str(),
";Purity;Tracks", 110, 0, 1.1);
213 p.
efficiency = tfs->
make<TH1D>((name+
"Efficiency").c_str(),
";Efficiency;Tracks", 110, 0, 1.1);
214 p.
deltaLength = tfs->
make<TH1D>((name+
"DeltaLength").c_str(),
";#DeltaL (cm);Tracks", 100, -50, +50);
215 p.
purityVsLength = tfs->
make<TH2D>((name+
"PurityVsLength").c_str(),
";True length (cm);Purity;Tracks", 300, 0, 1500, 110, 0, 1.1);
216 p.
efficiencyVsLength = tfs->
make<TH2D>((name+
"EfficiencyVsLength").c_str(),
";True length (cm);Efficiency;Tracks", 300, 0, 1500, 110, 0, 1.1);
217 p.
deltaLengthVsLength = tfs->
make<TH2D>((name+
"DeltaLengthVsLength").c_str(),
";True length (cm);#DeltaL (cm);Tracks", 300, 0, 1500, 100, -50, +50);
218 p.
effvpur = tfs->
make<TH2D>((name+
"EffvPurity").c_str(),
";Efficiency;Purity", 110, 0, 1.1, 110, 0, 1.1);
247 if(!plots.
pdgs.empty()){
248 std::set<int> newTrackIDs;
249 for(std::set<int>::iterator
it = trackIDs.begin();
it != trackIDs.end(); ++
it){
251 if(plots.
pdgs.find(pdg) != plots.
pdgs.end()) newTrackIDs.insert(*
it);
253 trackIDs = newTrackIDs;
256 std::map<int, int> parents;
258 for(std::set<int>::iterator
it = trackIDs.begin();
it != trackIDs.end(); ++
it){
260 if(!mother)
continue;
274 std::map<int, double> bestPur[3], bestEff[3], bestLength[3];
278 for(
size_t c = 0;
c < clscol.size(); ++
c){
287 std::map<int, double> purMap, effMap;
291 for(std::set<int>::iterator
it = trackIDs.begin();
it != trackIDs.end(); ++
it){
292 const double p = purMap[*
it];
293 const double e = effMap[*
it];
310 for(std::set<int>::iterator
it = trackIDs.begin();
it != trackIDs.end(); ++
it){
343 mf::LogWarning(
"RecoCheckAna") <<
"Particle has no trajectory points";
351 for(
unsigned int n = 0;
n < N-1; ++
n){
354 if(p0.Z() < 0 || p1.Z() < 0 ||
367 dist += (p1-
p0).
Mag();
377 nHits[0] = nHits[1] = 0;
378 nOnlyHits[0] = nOnlyHits[1] = 0;
388 std::set<int>
ids)
const 390 std::map<int, Counts> counts;
395 for(
unsigned int i = 0;
i < allhits.
size(); ++
i){
398 const std::vector<TrackIDE> trackIDs = bt->
HitToTrackIDE(hit);
402 for(
unsigned int j = 0;
j < trackIDs.size(); ++
j){
403 const int trackID = trackIDs[
j].trackID;
410 if(ids.find(trackID) != ids.end()){
415 ++counts[loneID].nOnlyHits[hit->
View()];
419 const int kMinHits = 2;
420 const int kMinOnlyHit = 1;
423 for(std::set<int>::iterator
it = ids.begin();
it != ids.end(); ++
it){
::xsd::cxx::tree::id< char, ncname > id
virtual geo::View_t View() const
kXorY for 3D clusters.
unsigned int NumberTrajectoryPoints() const
back track the reconstruction to the simulation
const TLorentzVector & Position(const int i=0) const
std::vector< Plots > fPlots
TH2D * efficiencyVsLength
fvar< T > fabs(const fvar< T > &x)
std::vector< TrackIDE > HitToTrackIDE(const rb::CellHit &hit, bool useBirksE=false) const
Convenience function. HitsToTrackIDE but for a single hit.
static const int kNoiseId
flag for noise id
TH2D * deltaLengthVsLength
enum geo::_plane_proj View_t
Enumerate the possible plane projections.
Vertical planes which measure X.
A collection of associated CellHits.
RecoCheckAna(fhicl::ParameterSet const &p)
A rb::Prong with full reconstructed trajectory.
std::set< int > allowedDaughters
To not count eg Michel electrons against purity.
DEFINE_ART_MODULE(TestTMapFile)
int NumberDaughters() const
Horizontal planes which measure Y.
art::PtrVector< rb::CellHit > AllCells() const
Get all cells from both views.
double HitCollectionEfficiency(const std::set< int > &trackIDs, const std::vector< rb::WeightedHit > &whits, const std::vector< const rb::CellHit * > &allhits, const geo::View_t &view, std::map< int, double > *effMap=0, bool energyEff=false, double *desiredEnergy=0, double *totalEnergy=0, int *desiredHits=0, int *totalHits=0) const
Returns the fraction of all energy in an event from a specific set of Geant4 track IDs that are repre...
int Daughter(const int i) const
virtual void reconfigure(fhicl::ParameterSet const &p)
std::set< int > FindVisibleProngs(const art::PtrVector< rb::CellHit > &allhits, std::set< int > ids) const
const std::vector< Plot > plots
T get(std::string const &key) const
virtual void analyze(art::Event const &e)
bool fCheckTracks
should we check the reconstruction of tracks?
std::string fTrackModuleLabel
label for module making the tracks
double DetHalfHeight() const
EDAnalyzer(Table< Config > const &config)
double HitCollectionPurity(const std::set< int > &trackIDs, const std::vector< rb::WeightedHit > &whits, std::map< int, double > *purMap=0, std::map< int, int > *parents=0, bool energyPur=false) const
Returns the fraction of hits in a collection that come from the specified Geant4 track ids...
bool fCheckShowers
should we check the reconstruction of showers?
const sim::Particle * TrackIDToParticle(int const &id) const
Returns a pointer to the sim::Particle object corresponding to the given TrackID. ...
T * make(ARGS...args) const
code to link reconstructed objects back to the MC truth information
double DetHalfWidth() const
bool getByLabel(std::string const &label, std::string const &productInstanceName, Handle< PROD > &result) const
std::string fClusterAssnLabel
label for module clusters, showers, and tracks are associated with
MaybeLogger_< ELseverityLevel::ELsev_warning, false > LogWarning
const std::set< int > GetTrackIDList() const
Get all G4 track ids present in the event.
assert(nhit_max >=nhit_nbins)
std::string fShowerModuleLabel
label for module making the showers
double TotalLengthInDetector(const sim::Particle *part, geo::View_t view) const
Total length in either 2D or 3D.
void CheckRecoClusters(std::vector< art::Ptr< rb::Cluster > > const &clscol, art::PtrVector< rb::CellHit > const &allhits, Plots &plots)
Encapsulate the geometry of one entire detector (near, far, ndos)