s2let_ridgelet_compute_wav - Compute ridgelets Compute ridgelet wavelet and scaling functions. Default usage : [ridgelet_wav, ridgelet_scal] = s2let_compute_wav(L, <options>) where L is harmonic band-limit for the reconstruction on the sphere, [ridgelet_wav, ridgelet_scal] are the reconstruct ridgelet wavelets and scaling functions on the sphere for each j. Options consist of parameter type and value pairs. Valid options include: 'B' = { Dilation factor; B > 1 (default = 2) } 'Spin' = { Spin number; Spin >= 0 (default = 0) } 'J_min' = { Minimum wavelet scale to consider; 0 <= J_min < log_B(L) (default=0) } 'Sampling' = { 'MW' [McEwen & Wiaux sampling (default)], 'MWSS' [McEwen & Wiaux symmetric sampling] } 'Reality' = { false [do not assume f real (default)], true [assume f real (improves performance)] } S2LET package to perform Wavelet transform on the Sphere. Copyright (C) 2012-2015-2014 Boris Leistedt, Martin Büttner & Jason McEwen See LICENSE.txt for license details
0001 function [ridgelet_wav, ridgelet_scal] = ... 0002 s2let_ridgelet_compute_wav(L, varargin) 0003 % s2let_ridgelet_compute_wav - Compute ridgelets 0004 % 0005 % Compute ridgelet wavelet and scaling functions. 0006 % 0007 % Default usage : 0008 % 0009 % [ridgelet_wav, ridgelet_scal] = s2let_compute_wav(L, <options>) 0010 % 0011 % where L is harmonic band-limit for the reconstruction on the sphere, 0012 % [ridgelet_wav, ridgelet_scal] are the reconstruct ridgelet wavelets and 0013 % scaling functions on the sphere for each j. 0014 % 0015 % Options consist of parameter type and value pairs. 0016 % Valid options include: 0017 % 0018 % 'B' = { Dilation factor; B > 1 (default = 2) } 0019 % 'Spin' = { Spin number; Spin >= 0 (default = 0) } 0020 % 'J_min' = { Minimum wavelet scale to consider; 0021 % 0 <= J_min < log_B(L) (default=0) } 0022 % 'Sampling' = { 'MW' [McEwen & Wiaux sampling (default)], 0023 % 'MWSS' [McEwen & Wiaux symmetric sampling] } 0024 % 'Reality' = { false [do not assume f real (default)], 0025 % true [assume f real (improves performance)] } 0026 % 0027 % S2LET package to perform Wavelet transform on the Sphere. 0028 % Copyright (C) 2012-2015-2014 Boris Leistedt, Martin Büttner & Jason McEwen 0029 % See LICENSE.txt for license details 0030 0031 % Parse arguments. 0032 p = inputParser; 0033 p.addRequired('L', @isnumeric); 0034 p.addParamValue('B', 2, @isnumeric); 0035 p.addParamValue('Spin', 0, @isnumeric); 0036 p.addParamValue('J_min', 0, @isnumeric); 0037 p.addParamValue('Sampling', 'MW', @ischar); 0038 p.addParamValue('Reality', false, @islogical); 0039 p.parse(L, varargin{:}); 0040 args = p.Results; 0041 0042 N = 1; 0043 s = args.Spin; 0044 0045 % Compute ring in harmonic space. 0046 ring_lm = zeros(args.L^2,1); 0047 for el = max([0 abs(args.Spin)]):args.L-1 0048 0049 % Spin zero setting. 0050 % elon2 = el./2.0; 0051 % logp = gammaln(2*elon2+1) - 2*elon2 * log(2) - 2 * gammaln(elon2+1); 0052 % p0 = real((-1).^elon2) .* exp(logp); 0053 0054 % Spin setting but computed via associated Legendre functions (slow). 0055 % P = legendre(el, 0); 0056 % p0 = P(abs(args.Spin)+1); 0057 0058 % Spin setting but computed explicitly. 0059 0060 logp2 = gammaln(el+s+1) - el * log(2) - gammaln((el+s)./2+1) - gammaln((el-s)./2+1); 0061 p0 = real((-1).^((el+s)./2)) .* exp(logp2); 0062 0063 ind = ssht_elm2ind(el, 0); 0064 ring_lm(ind) = 2*pi * sqrt((2*el+1)/(4*pi)) * p0; 0065 ring_lm(ind) = ring_lm(ind) .* ... 0066 (-1).^args.Spin .* sqrt(exp(gammaln(el-s+1) - gammaln(el+s+1))); 0067 0068 end 0069 0070 0071 %% Compute ridgets via wavelet transforms 0072 % Can only do this for scalar setting 0073 0074 % % Compute ring on the sphere. 0075 % ring = ssht_inverse(ring_lm, L, ... 0076 % 'Reality', args.Reality, ... 0077 % 'Method', args.Sampling, ... 0078 % 'Spin', args.Spin); 0079 % 0080 % % Compute ridgelets. 0081 % [ridgelet_wav, ridgelet_scal] = ... 0082 % s2let_transform_analysis_mw(ring, 'B', args.B, 'J_min', args.J_min, ... 0083 % 'N', N, 'Upsample', true, 'Spin', args.Spin, ... 0084 % 'Reality', args.Reality, 'Sampling', args.Sampling); 0085 0086 % [ridgelet_wav, ridgelet_scal] = ... 0087 % s2let_transform_axisym_analysis_mw(ring, 'L', L, 'B', args.B, 'J_min', args.J_min, ... 0088 % 'Upsample', true, ... 0089 % 'Reality', args.Reality); 0090 0091 0092 %% Compute ridgelets directly 0093 % Can then compute explicitly for spin functions. 0094 % But not that spin ridgelets are identical to scalar wavelets since can 0095 % always view as convolution of ring with scalar wavelets. 0096 0097 J = s2let_jmax(args.L, args.B); 0098 [kappa, kappa0] = s2let_transform_axisym_tiling(args.B, args.L, args.J_min); 0099 0100 ridgelet_scal_lm = zeros(args.L^2,1); 0101 for el = 0:args.L-1 0102 ind = ssht_elm2ind(el, 0); 0103 scal_lm = kappa0(el+1) .* sqrt((2*el+1)/(4*pi)); 0104 ridgelet_scal_lm(ind) = scal_lm .* sqrt(4*pi / (2*el+1)) .* ring_lm(ind); 0105 end 0106 ridgelet_scal = ssht_inverse(ridgelet_scal_lm, L, 'Spin', args.Spin, ... 0107 'Reality', args.Reality, 'Method', args.Sampling); 0108 0109 ridgelet_wav = cell(J-args.J_min+1,1); 0110 for j = args.J_min:J 0111 0112 ridgelet_wav_lm = zeros(args.L^2,1); 0113 for el = 0:args.L-1 0114 ind = ssht_elm2ind(el, 0); 0115 wav_lm = kappa(j+1, el+1) .* sqrt((2*el+1)/(8*pi.^2)); 0116 ridgelet_wav_lm(ind) = wav_lm .* sqrt(4*pi / (2*el+1)) .* ring_lm(ind); 0117 end 0118 ridgelet_wav{j-args.J_min+1} = ssht_inverse(ridgelet_wav_lm, L, ... 0119 'Spin', args.Spin, 'Reality', args.Reality, 'Method', args.Sampling); 0120 0121 end 0122 0123 0124 0125 0126 0127