transforms.h

Go to the documentation of this file.

#ifndef TRANSFORMS_H
#define TRANSFORMS_H

#include <itpp/base/vec.h>
#include <itpp/base/mat.h>
#include <itpp/base/matfunc.h>
#include <itpp/itexports.h>


namespace itpp
{




ITPP_EXPORT bool have_fourier_transforms();
ITPP_EXPORT void fft(const cvec &in, cvec &out);
ITPP_EXPORT cvec fft(const cvec &in);
ITPP_EXPORT cvec fft(const cvec &in, const int N);
ITPP_EXPORT void ifft(const cvec &in, cvec &out);
ITPP_EXPORT cvec ifft(const cvec &in);
ITPP_EXPORT cvec ifft(const cvec &in, const int N);

ITPP_EXPORT void fft_real(const vec& in, cvec &out);
ITPP_EXPORT cvec fft_real(const vec& in);
ITPP_EXPORT cvec fft_real(const vec &in, const int N);
ITPP_EXPORT void ifft_real(const cvec &in, vec &out);
ITPP_EXPORT vec ifft_real(const cvec &in);
ITPP_EXPORT vec ifft_real(const cvec &in, const int N);





ITPP_EXPORT bool have_cosine_transforms();
ITPP_EXPORT void dct(const vec &in, vec &out);
ITPP_EXPORT vec dct(const vec &in);
ITPP_EXPORT vec dct(const vec &in, const int N);
ITPP_EXPORT void idct(const vec &in, vec &out);
ITPP_EXPORT vec idct(const vec &in);
ITPP_EXPORT vec idct(const vec &in, const int N);



template <class T> Vec<T> dht(const Vec<T> &v);
template <class T> void dht(const Vec<T> &vin, Vec<T> &vout);
template <class T> void self_dht(Vec<T> &v);

template <class T> Vec<T> dwht(const Vec<T> &v);
template <class T> void dwht(const Vec<T> &vin, Vec<T> &vout);
template <class T> void self_dwht(Vec<T> &v);

template <class T> Mat<T> dht2(const Mat<T> &m);
template <class T> Mat<T> dwht2(const Mat<T> &m);

template <class T>
Vec<T> dht(const Vec<T> &v)
{
  Vec<T> ret(v.size());
  dht(v, ret);
  return ret;
}

template <class T>
void bitrv(Vec<T> &out)
{
  int N = out.size();
  int j = 0;
  int N1 = N - 1;
  for(int i = 0; i < N1; ++i) {
    if(i < j) {
      T temp = out[j];
      out[j] = out[i];
      out[i] = temp;
    }
    int K = N / 2;
    while(K <= j) {
      j -= K;
      K /= 2;
    }
    j += K;
  }
}

template <class T>
void dht(const Vec<T> &vin, Vec<T> &vout)
{
  int N = vin.size();
  int m = levels2bits(N);
  it_assert_debug((1 << m) == N, "dht(): The vector size must be a power of two");

  vout.set_size(N);

  // This step is separated because it copies vin to vout
  for(int ib = 0; ib < N; ib += 2) {
    vout(ib) = vin(ib) + vin(ib + 1);
    vout(ib + 1) = vin(ib) - vin(ib + 1);
  }
  N /= 2;

  int l = 2;
  for(int i = 1; i < m; ++i) {
    N /= 2;
    int ib = 0;
    for(int k = 0; k < N; ++k) {
      for(int j = 0; j < l; ++j) {
        T t = vout(ib + j);
        vout(ib + j) += vout(ib + j + l);
        vout(ib + j + l) = t - vout(ib + j + l);
      }
      ib += 2 * l;
    }
    l *= 2;
  }

  vout /= static_cast<T>(std::sqrt(static_cast<double>(vin.size())));
}

template <class T>
void self_dht(Vec<T> &v)
{
  int N = v.size();
  int m = levels2bits(N);
  it_assert_debug((1 << m) == N, "self_dht(): The vector size must be a power "
                  "of two");

  int l = 1;
  for(int i = 0; i < m; ++i) {
    N /= 2;
    int ib = 0;
    for(int k = 0; k < N; ++k) {
      for(int j = 0; j < l; ++j) {
        T t = v(ib + j);
        v(ib + j) += v(ib + j + l);
        v(ib + j + l) = t - v(ib + j + l);
      }
      ib += 2 * l;
    }
    l *= 2;
  }

  v /= static_cast<T>(std::sqrt(static_cast<double>(v.size())));
}

template <class T>
Vec<T> dwht(const Vec<T> &v)
{
  Vec<T> ret(v.size());
  dwht(v, ret);
  return ret;
}

template <class T>
void dwht(const Vec<T> &vin, Vec<T> &vout)
{
  dht(vin, vout);
  bitrv(vout);
}


template <class T>
void self_dwht(Vec<T> &v)
{
  self_dht(v);
  bitrv(v);
}

template <class T>
Mat<T> dht2(const Mat<T> &m)
{
  Mat<T> ret(m.rows(), m.cols());
  Vec<T> v;

  for(int i = 0; i < m.rows(); ++i) {
    v = m.get_row(i);
    self_dht(v);
    ret.set_row(i, v);
  }
  for(int i = 0; i < m.cols(); ++i) {
    v = ret.get_col(i);
    self_dht(v);
    ret.set_col(i, v);
  }

  return transpose(ret);
}

template <class T>
Mat<T> dwht2(const Mat<T> &m)
{
  Mat<T> ret(m.rows(), m.cols());
  Vec<T> v;

  for(int i = 0; i < m.rows(); ++i) {
    v = m.get_row(i);
    self_dwht(v);
    ret.set_row(i, v);
  }
  for(int i = 0; i < m.cols(); ++i) {
    v = ret.get_col(i);
    self_dwht(v);
    ret.set_col(i, v);
  }

  return transpose(ret);
}


// ----------------------------------------------------------------------
// Instantiations
// ----------------------------------------------------------------------

ITPP_EXPORT_TEMPLATE template ITPP_EXPORT vec dht(const vec &v);
ITPP_EXPORT_TEMPLATE template ITPP_EXPORT cvec dht(const cvec &v);
ITPP_EXPORT_TEMPLATE template ITPP_EXPORT void dht(const vec &vin, vec &vout);
ITPP_EXPORT_TEMPLATE template ITPP_EXPORT void dht(const cvec &vin, cvec &vout);

ITPP_EXPORT_TEMPLATE template ITPP_EXPORT void self_dht(vec &v);
ITPP_EXPORT_TEMPLATE template ITPP_EXPORT void self_dht(cvec &v);

ITPP_EXPORT_TEMPLATE template ITPP_EXPORT vec dwht(const vec &v);
ITPP_EXPORT_TEMPLATE template ITPP_EXPORT cvec dwht(const cvec &v);
ITPP_EXPORT_TEMPLATE template ITPP_EXPORT void dwht(const vec &vin, vec &vout);
ITPP_EXPORT_TEMPLATE template ITPP_EXPORT void dwht(const cvec &vin, cvec &vout);

ITPP_EXPORT_TEMPLATE template ITPP_EXPORT void self_dwht(vec &v);
ITPP_EXPORT_TEMPLATE template ITPP_EXPORT void self_dwht(cvec &v);

ITPP_EXPORT_TEMPLATE template ITPP_EXPORT mat  dht2(const mat &m);
ITPP_EXPORT_TEMPLATE template ITPP_EXPORT cmat dht2(const cmat &m);

ITPP_EXPORT_TEMPLATE template ITPP_EXPORT mat  dwht2(const mat &m);
ITPP_EXPORT_TEMPLATE template ITPP_EXPORT cmat dwht2(const cmat &m);


} // namespace itpp

#endif // #ifndef TRANSFORMS_H