GCNZlibMakerProtoDUNE_module.cc
Go to the documentation of this file.
1 ////////////////////////////////////////////////////////////////////////
2 // \file GCNZlibMakerProtoDUNE_module.cc
3 // \brief Analyzer module for creating CVN gzip file objects
4 // \author Jeremy Hewes - jhewes15@fnal.gov
5 // Saul Alonso-Monsalve - saul.alonso.monsalve@cern.ch
6 // - wrote the zlib code used in this module
7 ////////////////////////////////////////////////////////////////////////
8 
9 // C/C++ includes
10 #include <iostream>
11 #include <sstream>
12 #include "boost/filesystem.hpp"
13 
14 // Framework includes
19 
21 // Data products
24 
25 // CVN includes
28 
29 // Compression
30 #include "zlib.h"
31 
32 namespace fs = boost::filesystem;
33 
34 namespace cvn {
35 
37  public:
38 
39  explicit GCNZlibMakerProtoDUNE(fhicl::ParameterSet const& pset);
41 
42  void beginJob() override;
43  void analyze(const art::Event& evt) override;
44  void reconfigure(const fhicl::ParameterSet& pset);
45 
46  private:
47 
50  unsigned int fTopologyHitsCut;
51 
53 
55 
56  };
57 
58  //......................................................................
60  : EDAnalyzer(pset)
61  {
62  this->reconfigure(pset);
63  }
64 
65  //......................................................................
67  { }
68 
69  //......................................................................
71  {
72  fOutputDir = pset.get<std::string>("OutputDir", "");
73  fGraphLabel = pset.get<std::string>("GraphLabel");
74  fTopologyHitsCut = pset.get<unsigned int>("TopologyHitsCut");
75 
76  fLArG4ModuleLabel = pset.get<std::string>("LArG4ModuleLabel");
77  }
78 
79  //......................................................................
81  {
82  if (fOutputDir != "")
84 
85  else
86  out_dir = ".";
87 
88  // Throw an error if the specified output directory doesn't exist
89  if (!fs::exists(out_dir))
91  << "Output directory " << out_dir << " does not exist!" << std::endl;
92 
93  // std::cout << "Writing files to output directory " << out_dir << std::endl;
94  }
95 
96  //......................................................................
98  {
99 
100  std::cout << "GCNZlibMakerProtoDUNE: looking for graphs with label " << fGraphLabel << std::endl;
101 
102  // Get the graphs
103  std::vector<art::Ptr<cvn::GCNGraph>> graphs;
104  auto h_graphs = evt.getHandle<std::vector<cvn::GCNGraph>>(fGraphLabel);
105  if (h_graphs)
106  art::fill_ptr_vector(graphs, h_graphs);
107 
108  // If no graphs, quit
109  if (graphs.size() == 0) return;
110 
111  // MC information
113 
115  unsigned short beamParticleInteraction;
116 
117  // We only have one beam particle->called "primary", so find it. It should be the
118  // first particle-> but best not to assume
119  float beamParticleEnergy = 0;
120  TVector3 beamParticleVtx; // This is the interaction vertex
121  int beamParticlePDG = 0;
122 
123 // bool gotPrimary = false;
124  for(auto const m: piService->ParticleList()){
125  const simb::MCParticle* particle = m.second;
126  if(particle->Process().compare("primary")==0){
127  beamParticleEnergy = particle->E();
128  beamParticleVtx.SetXYZ(particle->EndX(),particle->EndY(),particle->EndZ());
129  beamParticlePDG = particle->PdgCode();
130  beamParticleInteraction = labels.GetProtoDUNEBeamInteractionType(*particle);
131 // gotPrimary = true;
132  break;
133  }
134  }
135 
136 // if(!gotPrimary) return;
137  unsigned int counter = 0;
138  for(auto const graph : graphs){
139  // If the graph has no nodes then give up
140  if(graph->GetNumberOfNodes() == 0) return;
141 
142 // std::cout << "GCNZlibMakerProtoDUNE: found graph with " << graph->GetNumberOfNodes() << " nodes" << std::endl;
143 
144 
145  // Now write the zlib file using this information
146  // We need to extract all of the information into a single vector to write
147  // into the compressed file format
148  const std::vector<float> vectorToWrite = graph->ConvertGraphToVector();
149 
150  ulong src_len = vectorToWrite.size() *sizeof(float);
151  ulong dest_len = compressBound(src_len); // calculate size of the compressed data
152  char* ostream = (char *) malloc(dest_len); // allocate memory for the compressed data
153 
154  int res = compress((Bytef *) ostream, &dest_len, (Bytef *) &vectorToWrite[0], src_len);
155 
156  // Buffer error
157  if (res == Z_BUF_ERROR)
158  std::cout << "Buffer too small!" << std::endl;
159  // Memory error
160  else if (res == Z_MEM_ERROR)
161  std::cout << "Not enough memory for compression!" << std::endl;
162  // Compression ok
163  else {
164 
165  // Create output files
166  std::stringstream image_file_name;
167  image_file_name << out_dir << "/gcn_event_" << evt.event() << "_" << counter << ".gz";
168  std::stringstream info_file_name;
169  info_file_name << out_dir << "/gcn_event_" << evt.event() << "_" << counter << ".info";
170 
171  std::ofstream image_file (image_file_name.str(), std::ofstream::binary);
172  std::ofstream info_file (info_file_name.str());
173 
174  if(image_file.is_open() && info_file.is_open()) {
175 
176  // Write the graph to the file and close it
177  image_file.write(ostream, dest_len);
178  image_file.close(); // close file
179 
180  // Write the auxillary information to the text file
181  info_file << beamParticleVtx.X() << std::endl;
182  info_file << beamParticleVtx.Y() << std::endl;
183  info_file << beamParticleVtx.Z() << std::endl;
184  info_file << beamParticleEnergy << std::endl;
185  info_file << beamParticleInteraction << std::endl; // Interaction type first
186  info_file << beamParticlePDG << std::endl;
187 
188  // Number of nodes and node features is needed for unpacking
189  info_file << graph->GetNumberOfNodes() << std::endl;
190  info_file << graph->GetNode(0).GetNumberOfFeatures() << std::endl;
191 
192  info_file.close(); // close file
193  }
194  else {
195 
196  if (image_file.is_open())
197  image_file.close();
198  else
200  << "Unable to open file " << image_file_name.str() << "!" << std::endl;
201 
202  if (info_file.is_open())
203  info_file.close();
204  else
206  << "Unable to open file " << info_file_name.str() << "!" << std::endl;
207  }
208  }
209  ++counter;
210  free(ostream); // free allocated memory
211  }
212  }
213 
215 } // namespace cvn
double E(const int i=0) const
Definition: MCParticle.h:233
void reconfigure(const fhicl::ParameterSet &pset)
int PdgCode() const
Definition: MCParticle.h:212
EventNumber_t event() const
Definition: DataViewImpl.cc:85
double EndZ() const
Definition: MCParticle.h:228
Handle< PROD > getHandle(SelectorBase const &) const
Definition: DataViewImpl.h:382
std::string string
Definition: nybbler.cc:12
void analyze(const art::Event &evt) override
static constexpr double fs
Definition: Units.h:100
def graph(desc, maker=maker)
Definition: apa.py:294
EDAnalyzer(fhicl::ParameterSet const &pset)
Definition: EDAnalyzer.h:25
std::string Process() const
Definition: MCParticle.h:215
Utility class for truth labels.
Particle class.
double EndY() const
Definition: MCParticle.h:227
bool exists(std::string path)
#define DEFINE_ART_MODULE(klass)
Definition: ModuleMacros.h:67
unsigned short GetProtoDUNEBeamInteractionType(const simb::MCParticle &particle) const
T get(std::string const &key) const
Definition: ParameterSet.h:271
unsigned long ulong
Definition: qglobal.h:352
GCNZlibMakerProtoDUNE(fhicl::ParameterSet const &pset)
const sim::ParticleList & ParticleList() const
cet::coded_exception< errors::ErrorCodes, ExceptionDetail::translate > Exception
Definition: Exception.h:66
GCNGraph for GCN.
TCEvent evt
Definition: DataStructs.cxx:7
void fill_ptr_vector(std::vector< Ptr< T >> &ptrs, H const &h)
Definition: Ptr.h:297
double EndX() const
Definition: MCParticle.h:226
QTextStream & endl(QTextStream &s)