Home > Matching_SYM_LSM > src > 2D_LSM_62 > 2D_LSM_62_sym_warp > LSM_62_sym_warp_estimate_f.m

LSM_62_sym_warp_estimate_f

PURPOSE ^

% estimate function f from g and h using xa and smooth.

SYNOPSIS ^

function [fe,fes,fgi,fhi] = LSM_62_sym_warp_estimate_f(g,h,vg,vh,Nf,xa,sigma,pt)

DESCRIPTION ^

% estimate function f from g and h using xa and smooth.

 [fe,fes,fgi,fhi] = LSM_62_sym_warp_estimate_f(g,h,vg,vh,Nf,xa,sigma,pt)

 g,h       two images
 vg, vh    variance function
 Nf        half size of reference image
 xa        current geometric and radiometic transformation
 sigma     stdv for smoothing
 pt        in [0,1,2]: controling output

 fe,fes    estimated and smoothed image
 fgi,fhi   interpolated and transformed images g and h

 wf 2018/11

CROSS-REFERENCE INFORMATION ^

This function calls: This function is called by:

SOURCE CODE ^

0001 function [fe,fes,fgi,fhi] = LSM_62_sym_warp_estimate_f(g,h,vg,vh,Nf,xa,sigma,pt)
0002 %% estimate function f from g and h using xa and smooth.
0003 %
0004 % [fe,fes,fgi,fhi] = LSM_62_sym_warp_estimate_f(g,h,vg,vh,Nf,xa,sigma,pt)
0005 %
0006 % g,h       two images
0007 % vg, vh    variance function
0008 % Nf        half size of reference image
0009 % xa        current geometric and radiometic transformation
0010 % sigma     stdv for smoothing
0011 % pt        in [0,1,2]: controling output
0012 %
0013 % fe,fes    estimated and smoothed image
0014 % fgi,fhi   interpolated and transformed images g and h
0015 %
0016 % wf 2018/11
0017 
0018 %% geometric transformation
0019 B  = [xa(1) xa(3) xa(5);xa(2) xa(4) xa(6);0 0 1];
0020 S  = [xa(7);xa(8)];
0021 Bi = inv(B);
0022 
0023 %% sizes
0024 Ng = (size(g,1)-1)/2;
0025 Nh = (size(h,1)-1)/2;
0026 
0027 Mf = 2*Nf+1;
0028 
0029 fe  = zeros(2*Nf+1);
0030 fes = zeros(2*Nf+1);
0031 fgi = zeros(2*Nf+1);
0032 fhi = zeros(2*Nf+1);
0033 
0034 %% estimate f
0035 
0036 
0037 % parallel
0038 % combined transformation f->g
0039 BC = [1 0 Ng+1;0 1 Ng+1; 0 0 1]* Bi*[1 0 -Nf-1;0 1 -Nf-1; 0 0 1];
0040 T = maketform('affine',inv(BC)');
0041 xdata = [1,Mf]; ydata = [1,Mf];
0042 fg_v  = imtransform(g',T,'bicubic','xdata',xdata,'ydata',ydata)';
0043 fg_v  = S(1) * fg_v/255 + S(2);
0044 
0045 % enforce [0,1]
0046 fg_v = max(0,fg_v);
0047 fg_v = min(1,fg_v);
0048 
0049 % weights
0050 wg_v  = 1./vg(round(fg_v(:,:)*255+1))/S(1)^2;
0051 
0052 % combined transformation f->h
0053 BC    = [1 0 Nh+1;0 1 Nh+1; 0 0 1] * B * [1 0 -Nf-1;0 1 -Nf-1; 0 0 1];
0054 T     = maketform('affine',inv(BC)');
0055 xdata = [1,Mf]; ydata = [1,Mf];
0056 fh_v  = imtransform(h',T,'bicubic','xdata',xdata,'ydata',ydata)';
0057 fh_v  = (fh_v/255-S(2))/S(1);
0058 
0059 % enforce [0,1]
0060 fh_v = max(0,fh_v);
0061 fh_v = min(1,fh_v);
0062 
0063 % weights
0064 wh_v  = 1./vh(round(fh_v(:,:)*255+1))*S(1)^2;
0065 
0066 % weighted average
0067 fe_v = (wg_v.*fg_v+wh_v.*fh_v)./(wg_v+wh_v);
0068 
0069 
0070 % % sequential
0071 % for i=-Nf:Nf
0072 %     for j=-Nf:Nf
0073 %
0074 %         % image g: dermine g-values of grid in f
0075 %         xgij = Bi * [i;j;1];
0076 %         xg = xgij(1) + Ng+1;
0077 %         yg = xgij(2) + Ng+1;
0078 %         if xg > 2*Ng || yg > 2*Ng || xg < 2  || yg < 2
0079 %             keyboard
0080 %         end
0081 %         fg = S(1)*LSM_f_cubic_interpolation(xg,yg,g)/255+S(2);
0082 %
0083 %         % enforce range [0,1]
0084 %         fg=min(fg,1);
0085 %         fg=max(fg,0);
0086 %
0087 %         wg = 1/vg(round(fg*255)+1)/S(1)^2;
0088 %         fgi(i+Nf+1,j+Nf+1) = fg;
0089 %
0090 %         % image h: dermine h-values of grid in f
0091 %
0092 %         xhij = B * [i;j;1];
0093 %         xh = xhij(1) + Nh+1;
0094 %         yh = xhij(2) + Nh+1;
0095 %         if xh > 2*Nh || yh > 2*Nh || xh < 2  || yh < 2
0096 %             keyboard
0097 %         end
0098 %         fh = (LSM_f_cubic_interpolation(xh,yh,h)/255-S(2))/S(1);
0099 %
0100 %         % enforce range [0,1]
0101 %         fh=min(fh,1);
0102 %         fh=max(fh,0);
0103 %
0104 %         wh = 1/vh(round(fh*255)+1)*S(1)^2;
0105 %         fhi(i+Nf+1,j+Nf+1) = fh;
0106 %
0107 %         % take weighted average:
0108 %         fe(i+Nf+1,j+Nf+1) = (wg*fg+wh*fh)/(wg+wh);
0109 %     end
0110 % end
0111 
0112 fgi = fg_v;
0113 fhi = fh_v;
0114 fe = fe_v;
0115 
0116 if sigma > 0
0117     % smooth: to be used for derivatives.
0118     fes = gaussFFT(fe,sigma,'G');
0119 else
0120     fes = fe;
0121 end
0122 
0123 if pt > 1
0124     figure(13)
0125     subplot(1,2,1)
0126     hist(fgi(:)-fe(:),2*Nf+1);
0127     title(['mean(g-f) = ',num2str(mean(fgi(:)-fe(:)))])
0128     subplot(1,2,2)
0129     hist(fhi(:)-fe(:),2*Nf+1);
0130     title(['mean(h-f) = ',num2str(mean(fhi(:)-fe(:)))])
0131     Cov_fg=cov(fgi(:),fhi(:));
0132     title('Histogram of residuals of restaured images')
0133     %S1_est=sqrt(tan(0.5*atan2(2*Cov_fg(1,2),Cov_fg(2,2)-Cov_fg(1,1))))-1
0134 end
0135 end
0136

Generated on Sun 19-Jul-2020 23:00:25 by m2html © 2005