s2let_transform_curvelet_synthesis_lmn2lm
Compute (spin) curvelet transform,
input in Wigner space,
output in harmonic space.
Default usage :
flm_rec = s2let_transform_curvelet_synthesis_lmn2lm(f_cur_lmn, f_scal_lm, <options>)
f_cur_lmn is the input curvelet contribution,
f_cur_lmn is the input scaling function contribution,
cur_lm is the curvelet kernels,
scal_l is the scaling function kernel.
Option :
'B' = { Dilation factor; B > 1 (default=2) }
'L' = { Harmonic band-limit; L > 1 (default=guessed from input) }
'Spin' = { Spin; (default=0) }
'J_min' = { Minimum curvelet scale to consider;
0 <= J_min < log_B(L) (default=0) }
'Upsample' = { false [multiresolution algorithm (default)],
true [full resolution curvelets] }
'Reality' = { false [do not assume corresponding signal f real (default)],
true [assume f real (improves performance)] }
'SpinLowered' = { true [Apply normalisation factors for spin-lowered
curvelets and scaling function.],
false [Apply the usual normalisation factors such
that the curvelets fulfil 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)]
-----------------------------------------------------------
S2LET package to perform Wavelet Transform on the Sphere.
Copyright (C) 2015 Boris Leistedt, Martin B¸ttner,
Jennifer Chan & Jason McEwen
See LICENSE.txt for license details
-----------------------------------------------------------0001 function flm_rec= s2let_transform_curvelet_synthesis_lmn2lm(f_cur_lmn, f_scal_lm, cur_lm, scal_l, varargin) 0002 0003 % s2let_transform_curvelet_synthesis_lmn2lm 0004 % Compute (spin) curvelet transform, 0005 % input in Wigner space, 0006 % output in harmonic space. 0007 % 0008 % Default usage : 0009 % 0010 % flm_rec = s2let_transform_curvelet_synthesis_lmn2lm(f_cur_lmn, f_scal_lm, <options>) 0011 % 0012 % f_cur_lmn is the input curvelet contribution, 0013 % f_cur_lmn is the input scaling function contribution, 0014 % cur_lm is the curvelet kernels, 0015 % scal_l is the scaling function kernel. 0016 % 0017 % Option : 0018 % 'B' = { Dilation factor; B > 1 (default=2) } 0019 % 'L' = { Harmonic band-limit; L > 1 (default=guessed from input) } 0020 % 'Spin' = { Spin; (default=0) } 0021 % 'J_min' = { Minimum curvelet scale to consider; 0022 % 0 <= J_min < log_B(L) (default=0) } 0023 % 'Upsample' = { false [multiresolution algorithm (default)], 0024 % true [full resolution curvelets] } 0025 % 'Reality' = { false [do not assume corresponding signal f real (default)], 0026 % true [assume f real (improves performance)] } 0027 % 'SpinLowered' = { true [Apply normalisation factors for spin-lowered 0028 % curvelets and scaling function.], 0029 % false [Apply the usual normalisation factors such 0030 % that the curvelets fulfil the admissibility 0031 % condition (default)]} 0032 % 'SpinLoweredFrom' = [integer; if the SpinLowered option is used, this 0033 % option indicates which spin number the curvelets 0034 % should be lowered from (default = 0)] 0035 % ----------------------------------------------------------- 0036 % S2LET package to perform Wavelet Transform on the Sphere. 0037 % Copyright (C) 2015 Boris Leistedt, Martin B¸ttner, 0038 % Jennifer Chan & Jason McEwen 0039 % See LICENSE.txt for license details 0040 % ----------------------------------------------------------- 0041 0042 sz = length(scal_l(:)); 0043 Lguessed = sqrt(sz); 0044 0045 p = inputParser; 0046 p.addRequired('f_cur_lmn', @iscell); 0047 p.addRequired('f_scal_lm', @isnumeric); 0048 p.addRequired('cur_lm', @iscell); 0049 p.addRequired('scal_l', @isnumeric); 0050 p.addParamValue('B', 2, @isnumeric); 0051 p.addParamValue('L', Lguessed, @isnumeric); 0052 p.addParamValue('J_min', 0, @isnumeric); 0053 p.addParamValue('Spin', 0, @isnumeric); 0054 p.addParamValue('Upsample', false, @islogical); 0055 p.addParamValue('Sampling', 'MW', @ischar); 0056 p.addParamValue('Reality', false, @islogical); 0057 p.addParamValue('SpinLowered', false, @islogical); 0058 p.addParamValue('SpinLoweredFrom', 0, @isnumeric); 0059 p.parse(f_cur_lmn, f_scal_lm, cur_lm, scal_l, varargin{:}); 0060 args = p.Results; 0061 0062 J = s2let_jmax(args.L, args.B); 0063 0064 % --------------- 0065 % Signal synthesis: 0066 % --------------- 0067 % Curvelet contribution: 0068 flm_rec = zeros(args.L^2,1); 0069 for j = args.J_min:J, 0070 band_limit = min([ s2let_bandlimit(j,args.J_min,args.B,args.L) args.L ]); 0071 % Nj = orientational band-limit at j-th scale 0072 Nj = band_limit; 0073 if (args.Reality == 0) 0074 for en = -Nj+1:Nj-1, 0075 for el = max(abs(args.Spin),abs(en)):band_limit-1, 0076 ind_ln = ssht_elm2ind(el, en); 0077 psi = (cur_lm{j-args.J_min+1}(ind_ln)); 0078 for m = -el:el, 0079 ind_lm = ssht_elm2ind(el, m); 0080 if (args.Upsample ~= 0) 0081 ind_lmn = so3_elmn2ind(el,m,en,args.L,Nj); 0082 else 0083 ind_lmn = so3_elmn2ind(el,m,en,band_limit,Nj); 0084 end 0085 flm_rec(ind_lm)= flm_rec(ind_lm)+ f_cur_lmn{j-args.J_min+1}(ind_lmn)* psi; 0086 end 0087 end 0088 end 0089 else % real setting: 0090 for el = 0:band_limit-1, 0091 for m = -el:el, 0092 ind_lm = ssht_elm2ind(el, m); 0093 % for n=0 terms: 0094 ind_lnzero = ssht_elm2ind(el, 0); 0095 if (args.Upsample ~= 0) 0096 ind_lmnzero = so3_elmn2ind(el,m,0,args.L,Nj,'Reality', args.Reality); 0097 else 0098 ind_lmnzero = so3_elmn2ind(el,m,0,band_limit,Nj,'Reality', args.Reality); 0099 end 0100 psizero= cur_lm{j-args.J_min+1}(ind_lnzero); 0101 flm_rec(ind_lm)= flm_rec(ind_lm)+ f_cur_lmn{j-args.J_min+1}(ind_lmnzero)* psizero; 0102 % for n ~= 0 terms: 0103 if (args.Upsample ~= 0) 0104 ind_lml = so3_elmn2ind(el,m,el,args.L,Nj,'Reality', args.Reality); 0105 ind_l_nm_l = so3_elmn2ind(el,-m,el,args.L,Nj,'Reality', args.Reality); 0106 else 0107 ind_lml = so3_elmn2ind(el,m,el,band_limit,Nj,'Reality', args.Reality); 0108 ind_l_nm_l = so3_elmn2ind(el,-m,el,band_limit,Nj,'Reality', args.Reality); 0109 end 0110 if (mod((m+el),2) == 1) 0111 sign = -1; 0112 else 0113 sign = 1; 0114 end 0115 for en = 1:Nj-1, 0116 if (el >= en) 0117 % contribution from positive n terms % i.e. Sum_over_lmn((cur_l_n)*f_cur_l_m_n)) 0118 ind_ln = ssht_elm2ind(el, en); 0119 psi = cur_lm{j-args.J_min+1}(ind_ln); 0120 flm_rec(ind_lm)= flm_rec(ind_lm)+ f_cur_lmn{j-args.J_min+1}(ind_lml)* psi; 0121 % contribution from negative n terms: 0122 % i.e. Sum_over_lmn((cur_l_neg-n)*f_cur_l_m_neg-n)) =Sum_over_lmn((cur_l_neg-n)*(-1^(m+n))*f_cur_l_neg-m_n)) 0123 ind_l_nn = ssht_elm2ind(el, -en); 0124 npsi = cur_lm{j-args.J_min+1}(ind_l_nn); 0125 flm_rec(ind_lm)= flm_rec(ind_lm)+ sign*conj(f_cur_lmn{j-args.J_min+1}(ind_l_nm_l))* npsi; 0126 end % end if (el>=en) loop 0127 end % end en-loop for Reality Option 0128 end % end em-loop for Reality Option 0129 end % end el-loop for Reality Option 0130 end % end if-loop for Reality Option 0131 end % end j-loop 0132 0133 % ----------------- 0134 % Add the scaling function contribution: 0135 % ----------------- 0136 if (args.Upsample == 0) 0137 band_limit = min([ s2let_bandlimit(args.J_min-1, args.J_min, args.B,args.L) args.L ]); 0138 else 0139 band_limit = args.L ; 0140 end 0141 lm_ind=0; 0142 if (args.Reality == 0) 0143 for el = abs(args.Spin): band_limit-1, 0144 phi = sqrt(4.*pi/(2.*el+1))*scal_l(el^2+el+1,1); 0145 for m = -el:el, 0146 lm_ind=ssht_elm2ind(el, m); 0147 flm_rec(lm_ind) = flm_rec(lm_ind)+ f_scal_lm(lm_ind)* phi; 0148 end 0149 end 0150 else % real setting: 0151 for el = 0 :band_limit-1, 0152 phi = sqrt(4.*pi/(2.*el+1))*scal_l(el^2+el+1,1); 0153 for m = -el:el, 0154 lm_ind=ssht_elm2ind(el, m); 0155 flm_rec(lm_ind) = flm_rec(lm_ind)+ f_scal_lm(lm_ind)* phi; 0156 end 0157 end 0158 end % end if-loop for Reality Option 0159 0160 end 0161 0162 0163 0164 0165 0166