$\newcommand{\dotp}[2]{\langle #1, #2 \rangle} \newcommand{\enscond}[2]{\lbrace #1, #2 \rbrace} \newcommand{\pd}[2]{ \frac{ \partial #1}{\partial #2} } \newcommand{\umin}[1]{\underset{#1}{\min}\;} \newcommand{\umax}[1]{\underset{#1}{\max}\;} \newcommand{\umin}[1]{\underset{#1}{\min}\;} \newcommand{\uargmin}[1]{\underset{#1}{argmin}\;} \newcommand{\norm}[1]{\|#1\|} \newcommand{\abs}[1]{\left|#1\right|} \newcommand{\choice}[1]{ \left\{ \begin{array}{l} #1 \end{array} \right. } \newcommand{\pa}[1]{\left(#1\right)} \newcommand{\diag}[1]{{diag}\left( #1 \right)} \newcommand{\qandq}{\quad\text{and}\quad} \newcommand{\qwhereq}{\quad\text{where}\quad} \newcommand{\qifq}{ \quad \text{if} \quad } \newcommand{\qarrq}{ \quad \Longrightarrow \quad } \newcommand{\ZZ}{\mathbb{Z}} \newcommand{\CC}{\mathbb{C}} \newcommand{\RR}{\mathbb{R}} \newcommand{\EE}{\mathbb{E}} \newcommand{\Zz}{\mathcal{Z}} \newcommand{\Ww}{\mathcal{W}} \newcommand{\Vv}{\mathcal{V}} \newcommand{\Nn}{\mathcal{N}} \newcommand{\NN}{\mathcal{N}} \newcommand{\Hh}{\mathcal{H}} \newcommand{\Bb}{\mathcal{B}} \newcommand{\Ee}{\mathcal{E}} \newcommand{\Cc}{\mathcal{C}} \newcommand{\Gg}{\mathcal{G}} \newcommand{\Ss}{\mathcal{S}} \newcommand{\Pp}{\mathcal{P}} \newcommand{\Ff}{\mathcal{F}} \newcommand{\Xx}{\mathcal{X}} \newcommand{\Mm}{\mathcal{M}} \newcommand{\Ii}{\mathcal{I}} \newcommand{\Dd}{\mathcal{D}} \newcommand{\Ll}{\mathcal{L}} \newcommand{\Tt}{\mathcal{T}} \newcommand{\si}{\sigma} \newcommand{\al}{\alpha} \newcommand{\la}{\lambda} \newcommand{\ga}{\gamma} \newcommand{\Ga}{\Gamma} \newcommand{\La}{\Lambda} \newcommand{\si}{\sigma} \newcommand{\Si}{\Sigma} \newcommand{\be}{\beta} \newcommand{\de}{\delta} \newcommand{\De}{\Delta} \renewcommand{\phi}{\varphi} \renewcommand{\th}{\theta} \newcommand{\om}{\omega} \newcommand{\Om}{\Omega} $
This tour explores multiscale computation on a 3D multiresolution sphere using a face-based haar transform.
%load_ext pymatbridge
Starting MATLAB on ZMQ socket ipc:///tmp/pymatbridge Send 'exit' command to kill the server .MATLAB started and connected!
%%matlab
addpath('toolbox_signal')
addpath('toolbox_general')
addpath('toolbox_graph')
addpath('solutions/meshwav_4_haar_sphere')
You need to download numerical_tours and install the IPython notebook to run the code.
You must also install the python-matlab-bridge.
One can define a function on a discrete 3D mesh that assigns a value to each vertex. One can then perform processing of the function according to the geometry of the surface. Here we use a simple sphere.
First compute a multiresolution sphere.
%%matlab
options.base_mesh = 'ico';
options.relaxation = 1;
options.keep_subdivision = 1;
J = 7;
[vertex,face] = compute_semiregular_sphere(J,options);
n = size(face{end},2);
--------------------------------------------------------------------------- RuntimeError Traceback (most recent call last) <ipython-input-8-2b9f76530b0b> in <module>() ----> 1 get_ipython().run_cell_magic('matlab', '', "options.base_mesh = 'ico';\noptions.relaxation = 1;\noptions.keep_subdivision = 1;\nJ = 7;\n[vertex,face] = compute_semiregular_sphere(J,options);\n% n = size(face{end},2);") /Users/gpeyre/anaconda/envs/py3k/lib/python3.3/site-packages/IPython/core/interactiveshell.py in run_cell_magic(self, magic_name, line, cell) 2160 magic_arg_s = self.var_expand(line, stack_depth) 2161 with self.builtin_trap: -> 2162 result = fn(magic_arg_s, cell) 2163 return result 2164 /Users/gpeyre/anaconda/envs/py3k/lib/python3.3/site-packages/pymatbridge/matlab_magic.py in matlab(self, line, cell, local_ns) /Users/gpeyre/anaconda/envs/py3k/lib/python3.3/site-packages/IPython/core/magic.py in <lambda>(f, *a, **k) 191 # but it's overkill for just that one bit of state. 192 def magic_deco(arg): --> 193 call = lambda f, *a, **k: f(*a, **k) 194 195 if callable(arg): /Users/gpeyre/anaconda/envs/py3k/lib/python3.3/site-packages/pymatbridge/matlab_magic.py in matlab(self, line, cell, local_ns) 215 e_s += "\n-----------------------" 216 e_s += "\nAre you sure Matlab is started?" --> 217 raise RuntimeError(e_s) 218 219 RuntimeError: There was an error running the code: options.base_mesh = 'ico'; options.relaxation = 1; options.keep_subdivision = 1; J = 7; [vertex,face] = compute_semiregular_sphere(J,options); % n = size(face{end},2); ----------------------- Are you sure Matlab is started?
Display two examples of sphere.
%%matlab
clf;
for j=[1 2 3 4]
subplot(2, 2, j);
plot_mesh(vertex{j}, face{j});
shading faceted;
end
Comput the center of each face.
%%matlab
x = [];
for i=1:3
v = vertex{end}(i,:);
x(i,:) = mean(v(face{end}));
end
Load an image.
%%matlab
name = 'lena';
%%matlab
M = rescale( load_image(name, 512) );
Display it.
%%matlab
clf;
imageplot(crop(M));
Load a function on the sphere. Use the center of each face to sample the function.
%%matlab
f = rescale( load_spherical_function(M, x, options) );
Display the function on the sphere.
%%matlab
vv = [125,-15];
%%matlab
options.face_vertex_color = f;
clf;
plot_mesh(vertex{end}, face{end}, options);
view(vv);
colormap gray(256);
lighting none;
One can compute low pass approximation by iteratively averaging over 4 neighboring triangles.
Perform one low pass filtering.
%%matlab
f1 = mean( reshape(f, [length(f)/4 4]), 2);
Display.
%%matlab
clf;
options.face_vertex_color = f1;
plot_mesh(vertex{end-1}, face{end-1}, options);
view(vv);
lighting none;
Exercise 1
Compute the successive low pass approximations.
%%matlab
exo1()
%%matlab
%% Insert your code here.
One can compute a wavelet transform by extracting, at each scale, 3 orthogonal wavelet coefficient to represent the orthogonal complement between the successive resolutions.
Precompute the local wavelet matrix, which contains the local vector and three orthognal detail directions.
%%matlab
U = randn(4);
U(:,1) = 1;
[U,R] = qr(U);
Initialize the forward transform.
%%matlab
fw = f;
nj = length(f);
Extract the low pass component and apply the matrix U
%%matlab
fj = fw(1:nj);
fj = reshape(fj, [nj/4 4]);
fj = fj*U;
Store back the coefficients.
%%matlab
fw(1:nj) = fj(:);
nj = nj/4;
Exercise 2
Compute the full wavelet transform, and check for orthogonality (conservation of energy). heck for orthogonality.
%%matlab
exo2()
%%matlab
%% Insert your code here.
Display the coefficients "in place".
%%matlab
clf;
options.face_vertex_color = clamp(fw,-2,2);
plot_mesh(vertex{end}, face{end}, options);
view(vv);
colormap gray(256);
lighting none;
Display the decay of the coefficients.
%%matlab
clf;
plot(fw);
axis([1 n -5 5]);
Exercise 3
Implement the backward spherical Haar transform (replace U by U' to perform the reconstruction), and check for perfect reconstruction.
%%matlab
exo3()
%%matlab
%% Insert your code here.
%%matlab
clf;
options.face_vertex_color = clamp(f1);
plot_mesh(vertex{end}, face{end}, options);
view(vv);
colormap gray(256);
lighting none;
Exercise 4
Perform Haar wavelet approximation with only 10% of the coefficients. orward transform hresholding ackward transform isplay
%%matlab
exo4()
%%matlab
%% Insert your code here.
Exercise 5
Compare with the traditional 2D Haar approximation of |M|.
%%matlab
exo5()
%%matlab
%% Insert your code here.
Exercise 6
Implement Spherical denoising using the Haar transform. Compare it with vertex-based lifting scheme denoising.
%%matlab
exo6()
%%matlab
%% Insert your code here.