| 1 | #include "simtoipr.h"
|
|---|
| 2 | #include "toimanager.h"
|
|---|
| 3 | #include "pexceptions.h"
|
|---|
| 4 | #include "ctimer.h"
|
|---|
| 5 |
|
|---|
| 6 | SimpleDeglitcher::SimpleDeglitcher(int wsz, double ns, int maxnpt)
|
|---|
| 7 | {
|
|---|
| 8 | if (wsz < 5) wsz = 5;
|
|---|
| 9 | if (wsz%2 == 0) wsz++;
|
|---|
| 10 | if (maxnpt > wsz) maxnpt = wsz;
|
|---|
| 11 |
|
|---|
| 12 | cout << "SimpleDeglitcher::SimpleDeglitcher() WSize= " << wsz
|
|---|
| 13 | << " nSig=" << ns << " maxNPt=" << maxnpt << endl;
|
|---|
| 14 |
|
|---|
| 15 | wsize = wsz;
|
|---|
| 16 | nsig = ns;
|
|---|
| 17 | maxpoints = maxnpt;
|
|---|
| 18 | totnscount = glnscount = glcount = 0;
|
|---|
| 19 | deglitchdone = false;
|
|---|
| 20 | }
|
|---|
| 21 |
|
|---|
| 22 | SimpleDeglitcher::~SimpleDeglitcher()
|
|---|
| 23 | {
|
|---|
| 24 | }
|
|---|
| 25 |
|
|---|
| 26 | void SimpleDeglitcher::PrintStatus(ostream & os)
|
|---|
| 27 | {
|
|---|
| 28 | os << "\n ------------------------------------------------------ \n"
|
|---|
| 29 | << " SimpleDeglitcher::PrintStatus() - WindowSize=" << WSize()
|
|---|
| 30 | << " NbSigmas=" << NbSigmas() << " MaxPoints=" << MaxPoints() << endl;
|
|---|
| 31 | TOIProcessor::PrintStatus(os);
|
|---|
| 32 | if (deglitchdone) os << " Deglitching performed " << endl;
|
|---|
| 33 | else os << " NO deglitching done " << endl;
|
|---|
| 34 | double nst = (TotalSampleCount() > 0) ? TotalSampleCount() : 1.;
|
|---|
| 35 | os << " TotalSampleCount=" << TotalSampleCount() << " GlitchCount= "
|
|---|
| 36 | << GlitchCount() << " GlitchSampleCount=" << GlitchSampleCount()
|
|---|
| 37 | << "( " << (double)GlitchSampleCount()*100./nst << " % )" << endl;
|
|---|
| 38 | os << " ------------------------------------------------------ " << endl;
|
|---|
| 39 | }
|
|---|
| 40 |
|
|---|
| 41 | void SimpleDeglitcher::init() {
|
|---|
| 42 | cout << "SimpleDeglitcher::init" << endl;
|
|---|
| 43 | declareInput("in");
|
|---|
| 44 | declareOutput("out");
|
|---|
| 45 | declareOutput("mean");
|
|---|
| 46 | declareOutput("sigma");
|
|---|
| 47 | declareOutput("incopie");
|
|---|
| 48 | name = "SimpleDeglitcher";
|
|---|
| 49 | // upExtra = 1; A quoi ca sert ?
|
|---|
| 50 | }
|
|---|
| 51 |
|
|---|
| 52 | void SimpleDeglitcher::run() {
|
|---|
| 53 |
|
|---|
| 54 | // TOIManager* mgr = TOIManager::getManager();
|
|---|
| 55 | int snb = getMinIn();
|
|---|
| 56 | int sne = getMaxIn();
|
|---|
| 57 |
|
|---|
| 58 | bool fgout = checkOutputTOIIndex(0);
|
|---|
| 59 | bool fgmean = checkOutputTOIIndex(1);
|
|---|
| 60 | bool fgsigma = checkOutputTOIIndex(2);
|
|---|
| 61 | bool fgincopie = checkOutputTOIIndex(3);
|
|---|
| 62 |
|
|---|
| 63 | if (!checkInputTOIIndex(0)) {
|
|---|
| 64 | cerr << " SimpleDeglitcher::run() - Input TOI (in) not connected! "
|
|---|
| 65 | << endl;
|
|---|
| 66 | throw ParmError("SimpleDeglitcher::run() Input TOI (in) not connected!");
|
|---|
| 67 | }
|
|---|
| 68 | if (!fgout && !fgmean && !fgsigma &&!fgincopie) {
|
|---|
| 69 | cerr << " SimpleDeglitcher::run() - No Output TOI connected! "
|
|---|
| 70 | << endl;
|
|---|
| 71 | throw ParmError("SimpleDeglitcher::run() No output TOI connected!");
|
|---|
| 72 | }
|
|---|
| 73 |
|
|---|
| 74 | if (!fgout) {
|
|---|
| 75 | cout << "Warning: SimpleDeglitcher::run() - No TOI connected to out \n"
|
|---|
| 76 | << " No deglitching would be performed !" << endl;
|
|---|
| 77 | }
|
|---|
| 78 |
|
|---|
| 79 | cout << " SimpleDeglitcher::run() SNRange=" << snb << " - " << sne << endl;
|
|---|
| 80 | try {
|
|---|
| 81 | Timer tm("SimpleDeglitcher::run()");
|
|---|
| 82 | // Le debut
|
|---|
| 83 | int wsz2 = wsize/2;
|
|---|
| 84 | Vector vin(wsize);
|
|---|
| 85 | int k;
|
|---|
| 86 | for(k=0; k<wsize; k++)
|
|---|
| 87 | vin(k) = getData(0, k+snb);
|
|---|
| 88 | double s = vin.Sum();
|
|---|
| 89 | double mean = s/wsize;
|
|---|
| 90 | double s2 = vin.SumX2();
|
|---|
| 91 | double sigma = sqrt(s2/wsize-mean*mean);
|
|---|
| 92 | int kgl = -1;
|
|---|
| 93 | int ii,lastput;
|
|---|
| 94 | bool fgglitch = false;
|
|---|
| 95 | double valcur,valsub,valadd;
|
|---|
| 96 | for(k=0;k<=sne-snb;k++) {
|
|---|
| 97 | totnscount++;
|
|---|
| 98 | // if (k%10000 == 0) cout << " DBG: K=" << k << endl;
|
|---|
| 99 | // Calcul mean-sigma
|
|---|
| 100 | if ((k>wsz2) && (k<sne-snb-wsz2)) {
|
|---|
| 101 | valsub = vin((k+wsz2)%wsize);
|
|---|
| 102 | vin((k+wsz2)%wsize) = valadd = getData(0, k+wsz2+snb);
|
|---|
| 103 | s += (valadd-valsub);
|
|---|
| 104 | s2 += (valadd*valadd-valsub*valsub);
|
|---|
| 105 | mean = s/wsize;
|
|---|
| 106 | sigma = sqrt(s2/wsize-mean*mean);
|
|---|
| 107 | }
|
|---|
| 108 |
|
|---|
| 109 | // On gere les sorties Mean et Sigma
|
|---|
| 110 | if (fgmean)
|
|---|
| 111 | putData(1, k+snb, mean, 0);
|
|---|
| 112 | if (fgsigma)
|
|---|
| 113 | putData(2, k+snb, sigma, 0);
|
|---|
| 114 |
|
|---|
| 115 | valcur = vin(k%wsize);
|
|---|
| 116 | if (fgincopie)
|
|---|
| 117 | putData(3, k+snb, valcur, 0);
|
|---|
| 118 |
|
|---|
| 119 | if (!fgout) continue; // Pas de sortie out (deglitche)
|
|---|
| 120 |
|
|---|
| 121 |
|
|---|
| 122 | // if (k<100) {
|
|---|
| 123 | // cout << "DBG-A-Deglitch[" << k << "] mean="
|
|---|
| 124 | // << mean << " sigma=" << sigma << " valcur="
|
|---|
| 125 | // << valcur << " Kgl=" << kgl ;
|
|---|
| 126 | // if (fgglitch) cout << " In Glitch" ;
|
|---|
| 127 | // cout << endl;
|
|---|
| 128 | // }
|
|---|
| 129 |
|
|---|
| 130 | if (valcur < mean+nsig*sigma) { // inferieur au seuil
|
|---|
| 131 | if (fgglitch) {
|
|---|
| 132 | if (k-kgl < maxpoints) { // On vient de detecter un glitch
|
|---|
| 133 | glcount++;
|
|---|
| 134 | for(ii=kgl; ii<k; ii++) {
|
|---|
| 135 | putData(0, ii+snb, mean, 7);
|
|---|
| 136 | glnscount++;
|
|---|
| 137 | }
|
|---|
| 138 | lastput = snb+k-1;
|
|---|
| 139 | }
|
|---|
| 140 | else {
|
|---|
| 141 | for(ii=kgl; ii<k; ii++) {
|
|---|
| 142 | putData(0, ii+snb, vin(ii%wsize), 0);
|
|---|
| 143 | }
|
|---|
| 144 | lastput = snb+k-1;
|
|---|
| 145 | }
|
|---|
| 146 | }
|
|---|
| 147 | putData(0, k+snb, valcur, 0);
|
|---|
| 148 | lastput = snb+k;
|
|---|
| 149 | kgl = -1; fgglitch = false;
|
|---|
| 150 | }
|
|---|
| 151 | else { // Superieur au seuil
|
|---|
| 152 | if (fgglitch) {
|
|---|
| 153 | if (k-kgl+1 >= maxpoints) { // serie de points > seuil
|
|---|
| 154 | for(ii=kgl; ii<=k; ii++) // -> Donc pas glitch
|
|---|
| 155 | putData(0, ii+snb, vin(ii%wsize), 0);
|
|---|
| 156 | lastput = snb+k;
|
|---|
| 157 | fgglitch = false;
|
|---|
| 158 | }
|
|---|
| 159 | }
|
|---|
| 160 | else {
|
|---|
| 161 | if (kgl < 0) { // debut possible de glitch
|
|---|
| 162 | fgglitch = true; kgl = k;
|
|---|
| 163 | }
|
|---|
| 164 | else { // On est toujours dans une serie > seuil
|
|---|
| 165 | putData(0, k+snb, valcur, 0);
|
|---|
| 166 | lastput = snb+k;
|
|---|
| 167 | }
|
|---|
| 168 | }
|
|---|
| 169 | }
|
|---|
| 170 |
|
|---|
| 171 | // if (k%5000 == 0) cout << " ---DBG2: K=" << k << " glcount="
|
|---|
| 172 | // << glcount << " LastPut= " << lastput << endl;
|
|---|
| 173 | } // Fin de Boucle sur les num-sample
|
|---|
| 174 |
|
|---|
| 175 | //DBG cout << " La fin lastput=" << lastput << " SNE=" << sne;
|
|---|
| 176 | //DBG for(k=lastput-snb+1; k<sne-snb; k++)
|
|---|
| 177 | //DBG putData(0, k+snb, vin(k%wsize), 0);
|
|---|
| 178 | // cout << " DBG3- OUT of try bloc ! " << endl;
|
|---|
| 179 | cout << " SimpleDeglitcher::run() - End of processing "
|
|---|
| 180 | << " TotalSampleCount=" << TotalSampleCount()
|
|---|
| 181 | << " GlitchCount= " << GlitchCount() << endl;
|
|---|
| 182 | if (fgout) deglitchdone = true;
|
|---|
| 183 | } // Bloc try
|
|---|
| 184 | catch (PException & exc) {
|
|---|
| 185 | cerr << "SimpleDeglitcher: Catched Exception " << (string)typeid(exc).name()
|
|---|
| 186 | << "\n .... Msg= " << exc.Msg() << endl;
|
|---|
| 187 | }
|
|---|
| 188 | }
|
|---|
| 189 |
|
|---|
| 190 |
|
|---|
| 191 |
|
|---|
| 192 | SimpleFilter::SimpleFilter(int wsz, FilterKind fk, double a, double s)
|
|---|
| 193 | {
|
|---|
| 194 | if (wsz < 3) wsz = 3;
|
|---|
| 195 | if (wsz%2 == 0) wsz++;
|
|---|
| 196 | cout << "SimpleFilter::SimpleFilter() wsz= " << wsz
|
|---|
| 197 | << " fk=" << fk << endl;
|
|---|
| 198 | wsize = wsz;
|
|---|
| 199 | totnscount = 0;
|
|---|
| 200 | coef = new double[wsz];
|
|---|
| 201 | for(int kk=0; kk<wsz; kk++)
|
|---|
| 202 | coef[kk] = 1./wsize;
|
|---|
| 203 | }
|
|---|
| 204 |
|
|---|
| 205 | SimpleFilter::~SimpleFilter()
|
|---|
| 206 | {
|
|---|
| 207 | delete[] coef;
|
|---|
| 208 | }
|
|---|
| 209 |
|
|---|
| 210 | void SimpleFilter::PrintStatus(ostream & os)
|
|---|
| 211 | {
|
|---|
| 212 | os << "\n ------------------------------------------------------ \n"
|
|---|
| 213 | << " SimpleFilter::PrintStatus() - WindowSize=" << WSize()
|
|---|
| 214 | << " FilterKind= " << Type() << endl;
|
|---|
| 215 | TOIProcessor::PrintStatus(os);
|
|---|
| 216 | os << " Coeff= " ;
|
|---|
| 217 | for(int k=0; k<wsize; k++) os << coef[k] << " " ;
|
|---|
| 218 | os << endl;
|
|---|
| 219 | os << " TotalSampleCount=" << TotalSampleCount() << endl;
|
|---|
| 220 | os << " ------------------------------------------------------ " << endl;
|
|---|
| 221 | }
|
|---|
| 222 |
|
|---|
| 223 | void SimpleFilter::init() {
|
|---|
| 224 | cout << "SimpleFilter::init" << endl;
|
|---|
| 225 | declareInput("in");
|
|---|
| 226 | declareOutput("out");
|
|---|
| 227 | declareOutput("incopie");
|
|---|
| 228 | name = "SimpleFilter";
|
|---|
| 229 | // upExtra = 1;
|
|---|
| 230 | }
|
|---|
| 231 |
|
|---|
| 232 | void SimpleFilter::run() {
|
|---|
| 233 | // TOIManager* mgr = TOIManager::getManager();
|
|---|
| 234 | int snb = getMinIn();
|
|---|
| 235 | int sne = getMaxIn();
|
|---|
| 236 |
|
|---|
| 237 | bool fgout = checkOutputTOIIndex(0);
|
|---|
| 238 | bool fgincopie = checkOutputTOIIndex(1);
|
|---|
| 239 |
|
|---|
| 240 | if (!checkInputTOIIndex(0)) {
|
|---|
| 241 | cerr << " SimpleFilter::run() - Input TOI (in) not connected! "
|
|---|
| 242 | << endl;
|
|---|
| 243 | throw ParmError("SimpleFilter::run() Input TOI (in) not connected!");
|
|---|
| 244 | }
|
|---|
| 245 | if (!fgout) {
|
|---|
| 246 | cerr << " SimpleFilter::run() - No Output TOI connected! "
|
|---|
| 247 | << endl;
|
|---|
| 248 | throw ParmError("SimpleFilter::run() No output TOI connected!");
|
|---|
| 249 | }
|
|---|
| 250 |
|
|---|
| 251 | cout << " SimpleFilter::run() SNRange=" << snb << " - " << sne << endl;
|
|---|
| 252 |
|
|---|
| 253 |
|
|---|
| 254 | try {
|
|---|
| 255 | Timer tm("SimpleFilter::run()");
|
|---|
| 256 | // Le debut
|
|---|
| 257 | int wsz2 = wsize/2;
|
|---|
| 258 | Vector vin(wsize);
|
|---|
| 259 | int k;
|
|---|
| 260 | for(k=0; k<wsize; k++)
|
|---|
| 261 | vin(k) = getData(0, k+snb);
|
|---|
| 262 | double mean = vin.Sum()/wsize;
|
|---|
| 263 | for(k=0; k<wsize/2; k++) vin(k) = mean;
|
|---|
| 264 | for(k=0;k<=sne-snb;k++) {
|
|---|
| 265 | totnscount++;
|
|---|
| 266 | // Calcul mean-sigma
|
|---|
| 267 | if ((k>wsz2) && (k<sne-snb-wsz2))
|
|---|
| 268 | vin((k+wsz2)%wsize) = getData(0, k+snb);
|
|---|
| 269 | // if (k == sne-snb-wsz2) {
|
|---|
| 270 | // }
|
|---|
| 271 | double sortie = 0;
|
|---|
| 272 | for(int ii=0; ii<wsize; ii++) {
|
|---|
| 273 | sortie += vin((ii+k)%wsize)*coef[ii];
|
|---|
| 274 | }
|
|---|
| 275 | putData(0,k+snb,sortie,0);
|
|---|
| 276 | if (fgincopie) putData(1,k+snb,vin(k%wsize),0);
|
|---|
| 277 | } // Boucle sur les num-sample
|
|---|
| 278 | cout << " SimpleFilter::run() - End of processing " << endl;
|
|---|
| 279 | } // Bloc try
|
|---|
| 280 |
|
|---|
| 281 | catch (PException & exc) {
|
|---|
| 282 | cerr << "SimpleFilter: Catched Exception " << (string)typeid(exc).name()
|
|---|
| 283 | << "\n .... Msg= " << exc.Msg() << endl;
|
|---|
| 284 | }
|
|---|
| 285 |
|
|---|
| 286 |
|
|---|
| 287 | }
|
|---|
![(please configure the [header_logo] section in trac.ini)](/PI/chrome/site/your_project_logo.png){kind=logo}
