path(path,'/home/ruoying/Software/ecomviz/')
path(path,'/home/ruoying/ROMS/matlab/tools/')
path(path,'/home/ruoying/Software/opnml_matlab5/FCAST_1.2/sinegrid.m')

run /home/ruoying/Software/opnml_matlab5/startup

load /home/ruoying/MERHAB/CTD/2003oc391.mat

t=cruise.temp;
s=cruise.sal;
z=cruise.z;
x=cruise.lon;
y=cruise.lat;
jd=cruise.time;

%== if only use Leg 1 in-situ data inx=83
inx=257;
t=t(:,1:inx); s=s(:,1:inx);

x=x(1:inx); y=y(1:inx);
jd=jd(1:inx,:); jd=julian(jd)-julian([2003 1 1 0 0 0]);

%== Form tn and sn by sigma interpolation in the vertical
for k=1:inx
 tt=t(:,k); ss=s(:,k); zz=(find(~isnan(tt)));
 tt=t(zz,k);ss=s(zz,k);

 z0=sinegrid(max(zz),0,11,1);
 tn(:,k)=transpose(interp1(zz,tt,-z0));
 sn(:,k)=transpose(interp1(zz,ss,-z0));
 zn(:,k)=z0'; 
end

x=ones(11,1)*x;
y=ones(11,1)*y;
jd=ones(11,1)*jd';

%==
ii=find(~isnan(tn) & ~isnan(sn));

tn=tn(ii); sn=sn(ii); x=x(ii); y=y(ii); zn=zn(ii); jd=jd(ii);
x=-1*x;  [x,y]=ll2xy(y,x);

%== output .mir3d file
fm='%9.4f %9.4f %9.4f %9.4f %9.4f %9.4f %4i \n';
fid1=fopen('merhab.t.mir3d_sub','w');
fid2=fopen('merhab.s.mir3d_sub','w');

fprintf(fid1,'%s \n','Merhab  temperature obs');
fprintf(fid1,'%s \n','2003');
fprintf(fid2,'%s \n','Merhab  salinity obs');
fprintf(fid2,'%s \n','2003');

for k=1:length(ii)
 k
 fprintf(fid1,fm,jd(k),x(k),y(k),zn(k),tn(k),1,1);
 fprintf(fid2,fm,jd(k),x(k),y(k),zn(k),sn(k),1,1);
end
fclose(fid1)
fclose(fid2)

