s2let_demo8 - curvelet analysis on the Earth topographic map Plot curvelet coefficients on multiple Mollweide projections. The function generates one plot of the scaling function contribution and a grid of plots for each orientation of each scale of the curvelet contributions. ----------------------------------------------------------- f_cur is cell array with all the curvelet coefficients. Its first index is the curvelet scale j, the second index is the orientation g, and each element is a function on the sphere in MW sampling. scal is the corresponding scaling function contribution (i.e. just a single function on the sphere). B is the curvelet parameter. L is the angular band-limit. J_min is the first curvelet scale in f_cur. Options consist of parameter type and value pairs. Valid options include: 'Spin' = { non-negative integers (default=0) } 'Reality' = { false [do not assume f real (default)], true [assume f real (improves performance)] } 'Upsample' = { false [multiresolution algorithm (default)], true [full resolution curvelets] } 'SpinLowered' = { true [Apply normalisation factors for spin-lowered curvelets and scaling function.], false [Apply the usual normalisation factors such that the curvelets fulfill the admissibility condition (default)]} 'SpinLoweredFrom' = [integer; if the SpinLowered option is used, this option indicates which spin number the curvelets should be lowered from (default = 0)] For ssht_plot_mollweide(f, L, <options>), f is the sampled function and L is the harmonic band-limit. Valid options include: 'ColourBar' = { false [do not add colour bar (default)], true [add colour bar] } 'Mode' = { 0 Plot amplitude, or modulus is f complex (default), 1 Plot real part, 2 Plot imaginaty part, 3 Plot modulus and arrows for real/img angle } 'Spin' = { non-negative integers (default=0) } ----------------------------------------------------------- S2LET package to perform Wavelet transform on the Sphere. Copyright (C) 2012-2016 Boris Leistedt, Jennifer Chan & Jason McEwen See LICENSE.txt for license details -----------------------------------------------------------
0001 % s2let_demo8 - curvelet analysis on the Earth topographic map 0002 % 0003 % Plot curvelet coefficients on multiple Mollweide projections. 0004 % The function generates one plot of the scaling function 0005 % contribution and a grid of plots for each orientation of 0006 % each scale of the curvelet contributions. 0007 % 0008 % ----------------------------------------------------------- 0009 % f_cur is cell array with all the curvelet coefficients. 0010 % Its first index is the curvelet scale j, the second 0011 % index is the orientation g, and each element is a 0012 % function on the sphere in MW sampling. 0013 % scal is the corresponding scaling function contribution 0014 % (i.e. just a single function on the sphere). 0015 % 0016 % B is the curvelet parameter. 0017 % L is the angular band-limit. 0018 % J_min is the first curvelet scale in f_cur. 0019 % 0020 % Options consist of parameter type and value pairs. 0021 % Valid options include: 0022 % 'Spin' = { non-negative integers (default=0) } 0023 % 'Reality' = { false [do not assume f real (default)], 0024 % true [assume f real (improves performance)] } 0025 % 'Upsample' = { false [multiresolution algorithm (default)], 0026 % true [full resolution curvelets] } 0027 % 'SpinLowered' = { true [Apply normalisation factors for 0028 % spin-lowered curvelets and 0029 % scaling function.], 0030 % false [Apply the usual normalisation factors such 0031 % that the curvelets fulfill the admissibility 0032 % condition (default)]} 0033 % 'SpinLoweredFrom' = [integer; if the SpinLowered option is used, this 0034 % option indicates which spin number the curvelets 0035 % should be lowered from (default = 0)] 0036 % 0037 % 0038 % For ssht_plot_mollweide(f, L, <options>), 0039 % f is the sampled function and L is the harmonic band-limit. 0040 % Valid options include: 0041 % 'ColourBar' = { false [do not add colour bar (default)], 0042 % true [add colour bar] } 0043 % 'Mode' = { 0 Plot amplitude, or modulus is f complex (default), 0044 % 1 Plot real part, 0045 % 2 Plot imaginaty part, 0046 % 3 Plot modulus and arrows for real/img angle } 0047 % 'Spin' = { non-negative integers (default=0) } 0048 % ----------------------------------------------------------- 0049 % S2LET package to perform Wavelet transform on the Sphere. 0050 % Copyright (C) 2012-2016 Boris Leistedt, Jennifer Chan & Jason McEwen 0051 % See LICENSE.txt for license details 0052 % ----------------------------------------------------------- 0053 0054 clear all; 0055 close all ; 0056 0057 load('EGM2008_Topography_flms_L0128'); 0058 % --------------- 0059 % band limit the signals: 0060 % --------------- 0061 L = 64; % (set L>=64 to see multi-resolution effects) 0062 flm_gen = flm(1:L^2,1); 0063 f_gen = ssht_inverse(flm_gen, L,'Reality', true); 0064 0065 % --------------- 0066 % Set reality flag for the transform (to improve performance if f is real): 0067 % --------------- 0068 reality = true; 0069 0070 % --------------- 0071 % Define curvelet parameters: 0072 % --------------- 0073 Spin = 0; % Spin value of wavelet 0074 B = 2; % B=2 for dyadic sampling 0075 J_min = 3; % Minimum scale probed by curvelets 0076 J =s2let_jmax(L, B); % Maximum scale probed by curvelets =ceil(log L/ log B); 0077 0078 % --------------- 0079 % Define Plotting parameters: 0080 % --------------- 0081 zoomfactor= 1.0; 0082 ny = 7; 0083 nx = 3; 0084 0085 maxfigs = nx*ny; 0086 pltroot = '../../../figs' ; 0087 configstr = ['_L',int2str(L),'_B',int2str(B),'_Jmin',int2str(J_min)]; 0088 0089 % ================================================ 0090 % FULL RESOLUTION PLOT (Upsample: true) 0091 % ================================================ 0092 % ----------------- 0093 % Signal analysis: 0094 % ----------------- 0095 [f_cur, f_scal] = s2let_transform_curvelet_analysis_px2cur(f_gen, ... 0096 'B', B, 'L', L, ... 0097 'J_min', J_min, ... 0098 'Spin', Spin, ... 0099 'Reality', reality, ... 0100 'Upsample', true , ... 0101 'SpinLowered', false, ... 0102 'SpinLoweredFrom', 0); 0103 0104 figure('Position',[20 20 1500 1400]) %100 100 1300 1000 0105 subplot(ny, nx, 1); 0106 % --- plot initial data --- % 0107 ssht_plot_mollweide(f_gen, L, 'Mode', 1); 0108 % 0109 title('Initial data','FontSize',8) 0110 campos([0 0 -1]); camup([0 1 0]); zoom(zoomfactor) 0111 v = caxis; 0112 temp = max(abs(v)); 0113 caxis([-temp temp]) 0114 % 0115 subplot(ny, nx, 2); 0116 % --- plot scaling function contributions --- % 0117 ssht_plot_mollweide(f_scal, L, 'Mode', 1); 0118 % 0119 title('Scaling fct','FontSize',8) 0120 campos([0 0 -1]); camup([0 1 0]); zoom(zoomfactor) 0121 v = caxis; 0122 temp = max(abs(v)); 0123 caxis([-temp temp]) 0124 % --- plot curvelet kernel contributions --- % 0125 ind = 2; 0126 for j = J_min:J 0127 band_limit = min([ s2let_bandlimit(j,J_min,B,L) L ]); 0128 % Nj is the orientational band-limit at j-th scale. 0129 Nj = band_limit; 0130 if (reality == false) 0131 enmax = 2*Nj-1; 0132 else 0133 enmax = Nj; 0134 end 0135 for en = 1: enmax 0136 ind = ind + 1; 0137 if ind <= maxfigs 0138 subplot(ny, nx, ind); 0139 % 0140 ssht_plot_mollweide(reshape(f_cur{j-J_min+1}(en,:,:), L, 2*L-1), L, 'Mode', 1); 0141 % 0142 campos([0 0 -1]); camup([0 1 0]); zoom(zoomfactor) 0143 v = caxis; 0144 temp = max(abs(v)); 0145 caxis([-temp temp]) 0146 title(['Curvelet scale j=',int2str(j)-J_min+1,', n=',int2str(en)],'FontSize', 8) 0147 end 0148 end 0149 end 0150 colormap(jet) 0151 fname = [pltroot,'/s2let_demo8_', configstr, '_curvelet_EarthTopo_FullRes.png'] 0152 print('-r200', '-dpng', fname) 0153 0154 0155 % ---------- 0156 % Plot reconstructed map (to compare with the initial map): 0157 % ---------- 0158 % ----------------- 0159 % Signal synthesis: 0160 % ----------------- 0161 f_rec = s2let_transform_curvelet_synthesis_cur2px(f_cur, f_scal, ... 0162 'B', B, 'L', L, ... 0163 'J_min', J_min, ... 0164 'Spin', Spin, ... 0165 'Reality', reality, ... 0166 'Upsample', true, ... 0167 'SpinLowered', false, ... 0168 'SpinLoweredFrom', 0); 0169 0170 figure('Position',[100 100 900 200]) 0171 subplot(2, 2, 1); 0172 % --- plot initial data --- % 0173 ssht_plot_mollweide(f_gen, L, 'Mode', 1); 0174 title('initial signal') 0175 hold on 0176 subplot(2, 2, 2); 0177 % --- plot reconstructed data --- % 0178 ssht_plot_mollweide(f_rec,L, 'Mode', 1); 0179 title('reconstructed signal') 0180 fname = [pltroot,'/s2let_demo8_', configstr, '_spin_curvelet_EarthTomo_FullRes_Int_Rec_signal.png'] 0181 print('-r200', '-dpng', fname) 0182 % Check error: 0183 check_error = max(abs(f_gen(:)-f_rec(:))) 0184 0185 0186 0187 % ================================================ 0188 % MULTI-RESOLUTION PLOT (Upsample: false) 0189 % ================================================ 0190 [f_cur, f_scal] = s2let_transform_curvelet_analysis_px2cur(f_gen, ... 0191 'B', B, 'L', L, ... 0192 'J_min', J_min, ... 0193 'Spin', Spin, ... 0194 'Reality', reality, ... 0195 'Upsample', false, ... 0196 'SpinLowered', false, ... 0197 'SpinLoweredFrom', 0); 0198 0199 figure('Position',[20 20 1500 1400]) %100 100 1300 1000 0200 subplot(ny, nx, 1); 0201 % --- plot initial data --- % 0202 ssht_plot_mollweide(f_gen, L, 'Mode', 1); 0203 % 0204 title('Initial data','FontSize',8) 0205 campos([0 0 -1]); camup([0 1 0]); zoom(zoomfactor) 0206 v = caxis; 0207 temp = max(abs(v)); 0208 caxis([-temp temp]) 0209 % 0210 subplot(ny, nx, 2); 0211 % --- plot scaling function contributions --- % 0212 bl = min([ s2let_bandlimit(J_min-1,J_min,B,L) L ]); 0213 ssht_plot_mollweide(f_scal, bl, 'Mode', 1); 0214 % 0215 title('Scaling fct','FontSize',8) 0216 campos([0 0 -1]); camup([0 1 0]); zoom(zoomfactor) 0217 v = caxis; 0218 temp = max(abs(v)); 0219 caxis([-temp temp]) 0220 % --- plot curvelet kernel contributions --- % 0221 ind = 2; 0222 for j = J_min:J 0223 band_limit = min([ s2let_bandlimit(j,J_min,B,L) L ]); 0224 Nj = band_limit; 0225 if (reality == false) 0226 enmax = 2*Nj-1; 0227 else 0228 enmax = Nj; 0229 end 0230 for en = 1: enmax 0231 ind = ind + 1; 0232 if ind <= maxfigs 0233 subplot(ny, nx, ind); 0234 % 0235 ssht_plot_mollweide(reshape(f_cur{j-J_min+1}(en,:), band_limit , 2*band_limit -1), band_limit , 'Mode', 1); 0236 % 0237 campos([0 0 -1]); camup([0 1 0]); zoom(zoomfactor) 0238 v = caxis; 0239 temp = max(abs(v)); 0240 caxis([-temp temp]) 0241 title(['Curvelet scale j=',int2str(j)-J_min+1,', n=',int2str(en)], 'FontSize', 8) 0242 end 0243 end 0244 end 0245 colormap(jet) 0246 fname = [pltroot,'/s2let_demo8_', configstr, '_curvelet_EarthTopo_MultiRes.png'] 0247 print('-r200', '-dpng', fname) 0248 0249 0250 % ---------- 0251 % Plot reconstructed map (to compare with the initial map): 0252 % ---------- 0253 % ----------------- 0254 % Signal synthesis: 0255 % ----------------- 0256 f_rec = s2let_transform_curvelet_synthesis_cur2px(f_cur, f_scal, ... 0257 'B', B, 'L', L, ... 0258 'J_min', J_min, ... 0259 'Spin', Spin, ... 0260 'Reality', reality, ... 0261 'Upsample', false, ... 0262 'SpinLowered', false, ... 0263 'SpinLoweredFrom', 0); 0264 0265 figure('Position',[100 100 900 200]) 0266 subplot(2, 2, 1); 0267 % --- plot initial data --- % 0268 ssht_plot_mollweide(f_gen,L, 'Mode', 1); 0269 title('initial signal') 0270 hold on 0271 subplot(2, 2, 2); 0272 % --- plot reconstructed data --- % 0273 ssht_plot_mollweide(f_rec,L, 'Mode', 1); 0274 title('reconstructed signal') 0275 % Check error: 0276 check_error = max(abs(f_gen(:)-f_rec(:))) 0277 0278 fname = [pltroot,'/s2let_demo8_', configstr, '_curvelet_EarthTopo_MultiRes_Int_Rec_signal.png'] 0279 print('-r200', '-dpng', fname) 0280