/[MITgcm]/MITgcm_contrib/gael/matlab_class/gcmfaces_IO/prep2nctiles.m
ViewVC logotype

Contents of /MITgcm_contrib/gael/matlab_class/gcmfaces_IO/prep2nctiles.m

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


Revision 1.7 - (show annotations) (download)
Wed Feb 4 14:08:34 2015 UTC (10 years, 5 months ago) by gforget
Branch: MAIN
CVS Tags: checkpoint65r, checkpoint65p, checkpoint65q, checkpoint65t
Changes since 1.6: +96 -2 lines
- treat the case of initial and final snapshots

1
2 choiceStruct=2;
3
4 %===============
5
6 if choiceStruct==1;
7 structIn=[];
8 %
9 structIn.vars=mygrid;
10 structIn.vars.RF=structIn.vars.RF(1:50);
11 %
12 structIn.descr={'C-grid parameters (see MITgcm documentation for details).'};
13 %
14 vars=[]; nv=length(vars)+1;
15 vars(nv).fldName='XC'; vars(nv).longName='longitude'; vars(nv).units='degrees_east'; nv=length(vars)+1;
16 vars(nv).fldName='YC'; vars(nv).longName='latitude'; vars(nv).units='degrees_north'; nv=length(vars)+1;
17 vars(nv).fldName='XG'; vars(nv).longName='longitude'; vars(nv).units='degrees_east'; nv=length(vars)+1;
18 vars(nv).fldName='YG'; vars(nv).longName='latitude'; vars(nv).units='degrees_north'; nv=length(vars)+1;
19 vars(nv).fldName='RAC'; vars(nv).longName='grid cell area'; vars(nv).units='m^2'; nv=length(vars)+1;
20 vars(nv).fldName='RAZ'; vars(nv).longName='grid cell area'; vars(nv).units='m^2'; nv=length(vars)+1;
21 vars(nv).fldName='DXC'; vars(nv).longName='grid spacing'; vars(nv).units='m'; nv=length(vars)+1;
22 vars(nv).fldName='DYC'; vars(nv).longName='grid spacing'; vars(nv).units='m'; nv=length(vars)+1;
23 vars(nv).fldName='DXG'; vars(nv).longName='grid spacing'; vars(nv).units='m'; nv=length(vars)+1;
24 vars(nv).fldName='DYG'; vars(nv).longName='grid spacing'; vars(nv).units='m'; nv=length(vars)+1;
25 vars(nv).fldName='hFacC'; vars(nv).longName='fractional thickness'; vars(nv).units='1'; nv=length(vars)+1;
26 vars(nv).fldName='hFacW'; vars(nv).longName='fractional thickness'; vars(nv).units='1'; nv=length(vars)+1;
27 vars(nv).fldName='hFacS'; vars(nv).longName='fractional thickness'; vars(nv).units='1'; nv=length(vars)+1;
28 vars(nv).fldName='Depth'; vars(nv).longName='sea floor depth'; vars(nv).units='m'; nv=length(vars)+1;
29 vars(nv).fldName='AngleCS'; vars(nv).longName='grid orientation (cosine)'; vars(nv).units='m'; nv=length(vars)+1;
30 vars(nv).fldName='AngleSN'; vars(nv).longName='grid orientation (sine)'; vars(nv).units='m'; nv=length(vars)+1;
31 vars(nv).fldName='RC'; vars(nv).longName='vertical coordinate'; vars(nv).units='m'; nv=length(vars)+1;
32 vars(nv).fldName='RF'; vars(nv).longName='vertical coordinate'; vars(nv).units='m'; nv=length(vars)+1;
33 vars(nv).fldName='DRC'; vars(nv).longName='grid spacing'; vars(nv).units='m'; nv=length(vars)+1;
34 vars(nv).fldName='DRF'; vars(nv).longName='grid spacing'; vars(nv).units='m'; nv=length(vars)+1;
35 %
36 structIn.defs=vars;
37 %
38 struct2nctiles('release1/','GRID',structIn,[90 90]);
39 end;
40
41 %===============
42
43 if choiceStruct==2;
44
45 listBudgs={'budgMo','budgHo','budgSo','budgMi','budgHi','budgSi'};
46
47 %directories and snapshot
48 dirIn='r4it11.c65i/';
49 dirMat='mat_budg3d/';
50 fldName='snapshot';
51 timName='0000000732'; suff='initial';
52 % timName='0000174564'; suff='final';
53 t0Name=num2str(3600*str2num(timName));
54
55 %load myparms
56 eval(['load ' dirIn dirMat 'diags_grid_parms.mat myparms;']);
57
58 %read THETA, SALT
59 diagName=['budg3d_snap_set1.' timName];
60 diagName=fullfile(dirIn,'diags',filesep,'BUDG',filesep,diagName);
61 THETA=rdmds2gcmfaces(diagName,'rec',1);
62 SALT=rdmds2gcmfaces(diagName,'rec',2);
63 ONE=1*(SALT>0);
64
65 %compute cell thickness
66 etanName=['budg2d_snap_set1.' timName];
67 etanName=fullfile(dirIn,'diags',filesep,'BUDG',filesep,etanName);
68 etan=rdmds2gcmfaces(etanName,'rec',1);
69 %compute time variable thickness
70 tmp1=mk3D(mygrid.DRF,mygrid.hFacC).*mygrid.hFacC;
71 tmp2=tmp1./mk3D(mygrid.Depth,tmp1);
72 tmp2=tmp2.*mk3D(etan,tmp1);
73 tmp3=mk3D(mygrid.RAC,tmp1);
74 vol=tmp3.*(tmp1+tmp2);
75
76 %load ice and snow volume/m2
77 SIheff=rdmds2gcmfaces(etanName,'rec',2);
78 SIsnow=rdmds2gcmfaces(etanName,'rec',3);
79
80 for ii=1:length(listBudgs);
81 budgName=listBudgs{ii};
82 [budgName ' -- ' fldName ' -- ' timName]
83
84 %directories and snapshot
85 dirOut=fullfile(dirIn,'nctiles_budg',filesep);
86 if ~isdir(dirOut); mkdir(dirOut); end;
87 dirOut=fullfile(dirOut,budgName,filesep);
88 if ~isdir(dirOut); mkdir(dirOut); end;
89
90 switch budgName;
91 case 'budgMo'; fld=ONE; fac=myparms.rhoconst;
92 case 'budgHo'; fld=THETA; fac=myparms.rcp;
93 case 'budgSo'; fld=SALT; fac=myparms.rhoconst;
94 case 'budgMi'; FACheff=myparms.rhoi; FACsnow=myparms.rhosn;
95 case 'budgHi'; FACheff=-myparms.flami*myparms.rhoi; FACsnow=-myparms.flami*myparms.rhosn;
96 case 'budgSi'; FACheff=myparms.SIsal0*myparms.rhoi; FACsnow=0;
97 end;
98
99 %compute extensive snapshot
100 structIn=[];
101 if budgName(end)=='o';
102 structIn.vars.snapshot=fac*fld.*vol;
103 else;
104 structIn.vars.snapshot=mygrid.RAC.*(FACheff*SIheff+FACsnow*SIsnow);
105 end;
106
107 %general description
108 specs=[];
109 if strcmp(budgName(5),'H');
110 structIn.descr={'Heat budget in extensive form (in Watt, on C-Grid)'};
111 specs.units='J';
112 specs.name='heat content';
113 elseif strcmp(budgName(5),'M');
114 structIn.descr={'Mass budget in extensive form (in kg/s, on C-Grid)'};
115 specs.units='kg';
116 specs.name='mass content';
117 elseif strcmp(budgName(5),'S');
118 structIn.descr={'Salt budget in extensive form (in g/s, on C-Grid)'};
119 specs.units='g';
120 specs.name='salt content';
121 end;
122
123 %variables description
124 vars=[]; nv=length(vars)+1;
125 vars(nv).fldName=fldName; vars(nv).units=specs.units;
126 if strcmp(budgName(6),'o'); vars(nv).longName=['ocean ' specs.name ' snapshot'];
127 elseif strcmp(budgName(6),'i'); vars(nv).longName=['ice+snow ' specs.name ' snapshot'];
128 end;
129 vars(nv).longName=[vars(nv).longName ' at ' suff ' time=' t0Name 's'];
130 %
131 structIn.defs=vars;
132
133 %create file
134 tic; struct2nctiles(dirIn,fldName,structIn,[90 90]); toc;
135 eval(['!mv ' pwd filesep dirIn filesep 'tmp_nctiles' filesep fldName ' ' pwd filesep dirOut suff]);
136
137 end;%for ii=1:length(listBudgs);
138
139 end;
140
141 %===============
142
143 if choiceStruct==3;
144
145 %budgName='budgHo';
146
147 if budgName(end)=='o'; listFlds={'tend','trU','trV','trWtop'};
148 else; listFlds={'tend','trU','trV','trWtop','trWbot'};
149 end;
150
151 for ii=1:length(listFlds);
152 fldName=listFlds{ii};
153 %fldName='trWtop';
154 [budgName ' -- ' fldName]
155
156 %directories
157 dirIn='r4it11.c65i/';
158 dirMat='mat_budg3d/';
159 dirOut=fullfile(dirIn,'nctiles_budg',filesep);
160 if ~isdir(dirOut); mkdir(dirOut); end;
161 dirOut=fullfile(dirOut,budgName,filesep);
162 if ~isdir(dirOut); mkdir(dirOut); end;
163
164 %load variable
165 eval(['load ' dirIn dirMat 'diags_grid_parms.mat myparms;']);
166 dirMat=[dirIn dirMat 'diags_set_' budgName '/'];
167 fileMat=[budgName '_*.mat'];
168 tic; structIn.vars=diags_read_from_mat(dirMat,fileMat,fldName); toc;
169
170 %time vectors
171 [listTimes]=diags_list_times({[dirIn 'diags/BUDG/']},{'budg2d_hflux_set1'});
172 structIn.vars.t0=3600*listTimes(1:end-2);
173 structIn.vars.t1=3600*listTimes(2:end-1);
174
175 %rename trWtop as trW if adequate
176 if strcmp(fldName,'trWtop')&(budgName(end)=='o');
177 structIn.vars=setfield(structIn.vars,'trW',structIn.vars.trWtop);
178 structIn.vars=rmfield(structIn.vars,'trWtop');
179 structIn.vars.listDiags={'trW'};
180 fldName='trW';
181 end;
182
183 %switch back to upward convention
184 if strcmp(fldName(1:3),'trW');
185 tmp1=getfield(structIn.vars,fldName);
186 structIn.vars=setfield(structIn.vars,fldName,-tmp1);
187 end;
188
189 %general description
190 tmp1=diags_read_from_mat(dirMat,fileMat,'specs',1);
191 specs=tmp1.specs;
192 if strcmp(specs.units,'W');
193 structIn.descr={'Heat budget in extensive form (in Watt, on C-Grid)'};
194 elseif strcmp(specs.units,'kg/s');
195 structIn.descr={'Mass budget in extensive form (in kg/s, on C-Grid)'};
196 elseif strcmp(specs.units,'g/s');
197 structIn.descr={'Salt budget in extensive form (in g/s, on C-Grid)'};
198 else;
199 error('unknown budget');
200 end;
201
202 %variables description
203 vars=[]; nv=length(vars)+1;
204 vars(nv).fldName=fldName; vars(nv).units=specs.units;
205 switch fldName;
206 case 'tend'; vars(nv).longName='tendency term';
207 case 'trU'; vars(nv).longName='horizontal transport (U)';
208 case 'trV'; vars(nv).longName='horizontal transport (V)';
209 case 'trW'; vars(nv).longName='upward vertical transport (W)';
210 case 'trWtop'; vars(nv).longName='upward vertical transport (W)';
211 case 'trWbot'; vars(nv).longName='upward vertical transport (W)';
212 end;
213 nv=length(vars)+1;
214 vars(nv).fldName='t0'; vars(nv).longName='initial time'; vars(nv).units='s'; nv=length(vars)+1;
215 vars(nv).fldName='t1'; vars(nv).longName='final time'; vars(nv).units='s'; nv=length(vars)+1;
216 %
217 structIn.defs=vars;
218
219 %create file
220 tic; struct2nctiles(dirIn,fldName,structIn,[90 90]); toc;
221 eval(['!mv ' pwd filesep dirIn filesep 'tmp_nctiles' filesep fldName ' ' pwd filesep dirOut]);
222
223 end;%for fldName=listFlds;
224
225 end;
226

  ViewVC Help
Powered by ViewVC 1.1.22