Skip to content

Instantly share code, notes, and snippets.

@bsdis
Created September 17, 2016 05:14
Show Gist options
  • Save bsdis/db16cf38f599cd7e83f1f84f5f8e2f03 to your computer and use it in GitHub Desktop.
Save bsdis/db16cf38f599cd7e83f1f84f5f8e2f03 to your computer and use it in GitHub Desktop.
numba error
import sys
import pyfftw
import numpy as np
import scipy.signal
import scipy.fftpack
import matplotlib.pyplot as plt
from timeit import Timer
import time
from numba import jit
def create_signal(fs = 10000, N = 20000):
# Timebase
t = np.linspace(1,N,num=N)
# Signal frequencies, f1 is musical 'A'
f1 = 440
f2 = 8 * f1
# Signal power, 20Vrms here, try making this smaller or bigger
a = 20.0 * np.sqrt(2)
# Noise power, equal to 0.1 V**2/Hz, make this smaller like 0.0001
noise = 0.0001 * fs / 2
# Construct the data
d = a * np.sin(2*np.pi*t*f1/fs) + 2 * a * np.sin(2*np.pi*t*f2/fs)
d += np.random.normal(scale = np.sqrt(noise), size = t.shape)
return d
@jit(nopython=True)
def runpyfftw(s, dbins, prf, win, ws, fft):
rdysig = np.tile(win,(dbins,1)) * scipy.signal.detrend(s)
for i in range(rdysig.shape[0]):
v = fft( rdysig[i,:].flatten() )
psd = 2 * (np.real(v) * np.imag( v ))/(prf * ws)
def main():
pyfftw.interfaces.cache.enable()
prf = 10000
secs = 2
fftsize = 10000
dbins = 800
a = pyfftw.empty_aligned(prf * secs, dtype='float32')
fft = pyfftw.builders.fft(a)
s = np.tile( create_signal(prf, prf * secs)[np.newaxis] , (dbins,1) )
win = np.hanning(prf * secs)
ws = np.sum(win**2)
t = time.process_time()
runpyfftw(s, dbins, prf, win, ws, fft)
print(time.process_time() - t)
Sign up for free to join this conversation on GitHub. Already have an account? Sign in to comment