I'm all new to wavelet analysis. I'm trying to get a working understanding of the continuous wavelet transform and its inverse. By "working understanding", I really mean "getting some code to execute for 1D signals and have an intuitive understanding of what the result means".
I'm using the ssqueezepy
library in python, by https://dsp.stackexchange.com/users/50076/overlordgolddragon
This library provide cwt
and icwt
routines for the continuous wavelet transform and its inverse ... but the inverse process does not seem to get the original signal back.
More precisely, using the following code, I expected to see the same signal several times:
import numpy as np
import itertools
import ssqueezepy
import matplotlib.pyplot as plt
wavelet = 'gmw'
nv = 16
n = 100
t = np.linspace(0.,1.,n)
x = np.sin(10*t*t)
x_mean = x.mean()
l1_norm_choices = [True,False]
one_int_choices = [True,False]
plt.plot(x,color='b')
for (l1_norm,one_int) in itertools.product(l1_norm_choices,one_int_choices):
cwtmatr, freqs = ssqueezepy.cwt(x,wavelet=wavelet,nv=nv,l1_norm=l1_norm)
x_inv = ssqueezepy.icwt(cwtmatr, wavelet=wavelet, nv=nv, one_int=one_int, l1_norm=l1_norm, x_mean=x_mean)
plt.plot(x_inv,color='r')
plt.show()
Instead, this is my result (original signal in blue, tentative reconstructions in red) :
Can someone help me
- understand why I don't retrieve the original signal with my code
- tweak the code so that I do retrieve the original signal
UPDATE: Here is a more thorough piece of code that searches for combinations of arguments that makes icwt the inverse of cwt. Found nothing.
import numpy as np
import itertools
import ssqueezepy
import matplotlib.pyplot as plt
import os
output_folder = os.path.dirname(__file__)
basename = os.path.splitext(os.path.basename(__file__))[0]
Do_plot = True
# Do_plot = False
plot_ext = [
'.png',
# '.eps',
# '.pdf',
]
eps = 1e-7
nv = 256 # Does this matter ?
n = 256
t = np.linspace(0.,1.,n)
x = np.sin(10*t*t+1) # or whatever signal
x_mean = x.mean()
wavelet_choices = list(ssqueezepy.wavs())
l1_norm_choices = [True,False]
one_int_choices = [True,False]
scales_choices = ['log','linear','log-piecewise']
padtype_choices = ['reflect', 'symmetric', 'replicate', 'wrap', 'zero' ]
# rpadded_choices = [True,False] # True fails an assertion check
rpadded_choices = [False]
i_plot = -1
for (l1_norm,one_int,scales,wavelet,padtype,rpadded) in itertools.product(l1_norm_choices,one_int_choices,scales_choices,wavelet_choices,padtype_choices,rpadded_choices):
cwtmatr, freqs = ssqueezepy.cwt(x,wavelet=wavelet,nv=nv,l1_norm=l1_norm,scales=scales,padtype=padtype,rpadded=rpadded)
x_inv = ssqueezepy.icwt(cwtmatr, wavelet=wavelet, nv=nv, one_int=one_int, l1_norm=l1_norm, x_mean=x_mean,scales=scales,padtype=padtype,rpadded=rpadded)
all_close = np.linalg.norm(x-x_inv) < eps * np.linalg.norm(x)
if (all_close):
print('Yay, found one ')
print(f'l1_norm={l1_norm}, one_int={one_int}, scales={scales}, wavelet={wavelet}, padtype={padtype}, rpadded={rpadded}')
if Do_plot:
plt.plot(x,color='b')
plt.plot(x_inv,color='r')
plt.title(f'l1_norm={l1_norm}, one_int={one_int}, scales={scales}, wavelet={wavelet}, padtype={padtype}, rpadded={rpadded}')
i_plot += 1
for ext in plot_ext:
plt.savefig(os.path.join(output_folder,f'{basename}_{i_plot}{ext}'))
plt.close()
```
help(icwt)
under References: "One integral inverse CWT. dsp.stackexchange.com/a/71274/50076", and "Inverse CWT derivation. dsp.stackexchange.com/a/71148/50076" $\endgroup$padtype
, I'll note it, thanks. Note it doesn't affectone_int=True
. $\endgroup$