All Classes Namespaces Files Functions Variables Typedefs Enumerations Enumerator Friends Macros Groups Pages
Public Member Functions | Private Member Functions | Private Attributes | List of all members
hit::RawHitFinder Class Reference
Inheritance diagram for hit::RawHitFinder:

Public Member Functions

 RawHitFinder (fhicl::ParameterSet const &pset)
 

Private Member Functions

void produce (art::Event &evt) override
 

Private Attributes

unsigned int fDataSize
 
art::InputTag fDigitModuleLabel
 
std::string fSpillName
 
std::string fCalDataModuleLabel
 
std::string fHitLabelName
 
double fMinSigInd
 
double fMinSigCol
 
double fIndWidth
 
double fColWidth
 
double fIndMinWidth
 
double fColMinWidth
 
int fMaxMultiHit
 
int fAreaMethod
 
std::vector< double > fAreaNorms
 
bool fUncompressWithPed
 
bool fSkipInd
 
double fIncludeMoreTail
 
int fColMinWindow
 
int fIndCutoff
 

Detailed Description

Definition at line 36 of file RawHitFinder_module.cc.

Constructor & Destructor Documentation

hit::RawHitFinder::RawHitFinder ( fhicl::ParameterSet const &  pset)
explicit

Definition at line 73 of file RawHitFinder_module.cc.

74  : EDProducer{pset}
75  {
76  fDigitModuleLabel = pset.get< art::InputTag >("DigitModuleLabel", "daq");
77  fCalDataModuleLabel = pset.get< std::string >("CalDataModuleLabel");
78  fMinSigInd = pset.get< double >("MinSigInd");
79  fMinSigCol = pset.get< double >("MinSigCol");
80  fIncludeMoreTail = pset.get< double >("IncludeMoreTail", 0.);
81  fIndWidth = pset.get< int >("IndWidth",20);
82  fColWidth = pset.get< double >("ColWidth");
83  fIndMinWidth = pset.get< double >("IndMinWidth");
84  fColMinWidth = pset.get< double >("ColMinWidth",0.);
85  fMaxMultiHit = pset.get< int >("MaxMultiHit");
86  fAreaMethod = pset.get< int >("AreaMethod");
87  fAreaNorms = pset.get< std::vector< double > >("AreaNorms");
88  fUncompressWithPed = pset.get< bool >("UncompressWithPed", true);
89  fSkipInd = pset.get< bool >("SkipInd", false);
90  fColMinWindow = pset.get< int >("ColMinWindow",0);
91  fIndCutoff = pset.get< int >("IndCutoff",20);
92  mf::LogInfo("RawHitFinder_module") << "fDigitModuleLabel: " << fDigitModuleLabel << std::endl;
93 
94  //LET HITCOLLECTIONCREATOR DECLARE THAT WE ARE GOING TO PRODUCE
95  //HITS AND ASSOCIATIONS TO RAW DIGITS BUT NOT ASSOCIATIONS TO WIRES
96  //(WITH NO PARTICULAR PRODUCT LABEL).
98  /*instance_name*/"",
99  /*doWireAssns*/false,
100  /*doRawDigitAssns*/true);
101  }
std::vector< double > fAreaNorms
static void declare_products(art::ProducesCollector &collector, std::string instance_name="", bool doWireAssns=true, bool doRawDigitAssns=true)
Declares the hit products we are going to fill.
Definition: HitCreator.cxx:248
std::string fCalDataModuleLabel
art::InputTag fDigitModuleLabel

Member Function Documentation

void hit::RawHitFinder::produce ( art::Event &  evt)
overrideprivate

Definition at line 104 of file RawHitFinder_module.cc.

105  {
106  //GET THE GEOMETRY.
107  art::ServiceHandle<geo::Geometry const> geom;
108 
109  //READ IN THE DIGIT LIST OBJECTS(S).
110  art::Handle< std::vector<raw::RawDigit> > digitVecHandle;
111 
112  bool retVal = evt.getByLabel(fDigitModuleLabel, digitVecHandle);
113  if(retVal==true)
114  mf::LogInfo("RawHitFinder_module") << "I got fDigitModuleLabel: " << fDigitModuleLabel << std::endl;
115  else
116  mf::LogWarning("RawHitFinder_module") << "Could not get fDigitModuleLabel: " << fDigitModuleLabel << std::endl;
117 
118  std::vector<float> holder; //HOLDS SIGNAL DATA.
119  std::vector<short> rawadc; //UNCOMPRESSED ADC VALUES.
120 
121  // ###############################################
122  // ### Making a ptr vector to put on the event ###
123  // ###############################################
124  // THIS CONTAINS THE HIT COLLECTION AND ITS ASSOCIATIONS TO WIRES AND RAW DIGITS.
125  recob::HitCollectionCreator hcol(evt, false /* doWireAssns */, true /* doRawDigitAssns */);
126 
127  std::vector<float> startTimes; //STORES TIME OF WINDOW START.
128  std::vector<float> maxTimes; //STORES TIME OF LOCAL MAXIMUM.
129  std::vector<float> endTimes; //STORES TIME OF WINDOW END.
130  std::vector<float> peakHeight; //STORES ADC COUNT AT THE MAXIMUM.
131  std::vector<float> hitrms; //STORES CHARGE WEIGHTED RMS OF TIME ACROSS THE HIT.
132  std::vector<double> charge; //STORES THE TOTAL CHARGE ASSOCIATED WITH THE HIT.
133 
134  uint32_t channel = 0; //CHANNEL NUMBER.
135  double threshold = 0; //MINIMUM SIGNAL SIZE FOR ID'ING A HIT.
136  double totSig = 0;
137  double myrms = 0;
138  double mynorm = 0;
140  std::stringstream numConv;
141 
142  hcol.reserve(digitVecHandle->size());
143  for(size_t rdIter = 0; rdIter < digitVecHandle->size(); ++rdIter){
144  holder.clear();
145 
146  //GET THE REFERENCE TO THE CURRENT raw::RawDigit.
147  art::Ptr<raw::RawDigit> digitVec(digitVecHandle, rdIter);
148  channel = digitVec->Channel();
149  fDataSize = digitVec->Samples();
150 
151  rawadc.resize(fDataSize);
152  holder.resize(fDataSize);
153 
154  //UNCOMPRESS THE DATA.
155  if (fUncompressWithPed){
156  int pedestal = (int)digitVec->GetPedestal();
157  raw::Uncompress(digitVec->ADCs(), rawadc, pedestal, digitVec->Compression());
158  }
159  else{
160  raw::Uncompress(digitVec->ADCs(), rawadc, digitVec->Compression());
161  }
162 
163  //GET THE LIST OF BAD CHANNELS.
164  lariov::ChannelStatusProvider const& channelStatus
165  = art::ServiceHandle<lariov::ChannelStatusService const>()->GetProvider();
166 
168  = channelStatus.BadChannels();
169 
170  for(unsigned int bin = 0; bin < fDataSize; ++bin){
171  holder[bin]=(rawadc[bin]-digitVec->GetPedestal());
172  }
173 
174  sigType = geom->SignalType(channel);
175 
176  peakHeight.clear();
177  endTimes.clear();
178  startTimes.clear();
179  maxTimes.clear();
180  charge.clear();
181  hitrms.clear();
182 
183  bool channelSwitch = false;
184 
185  for(auto it = BadChannels.begin(); it != BadChannels.end(); it++)
186  {
187  if(channel==*it)
188  {
189  channelSwitch = true;
190  break;
191  }
192  }
193 
194  if(channelSwitch==false)
195  {
196  // ###############################################
197  // ### Induction Planes ###
198  // ###############################################
199 
200  //THE INDUCTION PLANE METHOD HAS NOT YET BEEN MODIFIED AND TESTED FOR REAL DATA.
201  // Or for detectors without a grid plane
202  //
203  if(sigType == geo::kInduction && !fSkipInd){
204  threshold = fMinSigInd;
205  // std::cout<< "Threshold is " << threshold << std::endl;
206  // fitWidth = fIndWidth;
207  // minWidth = fIndMinWidth;
208  // continue;
209  float negthr=-1.0*threshold;
210  unsigned int bin =1;
211  float minadc=0;
212 
213  // find the dips
214  while (bin<(fDataSize-1)) { // loop over ticks
215  float thisadc = holder[bin]; float nextadc = holder[bin+1];
216  if (thisadc<negthr && nextadc < negthr) { // new region, require two ticks above threshold
217  // std::cout << "new region" << bin << " " << thisadc << std::endl;
218  // step back to find zero crossing
219  unsigned int place = bin;
220  while (thisadc<=0 && bin>0) {
221  // std::cout << bin << " " << thisadc << std::endl;
222  bin--;
223  thisadc=holder[bin];
224  }
225  float hittime = bin+thisadc/(thisadc-holder[bin+1]);
226  maxTimes.push_back(hittime);
227 
228  // step back more to find the hit start time
229  uint32_t stop;
230  if (fIndCutoff<(int)bin) {stop=bin-fIndCutoff;} else {stop=0;}
231  while (thisadc<threshold && bin>stop) {
232  // std::cout << bin << " " << thisadc << std::endl;
233  bin--;
234  thisadc=holder[bin];
235  }
236  if (bin>=2) bin-=2;
237  while (thisadc>threshold && bin>stop) {
238  // std::cout << bin << " " << thisadc << std::endl;
239  bin--;
240  thisadc=holder[bin];
241  }
242  startTimes.push_back(bin+1);
243  // now step forward from hit time to find end time, area of dip
244  bin=place;
245  thisadc=holder[bin];
246  minadc=thisadc;
247  bin++;
248  totSig = fabs(thisadc);
249  while (thisadc<negthr && bin<fDataSize) {
250  totSig += fabs(thisadc);
251  thisadc=holder[bin];
252  if (thisadc<minadc) minadc=thisadc;
253  bin++;
254  }
255  endTimes.push_back(bin-1);
256  peakHeight.push_back(-1.0*minadc);
257  charge.push_back(totSig);
258  hitrms.push_back(5.0);
259  // std::cout << "TOTAL SIGNAL INDUCTION " << totSig << " 5.0" << std::endl;
260  // std::cout << "filled end times " << bin-1 << "peak height vector size " << peakHeight.size() << std::endl;
261 
262  // don't look for a new hit until it returns to baseline
263  while (thisadc<0 && bin<fDataSize) {
264  // std::cout << bin << " " << thisadc << std::endl;
265  bin++;
266  if (bin == fDataSize) break;
267  thisadc=holder[bin];
268  }
269  } // end region
270  bin++;
271  }// loop over ticks
272  }
273 
274  // ###############################################
275  // ### Collection Plane ###
276  // ###############################################
277 
278  else if(sigType == geo::kCollection)
279  {
280  threshold = fMinSigCol;
281 
282  float madc = threshold;
283  int ibin = 0;
284  int start = 0;
285  int end = 0;
286  unsigned int bin = 0;
287 
288  while (bin<fDataSize)
289  {
290  float thisadc = holder[bin];
291  madc = threshold;
292  ibin = 0;
293 
294  if (thisadc>madc)
295  {
296  start = bin;
297 
298  if(thisadc>threshold && bin<fDataSize)
299  {
300  while (thisadc>threshold && bin<fDataSize)
301  {
302  if (thisadc>madc)
303  {
304  ibin=bin;
305  madc=thisadc;
306  }
307  bin++;
308  if (bin == fDataSize) break;
309  thisadc=holder[bin];
310  }
311  }
312  else
313  {
314  bin++;
315  }
316 
317  end = bin-1;
318 
319  if(start!=end)
320  {
321  maxTimes.push_back(ibin);
322  peakHeight.push_back(madc);
323  startTimes.push_back(start);
324  endTimes.push_back(end);
325 
326  totSig = 0;
327  myrms = 0;
328  mynorm = 0;
329 
330  int moreTail = std::ceil(fIncludeMoreTail*(end-start));
331  if (moreTail<fColMinWindow) moreTail=fColMinWindow;
332 
333  for(int i = start-moreTail; i <= end+moreTail; i++)
334  {
335  if(i<(int)(holder.size()) && i>=0)
336  {
337  float temp = ibin-i;
338  myrms += temp*temp*holder[i];
339 
340  totSig += holder[i];
341  }
342  }
343 
344  charge.push_back(totSig);
345  mynorm = totSig;
346  myrms/=mynorm;
347  hitrms.push_back(sqrt(myrms));
348 
349  //PRE CHANGES MADE 04/14/16. A BOOTH, DUNE 35T.
350  /*
351  int moreTail = std::ceil(fIncludeMoreTail*(end-start));
352 
353  for(int i = start-moreTail; i <= end+moreTail; i++)
354  {
355  totSig += holder[i];
356  float temp2 = holder[i]*holder[i];
357  mynorm += temp2;
358  float temp = ibin-i;
359  myrms += temp*temp*temp2;
360  }
361 
362  charge.push_back(totSig);
363  myrms/=mynorm;
364  if((end-start+2*moreTail+1)!=0)
365  {
366  myrms/=(float)(end-start+2*moreTail+1);
367  hitrms.push_back(sqrt(myrms));
368  }
369  else
370  {
371  hitrms.push_back(sqrt(myrms));
372  }*/
373  }
374  }
375  start = 0;
376  end = 0;
377  bin++;
378  }
379  }
380  }
381 
382  int numHits(0); //NUMBER OF CONSECUTIVE HITS BEING FITTED.
383  int hitIndex(0); //INDEX OF CURRENT HIT IN SEQUENCE.
384  double amplitude(0), position(0); //FIT PARAMETERS.
385  double start(0), end(0);
386  double amplitudeErr(0), positionErr(0); //FIT ERRORS.
387  double goodnessOfFit(0), chargeErr(0); //CHI2/NDF and error on charge.
388  double hrms(0);
389 
390  numHits = maxTimes.size();
391  for (int i = 0; i < numHits; ++i)
392  {
393  amplitude = peakHeight[i];
394  position = maxTimes[i];
395  start = startTimes[i];
396  end = endTimes[i];
397  hrms = hitrms[i];
398  amplitudeErr = -1;
399  positionErr = 1.0;
400  goodnessOfFit = -1;
401  chargeErr = -1;
402  totSig = charge[i];
403 
404 
405  std::vector<geo::WireID> wids = geom->ChannelToWire(channel);
406  geo::WireID wid = wids[0];
407 
408  if (start>=end)
409  {
410  mf::LogWarning("RawHitFinder_module") << "Hit start " << start << " is >= hit end " << end;
411  continue;
412  }
413 
415  *digitVec, //RAW DIGIT REFERENCE.
416  wid, //WIRE ID.
417  start, //START TICK.
418  end, //END TICK.
419  hrms, //RMS.
420  position, //PEAK_TIME.
421  positionErr, //SIGMA_PEAK_TIME.
422  amplitude, //PEAK_AMPLITUDE.
423  amplitudeErr, //SIGMA_PEAK_AMPLITUDE.
424  totSig, //HIT_INTEGRAL.
425  chargeErr, //HIT_SIGMA_INTEGRAL.
426  std::accumulate(holder.begin() + (int) start, holder.begin() + (int) end, 0.), //SUMMED CHARGE.
427  1, //MULTIPLICITY.
428  -1, //LOCAL_INDEX.
429  goodnessOfFit, //WIRE ID.
430  int(end-start+1) //DEGREES OF FREEDOM.
431  );
432  hcol.emplace_back(hit.move(), digitVec);
433 
434  ++hitIndex;
435  }
436  }
437 
438  hcol.put_into(evt);
439  }
std::set< raw::ChannelID_t > ChannelSet_t
Type of set of channel IDs.
process_name hit
Definition: cheaterreco.fcl:51
virtual ChannelSet_t BadChannels() const =0
Returns a copy of set of bad channel IDs for the current run.
constexpr details::BinObj< T > bin(T value)
Returns a wrapper to print the specified data in binary format.
Class managing the creation of a new recob::Hit object.
Definition: HitCreator.h:83
A class handling a collection of hits and its associations.
Definition: HitCreator.h:508
Signal from induction planes.
Definition: geo_types.h:145
enum geo::_plane_sigtype SigType_t
auto end(FixedBins< T, C > const &) noexcept
Definition: FixedBins.h:585
Class providing information about the quality of channels.
art::InputTag fDigitModuleLabel
TCEvent evt
Definition: DataStructs.cxx:8
void Uncompress(const std::vector< short > &adc, std::vector< short > &uncompressed, raw::Compress_t compress)
Uncompresses a raw data buffer.
Definition: raw.cxx:776
Signal from collection planes.
Definition: geo_types.h:146

Member Data Documentation

int hit::RawHitFinder::fAreaMethod
private

Definition at line 60 of file RawHitFinder_module.cc.

std::vector<double> hit::RawHitFinder::fAreaNorms
private

Definition at line 61 of file RawHitFinder_module.cc.

std::string hit::RawHitFinder::fCalDataModuleLabel
private

Definition at line 51 of file RawHitFinder_module.cc.

double hit::RawHitFinder::fColMinWidth
private

Definition at line 58 of file RawHitFinder_module.cc.

int hit::RawHitFinder::fColMinWindow
private

Definition at line 66 of file RawHitFinder_module.cc.

double hit::RawHitFinder::fColWidth
private

Definition at line 56 of file RawHitFinder_module.cc.

unsigned int hit::RawHitFinder::fDataSize
private

Definition at line 46 of file RawHitFinder_module.cc.

art::InputTag hit::RawHitFinder::fDigitModuleLabel
private

Definition at line 47 of file RawHitFinder_module.cc.

std::string hit::RawHitFinder::fHitLabelName
private

Definition at line 52 of file RawHitFinder_module.cc.

double hit::RawHitFinder::fIncludeMoreTail
private

Definition at line 64 of file RawHitFinder_module.cc.

int hit::RawHitFinder::fIndCutoff
private

Definition at line 67 of file RawHitFinder_module.cc.

double hit::RawHitFinder::fIndMinWidth
private

Definition at line 57 of file RawHitFinder_module.cc.

double hit::RawHitFinder::fIndWidth
private

Definition at line 55 of file RawHitFinder_module.cc.

int hit::RawHitFinder::fMaxMultiHit
private

Definition at line 59 of file RawHitFinder_module.cc.

double hit::RawHitFinder::fMinSigCol
private

Definition at line 54 of file RawHitFinder_module.cc.

double hit::RawHitFinder::fMinSigInd
private

Definition at line 53 of file RawHitFinder_module.cc.

bool hit::RawHitFinder::fSkipInd
private

Definition at line 63 of file RawHitFinder_module.cc.

std::string hit::RawHitFinder::fSpillName
private

Definition at line 48 of file RawHitFinder_module.cc.

bool hit::RawHitFinder::fUncompressWithPed
private

Definition at line 62 of file RawHitFinder_module.cc.


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