<< Chapter < Page | Chapter >> Page > |
You are free to use these programs or any derivative of them for any scientific purpose but please reference this book. Up-dated versionsof these programs and others can be found on our web page at: http: www-dsp.rice.edu/
function p = psa(h,kk)
% p = psa(h,kk) calculates samples of the scaling function% phi(t) = p by kk successive approximations from the
% scaling coefficients h. Initial iteration is a constant.% phi_k(t) is plotted at each iteration. csb 5/19/93
%if nargin==1, kk=11; end; % Default number of iterations
h2= h*2/sum(h); % normalize h(n)K = length(h2)-1; S = 128; % Sets sample density
p = [ones(1,3*S*K),0]/(3*K); % Sets initial iteration
P = p(1:K*S); % Store for later plottingaxis([0 K*S+2 -.5 1.4]);hu = upsam(h2,S); % upsample h(n) by S
for iter = 0:kk % Successive approx.p = dnsample(conv(hu,p)); % convolve and down-sample
plot(p); pause; % plot each iteration% P = [P;p(1:K*S)]; % store each iter. for plottingend
p = p(1:K*S); % only the supported partL = length(p);
x = ([1:L])/(S);
axis([0 3 -.5 1.4]);
plot(x,p); % Final plottitle('Scaling Function by Successive Approx.');
ylabel('Scaling Function');xlabel('x');
function p = pdyad(h,kk)
% p = pdyad(h,kk) calculates approx. (L-1)*2^(kk+2) samples of the% scaling function phi(t) = p by kk+3 dyadic expansions
% from the scaling coefficient vector h where L=length(h).% Also plots phi_k(t) at each iteration. csb 5/19/93
%if nargin==1, kk = 8; end % Default iterations
h2 = h*2/sum(h); % NormalizeN = length(h2); hr = h2(N:-1:1); hh = h2;
axis([0,N-1,-.5,1.4]);
MR = [hr,zeros(1,2*N-2)]; % Generater row for M0
MT = MR; M0 = [];
for k = 1:N-1 % Generate convolution andMR = [0, 0, MR(1:3*N-4)]; % downsample matrix from h(n)MT = [MT; MR];end
M0 = MT(:,N:2*N-1); % M0*p = p if p samples of phiMI = M0 - eye(N);
MJ = [MI(1:N-1,:);ones(1,N)];
pp = MJ\[zeros(N-1,1);1]; % Samples of phi at integers
p = pp(2:length(pp)-1).';x = [0:length(p)+1]*(N-1)/(length(p)+1); plot(x,[0,p,0]); pause
p = conv(h2,p); % value on half integersx = [0:length(p)+1]*(N-1)/(length(p)+1); plot(x,[0,p,0]); pause
y = conv(h2,dnsample(p)); % convolve and downsamplep = merge(y,p); % interleave values on Z and Z/2
x = [0:length(p)+1]*(N-1)/(length(p)+1); plot(x,[0,p,0]); pausefor k = 1:kk
hh = upsample(hh); % upsample coefficientsy = conv(hh,y); % calculate intermediate terms
p = merge(y,p); % insert new terms between oldx = [0:length(p)+1]*(N-1)/(length(p)+1); plot(x,[0,p,0]); pause;
endtitle('Scaling Function by Dyadic Expansion');
ylabel('Scaling Function');xlabel('x');
axis;
function [hf,ht] = pf(h,kk)% [hf,ht] = pf(h,kk) computes and plots hf, the Fourier transform% of the scaling function phi(t) using the freq domain
% infinite product formulation with kk iterations from the scaling% function coefficients h. Also calculates and plots ht = phi(t)
% using the inverse FFT csb 5/19/93if nargin==1, kk=8; end % Default iterations
L = 2^12; P = L; % Sets number of sample pointshp = fft(h,L); hf = hp; % Initializes iteration
plot(abs(hf));pause; % Plots first iterationfor k = 1:kk % Iterations
hp = [hp(1:2:L), hp(1:2:L)]; % Sample
hf = hf.*hp/sqrt(2); % Productplot(abs(hf(1:P/2)));pause; % Plot Phi(omega) each iteration
P=P/2; % Scales axis for plotend;
ht = real(ifft(hf)); % phi(t) from inverse FFTht = ht(1:8*2^kk); plot(ht(1:6*2^kk)); % Plot phi(t)
Notification Switch
Would you like to follow the 'Wavelets and wavelet transforms' conversation and receive update notifications?