Public Member Functions | Private Types | Private Member Functions | Static Private Member Functions | Private Attributes | List of all members
gar::CAFMaker Class Reference
Inheritance diagram for gar::CAFMaker:
art::EDAnalyzer art::detail::Analyzer art::detail::LegacyModule art::Observer art::ModuleBase

Public Member Functions

 CAFMaker (fhicl::ParameterSet const &p)
 
 CAFMaker (CAFMaker const &)=delete
 
 CAFMaker (CAFMaker &&)=delete
 
CAFMaker & operator= (CAFMaker const &)=delete
 
CAFMaker & operator= (CAFMaker &&)=delete
 
virtual void beginJob () override
 
void analyze (art::Event const &e) override
 
void endRun (const art::Run &run) override
 
- Public Member Functions inherited from art::EDAnalyzer
 EDAnalyzer (fhicl::ParameterSet const &pset)
 
template<typename Config >
 EDAnalyzer (Table< Config > const &config)
 
std::string workerType () const
 
- Public Member Functions inherited from art::detail::Analyzer
virtual ~Analyzer () noexcept
 
 Analyzer (fhicl::ParameterSet const &pset)
 
template<typename Config >
 Analyzer (Table< Config > const &config)
 
void doBeginJob (SharedResources const &resources)
 
void doEndJob ()
 
void doRespondToOpenInputFile (FileBlock const &fb)
 
void doRespondToCloseInputFile (FileBlock const &fb)
 
void doRespondToOpenOutputFiles (FileBlock const &fb)
 
void doRespondToCloseOutputFiles (FileBlock const &fb)
 
bool doBeginRun (RunPrincipal &rp, ModuleContext const &mc)
 
bool doEndRun (RunPrincipal &rp, ModuleContext const &mc)
 
bool doBeginSubRun (SubRunPrincipal &srp, ModuleContext const &mc)
 
bool doEndSubRun (SubRunPrincipal &srp, ModuleContext const &mc)
 
bool doEvent (EventPrincipal &ep, ModuleContext const &mc, std::atomic< std::size_t > &counts_run, std::atomic< std::size_t > &counts_passed, std::atomic< std::size_t > &counts_failed)
 
- Public Member Functions inherited from art::Observer
 ~Observer () noexcept
 
 Observer (Observer const &)=delete
 
 Observer (Observer &&)=delete
 
Observer & operator= (Observer const &)=delete
 
Observer & operator= (Observer &&)=delete
 
void registerProducts (ProductDescriptions &, ModuleDescription const &)
 
void fillDescriptions (ModuleDescription const &)
 
fhicl::ParameterSetID selectorConfig () const
 
- Public Member Functions inherited from art::ModuleBase
virtual ~ModuleBase () noexcept
 
 ModuleBase ()
 
ModuleDescription const & moduleDescription () const
 
void setModuleDescription (ModuleDescription const &)
 
std::array< std::vector< ProductInfo >, NumBranchTypes > const & getConsumables () const
 
void sortConsumables (std::string const &current_process_name)
 
template<typename T , BranchType BT>
ViewToken< T > consumesView (InputTag const &tag)
 
template<typename T , BranchType BT>
ViewToken< T > mayConsumeView (InputTag const &tag)
 

Private Types

typedef int TrkId
 

Private Member Functions

void FillGeneratorMonteCarloInfo (art::Event const &e, caf::StandardRecord &rec)
 
void FillRawInfo (art::Event const &e, caf::StandardRecord &rec)
 
void FillRecoInfo (art::Event const &e, caf::StandardRecord &rec)
 
void FillHighLevelRecoInfo (art::Event const &e, caf::StandardRecord &rec)
 
void processIonizationInfo (rec::TrackIoniz &ion, float ionizeTruncate, float &forwardIonVal, float &backwardIonVal)
 
float processOneDirection (std::vector< std::pair< float, float >> SigData, float ionizeTruncate)
 
float computeT (simb::MCTruth theMCTruth)
 
std::vector< std::pair< int, float > > processPIDInfo (float p)
 

Static Private Member Functions

static bool lessThan_byE (std::pair< float, float > a, std::pair< float, float > b)
 

Private Attributes

cheat::BackTrackerCore * BackTrack
 
float ItsInTulsa [3]
 
float xTPC
 
float rTPC
 
std::vector< std::string > fGeneratorLabels
 
std::vector< std::string > fGENIEGeneratorLabels
 
std::string fGeantLabel
 module label for geant4 simulated hits More...
 
std::string fInstanceLabelCalo
 Instance name for ECAL. More...
 
std::string fInstanceLabelMuID
 Instance name for MuID. More...
 
std::string fHitLabel
 module label for reco TPC hits rec::Hit More...
 
std::string fTPCClusterLabel
 module label for TPC Clusters rec::TPCCluster More...
 
std::string fTrackLabel
 module label for TPC Tracks rec:Track More...
 
std::string fTrackTrajectoryLabel
 module label for TPC Track Trajectories rec:TrackTrajectory More...
 
std::string fVertexLabel
 module label for vertexes rec:Vertex More...
 
std::string fVeeLabel
 module label for conversion/decay vertexes rec:Vee More...
 
std::string fRawCaloHitLabel
 module label for digitized calo hits raw::CaloRawDigit More...
 
std::string fRawMuIDHitLabel
 
std::string fCaloHitLabel
 module label for reco calo hits rec::CaloHit More...
 
std::string fMuIDHitLabel
 
std::string fClusterLabel
 module label for calo clusters rec::Cluster More...
 
std::string fClusterMuIDLabel
 module label for calo clusters rec::Cluster in MuID More...
 
std::string fPFLabel
 module label for reco particles rec::PFParticle More...
 
std::string fECALAssnLabel
 module label for track-clusters associations More...
 
std::string fPOTSummaryLabel
 
bool fWriteMCinfo
 Info from MCTruth, GTruth Default=true. More...
 
bool fWriteMCPTrajectory
 Write MCP Trajectory Default=true. More...
 
bool fWriteMCPTrajMomenta
 Write Momenta associated with MCP Trajectory Default=false. More...
 
bool fWriteMCCaloInfo
 Write MC info for calorimeter Default=true. More...
 
float fMatchMCPtoVertDist
 MCParticle to MC vertex match Default=roundoff. More...
 
bool fWriteHits
 Write info about TPC Hits Default=false. More...
 
bool fWriteTPCClusters
 Write TPCClusters info Default=true. More...
 
bool fWriteTracks
 Start/end X, P for tracks Default=true. More...
 
bool fWriteTrackTrajectories
 Point traj of reco tracks Default=false. More...
 
bool fWriteVertices
 Reco vertexes & their tracks Default=true. More...
 
bool fWriteVees
 Reco vees & their tracks Default=true. More...
 
bool fWriteMuID
 
bool fWriteCaloDigits
 Raw digits for calorimetry. Default=false. More...
 
bool fWriteCaloHits
 Write ECAL hits. Default=true. More...
 
bool fWriteCaloClusters
 Write ECAL clusters. Default=true. More...
 
bool fWriteMatchedTracks
 Write ECAL-track Assns Default=true. More...
 
float fIonizTruncate
 Default=1.00;. More...
 
bool fWriteCohInfo
 MC level t for coherent pi+. Default=false. More...
 
TTree * fTree
 
TH1 * fPOTHist
 
Float_t fTPC_X
 center of TPC stored as per-event & compressed by root More...
 
Float_t fTPC_Y
 
Float_t fTPC_Z
 
const geo::GeometryCore * fGeo
 pointer to the geometry More...
 
std::unordered_map< TrkId, Int_t > TrackIdToIndex
 
CLHEP::HepRandomEngine & fEngine
 random engine More...
 
std::unordered_map< int, TH2F * > m_pidinterp
 

Additional Inherited Members

- Public Types inherited from art::EDAnalyzer
using WorkerType = WorkerT< EDAnalyzer >
 
using ModuleType = EDAnalyzer
 
- Protected Member Functions inherited from art::Observer
std::string const & processName () const
 
bool wantAllEvents () const noexcept
 
bool wantEvent (ScheduleID id, Event const &e) const
 
Handle< TriggerResults > getTriggerResults (Event const &e) const
 
 Observer (fhicl::ParameterSet const &config)
 
 Observer (std::vector< std::string > const &select_paths, std::vector< std::string > const &reject_paths, fhicl::ParameterSet const &config)
 
- Protected Member Functions inherited from art::ModuleBase
ConsumesCollector & consumesCollector ()
 
template<typename T , BranchType = InEvent>
ProductToken< T > consumes (InputTag const &)
 
template<typename Element , BranchType = InEvent>
ViewToken< Element > consumesView (InputTag const &)
 
template<typename T , BranchType = InEvent>
void consumesMany ()
 
template<typename T , BranchType = InEvent>
ProductToken< T > mayConsume (InputTag const &)
 
template<typename Element , BranchType = InEvent>
ViewToken< Element > mayConsumeView (InputTag const &)
 
template<typename T , BranchType = InEvent>
void mayConsumeMany ()
 

Detailed Description

Definition at line 75 of file CAFMaker_module.cc.

Member Typedef Documentation

typedef int gar::CAFMaker::TrkId
private

Definition at line 184 of file CAFMaker_module.cc.

Constructor & Destructor Documentation

gar::CAFMaker::CAFMaker ( fhicl::ParameterSet const &  p)
explicit

Definition at line 199 of file CAFMaker_module.cc.

200  : EDAnalyzer(p),
201  fEngine(art::ServiceHandle<rndm::NuRandomService>()->createEngine(*this, p, "Seed")) {
202  fGeo = gar::providerFrom<geo::GeometryGAr>();
203 
204  bool usegenlabels =
205  p.get_if_present<std::vector<std::string> >("GeneratorLabels",fGeneratorLabels);
206  if (!usegenlabels) fGeneratorLabels.clear();
207 
208  bool usegeniegenlabels =
209  p.get_if_present<std::vector<std::string> >("GENIEGeneratorLabels",fGENIEGeneratorLabels);
210  if (!usegeniegenlabels) fGENIEGeneratorLabels.clear();
211 
212  //Sim Hits
213  fGeantLabel = p.get<std::string>("GEANTLabel","geant");
214  fInstanceLabelCalo = p.get<std::string>("InstanceLabelCalo","ECAL");
215  fInstanceLabelMuID = p.get<std::string>("InstanceLabelMuID","MuID");
216 
217  fHitLabel = p.get<std::string>("HitLabel","hit");
218  fTPCClusterLabel = p.get<std::string>("TPCClusterLabel","tpccluster");
219  fTrackLabel = p.get<std::string>("TrackLabel","track");
220  fTrackTrajectoryLabel = p.get<std::string>("TrackTrajectoryLabel","track");
221  fVertexLabel = p.get<std::string>("VertexLabel","vertex");
222  fVeeLabel = p.get<std::string>("VeeLabel","veefinder1");
223 
224  //Calorimetric related ECAL/MuID
225  fRawCaloHitLabel = p.get<std::string>("RawCaloHitLabel","daqsipm");
226  fRawMuIDHitLabel = p.get<std::string>("RawMuIDHitLabel","daqsipmmuid");
227  fCaloHitLabel = p.get<std::string>("CaloHitLabel","sipmhit");
228  fMuIDHitLabel = p.get<std::string>("MuIDHitLabel","sipmhit");
229 
230  fClusterLabel = p.get<std::string>("ClusterLabel","calocluster");
231  fClusterMuIDLabel = p.get<std::string>("MuIDClusterLabel","caloclustermuid");
232  fPFLabel = p.get<std::string>("PFLabel","pandora");
233  fECALAssnLabel = p.get<std::string>("ECALAssnLabel","trkecalassn");
234 
235  fPOTSummaryLabel = p.get<std::string>("POTSummaryLabel");
236 
237  // What to write
238  fWriteMCinfo = p.get<bool>("WriteMCinfo", true);
239  fWriteMCPTrajectory = p.get<bool>("WriteMCPTrajectory", true);
240  fWriteMCPTrajMomenta = p.get<bool>("WriteMCPTrajMomenta",false);
241  fWriteMCCaloInfo = p.get<bool>("WriteMCCaloInfo", true);
242  float MCPtoVertDefault = 10.0*std::numeric_limits<Float_t>::epsilon();
243  fMatchMCPtoVertDist = p.get<float>("MatchMCPtoVertDist",MCPtoVertDefault);
244 
245  fWriteHits = p.get<bool>("WriteHits", false);
246  fWriteTPCClusters = p.get<bool>("WriteTPCClusters", true);
247  fWriteTracks = p.get<bool>("WriteTracks", true);
248  fWriteTrackTrajectories = p.get<bool>("WriteTrackTrajectories", false);
249  fWriteVertices = p.get<bool>("WriteVertices", true);
250  fWriteVees = p.get<bool>("WriteVees", true);
251 
252  fWriteMuID = p.get<bool>("WriteMuID", true);
253  fWriteCaloDigits = p.get<bool>("WriteCaloDigits", false);
254  fWriteCaloHits = p.get<bool>("WriteCaloHits", true);
255  fWriteCaloClusters = p.get<bool>("WriteCaloClusters", true);
256  fWriteMatchedTracks = p.get<bool>("WriteMatchedTracks",true);
257 
258  fIonizTruncate = p.get<float>("IonizTruncate", 0.70);
259  fWriteCohInfo = p.get<bool> ("WriteCohInfo", false);
260 
261  if (usegenlabels) {
262  for (size_t i=0; i<fGeneratorLabels.size(); ++i) {
263  consumes<std::vector<simb::MCTruth> >(fGeneratorLabels.at(i));
264  }
265  } else {
266  consumesMany<std::vector<simb::MCTruth> >();
267  }
268 
269  if (usegeniegenlabels) {
270  for (size_t i=0; i<fGENIEGeneratorLabels.size(); ++i) {
271  consumes<std::vector<simb::GTruth> >(fGENIEGeneratorLabels.at(i));
272  }
273  } else {
274  consumesMany<std::vector<simb::GTruth> >();
275  }
276 
277  consumesMany<std::vector<sdp::GenieParticle> >();
278  //consumes<art::Assns<simb::MCTruth, simb::MCParticle> >(fGeantLabel);
279  consumes<std::vector<simb::MCParticle> >(fGeantLabel);
280 
281  //TPC related
282  consumes<std::vector<sdp::EnergyDeposit> >(fGeantLabel);
283  consumes<std::vector<rec::TPCCluster> >(fTPCClusterLabel);
284  consumes<art::Assns<rec::Track, rec::TPCCluster> >(fTPCClusterLabel);
285  consumes<std::vector<rec::Hit> >(fHitLabel);
286  consumes<std::vector<rec::Track> >(fTrackLabel);
287  consumes<std::vector<rec::TrackTrajectory> >(fTrackTrajectoryLabel);
288  consumes<std::vector<rec::Vertex> >(fVertexLabel);
289  consumes<art::Assns<rec::Track, rec::Vertex> >(fVertexLabel);
290  consumes<std::vector<rec::Vee> >(fVeeLabel);
291  consumes<art::Assns<rec::Track, rec::Vee> >(fVeeLabel);
292  consumes<std::vector<rec::TrackIoniz>>(fTrackLabel);
293  consumes<art::Assns<rec::TrackIoniz, rec::Track>>(fTrackLabel);
294 
295  //Calorimetry related
296  art::InputTag ecalgeanttag(fGeantLabel, fInstanceLabelCalo);
297  consumes<std::vector<gar::sdp::CaloDeposit> >(ecalgeanttag);
298 
299  art::InputTag ecalrawtag(fRawCaloHitLabel, fInstanceLabelCalo);
300  consumes<std::vector<raw::CaloRawDigit> >(ecalrawtag);
301 
302  art::InputTag ecalhittag(fCaloHitLabel, fInstanceLabelCalo);
303  consumes<std::vector<rec::CaloHit> >(ecalhittag);
304 
305  art::InputTag ecalclustertag(fClusterLabel, fInstanceLabelCalo);
306  consumes<std::vector<rec::Cluster> >(ecalclustertag);
307  consumes<art::Assns<rec::Cluster, rec::CaloHit>>(ecalclustertag);
308 
309  //Muon system related
310  if (fGeo->HasMuonDetector() && fWriteMuID) {
311  art::InputTag muidgeanttag(fGeantLabel, fInstanceLabelMuID);
312  consumes<std::vector<gar::sdp::CaloDeposit> >(muidgeanttag);
313 
314  art::InputTag muidrawtag(fRawMuIDHitLabel, fInstanceLabelMuID);
315  consumes<std::vector<raw::CaloRawDigit> >(muidrawtag);
316 
317  art::InputTag muidhittag(fCaloHitLabel, fInstanceLabelMuID);
318  consumes<std::vector<rec::CaloHit> >(muidhittag);
319 
320  art::InputTag muidclustertag(fClusterMuIDLabel, fInstanceLabelMuID);
321  consumes<std::vector<rec::Cluster> >(muidclustertag);
322  consumes<art::Assns<rec::Cluster, rec::CaloHit>>(muidclustertag);
323  }
324 
325  consumes<art::Assns<rec::Cluster, rec::Track>>(fECALAssnLabel);
326 
327  return;
328 } // end constructor
art::detail::EngineCreator::createEngine
base_engine_t & createEngine(seed_t seed)
Definition: EngineCreator.cc:20
gar::CAFMaker::fClusterLabel
std::string fClusterLabel
module label for calo clusters rec::Cluster
Definition: CAFMaker_module.cc:141
gar::CAFMaker::fWriteCaloClusters
bool fWriteCaloClusters
Write ECAL clusters. Default=true.
Definition: CAFMaker_module.cc:165
gar::CAFMaker::fRawCaloHitLabel
std::string fRawCaloHitLabel
module label for digitized calo hits raw::CaloRawDigit
Definition: CAFMaker_module.cc:135
gar::CAFMaker::fInstanceLabelCalo
std::string fInstanceLabelCalo
Instance name for ECAL.
Definition: CAFMaker_module.cc:125
gar::CAFMaker::fWriteMCinfo
bool fWriteMCinfo
Info from MCTruth, GTruth Default=true.
Definition: CAFMaker_module.cc:149
gar::CAFMaker::fTPCClusterLabel
std::string fTPCClusterLabel
module label for TPC Clusters rec::TPCCluster
Definition: CAFMaker_module.cc:129
string
std::string string
Definition: nybbler.cc:12
gar::CAFMaker::fWriteCaloDigits
bool fWriteCaloDigits
Raw digits for calorimetry. Default=false.
Definition: CAFMaker_module.cc:163
gar::CAFMaker::fWriteTPCClusters
bool fWriteTPCClusters
Write TPCClusters info Default=true.
Definition: CAFMaker_module.cc:156
gar::CAFMaker::fWriteMCCaloInfo
bool fWriteMCCaloInfo
Write MC info for calorimeter Default=true.
Definition: CAFMaker_module.cc:152
gar::CAFMaker::fWriteVertices
bool fWriteVertices
Reco vertexes & their tracks Default=true.
Definition: CAFMaker_module.cc:159
gar::CAFMaker::fWriteCaloHits
bool fWriteCaloHits
Write ECAL hits. Default=true.
Definition: CAFMaker_module.cc:164
art::EDAnalyzer::EDAnalyzer
EDAnalyzer(fhicl::ParameterSet const &pset)
Definition: EDAnalyzer.h:25
gar::CAFMaker::fWriteMCPTrajectory
bool fWriteMCPTrajectory
Write MCP Trajectory Default=true.
Definition: CAFMaker_module.cc:150
gar::CAFMaker::fRawMuIDHitLabel
std::string fRawMuIDHitLabel
Definition: CAFMaker_module.cc:136
gar::CAFMaker::fMuIDHitLabel
std::string fMuIDHitLabel
Definition: CAFMaker_module.cc:139
gar::CAFMaker::fWriteTracks
bool fWriteTracks
Start/end X, P for tracks Default=true.
Definition: CAFMaker_module.cc:157
gar::CAFMaker::fGeneratorLabels
std::vector< std::string > fGeneratorLabels
Definition: CAFMaker_module.cc:122
gar::CAFMaker::fTrackLabel
std::string fTrackLabel
module label for TPC Tracks rec:Track
Definition: CAFMaker_module.cc:130
gar::CAFMaker::fGeo
const geo::GeometryCore * fGeo
pointer to the geometry
Definition: CAFMaker_module.cc:182
gar::CAFMaker::fWriteVees
bool fWriteVees
Reco vees & their tracks Default=true.
Definition: CAFMaker_module.cc:160
gar::CAFMaker::fGeantLabel
std::string fGeantLabel
module label for geant4 simulated hits
Definition: CAFMaker_module.cc:124
gar::CAFMaker::fWriteMatchedTracks
bool fWriteMatchedTracks
Write ECAL-track Assns Default=true.
Definition: CAFMaker_module.cc:166
gar::CAFMaker::fWriteCohInfo
bool fWriteCohInfo
MC level t for coherent pi+. Default=false.
Definition: CAFMaker_module.cc:170
gar::CAFMaker::fPOTSummaryLabel
std::string fPOTSummaryLabel
Definition: CAFMaker_module.cc:146
test.p
p
Definition: test.py:223
gar::CAFMaker::fClusterMuIDLabel
std::string fClusterMuIDLabel
module label for calo clusters rec::Cluster in MuID
Definition: CAFMaker_module.cc:142
gar::CAFMaker::fPFLabel
std::string fPFLabel
module label for reco particles rec::PFParticle
Definition: CAFMaker_module.cc:143
gar::CAFMaker::fMatchMCPtoVertDist
float fMatchMCPtoVertDist
MCParticle to MC vertex match Default=roundoff.
Definition: CAFMaker_module.cc:153
gar::CAFMaker::fVertexLabel
std::string fVertexLabel
module label for vertexes rec:Vertex
Definition: CAFMaker_module.cc:132
gar::CAFMaker::fCaloHitLabel
std::string fCaloHitLabel
module label for reco calo hits rec::CaloHit
Definition: CAFMaker_module.cc:138
gar::CAFMaker::fTrackTrajectoryLabel
std::string fTrackTrajectoryLabel
module label for TPC Track Trajectories rec:TrackTrajectory
Definition: CAFMaker_module.cc:131
gar::CAFMaker::fGENIEGeneratorLabels
std::vector< std::string > fGENIEGeneratorLabels
Definition: CAFMaker_module.cc:123
gar::CAFMaker::fWriteTrackTrajectories
bool fWriteTrackTrajectories
Point traj of reco tracks Default=false.
Definition: CAFMaker_module.cc:158
gar::CAFMaker::fIonizTruncate
float fIonizTruncate
Default=1.00;.
Definition: CAFMaker_module.cc:169
gar::CAFMaker::fHitLabel
std::string fHitLabel
module label for reco TPC hits rec::Hit
Definition: CAFMaker_module.cc:128
gar::geo::GeometryCore::HasMuonDetector
bool HasMuonDetector() const
Definition: GeometryCore.h:999
gar::CAFMaker::fVeeLabel
std::string fVeeLabel
module label for conversion/decay vertexes rec:Vee
Definition: CAFMaker_module.cc:133
gar::CAFMaker::fWriteHits
bool fWriteHits
Write info about TPC Hits Default=false.
Definition: CAFMaker_module.cc:155
gar::CAFMaker::fWriteMCPTrajMomenta
bool fWriteMCPTrajMomenta
Write Momenta associated with MCP Trajectory Default=false.
Definition: CAFMaker_module.cc:151
gar::CAFMaker::fEngine
CLHEP::HepRandomEngine & fEngine
random engine
Definition: CAFMaker_module.cc:188
gar::CAFMaker::fInstanceLabelMuID
std::string fInstanceLabelMuID
Instance name for MuID.
Definition: CAFMaker_module.cc:126
gar::CAFMaker::fECALAssnLabel
std::string fECALAssnLabel
module label for track-clusters associations
Definition: CAFMaker_module.cc:144
gar::CAFMaker::CAFMaker ( CAFMaker const &  )
delete
gar::CAFMaker::CAFMaker ( CAFMaker &&  )
delete

Member Function Documentation

void gar::CAFMaker::analyze ( art::Event const &  e)
overridevirtual

Implements art::EDAnalyzer.

Definition at line 419 of file CAFMaker_module.cc.

419  {
420 
421  caf::StandardRecord rec;
422  caf::StandardRecord* prec = &rec;
423  fTree->SetBranchAddress("rec", &prec);
424 
425  rec.hdr.run = e.run();
426  rec.hdr.subrun = e.subRun();
427  rec.hdr.event = e.id().event();
428 
429  rec.hdr.TPC_X = fTPC_X;
430  rec.hdr.TPC_Y = fTPC_Y;
431  rec.hdr.TPC_Z = fTPC_Z;
432 
433 
434  // Need a non-constant backtracker instance, for now, in anayze ot beginJob
435  cheat::BackTrackerCore const* const_bt = gar::providerFrom<cheat::BackTracker>();
436  BackTrack = const_cast<cheat::BackTrackerCore*>(const_bt);
437 
438  //Fill generator and MC Information
439  if (fWriteMCinfo) FillGeneratorMonteCarloInfo(e, rec);
440 
441  //Fill Raw Information
442  if (fWriteCaloDigits) FillRawInfo(e, rec);
443 
444  //Fill Reco Information
445  FillRecoInfo(e, rec);
446 
447  //Fill HL Reco Information
448  FillHighLevelRecoInfo(e, rec);
449 
450  fTree->Fill();
451  return;
452 }
caf::SRHeader::event
int event
number of the event being processed
Definition: SRHeader.h:12
tracks.rec
rec
Definition: tracks.py:88
caf::StandardRecord::hdr
SRHeader hdr
global event info
Definition: StandardRecord.h:24
gar::CAFMaker::fWriteMCinfo
bool fWriteMCinfo
Info from MCTruth, GTruth Default=true.
Definition: CAFMaker_module.cc:149
caf::SRHeader::TPC_X
float TPC_X
center of TPC stored as per-event & compressed by root
Definition: SRHeader.h:16
gar::CAFMaker::fWriteCaloDigits
bool fWriteCaloDigits
Raw digits for calorimetry. Default=false.
Definition: CAFMaker_module.cc:163
gar::CAFMaker::FillRawInfo
void FillRawInfo(art::Event const &e, caf::StandardRecord &rec)
Definition: CAFMaker_module.cc:802
caf::SRHeader::run
int run
number of the run being processed
Definition: SRHeader.h:13
gar::CAFMaker::FillHighLevelRecoInfo
void FillHighLevelRecoInfo(art::Event const &e, caf::StandardRecord &rec)
Definition: CAFMaker_module.cc:944
e
const double e
Definition: gUpMuFluxGen.cxx:165
gar::CAFMaker::BackTrack
cheat::BackTrackerCore * BackTrack
Definition: CAFMaker_module.cc:95
gar::CAFMaker::FillRecoInfo
void FillRecoInfo(art::Event const &e, caf::StandardRecord &rec)
Definition: CAFMaker_module.cc:856
gar::CAFMaker::FillGeneratorMonteCarloInfo
void FillGeneratorMonteCarloInfo(art::Event const &e, caf::StandardRecord &rec)
Definition: CAFMaker_module.cc:467
caf::SRHeader::TPC_Y
float TPC_Y
Definition: SRHeader.h:17
caf::SRHeader::TPC_Z
float TPC_Z
Definition: SRHeader.h:18
gar::CAFMaker::fTPC_X
Float_t fTPC_X
center of TPC stored as per-event & compressed by root
Definition: CAFMaker_module.cc:176
caf::StandardRecord
The StandardRecord is the primary top-level object in the Common Analysis File trees.
Definition: StandardRecord.h:24
caf::SRHeader::subrun
int subrun
number of the sub-run being processed
Definition: SRHeader.h:14
void gar::CAFMaker::beginJob ( )
overridevirtual

Reimplemented from art::EDAnalyzer.

Definition at line 335 of file CAFMaker_module.cc.

335  {
336 
337  fTPC_X = ItsInTulsa[0] = fGeo->TPCXCent();
338  fTPC_Y = ItsInTulsa[1] = fGeo->TPCYCent();
339  fTPC_Z = ItsInTulsa[2] = fGeo->TPCZCent();
340 
341  xTPC = fGeo->TPCLength() / 2.;
342  rTPC = fGeo->TPCRadius();
343 
344  art::ServiceHandle<art::TFileService> tfs;
345  fTree = tfs->make<TTree>("recTree", "recTree");
346 
347  // Create the branch. We will update the address before we write the tree
348  caf::StandardRecord* rec = 0;
349  fTree->Branch("rec", &rec);
350 
351  fPOTHist = fPOTSummaryLabel.empty() ? 0 : tfs->make<TH1F>("TotalPOT", "", 1, 0, 1);
352 
353  if (fWriteMCPTrajectory) {
354  // Checking for parameter file validity only once to simplify
355  // later code in CAFMaker::process etc.
356  if (!fWriteMCinfo) {
357  throw cet::exception("CAFMaker")
358  << " fWriteMCPTrajectory, but !fWriteMCinfo."
359  << " Line " << __LINE__ << " in file " << __FILE__ << std::endl;
360  }
361  }
362 
363  if (fWriteMCCaloInfo) {
364  if (!fWriteMCinfo) {
365  throw cet::exception("CAFMaker")
366  << " fWriteMCCaloInfo, but !fWriteMCinfo."
367  << " Line " << __LINE__ << " in file " << __FILE__ << std::endl;
368  }
369  }
370 
371  // Reco'd verts & their track-ends
372  if (fWriteVertices) {
373  if (!fWriteTracks) {
374  throw cet::exception("CAFMaker")
375  << " fWriteVertices, but !fWriteTracks."
376  << " Line " << __LINE__ << " in file " << __FILE__ << std::endl;
377  }
378  }
379 
380  // Reco'd vees & their track-ends
381  if (fWriteVees) {
382  if (!fWriteTracks) {
383  throw cet::exception("CAFMaker")
384  << " fWriteVees, but !fWriteTracks."
385  << " Line " << __LINE__ << " in file " << __FILE__ << std::endl;
386  }
387  }
388 
389  if (fWriteMatchedTracks) {
390  if (!fWriteTracks || !fWriteCaloClusters) {
391  throw cet::exception("CAFMaker")
392  << " fWriteMatchedTracks, but (!fWriteTracks || !fWriteCaloClusters)."
393  << " Line " << __LINE__ << " in file " << __FILE__ << std::endl;
394  }
395  }
396 
397  std::string filename = "${DUNE_PARDATA_DIR}/MPD/dedxPID/dedxpidmatrices8kevcm.root";
398  TFile infile(filename.c_str(), "READ");
399 
400  m_pidinterp.clear();
401  char str[11];
402  for (int q = 0; q < 501; ++q) {
403  sprintf(str, "%d", q);
404  std::string s = "pidmatrix";
405  s.append(str);
406  // read the 500 histograms one by one; each histogram is a
407  // 6 by 6 matrix of probabilities for a given momentum value
408  m_pidinterp.insert( std::make_pair(q, (TH2F*) infile.Get(s.c_str())->Clone("pidinterp")) );
409  }
410 
411  return;
412 } // End of :CAFMaker::beginJob
gar::CAFMaker::fWriteCaloClusters
bool fWriteCaloClusters
Write ECAL clusters. Default=true.
Definition: CAFMaker_module.cc:165
tracks.rec
rec
Definition: tracks.py:88
gar::CAFMaker::fWriteMCinfo
bool fWriteMCinfo
Info from MCTruth, GTruth Default=true.
Definition: CAFMaker_module.cc:149
string
std::string string
Definition: nybbler.cc:12
gar::CAFMaker::fWriteMCCaloInfo
bool fWriteMCCaloInfo
Write MC info for calorimeter Default=true.
Definition: CAFMaker_module.cc:152
gar::CAFMaker::fWriteVertices
bool fWriteVertices
Reco vertexes & their tracks Default=true.
Definition: CAFMaker_module.cc:159
gar::geo::GeometryCore::TPCXCent
float TPCXCent() const
Returns the X location of the center of the TPC in cm.
Definition: GeometryCore.h:778
gar::CAFMaker::fWriteMCPTrajectory
bool fWriteMCPTrajectory
Write MCP Trajectory Default=true.
Definition: CAFMaker_module.cc:150
train.filename
string filename
Definition: train.py:213
gar::CAFMaker::fWriteTracks
bool fWriteTracks
Start/end X, P for tracks Default=true.
Definition: CAFMaker_module.cc:157
gar::CAFMaker::fGeo
const geo::GeometryCore * fGeo
pointer to the geometry
Definition: CAFMaker_module.cc:182
infile
string infile
Definition: INukeNucleonCorr.cxx:62
gar::CAFMaker::fWriteVees
bool fWriteVees
Reco vees & their tracks Default=true.
Definition: CAFMaker_module.cc:160
gar::geo::GeometryCore::TPCRadius
float TPCRadius() const
Returns the radius of the TPC (y or z direction)
Definition: GeometryCore.h:755
gar::CAFMaker::fWriteMatchedTracks
bool fWriteMatchedTracks
Write ECAL-track Assns Default=true.
Definition: CAFMaker_module.cc:166
gar::CAFMaker::fPOTSummaryLabel
std::string fPOTSummaryLabel
Definition: CAFMaker_module.cc:146
gar::geo::GeometryCore::TPCZCent
float TPCZCent() const
Returns the Z location of the center of the TPC in cm.
Definition: GeometryCore.h:792
gar::CAFMaker::fTPC_X
Float_t fTPC_X
center of TPC stored as per-event & compressed by root
Definition: CAFMaker_module.cc:176
gar::CAFMaker::m_pidinterp
std::unordered_map< int, TH2F * > m_pidinterp
Definition: CAFMaker_module.cc:189
caf::StandardRecord
The StandardRecord is the primary top-level object in the Common Analysis File trees.
Definition: StandardRecord.h:24
gar::geo::GeometryCore::TPCYCent
float TPCYCent() const
Returns the Y location of the center of the TPC in cm.
Definition: GeometryCore.h:785
gar::geo::GeometryCore::TPCLength
float TPCLength() const
Returns the length of the TPC (x direction)
Definition: GeometryCore.h:763
s
static QCString * s
Definition: config.cpp:1042
str
static QCString str
Definition: fortrancode.cpp:27098
fhicl::exception
cet::coded_exception< error, detail::translate > exception
Definition: exception.h:33
endl
QTextStream & endl(QTextStream &s)
Definition: qtextstream.cpp:2030
QCString::append
QCString & append(const char *s)
Definition: qcstring.cpp:383
float gar::CAFMaker::computeT ( simb::MCTruth  theMCTruth)
private

Definition at line 1482 of file CAFMaker_module.cc.

1482  {
1483  // Warning. You probably want the absolute value of t, not t.
1484  int nPart = theMCTruth.NParticles();
1485  enum { nu, mu, pi};
1486  float E[3], Px[3], Py[3], Pz[3];
1487  E[nu] = E[mu] = E[pi] = -1e42;
1488 
1489  for (int i=0; i<3;++i) {
1490  Px[i] = 0;
1491  Py[i] = 0;
1492  Pz[i] = 0;
1493  E[i] = 0;
1494  }
1495  // Find t from the MCParticles via the
1496  for (int iPart=0; iPart<nPart; iPart++) {
1497  simb::MCParticle Part = theMCTruth.GetParticle(iPart);
1498  int code = Part.PdgCode();
1499  int mom = Part.Mother();
1500 
1501  // get the neutrino
1502  if ( abs(code) == 12 || abs(code) == 14 || abs(code) == 16 ) {
1503  if (mom == -1) {
1504  E[nu] = Part.E(); Px[nu] = Part.Px(); Py[nu] = Part.Py(); Pz[nu] = Part.Pz();
1505  }
1506  }
1507 
1508  // get the lepton
1509  if ( abs(code) == 11 || abs(code) == 13 || abs(code) == 15 ) {
1510  if (mom == 0) {
1511  E[mu] = Part.E(); Px[mu] = Part.Px(); Py[mu] = Part.Py(); Pz[mu] = Part.Pz();
1512  }
1513  }
1514 
1515  // get the pion
1516  if ( code==111 || abs(code)==211 ) {
1517  if (mom == 1) {
1518  E[pi] = Part.E(); Px[pi] = Part.Px(); Py[pi] = Part.Py(); Pz[pi] = Part.Pz();
1519  }
1520  }
1521 
1522  // get outa here
1523  if ( E[nu]!=0 && E[mu]!=0 && E[pi]!=0) break;
1524 
1525  }
1526 
1527  // Compute t; reuse nu 4-vector to get first q, then t.
1528  E[nu] -= E[mu]; Px[nu] -= Px[mu]; Py[nu] -= Py[mu]; Pz[nu] -= Pz[mu];
1529  E[nu] -= E[pi]; Px[nu] -= Px[pi]; Py[nu] -= Py[pi]; Pz[nu] -= Pz[pi];
1530  float t = E[nu]*E[nu] -Px[nu]*Px[nu] -Py[nu]*Py[nu] -Pz[nu]*Pz[nu];
1531  return t;
1532 }
simb::MCParticle::E
double E(const int i=0) const
Definition: MCParticle.h:233
simb::MCParticle::PdgCode
int PdgCode() const
Definition: MCParticle.h:212
simb::MCParticle::Py
double Py(const int i=0) const
Definition: MCParticle.h:231
simb::MCParticle::Mother
int Mother() const
Definition: MCParticle.h:213
simb::MCParticle::Px
double Px(const int i=0) const
Definition: MCParticle.h:230
simb::MCTruth::NParticles
int NParticles() const
Definition: MCTruth.h:75
abs
T abs(T value)
Definition: sparse_vector_test.cc:30
reco_momentum_tuples.t
t
Definition: reco_momentum_tuples.py:25
code
CodeOutputInterface * code
Definition: tclscanner.cpp:1114
simb::MCTruth::GetParticle
const simb::MCParticle & GetParticle(int i) const
Definition: MCTruth.h:76
simb::MCParticle::Pz
double Pz(const int i=0) const
Definition: MCParticle.h:232
wirecell.units.pi
float pi
Definition: units.py:11
void gar::CAFMaker::endRun ( const art::Run &  run)
override

Definition at line 455 of file CAFMaker_module.cc.

456 {
457  if(fPOTHist){
458  art::InputTag itag(fPOTSummaryLabel);
459  auto potsum = run.getHandle<gar::sumdata::POTSummary>(itag);
460  fPOTHist->Fill(.5, potsum->totgoodpot);
461  }
462 }
art::DataViewImpl::getHandle
Handle< PROD > getHandle(SelectorBase const &) const
Definition: DataViewImpl.h:382
gar::CAFMaker::fPOTSummaryLabel
std::string fPOTSummaryLabel
Definition: CAFMaker_module.cc:146
void gar::CAFMaker::FillGeneratorMonteCarloInfo ( art::Event const &  e,
caf::StandardRecord &  rec 
)
private

Definition at line 467 of file CAFMaker_module.cc.

467  {
468 
469  // ============= Get art handles ==========================================
470  // Get handles for MCinfo, also good for MCPTrajectory
471  std::vector< art::Handle< std::vector<simb::MCTruth> > > mcthandlelist;
472  std::vector< art::Handle< std::vector<simb::GTruth> > > gthandlelist;
473 
474  if (fGeneratorLabels.size()<1) {
475  mcthandlelist = e.getMany<std::vector<simb::MCTruth> >(); // get them all (even if there are none)
476  } else {
477  mcthandlelist.resize(fGeneratorLabels.size());
478  for (size_t i=0; i< fGeneratorLabels.size(); ++i) {
479  // complain if we wanted a specific one but didn't find it
480  art::InputTag itag(fGeneratorLabels.at(i));
481  mcthandlelist.at(i) = e.getHandle< std::vector<simb::MCTruth> >(itag);
482  if (!mcthandlelist.at(i)) {
483  throw cet::exception("CAFMaker") << " No simb::MCTruth branch."
484  << " Line " << __LINE__ << " in file " << __FILE__ << std::endl;
485  }
486  }
487  }
488 
489  if (fGENIEGeneratorLabels.size()<1) {
490  gthandlelist = e.getMany<std::vector<simb::GTruth> >(); // get them all (even if there are none)
491  } else {
492  gthandlelist.resize(fGENIEGeneratorLabels.size());
493  for (size_t i=0; i< fGENIEGeneratorLabels.size(); ++i) {
494  gthandlelist.at(i) = e.getHandle<std::vector<simb::GTruth> >(fGENIEGeneratorLabels.at(i));
495  // complain if we wanted a specific one but didn't find it
496  if (!gthandlelist.at(i)) {
497  throw cet::exception("CAFMaker") << " No simb::GTruth branch."
498  << " Line " << __LINE__ << " in file " << __FILE__ << std::endl;
499  }
500  }
501  }
502 
503  // save MCTruth info
504  for (size_t imchl = 0; imchl < mcthandlelist.size(); ++imchl) {
505  for ( auto const& mct : (*mcthandlelist.at(imchl)) ) {
506  if (mct.NeutrinoSet()) {
507  simb::MCNeutrino nuw = mct.GetNeutrino();
508  caf::SRNeutrino srnu;
509  srnu.pdg = nuw.Nu().PdgCode();
510  srnu.ccnc = nuw.CCNC();
511  srnu.mode = nuw.Mode();
512  srnu.interactionType = nuw.InteractionType();
513  srnu.Q2 = nuw.QSqr();
514  srnu.W = nuw.W();
515  srnu.X = nuw.X();
516  srnu.Y = nuw.Y();
517  srnu.theta = nuw.Theta();
518  if (fWriteCohInfo) {
519  srnu.T = computeT(mct);
520  }
521  srnu.vtx = caf::SRVector3D(nuw.Nu().EndX(), nuw.Nu().EndY(), nuw.Nu().EndZ());
522  srnu.p = caf::SRVector3D(nuw.Nu().Px(), nuw.Nu().Py(), nuw.Nu().Pz());
523 
524  rec.mc.nu.push_back(srnu);
525  } // end MC info from MCTruth
526  }
527  }
528 
529  rec.mc.nnu = rec.mc.nu.size();
530 
531  unsigned int nuidx = 0;
532  // save GTruth info
533  for (size_t igthl = 0; igthl < gthandlelist.size(); ++igthl) {
534  for ( auto const& gt : (*gthandlelist.at(igthl)) ) {
535  if(nuidx >= rec.mc.nu.size()) abort();
536 
537  // We're assuming we visit the gtruths in the same order (should
538  // use an Assn?)
539  rec.mc.nu[nuidx].gint = gt.fGint;
540  rec.mc.nu[nuidx].tgtPDG = gt.ftgtPDG;
541  rec.mc.nu[nuidx].weight = gt.fweight;
542  rec.mc.nu[nuidx].gT = gt.fgT;
543 
544  ++nuidx;
545  }
546  }
547 
548  // Save the particle list from the GENIE event record
549 
550  auto gparthandlelist = e.getMany<std::vector<sdp::GenieParticle>>();
551 
552  for (size_t igphl = 0; igphl < gparthandlelist.size(); ++igphl) {
553  for ( auto const& gpart : (*gparthandlelist.at(igphl)) ) {
554  caf::SRGenieParticle srgpart;
555  srgpart.intidx = gpart.InteractionIndex();
556  srgpart.idx = gpart.Index();
557  srgpart.pdg = gpart.Pdg();
558  srgpart.status = gpart.Status();
559  srgpart.name = gpart.Name();
560  srgpart.firstmom = gpart.FirstMother();
561  srgpart.lastmom = gpart.LastMother();
562  srgpart.firstdaugh = gpart.FirstDaughter();
563  srgpart.lastdaugh = gpart.LastDaughter();
564  srgpart.p = caf::SRVector3D(gpart.Px(), gpart.Py(), gpart.Pz());
565  srgpart.E = gpart.E();
566  srgpart.mass = gpart.Mass();
567 
568  rec.mc.gpart.push_back(srgpart);
569  }
570  }
571 
572  rec.mc.ngpart = rec.mc.gpart.size();
573 
574  auto MCPHandle = e.getHandle< std::vector<simb::MCParticle> >(fGeantLabel);
575  if (!MCPHandle) {
576  throw cet::exception("CAFMaker") << " No simb::MCParticle branch."
577  << " Line " << __LINE__ << " in file " << __FILE__ << std::endl;
578  }
579 
580  // Save MCParticle info (post-GEANT; for pre-GEANT, maybe try MCTruth?)
581  // Helps to have a map from the TrackId from generator to index number
582  // into *MCPHandle, which is a vector<MCParticle>. Not all TrackId value
583  // appear in MCPHandle, although I think the missing ones are all from
584  // the GEANT rather than the GENIE stage of GENIEGen. (*MCPHandle).size()
585  // is the same as eg fMCPStartX.size() by the time this entry is written.
586  // TracIdToIndex isa class variable for accessibility by per-Track and
587  // per-Cluster code
588  Int_t index = 0;
589  for ( auto const& mcp : (*MCPHandle) ) {
590  int TrackId = mcp.TrackId();
591  TrackIdToIndex[TrackId] = index++;
592  }
593 
594  for ( auto const& mcp : (*MCPHandle) ) {
595  caf::SRParticle srpart;
596  srpart.trkid = mcp.TrackId();
597  srpart.pdg = mcp.PdgCode();
598 
599  // If mcp.Mother() == 0, particle is from initial vertex;
600  // Set MCMotherIndex to -1 in that case.
601  TrkId momTrkId = mcp.Mother();
602  Int_t momIndex = -1;
603  int momPDG = 0;
604  if (momTrkId>0) {
605  //Check if it exists!
606  if(TrackIdToIndex.find(momTrkId) != TrackIdToIndex.end()){
607  momIndex = TrackIdToIndex[momTrkId];
608  momPDG = (*MCPHandle).at(momIndex).PdgCode();
609  } else {
610  MF_LOG_DEBUG("CAFMaker_module")
611  << " mcp trkid " << mcp.TrackId()
612  << " pdg code " << mcp.PdgCode()
613  << " could not find mother trk id " << momTrkId
614  << " in the TrackIdToIndex map"
615  << " creating process is [ " << mcp.Process() << " ]";
616  }
617  }
618 
619  srpart.momidx = momIndex;
620  srpart.momtrkid = mcp.Mother(); //directly trackid not index
621  srpart.mompdg = momPDG;
622 
623  const TLorentzVector& position = mcp.Position(0);
624  const TLorentzVector& momentum = mcp.Momentum(0);
625  srpart.start.x = position.X();
626  srpart.start.y = position.Y();
627  srpart.start.z = position.Z();
628  srpart.time = mcp.T();
629  srpart.startp.x = momentum.Px();
630  srpart.startp.y = momentum.Py();
631  srpart.startp.z = momentum.Pz();
632 
633  const TLorentzVector& positionEnd = mcp.EndPosition();
634  const TLorentzVector& momentumEnd = mcp.EndMomentum();
635  srpart.end.x = positionEnd.X();
636  srpart.end.y = positionEnd.Y();
637  srpart.end.z = positionEnd.Z();
638  srpart.endp.x = momentumEnd.Px();
639  srpart.endp.y = momentumEnd.Py();
640  srpart.endp.z = momentumEnd.Pz();
641  srpart.proc = mcp.Process();
642  srpart.endproc = mcp.EndProcess();
643 
644  rec.mc.part.push_back(srpart);
645  }
646 
647  rec.mc.npart = rec.mc.part.size();
648 
649  // Get vertex for each MCParticle. In principle, the mct.GetParticle()
650  // method should produce all of them; filter out the ones with a non-stable
651  // status code and the GENIE pseudo-particles (Pid>= 2000000000) and you
652  // should have the contents of (*MCPHandle). But you don't. On the 1st
653  // event I (Leo Bellantoni) tested, with 40 vertices and 351 primary
654  // particles, particle 143, a proton, evidently decayed or something and is
655  // not in (*MCPHandle). So we use the following primitive earth technology.
656  size_t nMCParticles = (*MCPHandle).size();
657  size_t iMCParticle = 0;
658  for (; iMCParticle<nMCParticles; ++iMCParticle) {
659  // Assign noprimary to start with
660  rec.mc.part[iMCParticle].vtxidx = -1;
661  // Do the primaries first
662  if (rec.mc.part[iMCParticle].momidx!=-1) break;
663  Float_t trackX = rec.mc.part[iMCParticle].start.x;
664  Float_t trackY = rec.mc.part[iMCParticle].start.y;
665  Float_t trackZ = rec.mc.part[iMCParticle].start.z;
666  int vertexIndex = 0;
667  for (size_t imchl = 0; imchl < mcthandlelist.size(); ++imchl) {
668  for ( auto const& mct : (*mcthandlelist.at(imchl)) ) {
669  if (mct.NeutrinoSet()) {
670  simb::MCNeutrino nuw = mct.GetNeutrino();
671  Float_t vertX = nuw.Nu().EndX();
672  Float_t vertY = nuw.Nu().EndY();
673  Float_t vertZ = nuw.Nu().EndZ();
674  Float_t dist = std::hypot(trackX-vertX,trackY-vertY,trackZ-vertZ);
675  if ( dist <= fMatchMCPtoVertDist ) {
676  rec.mc.part[iMCParticle].vtxidx = vertexIndex;
677  goto foundMCvert; // break out of all inner loops
678  }
679  }
680  ++vertexIndex;
681  }
682  }
683  foundMCvert:
684  ;
685  }
686 
687  // Now the secondaries. As they are after the primaries, do not re-init iMCParticle
688  for (; iMCParticle<nMCParticles; ++iMCParticle) {
689  int momIndex = rec.mc.part[iMCParticle].momidx;
690  int lastMCParticle = iMCParticle;
691  while (momIndex != -1) {
692  lastMCParticle = momIndex;
693  momIndex = rec.mc.part[momIndex].momidx;
694  }
695  rec.mc.part[iMCParticle].vtxidx = rec.mc.part[lastMCParticle].vtxidx;
696  }
697 
698  if (fWriteMCPTrajectory) {
699  // It's in the MCParticle table
700  Int_t mcpIndex = 0;
701  for ( auto const& mcp : (*MCPHandle) ) {
702  const TDatabasePDG* databasePDG = TDatabasePDG::Instance();
703  const TParticlePDG* definition = databasePDG->GetParticle( mcp.PdgCode() );
704  //No charge don't store the trajectory
705  if (definition==nullptr || definition->Charge() == 0) continue;
706  //TrackID of the mcp to keep track to which mcp this trajectory is
707  int trackId = mcp.TrackId();
708  for(uint iTraj=0; iTraj < mcp.Trajectory().size(); iTraj++) {
709  float xTraj = mcp.Trajectory().X(iTraj);
710  float yTraj = mcp.Trajectory().Y(iTraj);
711  float zTraj = mcp.Trajectory().Z(iTraj);
712  float rTraj = std::hypot( yTraj - ItsInTulsa[1], zTraj - ItsInTulsa[2]);
713 
714  if (abs(xTraj - ItsInTulsa[0]) > xTPC) continue;
715  if (rTraj > rTPC) continue;
716 
717  caf::SRParticleTrajectoryPoint srtrajpt;
718  srtrajpt.x = xTraj;
719  srtrajpt.y = yTraj;
720  srtrajpt.z = zTraj;
721  srtrajpt.t = mcp.Trajectory().T(iTraj);
722  if (fWriteMCPTrajMomenta) {
723  srtrajpt.px = mcp.Trajectory().Px(iTraj);
724  srtrajpt.py = mcp.Trajectory().Py(iTraj);
725  srtrajpt.pz = mcp.Trajectory().Pz(iTraj);
726  }
727  srtrajpt.E = mcp.Trajectory().E(iTraj);
728 
729  rec.mc.part[mcpIndex].traj.push_back(srtrajpt);
730  }
731  rec.mc.part[mcpIndex].trkid = trackId;
732 
733  mcpIndex++;
734  }
735  }
736 
737  // Get handles for MCCaloInfo
738  art::Handle< std::vector<gar::sdp::CaloDeposit> > SimHitHandle;//ecal
739  art::Handle< std::vector<gar::sdp::CaloDeposit> > MuIDSimHitHandle;//ecal
740  art::InputTag ecalgeanttag(fGeantLabel, fInstanceLabelCalo);
741  art::InputTag muidgeanttag(fGeantLabel, fInstanceLabelMuID);
742  if (fWriteMCCaloInfo) {
743  SimHitHandle = e.getHandle< std::vector<gar::sdp::CaloDeposit> >(ecalgeanttag);
744  if (!SimHitHandle) {
745  throw cet::exception("CAFMaker") << " No gar::sdp::CaloDeposit branch."
746  << " Line " << __LINE__ << " in file " << __FILE__ << std::endl;
747  }
748  if(fWriteMuID) {
749  if (fGeo->HasMuonDetector())
750  {
751  MuIDSimHitHandle = e.getHandle< std::vector<gar::sdp::CaloDeposit> >(muidgeanttag);
752  if (!MuIDSimHitHandle) {
753  throw cet::exception("CAFMaker") << " No gar::sdp::CaloDeposit branch."
754  << " Line " << __LINE__ << " in file " << __FILE__ << std::endl;
755  }
756  }
757  }
758 
759  // Save simulation ecal hit info
760  for ( auto const& SimHit : (*SimHitHandle) ) {
761  caf::SRSimHit srsimhit;
762  srsimhit.x = SimHit.X();
763  srsimhit.y = SimHit.Y();
764  srsimhit.z = SimHit.Z();
765  srsimhit.t = SimHit.Time();
766  srsimhit.E = SimHit.Energy();
767  srsimhit.trkid = SimHit.TrackID();
768  srsimhit.cellid = SimHit.CellID();
769 
770  rec.mc.simhit.ecal.push_back(srsimhit);
771 
772  rec.mc.simhit.ecaltotE += SimHit.Energy();
773  }
774 
775  rec.mc.simhit.necal = rec.mc.simhit.ecal.size();
776 
777  // Save simulation muon system hit info
778  for ( auto const& SimHit : (*MuIDSimHitHandle) ) {
779  caf::SRSimHit srsimhit;
780  srsimhit.x = SimHit.X();
781  srsimhit.y = SimHit.Y();
782  srsimhit.z = SimHit.Z();
783  srsimhit.t = SimHit.Time();
784  srsimhit.E = SimHit.Energy();
785  srsimhit.trkid = SimHit.TrackID();
786  srsimhit.cellid = SimHit.CellID();
787 
788  rec.mc.simhit.muid.push_back(srsimhit);
789 
790  rec.mc.simhit.muidtotE += SimHit.Energy();
791  }
792 
793  rec.mc.simhit.nmuid = rec.mc.simhit.muid.size();
794  }
795 }
caf::SRTruth::gpart
std::vector< SRGenieParticle > gpart
Definition: SRTruth.h:27
caf::SRNeutrino::vtx
SRVector3D vtx
Definition: SRNeutrino.h:30
simb::MCParticle::PdgCode
int PdgCode() const
Definition: MCParticle.h:212
simb::MCNeutrino::CCNC
int CCNC() const
Definition: MCNeutrino.h:148
simb::MCNeutrino::QSqr
double QSqr() const
Definition: MCNeutrino.h:157
simb::MCNeutrino::Theta
double Theta() const
angle between incoming and outgoing leptons, in radians
Definition: MCNeutrino.cxx:63
caf::SRVector3D
A 3-vector with more efficient storage than TVector3.
Definition: SRVector3D.h:25
simb::MCParticle::Py
double Py(const int i=0) const
Definition: MCParticle.h:231
gar::CAFMaker::fInstanceLabelCalo
std::string fInstanceLabelCalo
Instance name for ECAL.
Definition: CAFMaker_module.cc:125
caf::SRSimHitBranch::necal
unsigned int necal
Definition: SRSimHitBranch.h:15
simb::MCParticle::EndZ
double EndZ() const
Definition: MCParticle.h:228
caf::SRSimHitBranch::muid
std::vector< SRSimHit > muid
Definition: SRSimHitBranch.h:20
caf::SRNeutrino::p
SRVector3D p
Definition: SRNeutrino.h:31
simb::MCNeutrino::Nu
const simb::MCParticle & Nu() const
Definition: MCNeutrino.h:146
simb::MCParticle::Px
double Px(const int i=0) const
Definition: MCParticle.h:230
gar::CAFMaker::fWriteMCCaloInfo
bool fWriteMCCaloInfo
Write MC info for calorimeter Default=true.
Definition: CAFMaker_module.cc:152
gar::CAFMaker::fWriteMCPTrajectory
bool fWriteMCPTrajectory
Write MCP Trajectory Default=true.
Definition: CAFMaker_module.cc:150
simb::MCParticle::EndY
double EndY() const
Definition: MCParticle.h:227
caf::SRParticle::proc
std::string proc
Definition: SRParticle.h:27
cet::hypot
std::enable_if_t< std::is_arithmetic_v< T >, T > hypot(T x, T y)
Definition: hypot.h:60
caf::SRParticle::start
SRVector3D start
Definition: SRParticle.h:22
gar::CAFMaker::fGeneratorLabels
std::vector< std::string > fGeneratorLabels
Definition: CAFMaker_module.cc:122
caf::SRParticle::end
SRVector3D end
Definition: SRParticle.h:25
simb::MCNeutrino::InteractionType
int InteractionType() const
Definition: MCNeutrino.h:150
abs
T abs(T value)
Definition: sparse_vector_test.cc:30
gar::CAFMaker::fGeo
const geo::GeometryCore * fGeo
pointer to the geometry
Definition: CAFMaker_module.cc:182
e
const double e
Definition: gUpMuFluxGen.cxx:165
gar::CAFMaker::TrackIdToIndex
std::unordered_map< TrkId, Int_t > TrackIdToIndex
Definition: CAFMaker_module.cc:185
simb::MCNeutrino::W
double W() const
Definition: MCNeutrino.h:154
simb::MCNeutrino::Y
double Y() const
Definition: MCNeutrino.h:156
gar::CAFMaker::fGeantLabel
std::string fGeantLabel
module label for geant4 simulated hits
Definition: CAFMaker_module.cc:124
gar::CAFMaker::fWriteCohInfo
bool fWriteCohInfo
MC level t for coherent pi+. Default=false.
Definition: CAFMaker_module.cc:170
caf::SRTruth::simhit
SRSimHitBranch simhit
Definition: SRTruth.h:33
simb::MCNeutrino::X
double X() const
Definition: MCNeutrino.h:155
caf::SRSimHitBranch::ecal
std::vector< SRSimHit > ecal
Definition: SRSimHitBranch.h:16
gar::CAFMaker::fMatchMCPtoVertDist
float fMatchMCPtoVertDist
MCParticle to MC vertex match Default=roundoff.
Definition: CAFMaker_module.cc:153
caf::SRTruth::nu
std::vector< SRNeutrino > nu
Definition: SRTruth.h:19
gar::CAFMaker::fGENIEGeneratorLabels
std::vector< std::string > fGENIEGeneratorLabels
Definition: CAFMaker_module.cc:123
caf::SRTruth::nnu
unsigned int nnu
Definition: SRTruth.h:18
caf::SRTruth::part
std::vector< SRParticle > part
Definition: SRTruth.h:30
caf::SRSimHit::cellid
size_t cellid
Definition: SRSimHit.h:17
larg4::dist
constexpr double dist(const TReal *x, const TReal *y, const unsigned int dimension)
Definition: OpFastScintillation.hh:563
caf::SRSimHitBranch::nmuid
unsigned int nmuid
Definition: SRSimHitBranch.h:19
gar::geo::GeometryCore::HasMuonDetector
bool HasMuonDetector() const
Definition: GeometryCore.h:999
caf::SRTruth::npart
unsigned int npart
Definition: SRTruth.h:29
MF_LOG_DEBUG
#define MF_LOG_DEBUG(id)
Definition: MessageLogger.h:268
simb::MCParticle::Pz
double Pz(const int i=0) const
Definition: MCParticle.h:232
caf::SRTruth::ngpart
unsigned int ngpart
Definition: SRTruth.h:26
gar::CAFMaker::fWriteMCPTrajMomenta
bool fWriteMCPTrajMomenta
Write Momenta associated with MCP Trajectory Default=false.
Definition: CAFMaker_module.cc:151
caf::SRParticle::startp
SRVector3D startp
Definition: SRParticle.h:24
caf::SRParticle::endproc
std::string endproc
Definition: SRParticle.h:28
reco_momentum_tuples.momentum
def momentum(x1, x2, x3, scale=1.)
Definition: reco_momentum_tuples.py:198
simb::MCParticle::EndX
double EndX() const
Definition: MCParticle.h:226
gar::CAFMaker::computeT
float computeT(simb::MCTruth theMCTruth)
Definition: CAFMaker_module.cc:1482
simb::MCNeutrino
Event generator information.
Definition: MCNeutrino.h:18
uint
unsigned uint
Definition: qglobal.h:351
simb::MCNeutrino::Mode
int Mode() const
Definition: MCNeutrino.h:149
gar::CAFMaker::fInstanceLabelMuID
std::string fInstanceLabelMuID
Instance name for MuID.
Definition: CAFMaker_module.cc:126
caf::SRParticle::endp
SRVector3D endp
Definition: SRParticle.h:26
fhicl::exception
cet::coded_exception< error, detail::translate > exception
Definition: exception.h:33
endl
QTextStream & endl(QTextStream &s)
Definition: qtextstream.cpp:2030
void gar::CAFMaker::FillHighLevelRecoInfo ( art::Event const &  e,
caf::StandardRecord &  rec 
)
private

Definition at line 944 of file CAFMaker_module.cc.

944  {
945 
946  // Get handle for TPCClusters
947  art::Handle< std::vector<rec::TPCCluster> > TPCClusterHandle;
948  if (fWriteTPCClusters) {
949  TPCClusterHandle = e.getHandle< std::vector<rec::TPCCluster> >(fTPCClusterLabel);
950  if (!TPCClusterHandle) {
951  throw cet::exception("CAFMaker") << " No rec::TPCCluster branch."
952  << " Line " << __LINE__ << " in file " << __FILE__ << std::endl;
953  }
954  }
955 
956  // Get handles for Tracks and their ionizations; also Assn's to TPCClusters, TrackIoniz
957  art::Handle< std::vector<rec::Track> > TrackHandle;
958  art::Handle< std::vector<rec::TrackIoniz> > TrackIonHandle;
959  art::Handle< std::vector<rec::TrackTrajectory> > TrackTrajHandle;
960  art::FindManyP<rec::TPCCluster>* findManyTPCClusters = NULL;
961  art::FindOneP<rec::TrackIoniz>* findIonization = NULL;
962  if(fWriteTracks) {
963  TrackHandle = e.getHandle< std::vector<rec::Track> >(fTrackLabel);
964  if (!TrackHandle) {
965  throw cet::exception("CAFMaker") << " No rec::Track branch."
966  << " Line " << __LINE__ << " in file " << __FILE__ << std::endl;
967  }
968  TrackIonHandle = e.getHandle< std::vector<rec::TrackIoniz> >(fTrackLabel);
969  if (!TrackIonHandle) {
970  throw cet::exception("CAFMaker") << " No rec::TrackIoniz branch."
971  << " Line " << __LINE__ << " in file " << __FILE__ << std::endl;
972  }
973 
974  findManyTPCClusters = new art::FindManyP<rec::TPCCluster>(TrackHandle,e,fTrackLabel);
975  findIonization = new art::FindOneP<rec::TrackIoniz>(TrackHandle,e,fTrackLabel);
976 
977  if(fWriteTrackTrajectories) {
978  TrackTrajHandle = e.getHandle< std::vector<rec::TrackTrajectory> >(fTrackTrajectoryLabel);
979  if (!TrackTrajHandle) {
980  throw cet::exception("CAFMaker") << " No rec::TrackTrajectory branch."
981  << " Line " << __LINE__ << " in file " << __FILE__ << std::endl;
982  }
983  }
984  }
985 
986  // Get handle for Vertices; also Assn's to Tracks
987  art::Handle< std::vector<rec::Vertex> > VertexHandle;
988  art::FindManyP<rec::Track, rec::TrackEnd>* findManyTrackEnd = NULL;
989  if (fWriteVertices) {
990  VertexHandle = e.getHandle< std::vector<rec::Vertex> >(fVertexLabel);
991  if (!VertexHandle) {
992  throw cet::exception("CAFMaker") << " No rec::Vertex branch."
993  << " Line " << __LINE__ << " in file " << __FILE__ << std::endl;
994  }
995  findManyTrackEnd = new art::FindManyP<rec::Track, rec::TrackEnd>(VertexHandle,e,fVertexLabel);
996  }
997 
998  // Get handle for Vees; also Assn's to Tracks
999  art::Handle< std::vector<rec::Vee> > VeeHandle;
1000  art::FindManyP<rec::Track, rec::TrackEnd>* findManyVeeTrackEnd = NULL;
1001  if (fWriteVees) {
1002  VeeHandle = e.getHandle< std::vector<rec::Vee> >(fVeeLabel);
1003  if (!VeeHandle) {
1004  throw cet::exception("CAFMaker") << " No rec::Vee branch."
1005  << " Line " << __LINE__ << " in file " << __FILE__ << std::endl;
1006  }
1007  findManyVeeTrackEnd = new art::FindManyP<rec::Track, rec::TrackEnd>(VeeHandle,e,fVeeLabel);
1008  }
1009 
1010  // Get handle for CaloClusters; also Assn for matching tracks
1011  art::InputTag ecalclustertag(fClusterLabel, fInstanceLabelCalo);
1012  art::InputTag muidclustertag(fClusterMuIDLabel, fInstanceLabelMuID);
1013  art::Handle< std::vector<rec::Cluster> > RecoClusterHandle;
1014  art::Handle< std::vector<rec::Cluster> > RecoClusterMuIDHandle;
1015  art::FindManyP<rec::Track, rec::TrackEnd>* findManyCALTrackEnd = NULL;
1016  // art::FindManyP<gar::rec::CaloHit>* findManyClusterRecoHit = NULL;
1017  // art::FindManyP<gar::rec::CaloHit>* findManyClusterMuIDHit = NULL;
1018 
1019  if (fWriteCaloClusters) {
1020  RecoClusterHandle = e.getHandle< std::vector<rec::Cluster> >(ecalclustertag);
1021  if (!RecoClusterHandle) {
1022  throw cet::exception("CAFMaker") << " No rec::Cluster branch."
1023  << " Line " << __LINE__ << " in file " << __FILE__ << std::endl;
1024  }
1025 
1026  if(fWriteMuID) {
1027  if (fGeo->HasMuonDetector())
1028  {
1029  RecoClusterMuIDHandle = e.getHandle< std::vector<rec::Cluster> >(muidclustertag);
1030  if (!RecoClusterMuIDHandle) {
1031  throw cet::exception("CAFMaker") << " No rec::Cluster (MuID) branch."
1032  << " Line " << __LINE__ << " in file " << __FILE__ << std::endl;
1033  }
1034  }
1035  }
1036 
1037  // findManyClusterRecoHit = new art::FindManyP<gar::rec::CaloHit>(RecoClusterHandle,e,ecalclustertag);
1038 
1039  // if (fGeo->HasMuonDetector() && fWriteMuID)
1040  // findManyClusterMuIDHit = new art::FindManyP<gar::rec::CaloHit>(RecoClusterMuIDHandle,e,muidclustertag);
1041 
1042  if (fWriteTracks)
1043  findManyCALTrackEnd = new art::FindManyP<rec::Track, rec::TrackEnd>(RecoClusterHandle,e,fECALAssnLabel);
1044  }
1045 
1046  // save clusters in the TPC. For some reason, can't get FindOneP<rec::Track> or
1047  // FindManyP<rec::Track> to work; seems the underlying Assn isn't found. Have
1048  // to FindManyP<TPCCluster> instead and iterate if (fWriteTracks). :(
1049  if (fWriteTPCClusters) {
1050  for ( auto const& TPCCluster : (*TPCClusterHandle) ) {
1051  caf::SRTPCCluster srclust;
1052  srclust.x = TPCCluster.Position()[0];
1053  srclust.y = TPCCluster.Position()[1];
1054  srclust.z = TPCCluster.Position()[2];
1055  srclust.sig = TPCCluster.Signal();
1056  srclust.rms = TPCCluster.RMS();
1057  const float* cov = TPCCluster.CovMatPacked();
1058  srclust.cov.xx = cov[0];
1059  srclust.cov.xy = cov[1];
1060  srclust.cov.xz = cov[2];
1061  srclust.cov.yy = cov[3];
1062  srclust.cov.yz = cov[4];
1063  srclust.cov.zz = cov[5];
1064 
1065  Int_t trackForThisTPCluster = -1;
1066  if (fWriteTracks) {
1067  size_t iTrack = 0;
1068  for ( auto const& track : (*TrackHandle) ) {
1069  for (size_t iCluster=0; iCluster<track.NHits(); iCluster++) {
1070  auto const& trackedCluster =
1071  *(findManyTPCClusters->at(iTrack).at(iCluster));
1072  if (TPCCluster==trackedCluster) {
1073  trackForThisTPCluster = track.getIDNumber();
1074  // No cluster is in 2 tracks (don't mess up dE/dx!)
1075  goto pushit; // break 2 loops
1076  }
1077  }
1078  iTrack++;
1079  }
1080  }
1081  pushit:
1082  srclust.trkid = trackForThisTPCluster;
1083 
1084  rec.clust.tpc.push_back(srclust);
1085  }
1086  rec.clust.ntpc = rec.clust.tpc.size();
1087  }
1088 
1089  // save per-track info
1090  if (fWriteTracks) {
1091  size_t iTrack = 0;
1092  for ( auto const& track : (*TrackHandle) ) {
1093  caf::SRTrack srtrk;
1094 
1095  // track is a rec::Track, not a rec::TrackPar
1096  srtrk.id = track.getIDNumber();
1097 
1098  srtrk.start = caf::SRVector3D(track.Vertex()[0],
1099  track.Vertex()[1],
1100  track.Vertex()[2]);
1101  srtrk.startp = caf::SRVector3D(track.Momentum_beg()*track.VtxDir()[0],
1102  track.Momentum_beg()*track.VtxDir()[1],
1103  track.Momentum_beg()*track.VtxDir()[2]);
1104  srtrk.startq = track.ChargeBeg();
1105 
1106  srtrk.end = caf::SRVector3D(track.End()[0],
1107  track.End()[1],
1108  track.End()[2]);
1109  srtrk.endp = caf::SRVector3D(track.Momentum_end()*track.EndDir()[0],
1110  track.Momentum_end()*track.EndDir()[1],
1111  track.Momentum_end()*track.EndDir()[2]);
1112  srtrk.endq = track.ChargeEnd();
1113 
1114  srtrk.fwd.len = track.LengthForward();
1115  srtrk.bak.len = track.LengthBackward();
1116  srtrk.fwd.chi2 = track.ChisqForward();
1117  srtrk.bak.chi2 = track.ChisqBackward();
1118  srtrk.NTPCClustersOnTrack = track.NHits();
1119 
1120  //Add the PID information based on Tom's parametrization
1121  TVector3 momF(track.Momentum_beg()*track.VtxDir()[0], track.Momentum_beg()*track.VtxDir()[1], track.Momentum_beg()*track.VtxDir()[2]);
1122  TVector3 momB(track.Momentum_end()*track.EndDir()[0], track.Momentum_end()*track.EndDir()[1], track.Momentum_end()*track.EndDir()[2]);
1123 
1124  //Reconstructed momentum forward and backward
1125  float pF = momF.Mag();
1126  float pB = momB.Mag();
1127  std::vector< std::pair<int, float> > pidF = processPIDInfo( pF );
1128  std::vector< std::pair<int, float> > pidB = processPIDInfo( pB );
1129 
1130  //Fill the pid and its probability
1131  for(size_t ipid = 0; ipid < pidF.size(); ipid++) {
1132  caf::SRTrackPID srpid;
1133  srpid.id = pidF.at(ipid).first;
1134  srpid.prob = pidF.at(ipid).second;
1135  srtrk.fwd.pid.push_back(srpid);
1136  }
1137  for(size_t ipid = 0; ipid < pidB.size(); ipid++) {
1138  caf::SRTrackPID srpid;
1139  srpid.id = pidB.at(ipid).first;
1140  srpid.prob = pidB.at(ipid).second;
1141  srtrk.bak.pid.push_back(srpid);
1142  }
1143 
1144  // Matching MCParticle info
1145  std::vector<std::pair<simb::MCParticle*,float>> trakt;
1146  trakt = BackTrack->TrackToMCParticles( const_cast<rec::Track*>(&track) );
1147  srtrk.mcidx = -1;
1148  if (trakt.size()>0 && TrackIdToIndex.size()!=0) {
1149  srtrk.mcidx = TrackIdToIndex[trakt[0].first->TrackId()];
1150  }
1151 
1152  if (srtrk.mcidx > -1) {
1153  srtrk.mcfrac = trakt[0].second;
1154  } else {
1155  srtrk.mcfrac = 0;
1156  }
1157 
1158  if (findIonization->isValid()) {
1159  // No calibration for now. Someday this should all be in reco
1160  rec::TrackIoniz ionization = *(findIonization->at(iTrack));
1161  float avgIonF, avgIonB;
1162  processIonizationInfo(ionization, fIonizTruncate, avgIonF, avgIonB);
1163  srtrk.fwd.avgion = avgIonF;
1164  srtrk.bak.avgion = avgIonB;
1165  } else {
1166  srtrk.fwd.avgion = 0;
1167  srtrk.bak.avgion = 0;
1168  }
1169  iTrack++;
1170 
1171  rec.trk.push_back(srtrk);
1172  } // end loop over TrackHandle
1173  }
1174 
1175  rec.ntrk = rec.trk.size();
1176 
1177  //TrackTrajectories
1178  if(fWriteTrackTrajectories) {
1179  size_t iTrackTraj = 0;
1180  for ( auto const& tracktraj : (*TrackTrajHandle) ) {
1181  caf::SRTrack& srtrk = rec.trk[iTrackTraj];
1182 
1183  std::vector<TVector3> temp = tracktraj.getFWDTrajectory();
1184  for(size_t i = 0; i < temp.size(); i++) {
1185  srtrk.fwd.traj.emplace_back(temp.at(i).X(),
1186  temp.at(i).Y(),
1187  temp.at(i).Z());
1188  }
1189 
1190  temp = tracktraj.getBAKTrajectory();
1191  for(size_t i = 0; i < temp.size(); i++) {
1192  srtrk.bak.traj.emplace_back(temp.at(i).X(),
1193  temp.at(i).Y(),
1194  temp.at(i).Z());
1195  }
1196  iTrackTraj++;
1197  }
1198  }
1199 
1200  // save Vertex and Track-Vertex association info
1201  if (fWriteVertices) {
1202  size_t iVertex = 0;
1203  for ( auto const& vertex : (*VertexHandle) ) {
1204  caf::SRVertex srvtx;
1205  srvtx.id = vertex.getIDNumber();
1206  srvtx.x = vertex.Position()[0];
1207  srvtx.y = vertex.Position()[1];
1208  srvtx.z = vertex.Position()[2];
1209  srvtx.t = vertex.Time();
1210 
1211  srvtx.ntrks = 0;
1212  if ( findManyTrackEnd->isValid() ) {
1213  srvtx.ntrks = findManyTrackEnd->at(iVertex).size();
1214  }
1215 
1216  int vertexCharge = 0;
1217  for (int iVertexedTrack=0; iVertexedTrack < srvtx.ntrks; ++iVertexedTrack) {
1218  // TODO can we just make a vector of tracks belonging to this vtx?
1219  caf::SRVertexAssn srassn;
1220  srassn.vtxid = vertex.getIDNumber();
1221 
1222  // Get this vertexed track.
1223  rec::Track track = *(findManyTrackEnd->at(iVertex).at(iVertexedTrack));
1224  srassn.trkid = track.getIDNumber();
1225 
1226  // Get the end of the track in the vertex. It isn't that odd for the end
1227  // of the track not used in the vertex to be closer to the vertex than the
1228  // one actually used; you might have a very short stub track in a 3 track
1229  // vertex with small opening angles and the other 2 tracks might pull the
1230  // vertex some distance towards the far end of the stub track
1231  rec::TrackEnd fee = *(findManyTrackEnd->data(iVertex).at(iVertexedTrack));
1232  // TrackEnd is defined in Track.h; 1 means use Beg values, 0 means use End
1233  srassn.trkend = fee;
1234 
1235  rec.assn.vtx.push_back(srassn);
1236 
1237  if (fee==rec::TrackEndBeg) {
1238  vertexCharge += track.ChargeBeg();
1239  } else {
1240  vertexCharge += track.ChargeEnd();
1241  }
1242  }
1243  srvtx.Q = vertexCharge;
1244  ++iVertex;
1245 
1246  rec.vtx.push_back(srvtx);
1247  } // end loop over VertexHandle
1248  }
1249 
1250  rec.nvtx = rec.vtx.size();
1251  rec.assn.nvtx = rec.assn.vtx.size();
1252 
1253  // save Vee and Track-Vee association info
1254  if (fWriteVees) {
1255  size_t iVee = 0;
1256  for ( auto const& vee : (*VeeHandle) ) {
1257  caf::SRVee srvee;
1258  srvee.id = vee.getIDNumber();
1259  srvee.x = vee.Position()[0];
1260  srvee.y = vee.Position()[1];
1261  srvee.z = vee.Position()[2];
1262  srvee.t = vee.Time();
1263  srvee.Kpipi.p.x = vee.FourMomentum(0).X();
1264  srvee.Kpipi.p.y = vee.FourMomentum(0).Y();
1265  srvee.Kpipi.p.z = vee.FourMomentum(0).Z();
1266  srvee.Kpipi.E = vee.FourMomentum(0).E();
1267  srvee.Kpipi.m = vee.FourMomentum(0).M();
1268  srvee.Lppi.p.x = vee.FourMomentum(1).X();
1269  srvee.Lppi.p.y = vee.FourMomentum(1).Y();
1270  srvee.Lppi.p.z = vee.FourMomentum(1).Z();
1271  srvee.Lppi.E = vee.FourMomentum(1).E();
1272  srvee.Lppi.m = vee.FourMomentum(1).M();
1273  srvee.Lpip.p.x = vee.FourMomentum(2).X();
1274  srvee.Lpip.p.y = vee.FourMomentum(2).Y();
1275  srvee.Lpip.p.z = vee.FourMomentum(2).Z();
1276  srvee.Lpip.E = vee.FourMomentum(2).E();
1277  srvee.Lpip.m = vee.FourMomentum(2).M();
1278 
1279  rec.vee.push_back(srvee);
1280 
1281  int nVeeTracks = 0;
1282  if ( findManyVeeTrackEnd->isValid() ) {
1283  nVeeTracks = findManyVeeTrackEnd->at(iVee).size();
1284  }
1285 
1286  for (int iVeeTrack=0; iVeeTrack<nVeeTracks; ++iVeeTrack) {
1287  caf::SRVeeAssn srassn;
1288  srassn.veeid = vee.getIDNumber();
1289 
1290  // Get this vertexed track.
1291  rec::Track track = *(findManyVeeTrackEnd->at(iVee).at(iVeeTrack));
1292  srassn.trkid = track.getIDNumber();
1293 
1294  rec::TrackEnd fee = *(findManyVeeTrackEnd->data(iVee).at(iVeeTrack));
1295  // TrackEnd is defined in Track.h; 1 means use Beg values, 0 means use End
1296  srassn.trkend = fee;
1297 
1298  rec.assn.vee.push_back(srassn);
1299  }
1300  ++iVee;
1301  } // end loop over VeeHandle
1302  }
1303 
1304  rec.nvee = rec.vee.size();
1305  rec.assn.nvee = rec.assn.vee.size();
1306 
1307  // save Cluster info
1308  if (fWriteCaloClusters) {
1309  for ( auto const& cluster : (*RecoClusterHandle) ) {
1310  caf::SRECalCluster srclust;
1311  srclust.id = cluster.getIDNumber();
1312  srclust.nhits = cluster.CalorimeterHits().size();
1313  srclust.E = cluster.Energy();
1314  srclust.t = cluster.Time();
1315  srclust.TimeDiffFirstLast = cluster.TimeDiffFirstLast();
1316  srclust.x = cluster.Position()[0];
1317  srclust.y = cluster.Position()[1];
1318  srclust.z = cluster.Position()[2];
1319  srclust.theta = cluster.ITheta();
1320  srclust.phi = cluster.IPhi();
1321  srclust.pid = cluster.ParticleID();
1322  // srclust.Shape = cluster.Shape();
1323  srclust.mainAxis.x = cluster.EigenVectors()[0];
1324  srclust.mainAxis.y = cluster.EigenVectors()[1];
1325  srclust.mainAxis.z = cluster.EigenVectors()[2];
1326 
1327  // Matching MCParticle info
1328  std::vector<std::pair<simb::MCParticle*,float>> trakt;
1329  trakt = BackTrack->ClusterToMCParticles( const_cast<rec::Cluster*>(&cluster) );
1330  srclust.mcidx = -1;
1331  if (trakt.size()>0 && TrackIdToIndex.size()!=0) {
1332  srclust.mcidx = TrackIdToIndex[trakt[0].first->TrackId()];
1333  }
1334  if (srclust.mcidx > -1) {
1335  srclust.mcfrac = trakt[0].second;
1336  } else {
1337  srclust.mcfrac = 0;
1338  }
1339 
1340  rec.clust.ecal.push_back(srclust);
1341  }
1342 
1343  if(fGeo->HasMuonDetector() && fWriteMuID) {
1344  for ( auto const& cluster : (*RecoClusterMuIDHandle) ) {
1345  caf::SRMuIDCluster srclust;
1346  srclust.id = cluster.getIDNumber();
1347  srclust.nhits = cluster.CalorimeterHits().size();
1348  srclust.E = cluster.Energy();
1349  srclust.t = cluster.Time();
1350  srclust.TimeDiffFirstLast = cluster.TimeDiffFirstLast();
1351  srclust.x = cluster.Position()[0];
1352  srclust.y = cluster.Position()[1];
1353  srclust.z = cluster.Position()[2];
1354  srclust.theta = cluster.ITheta();
1355  srclust.phi = cluster.IPhi();
1356  srclust.pid = cluster.ParticleID();
1357  // srclust.Shape = cluster.Shape();
1358  srclust.mainAxis.x = cluster.EigenVectors()[0];
1359  srclust.mainAxis.y = cluster.EigenVectors()[1];
1360  srclust.mainAxis.z = cluster.EigenVectors()[2];
1361 
1362  // Matching MCParticle info
1363  std::vector<std::pair<simb::MCParticle*,float>> trakt;
1364  trakt = BackTrack->ClusterToMCParticles( const_cast<rec::Cluster*>(&cluster) );
1365  srclust.mcidx = -1;
1366  if (trakt.size()>0 && TrackIdToIndex.size()!=0) {
1367  srclust.mcidx = TrackIdToIndex[trakt[0].first->TrackId()];
1368  }
1369  if (srclust.mcidx > -1) {
1370  srclust.mcfrac = trakt[0].second;
1371  } else {
1372  srclust.mcfrac = 0;
1373  }
1374 
1375  rec.clust.muid.push_back(srclust);
1376  }
1377  }
1378  }
1379 
1380  rec.clust.necal = rec.clust.ecal.size();
1381  rec.clust.nmuid = rec.clust.muid.size();
1382 
1383  // Write info for ECAL-matched tracks
1384  if (fWriteMatchedTracks) {
1385  size_t iCluster = 0;
1386  for ( auto const& cluster : (*RecoClusterHandle) ) {
1387  int nCALedTracks(0);
1388  if ( findManyCALTrackEnd->isValid() ) {
1389  nCALedTracks = findManyCALTrackEnd->at(iCluster).size();
1390  }
1391  for (int iCALedTrack=0; iCALedTrack<nCALedTracks; ++iCALedTrack) {
1392  caf::SRClusterAssn srassn;
1393  srassn.clustid = cluster.getIDNumber();
1394  rec::Track track = *(findManyCALTrackEnd->at(iCluster).at(iCALedTrack));
1395  srassn.trkid = track.getIDNumber();
1396 
1397  rec::TrackEnd fee = *(findManyCALTrackEnd->data(iCluster).at(iCALedTrack));
1398  srassn.trkend = fee; // The rec::TrackEnd (see Track.h) that extrapolated to cluster
1399  rec.assn.ecalclust.push_back(srassn);
1400  }
1401  iCluster++;
1402  }
1403  } // end branch on fWriteCaloInfo
1404 
1405  rec.assn.necalclust = rec.assn.ecalclust.size();
1406 
1407  return;
1408 } // end :CAFMaker::FillVectors
gar::CAFMaker::fClusterLabel
std::string fClusterLabel
module label for calo clusters rec::Cluster
Definition: CAFMaker_module.cc:141
caf::SRTrackTrajectory::pid
std::vector< SRTrackPID > pid
Definition: SRTrackTrajectory.h:20
caf::SRVee::y
float y
Definition: SRVee.h:16
gar::CAFMaker::fWriteCaloClusters
bool fWriteCaloClusters
Write ECAL clusters. Default=true.
Definition: CAFMaker_module.cc:165
gar::rec::TrackEnd
int TrackEnd
Definition: Track.h:32
caf::SRVector3D
A 3-vector with more efficient storage than TVector3.
Definition: SRVector3D.h:25
caf::SRMuIDCluster::mainAxis
SRVector3D mainAxis
Definition: SRMuIDCluster.h:23
gar::CAFMaker::fInstanceLabelCalo
std::string fInstanceLabelCalo
Instance name for ECAL.
Definition: CAFMaker_module.cc:125
caf::SRVee::id
size_t id
Definition: SRVee.h:15
caf::SRVertex::id
size_t id
Definition: SRVertex.h:14
gar::CAFMaker::fTPCClusterLabel
std::string fTPCClusterLabel
module label for TPC Clusters rec::TPCCluster
Definition: CAFMaker_module.cc:129
gar::CAFMaker::fWriteTPCClusters
bool fWriteTPCClusters
Write TPCClusters info Default=true.
Definition: CAFMaker_module.cc:156
gar::CAFMaker::fWriteVertices
bool fWriteVertices
Reco vertexes & their tracks Default=true.
Definition: CAFMaker_module.cc:159
caf::SRTrack::bak
SRTrackTrajectory bak
Definition: SRTrack.h:25
caf::SRTrack::NTPCClustersOnTrack
int NTPCClustersOnTrack
Definition: SRTrack.h:23
caf::SRVeeAssn::trkid
size_t trkid
Definition: SRVeeAssn.h:14
cluster
Cluster finding and building.
Definition: SmallClusterFilter_module.cc:34
caf::SRTrack::startp
SRVector3D startp
Definition: SRTrack.h:16
caf::SRAssnBranch::ecalclust
std::vector< SRClusterAssn > ecalclust
Definition: SRAssnBranch.h:24
caf::StandardRecord::clust
SRClusterBranch clust
Definition: StandardRecord.h:41
caf::SRTrack::startq
int startq
Definition: SRTrack.h:17
caf::SRVee::Lpip
SRVeeHypothesis Lpip
Definition: SRVee.h:19
caf::SRAssnBranch::vee
std::vector< SRVeeAssn > vee
Definition: SRAssnBranch.h:21
gar::CAFMaker::fWriteTracks
bool fWriteTracks
Start/end X, P for tracks Default=true.
Definition: CAFMaker_module.cc:157
gar::CAFMaker::fTrackLabel
std::string fTrackLabel
module label for TPC Tracks rec:Track
Definition: CAFMaker_module.cc:130
caf::SRVertex::t
size_t t
Definition: SRVertex.h:16
gar::CAFMaker::processIonizationInfo
void processIonizationInfo(rec::TrackIoniz &ion, float ionizeTruncate, float &forwardIonVal, float &backwardIonVal)
Definition: CAFMaker_module.cc:1416
caf::SRECalCluster::nhits
unsigned int nhits
Definition: SRECalCluster.h:16
gar::cheat::BackTrackerCore::ClusterToMCParticles
std::vector< std::pair< simb::MCParticle *, float > > ClusterToMCParticles(rec::Cluster *const c)
Definition: BackTrackerCore.cxx:933
caf::SRTrack::id
size_t id
Definition: SRTrack.h:14
gar::CAFMaker::fGeo
const geo::GeometryCore * fGeo
pointer to the geometry
Definition: CAFMaker_module.cc:182
e
const double e
Definition: gUpMuFluxGen.cxx:165
gar::CAFMaker::TrackIdToIndex
std::unordered_map< TrkId, Int_t > TrackIdToIndex
Definition: CAFMaker_module.cc:185
gar::CAFMaker::fWriteVees
bool fWriteVees
Reco vees & their tracks Default=true.
Definition: CAFMaker_module.cc:160
caf::SRVee::Kpipi
SRVeeHypothesis Kpipi
Definition: SRVee.h:19
caf::SRTrack::mcfrac
float mcfrac
Definition: SRTrack.h:28
caf::SRTrack::mcidx
int mcidx
Definition: SRTrack.h:27
submit_ppfxjobs.temp
temp
Definition: submit_ppfxjobs.py:110
caf::StandardRecord::vee
std::vector< SRVee > vee
Definition: StandardRecord.h:35
caf::SRTrackTrajectory::avgion
float avgion
average ionization
Definition: SRTrackTrajectory.h:18
gar::CAFMaker::BackTrack
cheat::BackTrackerCore * BackTrack
Definition: CAFMaker_module.cc:95
gar::CAFMaker::fWriteMatchedTracks
bool fWriteMatchedTracks
Write ECAL-track Assns Default=true.
Definition: CAFMaker_module.cc:166
caf::SRTrack::endq
int endq
Definition: SRTrack.h:21
caf::StandardRecord::trk
std::vector< SRTrack > trk
Definition: StandardRecord.h:29
caf::SRTrack::start
SRVector3D start
Track 3D start point.
Definition: SRTrack.h:18
caf::SRMuIDCluster::nhits
unsigned int nhits
Definition: SRMuIDCluster.h:16
caf::SRTrack::fwd
SRTrackTrajectory fwd
Definition: SRTrack.h:25
caf::SRVee::Lppi
SRVeeHypothesis Lppi
Definition: SRVee.h:19
caf::SRTrack::end
SRVector3D end
Track 3D end point.
Definition: SRTrack.h:19
gar::CAFMaker::fClusterMuIDLabel
std::string fClusterMuIDLabel
module label for calo clusters rec::Cluster in MuID
Definition: CAFMaker_module.cc:142
caf::SRAssnBranch::nvtx
unsigned int nvtx
Definition: SRAssnBranch.h:17
caf::SRVeeAssn::veeid
size_t veeid
Definition: SRVeeAssn.h:13
gar::CAFMaker::fVertexLabel
std::string fVertexLabel
module label for vertexes rec:Vertex
Definition: CAFMaker_module.cc:132
caf::SRTrackTrajectory::traj
std::vector< SRVector3D > traj
Definition: SRTrackTrajectory.h:22
caf::SRClusterBranch::muid
std::vector< SRMuIDCluster > muid
Definition: SRClusterBranch.h:24
gar::CAFMaker::processPIDInfo
std::vector< std::pair< int, float > > processPIDInfo(float p)
Definition: CAFMaker_module.cc:1539
gar::CAFMaker::fTrackTrajectoryLabel
std::string fTrackTrajectoryLabel
module label for TPC Track Trajectories rec:TrackTrajectory
Definition: CAFMaker_module.cc:131
caf::SRCovMx::xx
float xx
Definition: SRCovMx.h:11
caf::SRECalCluster::mainAxis
SRVector3D mainAxis
Definition: SRECalCluster.h:23
gar::rec::TrackEndBeg
TrackEnd const TrackEndBeg
Definition: Track.h:33
gar::CAFMaker::fWriteTrackTrajectories
bool fWriteTrackTrajectories
Point traj of reco tracks Default=false.
Definition: CAFMaker_module.cc:158
gar::CAFMaker::fIonizTruncate
float fIonizTruncate
Default=1.00;.
Definition: CAFMaker_module.cc:169
gar::cheat::BackTrackerCore::TrackToMCParticles
std::vector< std::pair< simb::MCParticle *, float > > TrackToMCParticles(rec::Track *const t)
Definition: BackTrackerCore.cxx:823
caf::SRAssnBranch::nvee
unsigned int nvee
Definition: SRAssnBranch.h:20
gar::geo::GeometryCore::HasMuonDetector
bool HasMuonDetector() const
Definition: GeometryCore.h:999
caf::SRTrack::endp
SRVector3D endp
Definition: SRTrack.h:20
caf::SRCovMx::xy
float xy
Definition: SRCovMx.h:11
gar::CAFMaker::fVeeLabel
std::string fVeeLabel
module label for conversion/decay vertexes rec:Vee
Definition: CAFMaker_module.cc:133
caf::SRAssnBranch::necalclust
unsigned int necalclust
Definition: SRAssnBranch.h:23
caf::SRCovMx::yz
float yz
Definition: SRCovMx.h:11
proxy::Track
TrackCollectionProxyElement< TrackCollProxy > Track
Proxy to an element of a proxy collection of recob::Track objects.
Definition: Track.h:1036
caf::SRVee::t
size_t t
Definition: SRVee.h:17
caf::SRCovMx::xz
float xz
Definition: SRCovMx.h:11
caf::SRVee::x
float x
Definition: SRVee.h:16
caf::SRCovMx::yy
float yy
Definition: SRCovMx.h:11
caf::SRCovMx::zz
float zz
Definition: SRCovMx.h:11
caf::SRClusterBranch::ecal
std::vector< SRECalCluster > ecal
Definition: SRClusterBranch.h:21
gar::CAFMaker::fInstanceLabelMuID
std::string fInstanceLabelMuID
Instance name for MuID.
Definition: CAFMaker_module.cc:126
caf::SRVee::z
float z
Definition: SRVee.h:16
caf::SRAssnBranch::vtx
std::vector< SRVertexAssn > vtx
Definition: SRAssnBranch.h:18
caf::StandardRecord::vtx
std::vector< SRVertex > vtx
Definition: StandardRecord.h:32
fhicl::exception
cet::coded_exception< error, detail::translate > exception
Definition: exception.h:33
endl
QTextStream & endl(QTextStream &s)
Definition: qtextstream.cpp:2030
gar::CAFMaker::fECALAssnLabel
std::string fECALAssnLabel
module label for track-clusters associations
Definition: CAFMaker_module.cc:144
caf::SRClusterBranch::tpc
std::vector< SRTPCCluster > tpc
Definition: SRClusterBranch.h:18
vertex
vertex reconstruction
Definition: AggregateVertexAna_module.cc:32
void gar::CAFMaker::FillRawInfo ( art::Event const &  e,
caf::StandardRecord &  rec 
)
private

Definition at line 802 of file CAFMaker_module.cc.

802  {
803 
804 
805  // Get handle for CaloDigits
806  art::InputTag ecalrawtag(fRawCaloHitLabel, fInstanceLabelCalo);
807  auto RawHitHandle = e.getHandle<std::vector<gar::raw::CaloRawDigit> >(ecalrawtag);
808  if (!RawHitHandle) {
809  throw cet::exception("CAFMaker") << " No :raw::CaloRawDigit branch."
810  << " Line " << __LINE__ << " in file " << __FILE__ << std::endl;
811  }
812 
813  art::InputTag muidrawtag(fRawMuIDHitLabel, fInstanceLabelMuID);
814  art::Handle< std::vector<gar::raw::CaloRawDigit> > MuIDRawHitHandle;
815  if (fGeo->HasMuonDetector())
816  {
817  MuIDRawHitHandle = e.getHandle< std::vector<gar::raw::CaloRawDigit> >(muidrawtag);
818  if (!MuIDRawHitHandle) {
819  throw cet::exception("CAFMaker") << " No :raw::CaloRawDigit branch."
820  << " Line " << __LINE__ << " in file " << __FILE__ << std::endl;
821  }
822  }
823 
824  // save ecal raw digits info
825  for ( auto const& DigiHit : (*RawHitHandle) ) {
826  caf::SRDigit srdig;
827  srdig.x = DigiHit.X();
828  srdig.y = DigiHit.Y();
829  srdig.z = DigiHit.Z();
830  srdig.t = (DigiHit.Time().first + DigiHit.Time().second) / 2.0 ;
831  srdig.adc = DigiHit.ADC().first;
832  srdig.cellid = DigiHit.CellID();
833 
834  rec.dig.ecal.push_back(srdig);
835  }
836  rec.dig.necal = rec.dig.ecal.size();
837 
838  // save muon system raw digits info
839  for ( auto const& DigiHit : (*MuIDRawHitHandle) ) {
840  caf::SRDigit srdig;
841  srdig.x = DigiHit.X();
842  srdig.y = DigiHit.Y();
843  srdig.z = DigiHit.Z();
844  srdig.t = (DigiHit.Time().first + DigiHit.Time().second) / 2.0 ;
845  srdig.adc = DigiHit.ADC().first;
846  srdig.cellid = DigiHit.CellID();
847  }
848  rec.dig.nmuid = rec.dig.muid.size();
849 }
gar::CAFMaker::fRawCaloHitLabel
std::string fRawCaloHitLabel
module label for digitized calo hits raw::CaloRawDigit
Definition: CAFMaker_module.cc:135
gar::CAFMaker::fInstanceLabelCalo
std::string fInstanceLabelCalo
Instance name for ECAL.
Definition: CAFMaker_module.cc:125
caf::StandardRecord::dig
SRDigitBranch dig
Definition: StandardRecord.h:37
caf::SRDigit::y
float y
Definition: SRDigit.h:13
caf::SRDigitBranch::muid
std::vector< SRDigit > muid
Definition: SRDigitBranch.h:19
caf::SRDigitBranch::nmuid
unsigned int nmuid
Definition: SRDigitBranch.h:18
caf::SRDigit::t
float t
Definition: SRDigit.h:14
gar::CAFMaker::fRawMuIDHitLabel
std::string fRawMuIDHitLabel
Definition: CAFMaker_module.cc:136
caf::SRDigit::z
float z
Definition: SRDigit.h:13
caf::SRDigitBranch::ecal
std::vector< SRDigit > ecal
Definition: SRDigitBranch.h:16
gar::CAFMaker::fGeo
const geo::GeometryCore * fGeo
pointer to the geometry
Definition: CAFMaker_module.cc:182
e
const double e
Definition: gUpMuFluxGen.cxx:165
caf::SRDigit::cellid
size_t cellid
Definition: SRDigit.h:16
caf::SRDigitBranch::necal
unsigned int necal
Definition: SRDigitBranch.h:15
gar::geo::GeometryCore::HasMuonDetector
bool HasMuonDetector() const
Definition: GeometryCore.h:999
gar::CAFMaker::fInstanceLabelMuID
std::string fInstanceLabelMuID
Instance name for MuID.
Definition: CAFMaker_module.cc:126
caf::SRDigit::adc
unsigned int adc
Definition: SRDigit.h:15
fhicl::exception
cet::coded_exception< error, detail::translate > exception
Definition: exception.h:33
endl
QTextStream & endl(QTextStream &s)
Definition: qtextstream.cpp:2030
caf::SRDigit::x
float x
Definition: SRDigit.h:13
void gar::CAFMaker::FillRecoInfo ( art::Event const &  e,
caf::StandardRecord &  rec 
)
private

Definition at line 856 of file CAFMaker_module.cc.

856  {
857 
858  // Get handle for TPC hit data
859  if (fWriteHits) {
860  auto HitHandle = e.getHandle< std::vector<rec::Hit> >(fHitLabel);
861  if (!HitHandle) {
862  throw cet::exception("CAFMaker") << " No rec::Hit branch."
863  << " Line " << __LINE__ << " in file " << __FILE__ << std::endl;
864  }
865 
866  // save hits in the TPC
867  for ( auto const& Hit : (*HitHandle) ) {
868  caf::SRTPCRecoHit srhit;
869  srhit.x = Hit.Position()[0];
870  srhit.y = Hit.Position()[1];
871  srhit.z = Hit.Position()[2];
872  srhit.sig = Hit.Signal();
873  srhit.rms = Hit.RMS();
874  srhit.chan = Hit.Channel();
875  rec.hit.tpc.push_back(srhit);
876  }
877  rec.hit.ntpc = rec.hit.tpc.size();
878  }
879 
880  // Get handle for CaloHits
881 
882  if (fWriteCaloHits) {
883 
884  art::InputTag ecalrecotag(fCaloHitLabel, fInstanceLabelCalo);
885  auto RecoHitHandle = e.getHandle< std::vector<rec::CaloHit> >(ecalrecotag);
886  if (!RecoHitHandle) {
887  throw cet::exception("CAFMaker") << " No rec::CaloHit branch."
888  << " Line " << __LINE__ << " in file " << __FILE__ << std::endl;
889  }
890 
891  // save reco'd Calorimetry hits
892  for ( auto const& Hit : (*RecoHitHandle) ) {
893  caf::SRRecoHit srhit;
894  srhit.id = Hit.getIDNumber();
895  srhit.x = Hit.Position()[0];
896  srhit.y = Hit.Position()[1];
897  srhit.z = Hit.Position()[2];
898  srhit.t = Hit.Time().first;
899  srhit.E = Hit.Energy();
900  srhit.cellid = Hit.CellID();
901  rec.hit.ecal.push_back(srhit);
902 
903  rec.hit.ecaltotE += srhit.E;
904  }
905  rec.hit.necal = rec.hit.ecal.size();
906 
907 
908  if(fWriteMuID) {
909  art::InputTag muirecotag(fMuIDHitLabel, fInstanceLabelMuID);
910  art::Handle< std::vector<rec::CaloHit> > MuIDRecoHitHandle;
911  if (fGeo->HasMuonDetector())
912  {
913  MuIDRecoHitHandle = e.getHandle< std::vector<rec::CaloHit> >(muirecotag);
914  if(!MuIDRecoHitHandle) {
915  throw cet::exception("CAFMaker") << " No rec::CaloHit branch."
916  << " Line " << __LINE__ << " in file " << __FILE__ << std::endl;
917  }
918  }
919  for ( auto const& Hit : (*MuIDRecoHitHandle) ) {
920  caf::SRRecoHit srhit;
921  srhit.id = Hit.getIDNumber();
922  srhit.x = Hit.Position()[0];
923  srhit.y = Hit.Position()[1];
924  srhit.z = Hit.Position()[2];
925  srhit.t = Hit.Time().first;
926  srhit.E = Hit.Energy();
927  srhit.cellid = Hit.CellID();
928  rec.hit.muid.push_back(srhit);
929 
930  rec.hit.muidtotE += srhit.E;
931  }
932  rec.hit.nmuid = rec.hit.muid.size();
933 
934  }
935 
936  }
937 }
gar::CAFMaker::fInstanceLabelCalo
std::string fInstanceLabelCalo
Instance name for ECAL.
Definition: CAFMaker_module.cc:125
gar::CAFMaker::fWriteCaloHits
bool fWriteCaloHits
Write ECAL hits. Default=true.
Definition: CAFMaker_module.cc:164
caf::SRRecoHitBranch::muid
std::vector< SRRecoHit > muid
Definition: SRRecoHitBranch.h:25
caf::SRRecoHitBranch::ecal
std::vector< SRRecoHit > ecal
Definition: SRRecoHitBranch.h:21
gar::CAFMaker::fMuIDHitLabel
std::string fMuIDHitLabel
Definition: CAFMaker_module.cc:139
gar::CAFMaker::fGeo
const geo::GeometryCore * fGeo
pointer to the geometry
Definition: CAFMaker_module.cc:182
e
const double e
Definition: gUpMuFluxGen.cxx:165
caf::SRRecoHitBranch::tpc
std::vector< SRTPCRecoHit > tpc
Definition: SRRecoHitBranch.h:17
caf::SRRecoHit::cellid
size_t cellid
Definition: SRRecoHit.h:17
caf::StandardRecord::hit
SRRecoHitBranch hit
Definition: StandardRecord.h:39
gar::CAFMaker::fCaloHitLabel
std::string fCaloHitLabel
module label for reco calo hits rec::CaloHit
Definition: CAFMaker_module.cc:138
caf::SRTPCRecoHit::chan
unsigned int chan
Definition: SRTPCRecoHit.h:16
gar::CAFMaker::fHitLabel
std::string fHitLabel
module label for reco TPC hits rec::Hit
Definition: CAFMaker_module.cc:128
gar::geo::GeometryCore::HasMuonDetector
bool HasMuonDetector() const
Definition: GeometryCore.h:999
gar::CAFMaker::fWriteHits
bool fWriteHits
Write info about TPC Hits Default=false.
Definition: CAFMaker_module.cc:155
gar::CAFMaker::fInstanceLabelMuID
std::string fInstanceLabelMuID
Instance name for MuID.
Definition: CAFMaker_module.cc:126
fhicl::exception
cet::coded_exception< error, detail::translate > exception
Definition: exception.h:33
endl
QTextStream & endl(QTextStream &s)
Definition: qtextstream.cpp:2030
static bool gar::CAFMaker::lessThan_byE ( std::pair< float, float >  a,
std::pair< float, float >  b 
)
inlinestaticprivate

Definition at line 107 of file CAFMaker_module.cc.

108  {return a.first < b.first;}
a
const double a
Definition: gUpMuFluxGen.cxx:164
b
static bool * b
Definition: config.cpp:1043
CAFMaker& gar::CAFMaker::operator= ( CAFMaker const &  )
delete
CAFMaker& gar::CAFMaker::operator= ( CAFMaker &&  )
delete
void gar::CAFMaker::processIonizationInfo ( rec::TrackIoniz &  ion,
float  ionizeTruncate,
float &  forwardIonVal,
float &  backwardIonVal 
)
private

Definition at line 1416 of file CAFMaker_module.cc.

1417  {
1418 
1419  // NO CALIBRATION SERVICE FOR NOW
1420 
1421  std::vector<std::pair<float,float>> SigData = ion.getFWD_dSigdXs();
1422  forwardIonVal = processOneDirection(SigData, ionizeTruncate);
1423 
1424  SigData = ion.getBAK_dSigdXs();
1425  backwardIonVal = processOneDirection(SigData, ionizeTruncate);
1426 
1427  return;
1428 }
gar::CAFMaker::processOneDirection
float processOneDirection(std::vector< std::pair< float, float >> SigData, float ionizeTruncate)
Definition: CAFMaker_module.cc:1432
float gar::CAFMaker::processOneDirection ( std::vector< std::pair< float, float >>  SigData,
float  ionizeTruncate 
)
private

Definition at line 1432 of file CAFMaker_module.cc.

1432  {
1433 
1434  std::vector<std::pair<float,float>> dEvsX; // Ionization vs distance along track
1435 
1436  // The first hit on the track never had its ionization info stored. Not a problem
1437  // really. Each pair is a hit and the step along the track that ends at the hit
1438  // For the last hit, just take the step from the n-1 hit; don't guess some distance
1439  // to (nonexistant!) n+1 hit. Using pointer arithmetic because you are a real K&R
1440  // C nerd! Except that C++ doesn't know you are such a nerd and if
1441  // SigData.size()==0, then SigData.end()-1 is 0xFFFFFFFFFFFFFFF8.
1442  if (SigData.size()==0) return 0.0;
1443  float distAlongTrack = 0;
1444  std::vector<std::pair<float,float>>::iterator littlebit = SigData.begin();
1445  for (; littlebit<(SigData.end()-1); ++littlebit) {
1446  float dE = std::get<0>(*littlebit);
1447  // tpctrackfit2_module.cc fills the TrackIoniz data product so that
1448  // this quantity is really dL > 0 not dX, a coordinate on the drift axis
1449  float dX = std::get<1>(*littlebit);
1450  distAlongTrack += dX; // But count full step to get hit position on track
1451  // Take dX to be 1/2 the previous + last segment
1452  dX += std::get<1>(*(littlebit+1));
1453  float dEdX = dE/(0.5*dX);
1454 
1455  std::pair pushme = std::make_pair(dEdX,distAlongTrack);
1456  dEvsX.push_back( pushme );
1457  }
1458 
1459  // Get the truncated mean; first sort then take mean
1460  std::sort(dEvsX.begin(),dEvsX.end(), lessThan_byE);
1461 
1462  // Get the dEdX vs length data, truncated.
1463  int goUpTo = ionizeTruncate * dEvsX.size() +0.5;
1464  if (goUpTo > (int)dEvsX.size()) goUpTo = dEvsX.size();
1465  int i = 1; float returnvalue = 0;
1466  littlebit = dEvsX.begin();
1467  for (; littlebit<dEvsX.end(); ++littlebit) {
1468  returnvalue += std::get<0>(*littlebit);
1469  ++i;
1470  if (i>goUpTo) break;
1471  }
1472  returnvalue /= goUpTo;
1473  return returnvalue;
1474 }
pionana::dEdX
double dEdX(double KE, const simb::MCParticle *part)
Definition: TruthAnalyzer_module.cc:67
gar::CAFMaker::lessThan_byE
static bool lessThan_byE(std::pair< float, float > a, std::pair< float, float > b)
Definition: CAFMaker_module.cc:107
std::vector< std::pair< int, float > > gar::CAFMaker::processPIDInfo ( float  p)
private

Definition at line 1539 of file CAFMaker_module.cc.

1539  {
1540 
1541  std::vector<std::string> recopnamelist = {"#pi", "#mu", "p", "K", "d", "e"};
1542  std::vector<int> pdg_charged = {211, 13, 2212, 321, 1000010020, 11};
1543  std::vector< std::pair<int, float> > pid;
1544  pid.resize(6);
1545 
1546  int qclosest = 0;
1547  float dist = 100000000.;
1548  CLHEP::RandFlat FlatRand(fEngine);
1549 
1550  for (int q = 0; q < 501; ++q) {
1551  //Check the title and the reco momentum take only the one that fits
1552  std::string fulltitle = m_pidinterp[q]->GetTitle();
1553  unsigned first = fulltitle.find("=");
1554  unsigned last = fulltitle.find("GeV");
1555  std::string substr = fulltitle.substr(first+1, last - first-1);
1556  float pidinterp_mom = std::atof(substr.c_str());
1557  //calculate the distance between the bin and mom, store the q the closest
1558  float disttemp = std::abs(pidinterp_mom - p);
1559 
1560  if( disttemp < dist ) {
1561  dist = disttemp;
1562  qclosest = q;
1563  }
1564  } // closes the "pidmatrix" loop
1565 
1566  //Compute all the probabities for each type of true to reco
1567  //loop over the columns (true pid)
1568  for (int pidm = 0; pidm < 6; ++pidm) {
1569 
1570  //loop over the columns (true pid)
1571  std::vector< std::pair<float, std::string> > v_prob;
1572 
1573  //loop over the rows (reco pid)
1574  for (int pidr = 0; pidr < 6; ++pidr) {
1575  std::string recoparticlename = m_pidinterp[qclosest]->GetYaxis()->GetBinLabel(pidr+1);
1576  float prob = m_pidinterp[qclosest]->GetBinContent(pidm+1, pidr+1);
1577  //Need to check random number value and prob value then associate the recopdg to the reco prob
1578  v_prob.push_back( std::make_pair(prob, recoparticlename) );
1579  }
1580 
1581  //Compute the pid from it
1582  if (v_prob.size() > 1) {
1583  //Order the vector of prob
1584  std::sort(v_prob.begin(), v_prob.end());
1585  //Throw a random number between 0 and 1
1586  float random_number = FlatRand.fire();
1587  //Make cumulative sum to get the range
1588  std::partial_sum(v_prob.begin(), v_prob.end(), v_prob.begin(), [](const P& _x, const P& _y){return P(_x.first + _y.first, _y.second);});
1589 
1590  for(size_t ivec = 0; ivec < v_prob.size()-1; ivec++) {
1591  if( random_number < v_prob.at(ivec+1).first && random_number >= v_prob.at(ivec).first ) {
1592  pid.push_back( std::make_pair(pdg_charged.at( std::distance( recopnamelist.begin(), std::find(recopnamelist.begin(), recopnamelist.end(), v_prob.at(ivec+1).second) ) ), v_prob.at(ivec+1).first) );
1593  }
1594  }
1595  } else {
1596  pid.push_back( std::make_pair(pdg_charged.at( std::distance( recopnamelist.begin(), std::find(recopnamelist.begin(), recopnamelist.end(), v_prob.at(0).second) ) ), v_prob.at(0).first) );
1597  }
1598  }
1599 
1600  //return a vector of pid and prob
1601  return pid;
1602 }
string
std::string string
Definition: nybbler.cc:12
gar::P
std::pair< float, std::string > P
Definition: anatree_module.cc:75
abs
T abs(T value)
Definition: sparse_vector_test.cc:30
test.p
p
Definition: test.py:223
protoana::distance
double distance(double x1, double y1, double z1, double x2, double y2, double z2)
Definition: truepionXsection_module.cc:217
gar::CAFMaker::m_pidinterp
std::unordered_map< int, TH2F * > m_pidinterp
Definition: CAFMaker_module.cc:189
larg4::dist
constexpr double dist(const TReal *x, const TReal *y, const unsigned int dimension)
Definition: OpFastScintillation.hh:563
gar::CAFMaker::fEngine
CLHEP::HepRandomEngine & fEngine
random engine
Definition: CAFMaker_module.cc:188

Member Data Documentation

cheat::BackTrackerCore* gar::CAFMaker::BackTrack
private

Definition at line 95 of file CAFMaker_module.cc.

std::string gar::CAFMaker::fCaloHitLabel
private

module label for reco calo hits rec::CaloHit

Definition at line 138 of file CAFMaker_module.cc.

std::string gar::CAFMaker::fClusterLabel
private

module label for calo clusters rec::Cluster

Definition at line 141 of file CAFMaker_module.cc.

std::string gar::CAFMaker::fClusterMuIDLabel
private

module label for calo clusters rec::Cluster in MuID

Definition at line 142 of file CAFMaker_module.cc.

std::string gar::CAFMaker::fECALAssnLabel
private

module label for track-clusters associations

Definition at line 144 of file CAFMaker_module.cc.

CLHEP::HepRandomEngine& gar::CAFMaker::fEngine
private

random engine

Definition at line 188 of file CAFMaker_module.cc.

std::string gar::CAFMaker::fGeantLabel
private

module label for geant4 simulated hits

Definition at line 124 of file CAFMaker_module.cc.

std::vector<std::string> gar::CAFMaker::fGeneratorLabels
private

Definition at line 122 of file CAFMaker_module.cc.

std::vector<std::string> gar::CAFMaker::fGENIEGeneratorLabels
private

Definition at line 123 of file CAFMaker_module.cc.

const geo::GeometryCore* gar::CAFMaker::fGeo
private

pointer to the geometry

Definition at line 182 of file CAFMaker_module.cc.

std::string gar::CAFMaker::fHitLabel
private

module label for reco TPC hits rec::Hit

Definition at line 128 of file CAFMaker_module.cc.

std::string gar::CAFMaker::fInstanceLabelCalo
private

Instance name for ECAL.

Definition at line 125 of file CAFMaker_module.cc.

std::string gar::CAFMaker::fInstanceLabelMuID
private

Instance name for MuID.

Definition at line 126 of file CAFMaker_module.cc.

float gar::CAFMaker::fIonizTruncate
private

Default=1.00;.

Definition at line 169 of file CAFMaker_module.cc.

float gar::CAFMaker::fMatchMCPtoVertDist
private

MCParticle to MC vertex match Default=roundoff.

Definition at line 153 of file CAFMaker_module.cc.

std::string gar::CAFMaker::fMuIDHitLabel
private

Definition at line 139 of file CAFMaker_module.cc.

std::string gar::CAFMaker::fPFLabel
private

module label for reco particles rec::PFParticle

Definition at line 143 of file CAFMaker_module.cc.

TH1* gar::CAFMaker::fPOTHist
private

Definition at line 174 of file CAFMaker_module.cc.

std::string gar::CAFMaker::fPOTSummaryLabel
private

Definition at line 146 of file CAFMaker_module.cc.

std::string gar::CAFMaker::fRawCaloHitLabel
private

module label for digitized calo hits raw::CaloRawDigit

Definition at line 135 of file CAFMaker_module.cc.

std::string gar::CAFMaker::fRawMuIDHitLabel
private

Definition at line 136 of file CAFMaker_module.cc.

Float_t gar::CAFMaker::fTPC_X
private

center of TPC stored as per-event & compressed by root

Definition at line 176 of file CAFMaker_module.cc.

Float_t gar::CAFMaker::fTPC_Y
private

Definition at line 177 of file CAFMaker_module.cc.

Float_t gar::CAFMaker::fTPC_Z
private

Definition at line 178 of file CAFMaker_module.cc.

std::string gar::CAFMaker::fTPCClusterLabel
private

module label for TPC Clusters rec::TPCCluster

Definition at line 129 of file CAFMaker_module.cc.

std::string gar::CAFMaker::fTrackLabel
private

module label for TPC Tracks rec:Track

Definition at line 130 of file CAFMaker_module.cc.

std::string gar::CAFMaker::fTrackTrajectoryLabel
private

module label for TPC Track Trajectories rec:TrackTrajectory

Definition at line 131 of file CAFMaker_module.cc.

TTree* gar::CAFMaker::fTree
private

Definition at line 173 of file CAFMaker_module.cc.

std::string gar::CAFMaker::fVeeLabel
private

module label for conversion/decay vertexes rec:Vee

Definition at line 133 of file CAFMaker_module.cc.

std::string gar::CAFMaker::fVertexLabel
private

module label for vertexes rec:Vertex

Definition at line 132 of file CAFMaker_module.cc.

bool gar::CAFMaker::fWriteCaloClusters
private

Write ECAL clusters. Default=true.

Definition at line 165 of file CAFMaker_module.cc.

bool gar::CAFMaker::fWriteCaloDigits
private

Raw digits for calorimetry. Default=false.

Definition at line 163 of file CAFMaker_module.cc.

bool gar::CAFMaker::fWriteCaloHits
private

Write ECAL hits. Default=true.

Definition at line 164 of file CAFMaker_module.cc.

bool gar::CAFMaker::fWriteCohInfo
private

MC level t for coherent pi+. Default=false.

Definition at line 170 of file CAFMaker_module.cc.

bool gar::CAFMaker::fWriteHits
private

Write info about TPC Hits Default=false.

Definition at line 155 of file CAFMaker_module.cc.

bool gar::CAFMaker::fWriteMatchedTracks
private

Write ECAL-track Assns Default=true.

Definition at line 166 of file CAFMaker_module.cc.

bool gar::CAFMaker::fWriteMCCaloInfo
private

Write MC info for calorimeter Default=true.

Definition at line 152 of file CAFMaker_module.cc.

bool gar::CAFMaker::fWriteMCinfo
private

Info from MCTruth, GTruth Default=true.

Definition at line 149 of file CAFMaker_module.cc.

bool gar::CAFMaker::fWriteMCPTrajectory
private

Write MCP Trajectory Default=true.

Definition at line 150 of file CAFMaker_module.cc.

bool gar::CAFMaker::fWriteMCPTrajMomenta
private

Write Momenta associated with MCP Trajectory Default=false.

Definition at line 151 of file CAFMaker_module.cc.

bool gar::CAFMaker::fWriteMuID
private

Definition at line 162 of file CAFMaker_module.cc.

bool gar::CAFMaker::fWriteTPCClusters
private

Write TPCClusters info Default=true.

Definition at line 156 of file CAFMaker_module.cc.

bool gar::CAFMaker::fWriteTracks
private

Start/end X, P for tracks Default=true.

Definition at line 157 of file CAFMaker_module.cc.

bool gar::CAFMaker::fWriteTrackTrajectories
private

Point traj of reco tracks Default=false.

Definition at line 158 of file CAFMaker_module.cc.

bool gar::CAFMaker::fWriteVees
private

Reco vees & their tracks Default=true.

Definition at line 160 of file CAFMaker_module.cc.

bool gar::CAFMaker::fWriteVertices
private

Reco vertexes & their tracks Default=true.

Definition at line 159 of file CAFMaker_module.cc.

float gar::CAFMaker::ItsInTulsa[3]
private

Definition at line 117 of file CAFMaker_module.cc.

std::unordered_map<int, TH2F*> gar::CAFMaker::m_pidinterp
private

Definition at line 189 of file CAFMaker_module.cc.

float gar::CAFMaker::rTPC
private

Definition at line 119 of file CAFMaker_module.cc.

std::unordered_map<TrkId, Int_t> gar::CAFMaker::TrackIdToIndex
private

Definition at line 185 of file CAFMaker_module.cc.

float gar::CAFMaker::xTPC
private

Definition at line 118 of file CAFMaker_module.cc.

The documentation for this class was generated from the following file:
Generated by   doxygen 1.8.11