/[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.7 - (hide 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 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.2 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 gforget 1.5 %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 gforget 1.2 %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 gforget 1.5 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 gforget 1.2 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 gforget 1.4 end;%for fldName=listFlds;
224    
225 gforget 1.2 end;
226    

  ViewVC Help
Powered by ViewVC 1.1.22