| 1 | #include <stdio.h>
|
|---|
| 2 | #include <iostream>
|
|---|
| 3 |
|
|---|
| 4 |
|
|---|
| 5 | // this function reads a (hopefully smoothed, and baseline corrected) DRS pipeline and produces
|
|---|
| 6 | // three vectors
|
|---|
| 7 | // 1. start of 'over threshold'
|
|---|
| 8 | // 2. end of 'over threshold'
|
|---|
| 9 | // 3. max within this 'over threshold' region
|
|---|
| 10 | //
|
|---|
| 11 | // it needs a threshold and the length of the falling edge
|
|---|
| 12 |
|
|---|
| 13 | #include "discriminator.h"
|
|---|
| 14 |
|
|---|
| 15 | vector<DiscOut> * discriminator(
|
|---|
| 16 | vector<float>& input, // vector of floats, the discriminator acts on
|
|---|
| 17 | float thr = 5., // threshold
|
|---|
| 18 | int fallingEdge = 100 ) // number of slices, after the maximum, which should be discarded...
|
|---|
| 19 | {
|
|---|
| 20 |
|
|---|
| 21 | vector<DiscOut> * result = new vector<DiscOut>;
|
|---|
| 22 | DiscOut disc;
|
|---|
| 23 | disc.begin = 0;
|
|---|
| 24 | disc.end = 0;
|
|---|
| 25 | disc.maxPos = 0;
|
|---|
| 26 | disc.maxVal = 0.;
|
|---|
| 27 | bool start_under_thr = false;
|
|---|
| 28 | bool over_thr_found = false;
|
|---|
| 29 |
|
|---|
| 30 |
|
|---|
| 31 | bool debug = false; // switch to true in order to generate some output.
|
|---|
| 32 |
|
|---|
| 33 | for ( int sl = 0; sl < input.size() ; sl++ ){
|
|---|
| 34 |
|
|---|
| 35 |
|
|---|
| 36 | if ( input[sl] > thr ) {
|
|---|
| 37 |
|
|---|
| 38 | if ( !over_thr_found ) {
|
|---|
| 39 | over_thr_found = true;
|
|---|
| 40 | disc.begin = sl;
|
|---|
| 41 | }
|
|---|
| 42 |
|
|---|
| 43 | if ( input[sl] > disc.maxVal) {
|
|---|
| 44 | disc.maxVal = input[sl];
|
|---|
| 45 | disc.maxPos = sl;
|
|---|
| 46 | }
|
|---|
| 47 | }
|
|---|
| 48 | if ( input[sl] < thr ) {
|
|---|
| 49 | if ( over_thr_found ) {
|
|---|
| 50 | over_thr_found = false;
|
|---|
| 51 | disc.end = sl;
|
|---|
| 52 |
|
|---|
| 53 | result->push_back(disc);
|
|---|
| 54 | disc.begin = 0;
|
|---|
| 55 | disc.end = 0;
|
|---|
| 56 | disc.maxPos = 0;
|
|---|
| 57 | disc.maxVal = 0.;
|
|---|
| 58 |
|
|---|
| 59 | }
|
|---|
| 60 | }
|
|---|
| 61 | } // end of for llop over all slices
|
|---|
| 62 |
|
|---|
| 63 | if (debug){ // output the vector
|
|---|
| 64 | for (int p=0; p<result->size(); p++ ){
|
|---|
| 65 | cout << p << ":\t";
|
|---|
| 66 | cout << result->at(p).begin << "\t";
|
|---|
| 67 | cout << result->at(p).end << "\t";
|
|---|
| 68 | cout << result->at(p).maxPos << "\t";
|
|---|
| 69 | cout << result->at(p).maxVal << endl;
|
|---|
| 70 | }
|
|---|
| 71 | }
|
|---|
| 72 |
|
|---|
| 73 | DiscOut last;
|
|---|
| 74 | vector<DiscOut> * realresult = new vector<DiscOut>;
|
|---|
| 75 |
|
|---|
| 76 | if (result->size() > 0) {
|
|---|
| 77 | while (!result->empty()){
|
|---|
| 78 |
|
|---|
| 79 | last.begin = result->back().begin;
|
|---|
| 80 | last.end = result->back().end;
|
|---|
| 81 | last.maxPos = result->back().maxPos;
|
|---|
| 82 | last.maxVal = result->back().maxVal;
|
|---|
| 83 | result->pop_back();
|
|---|
| 84 |
|
|---|
| 85 | if (result->empty()){
|
|---|
| 86 | if (last.maxPos > 0)
|
|---|
| 87 | realresult->push_back(last);
|
|---|
| 88 | break;
|
|---|
| 89 | }
|
|---|
| 90 | if (last.maxPos - result->back().maxPos > fallingEdge)
|
|---|
| 91 | realresult->push_back(last);
|
|---|
| 92 | }
|
|---|
| 93 | }
|
|---|
| 94 | /*
|
|---|
| 95 | vector<vector<DiscOut>::iterator> toBeDeleted;
|
|---|
| 96 |
|
|---|
| 97 | for (vector<DiscOut>::iterator it=result->begin(); it != result->end()-1; ++it ){
|
|---|
| 98 | if ( (*it).begin==0 ){
|
|---|
| 99 | toBeDeleted.push_back(it);
|
|---|
| 100 | }
|
|---|
| 101 | if ( ((*(it+1)).maxPos - (*it).maxPos) < fallingEdge ){
|
|---|
| 102 | toBeDeleted.push_back(it);
|
|---|
| 103 | }
|
|---|
| 104 | }
|
|---|
| 105 |
|
|---|
| 106 | for (int i=0; i<toBeDeleted.size(); i++) {
|
|---|
| 107 | result->erase(toBeDeleted[i]);
|
|---|
| 108 | if (debug)
|
|---|
| 109 | cout << "erased peak at:" << (*toBeDeleted[i]).maxPos << endl;
|
|---|
| 110 | }
|
|---|
| 111 | */
|
|---|
| 112 |
|
|---|
| 113 | delete result;
|
|---|
| 114 | return realresult;
|
|---|
| 115 |
|
|---|
| 116 | } // end of function - discriminator
|
|---|
