4
$\begingroup$

I have a problem in understanding how to create my own wavelets in Mathematica's wavelet framework.

There is a tutorial.

I have read the tutorial, but I don’t understand what exactly is needed to define my own wavelet type. I tried to implement a very simple wavelet (Haar wavelet) in the framework discussed in the tutorial to get a basic understanding, but it did not work out.

Does anyone has experience with this framework or directly sees how or what I need to do? I want to use this framework because it is the only way of using internal Mathematica functions with your own wavelets.

The tutorial says the following:

FranklinWavelet[lim_Integer]["WaveletQ"] := True
FranklinWavelet[lim_]["OrthogonalQ"] := True

a[k_] := Sqrt[3] HypergeometricPFQRegularized[{1/2, 1/2, 1}, {1 - k, 1 + k}, -2]

 Z[k_, n_] := 1/24 a[2 k - n - 2] + 1/4 a[2 k - n - 1] + 5/12 a[2 k - n] + 1/4 a[2 k - n + 1] + 1/24 a[2 k - n + 2]

FranklinWavelet[lim_Integer]["PrimalLowpass", prec_ : MachinePrecision] := 
 ParallelTable[{n, N[Total[Table[a[k] Z[k, n], {k, -15, 15}]], prec]}, {n,-lim, lim}]

wave = FranklinWavelet[14];

WaveletFilterCoefficients[wave]

{{-14, -9.77137*10^-6}, {-13, 0.0000212211}, {-12, 0.0000391428}, {-11, -0.0000867523}, {-10, -0.000158601}, {-9, 0.000361781}, {-8, 0.000652922}, {-7, -0.00155701}, {-6, -0.00274588}, {-5, 0.00706442}, {-4, 0.0120003}, {-3, -0.0367309}, {-2, -0.0488618}, {-1, 0.280931}, {0, 0.578163}, {1, 0.280931}, {2, -0.0488618}, {3, -0.0367309}, {4, 0.0120003}, {5, 0.00706442}, {6, -0.00274588}, {7, -0.00155701}, {8,0.000652922}, {9, 0.000361781}, {10, -0.000158601}, {11, -0.0000867523}, {12, 0.0000391428}, {13, 0.0000212211}, {14, -9.77137*10^-6}}

My code until now:

f[x_] := Piecewise[{{1, 0 <= x <= 1}}]

g[x_] := f[2 x] - f[2 x - 1];

MyWavelet[lim_Integer]["WaveletQ"] := True
MyWavelet[lim_]["OrthogonalQ"] := True

MyWavelet[lim_Integer]["PrimalLowpass", prec_ : MachinePrecision] := g[lim]

The Haar Wavelet in defined by it´s scaling function f(x) and has the wavelet function g(x). But obviously mathematica is waiting for some values that define the wavelet but how?

Improve this question
$\endgroup$
7
  • 1
    $\begingroup$ "I tried to implement..." - so why not include the code you wrote, which might help us diagnose the problem? $\endgroup$
    – J. M. can't deal with it ♦
    Mar 3 '16 at 11:14
  • $\begingroup$ To come up with an answer, most of us (and certainly me), would need to experiment with your code. But you don't give us any code. You are unlikely to get an answer unless you edit your question to include all the code needed to reproduce your problem. $\endgroup$
    – m_goldberg
    Mar 3 '16 at 11:20
  • $\begingroup$ I added some code. I hope u can open the topic again. ty $\endgroup$
    – Dan
    Mar 4 '16 at 12:59
  • $\begingroup$ It might be useful, mathematica.stackexchange.com/questions/33549/… $\endgroup$
    – s.s.o
    Mar 4 '16 at 14:34
  • 1
    $\begingroup$ @Dan, In every example on that page, the wavelet primallowpass returns a list, but yours just returns a number. $\endgroup$
    – Jason B.
    Mar 4 '16 at 14:59
3
$\begingroup$

You almost have it working, just a couple minor fixes and it seems to work. Fist, in order to use DiscreteWaveletTransform I had to set it to be an orthogonal wavelet. Then I corrected your definition of the PrimalLowpass to include a second optional argument of the precision.

f[x_] := Piecewise[{{1, 0 <= x <= 1}}]
g[x_] := f[2 x] - f[2 x - 1]

MyWavelet[]["WaveletQ"] := True
MyWavelet[]["OrthogonalQ"] := True
MyWavelet[]["BiorthogonalQ"] := False
MyWavelet[]["WaveletFunction"] := g[#1] &
MyWavelet[]["FourierFactor"] := 1/Sqrt[2 Pi]
MyWavelet[]["FourierTransform"] := Sinc[#1/2] &
MyWavelet[]["PrimalLowpass", 
  prec_: MachinePrecision] := {{0, 0.5}, {1, 0.5}}

Now you can use DiscreteWaveletTransform,

WaveletScalogram[
 DiscreteWaveletTransform[Table[Sin[x^2], {x, -6, 6, 0.01}], 
  MyWavelet[]]]

enter image description here

As suggested by Sektor, you can test the validity of the transform by comparing the back-transformed data to the original data. This transform is good for lists of 63 data points or less when using ContinuousWaveletTransform,

ListPlot[
 Table[
  data = RandomReal[5, n];
  {n, Max[
    data - InverseContinuousWaveletTransform@
      ContinuousWaveletTransform[data, MyWavelet[]]]}, {n, 10, 100}]
 ]

enter image description here

Why it fails after that? I haven't the foggiest, but there is no such problem with the discrete transform.

ListPlot[
 Table[
  data = RandomReal[5, n];
  {n, Max[
    data - InverseWaveletTransform@
      DiscreteWaveletTransform[data, MyWavelet[]]]}, {n, 10, 1000}]
 ]

enter image description here

Improve this answer
$\endgroup$
5
  • $\begingroup$ For completeness: alternative representations for f include UnitBox[x - 1/2] and BSplineBasis[0, 0, x]. $\endgroup$
    – J. M. can't deal with it ♦
    Mar 10 '16 at 15:31
  • $\begingroup$ And as always -- a good 'sign' that you are on the right track -- invert the DWT and compare the error between the input and the output. $\endgroup$
    – Sektor
    Mar 10 '16 at 16:05
  • $\begingroup$ perfect thank you :) $\endgroup$
    – Dan
    Mar 11 '16 at 16:51
  • $\begingroup$ So how to define a new continuous wavelet? $\endgroup$
    – Wang Yun
    Nov 4 '20 at 8:11
  • $\begingroup$ Dear all, how to define the new wavelet in the following post: mathematica.stackexchange.com/questions/255195/… $\endgroup$
    – 1_student
    Sep 7 '21 at 18:27
2
$\begingroup$

Until now I develop a partly solution of my own problem. With this solution someone can compute continuous wavelet transforms and inverse wavelet transforms as well as other examples below. But until now I dont fixed the problem with the low pass filter. As you can see I added the option for the filter coefficients but this is not working out. I think if this problem gets fixed someone can use every option of the framework with its own wavelet.

Some help to at least fix that would be realy nice.

f[x_] := Piecewise[{{1, 0 <= x <= 1}}]
g[x_] := f[2 x] - f[2 x - 1]

MyWavelet[]["WaveletQ"] := True
MyWavelet[]["OrthogonalQ"] := False
MyWavelet[]["BiorthogonalQ"] := False
MyWavelet[]["WaveletFunction"] := g[#1] &
MyWavelet[]["FourierFactor"] := 1/Sqrt[2 Pi]
MyWavelet[]["FourierTransform"] := Sinc[#1/2] &
MyWavelet[]["PrimalLowpass"] := {{0, 0.5}, {1, 0.5}}

With this options set someone can perform the following commands:

Plot[WaveletPsi[MyWavelet[], x], {x, -1, 2}, Exclusions -> None]

enter image description here

data = Table[Random[], {i, 7}]
CWTMyWavelet = ContinuousWaveletTransform[data, MyWavelet[]]
{0.607242, 0.866118, 0.898763, 0.802639, 0.445842,0.671422, 0.913055}

 InverseContinuousWaveletTransform[CWTMyWavelet]
{0.607242, 0.866118, 0.898763, 0.802639, 0.445842,0.671422, 0.913055}

WaveletScalogram[CWTMyWavelet, Ticks -> Full]

enter image description here

DiscreteWaveletTransform[data, MyWavelet]
DiscreteWaveletTransform::bdwave: The wavelet specification MyWavelet cannot be used to perform DiscreteWaveletTransform. >>
Improve this answer
$\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.