FillLSTME_module.cc
Go to the documentation of this file.
10 
11 #include "Geometry/Geometry.h"
12 #include "NovaDAQConventions/DAQConventions.h"
13 #include "RecoBase/Cluster.h"
14 #include "RecoBase/Vertex.h"
15 #include "SummaryData/SpillData.h"
16 #include "Utilities/AssociationUtil.h"
18 
19 
20 #include <memory>
21 
28 
29 namespace LSTME
30 {
31 
32 class FillLSTME : public art::EDProducer
33 {
34 
35 public:
36  explicit FillLSTME(const fhicl::ParameterSet &pset);
37 
38  FillLSTME(FillLSTME const &) = delete;
39  FillLSTME(FillLSTME &&) = delete;
40  FillLSTME & operator = (FillLSTME const &) = delete;
41  FillLSTME & operator = (FillLSTME &&) = delete;
42 
43  void produce(art::Event &evt) override;
44 
45  void reconfigure(const fhicl::ParameterSet &pset);
46 
47  void fillEnergies(
48  art::Event &evt,
49  std::unique_ptr<std::vector<LSTMEnergy>> &sliceEnergies,
50  std::unique_ptr<art::Assns<LSTMEnergy, rb::Cluster>> &sliceEnergyAssoc
51  );
52 
53 private:
55 
58 
63 
64  bool getIsRHC(const art::Event &evt);
65  Model& selectModel(const art::Event &evt);
66 
68  Model &model,
69  const art::Event &evt,
70  const art::Handle<std::vector<rb::Cluster>> &slices_h,
71  size_t sliceIdx,
72  const std::vector<art::Ptr<rb::Vertex>> &elastics,
73  const std::vector<bool> &sliceIsMC
74  );
75 
76 };
77 
78 
80  : modelFDFHC(util::EnvExpansion(pset.get<std::string>("modelFDFHC"))),
81  modelFDRHC(util::EnvExpansion(pset.get<std::string>("modelFDRHC"))),
82  modelNDFHC(util::EnvExpansion(pset.get<std::string>("modelNDFHC"))),
83  modelNDRHC(util::EnvExpansion(pset.get<std::string>("modelNDRHC")))
84 {
85  reconfigure(pset);
86 
87  produces<std::vector<LSTMEnergy>>();
88  produces<art::Assns<LSTMEnergy, rb::Cluster>>();
89 }
90 
92 {
93  beamLabel = pset.get<std::string>("beamLabel");
94  generatorLabel = pset.get<std::string>("generatorLabel");
95 
96  config.parse(pset);
97 }
98 
100 {
101  /* TODO: ask if this should be cached per file/run/job/etc ? */
102 
103  /*
104  * NOTE: Won't work for CRY files.
105  * Apparently need to change config file for them and set
106  * generatorLabel = beamLabel
107  */
109  if (! evt.isRealData()) {
110  evt.getByLabel(generatorLabel, spillPot);
111  }
112  else {
113  evt.getByLabel(beamLabel, spillPot);
114  }
115 
116  if (spillPot.failedToGet())
117  {
118  mf::LogError("FillLSTME") <<
119  "Spill Data not found, aborting without horn current information";
120  abort();
121  }
122 
123  return spillPot->isRHC;
124 }
125 
127 {
129  novadaq::cnv::DetId detId = geom->DetId();
130 
131  bool isRHC = getIsRHC(evt);
132 
133  if (detId == novadaq::cnv::kFARDET)
134  {
135  if (isRHC) {
136  return modelFDRHC;
137  }
138  else {
139  return modelFDFHC;
140  }
141  }
142 
143  if (detId == novadaq::cnv::kNEARDET)
144  {
145  if (isRHC) {
146  return modelNDRHC;
147  }
148  else {
149  return modelNDFHC;
150  }
151  }
152 
153  mf::LogError("FillLSTME") << "Unsupported detector type, aborting";
154  abort();
155 }
156 
158  Model &model,
159  const art::Event &evt,
160  const art::Handle<std::vector<rb::Cluster>> &slices_h,
161  size_t sliceIdx,
162  const std::vector<art::Ptr<rb::Vertex>> &elastics,
163  const std::vector<bool> &sliceIsMC
164 )
165 {
166  VarDictBuilder builder(
167  evt, slices_h, sliceIdx, elastics, config
168  );
169  VarDict varDict = builder.build();
170 
171  /* TODO: this must be removed once geometry is fixed */
173  TrackLengthCorrection trkLenCorr(
174  sliceIsMC[sliceIdx], (geom->DetId() == novadaq::cnv::kFARDET)
175  );
176  trkLenCorr.shift(varDict);
177 
178  return model.predict(varDict);
179 }
180 
182  art::Event &evt,
183  std::unique_ptr<std::vector<LSTMEnergy>> &sliceEnergies,
184  std::unique_ptr<art::Assns<LSTMEnergy, rb::Cluster>> &sliceEnergyAssoc
185 )
186 {
188  if (! evt.getByLabel(config.sliceLabel, slices_h)) {
189  mf::LogError("FillLSTME") << "Failed to get slice info" << std::endl;
190  return;
191  }
192 
193  art::FindManyP<rb::Vertex> elasticVxt_fmp(
194  slices_h, evt, config.elasticArmsLabel
195  );
196  if (! elasticVxt_fmp.isValid()) {
197  mf::LogError("FillLSTME") << "No elastics found" << std::endl;
198  return;
199  }
200 
201  Model &model = selectModel(evt);
202  std::vector<bool> sliceIsMC = getSliceIsMC(*slices_h);
203 
204  for (size_t sliceIdx = 0; sliceIdx < slices_h->size(); ++sliceIdx)
205  {
206  const art::Ptr<rb::Cluster> slicePtr(slices_h, sliceIdx);
207 
208  if (slicePtr->IsNoise()) {
209  continue;
210  }
211 
212  auto elastics = elasticVxt_fmp.at(sliceIdx);
213 
215  model, evt, slices_h, sliceIdx, elastics, sliceIsMC
216  );
217 
218  sliceEnergies->emplace_back(std::move(energy));
219 
221  *this, evt,
222  *(sliceEnergies.get()), slicePtr, *(sliceEnergyAssoc.get())
223  );
224  }
225 }
226 
228 {
229  std::unique_ptr<std::vector<LSTMEnergy>> sliceEnergies(
230  new std::vector<LSTMEnergy>
231  );
232  std::unique_ptr<art::Assns<LSTMEnergy, rb::Cluster>> sliceEnergyAssoc(
234  );
235 
236  fillEnergies(evt, sliceEnergies, sliceEnergyAssoc);
237 
238  evt.put(std::move(sliceEnergies));
239  evt.put(std::move(sliceEnergyAssoc));
240 }
241 
243 
244 }
bool isRHC
is the beam in antineutrino mode, aka RHC
Definition: SpillData.h:28
static bool CreateAssn(art::EDProducer const &prod, art::Event &evt, std::vector< T > &a, art::Ptr< U > b, art::Assns< T, U > &assn, size_t indx=UINT_MAX, std::string const &instance=std::string())
Create a 1 to 1 association between a new product and one already in the event.
Filter events based on their run/event numbers.
FillLSTME(const fhicl::ParameterSet &pset)
void fillEnergies(art::Event &evt, std::unique_ptr< std::vector< LSTMEnergy >> &sliceEnergies, std::unique_ptr< art::Assns< LSTMEnergy, rb::Cluster >> &sliceEnergyAssoc)
void reconfigure(const fhicl::ParameterSet &pset)
void produce(art::Event &evt) override
std::string EnvExpansion(const std::string &inString)
Function to expand environment variables.
Definition: EnvExpand.cxx:8
Definition: config.py:1
MaybeLogger_< ELseverityLevel::ELsev_error, false > LogError
bool isRealData() const
Definition: Event.h:83
DEFINE_ART_MODULE(TestTMapFile)
LSTMEnergy predictEnergy(Model &model, const art::Event &evt, const art::Handle< std::vector< rb::Cluster >> &slices_h, size_t sliceIdx, const std::vector< art::Ptr< rb::Vertex >> &elastics, const std::vector< bool > &sliceIsMC)
LSTMEnergy predict(const VarDict &varDict)
Definition: Model.cxx:27
VarDict build() const
ProductID put(std::unique_ptr< PROD > &&product)
Definition: Event.h:102
Far Detector at Ash River, MN.
T get(std::string const &key) const
Definition: ParameterSet.h:231
novadaq::cnv::DetId DetId() const
Prefer ds::DetectorService::DetId() instead.
Definition: GeometryBase.h:243
double energy
Definition: plottest35.C:25
Near Detector in the NuMI cavern.
void shift(VarDict &varDict) const
Model & selectModel(const art::Event &evt)
Vertex location in position and time.
std::string generatorLabel
::xsd::cxx::tree::string< char, simple_type > string
Definition: Database.h:154
FillLSTME & operator=(FillLSTME const &)=delete
bool getByLabel(std::string const &label, std::string const &productInstanceName, Handle< PROD > &result) const
Definition: DataViewImpl.h:344
void geom(int which=0)
Definition: geom.C:163
Definition: VarDict.h:7
std::vector< bool > getSliceIsMC(const std::vector< rb::Cluster > &slices)
Definition: utils.cxx:6
std::string beamLabel
bool IsNoise() const
Is the noise flag set?
Definition: Cluster.h:163
const XML_Char XML_Content * model
Definition: expat.h:151
Definition: fwd.h:28
Encapsulate the geometry of one entire detector (near, far, ndos)
bool getIsRHC(const art::Event &evt)
bool failedToGet() const
Definition: Handle.h:196