import numpy as np
def make_filter (x, noise_floor_db, noise_gain_db, filter_length, window_size):
assert noise_gain_db < 0
assert filter_length <= len (x)
assert filter_length <= window_size
filter_length /= 2
fft = np.fft.fft (x * (1.0 - np.cos (np.linspace (0.0, 2*np.pi, len (x)))) / 2.0)
t = len (fft)**2 / 4.0 * (10.0 ** (noise_floor_db/20.0))
g = 10.0 ** (noise_gain_db / 20.0)
strength = -t / np.log (1.0 - g)
flt = 1.0 - np.exp (-np.abs (fft)**2 / strength)
flt = nd.maximum_filter1d (flt, len (flt)/filter_length+1)
flt[-len (flt)/2:] = flt[len (flt)/2:0:-1]
impulse = np.fft.ifft (flt)
impulse[:filter_length] *= (np.cos (np.linspace (0, np.pi, filter_length))+1)/2
impulse = np.hstack ((
impulse[:filter_length],
np.zeros ((window_size - 2*filter_length), impulse.dtype),
impulse[filter_length:0:-1]))
return np.fft.fft (impulse)
def process_signal (input, noise_floor_db, noise_gain_db, filter_length):
window_size = 2 ** (filter_length-1).bit_length()
output = np.zeros_like (input)
output[:window_size/2] = input[:window_size/2] * (1 + np.cos (np.linspace (0.0, np.pi, window_size/2, endpoint = False))) / 2
work = np.zeros ((2*window_size,), input.dtype)
start = 0
while start + window_size < len (input):
work[:window_size] = input[start:start+window_size]
h = make_filter (work[:window_size], noise_floor_db, noise_gain_db, filter_length, 2*window_size)
work[:window_size] = input[start:start+window_size]
output[start:start+window_size] \
+= np.fft.ifft (np.fft.fft (work) * h)[:window_size].real \
* (1 - np.cos (np.linspace (0.0, 2*np.pi, window_size, endpoint = False)))/2
start += window_size/2
return output