%Queen Creek lower riser topo data and plottong

load gone.txt
load gtwo.txt
load gthree.txt
load gfour.txt

%refer to columns as:
%data(:,1) is point number
%data(:,2) is Easting
%data(:,3) is Northing
%data(:,4) is Elevation

%set parameters:
%number of x and y nodes for grid
%make these smaller if you get an error about excceding the max array size
nx = 50;
ny = nx;
%contour interval:
ci = 0.5;

%plot data group by group
figure(1)
clf
%group 1
subplot(2,2,1)
data = gone;
%Here is contour map:
e=data(:,2);
n=data(:,3);
h=data(:,4);
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,h,ee,nn');

v = floor(min(h)):ci:max(h);

[c,h] = contour(xi,yi,zi,v);

%label contours
%clabel(c,h);
colormap([0,0,0])
hold on

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

xlabel('easting (m)')    
ylabel('northing (m)')   
title('Group 1')
axis equal
axis([1000 1080 920 990])


%group 2
subplot(2,2,2)
data = gtwo;
%Here is contour map:
e=data(:,2);
n=data(:,3);
h=data(:,4);
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,h,ee,nn');

v = floor(min(h)):ci:max(h);

[c,h] = contour(xi,yi,zi,v);

%label contours
%clabel(c,h);
colormap([0,0,0])
hold on

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

xlabel('easting (m)')    
ylabel('northing (m)')   
title('Group 2')
axis equal
axis([1000 1080 920 990])


%group 3
subplot(2,2,3)
data = gthree;
%Here is contour map:
e=data(:,2);
n=data(:,3);
h=data(:,4);
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,h,ee,nn');

v = floor(min(h)):ci:max(h);

[c,h] = contour(xi,yi,zi,v);

%label contours
%clabel(c,h);
colormap([0,0,0])
hold on

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

xlabel('easting (m)')    
ylabel('northing (m)')   
title('Group 3')
axis equal
axis([1000 1080 920 990])


%group 4
subplot(2,2,4)
data = gfour;
%Here is contour map:
e=data(:,2);
n=data(:,3);
h=data(:,4);
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,h,ee,nn');

v = floor(min(h)):ci:max(h);

[c,h] = contour(xi,yi,zi,v);

%label contours
%clabel(c,h);
colormap([0,0,0])
hold on

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

xlabel('easting (m)')    
ylabel('northing (m)')   
title('Group 4')
axis equal
axis([1000 1080 920 990])

print -depsc grps.eps


%Now lump them all together
figure(2)
clf
data = [gone; gtwo;gthree;gfour];

%Here is contour map:
e=data(:,2);
n=data(:,3);
h=data(:,4);
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,h,ee,nn');

v = floor(min(h)):ci:max(h);

[c,h] = contour(xi,yi,zi,v);

%label contours
clabel(c,h, 'manual');
colormap([0,0,0])
hold on

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

xlabel('easting (m)')    
ylabel('northing (m)')   
title('All group data')
axis equal
axis([1000 1080 920 990])

%plot rebar control points  does not seem right
%rebar=[gtwo(length(gtwo),:);gtwo(length(gtwo)-1,:)];
%plot(rebar(:,2),rebar(:,3),'ko')

print -depsc all.eps

%oblique view  no mesh
figure(3)
clf
surfl(xi,yi,zi)
colormap(pink)
view(-78,46)
shading interp
xlabel('easting (m)')    
ylabel('northing (m)') 
zlabel('elevation (m)')
title('All group data')


%oblique view 
figure(4)
clf
surfl(xi,yi,zi)
colormap(pink)
view(-78,46)
xlabel('easting (m)')    
ylabel('northing (m)') 
zlabel('elevation (m)')
title('All group data')