Created on September 7, 2000 using Maple V.

Modified on December 8, 2001 using Maple 7.

This file contains programs used to plot the figures in the book:

Wavelets and Other Orthogonal Systems

Series: Studies in Advanced Mathematics  Volume: 36

(Second edition, ISBN 1584882271 )

Gilbert G.Walter and Xiaoping Shen

December 20, 2000

CRC Press, Boca Raton, FL, USA

Transparancy for figures in the book

Chapter 1. Orthogonal Systems

Figure 1.1 The scaling and mother wavelet for the Haar system

The scaling function

> phi:=t->piecewise(0<=t and t<1,1, t<0 or t>=1,0);

> plot(phi(t),t=-0.5..1.5,-1..1.5, tickmarks=[2,3]);

The wavelet function:

> psi:=t->phi(2*t)-phi(2*t-1);

> plot(psi(t),t=-0.5..1.5,-1..1.5,tickmarks=[2,3]);

Figure 1.2 The scaling function for the Shannon system

The scaling function

> shannon:=t->sin(Pi*t)/(Pi*t);

> plot(shannon(t),t=-4..4,-0.3..1.3,tickmarks=[3,2]);

The Fourier transform of the scaling function

> phihat:=w->piecewise(-Pi<=w and w<Pi,1, w<-Pi or w>=Pi,0);

> plot(phihat(w),w=-4..4,-0.5..1.3, tickmarks=[3,2]);

Figure 1.3 The wavelet for the Shannon system

> g:=t->-2*(sin(2*Pi*t)+cos(Pi*t))/(Pi*(2*t+1));

>

> plot(g(t),t=-4..4, -0.9..1.1,tickmarks=[3,2]);

Chapter 2. Tempered Distributions

Figure 2.1 A test function in the space S (The Hermite function h4(x))

> h0:=t->exp(-t^2/2);

> h2:=t->(2*t^2-1)*h0(t);

> evalf(Int(h2(t),t=-10..10));

> plot(2*h2(8*t),t=-1..1,-2.5..3, tickmarks=[2,3]);

Figure 2.2 Some approximations to the delta function in S'

> n:=4;delta:=t->piecewise(1/n<=t or t<-1/n,0, -1/n<=t and t<=0,n^2*t+n,0<t and t<1/n, -n^2*t+n);

> plot(delta(t),t=-3..3,-0.1..4.5,tickmarks=[3,3]);

> delta0:=t->piecewise(0<=t and t<1/n,n, t<0 or t>=1/n, 0);plot(delta0(t),t=-2..2,tickmarks=[3,3]);

> delta2:=t->sin(4*t)/(t);

> plot(delta2(t),t=-2..2,-1..4.5,tickmarks=[3,3]);

Chapter 3. Introduction of Wavelets

Figure 3.1 Typical functions in the subspaces V 0 of the multiresolution analysis for the Haar scaling function

> phi:=t->piecewise(0<=t and t<1,1, t<0 or t>=1, 0);

> f0:=t->phi(t)+.5*phi(t-1)-.7*phi(t-2)+.8*phi(t-3);

> plot (f0(t),t=-1..5,-0.8..1.2,tickmarks=[3,3]);

Figure 3.2 Typical functions in the subspaces V1 of the multiresolution analysis for the Haar scaling function

> f1:=t->f0(2*t)-.7*f0(2*t-4);

> plot(f1(t),t=-1..5,-0.8..1.2);

Figure 3.3 Typical functions in the subspaces V2 of the multiresolution analysis for the Haar scaling function

> f2:=t->f1(2*t)-.7*f1(2*t-4);

> plot(f2(t),t=-1..5,f2=-0.8..1.2);

Figure 3.4 The scaling function of the Franklin wavelet araising from the hat function

> pF:=t->piecewise(1<=t or t<-1,0, -1<=t and t<=0,t+1,0<t and t<1, -t+1);

> plot(pF(t),t=-3..3,-0.1..1.1,tickmarks=[3,2]);

> s1:=w->(1+w)^(-1/2);

> taylor(s1(w),w=0,4);

> phiF:=t->sqrt((3)/2)*((67/64)*pF(t)-(1/8)*pF(t-1)-(1/8)*pF(t+1)+(3/128)*pF(t-2)+(3/128)*pF(t+2));

> plot(phiF(t),t=-4..4,-0.5..1.5,tickmarks=[3,3]);

Figure 3.5 The Daubechies scaling function and mother wavelet of example 3

(Using Matlab Wavelet toolbox)

Figure 3.6 The Da ubechies mother wavelet of example 3

(Using Matlab Wavelet toolbox)

Figure 3.7 The scaling function and mother wavelet for a Meyer wavelet in frequency domain

> beta:=1/4;

> fphi11:= w->piecewise(0<w and w<Pi*(1-beta),1, Pi*(1-beta)<w and w<Pi*(1+beta),cos(w/(4*beta)-Pi*(1-beta)), w>Pi*(1+beta),0);

> fphi1:= w->piecewise(0<w, fphi11(w), w<0, fphi11(-w));

> plot(fphi1(w),w=-8..8,0..1.1,tickmarks=[3,2]);

> abfpsi1:=w->(fphi1(w+2*Pi)+fphi1(w-2*Pi))*fphi1(w/2);

> plot(abfpsi1(w),w=-10..10,0..1.1,tickmarks=[2,2]);

Figure 3.8 The Scaling function of Figure 3.7 in the time domain

Type one raised cosine wavelet (this function is the square root of raised cosine function used

by electrical enginners)

> beta:=1/4;

> phi1:=t->(sin(Pi*(1-beta)*t)+4*beta*t*cos(Pi*(1+beta)*t))/((Pi*t)*(1-(4*beta*t)^2));

> limit(phi1(t),t=0):limit(phi1(t),t=1):limit(phi1(t),t=-1):

>

> plot(phi1(t),t=-5..5,-0.24..1.2,tickmarks=[3,2]);

Figure 3.9 The mother wavelet of Figure 3.7 in the time domain

*Notice that here "psi1(t)" is, in fact, psi1(t+1/2). Same as in the definition of psi2(t).

> psi1:=t->(sin(Pi*(1+beta)*t)-4*beta*t*cos(Pi*(1-beta)*t))/((Pi*t)*((4*beta*t)^2-1))-(sin(2*Pi*(1-beta)*t)+8*beta*t*cos(2*Pi*(1+beta)*t))/((Pi*t)*((8*beta*t)^2-1));

> plot(psi1(t),t=-5..5,-0.9..1.3);

Figure 3.10 The Reproducing kernel q(x,t) for V 0 in the case of Haar wavelets

> phi:=t->piecewise(0<=t and t<1,1, t<0 or t>=1,0);

> q(s,t):=phi(s-trunc(t));

>

> plot3d(q(s,t),s=0..6,t=0..6);

Figure 3.11 The Daubechies scaling function and mother wavelet (N=4)

(Using MatLab Wavelet toolbox)

Figure 3.13 The system functions of some continuous filters

Low pass

> lowpass:=t->piecewise(t>-1 and t<1,10,0);

> lp:=t->int(lowpass(t-s)*exp(-4*s^2),s=-10..10);

> plot(lp(t),t=-3..3);

High pass

> highpass:=t->8.84-lp(t);

> plot(highpass(t),t=-5..5);

Band pass

> bandpass:=t->piecewise(abs(t)>2 and abs(t)<3,1,0);

> bp1:=t->lp(t-3)+lp(t+3);

> plot(bp1(t),t=-5..5);

Figure 3.14 The system function of discrete lowpass (halfband) filter

> plot(lp(t)/4.4,t=0..3, tickmarks=[3,3]);

Figure 3.15 The decomposition algorithm

Figure 3.16 The reconstruction algorithm

Chapter 4. Convergence and Summability of Fourier Series

Figure 4.1 Dirichlet kernel of Fourier series (n=6)

> Diri:=(n,x)->sin((n+1/2)*x)/(2*Pi*sin(x/2));

> plot(Diri(6,x),x=-5..5,-0.5..2.3,color=blue, tickmarks=[3,3]);

Showing periodicity

> plot(Diri(6,x),x=-10..10,-0.6..2.2,color=blue, tickmarks=[3,3]);

Figure 4.2 Fejer kernel of Fourier series (n=6)

> Fejer:=(n,x)->sin(n*x/2)^2/((2*Pi*n)*(sin(x/2))^2);

> plot(Fejer(6,x),x=-Pi..Pi,-0.05..1.12,ytickmarks=2);

Show the periodicity

> plot(Fejer(6,x),x=-3*Pi..3*Pi,-0.05..1.12,ytickmarks=2);

> 

Figure 4.3 The saw tooth function

> saw:=x->(Pi/2)*piecewise(x>0 and x<Pi,Pi-x,x<0 and x>-Pi,-Pi-x,0);

>

> plot(saw(x),x=-3.5..3.5,-6.5..6.5,tickmarks=[3,3],discont = true,color=blue);

Figure 4.4 Gibbs phenomenon for Fourier series; approximation to the saw tooth function using Dirichlet kernel.

> saw:=x->(Pi/2)*piecewise(x>0 and x<Pi,Pi-x,x<0 and x>-Pi,-Pi-x,0);

> plot(saw(x),x=-3.5..3.5,s=-6.5..6.5,discont = true,color=blue);

> b:=n->(2/Pi)*int(saw(x)*sin(n*x),x=0..Pi);

> s:=(n,x)->piecewise(abs(x)<Pi,sum(b(k)*sin(k*x), k=1..n),0);

> plot(s(80,x),x=-3.5..3.5,y=-6.5..6.5,discont = true,color=blue);

The approximation to the saw tooth function using Fejer kernel

The (truncated) Cesaro means of the Fourier series are defined by using Fejer kernel as in the following :

> sigma80:=x->(1/80)*sum(s(m,x),m=1..79);

> plot(sigma80(x),x=-3.5..3.5,-6.5..6.5,discont = true,color=blue);

Chapter 5. Wavelets and Tempered Distributions

Figure 5.1 A mother wavelet with point support. The vertical bars represent delta functions.

> with(plottools):
line1:= curve([[0,0], [0,0.5]],linestyle=1,color=red, thickness=3,x=-1..1.5):
line2:= curve([[1/2,0],[1/2,1]], linestyle=1, color=red,thickness=2,x=-1..1.5):line3:= curve([[1,0],[1,0.5]],linestyle=1,color=red, thickness=2,x=-1..1.5):

> point1:=curve([[1/2,1],[1/2,1.01]],color=blue, thickness=3,style=point):point2:=curve([[0,1/2],[0,0.51]],color=blue, thickness=3,style=point):point3:=curve([[1,1/2],[1,0.51]],color=blue, thickness=3,style=point):

> plots[display]({line1,line2,line3, point1,point2,point3},ytickmarks=2);

Figure 5.2 Continuous function and its impulse train

> f:=x->sin(x)+x*cos(x)+2;

> l:= plot(f(x), x=-4.5..4.7,color=blue):

> with(plottools):
l1:= curve([[-4,0], [-4,f(-4)]],linestyle=1,color=red, thickness=2):
l2:= curve([[-3,0], [-3,f(-3)]],linestyle=1,color=red, thickness=2):
l3:= curve([[-2,0], [-2,f(-2)]],linestyle=1,color=red, thickness=2):
l4:= curve([[-1,0], [-1,f(-1)]],linestyle=1,color=red, thickness=2):l5:= curve([[0,0], [0,f(0)]],linestyle=1,color=red, thickness=2):l6:= curve([[1,0], [1,f(1)]],linestyle=1,color=red,thickness=2):l7:= curve([[2,0], [2,f(2)]],linestyle=1,color=red, thickness=2):l8:= curve([[3,0], [3,f(3)]],linestyle=1,color=red,thickness=2):
l9:= curve([[4,0], [4,f(4)]],linestyle=1,color=red, thickness=2):

> p1:= curve([[-4,f(-4)], [-4,f(-4)+0.1]],linestyle=1,color=black,style=point, thickness=3):
p2:= curve([[-3,f(-3)], [-3,f(-3)+0.1]],linestyle=1,color=black,style=point, thickness=3):p3:= curve([[-2,f(-2)], [-2,f(-2)+0.1]],linestyle=1,color=black,style=point,thickness=3):
p4:= curve([[-1,f(-1)], [-1,f(-1)+0.1]],linestyle=1,color=black,style=point, thickness=3):p5:= curve([[0,f(0)], [0,f(0)+0.1]],linestyle=1,color=black,style=point, thickness=3):
p6:= curve([[1,f(1)], [1,f(1)+0.1]],linestyle=1,color=black, style=point,thickness=3):p7:= curve([[2,f(2)], [2,f(2)+0.1]],linestyle=1,color=black,style=point, thickness=3):
p8:= curve([[3,f(3)], [3,f(3)+0.1]],linestyle=1,color=black,style=point,thickness=3):p9:= curve([[4,f(4)], [4,f(4)+0.1]],linestyle=1,color=black,style=point,thickness=3):

> plots[display]({l,l1,l2,l3,l4,l5,l6,l7,l8,l9,p1,p2,p3,p4,p5,p6,p7,p8,p9});

Chapter 6. Orthogonal Polynomials

Figure 6.1 Some Legendre Polynomials (n=2, 3, and 6)

> with(orthopoly):

> plot({P(6,x),P(3,x),P(2,x)},x=-1..1,color=[black,blue,red],thickness=[1,2,3],tickmarks=[2,2],labels=[x,P],linestyle=[1,3,4]);

Figure 6.2 Some Laguerre Polynomials (alpha=1/2, n= 5, 7, and 10)

> plot({L(5,1/2,x),L(7,1/2,x),L(10,1/2,x)},x=0..9,color=[black,blue,red],thickness=[1,2,3],tickmarks=[2,2],labels=[x,L],linestyle=[1,3,4]);

Figure 6.3 Some Hermite Polynomials (modified by constant multiples, n= 4, 5, and 7).

> plot({2^(-4)*H(4,x),2^(-5)*H(5,x),2^(-7)*H(6,x)},x=-2..2,color=[black,blue,red],thickness=[1,2,3],tickmarks=[2,2],labels=[x,H],linestyle=[1,3,4]);

Chapter 7. Other Orthogonal Systems

Figure 7.1 The Haar mother wavelet

> HaarSc:=t->piecewise(0<=t and t<1,1, t<0 or t>=1,0):

> HaarWv:=t->HaarSc(2*t)-HaarSc(2*t-1):

> plot({HaarSc(t),HaarWv(t)},t=-1..2,-1.2..1.2,linestyle=[1,3],color=[blue,red],tickmarks=[3,3],thickness=[0,2]);

> plot({HaarWv(t)},t=-1..2,-1.2..1.2,tickmarks=[3,3],thickness=2);

Figure 7.2 One of the Rademacher functions

> r0:=t->HaarSc(t); r1:=t->HaarWv(t);

> r:=(t,y)->sum(HaarWv(2^(y-1)*t-n), n=0..2^(y-1));

> plot(r(t,2),t=-0.25..1.6,-1.2..1.2,color=red,thickness=1,tickmarks=[2,2]);

Figure 7.3 Two orthogonal Walsh functions

> w0:=t->HaarSc(t); w1:=t->HaarWv(t);

> w2:=t->HaarWv(2*t)+HaarWv(2*t-1);

> w3:=t->HaarWv(2*t)-HaarWv(2*t-1);

> w4:=t->HaarWv(2*t)+HaarWv(2*t-1);

> w5:=t->HaarWv(2*t)-HaarWv(2*t-1);

> w6:=t->HaarWv(2*t)+HaarWv(2*t-1);

>

> plot(w5(t),t=-0.1..1.1,-1.2..1.2,color=blue);

> plot(w6(t),t=-0.1..1.1,-1.2..1.2,color=red);

Figure 7.4 A bell used for a local cosine basis

> b:=x->piecewise(x<-0.25,0, x<0.25 and x>-0.25,2*(x+0.25), x>0.25 and x<0.75,1,x>0.75 and x<1.25, 2*(1.25-x), 0):

> b1:=x->sqrt(b(x)):

> plot(b(x),x=-0.5..1.5,y=-0.1..1.2,tickmarks=[4,2],labels=[x,b]);

> plot(b1(x),x=-0.5..1.5,y=-0.1..1.2,tickmarks=[4,2],labels=[x,b]);

Figure 7.5 Three elements in the local cosine basis with bell of Figure 7.4

> u:=(k,j,x)->sqrt(2)*b(x-k)*cos(0.5*(2*j+1)*Pi*(x-k));

>

> plot({u(-1,1,x),u(0,1,x),u(1,1,x)},x=-2.5..3.5,y=-1.5..1.5,color=[red,blue,black],labels=[x,u],linestyle=[1,3,4],thickness=[1,2,3]) ;

The local sine bases

> v:=(k,j,x)->sqrt(2)*b(x-k)*sin(0.5*(2*j+1)*Pi*(x-k));

> plot({v(-1,1,x),v(0,1,x),v(1,1,x)},x=-2.5..3.5,y=-1.5..1.5,color=[red,blue,black],labels=[x,v],linestyle=[1,3,4],thickness=[1,2,3]);

Figure 7.6 Two additional elements of the local cosine basis showing the bell

> u:=(k,j,x)->sqrt(2)*b1(x-k)*cos(0.5*(2*j+1)*Pi*(x-k));

> plot({u(1,5,x),sqrt(2)*b1(x-1),u(-1,10,x),sqrt(2)*b1(x+1)},x=-2.5..3.5,y=-1.5..1.5,color=[red,blue,black],labels=[x,u],linestyle=[1,1,1]) ;

Figure 7.7 A biorthogonal pair of scaling functions with the same MRA

The Matlab program

%Use Matlab Wavelet toolbox

%Plot two pairs of scaling functions and wavelet functions for
%compactly supported biorthogonal spline wavelets
%scaling regularity Nd=8, wavelet regularity Nr=6
%Reference: I. Daubechies, Ten lectures on wavelets, 271-280.

[phi1,psi1,phi2, psi2, XVAL] = WAVEFUN('bior6.8',5);
subplot(2,2,1) %pick the upper left of 4 subplots
plot(phi1);
subplot(2,2,2) %pick the upper right of 4 subplots
plot(psi1);
subplot(2,2,3) %pick the upper left of 4 subplots
plot(phi2);
subplot(2,2,4) %pick the upper right of 4 subplots
plot(psi2);

Graphs

Figure 7.8 An example of biorthogonal pairs of scaling functions and

wavelets with compact support

> chi(t):=piecewise(abs(t)<=1,1,0):

> plot(chi(t),t=-2..2, -0.1..1.1):

> BioSc1:=t->chi(t)*(1-abs(t)):

> plot(BioSc1(t), t=-1.5..1.5, -0.1..1.2,ytickmarks=2);

> BioSc1hat:=w->(sin(w/2)/(w/2))^2:

> plot(BioSc1hat(w),w=-20..20,-0.1..1.2,ytickmarks=2);

> a:=w->(1-(2/3)*sin(w/2)^2)^(-1):

> plot(a(w),w=-10..10,0.5..3.5,ytickmarks=5):

> BioSc2hat:=w->a(w)*(sin(w/2)/(w/2))^2:

> plot(BioSc2hat(w),w=-15..15,-0.1..1.5,ytickmarks=2):

> with(inttrans):

> BioSc2:=t->invfourier(BioSc2hat(w),w,t):

> pF:=t->piecewise(1<=t or t<-1,0, -1<=t and t<=0,t+1,0<t and t<1, -t+1):

> plot(pF(t),t=-3..3,F=-0.1..1.1);

> phiF:=t->((3)/2)*((9/8)*pF(t)-(19/64)*pF(t-1)-(19/64)*pF(t+1)+(1/16)*pF(t-2)+(1/16)*pF(t+2)):

>

> plot(phiF(t),t=-4..4,tickmarks=[3,3]);

(This is the biorthogonal function, used the Fourier transform of a(w) to get coefficients.)

Chapter 8. Pointwise Convergence of Wavelet Expansions

Figure 8.1 The delta sequences from Fourier series -- the Drichelet kernel, m=6.

> with(plots):

> CHI:=(x,y)->piecewise(abs(x-y)<Pi,1): plot3d(CHI(x,y),x=-4..4,y=-4..4,color=x):

> Diri_kernel:=(m,x,y)->piecewise(abs(x-y)>0.0000001,CHI(x,y)*sin((m+1/2)*(x-y))/(2*Pi*sin((x-y)/2)),2.069014260):limit(Diri_kernel(6,x,y),{x=y}):

> plot3d(Diri_kernel(6,x,y),x=-4..4,y=-4..4,grid=[40,40]);

> Fejer_kernel:=(m,x,y)->piecewise(abs(x-y)>0.000000001, CHI(x,y)*sin(m*(x-y)/2)^2/((2*Pi*m)*(sin((x-y)/2))^2),.9549296583);evalf(limit(Fejer_kernel(6,x,y),{x=y}));

> plot3d(Fejer_kernel(6,x,y),x=-3.2..3.2,y=-3.2..3.2);

Figure 8.2 The delta sequences from Fourier series -- the Fejer kernel, m=6.

> Fejer_kernel:=(m,x,y)->piecewise(abs(x-y)>0.000000001, CHI(x,y)*sin(m*(x-y)/2)^2/((2*Pi*m)*(sin((x-y)/2))^2),.9549296583);evalf(limit(Fejer_kernel(6,x,y),{x=y}));

> plot3d(Fejer_kernel(6,x,y),x=-3.2..3.2,y=-3.2..3.2);

Figure 8.3 The quasi positive delta sequence qm(x,y) for Daubechies wavelet 2phi(x), m=0,

(programmed using Mathcad).

Figure 8.4 the summability function rhor(x) for Daubechies wavelet 2 phi(x), m=0,

(programmed using Mathcad).

Figure 8.5 The positive delta sequence kr,m(x,y) for Daubechies wavelet 2 phi(x), m=0,

(programmed using Mathcad).

Chapter 9. A shannon Sampling Theorem in Wavelet Subspaces

Figure 9.1 The sampling function for the Daubechies wavelet 2 phi(t) with gamma=-1/sqrt(3).

(programed using MathCad)

Figure 9.2 The function h of Proposition 9.3

> h:=t->piecewise(t>0 and t<1,1-t,t<0 and t>-1,-1-t,t=0,0);

> plot(h(t),t=-1.2..1.2,discont=true,thickness=1,color=blue,tickmarks=[3,3]);

Figure 9.3 The partial sum of the Shannon series expansion

> s:=(t,n)->sum((1-k/(2^n))*sin(Pi*(t*(2^n)-k))/(Pi*((2^n)*t-k)),k=1..2^n)+sum((-1+k/(2^n))*sin(Pi*(t*(2^n)+k))/(Pi*((2^n)*t+k)),k=1..2^n);

> plot(s(t,6),t=-1..1, tickmarks=[3,2]);

Chapter 10. Extensions of Wavelet Sampling Theorems

Figure 10.1 An example of the scaling function of a Meyer wavelet at scale m=0 and sum of its 5 shifts in next scale m=1.

The type one raised cosine

> beta:=1/4:

> phi1:=t->(sin(Pi*(1-beta)*t)+4*beta*t*cos(Pi*(1+beta)*t))/((Pi*t)*(1-(4*beta*t)^2)):

> plot({phi1(2*t)+phi1(1/2)*(phi1(2*t-1)+phi1(2*t+1))+phi1(3/2)*(phi1(2*t-3)+phi1(2*t+3)),phi1(t)},t=-4..4,-0.25..1.3,linestyle=[1,4],color=[red,blue],thickness=[0,2]);

The type two raised cosine.

> phi2:=t->(sin(Pi*(1-beta)*t)+sin(Pi*(1+beta)*t))/((2*Pi*t)*(1+2*beta*t)):

> plot({phi2(2*t)+phi2(1/2)*(phi2(2*t-1)+phi2(2*t+1))+phi2(3/2)*(phi2(2*t-3)+phi2(2*t+3)),phi2(t)},t=-4..4,-0.4..1.3,linestyle=[1,4],color=[red,blue],thickness=[0,2]);

Figure 10.2 The summability function associated with Coiflet of degree 2 with r = 0.22

(programmed using MathCad).

Figure 10.3 The dual of the positive summability function in Figure 10.2

(programmed using MathCad)

Figure 10.4 A non negative function f(x)=exp(-x^2/2)Chi [-1/2,1/2] (x).

> f:=t->piecewise(abs(t)<0.5,exp(-t^2/2),abs(t)>=0,0):

> plot(f(t),t=-0.7..0.7, -0.1..1.2,tickmarks=[3,2]);

Figure 10.5 The positive hybrid sampling series of function in Figure 10.4 (m=4) using Coiflet of degree 2

Figure 10.6 The hybrid series in Figure 10.4

> 

Figure 10.7 An example of cardinal scaling function of Theorem 10.5 -- type two raised cosine wavelet

The scaling function

> beta:=1/4:

> phi2:=t->(sin(Pi*(1-beta)*t)+sin(Pi*(1+beta)*t))/((2*Pi*t)*(1+2*beta*t)):

> a:=limit(phi1(t),t=-1);b:=limit(phi1(t),t=0);

> c:=limit(phi1(t),t=1):

> evalf({a,b,c}):

> plot(phi2(t),t=-5..5, -0.35..1.2,tickmarks=[3,2]);

Figure 10.8 The mother wavelet for the scaling function in Figure 10.7

> psi2:=t->(sin(2*Pi*(1-beta)*t)+sin(2*Pi*(1+beta)*t))/((2*Pi*t)*(1+4*beta*t))-(sin(Pi*(1-beta)*t)+sin(Pi*(1+beta)*t))/((2*Pi*t)*(1-2*beta*t));

> a1:=limit(phi2(t),t=0);b1:=limit(phi2(t),t=1);c1:=limit(phi2(t),t=-1):

> plot(psi2(t),t=-5.5..5.5,-1..1.3,tickmarks=[3,3]);

The figure on the Book over

> plot({phi2(t),psi2(t-1)},t=-5.5..5.5,-1..1.3,resolution=2800,tickmarks=[3,3],color=[red,black],thickness=[3,3]);

plotsetup(ps,plotoptions=`color=rgb, width=6in,height=5in,noborder,landscape`);

Warning, plotoutput file set to postscript.out

> plot({phi2(t),psi2(t-1)}, t=-5.5..5.5,-1..1.3,resolution=2800,tickmarks=[3,3],color=[red,black],thickness=[3,3],axes=none);

Figure 10.9 The scaling function and wavelet for n=1.

Scaling function

> n:=1:

> delta:=t->piecewise(1/n<=t or t<-1/n,0, -1/n<=t and t<=0,n^2*t+n,0<t and t<1/n, -n^2*t+n):

> plot(delta(t-1),t=-1/n..3.5/n,d=-0.8..1.1,tickmarks=[4,2]);

wavelet

> psi1(t):=-.5*delta(2*t-3)+delta(2*t-4)-.5*delta(2*t-5):

> plot(psi1(t),t=-1..3.5,d=-0.7..1.1,tickmarks=[4,2]);

Figure 10.10 The two scaling functions and first wavelets for n=2.

Scaling Functions

> phi21:=t->piecewise(2<=t or t<0,0, 0<=t and t<=1,3*t^2-2*t^3,1<t and t<2,3*(2-t)^2-2*(2-t)^3):

> n:=1:

> plot(phi21(t),t=-1..3,tickmarks=[2,2]);

> phi22:=t->piecewise(2<=t or t<0,0, 0<=t and t<=1,-t^2+t^3,1<t and t<2,(2-t)^2-(2-t)^3):

> plot(phi22(t),t=-1..3,tickmarks=[2,2]);

First wavelet (use PHI(1/2))

> with(linalg):

> PHI:=t->matrix(2,2,[3*t^2-2*t^3,6*t-6*t^2,-t^2+t^3,-2*t+3*t^2]):

> Gamma:=matrix(2,2,[1,0,0,2]):

> C0:=evalm(PHI(1/2)&*inverse(Gamma)):

> C2:=matrix(2,2,[1/2,-3/4,1/8,-1/8]):

> D1:=-Gamma&*C2&*inverse(Gamma):

> D3:=-Gamma&*C0&*inverse(Gamma):

> evalm(D1):

> evalm(D3):

> Phi:=t->matrix(2,1,[phi21(t),phi22(t)]):

> psi:=evalm(D1&*Phi(2*t-1)+inverse(Gamma)&*Phi(2*t-2)+D3&*Phi(2*t-3)):

> simplify(psi):

> plot([psi[1,1],psi[2,1]],t=0..3,color=[red,blue],linestyle=[4,1],thickness=[2,3],tickmarks=[2,2]);

Chapter 11. Translation and Dilation Invariance in Orthogonal Systems

Preliminary

The Haar function and the shifted Haar

> phih:=t->piecewise(0<=t and t<1,1, t<0 or t>=1, 0);

> plot(phih(t),t=-1..2, -0.1..1.2);

> f(t):=phih(t-Pi):

> plot(f(t),t=3..5);

The integral for finding the coefficients

> I1:=int(phih(t*2^4-66)*4,t=Pi..Pi+1);

> evalf(I1);

> I2:=int(phih(t*2^4-50)*4,t=Pi..Pi+1);

> evalf(I2);

Figure 11.1 The approximation of the shifted scaling function by Haar series at scale m=1.

The approximation at m=1

> f2:=t->piecewise(3-t<=0 and t-7/2<0,.358, t-7/2>0 and t-4 < 0,.5,4-t <= 0 and t-9/2 < 0, .1416,t-3 < 0 or 9-2*t <= 0,0):

> plot([2*f2(t),f(t)],t=2..5,color=[red,blue],thickness=[1,2],tickmarks=[3,2],linestyle=[1,3]);

Chapter 12. Analytic Representations Via Orthogonal Series

Figure 12.1 The kernel Kr(t) given in Lemma 12.1 with m=2

The Poisson kernel

> Poisson:=(t,r)->(1/2)*(1-r^2)/(1-2*r*cos(t)+r^2);

The kernel K

> K:=(t,r)->(t^2)*diff(Poisson(t,r),t,t)/(2*Pi);

> plot(K(t,0.9),t=-2..2,-0.3..1.3,tickmarks=[3,2]);

> plot3d(K(t,r),t=-Pi..Pi, r=0..0.9);

Figure 12.2 The Hermite function h 2 (x) and the Hermite function of the second kind h~ 2 (x+i*0).

> with(orthopoly):

> h:=(n,t)->H(n,t)*exp(-t^2/2)/sqrt(2^n*n!*sqrt(Pi));

> h1(2,t):=sqrt(Pi)*(h(2,t)-sqrt(2)*h(0,t))*exp(t^2);

Hermite function of the first kind for n=2

> plot(h(2,t),t=-4..4,-0.6..0.7,tickmarks=[3,3]);

Real part of the Hermite function of the second kind for n=2. Notice the different scales in this and previous plot

> plot(h1(2,t),t=-2..2,-5..15,tickmarks=[3,3]);

The Hermite function of the first kind and the real part of the second kind for n=2

> plot({h1(2,t),h(2,t)},t=-3..3,-3.2..4,tickmarks=[3,3],color=[red,blue],linestyle=[0,3]);

Figure 12.3 A Legendre Polynomials (n=3)

> plot(P(3,t),t=-1..1, tickmarks=[2,3]);

Figure 12.4 The real and imaginary parts of the analytic representation at y=2 for the Legendre polynomial in Figure 12.3

> Q(x+0.2*I):=-1/2*int(P(3,t)/(x+0.2*I-t),t=-1..1):

> Q3_real:=x->Re(Q(x+0.2*I)):

> plot(Q3_real(x),x=-2.2..2.2,-0.35..0.3,tickmarks=[2,3]);

> Q3_im:=x->Im(Q(x+0.2*I)):

> plot(Q3_im(x),x=-2..2,-0.35..0.35,tickmarks=[3,3]);

Chapter 13. Orthogonal Series in Statistics

Figure 13.1 Histogram for data in Example 13.1

A random sample of size 50 based on a N(0, 1) distribution generated using the command:

stats[random, normald](50)

> with(stats):

> with(stats[statplots]):

> # The data was generated randomly and save to random_list,

> # so that the same histogram as in the book can be recovered.

> random_list:=[-1.382965471, .8186285636, -.3144064444, -.1696773949, -.3185691671e-2, -.1786961630, .4416301717, -.6050924001, 1.082710655, .5193987768, -.4177986377, -1.558096387, .5260012674, -.1603814437, -.4670074450, -.7404205817, -1.125369889, 1.471023262, 1.408348602, -2.006670749, .9375835459, -1.293794996, 2.064884399, -.9783131603e-1, -1.655509351, .2100036946e-1, -.2247647997, -1.143647391, .8806099291e-1, -1.087485888, 1.136724954, -.5025710218, 2.111905371, -.3581334786, .6895341061, .3094566959e-1, -.2763354912, .4201157131, 1.141442299, -.2137283818, .7210834223, -2.125973298, 1.492537181, -1.464467418, .5397258939, 1.282077475, .2741506152, .7247422778, -.3835794618, .5292857363e-2];

> histogram(random_list,tickmarks=[6,4],color=blue,area=count);

Figure 13.2 Smooth wavelet estimator for the same data as in Figure 13.1

(Using scaling function of Daubechies 4 wavelet)

The reproducing kernel (program by using MathCad) :

The wavelet density estimator (program by using MathCad) :

Figure 13.3 Daubechies wavelet scaling function (N=4).

(Program using MathCad)

Figuer 13.4 The summability function associated with Daubechies wavelet (N=4).

(Program by using MathCad)

Figuer 13.5 The positive kernel associated with Daubechies wavelet (N=4).

(Program by using MathCad)

Figure 13.6 Density estimate for the Old Faithful geyser data using reproducing kernel associated with Daubechies scaling function (N=4)

The density estimate for the Old Faithful geyser data using reproducing kernel (m=2)

(Program using MathCad)

The density estimate for the Old Faithful geyser data using reproducing (m=3).

(Program using MathCad)

Figure 13.7 Density estimate for the Old Faithful geyser data using positive kernel in Figure 13.4.

Density estimation of the Old Faithful geyser data m=2

(Program using MathCad)

Density estimation of the Old Faithful geyser data m=3.

(Program using MathCad)

(The end)

Please send your comments to:

ggw@csd.uwm.edushen@math.ohiou.edu

Thank you!