<div dir="ltr"><div><div><div>Hi Tim,<br><br></div>Brilliant!  Many thanks...  I think this is exactly what I need, I owe you a beer (or other beverage of your choice).<br><br></div>I'm now going to lock myself in the basement until I can work out an implementation of this for my use-case :)<br>
<br></div>/Carl<br></div><div class="gmail_extra"><br><br><div class="gmail_quote">On Tue, Sep 3, 2013 at 9:05 PM, Cera, Tim <span dir="ltr"><<a href="mailto:tim@cerazone.net" target="_blank">tim@cerazone.net</a>></span> wrote:<br>
<blockquote class="gmail_quote" style="margin:0 0 0 .8ex;border-left:1px #ccc solid;padding-left:1ex"><div class="im">> I am trying to take the rfft of a numpy array, like this:<br>
>>>> my_rfft = numpy.fft.rfft(my_numpy_array)<br>
><br>
> and replace the amplitudes that can be obtained with:<br>
><br>
>>>> my_amplitudes = numpy.abs(my_rfft)<br>
><br>
> with amplitudes from an arbitrary numpy array's rFFT, which is to then be<br>
> converted back using numpy.fft.irfft .  Alternately, some future plans will<br>
> involve having to modify individual array element amplitudes directly based<br>
> on other parameters.  I would think that modifying and re-synthesizing<br>
> signals using FFT is a fairly common use-case, but my attempts at Googling<br>
> example code have been fruitless.<br>
<br>
</div>I have FFT transform filter in my tidal analysis package.   See<br>
<a href="http://sourceforge.net/apps/mediawiki/tappy/index.php?title=CompareTidalFilters" target="_blank">http://sourceforge.net/apps/mediawiki/tappy/index.php?title=CompareTidalFilters</a><br>
for a comparison and short description.<br>
<br>
See my function below.  My earlier self made some poor variable name<br>
choices.  The 'low_bound' variable is actually where frequencies<br>
greater are set to zero ('factor[freq > low_bound] = 0.0'), then<br>
factor is ramped from 0 at 'low_bound' to 1 at 'high_bound'.  To<br>
filter out tidal signals if your water elevations are hourly then<br>
'low_bound' = 1/30.0 and 'high_bound' = 1/40.0.  Having this gradual<br>
change in the frequency domain rather than an abrupt change makes a<br>
better filter.<br>
<br>
def fft_lowpass(nelevation, low_bound, high_bound):<br>
     """ Performs a low pass filter on the nelevation series.<br>
     low_bound and high_bound specifies the boundary of the filter.<br>
     """<br>
     import numpy.fft as F<br>
     if len(nelevation) % 2:<br>
         result = F.rfft(nelevation, len(nelevation))<br>
     else:<br>
         result = F.rfft(nelevation)<br>
     freq = F.fftfreq(len(nelevation))[:len(nelevation)/2]<br>
     factor = np.ones_like(result)<br>
     factor[freq > low_bound] = 0.0<br>
<br>
     sl = np.logical_and(high_bound < freq, freq < low_bound)<br>
<br>
     a = factor[sl]<br>
     # Create float array of required length and reverse<br>
     a = np.arange(len(a) + 2).astype(float)[::-1]<br>
<br>
     # Ramp from 1 to 0 exclusive<br>
     a = (a/a[0])[1:-1]<br>
<br>
     # Insert ramp into factor<br>
     factor[sl] = a<br>
<br>
     result = result * factor<br>
     print 'result=', len(result)<br>
     relevation = F.irfft(result, len(nelevation))<br>
     print 'result=', len(relevation)<br>
     return relevation<br>
<br>
<br>
Kindest regards,<br>
Tim<br>
_______________________________________________<br>
NumPy-Discussion mailing list<br>
<a href="mailto:NumPy-Discussion@scipy.org">NumPy-Discussion@scipy.org</a><br>
<a href="http://mail.scipy.org/mailman/listinfo/numpy-discussion" target="_blank">http://mail.scipy.org/mailman/listinfo/numpy-discussion</a><br>
</blockquote></div><br></div>