miniprojet/tests/src/math.hpp

125 lines
3 KiB
C++

#pragma once
#include <vector>
#include <complex>
#include <opencv2/opencv.hpp>
#include <iterator>
#include <cmath>
namespace math {
using complex = std::complex<double>;
using signal = std::vector<double>;
using csignal = std::vector<complex>;
using contour = std::vector<cv::Point>;
constexpr double pi() {return std::atan(1)*4;}
csignal cont2sig(const contour& cont) {
csignal sig;
auto sig_it = sig.begin();
auto cont_it = cont.begin();
for (auto cont_it = cont.begin(); cont_it != cont.end(); ++cont_it) {
*(sig_it++) = complex((*cont_it).x, (*cont_it).y);
}
return sig;
};
complex mean(const csignal& sig) {
complex res = 0;
for (auto x: sig) {
res += x;
}
return complex(res.real()/sig.size(), res.imag()/sig.size());
};
//TODO: implémenter la fonction
csignal diff(const csignal& input, complex mean) {
return csignal();
}
//TODO implémenter la fft
csignal fft_rec(const csignal& input) {
int size = input.size();
if (size == 1) {
//TODO: que faire dans ce cas ?
return csignal();
} else {
csignal odd;
csignal even;
std::back_insert_iterator<csignal> odd_back_it(odd);
std::back_insert_iterator<csignal> even_back_it(even);
bool insert_in_even = false;
for (auto it = input.begin(); it != input.end(); ++it) {
if (insert_in_even) {
*even_back_it++ = *it;
insert_in_even = false;
} else {
*odd_back_it++ = *it;
insert_in_even = true;
}
}
csignal odd_fft = fft_rec(odd);
csignal even_fft = fft_rec(even);
csignal res;
res.reserve(size);
for (int k = 0; k<size/2; ++k) {
complex t = std::exp(complex(0, -2*pi()*k/size)) * odd[k];
res[k] = even[k] + t;
res[size/2+k] = even[k] - t;
}
return res;
}
}
csignal fft(const csignal& input) {
//TODO: s'assurer que le signal est bien formé (i.e. bonne taille)
return fft_rec(input);
};
contour coef2cont(const csignal& tfd, complex mean, int size, int cmax) {
contour cont;
auto tf_it = tfd.begin();
auto cont_it = cont.begin();
//TODO: trouver les bonnes valeurs
int kmin = 0;
int kmax = 2*cmax;
for (int m=0; m<tfd.size(); ++m) {
//TODO: retrouver la formule
complex sum;
for (int k=kmin; k<kmax; ++k) {
sum += tfd[k]*std::exp(complex(0, 2*pi()*k*m/tfd.size()));
}
complex zm = mean + sum;
*(cont_it++) = cv::Point(zm.real(), zm.imag());
}
return cont;
};
contour simplify_contour(const contour& cont, int cmax) {
contour res;
csignal z = cont2sig(cont);
complex zm = mean(z);
//TODO: fft(diff(z, zm));
csignal tfd = fft(z);
return coef2cont(tfd, zm, 0, cmax);
};
int max_cont(const std::vector<contour>& contours) {
int max = 0;
int id = 0;
for (int i=0; i<contours.size(); ++i) {
if (contours[i].size() > max) {
max = contours[i].size();
id = i;
}
}
return id;
};
}