src/complexmat.cpp

   1 #include "complexmat.hpp"
   2
   3 ComplexMat::ComplexMat() : cols(0), rows(0), n_channels(0) {}
   4 ComplexMat::ComplexMat(int _rows, int _cols, int _n_channels) : cols(_cols), rows(_rows), n_channels(_n_channels)
   5 {
   6     p_data.resize(n_channels*cols*rows);
   7 }
   8
   9
  10 //assuming that mat has 2 channels (real, img)
  11 ComplexMat::ComplexMat(const cv::Mat & mat) : cols(mat.cols), rows(mat.rows), n_channels(1)
  12 {
  13     p_data = convert(mat);
  14 }
  15
  16 void ComplexMat::create(int _rows, int _cols, int _n_channels)
  17 {
  18     rows = _rows;
  19     cols = _cols;
  20     n_channels = _n_channels;
  21     p_data.resize(n_channels*cols*rows);
  22 }
  23
  24 void ComplexMat::create(int _rows, int _cols, int _n_channels, int _n_scales)
  25 {
  26     rows = _rows;
  27     cols = _cols;
  28     n_channels = _n_channels;
  29     n_scales = _n_scales;
  30     p_data.resize(n_channels*cols*rows);
  31 }
  32 // cv::Mat API compatibility
  33 cv::Size ComplexMat::size() { return cv::Size(cols, rows); }
  34 int ComplexMat::channels() { return n_channels; }
  35 int ComplexMat::channels() const { return n_channels; }
  36
  37 //assuming that mat has 2 channels (real, imag)
  38 void ComplexMat::set_channel(int idx, const cv::Mat & mat)
  39 {
  40     assert(idx >= 0 && idx < n_channels);
  41     for (int i = 0; i < rows; ++i){
  42         const std::complex<float> *row = mat.ptr<std::complex<float>>(i);
  43         for (int j = 0; j < cols; ++j)
  44             p_data[idx*rows*cols+i*cols+j]=row[j];
  45     }
  46 }
  47
  48
  49 void ComplexMat::sqr_norm(float *sums_sqr_norms) const
  50 {
  51     int n_channels_per_scale = n_channels/n_scales;
  52     int scale_offset = n_channels_per_scale*rows*cols;
  53     float sum_sqr_norm;
  54     for (int scale = 0; scale < n_scales; ++scale) {
  55         sum_sqr_norm = 0;
  56         for (int i = 0; i < n_channels_per_scale; ++i)
  57             for (auto lhs = p_data.begin()+i*rows*cols+scale*scale_offset; lhs != p_data.begin()+(i+1)*rows*cols+scale*scale_offset; ++lhs)
  58                 sum_sqr_norm += lhs->real()*lhs->real() + lhs->imag()*lhs->imag();
  59         sums_sqr_norms[scale] = sum_sqr_norm/static_cast<float>(cols*rows);
  60     }
  61     return;
  62 }
  63
  64 ComplexMat ComplexMat::sqr_mag() const
  65 {
  66     return mat_const_operator( [](std::complex<float> & c) { c = c.real()*c.real() + c.imag()*c.imag(); } );
  67 }
  68
  69 ComplexMat ComplexMat::conj() const
  70 {
  71     return mat_const_operator( [](std::complex<float> & c) { c = std::complex<float>(c.real(), -c.imag()); } );
  72 }
  73
  74 ComplexMat ComplexMat::sum_over_channels() const
  75 {
  76     assert(p_data.size() > 1);
  77
  78     int n_channels_per_scale = n_channels/n_scales;
  79     int scale_offset = n_channels_per_scale*rows*cols;
  80
  81     ComplexMat result(this->rows, this->cols, n_scales);
  82     for (int scale = 0; scale < n_scales; ++scale) {
  83         std::copy(p_data.begin()+scale*scale_offset,p_data.begin()+rows*cols+scale*scale_offset, result.p_data.begin()+scale*rows*cols);
  84         for (int i = 1; i < n_channels_per_scale; ++i) {
  85             std::transform(result.p_data.begin()+scale*rows*cols, result.p_data.begin()+(scale+1)*rows*cols, p_data.begin()+i*rows*cols+scale*scale_offset,
  86                            result.p_data.begin()+scale*rows*cols, std::plus<std::complex<float>>());
  87         }
  88     }
  89     return result;
  90 }
  91
  92 //return 2 channels (real, imag) for first complex channel
  93 cv::Mat ComplexMat::to_cv_mat() const
  94 {
  95     assert(p_data.size() >= 1);
  96     return channel_to_cv_mat(0);
  97 }
  98 // return a vector of 2 channels (real, imag) per one complex channel
  99 std::vector<cv::Mat> ComplexMat::to_cv_mat_vector() const
 100 {
 101     std::vector<cv::Mat> result;
 102     result.reserve(n_channels);
 103
 104     for (int i = 0; i < n_channels; ++i)
 105         result.push_back(channel_to_cv_mat(i));
 106
 107     return result;
 108 }
 109
 110 std::complex<float>* ComplexMat::get_p_data() const
 111 {
 112     return p_data.data();
 113 }
 114
 115 //element-wise per channel multiplication, division and addition
 116 ComplexMat ComplexMat::operator*(const ComplexMat & rhs) const
 117 {
 118     return ComplexMat::mat_mat_operator( [](std::complex<float> & c_lhs, const std::complex<float> & c_rhs) { c_lhs *= c_rhs; }, rhs);
 119 }
 120 ComplexMat ComplexMat::operator/(const ComplexMat & rhs) const
 121 {
 122     return ComplexMat::mat_mat_operator( [](std::complex<float> & c_lhs, const std::complex<float> & c_rhs) { c_lhs /= c_rhs; }, rhs);
 123 }
 124 ComplexMat ComplexMat::operator+(const ComplexMat & rhs) const
 125 {
 126     return ComplexMat::mat_mat_operator( [](std::complex<float> & c_lhs, const std::complex<float> & c_rhs)  { c_lhs += c_rhs; }, rhs);
 127 }
 128
 129 //multiplying or adding constant
 130 ComplexMat ComplexMat::operator*(const float & rhs) const
 131 {
 132     return ComplexMat::mat_const_operator( [&rhs](std::complex<float> & c) { c *= rhs; });
 133 }
 134 ComplexMat ComplexMat::operator+(const float & rhs) const
 135 {
 136     return ComplexMat::mat_const_operator( [&rhs](std::complex<float> & c) { c += rhs; });
 137 }
 138
 139 //multiplying element-wise multichannel by one channel mats (rhs mat is with one channel)
 140 ComplexMat ComplexMat::mul(const ComplexMat & rhs) const
 141 {
 142     return ComplexMat::matn_mat1_operator( [](std::complex<float> & c_lhs, const std::complex<float> & c_rhs) { c_lhs *= c_rhs; }, rhs);
 143 }
 144
 145 //multiplying element-wise multichannel by one channel mats (rhs mat is with multiple channel)
 146 ComplexMat ComplexMat::mul2(const ComplexMat & rhs) const
 147 {
 148     return ComplexMat::matn_mat2_operator( [](std::complex<float> & c_lhs, const std::complex<float> & c_rhs) { c_lhs *= c_rhs; }, rhs);
 149 }
 150
 151
 152 //convert 2 channel mat (real, imag) to vector row-by-row
 153 std::vector<std::complex<float>> ComplexMat::convert(const cv::Mat & mat)
 154 {
 155     std::vector<std::complex<float>> result;
 156     result.reserve(mat.cols*mat.rows);
 157     for (int y = 0; y < mat.rows; ++y) {
 158         const float * row_ptr = mat.ptr<float>(y);
 159         for (int x = 0; x < 2*mat.cols; x += 2){
 160             result.push_back(std::complex<float>(row_ptr[x], row_ptr[x+1]));
 161         }
 162     }
 163     return result;
 164 }
 165
 166 ComplexMat ComplexMat::mat_mat_operator(void (*op)(std::complex<float> & c_lhs, const std::complex<float> & c_rhs), const ComplexMat & mat_rhs) const
 167 {
 168     assert(mat_rhs.n_channels == n_channels && mat_rhs.cols == cols && mat_rhs.rows == rows);
 169
 170     ComplexMat result = *this;
 171     for (int i = 0; i < n_channels; ++i) {
 172         auto lhs = result.p_data.begin()+i*rows*cols;
 173         auto rhs = mat_rhs.p_data.begin()+i*rows*cols;
 174         for ( ; lhs != result.p_data.begin()+(i+1)*rows*cols; ++lhs, ++rhs)
 175             op(*lhs, *rhs);
 176     }
 177
 178     return result;
 179 }
 180 ComplexMat ComplexMat::matn_mat1_operator(void (*op)(std::complex<float> & c_lhs, const std::complex<float> & c_rhs), const ComplexMat & mat_rhs) const
 181 {
 182     assert(mat_rhs.n_channels == 1 && mat_rhs.cols == cols && mat_rhs.rows == rows);
 183
 184     ComplexMat result = *this;
 185     for (int i = 0; i < n_channels; ++i) {
 186         auto lhs = result.p_data.begin()+i*rows*cols;
 187         auto rhs = mat_rhs.p_data.begin();
 188         for ( ; lhs != result.p_data.begin()+(i+1)*rows*cols; ++lhs, ++rhs)
 189             op(*lhs, *rhs);
 190     }
 191
 192     return result;
 193 }
 194 ComplexMat ComplexMat::matn_mat2_operator(void (*op)(std::complex<float> & c_lhs, const std::complex<float> & c_rhs), const ComplexMat & mat_rhs) const
 195 {
 196     assert(mat_rhs.n_channels == n_channels/n_scales && mat_rhs.cols == cols && mat_rhs.rows == rows);
 197
 198     int n_channels_per_scale = n_channels/n_scales;
 199     int scale_offset = n_channels_per_scale*rows*cols;
 200     ComplexMat result = *this;
 201     for (int i = 0; i < n_scales; ++i) {
 202         for (int j = 0; j < n_channels_per_scale; ++j) {
 203             auto lhs = result.p_data.begin()+(j*rows*cols)+(i*scale_offset);
 204             auto rhs = mat_rhs.p_data.begin()+(j*rows*cols);
 205             for ( ; lhs != result.p_data.begin()+((j+1)*rows*cols)+(i*scale_offset); ++lhs, ++rhs)
 206                 op(*lhs, *rhs);
 207         }
 208     }
 209
 210     return result;
 211 }
 212 ComplexMat ComplexMat::mat_const_operator(const std::function<void(std::complex<float> & c_rhs)> & op) const
 213 {
 214     ComplexMat result = *this;
 215     for (int i = 0; i < n_channels; ++i)
 216         for (auto lhs = result.p_data.begin()+i*rows*cols; lhs != result.p_data.begin()+(i+1)*rows*cols; ++lhs)
 217             op(*lhs);
 218     return result;
 219 }
 220
 221 cv::Mat ComplexMat::channel_to_cv_mat(int channel_id) const
 222 {
 223     cv::Mat result(rows, cols, CV_32FC2);
 224     for (int y = 0; y < rows; ++y) {
 225         std::complex<float> * row_ptr = result.ptr<std::complex<float>>(y);
 226         for (int x = 0; x < cols; ++x){
 227             row_ptr[x] = p_data[channel_id*rows*cols+y*cols+x];
 228         }
 229     }
 230     return result;
 231 }