% bathymetry
clear all, close all
nx=510;
LATC=readbin('LATC.bin',[nx 6 nx],1,'real*8');
LONC=readbin('LONC.bin',[nx 6 nx],1,'real*8');
ix=find(LONC<0); LONC(ix)=LONC(ix)+360;
load tbase12
BATH=LATC;
for i=1:nx, mydisp(i)
 for j=1:nx
  for k=1:6
   BATH(i,k,j)=map(closest(LATC(i,k,j),lat),closest(LONC(i,k,j),lon));
  end
 end
end
writebin('BATHY_510x6x510.bin',BATH,1)
clear BATH map

fix_bathy


%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%


% levitus_01 January temperature and salinity
% input files are from nemo:/data2/odap/levitus_01
clear all, close all
grid; DPT=dpt; nx=510;
LATC=readbin('LATC.bin',[nx 6 nx],1,'real*8');
LONC=readbin('LONC.bin',[nx 6 nx],1,'real*8');
ix=find(LONC<0); LONC(ix)=LONC(ix)+360;
lon=.5:359.5; lat=-89.5:89.5;
dpt=[0 10 20 30 50 75 100 125 150 200 250 300 400 500 ...
     600 700 800 900 1000 1100 1200 1300 1400 1500 ...
     1750 2000 2500 3000 3500 4000 4500 5000 5500];
t=readbin('temp01',[360 180 33],1);
s=readbin('sal01',[360 180 33],1);
in=find(t<-50);

% convert to potential
[Y X]=meshgrid(lat,lon);
for k=1:33
 p=pressure(dpt(k),Y);
 t(:,:,k)=insitutemp(s(:,:,k),t(:,:,k),p);
end
t(in)=nan; s(in)=nan;

% fill-in nans
for k=1:33, mydisp(k)
 t(:,:,k)=xpolate(t(:,:,k));
 s(:,:,k)=xpolate(s(:,:,k));
end
save LEV01_JAN_TS lat lon dpt t s



%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%

% monthly Levitus salinity climatology
clear all, close all
pn='/host/nemo/data2/odap/levitus_01/SALINITY/sal';
lon=.5:359.5; lat=-89.5:89.5;
dpt=[0 10 20 30 50 75 100 125 150 200 250 300 400 500 ...
     600 700 800 900 1000 1100 1200 1300 1400 1500 ...
     1750 2000 2500 3000 3500 4000 4500 5000 5500];
nx=length(lon); ny=length(lat); nz=length(dpt);
salt=zeros(nx,ny,12);
for m=1:12
  s=readbin([pn myint2str(m)],[nx ny nz]);
  tmp=(s(:,:,1)+s(:,:,2))/2;
  tmp(find(tmp<0))=nan;
  salt(:,:,m)=xpolate(tmp);
  mypcolor(lon,lat,salt(:,:,m)'); colorbar; pause(.01)
end
writebin('LEV01_SSS_360x180x12.bin',salt)


%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%

% interpolate to cube sphere
clear all, close all
grid; DPT=dpt; nx=510;
LATC=readbin('LATC.bin',[nx 6 nx],1,'real*8');
LONC=readbin('LONC.bin',[nx 6 nx],1,'real*8');
ix=find(LONC<0); LONC(ix)=LONC(ix)+360;
load LEV01_JAN_TS
t2=zeros(362,182,34);
t2(2:361,2:181,1:33)=t; t2(:,:,34) =t2(:,:,33);
t2(1,:,:)=t2(361,:,:);  t2(362,:,:)=t2(2,:,:);
t2(:,1,:)=t2(:,2,:);    t2(:,182,:)=t2(:,181,:);
s2=zeros(362,182,34);
s2(2:361,2:181,1:33)=s; s2(:,:,34) =s2(:,:,33);
s2(1,:,:)=s2(361,:,:);  s2(362,:,:)=s2(2,:,:);
s2(:,1,:)=s2(:,2,:);    s2(:,182,:)=s2(:,181,:);
lon2=-.5:360.5; lat2=-90.5:90.5; dpt=[dpt max(DPT)];
[Y X]=meshgrid(lat2,lon2);
LEVT=zeros(nx,6,nx,50); LEVS=LEVT;
for k=1:50
 ix=closest(DPT(k),dpt,0);
 tmp=(t2(:,:,ix(1))*abs(DPT(k)-dpt(ix(2)))+...
      t2(:,:,ix(2))*abs(DPT(k)-dpt(ix(1))))/...
     abs(dpt(ix(2))-dpt(ix(1)));
 LEVT(:,:,:,k)=interp2(Y,X,tmp,LATC,LONC);
 tmp=(s2(:,:,ix(1))*abs(DPT(k)-dpt(ix(2)))+...
      s2(:,:,ix(2))*abs(DPT(k)-dpt(ix(1))))/...
     abs(dpt(ix(2))-dpt(ix(1)));
 LEVS(:,:,:,k)=interp2(Y,X,tmp,LATC,LONC);
end
writebin('LEVT01_JAN_510x6x510.bin',LEVT,1)
writebin('LEVS01_JAN_510x6x510.bin',LEVS,1)


%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%

% convert to proper format for model input
clear all, close all
tx=85;ty=85;nt=216;cx=510;cy=510;

te=readbin('BATHY_510x6x510_filled.bin',[510 6 510],1);
te=permute(te,[1 3 2]);
nz=1; anti_foo
writebin('BATHY_18360x85_filled.bin',phi,1)

te=readbin('LEVT01_JAN_510x6x510.bin',[510 6 510 50],1);
te=permute(te,[1 3 2 4]);
nz=50; anti_foo
writebin('LEVT01_JAN_18360x85x50.bin',phi,1)

te=readbin('LEVS01_JAN_510x6x510.bin',[510 6 510 50],1);
te=permute(te,[1 3 2 4]);
nz=50; anti_foo
writebin('LEVS01_JAN_18360x85x50.bin',phi,1)



%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%

% convert to proper format for model input
clear all, close all
tx=85;ty=85;nt=216;cx=510;cy=510;

te=readbin('BATHY_510x6x510_filled.bin',[510 6 510],1);
te=permute(te,[1 3 2]);
nz=1; anti_foo
writebin('BATHY_18360x85_filled.bin',phi,1)

te=readbin('LEVT01_JAN_510x6x510.bin',[510 6 510 50],1);
te=permute(te,[1 3 2 4]);
nz=50; anti_foo
writebin('LEVT01_JAN_18360x85x50.bin',phi,1)

te=readbin('LEVS01_JAN_510x6x510.bin',[510 6 510 50],1);
te=permute(te,[1 3 2 4]);
nz=50; anti_foo
writebin('LEVS01_JAN_18360x85x50.bin',phi,1)


%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%

% convert to proper format for model input
clear all, close all
tx=51;ty=510;nt=60;cx=510;cy=510;

te=readbin('BATHY_510x6x510_filled.bin',[510 6 510],1);
te=permute(te,[1 3 2]);
nz=1; anti_foo
writebin('BATHY_3060x510_filled.bin',phi,1)

te=readbin('LEVT01_JAN_510x6x510.bin',[510 6 510 50],1);
te=permute(te,[1 3 2 4]);
nz=50; anti_foo
writebin('LEVT01_JAN_3060x510x50.bin',phi,1)

te=readbin('LEVS01_JAN_510x6x510.bin',[510 6 510 50],1);
te=permute(te,[1 3 2 4]);
nz=50; anti_foo
writebin('LEVS01_JAN_3060x510x50.bin',phi,1)