Download or view fourierTemp.frink in plain text format
//use englewoodarray.frink
a = new array
for line = lines["file:hourly.txt"]
a.push[eval[line]]
println[length[a]]
// Transform time-domain data into frequency-domain data
f = DFT[a]
g = new graphics
// fprime will contain just selected frequency components.
fprime = new array[]
len = length[f]
zeroUnit = 0 f@0
threshold = .1
vscale = 20
c = 0
str = ""
for y = f
{
x = (c + len/2) mod len // Draw low-freq components together in middle
// Draw phases
g.color[0,1,0]
phase = arctan[Im[y],Re[y]]
g.fillRectSides[x,0,x+1,-5 phase]
// Draw magnitudes
g.color[0,0,0,.7]
mag = abs[y]
if (c != 0) // Don't draw DC component because it's huge.
g.fillRectSides[x,0,x+1,-vscale mag]
if (mag > threshold) // Decide which components to keep.
{
daysphase = (phase/(2 pi) * len)
mag2 = 2 mag
println["Used component $c, mag is $mag, phase is $daysphase days"]
if c == 0
str = format[mag2/2, 1, 6]
if c > 0 && c < len/2
str = str + " +\n " + format[mag2,1,6] + " cos[2 pi $c/$len d - " + format[phase,1,6] + "]"
fprime.push[y]
} else
fprime.push[zeroUnit]
c = c + 1
}
// Draw threshold
g.color[0,0,1,.5]
g.line[0,-vscale threshold,len,-vscale threshold]
println[str]
g.show[]
g.write["fourierEnglewood.svg",800,600]
// Turn the frequency-domain data back into a time-domain waveform
t = InverseDFT[fprime]
// Plot the time-domain data that was generated from the truncated
// frequency domain.
g2 = new graphics
c = 0
scale = 3
for y = t
{
g2.fillRectCenter[c, -scale Re[y], 1, scale]
c = c+1
}
// Plot original data
g2.color[0,0,1,.7]
c = 0
for y = a
{
g2.fillRectCenter[c, -scale Re[y], 1, scale]
c = c+1
}
// Test that sin approximation gives the same values
g2.color[1,0,0,.5]
for d = 0 to len-1
g2.fillRectCenter[d, -scale eval[str], 1, scale]
g2.show[]
g2.write["recreatedEnglewood.svg", 800, 600]
Download or view fourierTemp.frink in plain text format
This is a program written in the programming language Frink.
For more information, view the Frink
Documentation or see More Sample Frink Programs.
Alan Eliasen was born 19966 days, 14 hours, 34 minutes ago.