I'm learning the fourier transform but I'm confused by a problem.
As an exercise, I've created a Python program to model the effect of a low pass RC filter on a digital waveform. The prog creates a waveform, converts it to the frequency domain with a fourier transform, applies the gain of an RC filter to each harmonic [ gain = 1/(1+jw/w0) ] then applies the inverse fourier transform to show the filtered waveform.
Image of output:
By choosing an extremely high cutoff frequency, I recovered the original waveform (shown in black) as expected. When choosing cutoff frequencies that were small multiples of the repetition rate of the waveform, I got waveforms with softened edges (blue and green) that were the expected shape except for being half the height of the original waveform which was unexpected.
To try to figure out why the waveforms were half the height, I tried increasing the cutoff frequency to see what would happen. At fairly high cutoff frequencies ( >100 times repetition rate ), I got crazy waveforms (red and orange) that were obviously wrong. Why?
The fourier transform that I'm using is the well-respected and well-tested 30 year old Fortran FFTPACK, that's wrapped into Python's "numpy" library, so it's a safe assumption that the tools are good but that I've misapplied it somehow. Any fourier transform experts here?
The short Python program below will work as-is in the default installation of recent Debian-derived distros, otherwise the "PIL" image library and the "numpy" maths library may need to be installed first. The relevant part of the program is the half-dozen lines in the "low_pass" function. Anybody know what's wrong with it?
#! /usr/bin/env python2 import numpy, Image, ImageDraw
# define input waveform as 010100 trapezoid lo = numpy.array( [0] * 400 ) hi = numpy.array(range(20) + [20] * 400 + range(19,-1,-1)) / 20.0 wave = numpy.concatenate( [lo,hi,lo,hi,lo,lo] )
# initialise PIL image width = len(wave) height = 2000 im = Image.new("RGB",(width,height),"white") draw = ImageDraw.Draw(im)
# Func for drawing a waveform to the image def draw_wave( waveform, hue, offset ): ypoints = (- waveform + 1.2 + offset / 4.0) * height / 5 xy = zip( numpy.arange(len(waveform)), ypoints) draw.line(xy, fill=hue, width=10)
# Func for low-pass RC filter of waveform def low_pass( waveform, freq, fmax, hue, offset ): fftwaveform = numpy.fft.fft(waveform) # fast fourier tr harmonics = numpy.arange(len(fftwaveform)) # 0,1,2,3,4... gains = 1.0 / (1.0 + 1j * (harmonics * freq) / fmax) # RC filter gains fftfiltered = fftwaveform * gains # filter the wave filtered = numpy.fft.ifft(fftfiltered).real # inverse fft draw_wave(filtered, hue, offset) # draw to image
# input_waveform, base_freq, RC_filter_freq, screen_colour, screen_offset low_pass(wave, 1, 1.0e10, "black", 0) low_pass( wave, 1, 8, "blue", 4) low_pass( wave, 1, 20, "green", 6) low_pass( wave, 1, 400, "red", 10) low_pass( wave, 1, 4000, "orange", 14)
smallim = im.resize( (width/4,height/4), Image.ANTIALIAS) smallim.save("image.png") smallim.show()