3
$\begingroup$

Following on from a previous question posted by @xslittlegrass and answered by @Sector and others, Extracting information from the result of ContinuousWaveletTransform I'd like to know if it is possible to recover the amplitude of the original Sin components from the wavelet power spectra. In the original post all four Sin components had amplitude 1 but lets say they were 1,2,3 and 4, can we recover that from the wavelet coefficients? The coefficients are larger than 1 (if you execute the following statement after the lines in the original post you see the values of the coefficients).

ListDensityPlot[Reverse[Abs[cwd[All, "Values"]]] , PlotRange -> All, AspectRatio -> 1/GoldenRatio, PlotLegends -> Automatic, PlotRangePadding -> None]

(note, using WaveletScalogram[cwd, ColorFunction -> "RustTones", PlotLegends -> Automatic] shows that the scalogram color scale is normalized to 1)

I was hoping it would be possible to do this based on reading Torrence and Compo if I understand correctly (http://paos.colorado.edu/research/wavelets/ see page 68 of their paper from the link, in the Reconstruction section). Here it is shown that the original data can be reconstructed using a different wavelet from the one used to produce the CWT; in particular a delta function will do.

So here's the naive question expanded into a couple of parts: Can I call these (suitably scaled) wavelet coefficients the instantaneous power spectrum in Fourier space? If I use the Morlet wavelet the Fourier period is 1.03*the wavelet scale, so only a small adjustment would need to be made to convert it to frequency space from wavelet space.

Also: @endolith from the signal-processing stackexchange illustrated nicely how to preserve the sign of the original data in the scalogram using a bipolar colormap: https://dsp.stackexchange.com/questions/7911/reading-the-wavelet-transform-plot (here's the actual scalogram: https://i.stack.imgur.com/8w6YK.jpg)

I guess we can multiply the CWT coefficients by the sign of the original data and overlay two images with different colormaps, one for the - values and one for the + values. Is there a simpler way to do that with CWT in mathematica?

@Sector has kindly implemented the complex Morlet wavelet in Mathematica at the following link (the Mathematica implementation is for the real Morlet wavelet only). Continuous wavelet transform with complex Morlet function which is where I originally started to ask this question as a comment, but realized it was really a new question. It might be that I am completely missing the point and that there's no correspondence whatsoever between wavelet coefficients and instantaneous frequency component amplitudes, or that by definition "instantaneous" makes "frequency" in this context meaningless.

$\endgroup$
6
  • $\begingroup$ I am really sorry for the long delay, but like everyone else around these parts I have other issues I have to attend to. If you have any questions - do not hesitate :) $\endgroup$
    – Sektor
    Nov 18, 2014 at 22:09
  • $\begingroup$ No need to apologies! Lit is great to get a solution and explanation of the details, and even better when it is a surprise,by hank you! $\endgroup$
    – DrBubbles
    Nov 21, 2014 at 10:57
  • $\begingroup$ I've been thinking about that problem for the last couple of days and I have few ideas - I will try to work on them and update the answer :) Hope it sheds some more light on the issues. What exactly are you trying to do BTW ? Time-frequency analysis ? Extracting magnitude of signal components ? $\endgroup$
    – Sektor
    Nov 21, 2014 at 10:59
  • $\begingroup$ I'm trying to see if I can extract more information from some experimental data that look like interference fringes. Two frequencies, where they meet they mix. Normally I look at just the regions where the wavelength is well defined. $\endgroup$
    – DrBubbles
    Nov 21, 2014 at 11:05
  • $\begingroup$ I will see what I can do :) $\endgroup$
    – Sektor
    Nov 21, 2014 at 11:07

1 Answer 1

3
$\begingroup$

To answer your first question - yes. The wavelet scalogram or wavelet energy density function is defined as $$S(t_{0}, s) = \left |W(t_{0}, s) \right |^{2}$$ where $W(t_{0}, s)$ is the wavelet transform at time $t_0$ and scale $s$. If you are interested in the region around time $t_0$ then $S(t_{0}, s)$ would be the instantaneous energy distribution at that time. Bear in mind that when working with wavelets you are usually analysing non-stationary signals, so you will have to go through all possible $t_0$ values to get the bigger picture. If you are after the energy distribution by scales, defined as: $$\left \langle E(s) \right \rangle = \int \left | W(t_{0}, s) \right |^{2} dt_{0}$$ then there's a known connection between $\left \langle E(s) \right \rangle$ and $P(\omega)$, $P(\omega) = \left | \hat{f}(\omega) \right |^2$ a.k.a the power spectrum: $$\left \langle E(s) \right \rangle \sim s \int P(\omega)\left | \hat{\psi _{0}}(s \omega) \right |^2 d\omega$$ So, you are smoothing your power spectrum using the mother wavelet and scale it.

The second part of your question - yes, you can easily create such scalograms in Mathematica. Suppose we have the test signal $$\begin{Bmatrix} sin (2 \pi 8 t) & &0\leq t\leq 0.4 \\ sin (2 \pi 8 t)+\sin (2 \pi 32 t) & & 0.4<t<0.6 \\ sin (2 \pi 8 t) & & 0.6\leq t\leq 1 \end{Bmatrix}$$ which looks like this

Mathematica graphics

We create a new function:

cv = Blend[
Partition[
 Riffle[Range[-1, 1, .201], {
       RGBColor[0.403921568627451, 0., 0.12156862745098039`],
       RGBColor[0.6980392156862745, 0.09411764705882353, 0.16862745098039217`],
       RGBColor[0.8392156862745098, 0.3764705882352941, 0.30196078431372547`], 
       RGBColor[0.9568627450980391, 0.6470588235294118, 0.5098039215686274], 
       RGBColor[0.9921568627450981, 0.8588235294117647, 0.7803921568627451], 
       RGBColor[0, 0, 0, 0.9],
       RGBColor[0.8196078431372549, 0.8980392156862745, 0.9411764705882353], 
       RGBColor[0.5725490196078431, 0.7725490196078432, 0.8705882352941177], 
       RGBColor[0.2627450980392157, 0.5764705882352941, 0.7647058823529411], 
       RGBColor[0.12941176470588234`, 0.4, 0.6745098039215687], 
       RGBColor[0.0196078431372549, 0.18823529411764706`, 0.3803921568627451]}], 2], #] &;

and then just execute (cwd is your ContinuousWaveletData object)

WaveletScalogram[cwd, All, Re, ColorFunction -> cv, ImageSize -> 500]

and you get

Mathematica graphics

It's that easy (and ugly! ^^")

$\endgroup$

Your Answer

By clicking “Post Your Answer”, you agree to our terms of service and acknowledge you have read our privacy policy.

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