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

Annotation 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.8 - (hide annotations) (download)
Fri Feb 12 16:19:16 2016 UTC (9 years, 5 months ago) by gforget
Branch: MAIN
CVS Tags: checkpoint65x, checkpoint65v, checkpoint65w, checkpoint65u, checkpoint66f, checkpoint66e, checkpoint66d, checkpoint66c, checkpoint66b, checkpoint66a, checkpoint66o, HEAD
Changes since 1.7: +5 -0 lines
- document implications of recent changes

1 gforget 1.1
2 gforget 1.7 choiceStruct=2;
3 gforget 1.1
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 gforget 1.7
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 gforget 1.1 end;
140    
141 gforget 1.2 %===============
142    
143     if choiceStruct==3;
144    
145 gforget 1.5 %budgName='budgHo';
146 gforget 1.4
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 gforget 1.2 [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 gforget 1.3 %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 gforget 1.2 %rename trWtop as trW if adequate
176 gforget 1.4 if strcmp(fldName,'trWtop')&(budgName(end)=='o');
177 gforget 1.8 %note: since geothermal heating was added we likely need to add a nr+1 level to trW
178     % using trWbot(:,:,nr) or output trWtop and trWbot themselves (as done for seaice)
179 gforget 1.2 structIn.vars=setfield(structIn.vars,'trW',structIn.vars.trWtop);
180     structIn.vars=rmfield(structIn.vars,'trWtop');
181     structIn.vars.listDiags={'trW'};
182     fldName='trW';
183     end;
184    
185 gforget 1.5 %switch back to upward convention
186     if strcmp(fldName(1:3),'trW');
187 gforget 1.8 %note: the following was needed in preparing the 03-Feb-2015 nctiles_budget (3d) files
188     % but is no longer needed since the trW sign convention in gcmfaces_diags was
189     % reversed to positive upward for consistency with MITgcm convention.
190 gforget 1.5 tmp1=getfield(structIn.vars,fldName);
191     structIn.vars=setfield(structIn.vars,fldName,-tmp1);
192     end;
193    
194 gforget 1.2 %general description
195     tmp1=diags_read_from_mat(dirMat,fileMat,'specs',1);
196     specs=tmp1.specs;
197     if strcmp(specs.units,'W');
198     structIn.descr={'Heat budget in extensive form (in Watt, on C-Grid)'};
199     elseif strcmp(specs.units,'kg/s');
200     structIn.descr={'Mass budget in extensive form (in kg/s, on C-Grid)'};
201     elseif strcmp(specs.units,'g/s');
202     structIn.descr={'Salt budget in extensive form (in g/s, on C-Grid)'};
203     else;
204     error('unknown budget');
205     end;
206    
207     %variables description
208     vars=[]; nv=length(vars)+1;
209     vars(nv).fldName=fldName; vars(nv).units=specs.units;
210     switch fldName;
211     case 'tend'; vars(nv).longName='tendency term';
212     case 'trU'; vars(nv).longName='horizontal transport (U)';
213     case 'trV'; vars(nv).longName='horizontal transport (V)';
214 gforget 1.5 case 'trW'; vars(nv).longName='upward vertical transport (W)';
215     case 'trWtop'; vars(nv).longName='upward vertical transport (W)';
216     case 'trWbot'; vars(nv).longName='upward vertical transport (W)';
217 gforget 1.2 end;
218     nv=length(vars)+1;
219     vars(nv).fldName='t0'; vars(nv).longName='initial time'; vars(nv).units='s'; nv=length(vars)+1;
220     vars(nv).fldName='t1'; vars(nv).longName='final time'; vars(nv).units='s'; nv=length(vars)+1;
221     %
222     structIn.defs=vars;
223    
224     %create file
225     tic; struct2nctiles(dirIn,fldName,structIn,[90 90]); toc;
226     eval(['!mv ' pwd filesep dirIn filesep 'tmp_nctiles' filesep fldName ' ' pwd filesep dirOut]);
227    
228 gforget 1.4 end;%for fldName=listFlds;
229    
230 gforget 1.2 end;
231    

  ViewVC Help
Powered by ViewVC 1.1.22