load keston.out

n=10;
x=keston;
cast=x(:,1);
yd0=x(:,2)-1;
lat=x(:,3);
lon=x(:,4);
z=-x(:,6);
temperature=x(:,7);
salinity=x(:,8);

[x,y]=ll2xy(lat,-lon);
cc=unique(cast);

YD=[];X=[];Y=[];Z=[];T=[];S=[];

for k=1:length(cc);
    j=find(cast==cc(k));
    zi=min(z(j))+[0:n]*(max(z(j)-min(z(j))))/n;
    ti=interp1(z(j),temperature(j),zi);
    si=interp1(z(j),salinity(j),zi);

    YD=[YD;yd0(j(1))+0*zi(:)];
    X=[X;x(j(1))+0*zi(:)];
    Y=[Y;y(j(1))+0*zi(:)];
    Z=[Z;zi(:)];
    T=[T;ti(:)];
    S=[S;si(:)];
end

    
    
fp=fopen('CTD.T.mir3d','w');
for k=1:length(YD)
    fprintf(fp,'%f %f %f %f %f %f %i\n',[YD(k),X(k),Y(k),Z(k),T(k),1.,1]);
end
fclose(fp);
    
fp=fopen('CTD.S.mir3d','w');
for k=1:length(YD)
    fprintf(fp,'%f %f %f %f %f %f %i\n',[YD(k),X(k),Y(k),Z(k),S(k),0.25,1]);
end
fclose(fp);
     