Ahhhhh, looking back I see what Pierre-André did, the penny drops. Sorry I miss-read.I'm only interested in this part;>zeropadded_fft_A = fft(seq_A, n=2**(ceil(log(len(seq_A),2))+1)) >zeropadded_fft_B = fft(seq_B, n=2**(ceil(log(len(seq_B),2))+1))>You could >remove the "+1" above to get faster results, but then that may >lead to >unwanted frequencies (coming from the fact that fft assumes >periodic signals, >read online about zero-padding). I did it in the following way and it seems to work. I've seen an increase in performance for sure. My Python sucks by the way :)def pad_to_power_of_two(audio): ''' Calculates next power of two and pads the audio sample array with zeros ''' # RIFF spec specifies: ckSize A 32-bit unsigned value identifying the size of ckData. max_shift = 32 shift = 1 power = len(audio) - 1 print "{0:b}".format(power) # Test if a power of two. Turn off the right-most bit (x & (x - 1)) and zero test while ((power + 1) & power) and (shift <= max_shift): power |= power >> shift shift *= 2 print "new power: " + "{0:b}".format(power) power = power + 1 print "{0:b}".format(power) ext = np.array([0x00] * (power - len(audio))) new_array = np.concatenate((audio, ext), axis=0) print "Next power of 2 greater than actual array length is " + str(power) print "Extending current array length (" + str(len(audio)) + ") by " + str(power - len(audio)) + " null bytes" print "New array length is " + str(len(new_array)) return new_array I will try using the suggested methods and see if I can boost the performance further.Thanks. > From: oscar.j.benja...@gmail.com > Date: Tue, 1 Sep 2015 16:14:41 +0100 > To: numpy-discussion@scipy.org > Subject: Re: [Numpy-discussion] Numpy FFT.FFT slow with certain samples > > On 1 September 2015 at 11:38, Joseph Codadeen <jdm...@hotmail.com> wrote: > > > >> And while you zero-pad, you can zero-pad to a sequence that is a power of > >> two, thus preventing awkward factorizations. > > > > Does numpy have an easy way to do this, i.e. for a given number, find the > > next highest number (within a range) that could be factored into small, > > prime numbers as Phil explained? It would help if it gave a list, > > prioritised by number of factors. > > Just use the next power of 2. Pure powers of 2 are the most efficient > for FFT algorithms so it potentially works out better than finding a > smaller but similarly composite size to pad to. Finding the next power > of 2 is easy to code and never a bad choice. > > To avoid the problems mentioned about zero-padding distorting the FFT > you can use the Bluestein transform as below. This pads up to a power > of two (greater than 2N-1) but does it in a special way that ensures > that it is still calculating the exact (up to floating point error) > same DFT: > > from numpy import array, exp, pi, arange, concatenate > from numpy.fft import fft, ifft > > def ceilpow2(N): > ''' > >>> ceilpow2(15) > 16 > >>> ceilpow2(16) > 16 > ''' > p = 1 > while p < N: > p *= 2 > return p > > def fftbs(x): > ''' > >>> data = [1, 2, 5, 2, 5, 2, 3] > >>> from numpy.fft import fft > >>> from numpy import allclose > >>> from numpy.random import randn > >>> for n in range(1, 1000): > ... data = randn(n) > ... assert allclose(fft(data), fftbs(data)) > ''' > N = len(x) > x = array(x) > > n = arange(N) > b = exp((1j*pi*n**2)/N) > a = x * b.conjugate() > > M = ceilpow2(N) * 2 > A = concatenate((a, [0] * (M - N))) > B = concatenate((b, [0] * (M - 2*N + 1), b[:0:-1])) > C = ifft(fft(A) * fft(B)) > c = C[:N] > return b.conjugate() * c > > if __name__ == "__main__": > import doctest > doctest.testmod() > > > -- > Oscar > _______________________________________________ > NumPy-Discussion mailing list > NumPy-Discussion@scipy.org > http://mail.scipy.org/mailman/listinfo/numpy-discussion
_______________________________________________ NumPy-Discussion mailing list NumPy-Discussion@scipy.org http://mail.scipy.org/mailman/listinfo/numpy-discussion