| 1 |
heimbach |
1.1 |
clear all, clf reset, cd /home/dimitri/mitgcm/setups/cap1deg |
| 2 |
|
|
nx=360; ny=360; |
| 3 |
|
|
yc=readbin('YC.bin',[nx ny]); xc=readbin('XC.bin',[nx ny]); |
| 4 |
|
|
yg=readbin('YG.bin',[nx+1 ny+1]); xg=readbin('XG.bin',[nx+1 ny+1]); |
| 5 |
|
|
[etopo2,lon,lat]=get_etopo2; |
| 6 |
|
|
|
| 7 |
|
|
% find median in each grid box |
| 8 |
|
|
bathy=0*xc; |
| 9 |
|
|
for i=1:nx, mydisp(xc(i,1)) |
| 10 |
|
|
ix=find(lon>=xg(i,1)&lon<=xg(i+1,1)); |
| 11 |
|
|
for j=1:ny |
| 12 |
|
|
iy=find(lat>=yg(i,j)&lat<=yg(i,j+1)); |
| 13 |
|
|
tmp=etopo2(ix,iy); in=find(tmp<0); |
| 14 |
|
|
if ~isempty(tmp) |
| 15 |
|
|
bathy(i,j)=median(tmp(:)); |
| 16 |
|
|
end, end, end |
| 17 |
|
|
|
| 18 |
|
|
% do not allow any numbers between 0 and 10 |
| 19 |
|
|
topog=bathy; topog(find(topog>-1.5))=0; |
| 20 |
|
|
topog(find(topog<0&topog>-10))=-10; |
| 21 |
|
|
|
| 22 |
|
|
% open a few shallow straits, just for fun |
| 23 |
|
|
topog(1,278)=-10; topog(12,285)=-10; topog(27:28,264)=-10; |
| 24 |
|
|
topog(43:44,226)=-10; topog(274:275,306)=-10; topog(276:278,316)=-10; |
| 25 |
|
|
topog(286,311)=-10; topog(285,312)=-10; topog(252,312)=-10; |
| 26 |
|
|
topog(253,313:314)=-10; topog(255,313)=-10; topog(74,320)=-10; |
| 27 |
|
|
topog(111,329)=-10; topog(355,258)=-284; |
| 28 |
|
|
|
| 29 |
|
|
% fill-in lakes |
| 30 |
|
|
tmp=logical(~topog); tmp=imfill(tmp,'holes'); |
| 31 |
|
|
topog(find(tmp==1))=0; topog(60:61,360)=0; topog(270:282,356:360)=0; |
| 32 |
|
|
|
| 33 |
|
|
% look at and save |
| 34 |
|
|
tmp=-topog; tmp(find(~tmp))=1; clf, mypcolor(log10(tmp)'), thincolorbar |
| 35 |
|
|
writebin('ETOPO2_360x360',topog); |