

%lon lat at shell point USGS sensor
lon_sp=-(82 + 0 + 15.95./3600);
lat_sp=26 + 31./60 + 26.6./3600;

k=16;  %sigma coord starts at bottom


cd ../uncompressed

time=[];

u=[];
v=[];

for f=10:34,
fn=['tmpnc' num2str(f) '.nc'];

% nc=netcdf(fn);  %opens file 
% time=nc{'ocean_time'}(:);
% lon=nc{'lon_u'}(:,:);
% lat=nc{'lat_u'}(:,:);

if f==10,
lon=nc_varget(fn,'lon_u');
lat=nc_varget(fn,'lat_u');

[i,j,dist_km]=near2d(lon,lat,lon_sp,lat_sp)

end

time=[time;nc_varget(fn,'ocean_time')];

%load velocity

u_=nc_varget(fn,'ubar');

v_=nc_varget(fn,'vbar');


u=[u ;u_(:,i,j)]; %j)  %,i);
v=[v ;v_(:,i,j)]; % j),i);

end

%close(nc)
cd ../analout

jd=datenum(1968,5,23)+time./(3600.*24);



figure
plot(jd,u,'.-',jd,v,'.-')
legend('u','v')
datetick
dat=load('shell001.vel');

u_all=dat(:,2:4:38);

uusgs=mean(u_all(:,1:5)')';


v_all=dat(:,3:4:39);

vusgs=mean(v_all(:,1:5)')';


%positive is ebb direction

dathd=load('shell001.dat');

jd_gmt=datenum(dathd(:,1),dathd(:,2),dathd(:,3),dathd(:,4),dathd(:,5),dat(:,6));


%inclusive days

 %  733611.7186799189    733616.9994107370  

inusgs=find(jd_gmt>=733611.7186799189 & jd_gmt<=733616.9994107370 );
inmod=find(jd>=733611.7186799189 & jd<=733616.9994107370 );

% figure
% subplot(211)
% plot(jd,u,'.-',jd,v,'.-')
% grid
% legend('u','v')
% title('model')
% datetick('x','keeplimits','keepticks')
% 
% subplot(212)
% 
% plot(jd_gmt(inmod),uusgs(inmod),...
%     jd_gmt(inmod),vusgs(inmod))
% legend('u','v')
% title('USGS')
% datetick('x','keeplimits','keepticks')


figure
subplot(211)
plot(jd,v,'.-',jd_gmt(inmod),-uusgs(inmod)./100)
grid
legend('model','usgs')
title('model')
datetick('x','keeplimits','keepticks')


figure
subplot(211)
plot(jd(inmod),v(inmod),'.-',jd_gmt(inusgs),sign(-uusgs(inusgs)).*sqrt...
    (uusgs(inusgs).^2+vusgs(inusgs).^2)./100)
grid
legend('model','usgs')
title('model')
datetick('x','keeplimits','keepticks')

