MedianFilter.h

Go to the documentation of this file.

#ifndef MEDIAN_FILTER_H
#define MEDIAN_FILTER_H

#include "ImageToImageOperator.h"
#include "PumaException.h"

namespace puma2 {

template <class T>
class MedianFilter : public ImageToImageOperator<T,T> {
  public:

    MedianFilter(int windowWidth, int windowHeight);

    MedianFilter(int windowSize);

    virtual ~MedianFilter();

          virtual void apply(const T& iImg, T& oImg);

    int getWindowWidth() const { return mWindowWidth; }

    int getWindowHeight() const { return mWindowHeight; }

  private:

    void computeHelpers();

    int mWindowWidth;

    int mWindowHeight;

    int mHalfWindowWidth;

    int mHalfWindowHeight;

    int mNumWindowElements;

    typename T::ElementType* mWindowElements;

    int mMedianIndex;
};


template <class T>
MedianFilter<T>::MedianFilter(int windowSize) {
  if (windowSize < 1 || windowSize % 2 == 0) {
    std::string msg("Illegal window size for median filter (has to be odd and positive)");
    throw PumaException(PumaException::intolerable, msg);
  }

  if (T::numberOfChannels() > 1) {
    //TODO throw warning
  }

  mWindowHeight = windowSize;
  mWindowWidth = windowSize;
  computeHelpers();
}


template <class T>
MedianFilter<T>::MedianFilter(int windowWidth, int windowHeight) {
  if (windowWidth < 1 || windowWidth % 2 == 0) {
    std::string msg("Illegal window width for median filter (has to be odd and positive)");
    throw PumaException(PumaException::intolerable, msg);
  }
  if (windowHeight < 1 || windowHeight % 2 == 0) {
    std::string msg("Illegal window height for median filter (has to be odd and positive)");
    throw PumaException(PumaException::intolerable, msg);
  }

  if (T::numberOfChannels() > 1) {
    //TODO throw warning
  }

  mWindowHeight = windowHeight;
  mWindowWidth = windowWidth;
  computeHelpers();
}

template <class T>
void MedianFilter<T>::computeHelpers() {
  mHalfWindowHeight = (mWindowHeight - 1) / 2;
  mHalfWindowWidth = (mWindowWidth - 1) / 2;
  mNumWindowElements = mWindowWidth * mWindowHeight;
  mWindowElements = new typename T::ElementType[mNumWindowElements];
  mMedianIndex = (mNumWindowElements - 1) / 2;
}

template <class T>
MedianFilter<T>::~MedianFilter() { 
  delete mWindowElements;
};


// NB: 
// - This is the first working implementation
// - much fuster algorithmes for computing the median exist
//   wich are even (almost) independent of the window size.
//   These algorithms will be implemented later

template <class T>
void MedianFilter<T>::apply(const T& iImg, T& oImg) {
  
  int minX = mHalfWindowWidth;
  int maxX = oImg.getWidth() - mHalfWindowWidth - 1;
  int minY = mHalfWindowHeight;
  int maxY = oImg.getHeight() - mHalfWindowHeight - 1;
  int numChannelsM1 = T::numberOfChannels() - 1;

  for (int y = minY; y <= maxY; y++) {
    for (int x = minX; x <= maxX; x++) {
      for (int n = 0; n <= numChannelsM1; n++) {
        int listIndex = 0;
        for (int wy = -mHalfWindowHeight; wy <= mHalfWindowHeight; wy++) {
          for (int wx = -mHalfWindowWidth; wx <= mHalfWindowWidth; wx++) {
            mWindowElements[listIndex++] = iImg.sample(x + wx, y + wy, n);
          }
        }
        // bubblesort mWindowElements
        for (int k = 1; k <= mNumWindowElements - 1; k++) {
          for (int l = 0; l < mNumWindowElements - k; l++) {
            if (mWindowElements[l] > mWindowElements[l + 1]) {
              typename T::ElementType temp = mWindowElements[l];
              mWindowElements[l] = mWindowElements[l + 1];
              mWindowElements[l + 1] = temp;
            }
          }
        }
        oImg.sample(x, y, n) = mWindowElements[mMedianIndex];
      }
    }
  }
}

}

#endif