0001
0002
0003
0004
0005
0006
0007
0008
0009
0010
0011
0012 function demo1
0013 close all
0014 fprintf('Welcome to SGWT demo #1\n');
0015
0016 npoints=500;
0017 fprintf('Creating Swiss Roll point cloud with %g points\n',npoints);
0018 dataparams=struct('n',npoints,'dataset',-1','noise',0,'state',0);
0019 r=create_synthetic_dataset(dataparams);
0020 x=rescale_center(r.x);
0021
0022
0023 fprintf('Computing edge weights and graph Laplacian\n');
0024
0025 d=distanz(x);
0026 s=.1;
0027 A=exp(-d.^2/(2*s^2));
0028 L=full(sgwt_laplacian(A));
0029
0030 fprintf('Measuring largest eigenvalue, lmax = ');
0031
0032
0033 Nscales=4;
0034 lmax=sgwt_rough_lmax(L);
0035 fprintf('%g\n',lmax);
0036
0037 fprintf('Designing transform in spectral domain\n');
0038 [g,gp,t]=sgwt_filter_design(lmax,Nscales);
0039 arange=[0 lmax];
0040
0041 figure
0042 sgwt_view_design(g,t,arange);
0043 ylim([0 3])
0044 set(gcf,'position',[0 780,600,300])
0045
0046 m=50;
0047 fprintf('Computing Chebyshev polynomials of order %g for fast transform \n',m);
0048 for k=1:numel(g)
0049 c{k}=sgwt_cheby_coeff(g{k},m,m+1,arange);
0050 end
0051
0052
0053 jcenter=32;
0054 fprintf('Computing forward transform of delta at vertex %g\n',jcenter);
0055 N=size(L,1);
0056 d=sgwt_delta(N,jcenter);
0057
0058 wpall=sgwt_cheby_op(d,L,c,arange);
0059
0060 fprintf('Displaying wavelets\n');
0061 msize=100;
0062 cp=[-1.4,-16.9,3.4];
0063
0064
0065
0066 ws=300;
0067 figure;
0068 xp=0; yp=ws+100;
0069 set(gcf,'position',[xp,yp,ws-10,ws+10]);
0070 scatter3(x(1,:),x(2,:),x(3,:),msize,d,'.');
0071 set(gcf,'Colormap',[.5 .5 .5;1 0 0]);
0072 clean_axes(cp);
0073 title(sprintf('Vertex %g',jcenter));
0074
0075
0076 for n=1:Nscales+1
0077 wp=wpall{n};
0078 figure
0079 xp=mod(n,3)*(ws+10);
0080 yp=(1-floor((n)/3))*(ws+100);
0081 set(gcf,'position',[xp,yp,ws-10,ws+10]);
0082 scatter3(x(1,:),x(2,:),x(3,:),msize,wp,'.');
0083 colormap jet
0084 caxis([-1 1]*max(abs(wp)));
0085 clean_axes(cp);
0086
0087 hcb=colorbar('location','north');
0088 cxt=get(hcb,'Xtick');
0089 cxt=[cxt(1),0,cxt(end)];
0090 set(hcb,'Xtick',cxt);
0091 cpos=get(hcb,'Position');
0092 cpos(4)=.02;
0093 set(hcb,'Position',cpos);
0094 set(hcb,'Position',[.25 .91 .6 .02]);
0095
0096 if n==1
0097 title('Scaling function');
0098 else
0099 title(sprintf('Wavelet at scale j=%g, t_j = %0.2f',n-1,t(end+1-(n-1))));
0100
0101 end
0102 end
0103
0104
0105 function clean_axes(cp)
0106 xlim([-1 1]);ylim([-1 1]);zlim([-1 1]);
0107 set(gca,'Xtick',[-1 0 1]);
0108 set(gca,'Ytick',[-1 0 1]);
0109 set(gca,'Ztick',[-1 0 1]);
0110 axis square
0111 set(gca,'CameraPosition',cp);