---

cmEstimator.h

Go to the documentation of this file.

#ifndef _cmEstimator_H
#define _cmEstimator_H

#include <numeric>
#include <algorithm>
#include <functional>
#include "General.h"
#include <iostream>
#include <vector>

template <typename TYPE>
class cmEstimator{

public:

  typedef typename std::vector<TYPE> contType;
  typedef typename contType::iterator iter_type;

  cmEstimator(TYPE outlier = TYPE(3),
              unsigned int iter = 1u):
  m_outlier(outlier),
  m_iter(iter){}

  virtual ~cmEstimator() {}

  TYPE outlierCut() const {return m_outlier;}

  contType calculate(contType& input);  

  int reject() const {return m_reject;}
  
private:

  void initMask(unsigned int contSize);
  
  TYPE mean(contType& input) const;

  TYPE variance(contType& input) const; 

  void subtractMean(contType& input, TYPE mean) const;

  void cmSubtraction(contType& input, TYPE slope) const;

  TYPE slope(contType& input) const;

  void setToZero(contType& input, TYPE s);

  TYPE m_outlier;
  std::vector<unsigned int> m_mask; 
  unsigned int m_iter;
  int m_reject;

};

template <typename TYPE> 
inline typename std::vector<TYPE> cmEstimator<TYPE>::calculate(contType& inputData){

  // initialization ....
  m_reject = 0;
  initMask(inputData.size());

  // mean 
  TYPE mean1 = mean(inputData);

  // subtract it
  contType firstPass;
  firstPass.insert(firstPass.begin(), inputData.begin(),inputData.end());
  subtractMean(firstPass,mean1);
  
  // calculate the slope..
  TYPE s1 = slope(firstPass);

  // common mode
  cmSubtraction(firstPass,s1);

  for (unsigned int it = 0u ; it < m_iter; ++it){

    contType tPass;
    tPass.insert(tPass.begin(), firstPass.begin(), firstPass.end());

    // set outliers to zero using variance
    setToZero(tPass, variance(tPass));

    // mean and subtract it and calc the slope
    TYPE mean2 = mean(tPass);
    subtractMean(tPass,mean2);
    TYPE s2 = slope(tPass);

    subtractMean(firstPass,mean2);  
    cmSubtraction(firstPass,s2);
  } // it

  return firstPass;
} 

template <typename TYPE> 
void cmEstimator<TYPE>::initMask(unsigned int contSize) {

   m_mask.clear();
  for (unsigned int i = 0; i< contSize; ++i){
     m_mask.push_back(0u);
  }
}

template <typename TYPE> 
inline TYPE cmEstimator<TYPE>::mean(contType& input) const{

  TYPE sum = TYPE(0);
  for (unsigned int i = 0; i< input.size(); ++i){
    if (m_mask[i] == 0) sum += input[i];
  } // iter
  return sum/TYPE(input.size()-m_reject);
}

template <typename TYPE> 
inline TYPE cmEstimator<TYPE>::variance(contType& input) const{
  TYPE sum2 = TYPE(0);
  for (unsigned int i = 0; i< input.size(); ++i){
    if (m_mask[i] == 0) sum2 += pow(input[i],2.);
  } // iter
  return sum2/TYPE(input.size()-m_reject);
}

template <typename TYPE> 
void cmEstimator<TYPE>::subtractMean(contType& input, TYPE mean) const{
  // subtract mean
  for (unsigned int i = 0; i< input.size(); ++i){
    if (m_mask[i] == 0){ 
       input[i] = input[i] - mean;  
    }
  } // iter
}

template <typename TYPE> 
inline TYPE cmEstimator<TYPE>::slope(contType& input) const{

  // estimate the slope
  TYPE s = TYPE(0); // zeros
  for (unsigned int i =0; i< input.size(); ++i){
    if (m_mask[i] == 0) s += (input[i] *(int(i)-16));
  } // iter

  //  return 3*s/8192;
  return s/(TYPE(2736));

}

template <typename TYPE> 
inline void cmEstimator<TYPE>::cmSubtraction(contType& input, TYPE s) const{

  // subtract common mode
  for ( unsigned int i = 0; i< input.size(); ++i){
    input[i] = input[i] - (s*(int(i)-16));
  } // i
}

template <typename TYPE> 
inline void cmEstimator<TYPE>::setToZero(contType& input, TYPE var){
 
  // mask off outliers.......
  for (  unsigned int i = 0; i< input.size(); ++i){
   if ((m_mask[i]  == 0)&&(pow(input[i],2.0)> var/m_outlier)){
     //input[i] = TYPE(0);
     m_mask[i] = 1u;
     if (i> 0) m_mask[i-1] = 1u;
     if (i< input.size()-1) m_mask[i+1] = 1u;
     ++m_reject;
   }
  } // iter   
}

#endif // _cmEstimator_H

---