123456789101112131415161718192021222324252627282930313233343536373839404142434445464748495051525354555657585960616263646566676869707172737475767778798081828384858687888990919293949596979899100101102103104105106107108109110111112113114115116117118119120121122123124125126127128129130131132133 |
- #include "simple/compress/fft.hpp"
- #include "simple/support/random.hpp"
- #include <cassert>
- #include <vector>
- #include <cmath>
- #include <algorithm>
- #include <numeric>
- #include <complex>
- #include <random>
- #include <iostream>
- #include <iomanip>
- using namespace simple::compress;
- const float tau = 2*std::acos(-1);
- template <typename Waveform>
- void DFT(std::size_t size, Waveform waveform)
- {
- std::vector<float> wave;
- wave.resize(size);
- std::iota(wave.begin(), wave.end(), 0);
- std::transform(wave.begin(), wave.end(), wave.begin(), [&](auto x) { return x/size; });
- std::transform(wave.begin(), wave.end(), wave.begin(), waveform);
- std::vector<std::complex<float>> control;
- control.resize(wave.size());
- for(std::size_t i = 0; i != wave.size(); ++i)
- {
- for(std::size_t j = 0; j != wave.size(); ++j)
- {
- control[i].real( control[i].real() + wave[j] * std::cos(tau*j*i/wave.size()) / control.size());
- control[i].imag( control[i].imag() + wave[j] * std::sin(tau*j*i/wave.size()) / control.size());
- }
- }
- std::vector<std::complex<float>> frequencies;
- frequencies.resize(wave.size());
- fft(std::multiplies<>{}, std::polar(1.f,tau/wave.size()), wave, frequencies);
- for(auto&& x : frequencies)
- x = std::conj(x) / float(frequencies.size());
- std::vector<std::complex<float>> aaanback;
- aaanback.resize(wave.size());
- fft(std::multiplies<>{}, std::polar(1.f,tau/wave.size()), frequencies, aaanback);
- for(std::size_t i = 0; i != aaanback.size(); ++i)
- {
- assert(std::abs(aaanback[i].real() - wave[i]) < 0.001);
- }
- for(std::size_t i = 0; i != frequencies.size(); ++i)
- {
- assert(std::abs(std::abs(frequencies[i]) - std::abs(control[i])) < 0.001);
- }
- #if defined SIMPLE_SUPPORT_DEBUG_HPP
- if(size <= 32)
- {
- std::cerr << std::fixed;
- std::cerr << std::setprecision(3);
- simple::support::print('\n');
- simple::support::print("IN : [ ");
- for(auto&& x : wave)
- simple::support::print(x, " ");
- simple::support::print("]\n");
- simple::support::print("DEC: [ ");
- for(auto&& x : aaanback)
- simple::support::print(x.real(), " ");
- simple::support::print("]\n");
- simple::support::print("ENC: [ ");
- for(auto&& x : frequencies)
- simple::support::print(std::abs(x), " ");
- simple::support::print("]\n");
- simple::support::print("CNR: [ ");
- for(auto&& x : control)
- simple::support::print(std::abs(x), " ");
- simple::support::print("]\n");
- }
- #endif
- }
- int main()
- {
- DFT(2, [&](auto x) { return std::sin(tau*x) + std::cos(tau*x*2); });
- DFT(2, [&](auto x) { return std::sin(tau*x) + std::cos(tau*x); });
- DFT(2, [&](auto x) { return std::sin(tau*x); });
- DFT(2, [&](auto x) { return std::cos(tau*x); });
- DFT(4, [&](auto x) { return std::sin(tau*x) + std::cos(tau*x*2); });
- DFT(4, [&](auto x) { return std::sin(tau*x) + std::cos(tau*x); });
- DFT(4, [&](auto x) { return std::sin(tau*x); });
- DFT(4, [&](auto x) { return std::cos(tau*x); });
- DFT(8, [&](auto x) { return std::sin(tau*x) + std::cos(tau*x*2); });
- DFT(8, [&](auto x) { return std::sin(tau*x) + std::cos(tau*x); });
- DFT(8, [&](auto x) { return std::sin(tau*x); });
- DFT(8, [&](auto x) { return std::cos(tau*x); });
- DFT(16, [&](auto x) { return std::sin(tau*x) + std::cos(tau*x*2); });
- DFT(16, [&](auto x) { return std::sin(tau*x) + std::cos(tau*x); });
- DFT(16, [&](auto x) { return std::sin(tau*x); });
- DFT(16, [&](auto x) { return std::cos(tau*x); });
- DFT(1024, [&](auto x) { return std::sin(tau*x) + std::cos(tau*x*2); });
- DFT(1024, [&](auto x) { return std::sin(tau*x) + std::cos(tau*x); });
- DFT(1024, [&](auto x) { return std::sin(tau*x); });
- DFT(1024, [&](auto x) { return std::cos(tau*x); });
- DFT(2048, [&](auto x) { return std::sin(tau*x) + std::cos(tau*x*2); });
- DFT(2048, [&](auto x) { return std::sin(tau*x) + std::cos(tau*x); });
- DFT(2048, [&](auto x) { return std::sin(tau*x); });
- DFT(2048, [&](auto x) { return std::cos(tau*x); });
- DFT(2048, [&](auto x) { return std::sin(tau*13*x) + std::cos(tau*345*x); });
- DFT(2048, [&](auto x) { return std::sin(tau*134*x) + std::sin(tau*345*x) + std::cos(tau*789*x); });
- auto seed = std::random_device{}();
- std::cout << "DFT random test seed: " << std::hex << std::showbase << seed << std::dec << std::endl;
- simple::support::random::engine::tiny<unsigned long long> random{seed};
- std::uniform_real_distribution<double> number{};
- DFT(2048, [&](auto) { return number(random); });
- return 0;
- }
|