2
$\begingroup$

I'm trying to use Morlet wavelet to reduce noise in my signal. I found formula for Morlet wavelet

$$\Psi(x) = \frac{1}{\sqrt{\pi\cdot \textrm{bandwidth}}} \cdot \exp\left(i \cdot 2\pi \cdot \textrm{centerFreq} \cdot x\right)\cdot\exp\left(\frac{x^2}{\textrm{bandwidth}}\right)$$

(sb told me I shoud use formula for complex Morlet wavelet)

I also got CWT formula

$$C = \int f(t) \cdot \frac{1}{\sqrt{\textrm{scale}}} \cdot \Psi\left(\frac{t-\textrm{shift}}{\textrm{scale}}\right)dt$$

I convert both formulas to java code

public double morletRe(double x)
{
return (1/Math.sqrt(Math.PI*fb))*Math.cos(2*Math.PI*fc*x)* Math.exp(x*x/fb); 
}

public double[] cwt(double[] data,double scale, double position)
{
double[] newData = new double[data.length];
    for (int t = 0; t<data.length; t++ )
    {
    newData[t] = data[t]*(morletRe(t-position)/scale)[0])/Math.sqrt(scale);
    }
return newData;
}

But in the results I get a strange chart which is completely not related to my data, is anyone can point me what I'm doing wrong? regards.

$\endgroup$
  • $\begingroup$ Why aren't you using a library that does the wavelet for you? $\endgroup$ – endolith Sep 3 '12 at 13:05
  • $\begingroup$ can you recommend me any? all what i found was jwave, and it doesn't have molet wavelet $\endgroup$ – user902383 Sep 3 '12 at 13:09
  • $\begingroup$ Where did you find that formula? exp(x^2 / bandwidth) grows very fast. Maybe it's exp(-x^2/bw)? Also, where are you calculating the integral? Does the function that calls cwt calculate a sum over the array that's returned? $\endgroup$ – Niki Estner Sep 3 '12 at 13:34
  • 1
    $\begingroup$ I'm not an expert. I think the general idea is: 1. Transform the signal to a wavelet representation. 2. Threshold the wavelet representation to reduce noise. 3. Transform the wavelet representation back to a signal. Since the Morlet wavelet is overcomplete, this means solving an overdetermined equation system (either using a direct method like QR decomposition if the signal is short enough, or an iterative method like conjugate gradient for larger data sizes). But honestly: I doubt you'll be able to piece that together using formulas you find on mathworks. I'm sure I couldn't. $\endgroup$ – Niki Estner Sep 4 '12 at 8:49
  • 1
    $\begingroup$ @endolith: Yes, overcompleteness means the you can remove some some of the transform output and still have a perfect reconstruction of the input data. And I think Gabor transform and STFT are overcomplete, too. (Gabor frames are actually an example for overcompleteness on wikipedia.) $\endgroup$ – Niki Estner Sep 5 '12 at 21:57
2
$\begingroup$

I'm not sure of the reliability of your source... I worked with wavelets on my master's thesis and the basic wavelet transform is a convolution of the wavelet an your signal. In most practical cases, it's best to go through Fast Fourier Transform (FFT), as it's mathematically equivalent.

With WT being the wavelet transform, and FFT() being the Fast Fourier Transform of what's inside the brackets, iFFT() is the inverse-FFT. You would compute the WT as such:

FFT(WT) = FFT(signal) * FFT(wavelet)
WT = iFFT(WT)

Now, for the wavelet formula itself, I believe there has to be a specific scaling factor to get a reliable result. Here's what we used:

$$\Psi = exp\Big(-\frac{ik_0x}{\lambda} - \frac{x^2}{2\lambda^2}\Big)$$

Where $\lambda$ is the wavelength (or scale) and $k_0 = 5.336$ (source: J.F. Kirby, Computers & Geosciences 31 (2005) 846–864) I did not include the normalization factor because you just have to divide your wavelet by it's highest value, and subtract its average, and you get the same result with less headhache. Off course, you can interchange $1/\lambda$ by $2\pi*frequency$ in the formula.

In order to qualify for the "Morlet" name, I think you have to use the same wavelength in both places in the formula. Else, you over-sample or under-sample the frequency part of your signal. $k_0$ is the best ratio between the two, as explained by Kirby in his article.

$\endgroup$

Your Answer

By clicking “Post Your Answer”, you agree to our terms of service, privacy policy and cookie policy

Not the answer you're looking for? Browse other questions tagged or ask your own question.