After a badly formulated question, I decided to make a new post searching for help.
The basic problem is the follows: I have a wavelet function $\psi(t)$ (real or complex) and would like to compute (a) the normalization factor to guarantee that the wavelet function has unitary energy and (b) the admissibility condition $C_g$ that is very important to the inverse of the continuous wavelet transform.
I wrote two functions for this and they work very well with Mexican hat wavelet (a real function) according to my tests based on the book "The Illustrated Wavelet Transform Handbook". For, on page 10 (Equation 2.10) the author found a normalization constant to asseverates unitary energy of this wavelet that is the same that I computed using my code below:
Note: I imported more that I used because my code has more functions
from sympy import exp, fourier_transform, sqrt, integrate, oo, pi, sin, cos
from sympy import Symbol
omega = Symbol('omega', real=True)
t = Symbol('t', real=True)
def normalization_constant_L2(wave):
'''
# This function computes the normalization of a wavelet
# considering the L2 normalization for that ||\psi(t)|| = 1.
# This function uses the Sympy Module.
# Important read the Sympy documentation if you are a student.
# Variables:
# t is the time variable (independent variable) of wavelets
# $omega" is the frequency angular computed by the Fourier Transform.
# INPUT: Wavelet Function (symbolic function) without normalization
# OUTPUTS:
# (a) Only normalization factor
# (b) Wavelet Function of Input with normalization factor
# REFERENCE::
# [BOOK] ADDISON, Paul S.
# The illustrated wavelet transform handbook: introductory theory and /
# applications in science, engineering, medicine, and finance.
# CRC press, 2017. (2a. edition). pages 07 - 12.
'''
# Main Code
# Fourier Transform of Wavelet Function
fourier_wave = fourier_transform(wave, t, omega)
# Computing of the Spectrum Fourier
# Remember the isometry of time and frequency domains (Plancherel)
spectrum = abs(fourier_wave)**2
# Energy as improper integral of Fourier spectrum on the real line
energy_wave = integrate(spectrum, (omega, -oo, oo))
# Normalization : (square root of the energy) **(-1)
normalization = sqrt((energy_wave)**(-1))
return normalization
I defined the Mexican Hat as made in the literature to make tests. Your definition is (this function is the second derivative of a Gaussian Function)
wave = (1-t**2)*exp(-t**2/2)
And my code above gives a normalization constant equal to $k = 2*\sqrt(3)/(3*pi**(1/4)) = 2/(\sqrt(3)*pi**0.25$ that is the correct value according to the reference above where we can read the same result to normalization constant (Eq. 10, p.12). (Question: Is possible to make Sympy simplify for me? Here was very simple but in hard cases the problem is impracticable).
In the sequence, I trying to compute the value of admissibility condition $C_g$ that is fundamental to calculate the inverse of the continuous wavelet transform. Here, I wrote this code:
By definition we have $$ C_g=\int_{0}^{\infty}\dfrac{|\widehat{\psi(f)}|^2}{f} df $$ where
- $f$ is the frequency variable and
- $\widehat{\psi(f)}$ is the Fourier Transform of the $\psi(t)$ function.
def admissibility_condition_L2(wave, normalization_wave):
'''
# This function computes the admissibility condition of the wavelet
# Zero values indicate that the inverse of the CWT is not performed.
# INPUT: (a) Wavelet Function (symbolic function) without normalization
# (b) Normalization parameter computed with function \
"normalization_constant_L2"
# OUTPUTS:
# (a) Admissibility Condition Cg with normalization constant
# (b) Admissibility Condition Cg without normalization constant
# REFERENCE::
# [BOOK] ADDISON, Paul S.
# The illustrated wavelet transform handbook: introductory theory and /
# applications in science, engineering, medicine, and finance.
# CRC press, 2017. (2a. edition). Equation 2.4 pag. 10.
'''
# Definition of Cg as the improper integral of the magnitude of square \
# wavelef coefficients split to frequencies f in the interval [0, oo]
admiss1 = simplify(integrate(abs(normalization_wave)**2, (omega, 0, oo)))
admiss2 = simplify(integrate(abs(wave)**2, (omega, 0, oo)))
return admiss1, admiss2
OBS: I defined two returns because I would like to study the difference that happens when I use a normalization constant or not.
Unfortunately, to Mexican hat, the result correct is 3.451 but my inputs are infinities. The integrals computed via code diverges.
Same suggestions?
I appreciate each contribution.
Luciano.
