LArSoft  v07_13_02
Liquid Argon Software toolkit - http://larsoft.org/
pmtana::AlgoCFD Class Reference

#include "AlgoCFD.h"

Inheritance diagram for pmtana::AlgoCFD:
pmtana::PMTPulseRecoBase

Public Member Functions

 AlgoCFD (const std::string name="CFD")
 Default constructor. More...
 
 AlgoCFD (const fhicl::ParameterSet &pset, const std::string name="CFD")
 Alternative ctor. More...
 
virtual ~AlgoCFD ()
 Default destructor. More...
 
void Reset ()
 Implementation of AlgoCFD::reset() method. More...
 
const std::string & Name () const
 Name getter. More...
 
const bool Status () const
 Status getter. More...
 
bool Reconstruct (const pmtana::Waveform_t &, const pmtana::PedestalMean_t &, const pmtana::PedestalSigma_t &)
 
const pulse_paramGetPulse (size_t index=0) const
 
const pulse_param_arrayGetPulses () const
 A getter for the whole array of pulse_param struct object. More...
 
size_t GetNPulse () const
 A getter for the number of reconstructed pulses from the input waveform. More...
 

Protected Member Functions

bool RecoPulse (const pmtana::Waveform_t &, const pmtana::PedestalMean_t &, const pmtana::PedestalSigma_t &)
 Implementation of AlgoCFD::reco() method. More...
 
const std::map< unsigned, double > LinearZeroPointX (const std::vector< double > &trace)
 
bool Integral (const std::vector< short > &wf, double &result, size_t begin=0, size_t end=0) const
 
bool Derivative (const std::vector< short > &wf, std::vector< int32_t > &diff, size_t begin=0, size_t end=0) const
 
size_t Max (const std::vector< short > &wf, double &result, size_t begin=0, size_t end=0) const
 
size_t Min (const std::vector< short > &wf, double &result, size_t begin=0, size_t end=0) const
 

Protected Attributes

pulse_param_array _pulse_v
 A container array of pulse_param struct objects to store (possibly multiple) reconstructed pulse(s). More...
 
pulse_param _pulse
 A subject pulse_param object to be filled with the last reconstructed pulse parameters. More...
 

Private Attributes

float _F
 
int _D
 
double _peak_thresh
 
double _start_thresh
 
double _end_thresh
 

Detailed Description

This class implements threshold algorithm to AlgoCFD class.

Definition at line 29 of file AlgoCFD.h.

Constructor & Destructor Documentation

pmtana::AlgoCFD::AlgoCFD ( const std::string  name = "CFD")

Default constructor.

Definition at line 18 of file AlgoCFD.cxx.

19  : PMTPulseRecoBase(name)
20  //*********************************************************************
21  {}
PMTPulseRecoBase(const std::string name="noname")
Default constructor with fhicl parameters.
pmtana::AlgoCFD::AlgoCFD ( const fhicl::ParameterSet pset,
const std::string  name = "CFD" 
)

Alternative ctor.

Definition at line 24 of file AlgoCFD.cxx.

References _D, _end_thresh, _F, _peak_thresh, _start_thresh, fhicl::ParameterSet::get(), and Reset().

27  : PMTPulseRecoBase(name)
28  //*********************************************************************
29  {
30 
31  _F = pset.get<float>("Fraction");
32  _D = pset.get<int> ("Delay");
33 
34  //_number_presample = pset.get<int> ("BaselinePreSample");
35  _peak_thresh = pset.get<double>("PeakThresh");
36  _start_thresh = pset.get<double>("StartThresh");
37  _end_thresh = pset.get<double>("EndThresh");
38 
39  Reset();
40 
41  }
void Reset()
Implementation of AlgoCFD::reset() method.
Definition: AlgoCFD.cxx:49
double _end_thresh
Definition: AlgoCFD.h:63
double _start_thresh
Definition: AlgoCFD.h:62
T get(std::string const &key) const
Definition: ParameterSet.h:231
double _peak_thresh
Definition: AlgoCFD.h:61
PMTPulseRecoBase(const std::string name="noname")
Default constructor with fhicl parameters.
pmtana::AlgoCFD::~AlgoCFD ( )
virtual

Default destructor.

Definition at line 44 of file AlgoCFD.cxx.

46  {}

Member Function Documentation

bool pmtana::PMTPulseRecoBase::Derivative ( const std::vector< short > &  wf,
std::vector< int32_t > &  diff,
size_t  begin = 0,
size_t  end = 0 
) const
protectedinherited

A method to compute derivative, which is a simple subtraction of previous ADC sample from each sample. The result is stored in the input "diff" reference vector which is int32_t type as a derivative could be negative.

Definition at line 122 of file PMTPulseRecoBase.cxx.

References evd::details::begin(), pmtana::CheckIndex(), and evd::details::end().

127  {
128 
129  if(CheckIndex(wf,begin,end)){
130 
131  diff.clear();
132  diff.reserve(end - begin);
133 
134  for(size_t index = begin ; index <= end ; ++index)
135 
136  diff.push_back(wf.at(index+1) - wf.at(index));
137 
138  return true;
139  }
140 
141  return false;
142 
143  }
bool CheckIndex(const std::vector< short > &wf, const size_t &begin, size_t &end)
std::vector< evd::details::RawDigitInfo_t >::const_iterator begin(RawDigitCacheDataClass const &cache)
std::vector< evd::details::RawDigitInfo_t >::const_iterator end(RawDigitCacheDataClass const &cache)
size_t pmtana::PMTPulseRecoBase::GetNPulse ( ) const
inlineinherited

A getter for the number of reconstructed pulses from the input waveform.

Definition at line 106 of file PMTPulseRecoBase.h.

Referenced by opdet::LEDCalibrationAna::analyze().

106 {return _pulse_v.size();};
pulse_param_array _pulse_v
A container array of pulse_param struct objects to store (possibly multiple) reconstructed pulse(s)...
const pulse_param & pmtana::PMTPulseRecoBase::GetPulse ( size_t  index = 0) const
inherited

A getter for the pulse_param struct object. Reconstruction algorithm may have more than one pulse reconstructed from an input waveform. Note you must, accordingly, provide an index key to specify which pulse_param object to be retrieved.

Definition at line 73 of file PMTPulseRecoBase.cxx.

References pmtana::PMTPulseRecoBase::_pulse_v.

Referenced by opdet::LEDCalibrationAna::analyze().

75  {
76 
77  if(index >= _pulse_v.size()) {
78 
79  std::cerr << "\033[93m"
80  << "Invalid pulse index: " << index
81  << "\033[00m"
82  << std::endl;
83 
84  throw std::exception();
85  }
86 
87  else return _pulse_v.at(index);
88 
89  }
pulse_param_array _pulse_v
A container array of pulse_param struct objects to store (possibly multiple) reconstructed pulse(s)...
cet::coded_exception< error, detail::translate > exception
Definition: exception.h:33
const pulse_param_array & pmtana::PMTPulseRecoBase::GetPulses ( ) const
inherited

A getter for the whole array of pulse_param struct object.

Definition at line 92 of file PMTPulseRecoBase.cxx.

References pmtana::PMTPulseRecoBase::_pulse_v.

Referenced by opdet::RunHitFinder().

94  {
95  return _pulse_v;
96  }
pulse_param_array _pulse_v
A container array of pulse_param struct objects to store (possibly multiple) reconstructed pulse(s)...
bool pmtana::PMTPulseRecoBase::Integral ( const std::vector< short > &  wf,
double &  result,
size_t  begin = 0,
size_t  end = 0 
) const
protectedinherited

A method to integrate an waveform from index "begin" to the "end". The result is filled in "result" reference. If the "end" is default (=0), then "end" is set to the last index of the waveform.

Definition at line 99 of file PMTPulseRecoBase.cxx.

References evd::details::begin(), pmtana::CheckIndex(), and evd::details::end().

Referenced by pmtana::AlgoFixedWindow::RecoPulse().

104  {
105 
106  if(!CheckIndex(wf,begin,end)) return false;
107 
108  std::vector<short>::const_iterator begin_iter(wf.begin());
109 
110  std::vector<short>::const_iterator end_iter(wf.begin());
111 
112  begin_iter = begin_iter + begin;
113 
114  end_iter = end_iter + end + 1;
115 
116  result = (double)(std::accumulate(begin_iter, end_iter, 0));
117 
118  return true;
119  }
intermediate_table::const_iterator const_iterator
bool CheckIndex(const std::vector< short > &wf, const size_t &begin, size_t &end)
std::vector< evd::details::RawDigitInfo_t >::const_iterator begin(RawDigitCacheDataClass const &cache)
std::vector< evd::details::RawDigitInfo_t >::const_iterator end(RawDigitCacheDataClass const &cache)
const std::map< unsigned, double > pmtana::AlgoCFD::LinearZeroPointX ( const std::vector< double > &  trace)
protected

Definition at line 240 of file AlgoCFD.cxx.

References pmtana::sign().

Referenced by RecoPulse().

240  {
241 
242  std::map<unsigned,double> crossing;
243 
244  //step through the trace and find where slope is POSITIVE across zero
245  for ( unsigned i = 0; i < trace.size() - 1; ++i) {
246 
247  auto si = ::pmtana::sign(trace.at(i));
248  auto sf = ::pmtana::sign(trace.at(i+1));
249 
250  if ( si == sf ) //no sign flip, no zero cross
251  continue;
252 
253  if ( sf < si ) //this is a negative slope, continue
254  continue;
255 
256  //calculate the crossing X based on linear interpolation bt two pts
257 
258  crossing[i] = (double) i - trace.at(i) * ( 1.0 / ( trace.at(i+1) - trace.at(i) ) );
259 
260  }
261 
262 
263  return crossing;
264 
265  }
int sign(double val)
Definition: UtilFunc.cxx:106
size_t pmtana::PMTPulseRecoBase::Max ( const std::vector< short > &  wf,
double &  result,
size_t  begin = 0,
size_t  end = 0 
) const
protectedinherited

A method to return the maximum value of ADC sample within the index from "begin" to "end". If the "end" is default (=0), then "end" is set to the last index of the waveform.

Definition at line 146 of file PMTPulseRecoBase.cxx.

References evd::details::begin(), pmtana::CheckIndex(), and evd::details::end().

Referenced by pmtana::AlgoFixedWindow::RecoPulse().

151  {
152 
153  size_t target_index = wf.size() + 1;
154 
155  result = 0;
156 
157  if(CheckIndex(wf,begin,end)) {
158 
159  for(size_t index = begin; index <= end; ++index)
160 
161  if( result < wf.at(index)) { target_index = index; result = (double)(wf.at(index)); }
162 
163  }
164 
165  return target_index;
166 
167  }
bool CheckIndex(const std::vector< short > &wf, const size_t &begin, size_t &end)
std::vector< evd::details::RawDigitInfo_t >::const_iterator begin(RawDigitCacheDataClass const &cache)
std::vector< evd::details::RawDigitInfo_t >::const_iterator end(RawDigitCacheDataClass const &cache)
size_t pmtana::PMTPulseRecoBase::Min ( const std::vector< short > &  wf,
double &  result,
size_t  begin = 0,
size_t  end = 0 
) const
protectedinherited

A method to return the minimum value of ADC sample within the index from "begin" to "end". If the "end" is default (=0), then "end" is set to the last index of the waveform.

Definition at line 170 of file PMTPulseRecoBase.cxx.

References evd::details::begin(), pmtana::CheckIndex(), and evd::details::end().

175  {
176 
177  size_t target_index = wf.size() + 1;
178 
179  result = 4096;
180 
181  if(CheckIndex(wf,begin,end)) {
182 
183  for(size_t index = begin; index <= end; ++index)
184 
185  if( result > wf.at(index)) { target_index = index; result = (double)(wf.at(index)); }
186 
187  }
188 
189  return target_index;
190 
191  }
bool CheckIndex(const std::vector< short > &wf, const size_t &begin, size_t &end)
std::vector< evd::details::RawDigitInfo_t >::const_iterator begin(RawDigitCacheDataClass const &cache)
std::vector< evd::details::RawDigitInfo_t >::const_iterator end(RawDigitCacheDataClass const &cache)
const std::string & pmtana::PMTPulseRecoBase::Name ( ) const
inherited

Name getter.

Definition at line 21 of file PMTPulseRecoBase.cxx.

References pmtana::PMTPulseRecoBase::_name.

23  { return _name; }
std::string _name
Unique name.
bool pmtana::PMTPulseRecoBase::Reconstruct ( const pmtana::Waveform_t wf,
const pmtana::PedestalMean_t mean_v,
const pmtana::PedestalSigma_t sigma_v 
)
inherited

A core method: this executes the algorithm and stores reconstructed parameters in the pulse_param struct object.

Definition at line 36 of file PMTPulseRecoBase.cxx.

References pmtana::PMTPulseRecoBase::_status, and pmtana::PMTPulseRecoBase::RecoPulse().

40  {
41  _status = this->RecoPulse(wf,mean_v,sigma_v);
42  return _status;
43  }
bool _status
Status after pulse reconstruction.
virtual bool RecoPulse(const pmtana::Waveform_t &, const pmtana::PedestalMean_t &, const pmtana::PedestalSigma_t &)=0
bool pmtana::AlgoCFD::RecoPulse ( const pmtana::Waveform_t wf,
const pmtana::PedestalMean_t mean_v,
const pmtana::PedestalSigma_t sigma_v 
)
protectedvirtual

Implementation of AlgoCFD::reco() method.

Implements pmtana::PMTPulseRecoBase.

Definition at line 56 of file AlgoCFD.cxx.

References _D, _end_thresh, _F, _peak_thresh, pmtana::PMTPulseRecoBase::_pulse, pmtana::PMTPulseRecoBase::_pulse_v, _start_thresh, pmtana::pulse_param::area, evd::details::begin(), geo::vect::cross(), LinearZeroPointX(), pmtana::pulse_param::peak, pmtana::pulse_param::ped_mean, Reset(), pmtana::pulse_param::reset_param(), pmtana::pulse_param::t_cfdcross, pmtana::pulse_param::t_end, pmtana::pulse_param::t_max, and pmtana::pulse_param::t_start.

60  {
61 
62  Reset();
63 
64  std::vector<double> cfd; cfd.reserve(wf.size());
65 
66  // follow cfd procedure: invert waveform, multiply by constant fraction
67  // add to delayed waveform.
68  for (unsigned int k = 0; k < wf.size(); ++k) {
69 
70  auto delayed = -1.0 * _F * ( (float) wf.at(k) - mean_v.at(k) );
71 
72  if ((int)k < _D)
73 
74  cfd.push_back( delayed );
75 
76  else
77 
78  cfd.push_back(delayed + ( (float) wf.at(k - _D) - mean_v.at(k) ) );
79  }
80 
81 
82  // Get the zero point crossings, how can I tell which are meaningful?
83  // go to each crossing, see if waveform is above pedestal (high above pedestal)
84 
85  auto crossings = LinearZeroPointX(cfd);
86 
87  // lambda criteria to determine if inside pulse
88 
89  auto in_peak = [&wf,&sigma_v,&mean_v](int i, float thresh) -> bool
90  { return wf.at(i) > sigma_v.at(i) * thresh + mean_v.at(i); };
91 
92  // loop over CFD crossings
93  for(const auto& cross : crossings) {
94 
95  if( in_peak( cross.first, _peak_thresh) ) {
97 
98  int i = cross.first;
99 
100  //backwards
101  while ( in_peak(i, _start_thresh) ){
102  i--;
103  if ( i < 0 ) { i = 0; break; }
104  }
105  _pulse.t_start = i;
106 
107  //walk a little further backwards to see if we can get 5 low RMS
108  // while ( !in_peak(i,_start_thresh) ) {
109  // if (i == ( _pulse.t_start - _number_presample ) ) break;
110  // i--;
111  // if ( i < 0 ) { i = 0; break; }
112  // }
113 
114  // auto before_mean = double{0.0};
115 
116  // if ( _pulse.t_start - i > 0 )
117  // before_mean = std::accumulate(std::begin(mean_v) + i,
118  // std::begin(mean_v) + _pulse.t_start, 0.0) / ((double) (_pulse.t_start - i));
119 
120  i = _pulse.t_start + 1;
121 
122  //forwards
123  while ( in_peak(i,_end_thresh) ) {
124  i++;
125  if ( i > (int)(wf.size()) - 1 ) { i = (int)(wf.size()) - 1; break; }
126  }
127 
128  _pulse.t_end = i;
129 
130  // //walk a little further forwards to see if we can get 5 low RMS
131  // while ( !in_peak(i,_end_thresh) ) {
132  // if (i == ( _pulse.t_end + _number_presample ) ) break;
133  // i++;
134  // if ( i > wf.size() - 1 ) { i = wf.size() - 1; break; }
135  // }
136 
137  // auto after_mean = double{0.0};
138 
139  // if( i - _pulse.t_end > 0)
140  // after_mean = std::accumulate(std::begin(mean_v) + _pulse.t_end + 1,
141  // std::begin(mean_v) + i + 1, 0.0) / ((double) (i - _pulse.t_end));
142 
143 
144  //how to decide before or after? set before for now
145  //if ( wf.size() < 1500 ) //it's cosmic discriminator
146  //before_mean = mean_v.front();
147 
148  // if( after_mean <= 0 and before_mean <= 0 ) {
149  // std::cerr << "\033[93m<<" << __FUNCTION__ << ">>\033[00m Could not find good pedestal for CDF"
150  // << " both before_mean and after_mean are zero or less? Ignoring this crossing." << std::endl;
151  // continue;
152  // }
153 
154  //x
155 
156  auto start_ped = mean_v.at(_pulse.t_start);
157  auto end_ped = mean_v.at(_pulse.t_end);
158 
159  //just take the "smaller one"
160  _pulse.ped_mean = start_ped <= end_ped ? start_ped : end_ped;
161 
162  if(wf.size() < 50) _pulse.ped_mean = mean_v.front(); //is COSMIC DISCRIMINATOR
163 
164  auto it = std::max_element(std::begin(wf) + _pulse.t_start, std::begin(wf) + _pulse.t_end);
165 
166  _pulse.t_max = it - std::begin(wf);
167  _pulse.peak = *it - _pulse.ped_mean;
168  _pulse.t_cfdcross = cross.second;
169 
170  for(auto k = _pulse.t_start; k <= _pulse.t_end; ++k) {
171  auto a = wf.at(k) - _pulse.ped_mean;
172  if ( a > 0 ) _pulse.area += a;
173  }
174 
175  _pulse_v.push_back(_pulse);
176  }
177 
178  }
179 
180  // Vic:
181  // Very close in time pulses have multiple CFD
182  // crossing points. Should we check that pulses now have
183  // some multiplicity? No lets just delete them.
184 
185  auto pulses_copy = _pulse_v;
186  _pulse_v.clear();
187 
188  std::unordered_map<unsigned,pulse_param> delta;
189 
190  //unsigned width = 0;
191  for( const auto& p : pulses_copy ) {
192 
193  if ( delta.count(p.t_start) ) {
194  if ( (p.t_end - p.t_start) > (delta[p.t_start].t_end - delta[p.t_start].t_start) )
195  delta[p.t_start] = p;
196  else
197  continue;
198  }
199  else {
200  delta[p.t_start] = p;
201  }
202  }
203 
204  for(const auto & p : delta)
205  _pulse_v.push_back(p.second);
206 
207 
208  //do the same now ensure t_final's are all unique
209  //width = 0;
210 
211  pulses_copy.clear();
212  pulses_copy = _pulse_v;
213 
214  _pulse_v.clear();
215  delta.clear();
216 
217  for( const auto& p : pulses_copy ) {
218 
219  if ( delta.count(p.t_end) ) {
220  if ( (p.t_end - p.t_start) > (delta[p.t_end].t_end - delta[p.t_end].t_start) )
221  delta[p.t_end] = p;
222  else
223  continue;
224  }
225  else {
226  delta[p.t_end] = p;
227  }
228  }
229 
230  for(const auto & p : delta)
231  _pulse_v.push_back(p.second);
232 
233  //there should be no overlapping pulses now...
234 
235  return true;
236 
237  }
const std::map< unsigned, double > LinearZeroPointX(const std::vector< double > &trace)
Definition: AlgoCFD.cxx:240
void Reset()
Implementation of AlgoCFD::reset() method.
Definition: AlgoCFD.cxx:49
double _end_thresh
Definition: AlgoCFD.h:63
pulse_param _pulse
A subject pulse_param object to be filled with the last reconstructed pulse parameters.
double _start_thresh
Definition: AlgoCFD.h:62
std::vector< evd::details::RawDigitInfo_t >::const_iterator begin(RawDigitCacheDataClass const &cache)
double _peak_thresh
Definition: AlgoCFD.h:61
Vector cross(Vector const &a, Vector const &b)
Return cross product of two vectors.
pulse_param_array _pulse_v
A container array of pulse_param struct objects to store (possibly multiple) reconstructed pulse(s)...
void pmtana::AlgoCFD::Reset ( )
virtual

Implementation of AlgoCFD::reset() method.

Reimplemented from pmtana::PMTPulseRecoBase.

Definition at line 49 of file AlgoCFD.cxx.

References pmtana::PMTPulseRecoBase::Reset().

Referenced by AlgoCFD(), and RecoPulse().

51  {
53  }
virtual void Reset()
A method to be called event-wise to reset parameters.
const bool pmtana::PMTPulseRecoBase::Status ( ) const
inherited

Status getter.

Definition at line 26 of file PMTPulseRecoBase.cxx.

References pmtana::PMTPulseRecoBase::_status.

28  { return _status; }
bool _status
Status after pulse reconstruction.

Member Data Documentation

int pmtana::AlgoCFD::_D
private

Definition at line 58 of file AlgoCFD.h.

Referenced by AlgoCFD(), and RecoPulse().

double pmtana::AlgoCFD::_end_thresh
private

Definition at line 63 of file AlgoCFD.h.

Referenced by AlgoCFD(), and RecoPulse().

float pmtana::AlgoCFD::_F
private

Definition at line 57 of file AlgoCFD.h.

Referenced by AlgoCFD(), and RecoPulse().

double pmtana::AlgoCFD::_peak_thresh
private

Definition at line 61 of file AlgoCFD.h.

Referenced by AlgoCFD(), and RecoPulse().

pulse_param pmtana::PMTPulseRecoBase::_pulse
protectedinherited
double pmtana::AlgoCFD::_start_thresh
private

Definition at line 62 of file AlgoCFD.h.

Referenced by AlgoCFD(), and RecoPulse().


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