1 |
dimitri |
1.1 |
% function create_bc_grid(dir_in, nx, ny, nz, cbdry) |
2 |
|
|
% |
3 |
|
|
% Input: dir_in, directory that has MITgcm grid files |
4 |
|
|
% (following files are expected, |
5 |
|
|
% YC.data. XC.data, hFacC.data |
6 |
|
|
% RF.data, AngleCS.data, AngleSN.data ) |
7 |
|
|
% nx, first dimension of MITgcm grid, (west-east) |
8 |
|
|
% ny, second dimension of MITgcm grid (south-north) |
9 |
|
|
% nz, dimension in z (or r) direction |
10 |
|
|
% cbdry, character for bdry to generate: s, n, w, e, (south, north, west, east) |
11 |
|
|
% |
12 |
|
|
% Output: No output |
13 |
|
|
% The script will generate a grid file named |
14 |
|
|
% bdry_grid_[snwe].mat |
15 |
|
|
% |
16 |
|
|
% |
17 |
|
|
% XC WANG 25jan2013 |
18 |
|
|
% Modified 01feb2013 |
19 |
|
|
% |
20 |
|
|
% |
21 |
|
|
function create_bc_grid(dir_in,nx,ny,nz,cbdry) |
22 |
|
|
|
23 |
|
|
% addpath /data17/home/xiao/Mat_Tools/MITgcm_Tool |
24 |
|
|
|
25 |
|
|
% Default location of boundary |
26 |
|
|
ny_north = ny -1 ; ny_south=1; |
27 |
|
|
nx_west = 1; nx_east=nx-1 ; |
28 |
|
|
|
29 |
|
|
%dir_in=['/data17/XCW_output/Grid/data/']; |
30 |
|
|
% Three boundaries, South, North and West |
31 |
|
|
%ibc = 2; |
32 |
|
|
%c_bdry3 = {'s', 'n', 'w', 'e'} ; |
33 |
|
|
%c_bdry = c_bdry3{ibc} ; |
34 |
|
|
%%bdry_file=['bdry_grid_',c_bdry, '.mat.test']; |
35 |
|
|
%bdry_file=['bdry_grid_',c_bdry, '.mat']; |
36 |
|
|
|
37 |
|
|
ufile = ['YC.data']; |
38 |
|
|
ufile1 = [dir_in, ufile]; |
39 |
|
|
lat = readbin(ufile1, [nx ny]) ; |
40 |
|
|
|
41 |
|
|
ufile = ['XC.data']; |
42 |
|
|
ufile1 = [dir_in, ufile]; |
43 |
|
|
lon = readbin(ufile1, [nx ny]) ; |
44 |
|
|
|
45 |
|
|
%ufile = ['Depth.data']; |
46 |
|
|
%ufile1 = [dir_in, ufile]; |
47 |
|
|
%depth = readbin(ufile1, [nx ny]) ; |
48 |
|
|
|
49 |
|
|
rf=-readbin([dir_in, 'RF.data'], nz+1) ; |
50 |
|
|
thk=diff(rf) ; |
51 |
|
|
|
52 |
|
|
hfac =readbin([dir_in, 'hFacC.data'], [nx ny nz]) ; |
53 |
|
|
%obmask = squeeze(hfac(:, ny_north, :)) ; |
54 |
|
|
thk3d = zeros(size(hfac)) ; |
55 |
|
|
for k =1:nz |
56 |
|
|
temp1 = squeeze(hfac(:,:, k)) ; |
57 |
|
|
thk3d(:,:, k) = thk(k) * temp1 ; |
58 |
|
|
end |
59 |
|
|
depth = sum(thk3d, 3) ; |
60 |
|
|
%depth_bc = depth2d(:, ny_north) ; |
61 |
|
|
|
62 |
|
|
anglefile = ['AngleCS.data'] ; |
63 |
|
|
tfile1 = [dir_in, anglefile]; |
64 |
|
|
AngleCS = readbin(tfile1, [nx ny]) ; |
65 |
|
|
|
66 |
|
|
anglefile = ['AngleSN.data'] ; |
67 |
|
|
tfile1 = [dir_in, anglefile]; |
68 |
|
|
AngleSN = readbin(tfile1, [nx ny]) ; |
69 |
|
|
|
70 |
|
|
%%bdry_file=['bdry_grid_',c_bdry, '.mat.test']; |
71 |
|
|
bdry_file=['bdry_grid_',cbdry, '.mat']; |
72 |
|
|
switch cbdry |
73 |
|
|
case 's' |
74 |
|
|
lat_bc=lat(:, ny_south); |
75 |
|
|
lon_bc=lon(:, ny_south) ; |
76 |
|
|
depth_bc=depth(:, ny_south) ; |
77 |
|
|
AngleCS_bc = AngleCS(:,ny_south); |
78 |
|
|
AngleSN_bc = AngleSN(:,ny_south); |
79 |
|
|
%bdry_file=['bdry_grid_',c_bdry, '.mat']; |
80 |
|
|
case 'n' |
81 |
|
|
lat_bc=lat(:, ny_north); |
82 |
|
|
lon_bc=lon(:, ny_north) ; |
83 |
|
|
depth_bc=depth(:, ny_north) ; |
84 |
|
|
AngleCS_bc = AngleCS(:,ny_north); |
85 |
|
|
AngleSN_bc = AngleSN(:,ny_north); |
86 |
|
|
case 'w' |
87 |
|
|
lat_bc=lat(nx_west,:); |
88 |
|
|
lon_bc=lon(nx_west,:) ; |
89 |
|
|
depth_bc=depth(nx_west,:) ; |
90 |
|
|
AngleCS_bc = AngleCS(nx_west,:); |
91 |
|
|
AngleSN_bc = AngleSN(nx_west,:); |
92 |
|
|
case 'e' |
93 |
|
|
lat_bc=lat(nx_east,:); |
94 |
|
|
lon_bc=lon(nx_east,:) ; |
95 |
|
|
depth_bc=depth(nx_east,:) ; |
96 |
|
|
AngleCS_bc = AngleCS(nx_east,:); |
97 |
|
|
AngleSN_bc = AngleSN(nx_east,:); |
98 |
|
|
otherwise |
99 |
|
|
display('No such selection ', cbdry) |
100 |
|
|
return |
101 |
|
|
end |
102 |
|
|
|
103 |
|
|
%rad2deg = 180/pi; |
104 |
|
|
%angle = atan2(AngleSN, AngleCS); |
105 |
|
|
%angle = angle*rad2deg ; |
106 |
|
|
|
107 |
|
|
save(bdry_file, 'lat_bc', 'lon_bc', 'depth_bc', 'AngleCS_bc', 'AngleSN_bc') ; |
108 |
|
|
|
109 |
|
|
return |