[Numpy-discussion] (no subject)
carl.canuck.official at gmail.com
Tue Sep 3 22:04:13 EDT 2013
Brilliant! Many thanks... I think this is exactly what I need, I owe you
a beer (or other beverage of your choice).
I'm now going to lock myself in the basement until I can work out an
implementation of this for my use-case :)
On Tue, Sep 3, 2013 at 9:05 PM, Cera, Tim <tim at cerazone.net> wrote:
> > I am trying to take the rfft of a numpy array, like this:
> >>>> my_rfft = numpy.fft.rfft(my_numpy_array)
> > and replace the amplitudes that can be obtained with:
> >>>> my_amplitudes = numpy.abs(my_rfft)
> > with amplitudes from an arbitrary numpy array's rFFT, which is to then be
> > converted back using numpy.fft.irfft . Alternately, some future plans
> > involve having to modify individual array element amplitudes directly
> > on other parameters. I would think that modifying and re-synthesizing
> > signals using FFT is a fairly common use-case, but my attempts at
> > example code have been fruitless.
> I have FFT transform filter in my tidal analysis package. See
> for a comparison and short description.
> See my function below. My earlier self made some poor variable name
> choices. The 'low_bound' variable is actually where frequencies
> greater are set to zero ('factor[freq > low_bound] = 0.0'), then
> factor is ramped from 0 at 'low_bound' to 1 at 'high_bound'. To
> filter out tidal signals if your water elevations are hourly then
> 'low_bound' = 1/30.0 and 'high_bound' = 1/40.0. Having this gradual
> change in the frequency domain rather than an abrupt change makes a
> better filter.
> def fft_lowpass(nelevation, low_bound, high_bound):
> """ Performs a low pass filter on the nelevation series.
> low_bound and high_bound specifies the boundary of the filter.
> import numpy.fft as F
> if len(nelevation) % 2:
> result = F.rfft(nelevation, len(nelevation))
> result = F.rfft(nelevation)
> freq = F.fftfreq(len(nelevation))[:len(nelevation)/2]
> factor = np.ones_like(result)
> factor[freq > low_bound] = 0.0
> sl = np.logical_and(high_bound < freq, freq < low_bound)
> a = factor[sl]
> # Create float array of required length and reverse
> a = np.arange(len(a) + 2).astype(float)[::-1]
> # Ramp from 1 to 0 exclusive
> a = (a/a)[1:-1]
> # Insert ramp into factor
> factor[sl] = a
> result = result * factor
> print 'result=', len(result)
> relevation = F.irfft(result, len(nelevation))
> print 'result=', len(relevation)
> return relevation
> Kindest regards,
> NumPy-Discussion mailing list
> NumPy-Discussion at scipy.org
-------------- next part --------------
An HTML attachment was scrubbed...
More information about the NumPy-Discussion