5
$\begingroup$

I'm trying to code the Haar wavelet approximation of a function.
I'm a complete newbie.

Edit
Below is a pseudo-code of my algorithm for $f: [0,10] \to \mathbb R$, $f(x)=x$.
I picked it here and there on google.

1) Haar wave formula of order $N$ for $x \in [0, \infty)$ $$ \psi_{N,j}(x) = \frac{1}{\sqrt{N}} \cdot\begin{cases} 1 & j =0 \\ +2^{p/2} & a \le x - \left \lfloor \frac{x-a}{L}\right\rfloor L< b \\ -2^{p/2} & b \le x - \left \lfloor \frac{x-a}{L} \right\rfloor L < c\\ 0 & \text{otherwise} \end{cases} $$ where \begin{align} a &= \frac{q-1}{2^p}, b = \frac{q-0.5}{2^p}, c = \frac{q}{2^p}, L = c-a, \\ p &= \left\lfloor \log_2(j+1) \right\rfloor, q = j +1 - p. \end{align}

2) Haar series coefficients $$ \hat \psi_{N,j} = \int_0^{10} f(t) \psi_{N,j}(t) dt $$

3) Haar series $$ f(x) \approx \sum_{j=0}^N \hat \psi_{N,j} \cdot \psi_{N,j}(x) $$

4) For $N=4$, this gives the graph enter image description here

$\endgroup$
  • 1
    $\begingroup$ Are you referencing any specific paper or book? Also, it would help if you could write out your algorithm in very condensed pseudocode, since most people here will not be familiar with Python lambda functions. $\endgroup$ – Phonon Oct 22 '13 at 23:21
  • 1
    $\begingroup$ @Phonon I picked the algorithm here and there on google. The pseudocode idea is very good, I added it. Thank you for helping me clarify the question. $\endgroup$ – user5746 Oct 23 '13 at 0:06
  • $\begingroup$ @NicolasEssis-Breton Can you add your codes? $\endgroup$ – user1772257 Oct 24 '13 at 13:27
  • $\begingroup$ @user1772257 Jan on scicompSE shared his code, which helped me fix my code. You can find it here. I add it below for completeness. $\endgroup$ – user5746 Oct 25 '13 at 14:02
2
$\begingroup$

Thanks to Jan's code. I made my implementation work.
The code below compare: Haar vs Fourier vs Chebyshev.
enter image description here

from __future__ import division
from mpmath import *

# --------- Haar wavelet approximation of a function
# algorithm from : http://fourier.eng.hmc.edu/e161/lectures/wavelets/node5.html
# implementation only handle [0,1] for the moment: scaling and wavelet fcts need to be periodice

phi = lambda x : (0 <= x < 1) #scaling fct
psi = lambda x : (0 <= x < .5) - (.5 <= x < 1) #wavelet fct
phi_j_k = lambda x, j, k : 2**(j/2) * phi(2**j * x - k)
psi_j_k = lambda x, j, k : 2**(j/2) * psi(2**j * x - k)

def haar(f, interval, level):
    c0 = quadgl(  lambda t : f(t) * phi_j_k(t, 0, 0), interval  )

    coef = []
    for j in xrange(0, level):
        for k in xrange(0, 2**j):
                djk = quadgl(  lambda t: f(t) * psi_j_k(t, j, k), interval  )
                coef.append( (j, k, djk) )

    return c0, coef

def haarval(haar_coef, x):
    c0, coef = haar_coef
    s = c0 * phi_j_k(x, 0, 0)
    for j, k ,djk in coef:
            s += djk * psi_j_k(x, j, k)
    return s

# --------- to plot an Haar wave
interval = [0, 1]
plot([lambda x : phi_j_k(x,1,1)],interval)

# ---------- main
# below is code to compate : Haar vs Fourier vs Chebyshev

nb_coeff = 5
interval = [0, 1] # haar only handle [0,1] for the moment: scaling and wavelet fcts need to be periodice

fct = lambda x : x

haar_coef = haar(fct, interval, nb_coeff)
haar_series_apx = lambda x : haarval(haar_coef, x)

fourier_coef = fourier(fct, interval, nb_coeff)
fourier_series_apx = lambda x: fourierval(fourier_coef, interval, x)
chebyshev_coef = chebyfit(fct, interval, nb_coeff)
chebyshev_series_apx = lambda x : polyval(chebyshev_coef, x)


print 'fourier %d chebyshev %d haar %d' % ( len(fourier_coef[0]) + len(fourier_coef[1]),len(chebyshev_coef), 1 + len(haar_coef[1]))
print 'error:'
print 'fourier', quadgl(  lambda x : abs( fct(x) - fourier_series_apx(x) ), interval  )
print 'chebyshev', quadgl(  lambda x : abs( fct(x) - chebyshev_series_apx(x) ), interval  )
print 'haar', quadgl(  lambda x : abs( fct(x) - haar_series_apx(x) ), interval  )

plot([fct, fourier_series_apx, chebyshev_series_apx, haar_series_apx], interval) 
$\endgroup$

Your Answer

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