function contourmap (data, xblocks, yblocks, interpmethod)
%
%contourmap
%
%This will assign contours to the data
%EXAMPLE:  Here is contour map:
%figure(2)
%clf
%estep = abs((min(e)-max(e))/nx);
%nstep = abs((min(n)-max(n))/ny);
%ee = [min(e):estep:max(e)];
%nn = [min(n):nstep:max(n)];
%[xi,yi,zi] = griddata(e,n,hh,ee,nn');
%v = floor(min(hh)):ci:max(hh);
%[c,h] = contour(xi,yi,zi,v);
%label contours
%clabel(c,h);
%colormap(jet)
%hold on
%show data:
%plot(e,n,'o')
%
%show grid:
%plot(xi,yi,'+')
%
%xlabel('easting')    
%ylabel('northing')   
%title('Grid and observation points')
%axis equal
%
%print -depsc vmptsandgrid.eps
%
%Wendy Bohon, March 10, 2009

%this will almost always be linear, but putting it as an option allows the
%user to change it if they want to
interpmethod = 'linear'


%This plots the data and gives an unevenly spaced grid.
figure(2)
clf

%this brings up "help" and gives lots of instructions on griddata
help griddata

%tells griddata what the x and y grid looks like on which 
%it will operate-it's the easting/northing stepsize.
%It needs to be positive, so if it ends up being negative 
%you need to take the absolute value, which is the "abs" command
numsteps=100
estep = abs((min(e)-max(e))/numsteps);
nstep = abs((min(n)-max(n))/numsteps);
ee = [min(e):estep:max(e)];
nn = [min(n):nstep:max(n)];

%This grids the data.
%the 'n means transpose N (switch it from a row vector to a column vector)
% This was orginial  [xi,yi,zi] = griddata(e,n,hh,ee,nn');
[xi,yi,zi] = griddata(e,n,h,ee,nn');

%this was orginial [c,h] = contour(xi,yi,zi,v);
[c,h] = contour(xi,yi,zi);

%this labels the contours
%label contours
clabel(c,h);
colormap(jet)

hold on

%show data:
plot(e,n,'.')


xlabel('easting')    
ylabel('northing')   
title('Contour Map')
axis equal

print('contour.jpg','-djpeg');

end