LArSoft  v09_90_00
Liquid Argon Software toolkit - https://larsoft.org/
EventWritingAlgorithm.cc
Go to the documentation of this file.
1 
9 #include "Pandora/AlgorithmHeaders.h"
10 
11 #include "Persistency/BinaryFileWriter.h"
12 #include "Persistency/XmlFileWriter.h"
13 
16 
19 
21 
22 using namespace pandora;
23 
24 namespace lar_content
25 {
26 
27 EventWritingAlgorithm::EventWritingAlgorithm() :
28  m_geometryFileType(UNKNOWN_FILE_TYPE),
29  m_eventFileType(UNKNOWN_FILE_TYPE),
30  m_pEventFileWriter(nullptr),
31  m_pGeometryFileWriter(nullptr),
32  m_shouldWriteGeometry(false),
33  m_writtenGeometry(false),
34  m_shouldWriteEvents(true),
35  m_shouldWriteMCRelationships(true),
36  m_shouldWriteTrackRelationships(true),
37  m_shouldOverwriteEventFile(false),
38  m_shouldOverwriteGeometryFile(false),
39  m_useLArCaloHits(true),
40  m_larCaloHitVersion(1),
41  m_useLArMCParticles(true),
42  m_shouldFilterByNuanceCode(false),
43  m_filterNuanceCode(0),
44  m_shouldFilterByMCParticles(false),
45  m_neutrinoInducedOnly(false),
46  m_matchingMinPrimaryHits(15),
47  m_nNonNeutrons(0),
48  m_nMuons(0),
49  m_nElectrons(0),
50  m_nProtons(0),
51  m_nPhotons(0),
52  m_nChargedPions(0),
53  m_shouldFilterByNeutrinoVertex(false),
54  m_detectorHalfLengthX(-1.f),
55  m_detectorHalfLengthY(-1.f),
56  m_detectorHalfLengthZ(-1.f),
57  m_coordinateOffsetX(0.f),
58  m_coordinateOffsetY(0.f),
59  m_coordinateOffsetZ(0.f),
60  m_selectedBorderX(-1.f),
61  m_selectedBorderY(-1.f),
62  m_selectedBorderZ(-1.f)
63 {
64 }
65 
66 //------------------------------------------------------------------------------------------------------------------------------------------
67 
69 {
70  delete m_pEventFileWriter;
71  delete m_pGeometryFileWriter;
72 }
73 
74 //------------------------------------------------------------------------------------------------------------------------------------------
75 
77 {
79  {
80  const FileMode fileMode(m_shouldOverwriteGeometryFile ? OVERWRITE : APPEND);
81 
82  if (BINARY == m_geometryFileType)
83  {
84  m_pGeometryFileWriter = new BinaryFileWriter(this->GetPandora(), m_geometryFileName, fileMode);
85  }
86  else if (XML == m_geometryFileType)
87  {
88  m_pGeometryFileWriter = new XmlFileWriter(this->GetPandora(), m_geometryFileName, fileMode);
89  }
90  else
91  {
92  return STATUS_CODE_FAILURE;
93  }
94  }
95 
97  {
98  const FileMode fileMode(m_shouldOverwriteEventFile ? OVERWRITE : APPEND);
99 
100  if (BINARY == m_eventFileType)
101  {
102  m_pEventFileWriter = new BinaryFileWriter(this->GetPandora(), m_eventFileName, fileMode);
103  }
104  else if (XML == m_eventFileType)
105  {
106  m_pEventFileWriter = new XmlFileWriter(this->GetPandora(), m_eventFileName, fileMode);
107  }
108  else
109  {
110  return STATUS_CODE_FAILURE;
111  }
112 
113  if (m_useLArCaloHits)
115 
117  m_pEventFileWriter->SetFactory(new LArMCParticleFactory);
118  }
119 
120  return STATUS_CODE_SUCCESS;
121 }
122 
123 //------------------------------------------------------------------------------------------------------------------------------------------
124 
126 {
127  // ATTN Should complete geometry creation in LArSoft begin job, but some channel status service functionality unavailable at that point
129  {
130  PANDORA_RETURN_RESULT_IF(STATUS_CODE_SUCCESS, !=, m_pGeometryFileWriter->WriteGeometry());
131  m_writtenGeometry = true;
132  }
133 
134  bool matchNuanceCode(!m_shouldFilterByNuanceCode || this->PassNuanceCodeFilter());
135  bool matchParticles(!m_shouldFilterByMCParticles || this->PassMCParticleFilter());
136  bool matchNeutrinoVertexPosition(!m_shouldFilterByNeutrinoVertex || this->PassNeutrinoVertexFilter());
137 
138  if (matchNuanceCode && matchParticles && matchNeutrinoVertexPosition && m_pEventFileWriter && m_shouldWriteEvents)
139  {
140  const CaloHitList *pCaloHitList = nullptr;
141  PANDORA_RETURN_RESULT_IF(STATUS_CODE_SUCCESS, !=, PandoraContentApi::GetCurrentList(*this, pCaloHitList));
142 
143  const TrackList *pTrackList = nullptr;
144  PANDORA_RETURN_RESULT_IF(STATUS_CODE_SUCCESS, !=, PandoraContentApi::GetCurrentList(*this, pTrackList));
145 
146  const MCParticleList *pMCParticleList = nullptr;
147  PANDORA_RETURN_RESULT_IF(STATUS_CODE_SUCCESS, !=, PandoraContentApi::GetCurrentList(*this, pMCParticleList));
148 
149  PANDORA_RETURN_RESULT_IF(STATUS_CODE_SUCCESS, !=,
150  m_pEventFileWriter->WriteEvent(*pCaloHitList, *pTrackList, *pMCParticleList, m_shouldWriteMCRelationships, m_shouldWriteTrackRelationships));
151  }
152 
153  return STATUS_CODE_SUCCESS;
154 }
155 
156 //------------------------------------------------------------------------------------------------------------------------------------------
157 
159 {
160  const MCParticleList *pMCParticleList = nullptr;
161  PANDORA_RETURN_RESULT_IF(STATUS_CODE_SUCCESS, !=, PandoraContentApi::GetCurrentList(*this, pMCParticleList));
162 
163  MCParticleVector mcNeutrinoList;
164  LArMCParticleHelper::GetTrueNeutrinos(pMCParticleList, mcNeutrinoList);
165 
166  for (const MCParticle *const pMCNeutrino : mcNeutrinoList)
167  {
168  const LArMCParticle *const pLArMCNeutrino = dynamic_cast<const LArMCParticle *>(pMCNeutrino);
169 
170  if (pLArMCNeutrino && (pLArMCNeutrino->GetNuanceCode() == m_filterNuanceCode))
171  return true;
172  }
173 
174  return false;
175 }
176 
177 //------------------------------------------------------------------------------------------------------------------------------------------
178 
180 {
181  const MCParticleList *pMCParticleList = nullptr;
182  PANDORA_RETURN_RESULT_IF(STATUS_CODE_SUCCESS, !=, PandoraContentApi::GetCurrentList(*this, pMCParticleList));
183 
184  const CaloHitList *pCaloHitList = nullptr;
185  PANDORA_RETURN_RESULT_IF(STATUS_CODE_SUCCESS, !=, PandoraContentApi::GetCurrentList(*this, pCaloHitList));
186 
188  LArMCParticleHelper::MCContributionMap mcParticlesToGoodHitsMap;
190  pMCParticleList, pCaloHitList, parameters, LArMCParticleHelper::IsBeamNeutrinoFinalState, mcParticlesToGoodHitsMap);
191 
193  {
195  pMCParticleList, pCaloHitList, parameters, LArMCParticleHelper::IsBeamParticle, mcParticlesToGoodHitsMap);
197  pMCParticleList, pCaloHitList, parameters, LArMCParticleHelper::IsCosmicRay, mcParticlesToGoodHitsMap);
198  }
199 
200  unsigned int nNonNeutrons(0), nMuons(0), nElectrons(0), nProtons(0), nPhotons(0), nChargedPions(0);
201 
202  MCParticleVector mcPrimaryVector;
203  for (const auto &mapEntry : mcParticlesToGoodHitsMap)
204  mcPrimaryVector.push_back(mapEntry.first);
205  std::sort(mcPrimaryVector.begin(), mcPrimaryVector.end(), LArMCParticleHelper::SortByMomentum);
206 
207  for (const MCParticle *const pMCPrimary : mcPrimaryVector)
208  {
209  const unsigned int particleId(std::abs(pMCPrimary->GetParticleId()));
210  if (NEUTRON != particleId)
211  ++nNonNeutrons;
212 
213  if (MU_MINUS == particleId)
214  ++nMuons;
215  else if (E_MINUS == particleId)
216  ++nElectrons;
217  else if (PROTON == particleId)
218  ++nProtons;
219  else if (PHOTON == particleId)
220  ++nPhotons;
221  else if (PI_PLUS == particleId)
222  ++nChargedPions;
223  }
224 
225  if ((nNonNeutrons == m_nNonNeutrons) && (nMuons == m_nMuons) && (nElectrons == m_nElectrons) && (nProtons == m_nProtons) &&
226  (nPhotons == m_nPhotons) && (nChargedPions == m_nChargedPions))
227  {
228  return true;
229  }
230 
231  return false;
232 }
233 
234 //------------------------------------------------------------------------------------------------------------------------------------------
235 
237 {
238  const MCParticleList *pMCParticleList = nullptr;
239  PANDORA_RETURN_RESULT_IF(STATUS_CODE_SUCCESS, !=, PandoraContentApi::GetCurrentList(*this, pMCParticleList));
240 
241  MCParticleVector mcNeutrinoVector;
242  LArMCParticleHelper::GetTrueNeutrinos(pMCParticleList, mcNeutrinoVector);
243 
244  for (const MCParticle *const pMCNeutrino : mcNeutrinoVector)
245  {
246  const CartesianVector &neutrinoInteractionVertex(pMCNeutrino->GetEndpoint());
247 
248  if ((neutrinoInteractionVertex.GetX() < (m_detectorHalfLengthX - m_coordinateOffsetX - m_selectedBorderX)) &&
249  (neutrinoInteractionVertex.GetX() > (-m_coordinateOffsetX + m_selectedBorderX)) &&
250  (neutrinoInteractionVertex.GetY() < (m_detectorHalfLengthY - m_coordinateOffsetY - m_selectedBorderY)) &&
251  (neutrinoInteractionVertex.GetY() > (-m_coordinateOffsetY + m_selectedBorderY)) &&
252  (neutrinoInteractionVertex.GetZ() < (m_detectorHalfLengthZ - m_coordinateOffsetZ - m_selectedBorderZ)) &&
253  (neutrinoInteractionVertex.GetZ() > (-m_coordinateOffsetZ + m_selectedBorderZ)))
254  {
255  return true;
256  }
257  }
258 
259  return false;
260 }
261 
262 //------------------------------------------------------------------------------------------------------------------------------------------
263 
264 StatusCode EventWritingAlgorithm::ReadSettings(const TiXmlHandle xmlHandle)
265 {
266  PANDORA_RETURN_RESULT_IF_AND_IF(
267  STATUS_CODE_SUCCESS, STATUS_CODE_NOT_FOUND, !=, XmlHelper::ReadValue(xmlHandle, "ShouldWriteGeometry", m_shouldWriteGeometry));
268 
270  {
271  PANDORA_RETURN_RESULT_IF(STATUS_CODE_SUCCESS, !=, XmlHelper::ReadValue(xmlHandle, "GeometryFileName", m_geometryFileName));
272 
273  std::string fileExtension(m_geometryFileName.substr(m_geometryFileName.find_last_of(".")));
274  std::transform(fileExtension.begin(), fileExtension.end(), fileExtension.begin(), ::tolower);
275 
276  if (std::string(".xml") == fileExtension)
277  {
278  m_geometryFileType = XML;
279  }
280  else if (std::string(".pndr") == fileExtension)
281  {
282  m_geometryFileType = BINARY;
283  }
284  else
285  {
286  std::cout << "EventReadingAlgorithm: Unknown geometry file type specified " << std::endl;
287  return STATUS_CODE_INVALID_PARAMETER;
288  }
289  }
290 
291  PANDORA_RETURN_RESULT_IF_AND_IF(
292  STATUS_CODE_SUCCESS, STATUS_CODE_NOT_FOUND, !=, XmlHelper::ReadValue(xmlHandle, "ShouldWriteEvents", m_shouldWriteEvents));
293 
295  {
296  PANDORA_RETURN_RESULT_IF(STATUS_CODE_SUCCESS, !=, XmlHelper::ReadValue(xmlHandle, "EventFileName", m_eventFileName));
297 
298  std::string fileExtension(m_eventFileName.substr(m_eventFileName.find_last_of(".")));
299  std::transform(fileExtension.begin(), fileExtension.end(), fileExtension.begin(), ::tolower);
300 
301  if (std::string(".xml") == fileExtension)
302  {
303  m_eventFileType = XML;
304  }
305  else if (std::string(".pndr") == fileExtension)
306  {
307  m_eventFileType = BINARY;
308  }
309  else
310  {
311  std::cout << "EventReadingAlgorithm: Unknown event file type specified " << std::endl;
312  return STATUS_CODE_INVALID_PARAMETER;
313  }
314  }
315 
316  PANDORA_RETURN_RESULT_IF_AND_IF(STATUS_CODE_SUCCESS, STATUS_CODE_NOT_FOUND, !=,
317  XmlHelper::ReadValue(xmlHandle, "ShouldWriteMCRelationships", m_shouldWriteMCRelationships));
318 
319  PANDORA_RETURN_RESULT_IF_AND_IF(STATUS_CODE_SUCCESS, STATUS_CODE_NOT_FOUND, !=,
320  XmlHelper::ReadValue(xmlHandle, "ShouldWriteTrackRelationships", m_shouldWriteTrackRelationships));
321 
322  PANDORA_RETURN_RESULT_IF_AND_IF(STATUS_CODE_SUCCESS, STATUS_CODE_NOT_FOUND, !=,
323  XmlHelper::ReadValue(xmlHandle, "ShouldOverwriteEventFile", m_shouldOverwriteEventFile));
324 
325  PANDORA_RETURN_RESULT_IF_AND_IF(STATUS_CODE_SUCCESS, STATUS_CODE_NOT_FOUND, !=,
326  XmlHelper::ReadValue(xmlHandle, "ShouldOverwriteGeometryFile", m_shouldOverwriteGeometryFile));
327 
328  PANDORA_RETURN_RESULT_IF_AND_IF(STATUS_CODE_SUCCESS, STATUS_CODE_NOT_FOUND, !=, XmlHelper::ReadValue(xmlHandle, "UseLArCaloHits", m_useLArCaloHits));
329 
330  PANDORA_RETURN_RESULT_IF_AND_IF(
331  STATUS_CODE_SUCCESS, STATUS_CODE_NOT_FOUND, !=, XmlHelper::ReadValue(xmlHandle, "LArCaloHitVersion", m_larCaloHitVersion));
332 
333  PANDORA_RETURN_RESULT_IF_AND_IF(
334  STATUS_CODE_SUCCESS, STATUS_CODE_NOT_FOUND, !=, XmlHelper::ReadValue(xmlHandle, "UseLArMCParticles", m_useLArMCParticles));
335 
336  PANDORA_RETURN_RESULT_IF_AND_IF(STATUS_CODE_SUCCESS, STATUS_CODE_NOT_FOUND, !=,
337  XmlHelper::ReadValue(xmlHandle, "ShouldFilterByNuanceCode", m_shouldFilterByNuanceCode));
338 
339  PANDORA_RETURN_RESULT_IF_AND_IF(STATUS_CODE_SUCCESS, STATUS_CODE_NOT_FOUND, !=,
340  XmlHelper::ReadValue(xmlHandle, "ShouldFilterByMCParticles", m_shouldFilterByMCParticles));
341 
342  PANDORA_RETURN_RESULT_IF_AND_IF(STATUS_CODE_SUCCESS, STATUS_CODE_NOT_FOUND, !=,
343  XmlHelper::ReadValue(xmlHandle, "ShouldFilterByNeutrinoVertex", m_shouldFilterByNeutrinoVertex));
344 
346  {
347  PANDORA_RETURN_RESULT_IF(STATUS_CODE_SUCCESS, !=, XmlHelper::ReadValue(xmlHandle, "FilterNuanceCode", m_filterNuanceCode));
348  }
349 
351  {
352  PANDORA_RETURN_RESULT_IF(STATUS_CODE_SUCCESS, !=, XmlHelper::ReadValue(xmlHandle, "NeutrinoInducedOnly", m_neutrinoInducedOnly));
353  PANDORA_RETURN_RESULT_IF(STATUS_CODE_SUCCESS, !=, XmlHelper::ReadValue(xmlHandle, "MatchingMinPrimaryHits", m_matchingMinPrimaryHits));
354  PANDORA_RETURN_RESULT_IF(STATUS_CODE_SUCCESS, !=, XmlHelper::ReadValue(xmlHandle, "NNonNeutrons", m_nNonNeutrons));
355  PANDORA_RETURN_RESULT_IF(STATUS_CODE_SUCCESS, !=, XmlHelper::ReadValue(xmlHandle, "NMuons", m_nMuons));
356  PANDORA_RETURN_RESULT_IF(STATUS_CODE_SUCCESS, !=, XmlHelper::ReadValue(xmlHandle, "NElectrons", m_nElectrons));
357  PANDORA_RETURN_RESULT_IF(STATUS_CODE_SUCCESS, !=, XmlHelper::ReadValue(xmlHandle, "NProtons", m_nProtons));
358  PANDORA_RETURN_RESULT_IF(STATUS_CODE_SUCCESS, !=, XmlHelper::ReadValue(xmlHandle, "NPhotons", m_nPhotons));
359  PANDORA_RETURN_RESULT_IF(STATUS_CODE_SUCCESS, !=, XmlHelper::ReadValue(xmlHandle, "NChargedPions", m_nChargedPions));
360  }
361 
363  {
364  PANDORA_RETURN_RESULT_IF(STATUS_CODE_SUCCESS, !=, XmlHelper::ReadValue(xmlHandle, "DetectorHalfLengthX", m_detectorHalfLengthX));
365  PANDORA_RETURN_RESULT_IF(STATUS_CODE_SUCCESS, !=, XmlHelper::ReadValue(xmlHandle, "DetectorHalfLengthY", m_detectorHalfLengthY));
366  PANDORA_RETURN_RESULT_IF(STATUS_CODE_SUCCESS, !=, XmlHelper::ReadValue(xmlHandle, "DetectorHalfLengthZ", m_detectorHalfLengthZ));
367  PANDORA_RETURN_RESULT_IF(STATUS_CODE_SUCCESS, !=, XmlHelper::ReadValue(xmlHandle, "CoordinateOffsetX", m_coordinateOffsetX));
368  PANDORA_RETURN_RESULT_IF(STATUS_CODE_SUCCESS, !=, XmlHelper::ReadValue(xmlHandle, "CoordinateOffsetY", m_coordinateOffsetY));
369  PANDORA_RETURN_RESULT_IF(STATUS_CODE_SUCCESS, !=, XmlHelper::ReadValue(xmlHandle, "CoordinateOffsetZ", m_coordinateOffsetZ));
370  PANDORA_RETURN_RESULT_IF(STATUS_CODE_SUCCESS, !=, XmlHelper::ReadValue(xmlHandle, "SelectedBorderX", m_selectedBorderX));
371  PANDORA_RETURN_RESULT_IF(STATUS_CODE_SUCCESS, !=, XmlHelper::ReadValue(xmlHandle, "SelectedBorderY", m_selectedBorderY));
372  PANDORA_RETURN_RESULT_IF(STATUS_CODE_SUCCESS, !=, XmlHelper::ReadValue(xmlHandle, "SelectedBorderZ", m_selectedBorderZ));
373  }
374 
375  return STATUS_CODE_SUCCESS;
376 }
377 
378 } // namespace lar_content
float m_detectorHalfLengthZ
Half length of detector in z dimension.
unsigned int m_matchingMinPrimaryHits
The minimum number of mc primary hits used in matching scheme.
bool m_shouldWriteMCRelationships
Whether to write mc relationship information to the events file.
pandora::FileWriter * m_pEventFileWriter
Address of the event file writer.
Header file for the event writing algorithm class.
std::unordered_map< const pandora::MCParticle *, pandora::CaloHitList > MCContributionMap
bool m_shouldFilterByMCParticles
Whether to filter output by mc particle constituents.
bool m_shouldWriteEvents
Whether to write events to a specified file.
Header file for the lar calo hit class.
bool PassNeutrinoVertexFilter() const
Whether current event passes neutrino vertex position filter (e.g. fiducial volume cut) ...
bool m_useLArCaloHits
Whether to write lar calo hits, or standard pandora calo hits.
float m_coordinateOffsetY
Origin offset (from detector corner) in y dimension.
int GetNuanceCode() const
Get the nuance code.
bool m_neutrinoInducedOnly
Whether to consider only mc particles that were neutrino induced.
unsigned int m_larCaloHitVersion
LArCaloHit version for LArCaloHitFactory.
unsigned int m_nPhotons
The requested number of mc primaries that are photons.
unsigned int m_nProtons
The requested number of mc primaries that are protons.
constexpr auto abs(T v)
Returns the absolute value of the argument.
float m_selectedBorderY
Required distance from detector edge in y dimension.
float m_detectorHalfLengthY
Half length of detector in y dimension.
float m_selectedBorderZ
Required distance from detector edge in z dimension.
pandora::StatusCode ReadSettings(const pandora::TiXmlHandle xmlHandle)
std::string m_eventFileName
Name of the output event file.
Header file for the lar monitoring helper helper class.
pandora::FileType m_geometryFileType
The geometry file type.
unsigned int m_nChargedPions
The requested number of mc primaries that are charged pions.
TFile f
Definition: plotHisto.C:6
bool PassNuanceCodeFilter() const
Whether current event passes nuance code filter.
bool m_shouldWriteTrackRelationships
Whether to write track relationship information to the events file.
bool m_useLArMCParticles
Whether to write lar mc particles, or standard pandora mc particles.
static bool SortByMomentum(const pandora::MCParticle *const pLhs, const pandora::MCParticle *const pRhs)
Sort mc particles by their momentum.
LAr mc particle class.
Definition: LArMCParticle.h:94
static bool IsCosmicRay(const pandora::MCParticle *const pMCParticle)
Return true if passed a primary cosmic ray MCParticle.
pandora::FileWriter * m_pGeometryFileWriter
Address of the geometry file writer.
static void SelectReconstructableMCParticles(const pandora::MCParticleList *pMCParticleList, const pandora::CaloHitList *pCaloHitList, const PrimaryParameters &parameters, std::function< bool(const pandora::MCParticle *const)> fCriteria, MCContributionMap &selectedMCParticlesToHitsMap)
Select target, reconstructable mc particles that match given criteria.
unsigned int m_nNonNeutrons
The requested number of mc primaries that are not neutrons.
Header file for the lar monte carlo particle helper helper class.
static bool IsBeamParticle(const pandora::MCParticle *const pMCParticle)
Returns true if passed a primary beam MCParticle.
float m_selectedBorderX
Required distance from detector edge in x dimension.
bool m_shouldFilterByNeutrinoVertex
Whether to filter output by neutrino vertex position (e.g. fiducial volume cut)
std::vector< art::Ptr< simb::MCParticle > > MCParticleVector
bool PassMCParticleFilter() const
Whether current event passes mc particle constituent filter.
int m_filterNuanceCode
The filter nuance code (required if specify filter by nuance code)
float m_detectorHalfLengthX
Half length of detector in x dimension.
Header file for the lar mc particle class.
unsigned int m_nElectrons
The requested number of mc primaries that are electrons.
LArMCParticleFactory responsible for object creation.
std::string m_geometryFileName
Name of the output geometry file.
static void GetTrueNeutrinos(const pandora::MCParticleList *const pMCParticleList, pandora::MCParticleVector &trueNeutrinos)
Get neutrino MC particles from an input MC particle list.
pandora::FileType m_eventFileType
The event file type.
bool m_shouldFilterByNuanceCode
Whether to filter output by nuance code.
bool m_shouldOverwriteEventFile
Whether to overwrite existing event file with specified name, or append.
unsigned int m_nMuons
The requested number of mc primaries that are muons.
LArCaloHitFactory responsible for object creation.
Definition: LArCaloHit.h:110
bool m_shouldWriteGeometry
Whether to write geometry to a specified file.
bool m_shouldOverwriteGeometryFile
Whether to overwrite existing geometry file with specified name, or append.
float m_coordinateOffsetZ
Origin offset (from detector corner) in z dimension.
static bool IsBeamNeutrinoFinalState(const pandora::MCParticle *const pMCParticle)
Returns true if passed a primary neutrino final state MCParticle.
bool m_writtenGeometry
Whether geometry has been written.
float m_coordinateOffsetX
Origin offset (from detector corner) in x dimension.