[Z,R] = arcgridread('montbays15asc.asc');
 Zsave = Z;
 
 %remove the rows and columns with all NaN
Z( all( isnan( Z ), 2 ), : ) = []; 
Z( :, all( isnan( Z ), 1 ) ) = []; 

%Need to re-map R array

%Build Haxby Color Map
ncolors=11;
c=[ 37    57   175;    40   127   251;    50   190   255;   106   235   255;
    138   236   174;   205   255   162;   240   236   121;   255   189    87;
    255   161    68;   255   186   133;   255   255   255];
pp=1:(m-1)/(ncolors-1):m;
r=interp1(pp,c(:,1),1:m);
g=interp1(pp,c(:,2),1:m);
b=interp1(pp,c(:,3),1:m);

mapshow(Z,R,'DisplayType','surface');
xlabel('x (easting in meters)'); 
ylabel('y (northing in meters)')
demcmap(Z)

% View the terrain in 3-D
axis normal; 
view(2); 
axis equal; 
grid on
zlabel('elevation in meters')