import std.stdio;
extern(C) __gshared string[] rt_options =
[ "gcopt=minPoolSize:2K maxPoolSize:1M incPoolSize:4K" ];
auto fastAbs(float x) @system
{
uint p = *(cast(uint*) &x);
p &=0x7fffffff;
return *(cast(float*) &p);
}
auto renormalize(Range)(Range r)
{
import std.algorithm;
auto maximalAbsoluteValue = r
.map!(a => fastAbs(a))
.fold!max;
return r
.map!(a => a / maximalAbsoluteValue);
}
void main()
{
import std.algorithm;
import std.complex;
import std.math;
import std.numeric;
import std.random;
import std.range;
import ppmformats;
auto rng = Random(unpredictableSeed);
auto signalWithNoise(float x)
{
return (3.3 * sin(2*x)) + uniform(-0.3, 0.3, rng);
}
alias squareWave = delegate(float x) {
return
((4 * 10) / (PI)) * (sin(x) + (sin(3 * x) / (3 * x)) + (sin(5 * x) /(5 * x)) + (sin(7 * x) / (7 * x)) + (sin(9 * x) / (9 * x)) + (sin(11 * x) / (11 * x))) + (sin(13 * x) / (13 * x)) + (sin(15 * x) / (15 * x)) + (sin(17 * x) / (17 * x));
};
auto signal = iota(0, 2_048)
.map!squareWave;
auto normalizedFFT(Range)(Range r)
{
return
r
.fft
.map!(a => 20.0f * log10(abs(a ^^ 2)))
.renormalize;
}
auto pd = signal
.slide(1_024)
.map!(a => normalizedFFT(a));
auto img = new P6Image(1024, 1024);
foreach (i, e; pd.enumerate)
{
foreach (j, s; e.enumerate)
{
img[i, j] = new RGBColor(
cast(int) (0.8 * 255.0 * abs(s)),
0,
255 - cast(int) (255 * abs(s))
);
}
}
img.save(`fft.ppm`);
}