LArSoft  v06_85_00
Liquid Argon Software toolkit - http://larsoft.org/
vertex::HarrisVertexFinder Class Reference
Inheritance diagram for vertex::HarrisVertexFinder:
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

 HarrisVertexFinder (fhicl::ParameterSet const &pset)
 
virtual ~HarrisVertexFinder ()
 
void beginJob ()
 
void produce (art::Event &evt)
 
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< T > consumes (InputTag const &)
 
template<typename T , art::BranchType BT>
art::ProductToken< T > consumes (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< T > consumesView (InputTag const &it)
 
template<typename T , BranchType = InEvent>
ProductToken< T > mayConsume (InputTag const &)
 
template<typename T , art::BranchType BT>
art::ProductToken< T > mayConsume (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< T > mayConsumeView (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 Member Functions

double Gaussian (int x, int y, double sigma)
 
double GaussianDerivativeX (int x, int y)
 
double GaussianDerivativeY (int x, int y)
 
void SaveBMPFile (const char *fileName, unsigned char *pix, int dx, int dy)
 

Private Attributes

std::string fDBScanModuleLabel
 
int fTimeBins
 
int fMaxCorners
 
double fGsigma
 
int fWindow
 
double fThreshold
 
int fSaveVertexMap
 
TH2F * fNoVertices
 

Detailed Description

Definition at line 64 of file HarrisVertexFinder_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

vertex::HarrisVertexFinder::HarrisVertexFinder ( fhicl::ParameterSet const &  pset)
explicit

Definition at line 95 of file HarrisVertexFinder_module.cc.

96  : fDBScanModuleLabel (pset.get< std::string >("DBScanModuleLabel"))
97  , fTimeBins (pset.get< int >("TimeBins") )
98  , fMaxCorners (pset.get< int >("MaxCorners") )
99  , fGsigma (pset.get< double >("Gsigma") )
100  , fWindow (pset.get< int >("Window") )
101  , fThreshold (pset.get< double >("Threshold") )
102  , fSaveVertexMap (pset.get< int >("SaveVertexMap") )
103 {
104  produces< std::vector<recob::EndPoint2D> >();
105  produces< art::Assns<recob::EndPoint2D, recob::Hit> >();
106 }
vertex::HarrisVertexFinder::~HarrisVertexFinder ( )
virtual

Definition at line 109 of file HarrisVertexFinder_module.cc.

110 {
111 }

Member Function Documentation

void vertex::HarrisVertexFinder::beginJob ( )
virtual

Reimplemented from art::EDProducer.

Definition at line 113 of file HarrisVertexFinder_module.cc.

References fNoVertices, and art::TFileDirectory::make().

113  {
114  // get access to the TFile service
116 
117  fNoVertices= tfs->make<TH2F>("fNoVertices", ";Event No; No of vertices", 100,0, 100, 30, 0, 30);
118 
119 }
T * make(ARGS...args) const
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 147 of file Consumer.h.

References art::InputTag::instance(), art::InputTag::label(), and art::InputTag::process().

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

Definition at line 162 of file Consumer.h.

163 {
164  if (!moduleContext_)
165  return;
166 
167  consumables_[BT].emplace_back(ConsumableType::Many, TypeID{typeid(T)});
168 }
ConsumableProducts consumables_
Definition: Consumer.h:138
bool moduleContext_
Definition: Consumer.h:136
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 172 of file Consumer.h.

References art::InputTag::instance(), art::InputTag::label(), and art::InputTag::process().

173 {
174  if (!moduleContext_)
175  return ViewToken<T>::invalid();
176 
177  consumables_[BT].emplace_back(ConsumableType::ViewElement,
178  TypeID{typeid(T)},
179  it.label(),
180  it.instance(),
181  it.process());
182  return ViewToken<T>{it};
183 }
static ViewToken< Element > invalid()
Definition: ProductToken.h:75
ConsumableProducts consumables_
Definition: Consumer.h:138
bool moduleContext_
Definition: Consumer.h:136
EngineCreator::base_engine_t & EngineCreator::createEngine ( seed_t  seed,
std::string const &  kind_of_engine_to_make 
)
inherited

Definition at line 32 of file EngineCreator.cc.

References art::EngineCreator::rng().

34 {
35  return rng()->createEngine(
36  placeholder_schedule_id(), seed, kind_of_engine_to_make);
37 }
long seed
Definition: chem4.cc:68
static art::ServiceHandle< art::RandomNumberGenerator > & rng()
EngineCreator::base_engine_t & EngineCreator::createEngine ( seed_t  seed,
std::string const &  kind_of_engine_to_make,
label_t const &  engine_label 
)
inherited

Definition at line 40 of file EngineCreator.cc.

References art::EngineCreator::rng().

43 {
44  return rng()->createEngine(
45  placeholder_schedule_id(), seed, kind_of_engine_to_make, engine_label);
46 }
long seed
Definition: chem4.cc:68
static art::ServiceHandle< art::RandomNumberGenerator > & rng()
CurrentProcessingContext const * art::EDProducer::currentContext ( ) const
protectedinherited

Definition at line 120 of file EDProducer.cc.

References art::EDProducer::current_context_.

121  {
122  return current_context_.get();
123  }
CPC_exempt_ptr current_context_
Definition: EDProducer.h:116
double vertex::HarrisVertexFinder::Gaussian ( int  x,
int  y,
double  sigma 
)
private

Definition at line 122 of file HarrisVertexFinder_module.cc.

References fhicl::detail::atom::value().

Referenced by produce().

123 {
124  double Norm = 1./std::sqrt(2*TMath::Pi()*pow(sigma,2));
125  double value = Norm*exp(-(pow(x,2)+pow(y,2))/(2*pow(sigma,2)));
126  return value;
127 }
Float_t x
Definition: compare.C:6
Float_t y
Definition: compare.C:6
std::string value(boost::any const &)
double vertex::HarrisVertexFinder::GaussianDerivativeX ( int  x,
int  y 
)
private

Definition at line 130 of file HarrisVertexFinder_module.cc.

References fGsigma, fhicl::detail::atom::value(), and x.

Referenced by produce().

131 {
132  double Norm = 1./(std::sqrt(2*TMath::Pi())*pow(fGsigma,3));
133  double value = Norm*(-x)*exp(-(pow(x,2)+pow(y,2))/(2*pow(fGsigma,2)));
134  return value;
135 }
Float_t x
Definition: compare.C:6
Float_t y
Definition: compare.C:6
std::string value(boost::any const &)
double vertex::HarrisVertexFinder::GaussianDerivativeY ( int  x,
int  y 
)
private

Definition at line 138 of file HarrisVertexFinder_module.cc.

References fGsigma, fhicl::detail::atom::value(), and y.

Referenced by produce().

139 {
140  double Norm = 1./(std::sqrt(2*TMath::Pi())*pow(fGsigma,3));
141  double value = Norm*(-y)*exp(-(pow(x,2)+pow(y,2))/(2*pow(fGsigma,2)));
142  return value;
143 }
Float_t x
Definition: compare.C:6
Float_t y
Definition: compare.C:6
std::string value(boost::any const &)
EngineCreator::seed_t EngineCreator::get_seed_value ( fhicl::ParameterSet const &  pset,
char const  key[] = "seed",
seed_t const  implicit_seed = -1 
)
inherited

Definition at line 49 of file EngineCreator.cc.

References fhicl::ParameterSet::get().

Referenced by art::MixFilter< T >::initEngine_().

52 {
53  auto const& explicit_seeds = pset.get<std::vector<int>>(key, {});
54  return explicit_seeds.empty() ? implicit_seed : explicit_seeds.front();
55 }
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_.

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 B, and 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  }
Int_t B
Definition: plot.C:25
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 190 of file Consumer.h.

References art::InputTag::instance(), art::InputTag::label(), and art::InputTag::process().

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

Definition at line 205 of file Consumer.h.

206 {
207  if (!moduleContext_)
208  return;
209 
210  consumables_[BT].emplace_back(ConsumableType::Many, TypeID{typeid(T)});
211 }
ConsumableProducts consumables_
Definition: Consumer.h:138
bool moduleContext_
Definition: Consumer.h:136
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 215 of file Consumer.h.

References art::InputTag::instance(), art::InputTag::label(), and art::InputTag::process().

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

Definition at line 40 of file ProducerBase.h.

References art::ProducerBase::getProductID().

41  {
42  return true;
43  }
void art::Consumer::prepareForJob ( fhicl::ParameterSet const &  pset)
protectedinherited

Definition at line 89 of file Consumer.cc.

References fhicl::ParameterSet::get_if_present().

Referenced by art::EDProducer::doBeginJob(), art::EDFilter::doBeginJob(), and art::EDAnalyzer::doBeginJob().

90 {
91  if (!moduleContext_)
92  return;
93 
94  pset.get_if_present("errorOnMissingConsumes", requireConsumes_);
95  for (auto& consumablesPerBranch : consumables_) {
96  cet::sort_all(consumablesPerBranch);
97  }
98 }
bool requireConsumes_
Definition: Consumer.h:137
ConsumableProducts consumables_
Definition: Consumer.h:138
bool moduleContext_
Definition: Consumer.h:136
void vertex::HarrisVertexFinder::produce ( art::Event evt)
virtual
Todo:
this is completely specific to ArgoNeuT!!!!

Implements art::EDProducer.

Definition at line 146 of file HarrisVertexFinder_module.cc.

References art::PtrVector< T >::clear(), util::CreateAssn(), geo::GeometryCore::Cryostat(), art::EventID::event(), fDBScanModuleLabel, fGsigma, fMaxCorners, fNoVertices, fSaveVertexMap, fThreshold, fTimeBins, fWindow, Gaussian(), GaussianDerivativeX(), GaussianDerivativeY(), art::DataViewImpl::getByLabel(), art::Event::id(), geo::GeometryCore::IteratePlaneIDs(), max, min, n, geo::PlaneGeo::Nwires(), geo::TPCGeo::Plane(), art::PtrVector< T >::push_back(), art::Event::put(), SaveBMPFile(), art::PtrVector< T >::size(), geo::CryostatGeo::TPC(), geo::GeometryCore::View(), w, x, and y.

147 {
148 
150 
151  art::Handle< std::vector<recob::Cluster> > clusterListHandle;
152  evt.getByLabel(fDBScanModuleLabel, clusterListHandle);
153 
154  //Point to a collection of vertices to output and the associations with the hits
155  std::unique_ptr<std::vector<recob::EndPoint2D> > vtxcol(new std::vector<recob::EndPoint2D>);
156  std::unique_ptr< art::Assns<recob::EndPoint2D, recob::Hit> > assn(new art::Assns<recob::EndPoint2D, recob::Hit>);
157 
158  std::vector< art::Ptr<recob::Hit> > cHits;
160 
162  for(unsigned int ii = 0; ii < clusterListHandle->size(); ++ii){
163  art::Ptr<recob::Cluster> cluster(clusterListHandle, ii);
164  clusIn.push_back(cluster);
165  }
166 
167  art::FindManyP<recob::Hit> fmh(clusterListHandle, evt, fDBScanModuleLabel);
168 
169  int flag = 0;
170  int n = 0; //index of window cell. There are 49 cells in the 7X7 Gaussian and Gaussian derivative windows
171  unsigned int numberwires = 0;
172  double MatrixAAsum = 0.;
173  double MatrixBBsum = 0.;
174  double MatrixCCsum = 0.;
175  std::vector<double> Cornerness2;
176  //gaussian window definitions. The cell weights are calculated here to help the algorithm's speed
177  double w[49] = {0.};
178  double wx[49] = {0.};
179  double wy[49] = {0.};
180  int ctr = 0;
181  for(int i = -3; i < 4; ++i){
182  for(int j = 3; j > -4; --j){
183  w[ctr] = Gaussian(i, j, fGsigma);
184  wx[ctr] = GaussianDerivativeX(i,j);
185  wy[ctr] = GaussianDerivativeY(i,j);
186  ++ctr;
187  }
188  }
189 
190  const unsigned int numbertimesamples
191  = lar::providerFrom<detinfo::DetectorPropertiesService>()->ReadOutWindowSize();
192 
193  const float BinsPerTick = fTimeBins / numbertimesamples;
194  const float TicksPerBin = numbertimesamples / fTimeBins;
195 
196  for(auto pid : geom->IteratePlaneIDs()){
198  geo::View_t view = geom->View(pid);
199  hit.clear();
200  cHits.clear();
201  for(size_t ci = 0; ci < clusIn.size(); ++ci) {
202  cHits = fmh.at(ci);
203  if(cHits.size() > 0)
204  for(unsigned int i = 0; i < cHits.size(); ++i) hit.push_back(cHits[i]);
205 
206  }
207  if(hit.size() == 0) continue;
208 
209  numberwires = geom->Cryostat(pid.Cryostat).TPC(pid.TPC).Plane(pid.Plane).Nwires();
210  std::vector< std::vector<double> > MatrixAsum(numberwires);
211  std::vector< std::vector<double> > MatrixBsum(numberwires);
212  std::vector< std::vector<double> > hit_map(numberwires);//the map of hits
213  std::vector< std::vector<int> > hit_loc(numberwires);//the index of the hit that corresponds to the potential corner
214  std::vector< std::vector<double> > Cornerness(numberwires);//the "weight" of a corner
215 
216  for(unsigned int wi = 0; wi < numberwires; ++wi){
217  MatrixAsum[wi].resize(fTimeBins);
218  MatrixBsum[wi].resize(fTimeBins);
219  hit_map[wi].resize(fTimeBins);
220  hit_loc[wi].resize(fTimeBins);
221  Cornerness[wi].resize(fTimeBins);
222  for(int timebin = 0; timebin < fTimeBins; ++timebin){
223  hit_map[wi][timebin] = 0.;
224  hit_loc[wi][timebin] = -1;
225  Cornerness[wi][timebin] = 0.;
226  MatrixAsum[wi][timebin] = 0.;
227  MatrixBsum[wi][timebin] = 0.;
228  }
229  }
230 
231  for(unsigned int i = 0; i < hit.size(); ++i){
232  unsigned int wire = hit[i]->WireID().Wire;
233  // pixelization using a Gaussian (3 sigmas around peak)
234  const float center = hit[i]->PeakTime(), sigma = hit[i]->RMS();
235  const int iFirstBin = int((center - 3*sigma) / TicksPerBin),
236  iLastBin = int((center + 3*sigma) / TicksPerBin);
237  for (int iBin = iFirstBin; iBin <= iLastBin; ++iBin) {
238  const float bin_center = iBin * TicksPerBin;
239  hit_map[wire][iBin] += Gaussian(bin_center, center, sigma);
240  }
241  }
242 
243  // Gaussian derivative convolution
244  for(unsigned int wire = 1;wire < numberwires-1; ++wire){
245  for(int timebin = 1;timebin < fTimeBins-1; ++timebin){
246  MatrixAsum[wire][timebin] = 0.;
247  MatrixBsum[wire][timebin] = 0.;
248  n = 0;
249  for(int i = -3; i <= 3; ++i) {
250  unsigned int windex = ((int)wire < -i)? 0: wire + i;
251  if (windex >= numberwires) windex = numberwires - 1;
252  for(int j = -3; j <= 3; ++j) {
253  unsigned int tindex = (unsigned int) std::min(std::max(timebin +j, 0), fTimeBins - 1);
254 
255  MatrixAsum[wire][timebin] += wx[n]*hit_map[windex][tindex];
256  MatrixBsum[wire][timebin] += wy[n]*hit_map[windex][tindex];
257  ++n;
258  }// end loop over j
259  }// end loop over i
260  }// end loop over time bins
261  }// end loop over wires
262 
263  //calculate the cornerness of each pixel while making sure not to fall off the hit map.
264  for(unsigned int wire = 1; wire < numberwires-1; ++wire){
265  for(int timebin = 1; timebin < fTimeBins-1; ++timebin){
266  MatrixAAsum = 0.;
267  MatrixBBsum = 0.;
268  MatrixCCsum = 0.;
269  //Gaussian smoothing convolution
270  n = 0;
271 
272  for(int i = -3; i <= 3; i++) {
273  unsigned int windex = ((int)wire < -i)? 0: wire + i;
274  if (windex >= numberwires) windex = numberwires - 1;
275  for(int j = -3; j <= 3; j++) {
276  unsigned int tindex = (unsigned int) std::min(std::max(timebin +j, 0), fTimeBins - 1);
277  MatrixAAsum += w[n]*pow(MatrixAsum[windex][tindex],2);
278  MatrixBBsum += w[n]*pow(MatrixBsum[windex][tindex],2);
279  MatrixCCsum += w[n]*MatrixAsum[windex][tindex]*MatrixBsum[windex][tindex];
280  ++n;
281  }// end loop over j
282  }// end loop over i
283 
284  if((MatrixAAsum+MatrixBBsum) > 0)
285  Cornerness[wire][timebin] = (MatrixAAsum*MatrixBBsum-pow(MatrixCCsum,2))/(MatrixAAsum+MatrixBBsum);
286  else
287  Cornerness[wire][timebin] = 0;
288 
289  if(Cornerness[wire][timebin]>0){
290  for(unsigned int i = 0;i < hit.size(); ++i){
291  unsigned int wire2 = hit[i]->WireID().Wire;
292  // make sure the vertex candidate coincides with an actual hit
293  // (time within one RMS from the center of the hit)
294  const float time = timebin / TicksPerBin;
295  if(wire == wire2 && std::abs(hit[i]->TimeDistanceAsRMS(time)) < 1.){
296  //this index keeps track of the hit number
297  hit_loc[wire][timebin] = i;
298  Cornerness2.push_back(Cornerness[wire][timebin]);
299  break;
300  }
301  }// end loop over hits
302  }// end if cornerness > 0
303  }// end loop over time bins
304  }// end loop over wires
305 
306  std::sort(Cornerness2.rbegin(),Cornerness2.rend());
307 
308  for(int vertexnum = 0; vertexnum < fMaxCorners; ++vertexnum){
309  flag = 0;
310  for(unsigned int wire = 0;wire < numberwires && flag == 0; ++wire){
311  for(int timebin = 0;timebin < fTimeBins && flag == 0; ++timebin){
312  if(Cornerness2.size() > (unsigned int)vertexnum)
313  if(Cornerness[wire][timebin] == Cornerness2[vertexnum]
314  && Cornerness[wire][timebin] > 0.
315  && hit_loc[wire][timebin]>-1){
316  ++flag;
317  //thresholding
318  if(Cornerness2.size())
319  if(Cornerness[wire][timebin] < (fThreshold*Cornerness2[0]))
320  vertexnum = fMaxCorners;
321 
322  vHits.push_back(hit[hit_loc[wire][timebin]]);
323  // get the total charge from the associated hits
324  double totalQ = 0.;
325  for(size_t vh = 0; vh < vHits.size(); ++vh) totalQ += vHits[vh]->Integral();
326 
327  recob::EndPoint2D vertex(hit[hit_loc[wire][timebin]]->PeakTime(),
328  hit[hit_loc[wire][timebin]]->WireID(), //for update to EndPoint2D ... WK 4/22/13
329  Cornerness[wire][timebin],
330  vtxcol->size(),
331  view,
332  totalQ);
333  vtxcol->push_back(vertex);
334 
335  util::CreateAssn(*this, evt, *(vtxcol.get()), vHits, *(assn.get()));
336 
337  vHits.clear();
338  // non-maximal suppression on a square window. The wire coordinate units are
339  // converted to time ticks so that the window is truly square.
340  // Note that there are 1/0.0743=13.46 time samples per 4.0 mm (wire pitch in ArgoNeuT),
341  // assuming a 1.5 mm/us drift velocity for a 500 V/cm E-field
342 
344 
345  for(unsigned int wireout=wire-(int)((fWindow*BinsPerTick*.0743)+.5);
346  wireout <= wire+(int)((fWindow*BinsPerTick*.0743)+.5) ; wireout++)
347  for(int timebinout=timebin-fWindow;timebinout <= timebin+fWindow; timebinout++)
348  if(std::sqrt(pow(wire-wireout,2)+pow(timebin-timebinout,2))<fWindow)//circular window
349  Cornerness[wireout][timebinout]=0;
350  }
351  }// end loop over time bins
352  }// end loop over wires
353  }// end loop over vertices
354  Cornerness2.clear();
355  hit.clear();
356 
357  if(pid.Plane == (unsigned int)fSaveVertexMap){
358  unsigned char *outPix = new unsigned char [fTimeBins*numberwires];
359  //finds the maximum cell in the map for image scaling
360  int cell, pix=0, maxCell=0;
361  for (int y = 0; y < fTimeBins; ++y){
362  for (unsigned int x = 0; x < numberwires; ++x){
363  cell = (int)(hit_map[x][y]*1000);
364  if (cell > maxCell){
365  maxCell = cell;
366  }
367  }
368  }
369 
370  for (unsigned int y = 0; y < (unsigned int)fTimeBins; ++y){
371  for (unsigned int x = 0; x < numberwires; ++x){
372  //scales the pixel weights based on the maximum cell value
373  if(maxCell>0)
374  pix = (int)((1500000*hit_map[x][y])/maxCell);
375  outPix[y*numberwires + x] = pix;
376  }
377  }
378  SaveBMPFile("harrisvertexmap.bmp", outPix, numberwires, fTimeBins);
379  delete [] outPix;
380  }
381  }// end loop over planes
382 
383  mf::LogInfo("HarrisVertexFinder") << "Size of vtxcol= " << vtxcol->size();
384  fNoVertices->Fill(evt.id().event(),vtxcol->size());
385 
386  evt.put(std::move(vtxcol));
387 }
Float_t x
Definition: compare.C:6
enum geo::_plane_proj View_t
Enumerate the possible plane projections.
MaybeLogger_< ELseverityLevel::ELsev_info, false > LogInfo
Float_t y
Definition: compare.C:6
void SaveBMPFile(const char *fileName, unsigned char *pix, int dx, int dy)
Cluster finding and building.
IteratorBox< plane_id_iterator,&GeometryCore::begin_plane_id,&GeometryCore::end_plane_id > IteratePlaneIDs() const
Enables ranged-for loops on all plane IDs of the detector.
ProductID put(std::unique_ptr< PROD > &&product)
Definition: Event.h:102
double Gaussian(int x, int y, double sigma)
Int_t max
Definition: plot.C:27
void push_back(Ptr< U > const &p)
Definition: PtrVector.h:441
CryostatGeo const & Cryostat(geo::CryostatID const &cryoid) const
Returns the specified cryostat.
bool CreateAssn(PRODUCER const &prod, art::Event &evt, std::vector< T > const &a, art::Ptr< U > const &b, art::Assns< U, T > &assn, std::string a_instance, size_t indx=UINT_MAX)
Creates a single one-to-one association.
View_t View(geo::PlaneID const &pid) const
Returns the view (wire orientation) on the channels of specified TPC plane.
size_type size() const
Definition: PtrVector.h:308
Detector simulation of raw signals on wires.
bool getByLabel(std::string const &label, std::string const &productInstanceName, Handle< PROD > &result) const
Definition: DataViewImpl.h:344
const TPCGeo & TPC(unsigned int itpc) const
Return the itpc&#39;th TPC in the cryostat.
Int_t min
Definition: plot.C:26
unsigned int Nwires() const
Number of wires in this plane.
Definition: PlaneGeo.h:250
EventNumber_t event() const
Definition: EventID.h:117
Char_t n[5]
PlaneGeo const & Plane(geo::View_t view) const
Return the plane in the tpc with View_t view.
Definition: TPCGeo.cxx:298
void clear()
Definition: PtrVector.h:537
Float_t w
Definition: plot.C:23
EventID id() const
Definition: Event.h:56
vertex reconstruction
void vertex::HarrisVertexFinder::SaveBMPFile ( const char *  fileName,
unsigned char *  pix,
int  dx,
int  dy 
)
private

Definition at line 391 of file HarrisVertexFinder_module.cc.

References DEFINE_ART_MODULE.

Referenced by produce().

392 {
393  std::ofstream bmpFile(fileName, std::ios::binary);
394  bmpFile.write("B", 1);
395  bmpFile.write("M", 1);
396  int bitsOffset = 54 +256*4;
397  int size = bitsOffset + dx*dy; //header plus 256 entry LUT plus pixels
398  bmpFile.write((const char *)&size, 4);
399  int reserved = 0;
400  bmpFile.write((const char *)&reserved, 4);
401  bmpFile.write((const char *)&bitsOffset, 4);
402  int bmiSize = 40;
403  bmpFile.write((const char *)&bmiSize, 4);
404  bmpFile.write((const char *)&dx, 4);
405  bmpFile.write((const char *)&dy, 4);
406  short planes = 1;
407  bmpFile.write((const char *)&planes, 2);
408  short bitCount = 8;
409  bmpFile.write((const char *)&bitCount, 2);
410  int i, temp = 0;
411  for (i=0; i<6; i++)
412  bmpFile.write((const char *)&temp, 4); // zero out optional color info
413  // write a linear LUT
414  char lutEntry[4]; // blue,green,red
415  lutEntry[3] = 0; // reserved part
416  for (i=0; i<256; i++)
417  {
418  lutEntry[0] = lutEntry[1] = lutEntry[2] = i;
419  bmpFile.write(lutEntry, sizeof lutEntry);
420  }
421  // write the actual pixels
422  bmpFile.write((const char *)pix, dx*dy);
423 }
void art::Consumer::showMissingConsumes ( ) const
protectedinherited

Definition at line 125 of file Consumer.cc.

Referenced by art::EDProducer::doEndJob(), art::EDFilter::doEndJob(), art::EDAnalyzer::doEndJob(), and art::RootOutput::endJob().

126 {
127  if (!moduleContext_)
128  return;
129 
130  // If none of the branches have missing consumes statements, exit early.
131  if (std::all_of(cbegin(missingConsumes_),
132  cend(missingConsumes_),
133  [](auto const& perBranch) { return perBranch.empty(); }))
134  return;
135 
136  constexpr cet::HorizontalRule rule{60};
137  mf::LogPrint log{"MTdiagnostics"};
138  log << '\n'
139  << rule('=') << '\n'
140  << "The following consumes (or mayConsume) statements are missing from\n"
141  << module_context(moduleDescription_) << '\n'
142  << rule('-') << '\n';
143 
144  cet::for_all_with_index(
145  missingConsumes_, [&log](std::size_t const i, auto const& perBranch) {
146  for (auto const& pi : perBranch) {
147  log << " "
148  << assemble_consumes_statement(static_cast<BranchType>(i), pi)
149  << '\n';
150  }
151  });
152  log << rule('=');
153 }
cet::exempt_ptr< ModuleDescription const > moduleDescription_
Definition: Consumer.h:140
constexpr T pi()
Returns the constant pi (up to 35 decimal digits of precision)
bool moduleContext_
Definition: Consumer.h:136
ConsumableProductSets missingConsumes_
Definition: Consumer.h:139
void art::Consumer::validateConsumedProduct ( BranchType const  bt,
ProductInfo const &  pi 
)
protectedinherited

Definition at line 101 of file Consumer.cc.

References art::errors::ProductRegistrationFailure.

103 {
104  // Early exits if consumes tracking has been disabled or if the
105  // consumed product is an allowed consumable.
106  if (!moduleContext_)
107  return;
108 
109  if (cet::binary_search_all(consumables_[bt], pi))
110  return;
111 
112  if (requireConsumes_) {
114  "Consumer: an error occurred during validation of a "
115  "retrieved product\n\n")
116  << "The following consumes (or mayConsume) statement is missing from\n"
117  << module_context(moduleDescription_) << ":\n\n"
118  << " " << assemble_consumes_statement(bt, pi) << "\n\n";
119  }
120 
121  missingConsumes_[bt].insert(pi);
122 }
cet::exempt_ptr< ModuleDescription const > moduleDescription_
Definition: Consumer.h:140
bool requireConsumes_
Definition: Consumer.h:137
cet::coded_exception< errors::ErrorCodes, ExceptionDetail::translate > Exception
Definition: Exception.h:66
constexpr T pi()
Returns the constant pi (up to 35 decimal digits of precision)
ConsumableProducts consumables_
Definition: Consumer.h:138
bool moduleContext_
Definition: Consumer.h:136
ConsumableProductSets missingConsumes_
Definition: Consumer.h:139

Member Data Documentation

std::string vertex::HarrisVertexFinder::fDBScanModuleLabel
private

Definition at line 81 of file HarrisVertexFinder_module.cc.

Referenced by produce().

double vertex::HarrisVertexFinder::fGsigma
private

Definition at line 84 of file HarrisVertexFinder_module.cc.

Referenced by GaussianDerivativeX(), GaussianDerivativeY(), and produce().

int vertex::HarrisVertexFinder::fMaxCorners
private

Definition at line 83 of file HarrisVertexFinder_module.cc.

Referenced by produce().

TH2F* vertex::HarrisVertexFinder::fNoVertices
private

Definition at line 88 of file HarrisVertexFinder_module.cc.

Referenced by beginJob(), and produce().

int vertex::HarrisVertexFinder::fSaveVertexMap
private

Definition at line 87 of file HarrisVertexFinder_module.cc.

Referenced by produce().

double vertex::HarrisVertexFinder::fThreshold
private

Definition at line 86 of file HarrisVertexFinder_module.cc.

Referenced by produce().

int vertex::HarrisVertexFinder::fTimeBins
private

Definition at line 82 of file HarrisVertexFinder_module.cc.

Referenced by produce().

int vertex::HarrisVertexFinder::fWindow
private

Definition at line 85 of file HarrisVertexFinder_module.cc.

Referenced by produce().


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