k = 10
ws = 2**k
tf = np.zeros(k)
data = np.zeros(n)
rem = snd.copy()
for i in range(k - 2, -1, -1):
    sz = 2 ** (i + 3)
    hsz = sz // 2
    freq = []
    rng = (0.5 + np.arange(sz)) / sz
    wdw = np.sin(np.pi * rng) ** 2
    rand = np.random.uniform(2 * np.pi)
    for j in range(n // hsz - 1):
        f = 2 * np.dot(
            rem[j * hsz : (j + 2) * hsz], np.exp(-2j * np.pi * rng)
        ) / sz
        data[j * hsz : (j + 2) * hsz] += np.real(
            f
            * wdw
            * np.exp(2j * np.pi * rng)
        )
    rem -= data
pp.figure(figsize=(10,6))
pp.plot(snd[:2000])
pp.plot(data[:2000])
wf.write("out.wav", sr, data.astype(np.float32))