/[MITgcm]/MITgcm/verification/tutorial_global_oce_latlon/diags_matlab/mit_loadgrid.m
ViewVC logotype

Annotation of /MITgcm/verification/tutorial_global_oce_latlon/diags_matlab/mit_loadgrid.m

Parent Directory Parent Directory | Revision Log Revision Log | View Revision Graph Revision Graph


Revision 1.3 - (hide annotations) (download)
Sat Aug 12 20:25:13 2006 UTC (17 years, 9 months ago) by jmc
Branch: MAIN
CVS Tags: checkpoint64y, checkpoint64x, checkpoint64z, checkpoint64q, checkpoint64p, checkpoint64s, checkpoint64r, checkpoint64u, checkpoint64t, checkpoint64w, checkpoint64v, checkpoint64i, checkpoint64h, checkpoint64k, checkpoint64j, checkpoint64m, checkpoint64l, checkpoint64o, checkpoint64n, checkpoint64a, checkpoint64c, checkpoint64b, checkpoint64e, checkpoint64d, checkpoint64g, checkpoint64f, checkpoint58u_post, checkpoint58w_post, checkpoint63p, checkpoint63q, checkpoint63r, checkpoint63s, checkpoint63l, checkpoint63m, checkpoint63n, checkpoint63o, checkpoint63h, checkpoint63i, checkpoint63j, checkpoint63k, checkpoint63d, checkpoint63e, checkpoint63f, checkpoint63g, checkpoint63a, checkpoint63b, checkpoint63c, checkpoint64, checkpoint65, checkpoint60, checkpoint61, checkpoint62, checkpoint63, checkpoint58r_post, checkpoint66g, checkpoint66f, checkpoint66e, checkpoint66d, checkpoint66c, checkpoint66b, checkpoint66a, checkpoint66o, checkpoint66n, checkpoint66m, checkpoint66l, checkpoint66k, checkpoint66j, checkpoint66i, checkpoint66h, checkpoint58n_post, checkpoint58x_post, checkpoint58t_post, checkpoint65z, checkpoint65x, checkpoint65y, checkpoint65r, checkpoint65s, checkpoint65p, checkpoint65q, checkpoint65v, checkpoint65w, checkpoint65t, checkpoint65u, checkpoint65j, checkpoint65k, checkpoint65h, checkpoint65i, checkpoint65n, checkpoint65o, checkpoint65l, checkpoint65m, checkpoint65b, checkpoint65c, checkpoint65a, checkpoint65f, checkpoint65g, checkpoint65d, checkpoint65e, checkpoint58q_post, checkpoint59q, checkpoint59p, checkpoint59r, checkpoint59e, checkpoint59d, checkpoint59g, checkpoint59f, checkpoint59a, checkpoint59c, checkpoint59b, checkpoint59m, checkpoint59l, checkpoint59o, checkpoint59n, checkpoint59i, checkpoint59h, checkpoint59k, checkpoint59j, checkpoint59, checkpoint58o_post, checkpoint62c, checkpoint62b, checkpoint62a, checkpoint62g, checkpoint62f, checkpoint62e, checkpoint62d, checkpoint62k, checkpoint62j, checkpoint62i, checkpoint62h, checkpoint62o, checkpoint62n, checkpoint62m, checkpoint62l, checkpoint62s, checkpoint62r, checkpoint62q, checkpoint62p, checkpoint62w, checkpoint62v, checkpoint62u, checkpoint62t, checkpoint62z, checkpoint62y, checkpoint62x, checkpoint58y_post, checkpoint58v_post, checkpoint58s_post, checkpoint61f, checkpoint61g, checkpoint61d, checkpoint61e, checkpoint61b, checkpoint61c, checkpoint58p_post, checkpoint61a, checkpoint61n, checkpoint61o, checkpoint61l, checkpoint61m, checkpoint61j, checkpoint61k, checkpoint61h, checkpoint61i, checkpoint61v, checkpoint61w, checkpoint61t, checkpoint61u, checkpoint61r, checkpoint61s, checkpoint61p, checkpoint61q, checkpoint61z, checkpoint61x, checkpoint61y, HEAD
Changes since 1.2: +1 -1 lines
accidentally removed ; put them back.

1 jmc 1.1 function grid = mit_loadgrid(varargin);
2     % function grid = mit_loadgrid(dname);
3     % load the geometry of an arbitrary run with the mitgcm
4    
5     % Aug 15, 2002: fixed a bug (?): yc' -> yc
6     %
7 jmc 1.3 % $Header: $
8 jmc 1.1 % $Name: $
9    
10     if nargin == 0
11     dname = '.';
12     else
13     dname = varargin{1};
14     end
15     % tracer time step
16     deltattracer = mit_getparm(fullfile(dname,'data'),'deltaTtracer');
17     if isempty(deltattracer)
18     error('deltaTtracer is empty')
19     end
20     gravity = mit_getparm(fullfile(dname,'data'),'gravity');
21     if isempty(gravity);
22     disp('assuming gravity = 9.81')
23     gravity = 9.81;
24     end
25     rhoNil = mit_getparm(fullfile(dname,'data'),'rhonil');
26     if isempty(rhoNil);
27     rhoNil = mit_getparm(fullfile(dname,'data'),'rhoNil');
28     end
29     if isempty(rhoNil);
30     rhoNil = mit_getparm(fullfile(dname,'data'),'rhoConst');
31     end
32     if isempty(rhoNil);
33     rhoNil = mit_getparm(fullfile(dname,'data'),'rhoconst');
34     end
35     if isempty(rhoNil);
36     disp('assuming rhoNil = 1035')
37     rhoNil = 1035;
38     end
39     % determine vertical coordinates
40     br = mit_getparm(fullfile(dname,'data'),'buoyancyRelation');
41     if isempty(br)
42     disp('assuming buoyancyRelation = OCEANIC')
43     br = 'OCEANIC';
44     end
45     if ~strcmpi(br,'OCEANIC');
46     pfac = 1/rhoNil/gravity;
47     else
48     pfac = 1;
49     end
50     % create masks for plotting
51     hfac=rdmds(fullfile(dname,'hFacC')); hfacc=change(hfac,'==',0,NaN);
52     hfac=rdmds(fullfile(dname,'hFacW')); hfacw=change(hfac,'==',0,NaN);
53     hfac=rdmds(fullfile(dname,'hFacS')); hfacs=change(hfac,'==',0,NaN);
54     clear hfac;
55     cmask=change(hfacc,'>',0,1);
56     umask=change(hfacw,'>',0,1);
57     vmask=change(hfacs,'>',0,1);
58    
59     [nx ny nz] = size(cmask);
60    
61     % find the index of the deepest wet tracer-cell
62     klowc = sum(change(cmask,'==',NaN,0),3);
63    
64     % create grid parameters
65     dz = mit_getparm(fullfile(dname,'data'),'delR');
66     if isempty(dz);
67     disp('trying delZ')
68     dz = mit_getparm(fullfile(dname,'data'),'delZ');
69     if isempty(dz)
70     error('vertical grid cannot be established')
71     end
72     end
73     if size(dz,1) == 1;
74     dz = dz';
75     end
76     if length(dz) ~= nz
77     error('dz could not be retrieved correctly from the data file')
78     end
79     dz = dz*pfac;
80    
81     if ~strcmp(br,'OCEANIC');
82     dz = dz(end:-1:1);
83     end
84     zgpsi = [0;cumsum(dz)];
85     zc = .5*(zgpsi(1:nz)+zgpsi(2:nz+1));
86     zg=zgpsi(1:nz);
87    
88     xc = rdmds(fullfile(dname,'XC'));
89     yc = rdmds(fullfile(dname,'YC'));
90     xg = rdmds(fullfile(dname,'XG'));
91     yg = rdmds(fullfile(dname,'YG'));
92     dxc = rdmds(fullfile(dname,'DXC'));
93     dyc = rdmds(fullfile(dname,'DYC'));
94     dxg = rdmds(fullfile(dname,'DXG'));
95     dyg = rdmds(fullfile(dname,'DYG'));
96    
97     lonc = xc(:,1);
98     latc = yc(1,:)';
99     long = xg(:,1);
100     latg = yg(1,:)';
101    
102     % depth
103     depth = rdmds(fullfile(dname,'Depth'))*pfac;
104    
105     % current directory
106     if strcmp(dname,'.') | strcmp(dname,'./')
107     [pathstr,dirname,ext]=fileparts(pwd);
108     dirname = [dirname ext];
109     else
110     [pathstr,dirname,ext]=fileparts(fullfile(pwd,dname));
111     dirname = [dirname ext];
112     end
113    
114     ius = findstr(dirname,'_');
115     if ~isempty(ius)
116     for k=1:length(ius)
117     dirname = [dirname(1:ius(k)-1) '\' dirname(ius(k):end)];
118     ius = ius+1;
119     end
120     end
121     %
122     if ~strcmp(br,'OCEANIC');
123     hfacc = hfacc(:,:,end:-1:1);
124     hfacs = hfacs(:,:,end:-1:1);
125     hfacw = hfacw(:,:,end:-1:1);
126     cmask = cmask(:,:,end:-1:1);
127     umask = umask(:,:,end:-1:1);
128     vmask = vmask(:,:,end:-1:1);
129     klowc = nz - klowc + 1;
130     klowc(find(klowc==nz+1)) = 0;
131     end
132    
133     % area and volume of grid cells has to be done after flipping the masks
134     % vertically
135     rac = rdmds(fullfile(dname,'RAC')).*cmask(:,:,1);
136     ras = rdmds(fullfile(dname,'RAS')).*vmask(:,:,1);
137     raw = rdmds(fullfile(dname,'RAW')).*umask(:,:,1);
138     rac3d = repmat(rac,[1 1 nz]).*hfacc;
139     volc = rac3d.*permute(repmat(dz,[1 nx ny]),[2 3 1]);
140     % create the structure
141     % Aug 15, 2002: fixed a bug (?): yc' -> yc
142     grid = struct('dname',dirname, ...
143     'deltattracer',deltattracer, ...
144     'gravity',gravity','rhonil',rhoNil, ...
145     'buoyancy',br,'pfac',pfac, ...
146     'nx',nx,'ny',ny,'nz',nz, ...
147     'xc',xc,'yc',yc,'xg',xg,'yg',yg, ...
148     'lonc',lonc,'latc',latc,'long',long,'latg',latg, ...
149     'dxc',dxc,'dyc',dyc,'dxg',dxg,'dyg',dyg, ...
150     'rac',rac,'ras',ras,'raw',raw, ...
151     'dz',dz,'zc',zc,'zg',zg,'zgpsi',zgpsi, ...
152     'depth',depth,'hfacc',hfacc,'hfacs',hfacs,'hfacw',hfacw, ...
153     'cmask',cmask,'umask',umask,'vmask',vmask,'volc',volc, ...
154     'klowc',klowc);
155    
156     return

  ViewVC Help
Powered by ViewVC 1.1.22