import numpy as np
import numpy.linalg as la
import matplotlib.pyplot as pt
wave = np.load("cs357.npy")
pt.plot(wave)
from simple_audio import play
play(wave)
Now build the Fourier Vandermonde matrix. This is also called the Discrete Fourier Transform (DFT) matrix.
n = 255
assert n % 2 == 1
cos_k = np.arange(0, n//2 + 1, dtype=np.float64)
sin_k = np.arange(1, n//2 + 1, dtype=np.float64)
x = np.linspace(0, 2*np.pi, n, endpoint=False)
DFT = np.zeros((n,n))
DFT[:, ::2] = np.cos(cos_k*x[:, np.newaxis])
DFT[:, 1::2] = np.sin(sin_k*x[:, np.newaxis])
IDFT = la.inv(DFT)
And transform the signal:
n_chunks = len(wave)//n
chunked = wave[:n_chunks*n].reshape(-1, n)
coeffs = np.dot(IDFT, chunked.T).T
pt.figure(figsize=(10,8))
pt.imshow(np.log10(1e-5+np.abs(coeffs)).T)
pt.colorbar()
modified = coeffs.copy()
#modified[:,5:-5] *= 1e-2
#modified *= 1e-2
#modified[:,5:-5] /= 1e-2
nshift = 3
modified.fill(0)
modified[:, nshift:] = coeffs[:, 0:-nshift]
pt.figure(figsize=(8,5))
pt.imshow(np.log10(1e-5+abs(modified)).T)
pt.colorbar()
result = np.dot(DFT, modified.T).T
result = result.reshape(-1)
new_wave = np.zeros_like(wave)
new_wave[:len(result)] = result
new_wave = new_wave/np.max(np.abs(new_wave))
play(new_wave)
play(wave)
pt.plot(new_wave)