20 srand(static_cast<unsigned int>(
time(0)));
63 for(
size_t i=0; i< wf.size(); ++i) {
68 if( i < _sample_size || i >= (wf.size() -
_sample_size) )
continue;
82 for(
size_t i=(wf.size() -
_sample_size); i<wf.size(); ++i) {
84 mean_v[i] = mean_v [wf.size() - _sample_size -1];
85 sigma_v[i] = sigma_v[wf.size() - _sample_size -1];
89 float best_sigma = 1.1e9;
90 size_t best_sigma_index = 0;
91 size_t num_good_adc = 0;
93 for(
size_t i=0; i<sigma_v.size(); ++i) {
95 auto const&
mean = mean_v[i];
99 auto const& sigma = sigma_v[i];
100 if(sigma < best_sigma) {
102 best_sigma_index = i;
109 if( num_good_adc < 1 ) {
110 std::cerr <<
"\033[93m<<" << __FUNCTION__ <<
">>\033[00m Could not find good pedestal at all..." <<
std::endl;
115 if(best_sigma >
_max_sigma || num_good_adc < 3) {
116 for(
size_t i=0; i<mean_v.size(); ++i) {
117 mean_v[i] = mean_v.at ( best_sigma_index );
118 sigma_v[i] = sigma_v.at ( best_sigma_index );
127 unsigned nbins = 1000;
142 int last_good_index = -1;
144 for(
size_t i=0; i < wf.size(); ++i) {
146 auto const mean = mean_v[i];
147 auto const sigma = sigma_v[i];
154 if(last_good_index<0) {
155 last_good_index = (
int)i;
159 if( ( last_good_index + 1 ) < (
int) i ) {
161 auto diff = fabs(mean_v.at(last_good_index) -
mean);
163 if ( diff > diff_cutoff) {
165 float slope = (
mean - mean_v.at(last_good_index)) / (
float(i - last_good_index));
167 for(
size_t j = last_good_index + 1; j < i; ++j) {
168 mean_v.at(j) = slope * (
float(j - last_good_index) ) + mean_v.at(last_good_index);
169 sigma_v.at(j) = mode_sigma;
178 auto diff_pre = fabs(presample_mean - mode_mean);
179 auto diff_post = fabs(postsample_mean - mode_mean);
181 auto constant_mean = diff_pre <= diff_post ? presample_mean : postsample_mean;
183 for(
size_t j = last_good_index + 1; j < i; ++j) {
185 mean_v.at(j) = constant_mean;
186 sigma_v.at(j) = mode_sigma;
190 last_good_index = (
int)i;
200 if(sigma_v.front() > mode_sigma) {
202 int first_index = -1;
203 int second_index = -1;
205 for(
size_t i=0; i < wf.size(); ++i) {
206 if( sigma_v.at(i) < mode_sigma ) {
207 if( first_index < 0 ) first_index = (
int)i;
208 else if( second_index < 0 ) {
209 second_index = (
int)i;
215 if(first_index < 0 || second_index < 0) {
216 std::cerr <<
"\033[93m<<" << __FUNCTION__ <<
">>\033[00m Could not find good pedestal for CDF" 218 <<
"first_index: " << first_index <<
"\n" 219 <<
"second_index: " << second_index <<
"\n" 220 <<
"If one of these is less than zero, this means there is pulse\n" 221 <<
"on first sample and baseline never went back down. Returning false here.";
226 auto diff = fabs(mean_v.at(first_index) - mean_v.at(second_index));
228 if ( diff > diff_cutoff) {
230 float slope = (mean_v.at(second_index) - mean_v.at(first_index)) / (
float(second_index - first_index));
232 for(
int j=0; j < first_index; ++j) {
233 mean_v.at(j) = mean_v.at(first_index) - slope * (first_index - j);
241 for(
int j=0; j < first_index; ++j) {
243 mean_v.at(j) = postsample_mean;
244 sigma_v.at(j) = mode_sigma;
251 if(sigma_v.back() > mode_sigma) {
253 int first_index = -1;
254 int second_index = -1;
256 for(
int i = wf.size()-1; i >= 0; --i) {
257 if(sigma_v.at(i) < mode_sigma) {
258 if( second_index < 0 ) second_index = (
int)i;
259 else if( first_index < 0 ) {
260 first_index = (
int)i;
266 if(first_index < 0 || second_index < 0) {
267 std::cerr <<
"\033[93m<<" << __FUNCTION__ <<
">>\033[00m Could not find good pedestal for CDF" 269 <<
"first_index: " << first_index <<
"\n" 270 <<
"second_index: " << second_index <<
"\n" 271 <<
"If one of these is less than zero, this means there is pulse\n" 272 <<
"on the last sample and baseline never went back down. Returning false here.";
277 auto diff = fabs(mean_v.at(second_index) - mean_v.at(first_index) );
279 if ( diff > diff_cutoff) {
281 float slope = (mean_v.at(second_index) - mean_v.at(first_index)) / (
float(second_index - first_index));
282 for(
int j = second_index+1; j <
int(wf.size()); ++j) {
283 mean_v.at(j) = mean_v.at(second_index) + slope * (j-second_index);
292 for(
int j = second_index+1; j <
int(wf.size()); ++j) {
294 mean_v.at(j) = presample_mean;
295 sigma_v.at(j) = mode_sigma;
double std(const std::vector< short > &wf, const double ped_mean, size_t start, size_t nsample)
double BinnedMaxOccurrence(const PedestalMean_t &mean_v, const size_t nbins)
std::vector< double > PedestalSigma_t
bool ComputePedestal(const pmtana::Waveform_t &wf, pmtana::PedestalMean_t &mean_v, pmtana::PedestalSigma_t &sigma_v)
Method to compute a pedestal of the input waveform using "nsample" ADC samples from "start" index...
double edge_aware_mean(const std::vector< short > &wf, int start, int end)
Class definition file of PedAlgoRollingMean.
T get(std::string const &key) const
std::vector< short > Waveform_t
double mean(const std::vector< short > &wf, size_t start, size_t nsample)
PedAlgoRollingMean(const std::string name="PedRollingMean")
Default constructor.
std::vector< double > PedestalMean_t
QTextStream & endl(QTextStream &s)