# High frequencies disappear when applying discrete wavelet transform

Trying to decompose and reconstruct a signal using a to some extent self-made implementation of DWT for some reason fails. The result looks highpass filtered and/or shifted. I wanted to write the code in Python not using the PyWavelets package (except for getting the wavelet object).

import numpy as np
from scipy import interpolate, signal
import pywt
import matplotlib.pyplot as plt

# Test signal and wavelet
sig = np.zeros((10*512))
for i in range(256):
sig += i * np.cos(2*np.pi*i*np.linspace(0, 10, 10*512)+i**2*2*np.pi/((i+1)/1024))
w = pywt.DiscreteContinuousWavelet('db6')


I convolved the test signal with the decomposition filters from the wavelet object and used scipy's signal.decimate(...). As far as I understood, the decomposition filters do not lowpass filter the signal (or do they and is this where I went wrong?). However, a simple downsampling instead of decimation did not yield better results.

# Downsampling by decimation, decomposition
coefs_own = []
cA = np.copy(sig)
for i in range(1):
cA_last = cA
cA = np.convolve(cA_last, w.dec_lo, 'same')
cD = np.convolve(cA_last, w.dec_hi, 'same')
cA = signal.decimate(cA, 2)
cD = signal.decimate(cD, 2)
coefs_own.insert(0, cD)
coefs_own.insert(0, cA)


And this is my reconstruction using linear interpolation and a convolution with the recomposition filters provided by the wavelet object.

# Helper function
def upsample(data):
x = np.array(range(data.shape[0]))
x_double = np.linspace(x[0], x[-1], x.shape[0]*2)
f = interpolate.interp1d(x, data, kind='linear')
data_upsampled = f(x_double)
return data_upsampled

# Upsampling by interpolation, reconstruction
cA = coefs_own[0]
for i in range(1):
cA_up = upsample(cA)
cD_up = upsample(coefs_own[i+1])
cA = np.convolve(cA_up, w.rec_lo, 'same') + np.convolve(cD_up, w.rec_hi, 'same')
sig_rec_own = cA


For comparison, see the multilevel decomposition by PyWavelets, which does a perfect reconstruction.

# PyWavelets
coefs = pywt.wavedec(sig, w, mode='zero', level=1)
sig_rec_pywt = pywt.waverec(coefs, w, mode='zero')

# Plot both decomposition versions
fig, axes = plt.subplots(2, 1)
axes[0].plot(coefs_own[0])
axes[0].plot(coefs_own[1])
axes[0].set_title('Own')
axes[0].legend(['Approx.', 'Detail'])
axes[1].plot(coefs[0])
axes[1].plot(coefs[1])
axes[1].set_title('PyWavelets')
axes[1].legend(['Approx.', 'Detail'])

# Plot both reconstructions
fig, axes = plt.subplots(3, 1)
axes[0].plot(sig)
axes[0].set_title('Original')
axes[1].plot(sig_rec_own)
axes[2].plot(sig_rec_pywt)
axes[2].set_title('PyWavelets reconstruction')

plt.show()


Does anybody have an idea, why the decomposition and the reconstruction both obviously are erroneous? What did I miss in my implementation?

• Why the linear interpolation? Usually the wavelet filters are designed that you fill with zeros while upsampling. – LutzL Apr 29 at 15:59
• Thanks, you are right. It must be interpolated by inserting zeros. I would mark your answer as accepted, but I cannot find any button for that. As you see, I am newly registered here so I might have overlooked the accept-button? – Sergeant Salty May 2 at 13:42
• I put a longer version of the comment as answer, as of now you can only vote on answers, later you can give kudos to comments. – LutzL May 2 at 13:55

Upsampling with interpolation is equal to applying an additional filter $$[0.5,1,0.5]$$ which is a low-pass, smoothening filter. Esp. on the high-frequency or detail component this will lead to large errors.