%  intdemo2.m  --- Chapter 18 : example M-file illustrating (in 2D)

%

%  Piecewise linear and cubic interpolations vs global LS appxomations.

%

help intdemo2; 

disp('                     Press any key to continue');pause 



xt = 0:0.01:1;  yt = 0:0.01:1 ;  

xplot = xt;  yplot = yt ; %<<=== 1D vector

[xt yt]=meshgrid(xt,yt) ; %<<::: 2D Points for Plots ::::::::::::::

x=xt;  y=yt ;             %<<==== Set up Plot Points 

 R = sqrt ( (x-1/2).^2 + (y-1/4).^2 ) ;

zt = R.*log(R+eps^2) ;   %% <<==== Main Function for Plots

                         %% True values at lots of points for plot

v1 = min(min(zt)) ;

v2 = max(max(zt)) ;

V = [v1 v1/2 v1/4 (v1+v2)/2 v2/4 v2/2 v2]



x = 0:0.25:1  ;  y = 0:0.25:1 ; % <<<=== Interpolation points given

[X Y]=meshgrid(x,y) ;

x = X;  y = Y ; clear X; clear Y;

 R = sqrt ( (x-1/2).^2 + (y-1/4).^2 ) ;

 z = R.*log(R+eps) ;



%_____________% Global

% First  test %        % Test points xi yi for error tests

%-------------%



xi = 0.1:.15:0.9 ;  yi = xi ; [n m] = size(xi) ; m = m*n;

[u v]=meshgrid(xi,yi) ;

xi = u; yi = v;    R = sqrt ( (xi-1/2).^2 + (yi-1/4).^2 );

yi_t = R.*log(R+eps) ;  % True values for error checks

 

%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%

% Part I  - Least squares fitting %

%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%

  CL = polyfit2(x,y,z, 'linear') ;

  ZL = polyval2(CL, xi,yi, 'linear');

  CC = polyfit2(x,y,z, 'cubic') ;

  ZC = polyval2(CC, xi,yi, 'cubic');

Error_L = sqrt( norm( ZL - yi_t, 'fro')^2 / m)  ; clear ZL;

Error_C = sqrt( norm( ZC - yi_t, 'fro')^2 / m)  ; clear ZC;



%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%

% Part II - Piecewise fitting %

%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%

li = interp2(x,y,z, xi,yi);

ci = interp2(x,y,z, xi,yi, 'cubic'); clear xi; clear yi;

error_L = sqrt( norm( li - yi_t, 'fro')^2 / m)  ;

error_C = sqrt( norm( ci - yi_t, 'fro')^2 / m)  ;





%_____________%         Piecewise 

% Second test %         % Re-select Test points xi yi for error tests

%-------------%



xi = 0.25:.05:0.65 ;  yi = 0.1:0.05:0.4;  [n m] = size(xi) ; m = m*n;

[u v]=meshgrid(xi,yi) ;

xi = u; yi = v;  R = sqrt ( (xi-1/2).^2 + (yi-1/4).^2 ) ;

yi_t = R.*log(R+eps);  % True values for error checks



li = interp2(x,y,z, xi,yi);

ci = interp2(x,y,z, xi,yi, 'cubic'); clear xi; clear yi;

error_L2 = sqrt( norm( li - yi_t, 'fro')^2 / m)  ;

error_C2 = sqrt( norm( ci - yi_t, 'fro')^2 / m)  ;



lt = interp2(x,y,z, xt,yt) ;

ct = interp2(x,y,z, xt,yt, 'cubic') ;



%________________________%

% Plot at lots of points %

%------------------------%

figure(1);clg

         contour(zt,V, xplot,yplot,'-y'); hold on;

         contour(lt,V,xplot,yplot, ':r');

    c2 = contour(ct,V,xplot,yplot, ':g');   %% clabel(c2) ;

ylabel('Function of SIN(x^2+y^3)');

xlabel('File intdemo2.m - example 2D (yellow : true)');

title('2D interpolations - Red Bi-linear vs Green Bi-Cubic');



figure(2);clg, meshz(lt); title('Viewing 3D graph - meshz');  

figure(3);clg, surf(ct);  title('Viewing 3D graph - surf'); 



disp('May try waterfall at this stage')



disp('%%%%%%%%%%%% RMS Error Comparisons %%%%%%%%%%%')

fprintf('AT set 1 of test points, errors are \n')

fprintf('   Linear least squares = %9.3f \n',Error_L)

fprintf('   Cubic  least squares = %9.3f \n',Error_C)

fprintf('   Linear interpolation = %9.3f (smallest = best)\n',error_L)

fprintf('   Cubic  interpolation = %9.3f \n',error_C)

fprintf('AT set 2 of test points, errors are \n')

fprintf('   Linear interpolation = %9.3f \n',error_L2)

fprintf('   Cubic  interpolation = %9.3f \n',error_C2)

disp('........................ End of intdemo2.m...............')



