Public Types | Public Member Functions | Static Public Member Functions | Protected Member Functions | Private Attributes | List of all members
numusand::FillSandbox Class Reference
Inheritance diagram for numusand::FillSandbox:
art::EDProducer art::ProducerBase art::Consumer art::EngineCreator art::ProductRegistryHelper

Public Types

using ModuleType = EDProducer
 
using WorkerType = WorkerT< EDProducer >
 
template<typename UserConfig , typename KeysToIgnore = void>
using Table = ProducerBase::Table< UserConfig, KeysToIgnore >
 

Public Member Functions

 FillSandbox (fhicl::ParameterSet const &pset)
 
virtual ~FillSandbox ()
 
void produce (art::Event &evt)
 
void beginRun (art::Run &run)
 
template<typename PROD , BranchType B = InEvent>
ProductID getProductID (std::string const &instanceName={}) const
 
template<typename PROD , BranchType B>
ProductID getProductID (ModuleDescription const &moduleDescription, std::string const &instanceName) const
 
bool modifiesEvent () const
 
template<typename T , BranchType = InEvent>
ProductToken< Tconsumes (InputTag const &)
 
template<typename T , art::BranchType BT>
art::ProductToken< Tconsumes (InputTag const &it)
 
template<typename T , BranchType = InEvent>
void consumesMany ()
 
template<typename Element , BranchType = InEvent>
ViewToken< Element > consumesView (InputTag const &)
 
template<typename T , art::BranchType BT>
art::ViewToken< TconsumesView (InputTag const &it)
 
template<typename T , BranchType = InEvent>
ProductToken< TmayConsume (InputTag const &)
 
template<typename T , art::BranchType BT>
art::ProductToken< TmayConsume (InputTag const &it)
 
template<typename T , BranchType = InEvent>
void mayConsumeMany ()
 
template<typename Element , BranchType = InEvent>
ViewToken< Element > mayConsumeView (InputTag const &)
 
template<typename T , art::BranchType BT>
art::ViewToken< TmayConsumeView (InputTag const &it)
 
base_engine_tcreateEngine (seed_t seed)
 
base_engine_tcreateEngine (seed_t seed, std::string const &kind_of_engine_to_make)
 
base_engine_tcreateEngine (seed_t seed, std::string const &kind_of_engine_to_make, label_t const &engine_label)
 
seed_t get_seed_value (fhicl::ParameterSet const &pset, char const key[]="seed", seed_t const implicit_seed=-1)
 

Static Public Member Functions

static cet::exempt_ptr< Consumernon_module_context ()
 

Protected Member Functions

CurrentProcessingContext const * currentContext () const
 
void validateConsumedProduct (BranchType const bt, ProductInfo const &pi)
 
void prepareForJob (fhicl::ParameterSet const &pset)
 
void showMissingConsumes () const
 

Private Attributes

numusand::NumuSandFxsfNumuSandFxs
 
std::string fSlicerLabel
 
std::string fTrackLabel
 
std::string fEnergyLabel
 
std::string fPIDLabel
 
std::string fVtxLabel
 

Detailed Description

Definition at line 37 of file FillSandbox_module.cc.

Member Typedef Documentation

using art::EDProducer::ModuleType = EDProducer
inherited

Definition at line 34 of file EDProducer.h.

template<typename UserConfig , typename KeysToIgnore = void>
using art::EDProducer::Table = ProducerBase::Table<UserConfig, KeysToIgnore>
inherited

Definition at line 43 of file EDProducer.h.

using art::EDProducer::WorkerType = WorkerT<EDProducer>
inherited

Definition at line 35 of file EDProducer.h.

Constructor & Destructor Documentation

numusand::FillSandbox::FillSandbox ( fhicl::ParameterSet const &  pset)
explicit

Definition at line 62 of file FillSandbox_module.cc.

62  :
63  fSlicerLabel (pset.get< std::string >("SlicerLabel")),
64  fTrackLabel (pset.get< std::string >("TrackLabel")),
65  fEnergyLabel (pset.get< std::string >("EnergyLabel")),
66  fPIDLabel (pset.get< std::string >("PIDLabel")),
67  fVtxLabel (pset.get< std::string >("VtxLabel"))
68  {
69  produces< std::vector<numusand::NumuSandObj> >();
70  produces< art::Assns<numusand::NumuSandObj, rb::Cluster> >();
71  }
enum BeamMode string
numusand::FillSandbox::~FillSandbox ( )
virtual

Definition at line 75 of file FillSandbox_module.cc.

76  {
77  }

Member Function Documentation

void numusand::FillSandbox::beginRun ( art::Run run)
virtual

Reimplemented from art::EDProducer.

Definition at line 80 of file FillSandbox_module.cc.

References fNumuSandFxs.

81  {
82  fNumuSandFxs = new NumuSandFxs();
83  } //end beginRun
numusand::NumuSandFxs * fNumuSandFxs
template<typename T , BranchType = InEvent>
ProductToken<T> art::Consumer::consumes ( InputTag const &  )
inherited
template<typename T , art::BranchType BT>
art::ProductToken<T> art::Consumer::consumes ( InputTag const &  it)
inherited

Definition at line 146 of file Consumer.h.

References art::InputTag::instance(), PandAna.reco_validation.prod5_pid_validation::invalid, art::InputTag::label(), art::InputTag::process(), and T.

147 {
148  if (!moduleContext_)
149  return ProductToken<T>::invalid();
150 
151  consumables_[BT].emplace_back(ConsumableType::Product,
152  TypeID{typeid(T)},
153  it.label(),
154  it.instance(),
155  it.process());
156  return ProductToken<T>{it};
157 }
set< int >::iterator it
static ProductToken< T > invalid()
Definition: ProductToken.h:47
ConsumableProducts consumables_
Definition: Consumer.h:137
double T
Definition: Xdiff_gwt.C:5
bool moduleContext_
Definition: Consumer.h:135
template<typename T , art::BranchType BT>
void art::Consumer::consumesMany ( )
inherited

Definition at line 161 of file Consumer.h.

References T.

162 {
163  if (!moduleContext_)
164  return;
165 
166  consumables_[BT].emplace_back(ConsumableType::Many, TypeID{typeid(T)});
167 }
ConsumableProducts consumables_
Definition: Consumer.h:137
double T
Definition: Xdiff_gwt.C:5
bool moduleContext_
Definition: Consumer.h:135
template<typename Element , BranchType = InEvent>
ViewToken<Element> art::Consumer::consumesView ( InputTag const &  )
inherited
template<typename T , art::BranchType BT>
art::ViewToken<T> art::Consumer::consumesView ( InputTag const &  it)
inherited

Definition at line 171 of file Consumer.h.

References art::InputTag::instance(), PandAna.reco_validation.prod5_pid_validation::invalid, art::InputTag::label(), art::InputTag::process(), and T.

172 {
173  if (!moduleContext_)
174  return ViewToken<T>::invalid();
175 
176  consumables_[BT].emplace_back(ConsumableType::ViewElement,
177  TypeID{typeid(T)},
178  it.label(),
179  it.instance(),
180  it.process());
181  return ViewToken<T>{it};
182 }
set< int >::iterator it
static ViewToken< Element > invalid()
Definition: ProductToken.h:75
ConsumableProducts consumables_
Definition: Consumer.h:137
double T
Definition: Xdiff_gwt.C:5
bool moduleContext_
Definition: Consumer.h:135
base_engine_t& art::EngineCreator::createEngine ( seed_t  seed)
inherited
base_engine_t& art::EngineCreator::createEngine ( seed_t  seed,
std::string const &  kind_of_engine_to_make 
)
inherited
base_engine_t& art::EngineCreator::createEngine ( seed_t  seed,
std::string const &  kind_of_engine_to_make,
label_t const &  engine_label 
)
inherited
CurrentProcessingContext const* art::EDProducer::currentContext ( ) const
protectedinherited
seed_t art::EngineCreator::get_seed_value ( fhicl::ParameterSet const &  pset,
char const  key[] = "seed",
seed_t const  implicit_seed = -1 
)
inherited
template<typename PROD , BranchType B>
ProductID art::EDProducer::getProductID ( std::string const &  instanceName = {}) const
inlineinherited

Definition at line 123 of file EDProducer.h.

References art::EDProducer::moduleDescription_.

Referenced by skim::NueSkimmer::CopyMichelSlice(), and skim::NueSkimmer::CopyMichelTrack().

124  {
125  return ProducerBase::getProductID<PROD, B>(moduleDescription_,
126  instanceName);
127  }
ModuleDescription moduleDescription_
Definition: EDProducer.h:115
template<typename PROD , BranchType B>
ProductID art::ProducerBase::getProductID ( ModuleDescription const &  moduleDescription,
std::string const &  instanceName 
) const
inherited

Definition at line 56 of file ProducerBase.h.

References art::ModuleDescription::moduleLabel().

Referenced by art::ProducerBase::modifiesEvent().

58  {
59  auto const& pd =
60  get_ProductDescription<PROD>(B, md.moduleLabel(), instanceName);
61  return pd.productID();
62  }
template<typename T , BranchType = InEvent>
ProductToken<T> art::Consumer::mayConsume ( InputTag const &  )
inherited
template<typename T , art::BranchType BT>
art::ProductToken<T> art::Consumer::mayConsume ( InputTag const &  it)
inherited

Definition at line 189 of file Consumer.h.

References art::InputTag::instance(), PandAna.reco_validation.prod5_pid_validation::invalid, art::InputTag::label(), art::InputTag::process(), and T.

190 {
191  if (!moduleContext_)
192  return ProductToken<T>::invalid();
193 
194  consumables_[BT].emplace_back(ConsumableType::Product,
195  TypeID{typeid(T)},
196  it.label(),
197  it.instance(),
198  it.process());
199  return ProductToken<T>{it};
200 }
set< int >::iterator it
static ProductToken< T > invalid()
Definition: ProductToken.h:47
ConsumableProducts consumables_
Definition: Consumer.h:137
double T
Definition: Xdiff_gwt.C:5
bool moduleContext_
Definition: Consumer.h:135
template<typename T , art::BranchType BT>
void art::Consumer::mayConsumeMany ( )
inherited

Definition at line 204 of file Consumer.h.

References T.

205 {
206  if (!moduleContext_)
207  return;
208 
209  consumables_[BT].emplace_back(ConsumableType::Many, TypeID{typeid(T)});
210 }
ConsumableProducts consumables_
Definition: Consumer.h:137
double T
Definition: Xdiff_gwt.C:5
bool moduleContext_
Definition: Consumer.h:135
template<typename Element , BranchType = InEvent>
ViewToken<Element> art::Consumer::mayConsumeView ( InputTag const &  )
inherited
template<typename T , art::BranchType BT>
art::ViewToken<T> art::Consumer::mayConsumeView ( InputTag const &  it)
inherited

Definition at line 214 of file Consumer.h.

References art::InputTag::instance(), PandAna.reco_validation.prod5_pid_validation::invalid, art::InputTag::label(), art::InputTag::process(), and T.

215 {
216  if (!moduleContext_)
217  return ViewToken<T>::invalid();
218 
219  consumables_[BT].emplace_back(ConsumableType::ViewElement,
220  TypeID{typeid(T)},
221  it.label(),
222  it.instance(),
223  it.process());
224  return ViewToken<T>{it};
225 }
set< int >::iterator it
static ViewToken< Element > invalid()
Definition: ProductToken.h:75
ConsumableProducts consumables_
Definition: Consumer.h:137
double T
Definition: Xdiff_gwt.C:5
bool moduleContext_
Definition: Consumer.h:135
bool art::ProducerBase::modifiesEvent ( ) const
inlineinherited

Definition at line 40 of file ProducerBase.h.

References art::ProducerBase::getProductID(), and string.

41  {
42  return true;
43  }
static cet::exempt_ptr<Consumer> art::Consumer::non_module_context ( )
staticinherited
void art::Consumer::prepareForJob ( fhicl::ParameterSet const &  pset)
protectedinherited
void numusand::FillSandbox::produce ( art::Event evt)
virtual

ADD HERE VARS FOR QE PID?

Implements art::EDProducer.

Definition at line 86 of file FillSandbox_module.cc.

References abs(), rb::Cluster::Add(), rb::Cluster::AllCells(), rb::Cluster::CalorimetricEnergy(), rb::CellHit::Cell(), rb::Cluster::Cell(), om::cout, util::CreateAssn(), d, DEFINE_ART_MODULE(), allTimeWatchdog::endl, fEnergyLabel, fNumuSandFxs, fPIDLabel, fSlicerLabel, fTrackLabel, fVtxLabel, geom(), art::Ptr< T >::get(), numusand::NumuSandFxs::getActivityVtx(), numusand::NumuSandFxs::getAveTrackdEdx(), numusand::NumuSandFxs::getAveTrackdEdxLast4Cells(), numusand::NumuSandFxs::getAveTrackdEdxLast6Cells(), numusand::NumuSandFxs::getAveTrackdEdxLast8Cells(), art::DataViewImpl::getByLabel(), numusand::NumuSandFxs::getMissingE(), numusand::NumuSandFxs::getScatt(), rb::RecoHit::GeV(), make_syst_table_plots::h, cheat::BackTracker::HaveTruthInfo(), cheat::BackTracker::HitsToParticle(), cheat::BackTracker::HitsToTrackIDE(), MECModelEnuComparisons::i, rb::Prong::Is3D(), rb::RecoHit::IsCalibrated(), rb::Cluster::IsNoise(), calib::j, geo::kX, geo::kY, cafExposure::match, std::min(), geo::LiveGeometry::NBadDCMBC(), rb::Cluster::NCell(), geo::PlaneGeo::Ncells(), geo::LiveGeometry::NDropouts(), submit_concat_project::parts, genie::utils::res::PdgCode(), rb::CellHit::PE(), rb::CellHit::Plane(), geo::GeometryBase::Plane(), art::Event::put(), rb::Cluster::RecoHit(), numusand::NumuSandObj::SetActVtx(), numusand::NumuSandObj::SetAvedEdxTrack1(), numusand::NumuSandObj::SetAvedEdxTrack1Last4Cells(), numusand::NumuSandObj::SetAvedEdxTrack1Last6Cells(), numusand::NumuSandObj::SetAvedEdxTrack1Last8Cells(), numusand::NumuSandObj::SetAvedEdxTrack2(), numusand::NumuSandObj::SetAvedEdxTrack2Last4Cells(), numusand::NumuSandObj::SetAvedEdxTrack2Last6Cells(), numusand::NumuSandObj::SetAvedEdxTrack2Last8Cells(), numusand::NumuSandObj::SetHadCellsFromEdge(), numusand::NumuSandObj::SetMissingE(), numusand::NumuSandObj::SetNBad(), numusand::NumuSandObj::SetNDropped(), numusand::NumuSandObj::SetNHadHits(), numusand::NumuSandObj::SetNMuTrks(), numusand::NumuSandObj::SetNProtons(), numusand::NumuSandObj::SetOffTrkFra(), numusand::NumuSandObj::SetPDGBest(), numusand::NumuSandObj::SetPDGSecondBest(), numusand::NumuSandObj::SetPiMuDecay(), numusand::NumuSandObj::SetScattAngleTrack1(), numusand::NumuSandObj::SetScattAngleTrack2(), numusand::NumuSandObj::SetVtxE20(), numusand::NumuSandObj::SetVtxE40(), numusand::NumuSandObj::SetVtxE60(), numusand::NumuSandObj::SetXMax2(), numusand::NumuSandObj::SetXMaxt(), numusand::NumuSandObj::SetXMin2(), numusand::NumuSandObj::SetXMint(), numusand::NumuSandObj::SetYMax2(), numusand::NumuSandObj::SetYMaxt(), numusand::NumuSandObj::SetYMin2(), numusand::NumuSandObj::SetYMint(), numusand::NumuSandObj::SetZMax2(), numusand::NumuSandObj::SetZMaxt(), numusand::NumuSandObj::SetZMin2(), numusand::NumuSandObj::SetZMint(), art::PtrVector< T >::size(), util::sqr(), std::sqrt(), rb::Cluster::TotalGeV(), rb::Track::TotalLength(), rb::PID::Value(), rb::CellHit::View(), rb::RecoHit::X(), submit_syst::x, xmax, make_mec_shifts_plots::xmin, rb::RecoHit::Y(), submit_syst::y, ymax, ymin, test::z, rb::RecoHit::Z(), and make_true_q0q3_plots::zmax.

87  {
88 
90  evt.getByLabel(fSlicerLabel,slicevec);
91 
92  if(slicevec->empty()) {
93  mf::LogWarning ("No Slices")<<"No Slices in the input file";
94  return;
95  }
96 
100 
101  art::FindManyP<rb::Track> trackAssnList(slicevec, evt, fTrackLabel);
102 
103  // set up associations between slices and energy
104  art::FindManyP<numue::NumuE> eAss(slicevec,evt,fEnergyLabel);
105 
106  std::unique_ptr< std::vector<numusand::NumuSandObj> > outputObjects(new std::vector<numusand::NumuSandObj>);
107  std::unique_ptr< art::Assns<numusand::NumuSandObj, rb::Cluster> > assoc(new art::Assns<numusand::NumuSandObj, rb::Cluster>);
108 
109  for(size_t i = 0; i < slicevec->size(); ++i){
110 
111  art::Ptr<rb::Cluster> slice(slicevec,i);
112 
113  //get the total no. of cell hits
114  art::PtrVector<rb::CellHit> allslchits = slice->AllCells();
115 
116  if(slice->IsNoise()) continue;
117 
118  numusand::NumuSandObj sandbox; // initialization in constructor
119 
120  float missE;
121 
122  fNumuSandFxs->getMissingE(slice, missE);
123 
124  sandbox.SetMissingE(missE);
125 
126  sandbox.SetNDropped(livegeom->NDropouts());
127  sandbox.SetNBad(livegeom->NBadDCMBC());
128 
129  float xmin = 9999, ymin = 9999, zmin = 9999, xmax = -9999, ymax = -9999, zmax = -9999;
130  float xmin2 = 9999, ymin2 = 9999, zmin2 = 9999, xmint = 9999, ymint = 9999, zmint = 9999;
131  float xmax2 = -9999, ymax2 = -9999, zmax2 = -9999, xmaxt = -9999, ymaxt = -9999, zmaxt = -9999;
132  int ntx = 0, nty = 0;
133 
134 
135  // Temporary code. Not a great place to put these and not much time to
136  // re-order anything. Calculate vtx E within 20/40/60 cm to origin.
137  art::FindManyP<rb::Vertex> fmElastic(slicevec, evt,fVtxLabel);
138  if(fmElastic.isValid()){
139  std::vector< art::Ptr<rb::Vertex> > elastics = fmElastic.at(i);
140  if (elastics.size() == 1){
141  // Now, we can get the vtx, and loop over cells to get vars
142  TVector3 vtx = elastics[0]->GetXYZ();
143  float vtxE20 = 0;
144  float vtxE40 = 0;
145  float vtxE60 = 0;
146 
147  for (unsigned int hIdx = 0; hIdx < slice->NCell(); hIdx++){
148  rb::RecoHit rhit = slice->RecoHit(hIdx);
149  if (!rhit.IsCalibrated()) continue;
150  bool isX = slice->Cell(hIdx)->View()==geo::kX;
151  // Calcualte dist from hit to vtx.
152  // Only calculate using coordinate cell measures! (DocDB 14330)
153  double d =
154  isX ? util::sqr(rhit.X()-vtx[0]) :util::sqr(rhit.Y()-vtx[1]);
155  d += util::sqr(rhit.Z()-vtx[2]);
156  d = sqrt(d);
157  if (d < 20){
158  vtxE20 += rhit.GeV();
159  vtxE40 += rhit.GeV();
160  vtxE60 += rhit.GeV();
161  }
162  else if (d < 40){
163  vtxE40 += rhit.GeV();
164  vtxE60 += rhit.GeV();
165  }
166  else if (d < 60)
167  vtxE60 += rhit.GeV();
168  } // End loop over NCells
169  sandbox.SetVtxE20(vtxE20);
170  sandbox.SetVtxE40(vtxE40);
171  sandbox.SetVtxE60(vtxE60);
172  } // End NVtx == 1
173  } // End isValid
174 
175 
176 
177  for(unsigned int j = 0; j < slice->NCell(); ++j) {
178  const rb::CellHit* h = slice->Cell(j).get();
179  double xyz[3];
180  geom->Plane(h->Plane())->Cell(h->Cell())->GetCenter(xyz);
181  const double x = xyz[0];
182  const double y = xyz[1];
183  const double z = xyz[2];
184 
185  if(h->View()==geo::kX) {
186  if(x < xmin) { xmin2 = xmin; xmin = x;}
187  else if(x < xmin2) xmin2 = x;
188  if(z < zmin) { zmin2 = zmin; zmin = z;}
189  else if(z < zmin2) zmin2 = z;
190  if(x > xmax) { xmax2 = xmax; xmax = x;}
191  else if(x > xmax2) xmax2 = x;
192  if(z > zmax) { zmax2 = zmax; zmax = z;}
193  else if (z > zmax2) zmax2 = zmax;
194 
195  if(h->PE() > 100) { // above threshold only TODO: add tracks and recohits, use pecorr
196  ntx++;
197  if(x < xmint) xmint = x;
198  if(z < zmint) zmint = z;
199 
200  if(x > xmaxt) xmaxt = x;
201  if(z > zmaxt) zmaxt = z;
202  }
203  }
204 
205  else if(h->View()==geo::kY) {
206  if(y < ymin) { ymin2 = ymin; ymin = y;}
207  else if (y < ymin2) ymin2 = y;
208  if(z < zmin) { zmin2 = zmin; zmin = z;}
209  else if (z < zmin2) zmin2 = z;
210 
211  if(y > ymax) { ymax2 = ymax; ymax = y;}
212  else if(y > ymax2) ymax2 = y;
213  if(z > zmax) { zmax2 = zmax; zmax = z;}
214  else if(z > zmax2) zmax2 = z;
215 
216  if(h->PE() > 100) { // above threshold only TODO: add tracks and recohits, use pecorr
217  nty++;
218  if(y < ymint) ymint = y;
219  if(z < zmint) zmint = z;
220 
221  if(y > ymaxt) ymaxt = y;
222  if(z > zmaxt) zmaxt = z;
223  }
224  }
225 
226  }
227  sandbox.SetXMin2(xmin2);
228  sandbox.SetYMin2(ymin2);
229  sandbox.SetZMin2(zmin2);
230  sandbox.SetXMax2(xmax2);
231  sandbox.SetYMax2(ymax2);
232  sandbox.SetZMax2(zmax2);
233 
234  if(ntx > 0 && nty > 0) {
235  sandbox.SetXMint(xmint);
236  sandbox.SetYMint(ymint);
237  sandbox.SetZMint(zmint);
238  sandbox.SetXMaxt(xmaxt);
239  sandbox.SetYMaxt(ymaxt);
240  sandbox.SetZMaxt(zmaxt);
241  }
242 
243  else { // -5 default can be confused as a real value; +=9999 too far from real values
244  sandbox.SetXMint(-999);
245  sandbox.SetYMint(-999);
246  sandbox.SetZMint(-999);
247  sandbox.SetXMaxt(-999);
248  sandbox.SetYMaxt(-999);
249  sandbox.SetZMaxt(-999);
250  }
251 
252 
253  if(!trackAssnList.isValid()) {
254  std::cout << "NumuSandbox: No Kalman Tracks!" << std::endl;
255  outputObjects->push_back(sandbox);
256  util::CreateAssn(*this,evt,*(outputObjects.get()),slice,*(assoc.get()));
257  continue;
258  }
259 
260  const std::vector< art::Ptr<rb::Track> > sliceTracks = trackAssnList.at(i);
261 
262  art::FindOneP<remid::ReMId> remidVec(sliceTracks, evt, fPIDLabel);
263  if(!remidVec.isValid()) {
264  // std::cout << "NumuSandbox: No ReMId information!" << std::endl;
265  outputObjects->push_back(sandbox);
266  util::CreateAssn(*this,evt,*(outputObjects.get()),slice,*(assoc.get()));
267  continue;
268  }
269 
270 
271  int bestIdx = -1; // for index in list of the highest PID track
272  int pimuIdx = -5; // for index in list of a pimu decay event
273  int Idx2 = -1;
274  int idx3d = -1;
275  double bestpimuE = -999.; // keep track of how much energy due to a non-primary muon was in a track
276  // this is for if a pi->mu muon contributes light to more than once track
277  // we only want to associate with the track that saw 'most' of it.
278  // note we assume all non-primary muons are from pi->mu decay.
279  double bestval = -999., val2 = -999;
280  int nmutrks = 0;
281 
282  // ----------------------------------
283  // WORKING AREA
284 
285  float avededxtrk1 = 0.0;
286  float avededxtrk2 = 0.0;
287  float avededxtrk1Last4Cells = 0.0;
288  float avededxtrk2Last4Cells = 0.0;
289  float avededxtrk1Last6Cells = 0.0;
290  float avededxtrk2Last6Cells = 0.0;
291  float avededxtrk1Last8Cells = 0.0;
292  float avededxtrk2Last8Cells = 0.0;
293  float scattangletrk1 = 0.0;
294  float scattangletrk2 = 0.0;
295  float actvtxx = 0.0;
296  float slccalE = 0.0;
297 
298  if((sliceTracks.size() == 2)){
299  // Loop over tracks in slice
300  for(size_t j = 0; j < sliceTracks.size(); ++j){
301  const art::Ptr<rb::Track> track = sliceTracks[j];
302  art::Ptr<remid::ReMId> trackRemid = remidVec.at(j);
303 
304  if(!track->Is3D()) continue;
305 
306  if(trackRemid->Value() > bestval) { // best track so far
307  val2 = bestval;
308  Idx2 = bestIdx;
309  bestval = trackRemid->Value();
310  bestIdx = j;
311  }
312  else if(trackRemid->Value() > val2) { // not best, but better second best track
313  val2 = trackRemid->Value();
314  Idx2 = j;
315  }
316 
317  } // end loop over tracks
318 
319  // Making sure is OK
320  if((!remidVec.isValid()) ||
321  (bestval < 0.0) || (val2 < 0.0)){
322  continue;
323  }
324 
325  // Get info on muon like track
326  //if ((sliceTracks[bestval]->NCell()>0) && (bestval>0.75)){
327  // get info about this track
328  const art::Ptr<rb::Track> track1 = sliceTracks[bestIdx];
329  //
330  // Fill here with dEdx info and other variables
331  //
332  fNumuSandFxs->getAveTrackdEdx(track1,avededxtrk1);
333  sandbox.SetAvedEdxTrack1(avededxtrk1);
334 
335  fNumuSandFxs->getAveTrackdEdxLast4Cells(track1,avededxtrk1Last4Cells);
336  sandbox.SetAvedEdxTrack1Last4Cells(avededxtrk1Last4Cells);
337 
338  fNumuSandFxs->getAveTrackdEdxLast6Cells(track1,avededxtrk1Last6Cells);
339  sandbox.SetAvedEdxTrack1Last6Cells(avededxtrk1Last6Cells);
340 
341  fNumuSandFxs->getAveTrackdEdxLast8Cells(track1,avededxtrk1Last8Cells);
342  sandbox.SetAvedEdxTrack1Last8Cells(avededxtrk1Last8Cells);
343 
344  scattangletrk1 = fNumuSandFxs->getScatt(track1);
345  if(track1->TotalLength() > 0) sandbox.SetScattAngleTrack1(scattangletrk1/track1->TotalLength());
346 
347  //}
348  //else if((sliceTracks[val2]->NCell()>0)){
349  // get info about this track
350  const art::Ptr<rb::Track> track2 = sliceTracks[Idx2];
351  //
352  // Fill here with dEdx info and other variables
353  //
354  fNumuSandFxs->getAveTrackdEdx(track2,avededxtrk2);
355  sandbox.SetAvedEdxTrack2(avededxtrk2);
356 
357  fNumuSandFxs->getAveTrackdEdxLast4Cells(track2,avededxtrk2Last4Cells);
358  sandbox.SetAvedEdxTrack2Last4Cells(avededxtrk2Last4Cells);
359 
360  fNumuSandFxs->getAveTrackdEdxLast6Cells(track2,avededxtrk2Last6Cells);
361  sandbox.SetAvedEdxTrack2Last6Cells(avededxtrk2Last6Cells);
362 
363  fNumuSandFxs->getAveTrackdEdxLast8Cells(track2,avededxtrk2Last8Cells);
364  sandbox.SetAvedEdxTrack2Last8Cells(avededxtrk2Last8Cells);
365 
366  scattangletrk2 = fNumuSandFxs->getScatt(track2);
367  if(track2->TotalLength() > 0) sandbox.SetScattAngleTrack2(scattangletrk2/track2->TotalLength());
368  //}
369 
370  ///////////////////////////////////
371  ///--------------------------------
372  /// ADD HERE VARS FOR QE PID?
373  ///-------------------------------
374  ///////////////////////////////////
375 
376  // WORKING AREA
377  //art::PtrVector<rb::CellHit> allslchits = slice->AllCells();
378  art::PtrVector<rb::CellHit> alltrk1hits = track1->AllCells();
379  art::PtrVector<rb::CellHit> alltrk2hits = track2->AllCells();
380 
381  fNumuSandFxs->getActivityVtx(alltrk1hits,alltrk2hits,allslchits,actvtxx);
382  slccalE = slice->CalorimetricEnergy();
383 
384  if(slice->CalorimetricEnergy() <= 0){
385  //continue;
386  std::cout<<"Calorimetric energy is zero!!"<<std::endl;
387  }
388  else{
389  sandbox.SetActVtx(actvtxx/slccalE);
390  //std::cout<<"TMVAvars[0]: "<< actvtxx/slccalE <<std::endl;
391  //TMVAvars[0] = actvtx/slccalE;
392  }
393 
394  //----
395  // Another working area
396  // from QEEventFinder/QeFinder_module.cc
397  // get numu energy estimates
398  float trkCCE = -5.0;
399  float RecoMuonE = -5.0;
400  float Reco2ndTrkE = -5.0;
401 
402  if(eAss.isValid()){
403  std::vector<art::Ptr<numue::NumuE> > numuEnergy = eAss.at(i);
404  if(numuEnergy.size() > 0){
405  // get the energy variables
406  trkCCE = numuEnergy[0]->TrkCCE();
407  }//end if valid size
408 
409  art::FindOneP<rb::Energy> trackEnergy(sliceTracks, evt, fEnergyLabel);
410  if(trackEnergy.isValid()){
411  RecoMuonE = trackEnergy.at(bestIdx)->E();
412  }
413  //Reco2ndTrkE = track2->TotalGeV(); // sloppy approach
414  Reco2ndTrkE = track2->CalorimetricEnergy(); // sloppy approach
415 
416  float tmpHadE = trkCCE - RecoMuonE;
417 
418  float offTrkEnFrac = (tmpHadE - Reco2ndTrkE)/trkCCE;
419 
420  sandbox.SetOffTrkFra(offTrkEnFrac);
421 
422  //std::cout<<"TMVAvars[0]: "<< offTrkEnFrac <<std::endl;
423 
424  }// eAss is valid
425  else{
426  std::cout<<"Event Association is not Valid!"<<std::endl;
427  }
428 
429  }//end slice tracks = 2
430 
431  // ----------------------------------
432 
433  // Loop over tracks in slice
434  for(size_t j = 0; j < sliceTracks.size(); ++j){
435  const art::Ptr<rb::Track> track = sliceTracks[j];
436  art::Ptr<remid::ReMId> trackRemid = remidVec.at(j);
437 
438  if(!track->Is3D()) {
439  continue;
440  }
441 
442  idx3d++;
443 
444  if(trackRemid->Value() > bestval) { // best track so far
445  val2 = bestval;
446  Idx2 = bestIdx;
447  bestval = trackRemid->Value();
448  bestIdx = j;
449  }
450  else if(trackRemid->Value() > val2) { // not best, but better second best track
451  val2 = trackRemid->Value();
452  Idx2 = j;
453  }
454 
455  if(trackRemid->Value() > 0.6) nmutrks++;
456 
457  if(bt->HaveTruthInfo()) {
458  art::PtrVector<rb::CellHit> trkhits = track->AllCells();
459  const std::vector<const sim::Particle*> parts = bt->HitsToParticle(trkhits);
460  const std::vector< cheat::TrackIDE > trids = bt->HitsToTrackIDE(trkhits);
461 
462  // check particle that contributes most to track for pi->mu decay
463  if( !parts.empty() && !trids.empty() ) {
464  if(abs(parts[0]->PdgCode())==13 && parts[0]->Process()=="Decay") {
465  double tempE = trids[0].energyFrac * track->TotalGeV();
466  if(tempE > bestpimuE) {
467  bestpimuE = tempE;
468  pimuIdx = idx3d+1;
469  }
470  }
471  }
472 
473  // check particle that contributes second most; track can contain both pi and mu
474  if(parts.size() > 1 && trids.size() > 1) {
475  if(abs(parts[1]->PdgCode())==13 && parts[1]->Process()=="Decay") {
476  double tempE = trids[1].energyFrac * track->TotalGeV();
477  if(tempE > bestpimuE) {
478  bestpimuE = tempE;
479  pimuIdx = idx3d+1;
480  }
481  }
482  }
483  }
484  } //end loop over tracks
485 
486  // ===============================================================
487  //New: Get info on muon like track
488  //if( (bestval > 0.0) && (val2 > 0.0) ){
489  if((bestval > 0.0)){
490 
491  const art::Ptr<rb::Track> mutrack = sliceTracks[bestIdx];
492  art::PtrVector<rb::CellHit> mutrkhits = mutrack->AllCells();
493 
494  //Getting all the hits in the slice that are not on our best track
495  ////rb::Cluster vertexHits = MakeVertexCluster(allslchits,mutrkhits);
496  rb::Cluster vertexCluster;
497  int nhadhits = -5;
498  for (unsigned int sliceHit = 0; sliceHit < allslchits.size(); ++sliceHit){
499 
500  bool match = 0;
501 
502  for (unsigned int trackHit = 0; trackHit < mutrkhits.size(); ++trackHit){
503  match = (allslchits[sliceHit] == mutrkhits[trackHit]);
504  if (match) break;
505  }//End of loop over hits in track
506 
507  if (!match) vertexCluster.Add(allslchits[sliceHit]);
508 
509  }//End of loop over hits in the slice
510 
511  nhadhits = vertexCluster.NCell();
512  //std::cout<<" ==== No. hits in had system: "<<vertexCluster.NCell()<<std::endl;
513 
514  sandbox.SetNHadHits(nhadhits);
515 
516  // New
517  unsigned int minCellsFromEdge = 99999999;
518  //int cellsFromEdge = -1;
519 
520  // TO DO: case when nhadhits is zero?
521  // Loop over vertex cluster
522  for(unsigned int hitIdx = 0; hitIdx < vertexCluster.NCell(); ++hitIdx){
523  const art::Ptr<rb::CellHit>& chit = vertexCluster.Cell(hitIdx);
524  const int planeNum = chit->Plane();
525  const unsigned int cellsFromEdge = std::min((unsigned int)chit->Cell(), geom->Plane(planeNum)->Ncells() - 1 - chit->Cell());
526 
527  if(cellsFromEdge < minCellsFromEdge){
528  minCellsFromEdge = cellsFromEdge;
529  } //end if
530 
531  } //end for
532  //std::cout<<" ==== Had ncells from edge: "<< minCellsFromEdge <<std::endl;
533 
534  sandbox.SetHadCellsFromEdge(minCellsFromEdge);
535 
536  }//end if best val
537  // ===============================================================
538 
539  sandbox.SetNMuTrks(nmutrks);
540 
541  if(bt->HaveTruthInfo()) { // no truth info case should correspond to unfilled
542 
543  if(bestIdx!=-1) {
544 
545  const art::Ptr<rb::Track> track1 = sliceTracks[bestIdx];
546  const std::vector<const sim::Particle*> parts1 = bt->HitsToParticle(track1->AllCells());
547  const std::vector< cheat::TrackIDE > trids1 = bt->HitsToTrackIDE(track1->AllCells());
548 
549  if(!parts1.empty()){
550  sandbox.SetPDGBest(parts1[0]->PdgCode());
551  }
552  }
553 
554  if(Idx2!=-1) {
555 
556  const art::Ptr<rb::Track> track2 = sliceTracks[Idx2];
557  const std::vector<const sim::Particle*> parts2 = bt->HitsToParticle(track2->AllCells());
558  const std::vector< cheat::TrackIDE > trids2 = bt->HitsToTrackIDE(track2->AllCells());
559 
560  if(!parts2.empty()){
561  sandbox.SetPDGSecondBest(parts2[0]->PdgCode());
562  }
563  }
564 
565  int nprotons = 0;
566  if(pimuIdx == -5) pimuIdx = 0; // differentiate unfilled with not found
567  art::PtrVector<rb::CellHit> slchits = slice->AllCells();
568  const std::vector<const sim::Particle*> parts = bt->HitsToParticle(slchits);
569  const std::vector< cheat::TrackIDE > trids = bt->HitsToTrackIDE(slchits);
570  if(!parts.empty() && !trids.empty()) {
571  for(unsigned int i = 0; i < parts.size(); i++) {
572  if(trids.size() > i) {
573  if(parts[i]->PdgCode()==2212 && trids[i].energyFrac * slice->TotalGeV() > 0.05) nprotons++;
574  }
575  }
576  }
577  sandbox.SetNProtons(nprotons);
578  }
579 
580  sandbox.SetPiMuDecay(pimuIdx);
581 
582  outputObjects->push_back(sandbox);
583  util::CreateAssn(*this,evt,*(outputObjects.get()),slice,*(assoc.get()));
584  } //end loop over slices
585 
586  evt.put(std::move(outputObjects));
587  evt.put(std::move(assoc));
588 
589  } //end produce
void SetNHadHits(int nhadhits)
void SetXMaxt(float dxmaxt)
void SetPDGSecondBest(int pdgsecondbest)
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.
unsigned int NCell(geo::View_t view) const
Number of cells in view view.
Definition: Cluster.cxx:134
std::map< std::string, double > xmax
void SetZMaxt(float dzmaxt)
void SetXMint(float dxmint)
void SetAvedEdxTrack1Last4Cells(float avededxtrk1Last4Cells)
unsigned short Plane() const
Definition: CellHit.h:39
void SetYMax2(float dymax2)
Definition: NumuSandObj.cxx:92
void SetAvedEdxTrack2Last4Cells(float avededxtrk2Last4Cells)
void SetOffTrkFra(float offtrk)
proxy for off track energy
geo::View_t View() const
Definition: CellHit.h:41
void getAveTrackdEdx(const art::Ptr< rb::Track > track, float &avededx)
T sqrt(T number)
Definition: d0nt_math.hpp:156
void getAveTrackdEdxLast4Cells(const art::Ptr< rb::Track > track, float &avededxlast4)
Vertical planes which measure X.
Definition: PlaneGeo.h:28
unsigned int Ncells() const
Number of cells in this plane.
Definition: PlaneGeo.h:43
float Z() const
Definition: RecoHit.h:38
Definition: event.h:19
A collection of associated CellHits.
Definition: Cluster.h:47
float getScatt(const art::Ptr< rb::Track > track)
void abs(TH1 *hist)
void SetAvedEdxTrack2Last6Cells(float avededxtrk2Last6Cells)
void SetScattAngleTrack1(float scattangletrk1)
void SetAvedEdxTrack2Last8Cells(float avededxtrk2Last8Cells)
const PlaneGeo * Plane(unsigned int i) const
void SetPDGBest(int pdgbest)
T sqr(T x)
More efficient square function than pow(x,2)
Definition: MathUtil.h:23
double Value() const
Definition: PID.h:22
Horizontal planes which measure Y.
Definition: PlaneGeo.h:29
void SetXMin2(float dxmin2)
minimum position of hit second farthest from edge in that dimension
Definition: NumuSandObj.cxx:68
art::PtrVector< rb::CellHit > AllCells() const
Get all cells from both views.
Definition: Cluster.cxx:180
void SetNMuTrks(int nmutrks)
Calibrated quantities relying on position in the orthogonal view. To generate a rb::CellHit from a rb...
Definition: RecoHit.h:19
void SetPiMuDecay(int pitomu)
Definition: NumuSandObj.cxx:56
virtual void Add(const art::Ptr< rb::CellHit > &cell, double weight=1)
Definition: Cluster.cxx:84
void SetVtxE20(float vtxE20)
Double_t ymax
Definition: plot.C:25
ProductID put(std::unique_ptr< PROD > &&product)
Definition: Event.h:102
virtual double TotalLength() const
Length (cm) of all the track segments.
Definition: Track.cxx:213
unsigned short Cell() const
Definition: CellHit.h:40
int PdgCode(Resonance_t res, int Q)
(resonance id, charge) -> PDG code
void getActivityVtx(const art::PtrVector< rb::CellHit > track1Hits, const art::PtrVector< rb::CellHit > track2Hits, const art::PtrVector< rb::CellHit > sliceHits, float &actvtx)
void SetYMaxt(float dymaxt)
minimum position of hit closest to edge in that dimension that has more than 100 photoelectons (not p...
double CalorimetricEnergy(EEnergyCalcScheme escheme=kRecomputeEnergy) const
Simple estimate of neutrino energy.
Definition: Cluster.cxx:439
void SetHadCellsFromEdge(int hadcellsfromedge)
void SetYMint(float dymint)
minimum position of hit closest to edge in that dimension that has more than 100 photoelectons (not p...
void SetAvedEdxTrack1Last6Cells(float avededxtrk1Last6Cells)
void SetAvedEdxTrack1(float avededxtrk1)
float PE() const
Definition: CellHit.h:42
Float_t d
Definition: plot.C:236
void getAveTrackdEdxLast8Cells(const art::Ptr< rb::Track > track, float &avededxlast8)
void SetZMint(float dzmint)
bool IsCalibrated() const
You MUST check here before accessing PECorr, MIP or GeV.
Definition: RecoHit.cxx:35
const double j
Definition: BetheBloch.cxx:29
void SetAvedEdxTrack1Last8Cells(float avededxtrk1Last8Cells)
void SetAvedEdxTrack2(float avededxtrk2)
z
Definition: test.py:28
void SetZMax2(float dzmax2)
Definition: NumuSandObj.cxx:98
size_type size() const
Definition: PtrVector.h:308
void SetMissingE(float misse)
total E of true particles leaving detector (>1 MeV, not neutron)
Definition: NumuSandObj.cxx:62
OStream cout
Definition: OStream.cxx:6
void SetYMin2(float dymin2)
Definition: NumuSandObj.cxx:74
numusand::NumuSandFxs * fNumuSandFxs
art::Ptr< rb::CellHit > Cell(geo::View_t view, unsigned int viewIdx) const
Get the ith cell from view view.
Definition: Cluster.cxx:145
std::vector< TrackIDE > HitsToTrackIDE(const std::vector< const rb::CellHit * > &hits, bool useBirksE=false) const
Returns vector of TrackIDE structs contributing to the given collection of hits.
void SetNDropped(int ndropped)
A rawdata::RawDigit with channel information decoded.
Definition: CellHit.h:27
float GeV() const
Definition: RecoHit.cxx:69
void getAveTrackdEdxLast6Cells(const art::Ptr< rb::Track > track, float &avededxlast6)
void SetXMax2(float dxmax2)
maximum position of hit second farthest from edge in that dimension
Definition: NumuSandObj.cxx:86
bool getByLabel(std::string const &label, std::string const &productInstanceName, Handle< PROD > &result) const
Definition: DataViewImpl.h:344
float X() const
Definition: RecoHit.h:36
MaybeLogger_< ELseverityLevel::ELsev_warning, false > LogWarning
void geom(int which=0)
Definition: geom.C:163
bool HaveTruthInfo() const
Is this a file with truth info in? (Is BackTracker going to be any use to you?)
Definition: BackTracker.h:133
double TotalGeV(EEnergyCalcScheme escheme=kRecomputeEnergy) const
Simple sum of the estimated GeV of all the hits.
Definition: Cluster.cxx:378
void SetActVtx(float startact)
proxy for hadE[all hits not in main track] - calE[2nd reco trk]
Double_t ymin
Definition: plot.C:24
void SetScattAngleTrack2(float scattangletrk2)
bool getMissingE(const art::Ptr< rb::Cluster > slice, float &missE1)
void SetZMin2(float dzmin2)
Definition: NumuSandObj.cxx:80
float Y() const
Definition: RecoHit.h:37
void SetNBad(int nbad)
T min(const caf::Proxy< T > &a, T b)
virtual bool Is3D() const
Definition: Prong.h:71
void SetVtxE60(float vtxE60)
void SetNProtons(int nprotons)
void SetVtxE40(float vtxE40)
std::vector< const sim::Particle * > HitsToParticle(const std::vector< const rb::CellHit * > &hits) const
Returns vector of sim::Particle objects contributing to the given collection of hits.
void art::Consumer::showMissingConsumes ( ) const
protectedinherited

Referenced by art::RootOutput::endJob().

void art::Consumer::validateConsumedProduct ( BranchType const  bt,
ProductInfo const &  pi 
)
protectedinherited

Member Data Documentation

std::string numusand::FillSandbox::fEnergyLabel
private

Definition at line 50 of file FillSandbox_module.cc.

Referenced by produce().

numusand::NumuSandFxs* numusand::FillSandbox::fNumuSandFxs
private

Definition at line 47 of file FillSandbox_module.cc.

Referenced by beginRun(), and produce().

std::string numusand::FillSandbox::fPIDLabel
private

Definition at line 51 of file FillSandbox_module.cc.

Referenced by produce().

std::string numusand::FillSandbox::fSlicerLabel
private

Definition at line 48 of file FillSandbox_module.cc.

Referenced by produce().

std::string numusand::FillSandbox::fTrackLabel
private

Definition at line 49 of file FillSandbox_module.cc.

Referenced by produce().

std::string numusand::FillSandbox::fVtxLabel
private

Definition at line 52 of file FillSandbox_module.cc.

Referenced by produce().


The documentation for this class was generated from the following file: