SICOPOLIS V3.0
 All Classes Files Functions Variables Macros
output_nc.F90
Go to the documentation of this file.
1 !+++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++
2 !
3 ! Subroutine : o u t p u t _ n c
4 !
5 !> @file
6 !!
7 !! Writing of time-slice files in NetCDF format.
8 !!
9 !! @section Copyright
10 !!
11 !! Copyright 2009-2013 Ralf Greve, Thomas Goelles
12 !!
13 !! @section License
14 !!
15 !! This file is part of SICOPOLIS.
16 !!
17 !! SICOPOLIS is free software: you can redistribute it and/or modify
18 !! it under the terms of the GNU General Public License as published by
19 !! the Free Software Foundation, either version 3 of the License, or
20 !! (at your option) any later version.
21 !!
22 !! SICOPOLIS is distributed in the hope that it will be useful,
23 !! but WITHOUT ANY WARRANTY; without even the implied warranty of
24 !! MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the
25 !! GNU General Public License for more details.
26 !!
27 !! You should have received a copy of the GNU General Public License
28 !! along with SICOPOLIS. If not, see <http://www.gnu.org/licenses/>.
29 !<
30 !+++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++
31 
32 !-------------------------------------------------------------------------------
33 !> Writing of time-slice files in NetCDF format.
34 !<------------------------------------------------------------------------------
35 subroutine output_nc(runname, time, delta_ts, glac_index, z_sl, &
36  flag_3d_output, ndat2d, ndat3d)
37 
38 use sico_types
40 use netcdf
41 
42 implicit none
43 
44 real(dp), intent (in) :: time, delta_ts, glac_index, z_sl
45 character (len=100), intent (in) :: runname
46 logical, intent (in) :: flag_3d_output
47 
48 integer(i4b), intent (inout) :: ndat2d, ndat3d
49 
50 integer(i4b) :: i, j, kc, kt, kr
51 integer(i4b) :: ndat, ndat_help, ndat_1000s, ndat_100s, ndat_10s, ndat_1s
52 real(dp) :: v_tot, a_grounded, a_floating
53 real(sp) :: time_conv, delta_ts_conv, glac_index_conv, z_sl_conv, &
54  v_tot_conv, a_grounded_conv, a_floating_conv, &
55  h_r_conv, xi_conv(0:imax), eta_conv(0:jmax), &
56  sigma_level_c_conv(0:kcmax), sigma_level_t_conv(0:ktmax), &
57  sigma_level_r_conv(0:krmax), &
58  lambda_conv(0:jmax,0:imax), phi_conv(0:jmax,0:imax), &
59  lon_conv(0:jmax,0:imax), lat_conv(0:jmax,0:imax), &
60  temp_s_conv(0:jmax,0:imax), as_perp_conv(0:jmax,0:imax), &
61  zs_conv(0:jmax,0:imax), zm_conv(0:jmax,0:imax), &
62  zb_conv(0:jmax,0:imax), zl_conv(0:jmax,0:imax), &
63  h_c_conv(0:jmax,0:imax), h_t_conv(0:jmax,0:imax), &
64  h_conv(0:jmax,0:imax), &
65  q_bm_conv(0:jmax,0:imax), q_tld_conv(0:jmax,0:imax), &
66  am_perp_conv(0:jmax,0:imax), &
67  qx_conv(0:jmax,0:imax), qy_conv(0:jmax,0:imax), &
68  temp_r_conv(0:krmax,0:jmax,0:imax), &
69  vx_c_conv(0:kcmax,0:jmax,0:imax), &
70  vy_c_conv(0:kcmax,0:jmax,0:imax), &
71  vz_c_conv(0:kcmax,0:jmax,0:imax), &
72  vx_t_conv(0:ktmax,0:jmax,0:imax), &
73  vy_t_conv(0:ktmax,0:jmax,0:imax), &
74  vz_t_conv(0:ktmax,0:jmax,0:imax), &
75  temp_c_conv(0:kcmax,0:jmax,0:imax), &
76  age_c_conv(0:kcmax,0:jmax,0:imax), &
77  omega_t_conv(0:ktmax,0:jmax,0:imax), &
78  age_t_conv(0:ktmax,0:jmax,0:imax), &
79  dzs_dtau_conv(0:jmax,0:imax), dzm_dtau_conv(0:jmax,0:imax), &
80  dzb_dtau_conv(0:jmax,0:imax), dzl_dtau_conv(0:jmax,0:imax), &
81  dh_c_dtau_conv(0:jmax,0:imax), dh_t_dtau_conv(0:jmax,0:imax), &
82  dh_dtau_conv(0:jmax,0:imax), &
83  vx_b_g_conv(0:jmax,0:imax), vy_b_g_conv(0:jmax,0:imax), &
84  vz_b_conv(0:jmax,0:imax), &
85  vx_s_g_conv(0:jmax,0:imax), vy_s_g_conv(0:jmax,0:imax), &
86  vz_s_conv(0:jmax,0:imax), &
87  temp_b_conv(0:jmax,0:imax), temph_b_conv(0:jmax,0:imax), &
88  p_b_w_conv(0:jmax,0:imax), h_w_conv(0:jmax,0:imax), &
89  q_gl_g_conv(0:jmax,0:imax)
90 real(sp) :: lon0, lat0
91 character (len=256) :: filename
92 character :: ch_1000s, ch_100s, ch_10s, ch_1s
93 
94 integer(i4b) :: ncid, ncv
95 ! ncid: ID of the output file
96 ! ncv: Variable ID
97 integer(i4b) :: ncd, nc1d, nc2d(2), nc3d(3)
98 ! ncd: Dimension ID
99 ! nc1d: Dimension of a 1-d array
100 ! nc2d: Vector with the dimensions of a 2-d array
101 ! nc3d: Vector with the dimensions of a 3-d array
102 integer(i4b) :: nc3flag(3), nc4flag(4)
103 ! nc3flag: Vector with the 3 possible values of a flag variable
104 ! nc4flag: Vector with the 4 possible values of a flag variable
105 integer(i4b) :: nc1cor(1), nc2cor(2), nc3cor(3)
106 ! nc1cor(1): Corner of a 1-d array
107 ! nc2cor(2): Corner of a 2-d array
108 ! nc3cor(3): Corner of a 3-d array
109 integer(i4b) :: nc1cnt(1), nc2cnt(2), nc3cnt(3)
110 ! nc1cnt(1): Count of a 1-d array
111 ! nc2cnt(2): Count of a 2-d array
112 ! nc3cnt(3): Count of a 3-d array
113 character (len= 16) :: ch_date, ch_time, ch_zone
114 character (len=256) :: buffer
115 
116 !-------- Create consecutively numbered file names --------
117 
118 if (flag_3d_output) then
119  ndat = ndat3d
120 else
121  ndat = ndat2d
122 end if
123 
124 if (ndat > 9999) stop ' output_nc: Too many time-slice files!'
125 
126 ndat_help = ndat
127 ndat_1000s = ndat_help/1000
128 ndat_help = ndat_help-ndat_1000s*1000
129 ndat_100s = ndat_help/100
130 ndat_help = ndat_help-ndat_100s*100
131 ndat_10s = ndat_help/10
132 ndat_help = ndat_help-ndat_10s*10
133 ndat_1s = ndat_help
134 
135 ch_1000s = char(ndat_1000s+ichar('0'))
136 ch_100s = char(ndat_100s +ichar('0'))
137 ch_10s = char(ndat_10s +ichar('0'))
138 ch_1s = char(ndat_1s +ichar('0'))
139 
140 if (flag_3d_output) then
141  filename = trim(runname)//ch_1000s//ch_100s//ch_10s//ch_1s//'.nc'
142 else
143  filename = trim(runname)//'_2d_'//ch_1000s//ch_100s//ch_10s//ch_1s//'.nc'
144 end if
145 
146 !-------- NetCDF initialization --------
147 
148 ! ------ Open NetCDF file
149 
150 buffer = outpath//'/'//trim(filename)
151 call check( nf90_create(trim(buffer), nf90_noclobber, ncid) )
152 
153 ! ------ Global attributes
154 
155 buffer = 'Time-slice output no. '//ch_1000s//ch_100s//ch_10s//ch_1s//' '// &
156  'of simulation '//trim(runname)
157 call check( nf90_put_att(ncid, nf90_global, 'title', trim(buffer)) )
158 
159 buffer = 'Institute of Low Temperature Science, Hokkaido University, '// &
160  'Sapporo, Japan'
161 call check( nf90_put_att(ncid, nf90_global, 'institution', trim(buffer)) )
162 
163 buffer = 'SICOPOLIS Version '//version
164 call check( nf90_put_att(ncid, nf90_global, 'source', trim(buffer)) )
165 
166 call date_and_time(ch_date, ch_time, ch_zone)
167 buffer = ch_date(1:4)//'-'//ch_date(5:6)//'-'//ch_date(7:8)//' '// &
168  ch_time(1:2)//':'//ch_time(3:4)//':'//ch_time(5:6)//' '// &
169  ch_zone(1:3)//':'//ch_zone(4:5)//' - Data produced'
170 call check( nf90_put_att(ncid, nf90_global, 'history', trim(buffer)) )
171 
172 buffer = 'http://sicopolis.greveweb.net/'
173 call check( nf90_put_att(ncid, nf90_global, 'references', trim(buffer)) )
174 
175 ! ------ Definition of the dimensions
176 
177 call check( nf90_def_dim(ncid, 'IMAX', imax+1, ncd) )
178 call check( nf90_def_dim(ncid, 'JMAX', jmax+1, ncd) )
179 call check( nf90_def_dim(ncid, 'KCMAX', kcmax+1, ncd) )
180 call check( nf90_def_dim(ncid, 'KTMAX', ktmax+1, ncd) )
181 call check( nf90_def_dim(ncid, 'KRMAX', krmax+1, ncd) )
182 
183 ! ------ Definition of the variables
184 
185 call check( nf90_def_var(ncid, 'crs', nf90_short, ncv) )
186 #if ( (GRID == 0) || (GRID == 1) )
187 buffer = 'polar_stereographic'
188 call check( nf90_put_att(ncid, ncv, 'grid_mapping_name', trim(buffer)) )
189 #elif (GRID == 2)
190 buffer = 'latitude_longitude'
191 call check( nf90_put_att(ncid, ncv, 'grid_mapping_name', trim(buffer)) )
192 #endif
193 #if ( defined(ANT) \
194  || defined(asf) \
195  || defined(grl) \
196  || defined(scand) \
197  || defined(nhem) \
198  || defined(tibet) \
199  || defined(emtp2sge) \
200  || defined(heino) ) /* terrestrial ice sheet */
201 buffer = 'WGS84'
202 call check( nf90_put_att(ncid, ncv, 'ellipsoid', trim(buffer)) )
203 #elif ( defined(NMARS) \
204  || defined(smars) ) /* martian ice sheet */
205 buffer = 'Mars_ellipsoid'
206 call check( nf90_put_att(ncid, ncv, 'ellipsoid', trim(buffer)) )
207 #else
208 stop ' output_nc: No valid domain (ANT, GRL etc.) specified!'
209 #endif
210 #if ( (GRID == 0) || (GRID == 1) )
211 call check( nf90_put_att(ncid, ncv, 'false_easting', 0.0) )
212 call check( nf90_put_att(ncid, ncv, 'false_northing', 0.0) )
213 lon0 = lambda0 *pi_180_inv
214 lon0 = modulo(lon0+180.0_sp, 360.0_sp)-180.0_sp
215 lon0 = nint(lon0*1.0e+04_sp)*1.0e-04_sp
216 lat0 = phi0 *pi_180_inv
217 if (lat0 > 90.0_sp) lat0 = 90.0_sp
218 if (lat0 < -90.0_sp) lat0 = -90.0_sp
219 lat0 = nint(lat0*1.0e+04_sp)*1.0e-04_sp
220  ! reference longitude and standard parallel in deg rounded to 4 digits
221 if (lat0 >= 0.0_sp) then
222  call check( nf90_put_att(ncid, ncv, &
223  'latitude_of_projection_origin', 90.0) )
224 else
225  call check( nf90_put_att(ncid, ncv, &
226  'latitude_of_projection_origin', -90.0) )
227 end if
228 call check( nf90_put_att(ncid, ncv, &
229  'straight_vertical_longitude_from_pole', lon0) )
230 call check( nf90_put_att(ncid, ncv, &
231  'standard_parallel', lat0) )
232 #endif
233 
234 call check( nf90_def_var(ncid, 'time', nf90_float, ncv) )
235 buffer = 's'
236 call check( nf90_put_att(ncid, ncv, 'units', trim(buffer)) )
237 buffer = 'time'
238 call check( nf90_put_att(ncid, ncv, 'standard_name', trim(buffer)) )
239 buffer = 'Time'
240 call check( nf90_put_att(ncid, ncv, 'long_name', trim(buffer)) )
241 
242 if (forcing_flag == 1) then
243 
244  call check( nf90_def_var(ncid, 'delta_ts', nf90_float, ncv) )
245  buffer = 'degC'
246  call check( nf90_put_att(ncid, ncv, 'units', trim(buffer)) )
247  buffer = 'surface_temperature_anomaly'
248  call check( nf90_put_att(ncid, ncv, 'standard_name', trim(buffer)) )
249  buffer = 'Surface temperature anomaly'
250  call check( nf90_put_att(ncid, ncv, 'long_name', trim(buffer)) )
251 
252 else if (forcing_flag == 2) then
253 
254  call check( nf90_def_var(ncid, 'glac_index', nf90_float, ncv) )
255  buffer = '1'
256  call check( nf90_put_att(ncid, ncv, 'units', trim(buffer)) )
257  buffer = 'glacial_index'
258  call check( nf90_put_att(ncid, ncv, 'standard_name', trim(buffer)) )
259  buffer = 'Glacial index'
260  call check( nf90_put_att(ncid, ncv, 'long_name', trim(buffer)) )
261 
262 else if (forcing_flag == 3) then
263 
264  call check( nf90_def_var(ncid, 'glac_index', nf90_float, ncv) )
265  buffer = '1'
266  call check( nf90_put_att(ncid, ncv, 'units', trim(buffer)) )
267  buffer = 'glacial_index'
268  call check( nf90_put_att(ncid, ncv, 'standard_name', trim(buffer)) )
269  buffer = 'Glacial index'
270  call check( nf90_put_att(ncid, ncv, 'long_name', trim(buffer)) )
271  buffer = 'This variable will be assigned a dummy value only!'
272  call check( nf90_put_att(ncid, ncv, 'comment', trim(buffer)) )
273 
274 end if
275 
276 call check( nf90_def_var(ncid, 'z_sl', nf90_float, ncv) )
277 buffer = 'm'
278 call check( nf90_put_att(ncid, ncv, 'units', trim(buffer)) )
279 buffer = 'global_average_sea_level_change'
280 call check( nf90_put_att(ncid, ncv, 'standard_name', trim(buffer)) )
281 buffer = 'Sea level'
282 call check( nf90_put_att(ncid, ncv, 'long_name', trim(buffer)) )
283 
284 call check( nf90_inq_dimid(ncid, 'IMAX', nc1d) )
285 call check( nf90_def_var(ncid, 'xi', nf90_float, nc1d, ncv) )
286 buffer = 'm'
287 call check( nf90_put_att(ncid, ncv, 'units', trim(buffer)) )
288 buffer = 'projection_x_coordinate'
289 call check( nf90_put_att(ncid, ncv, 'standard_name', trim(buffer)) )
290 buffer = 'x-coordinate of the grid point i'
291 call check( nf90_put_att(ncid, ncv, 'long_name', trim(buffer)) )
292 
293 call check( nf90_inq_dimid(ncid, 'JMAX', nc1d) )
294 call check( nf90_def_var(ncid, 'eta', nf90_float, nc1d, ncv) )
295 buffer = 'm'
296 call check( nf90_put_att(ncid, ncv, 'units', trim(buffer)) )
297 buffer = 'projection_y_coordinate'
298 call check( nf90_put_att(ncid, ncv, 'standard_name', trim(buffer)) )
299 buffer = 'y-coordinate of the grid point j'
300 call check( nf90_put_att(ncid, ncv, 'long_name', trim(buffer)) )
301 
302 call check( nf90_inq_dimid(ncid, 'KCMAX', nc1d) )
303 call check( nf90_def_var(ncid, 'sigma_level_c', nf90_float, nc1d, ncv) )
304 buffer = 'up'
305 call check( nf90_put_att(ncid, ncv, 'positive', trim(buffer)) )
306 buffer = 'land_ice_cold_layer_sigma_coordinate'
307 call check( nf90_put_att(ncid, ncv, 'standard_name', trim(buffer)) )
308 buffer = 'sigma-coordinate of the grid point kc'
309 call check( nf90_put_att(ncid, ncv, 'long_name', trim(buffer)) )
310 
311 call check( nf90_inq_dimid(ncid, 'KTMAX', nc1d) )
312 call check( nf90_def_var(ncid, 'sigma_level_t', nf90_float, nc1d, ncv) )
313 buffer = 'up'
314 call check( nf90_put_att(ncid, ncv, 'positive', trim(buffer)) )
315 buffer = 'land_ice_temperate_layer_sigma_coordinate'
316 call check( nf90_put_att(ncid, ncv, 'standard_name', trim(buffer)) )
317 buffer = 'sigma-coordinate of the grid point kt'
318 call check( nf90_put_att(ncid, ncv, 'long_name', trim(buffer)) )
319 
320 call check( nf90_inq_dimid(ncid, 'KRMAX', nc1d) )
321 call check( nf90_def_var(ncid, 'sigma_level_r', nf90_float, nc1d, ncv) )
322 buffer = 'up'
323 call check( nf90_put_att(ncid, ncv, 'positive', trim(buffer)) )
324 buffer = 'lithosphere_layer_sigma_coordinate'
325 call check( nf90_put_att(ncid, ncv, 'standard_name', trim(buffer)) )
326 buffer = 'sigma-coordinate of the grid point kr'
327 call check( nf90_put_att(ncid, ncv, 'long_name', trim(buffer)) )
328 
329 call check( nf90_inq_dimid(ncid, 'JMAX', nc2d(1)) )
330 call check( nf90_inq_dimid(ncid, 'IMAX', nc2d(2)) )
331 call check( nf90_def_var(ncid, 'lon', nf90_float, nc2d, ncv) )
332 buffer = 'degrees_E'
333 call check( nf90_put_att(ncid, ncv, 'units', trim(buffer)) )
334 buffer = 'longitude'
335 call check( nf90_put_att(ncid, ncv, 'standard_name', trim(buffer)) )
336 buffer = 'Geographical longitude'
337 call check( nf90_put_att(ncid, ncv, 'long_name', trim(buffer)) )
338 
339 call check( nf90_inq_dimid(ncid, 'JMAX', nc2d(1)) )
340 call check( nf90_inq_dimid(ncid, 'IMAX', nc2d(2)) )
341 call check( nf90_def_var(ncid, 'lat', nf90_float, nc2d, ncv) )
342 buffer = 'degrees_N'
343 call check( nf90_put_att(ncid, ncv, 'units', trim(buffer)) )
344 buffer = 'latitude'
345 call check( nf90_put_att(ncid, ncv, 'standard_name', trim(buffer)) )
346 buffer = 'Geographical latitude'
347 call check( nf90_put_att(ncid, ncv, 'long_name', trim(buffer)) )
348 
349 call check( nf90_inq_dimid(ncid, 'JMAX', nc2d(1)) )
350 call check( nf90_inq_dimid(ncid, 'IMAX', nc2d(2)) )
351 call check( nf90_def_var(ncid, 'lambda', nf90_float, nc2d, ncv) )
352 buffer = 'rad'
353 call check( nf90_put_att(ncid, ncv, 'units', trim(buffer)) )
354 buffer = 'longitude'
355 call check( nf90_put_att(ncid, ncv, 'standard_name', trim(buffer)) )
356 buffer = 'Geographical longitude'
357 call check( nf90_put_att(ncid, ncv, 'long_name', trim(buffer)) )
358 call check( nf90_put_att(ncid, ncv, 'grid_mapping', 'crs') )
359 call check( nf90_put_att(ncid, ncv, 'coordinates', 'lat lon') )
360 
361 call check( nf90_inq_dimid(ncid, 'JMAX', nc2d(1)) )
362 call check( nf90_inq_dimid(ncid, 'IMAX', nc2d(2)) )
363 call check( nf90_def_var(ncid, 'phi', nf90_float, nc2d, ncv) )
364 buffer = 'rad'
365 call check( nf90_put_att(ncid, ncv, 'units', trim(buffer)) )
366 buffer = 'latitude'
367 call check( nf90_put_att(ncid, ncv, 'standard_name', trim(buffer)) )
368 buffer = 'Geographical latitude'
369 call check( nf90_put_att(ncid, ncv, 'long_name', trim(buffer)) )
370 call check( nf90_put_att(ncid, ncv, 'grid_mapping', 'crs') )
371 call check( nf90_put_att(ncid, ncv, 'coordinates', 'lat lon') )
372 
373 call check( nf90_inq_dimid(ncid, 'JMAX', nc2d(1)) )
374 call check( nf90_inq_dimid(ncid, 'IMAX', nc2d(2)) )
375 call check( nf90_def_var(ncid, 'temp_s', nf90_float, nc2d, ncv) )
376 buffer = 'degC'
377 call check( nf90_put_att(ncid, ncv, 'units', trim(buffer)) )
378 buffer = 'surface_temperature'
379 call check( nf90_put_att(ncid, ncv, 'standard_name', trim(buffer)) )
380 buffer = 'Temperature at the ice surface'
381 call check( nf90_put_att(ncid, ncv, 'long_name', trim(buffer)) )
382 call check( nf90_put_att(ncid, ncv, 'grid_mapping', 'crs') )
383 call check( nf90_put_att(ncid, ncv, 'coordinates', 'lat lon') )
384 
385 call check( nf90_inq_dimid(ncid, 'JMAX', nc2d(1)) )
386 call check( nf90_inq_dimid(ncid, 'IMAX', nc2d(2)) )
387 call check( nf90_def_var(ncid, 'as_perp', nf90_float, nc2d, ncv) )
388 buffer = 'm s-1'
389 call check( nf90_put_att(ncid, ncv, 'units', trim(buffer)) )
390 buffer = 'land_ice_surface_mass_balance'
391 call check( nf90_put_att(ncid, ncv, 'standard_name', trim(buffer)) )
392 buffer = 'Mass balance at the ice surface'
393 call check( nf90_put_att(ncid, ncv, 'long_name', trim(buffer)) )
394 call check( nf90_put_att(ncid, ncv, 'grid_mapping', 'crs') )
395 call check( nf90_put_att(ncid, ncv, 'coordinates', 'lat lon') )
396 
397 call check( nf90_inq_dimid(ncid, 'JMAX', nc2d(1)) )
398 call check( nf90_inq_dimid(ncid, 'IMAX', nc2d(2)) )
399 call check( nf90_def_var(ncid, 'maske', nf90_short, nc2d, ncv) )
400 buffer = 'ice_land_sea_mask'
401 call check( nf90_put_att(ncid, ncv, 'standard_name', trim(buffer)) )
402 buffer = 'Ice-land-sea mask'
403 call check( nf90_put_att(ncid, ncv, 'long_name', trim(buffer)) )
404 nc4flag(1) = 0
405 nc4flag(2) = 1
406 nc4flag(3) = 2
407 nc4flag(4) = 3
408 call check( nf90_put_att(ncid, ncv, 'flag_values', nc4flag) )
409 buffer = 'glaciated_land '// &
410  'ice_free_land '// &
411  'sea '// &
412  'floating_ice'
413 call check( nf90_put_att(ncid, ncv, 'flag_meanings', trim(buffer)) )
414 call check( nf90_put_att(ncid, ncv, 'grid_mapping', 'crs') )
415 call check( nf90_put_att(ncid, ncv, 'coordinates', 'lat lon') )
416 
417 call check( nf90_inq_dimid(ncid, 'JMAX', nc2d(1)) )
418 call check( nf90_inq_dimid(ncid, 'IMAX', nc2d(2)) )
419 call check( nf90_def_var(ncid, 'n_cts', nf90_short, nc2d, ncv) )
420 buffer = 'polythermal_condition_mask'
421 call check( nf90_put_att(ncid, ncv, 'standard_name', trim(buffer)) )
422 buffer = 'Mask for polythermal conditions'
423 call check( nf90_put_att(ncid, ncv, 'long_name', trim(buffer)) )
424 nc3flag(1) = -1
425 nc3flag(2) = 0
426 nc3flag(3) = 1
427 call check( nf90_put_att(ncid, ncv, 'flag_values', nc3flag) )
428 buffer = 'cold_base '// &
429  'temperate_base_with_cold_ice_above '// &
430  'temperate_base_with_temperate_ice_above'
431 call check( nf90_put_att(ncid, ncv, 'flag_meanings', trim(buffer)) )
432 call check( nf90_put_att(ncid, ncv, 'grid_mapping', 'crs') )
433 call check( nf90_put_att(ncid, ncv, 'coordinates', 'lat lon') )
434 
435 call check( nf90_inq_dimid(ncid, 'JMAX', nc2d(1)) )
436 call check( nf90_inq_dimid(ncid, 'IMAX', nc2d(2)) )
437 call check( nf90_def_var(ncid, 'zs', nf90_float, nc2d, ncv) )
438 buffer = 'm'
439 call check( nf90_put_att(ncid, ncv, 'units', trim(buffer)) )
440 buffer = 'surface_altitude'
441 call check( nf90_put_att(ncid, ncv, 'standard_name', trim(buffer)) )
442 buffer = 'Topography of the free surface'
443 call check( nf90_put_att(ncid, ncv, 'long_name', trim(buffer)) )
444 call check( nf90_put_att(ncid, ncv, 'grid_mapping', 'crs') )
445 call check( nf90_put_att(ncid, ncv, 'coordinates', 'lat lon') )
446 
447 call check( nf90_inq_dimid(ncid, 'JMAX', nc2d(1)) )
448 call check( nf90_inq_dimid(ncid, 'IMAX', nc2d(2)) )
449 call check( nf90_def_var(ncid, 'zm', nf90_float, nc2d, ncv) )
450 buffer = 'm'
451 call check( nf90_put_att(ncid, ncv, 'units', trim(buffer)) )
452 buffer = 'cts_altitude'
453 call check( nf90_put_att(ncid, ncv, 'standard_name', trim(buffer)) )
454 buffer = 'Topography of the CTS (cold-temperate transition surface)'
455 call check( nf90_put_att(ncid, ncv, 'long_name', trim(buffer)) )
456 call check( nf90_put_att(ncid, ncv, 'grid_mapping', 'crs') )
457 call check( nf90_put_att(ncid, ncv, 'coordinates', 'lat lon') )
458 
459 call check( nf90_inq_dimid(ncid, 'JMAX', nc2d(1)) )
460 call check( nf90_inq_dimid(ncid, 'IMAX', nc2d(2)) )
461 call check( nf90_def_var(ncid, 'zb', nf90_float, nc2d, ncv) )
462 buffer = 'm'
463 call check( nf90_put_att(ncid, ncv, 'units', trim(buffer)) )
464 buffer = 'ice_base_altitude'
465 call check( nf90_put_att(ncid, ncv, 'standard_name', trim(buffer)) )
466 buffer = 'Topography of the ice base'
467 call check( nf90_put_att(ncid, ncv, 'long_name', trim(buffer)) )
468 call check( nf90_put_att(ncid, ncv, 'grid_mapping', 'crs') )
469 call check( nf90_put_att(ncid, ncv, 'coordinates', 'lat lon') )
470 
471 call check( nf90_inq_dimid(ncid, 'JMAX', nc2d(1)) )
472 call check( nf90_inq_dimid(ncid, 'IMAX', nc2d(2)) )
473 call check( nf90_def_var(ncid, 'zl', nf90_float, nc2d, ncv) )
474 buffer = 'm'
475 call check( nf90_put_att(ncid, ncv, 'units', trim(buffer)) )
476 buffer = 'bedrock_altitude'
477 call check( nf90_put_att(ncid, ncv, 'standard_name', trim(buffer)) )
478 buffer = 'Topography of the lithosphere surface'
479 call check( nf90_put_att(ncid, ncv, 'long_name', trim(buffer)) )
480 call check( nf90_put_att(ncid, ncv, 'grid_mapping', 'crs') )
481 call check( nf90_put_att(ncid, ncv, 'coordinates', 'lat lon') )
482 
483 call check( nf90_inq_dimid(ncid, 'JMAX', nc2d(1)) )
484 call check( nf90_inq_dimid(ncid, 'IMAX', nc2d(2)) )
485 call check( nf90_def_var(ncid, 'H_c', nf90_float, nc2d, ncv) )
486 buffer = 'm'
487 call check( nf90_put_att(ncid, ncv, 'units', trim(buffer)) )
488 buffer = 'land_ice_cold_layer_thickness'
489 call check( nf90_put_att(ncid, ncv, 'standard_name', trim(buffer)) )
490 buffer = 'Thickness of the cold ice layer'
491 call check( nf90_put_att(ncid, ncv, 'long_name', trim(buffer)) )
492 call check( nf90_put_att(ncid, ncv, 'grid_mapping', 'crs') )
493 call check( nf90_put_att(ncid, ncv, 'coordinates', 'lat lon') )
494 
495 call check( nf90_inq_dimid(ncid, 'JMAX', nc2d(1)) )
496 call check( nf90_inq_dimid(ncid, 'IMAX', nc2d(2)) )
497 call check( nf90_def_var(ncid, 'H_t', nf90_float, nc2d, ncv) )
498 buffer = 'm'
499 call check( nf90_put_att(ncid, ncv, 'units', trim(buffer)) )
500 buffer = 'land_ice_temperate_layer_thickness'
501 call check( nf90_put_att(ncid, ncv, 'standard_name', trim(buffer)) )
502 buffer = 'Thickness of the temperate ice layer'
503 call check( nf90_put_att(ncid, ncv, 'long_name', trim(buffer)) )
504 call check( nf90_put_att(ncid, ncv, 'grid_mapping', 'crs') )
505 call check( nf90_put_att(ncid, ncv, 'coordinates', 'lat lon') )
506 
507 call check( nf90_inq_dimid(ncid, 'JMAX', nc2d(1)) )
508 call check( nf90_inq_dimid(ncid, 'IMAX', nc2d(2)) )
509 call check( nf90_def_var(ncid, 'H', nf90_float, nc2d, ncv) )
510 buffer = 'm'
511 call check( nf90_put_att(ncid, ncv, 'units', trim(buffer)) )
512 buffer = 'land_ice_thickness'
513 call check( nf90_put_att(ncid, ncv, 'standard_name', trim(buffer)) )
514 buffer = 'Ice thickness'
515 call check( nf90_put_att(ncid, ncv, 'long_name', trim(buffer)) )
516 call check( nf90_put_att(ncid, ncv, 'grid_mapping', 'crs') )
517 call check( nf90_put_att(ncid, ncv, 'coordinates', 'lat lon') )
518 
519 call check( nf90_def_var(ncid, 'H_R', nf90_float, ncv) )
520 buffer = 'm'
521 call check( nf90_put_att(ncid, ncv, 'units', trim(buffer)) )
522 buffer = 'lithosphere_layer_thickness'
523 call check( nf90_put_att(ncid, ncv, 'standard_name', trim(buffer)) )
524 buffer = 'Thickness of the lithosphere layer'
525 call check( nf90_put_att(ncid, ncv, 'long_name', trim(buffer)) )
526 
527 if (flag_3d_output) then
528 
529  call check( nf90_inq_dimid(ncid, 'KCMAX', nc3d(1)) )
530  call check( nf90_inq_dimid(ncid, 'JMAX', nc3d(2)) )
531  call check( nf90_inq_dimid(ncid, 'IMAX', nc3d(3)) )
532  call check( nf90_def_var(ncid, 'vx_c', nf90_float, nc3d, ncv) )
533  buffer = 'm s-1'
534  call check( nf90_put_att(ncid, ncv, 'units', trim(buffer)) )
535  buffer = 'land_ice_cold_layer_x_velocity'
536  call check( nf90_put_att(ncid, ncv, 'standard_name', trim(buffer)) )
537  buffer = 'Horizontal velocity vx in the cold ice layer'
538  call check( nf90_put_att(ncid, ncv, 'long_name', trim(buffer)) )
539  buffer = 'Staggered grid variable, defined at (kc,j,i+1/2)'
540  call check( nf90_put_att(ncid, ncv, 'comment', trim(buffer)) )
541  call check( nf90_put_att(ncid, ncv, 'grid_mapping', 'crs') )
542  call check( nf90_put_att(ncid, ncv, 'coordinates', 'lat lon') )
543 
544  call check( nf90_inq_dimid(ncid, 'KCMAX', nc3d(1)) )
545  call check( nf90_inq_dimid(ncid, 'JMAX', nc3d(2)) )
546  call check( nf90_inq_dimid(ncid, 'IMAX', nc3d(3)) )
547  call check( nf90_def_var(ncid, 'vy_c', nf90_float, nc3d, ncv) )
548  buffer = 'm s-1'
549  call check( nf90_put_att(ncid, ncv, 'units', trim(buffer)) )
550  buffer = 'land_ice_cold_layer_y_velocity'
551  call check( nf90_put_att(ncid, ncv, 'standard_name', trim(buffer)) )
552  buffer = 'Horizontal velocity vy in the cold ice layer'
553  call check( nf90_put_att(ncid, ncv, 'long_name', trim(buffer)) )
554  buffer = 'Staggered grid variable, defined at (kc,j+1/2,i)'
555  call check( nf90_put_att(ncid, ncv, 'comment', trim(buffer)) )
556  call check( nf90_put_att(ncid, ncv, 'grid_mapping', 'crs') )
557  call check( nf90_put_att(ncid, ncv, 'coordinates', 'lat lon') )
558 
559  call check( nf90_inq_dimid(ncid, 'KCMAX', nc3d(1)) )
560  call check( nf90_inq_dimid(ncid, 'JMAX', nc3d(2)) )
561  call check( nf90_inq_dimid(ncid, 'IMAX', nc3d(3)) )
562  call check( nf90_def_var(ncid, 'vz_c', nf90_float, nc3d, ncv) )
563  buffer = 'm s-1'
564  call check( nf90_put_att(ncid, ncv, 'units', trim(buffer)) )
565  buffer = 'land_ice_cold_layer_z_velocity'
566  call check( nf90_put_att(ncid, ncv, 'standard_name', trim(buffer)) )
567  buffer = 'Vertical velocity vz in the cold ice layer'
568  call check( nf90_put_att(ncid, ncv, 'long_name', trim(buffer)) )
569  buffer = 'Staggered grid variable, defined at (kc+1/2,j,i)'
570  call check( nf90_put_att(ncid, ncv, 'comment', trim(buffer)) )
571  call check( nf90_put_att(ncid, ncv, 'grid_mapping', 'crs') )
572  call check( nf90_put_att(ncid, ncv, 'coordinates', 'lat lon') )
573 
574  call check( nf90_inq_dimid(ncid, 'KTMAX', nc3d(1)) )
575  call check( nf90_inq_dimid(ncid, 'JMAX', nc3d(2)) )
576  call check( nf90_inq_dimid(ncid, 'IMAX', nc3d(3)) )
577  call check( nf90_def_var(ncid, 'vx_t', nf90_float, nc3d, ncv) )
578  buffer = 'm s-1'
579  call check( nf90_put_att(ncid, ncv, 'units', trim(buffer)) )
580  buffer = 'land_ice_temperate_layer_x_velocity'
581  call check( nf90_put_att(ncid, ncv, 'standard_name', trim(buffer)) )
582  buffer = 'Horizontal velocity vx in the temperate ice layer'
583  call check( nf90_put_att(ncid, ncv, 'long_name', trim(buffer)) )
584  buffer = 'Staggered grid variable, defined at (kt,j,i+1/2)'
585  call check( nf90_put_att(ncid, ncv, 'comment', trim(buffer)) )
586  call check( nf90_put_att(ncid, ncv, 'grid_mapping', 'crs') )
587  call check( nf90_put_att(ncid, ncv, 'coordinates', 'lat lon') )
588 
589  call check( nf90_inq_dimid(ncid, 'KTMAX', nc3d(1)) )
590  call check( nf90_inq_dimid(ncid, 'JMAX', nc3d(2)) )
591  call check( nf90_inq_dimid(ncid, 'IMAX', nc3d(3)) )
592  call check( nf90_def_var(ncid, 'vy_t', nf90_float, nc3d, ncv) )
593  buffer = 'm s-1'
594  call check( nf90_put_att(ncid, ncv, 'units', trim(buffer)) )
595  buffer = 'land_ice_temperate_layer_y_velocity'
596  call check( nf90_put_att(ncid, ncv, 'standard_name', trim(buffer)) )
597  buffer = 'Horizontal velocity vy in the temperate ice layer'
598  call check( nf90_put_att(ncid, ncv, 'long_name', trim(buffer)) )
599  buffer = 'Staggered grid variable, defined at (kt,j+1/2,i)'
600  call check( nf90_put_att(ncid, ncv, 'comment', trim(buffer)) )
601  call check( nf90_put_att(ncid, ncv, 'grid_mapping', 'crs') )
602  call check( nf90_put_att(ncid, ncv, 'coordinates', 'lat lon') )
603 
604  call check( nf90_inq_dimid(ncid, 'KTMAX', nc3d(1)) )
605  call check( nf90_inq_dimid(ncid, 'JMAX', nc3d(2)) )
606  call check( nf90_inq_dimid(ncid, 'IMAX', nc3d(3)) )
607  call check( nf90_def_var(ncid, 'vz_t', nf90_float, nc3d, ncv) )
608  buffer = 'm s-1'
609  call check( nf90_put_att(ncid, ncv, 'units', trim(buffer)) )
610  buffer = 'land_ice_temperate_layer_z_velocity'
611  call check( nf90_put_att(ncid, ncv, 'standard_name', trim(buffer)) )
612  buffer = 'Vertical velocity vz in the temperate ice layer'
613  call check( nf90_put_att(ncid, ncv, 'long_name', trim(buffer)) )
614  buffer = 'Staggered grid variable, defined at (kt+1/2,j,i)'
615  call check( nf90_put_att(ncid, ncv, 'comment', trim(buffer)) )
616  call check( nf90_put_att(ncid, ncv, 'grid_mapping', 'crs') )
617  call check( nf90_put_att(ncid, ncv, 'coordinates', 'lat lon') )
618 
619  call check( nf90_inq_dimid(ncid, 'KCMAX', nc3d(1)) )
620  call check( nf90_inq_dimid(ncid, 'JMAX', nc3d(2)) )
621  call check( nf90_inq_dimid(ncid, 'IMAX', nc3d(3)) )
622  call check( nf90_def_var(ncid, 'temp_c', nf90_float, nc3d, ncv) )
623  buffer = 'degC'
624  call check( nf90_put_att(ncid, ncv, 'units', trim(buffer)) )
625  buffer = 'land_ice_cold_layer_temperature'
626  call check( nf90_put_att(ncid, ncv, 'standard_name', trim(buffer)) )
627  buffer = 'Temperature in the cold ice layer'
628  call check( nf90_put_att(ncid, ncv, 'long_name', trim(buffer)) )
629  call check( nf90_put_att(ncid, ncv, 'grid_mapping', 'crs') )
630  call check( nf90_put_att(ncid, ncv, 'coordinates', 'lat lon') )
631 
632  call check( nf90_inq_dimid(ncid, 'KTMAX', nc3d(1)) )
633  call check( nf90_inq_dimid(ncid, 'JMAX', nc3d(2)) )
634  call check( nf90_inq_dimid(ncid, 'IMAX', nc3d(3)) )
635  call check( nf90_def_var(ncid, 'omega_t', nf90_float, nc3d, ncv) )
636  buffer = '1'
637  call check( nf90_put_att(ncid, ncv, 'units', trim(buffer)) )
638  buffer = 'land_ice_temperate_layer_water_content'
639  call check( nf90_put_att(ncid, ncv, 'standard_name', trim(buffer)) )
640  buffer = 'Water content in the temperate ice layer'
641  call check( nf90_put_att(ncid, ncv, 'long_name', trim(buffer)) )
642  call check( nf90_put_att(ncid, ncv, 'grid_mapping', 'crs') )
643  call check( nf90_put_att(ncid, ncv, 'coordinates', 'lat lon') )
644 
645  call check( nf90_inq_dimid(ncid, 'KRMAX', nc3d(1)) )
646  call check( nf90_inq_dimid(ncid, 'JMAX', nc3d(2)) )
647  call check( nf90_inq_dimid(ncid, 'IMAX', nc3d(3)) )
648  call check( nf90_def_var(ncid, 'temp_r', nf90_float, nc3d, ncv) )
649  buffer = 'degC'
650  call check( nf90_put_att(ncid, ncv, 'units', trim(buffer)) )
651  buffer = 'lithosphere_layer_temperature'
652  call check( nf90_put_att(ncid, ncv, 'standard_name', trim(buffer)) )
653  buffer = 'Temperature in the lithosphere layer'
654  call check( nf90_put_att(ncid, ncv, 'long_name', trim(buffer)) )
655  call check( nf90_put_att(ncid, ncv, 'grid_mapping', 'crs') )
656  call check( nf90_put_att(ncid, ncv, 'coordinates', 'lat lon') )
657 
658 end if
659 
660 call check( nf90_inq_dimid(ncid, 'JMAX', nc2d(1)) )
661 call check( nf90_inq_dimid(ncid, 'IMAX', nc2d(2)) )
662 call check( nf90_def_var(ncid, 'Q_bm', nf90_float, nc2d, ncv) )
663 buffer = 'm s-1'
664 call check( nf90_put_att(ncid, ncv, 'units', trim(buffer)) )
665 buffer = 'land_ice_basal_melt_rate'
666 call check( nf90_put_att(ncid, ncv, 'standard_name', trim(buffer)) )
667 buffer = 'Basal melting rate'
668 call check( nf90_put_att(ncid, ncv, 'long_name', trim(buffer)) )
669 call check( nf90_put_att(ncid, ncv, 'grid_mapping', 'crs') )
670 call check( nf90_put_att(ncid, ncv, 'coordinates', 'lat lon') )
671 
672 call check( nf90_inq_dimid(ncid, 'JMAX', nc2d(1)) )
673 call check( nf90_inq_dimid(ncid, 'IMAX', nc2d(2)) )
674 call check( nf90_def_var(ncid, 'Q_tld', nf90_float, nc2d, ncv) )
675 buffer = 'm s-1'
676 call check( nf90_put_att(ncid, ncv, 'units', trim(buffer)) )
677 buffer = 'land_ice_temperate_layer_water_drainage'
678 call check( nf90_put_att(ncid, ncv, 'standard_name', trim(buffer)) )
679 buffer = 'Water drainage from the temperate layer'
680 call check( nf90_put_att(ncid, ncv, 'long_name', trim(buffer)) )
681 call check( nf90_put_att(ncid, ncv, 'grid_mapping', 'crs') )
682 call check( nf90_put_att(ncid, ncv, 'coordinates', 'lat lon') )
683 
684 call check( nf90_inq_dimid(ncid, 'JMAX', nc2d(1)) )
685 call check( nf90_inq_dimid(ncid, 'IMAX', nc2d(2)) )
686 call check( nf90_def_var(ncid, 'am_perp', nf90_float, nc2d, ncv) )
687 buffer = 'm s-1'
688 call check( nf90_put_att(ncid, ncv, 'units', trim(buffer)) )
689 buffer = 'land_ice_volume_flux_across_cts'
690 call check( nf90_put_att(ncid, ncv, 'standard_name', trim(buffer)) )
691 buffer = 'Volume flux across the CTS (cold-temperate transition surface)'
692 call check( nf90_put_att(ncid, ncv, 'long_name', trim(buffer)) )
693 call check( nf90_put_att(ncid, ncv, 'grid_mapping', 'crs') )
694 call check( nf90_put_att(ncid, ncv, 'coordinates', 'lat lon') )
695 
696 call check( nf90_inq_dimid(ncid, 'JMAX', nc2d(1)) )
697 call check( nf90_inq_dimid(ncid, 'IMAX', nc2d(2)) )
698 call check( nf90_def_var(ncid, 'qx', nf90_float, nc2d, ncv) )
699 buffer = 'm2 s-1'
700 call check( nf90_put_att(ncid, ncv, 'units', trim(buffer)) )
701 buffer = 'land_ice_vertical_integral_x_velocity'
702 call check( nf90_put_att(ncid, ncv, 'standard_name', trim(buffer)) )
703 buffer = 'Horizontal volume flux qx'
704 call check( nf90_put_att(ncid, ncv, 'long_name', trim(buffer)) )
705 buffer = 'Staggered grid variable, defined at (j,i+1/2)'
706 call check( nf90_put_att(ncid, ncv, 'comment', trim(buffer)) )
707 call check( nf90_put_att(ncid, ncv, 'grid_mapping', 'crs') )
708 call check( nf90_put_att(ncid, ncv, 'coordinates', 'lat lon') )
709 
710 call check( nf90_inq_dimid(ncid, 'JMAX', nc2d(1)) )
711 call check( nf90_inq_dimid(ncid, 'IMAX', nc2d(2)) )
712 call check( nf90_def_var(ncid, 'qy', nf90_float, nc2d, ncv) )
713 buffer = 'm2 s-1'
714 call check( nf90_put_att(ncid, ncv, 'units', trim(buffer)) )
715 buffer = 'land_ice_vertical_integral_y_velocity'
716 call check( nf90_put_att(ncid, ncv, 'standard_name', trim(buffer)) )
717 buffer = 'Horizontal volume flux qy'
718 call check( nf90_put_att(ncid, ncv, 'long_name', trim(buffer)) )
719 buffer = 'Staggered grid variable, defined at (j+1/2,i)'
720 call check( nf90_put_att(ncid, ncv, 'comment', trim(buffer)) )
721 call check( nf90_put_att(ncid, ncv, 'grid_mapping', 'crs') )
722 call check( nf90_put_att(ncid, ncv, 'coordinates', 'lat lon') )
723 
724 if (flag_3d_output) then
725 
726  call check( nf90_inq_dimid(ncid, 'KCMAX', nc3d(1)) )
727  call check( nf90_inq_dimid(ncid, 'JMAX', nc3d(2)) )
728  call check( nf90_inq_dimid(ncid, 'IMAX', nc3d(3)) )
729  call check( nf90_def_var(ncid, 'age_c', nf90_float, nc3d, ncv) )
730  buffer = 's'
731  call check( nf90_put_att(ncid, ncv, 'units', trim(buffer)) )
732  buffer = 'land_ice_cold_layer_age'
733  call check( nf90_put_att(ncid, ncv, 'standard_name', trim(buffer)) )
734  buffer = 'Age in the cold ice layer'
735  call check( nf90_put_att(ncid, ncv, 'long_name', trim(buffer)) )
736  call check( nf90_put_att(ncid, ncv, 'grid_mapping', 'crs') )
737  call check( nf90_put_att(ncid, ncv, 'coordinates', 'lat lon') )
738 
739  call check( nf90_inq_dimid(ncid, 'KTMAX', nc3d(1)) )
740  call check( nf90_inq_dimid(ncid, 'JMAX', nc3d(2)) )
741  call check( nf90_inq_dimid(ncid, 'IMAX', nc3d(3)) )
742  call check( nf90_def_var(ncid, 'age_t', nf90_float, nc3d, ncv) )
743  buffer = 's'
744  call check( nf90_put_att(ncid, ncv, 'units', trim(buffer)) )
745  buffer = 'land_ice_temperate_layer_age'
746  call check( nf90_put_att(ncid, ncv, 'standard_name', trim(buffer)) )
747  buffer = 'Age in the temperate ice layer'
748  call check( nf90_put_att(ncid, ncv, 'long_name', trim(buffer)) )
749  call check( nf90_put_att(ncid, ncv, 'grid_mapping', 'crs') )
750  call check( nf90_put_att(ncid, ncv, 'coordinates', 'lat lon') )
751 
752 end if
753 
754 call check( nf90_inq_dimid(ncid, 'JMAX', nc2d(1)) )
755 call check( nf90_inq_dimid(ncid, 'IMAX', nc2d(2)) )
756 call check( nf90_def_var(ncid, 'dzs_dtau', nf90_float, nc2d, ncv) )
757 buffer = 'm s-1'
758 call check( nf90_put_att(ncid, ncv, 'units', trim(buffer)) )
759 buffer = 'tendency_of_surface_altitude'
760 call check( nf90_put_att(ncid, ncv, 'standard_name', trim(buffer)) )
761 buffer = 'Rate of change of the topography of the free surface'
762 call check( nf90_put_att(ncid, ncv, 'long_name', trim(buffer)) )
763 call check( nf90_put_att(ncid, ncv, 'grid_mapping', 'crs') )
764 call check( nf90_put_att(ncid, ncv, 'coordinates', 'lat lon') )
765 
766 call check( nf90_inq_dimid(ncid, 'JMAX', nc2d(1)) )
767 call check( nf90_inq_dimid(ncid, 'IMAX', nc2d(2)) )
768 call check( nf90_def_var(ncid, 'dzm_dtau', nf90_float, nc2d, ncv) )
769 buffer = 'm s-1'
770 call check( nf90_put_att(ncid, ncv, 'units', trim(buffer)) )
771 buffer = 'tendency_of_cts_altitude'
772 call check( nf90_put_att(ncid, ncv, 'standard_name', trim(buffer)) )
773 buffer = 'Rate of change of the topography of the CTS '// &
774  '(cold-temperate transition surface)'
775 call check( nf90_put_att(ncid, ncv, 'long_name', trim(buffer)) )
776 call check( nf90_put_att(ncid, ncv, 'grid_mapping', 'crs') )
777 call check( nf90_put_att(ncid, ncv, 'coordinates', 'lat lon') )
778 
779 call check( nf90_inq_dimid(ncid, 'JMAX', nc2d(1)) )
780 call check( nf90_inq_dimid(ncid, 'IMAX', nc2d(2)) )
781 call check( nf90_def_var(ncid, 'dzb_dtau', nf90_float, nc2d, ncv) )
782 buffer = 'm s-1'
783 call check( nf90_put_att(ncid, ncv, 'units', trim(buffer)) )
784 buffer = 'tendency_of_ice_base_altitude'
785 call check( nf90_put_att(ncid, ncv, 'standard_name', trim(buffer)) )
786 buffer = 'Rate of change of the topography of the ice base'
787 call check( nf90_put_att(ncid, ncv, 'long_name', trim(buffer)) )
788 call check( nf90_put_att(ncid, ncv, 'grid_mapping', 'crs') )
789 call check( nf90_put_att(ncid, ncv, 'coordinates', 'lat lon') )
790 
791 call check( nf90_inq_dimid(ncid, 'JMAX', nc2d(1)) )
792 call check( nf90_inq_dimid(ncid, 'IMAX', nc2d(2)) )
793 call check( nf90_def_var(ncid, 'dzl_dtau', nf90_float, nc2d, ncv) )
794 buffer = 'm s-1'
795 call check( nf90_put_att(ncid, ncv, 'units', trim(buffer)) )
796 buffer = 'tendency_of_bedrock_altitude'
797 call check( nf90_put_att(ncid, ncv, 'standard_name', trim(buffer)) )
798 buffer = 'Rate of change of the topography of the lithosphere surface'
799 call check( nf90_put_att(ncid, ncv, 'long_name', trim(buffer)) )
800 call check( nf90_put_att(ncid, ncv, 'grid_mapping', 'crs') )
801 call check( nf90_put_att(ncid, ncv, 'coordinates', 'lat lon') )
802 
803 call check( nf90_inq_dimid(ncid, 'JMAX', nc2d(1)) )
804 call check( nf90_inq_dimid(ncid, 'IMAX', nc2d(2)) )
805 call check( nf90_def_var(ncid, 'dH_c_dtau', nf90_float, nc2d, ncv) )
806 buffer = 'm s-1'
807 call check( nf90_put_att(ncid, ncv, 'units', trim(buffer)) )
808 buffer = 'tendency_of_land_ice_cold_layer_thickness'
809 call check( nf90_put_att(ncid, ncv, 'standard_name', trim(buffer)) )
810 buffer = 'Rate of change of the thickness of the cold ice layer'
811 call check( nf90_put_att(ncid, ncv, 'long_name', trim(buffer)) )
812 call check( nf90_put_att(ncid, ncv, 'grid_mapping', 'crs') )
813 call check( nf90_put_att(ncid, ncv, 'coordinates', 'lat lon') )
814 
815 call check( nf90_inq_dimid(ncid, 'JMAX', nc2d(1)) )
816 call check( nf90_inq_dimid(ncid, 'IMAX', nc2d(2)) )
817 call check( nf90_def_var(ncid, 'dH_t_dtau', nf90_float, nc2d, ncv) )
818 buffer = 'm s-1'
819 call check( nf90_put_att(ncid, ncv, 'units', trim(buffer)) )
820 buffer = 'tendency_of_land_ice_temperate_layer_thickness'
821 call check( nf90_put_att(ncid, ncv, 'standard_name', trim(buffer)) )
822 buffer = 'Rate of change of the thickness of the temperate ice layer'
823 call check( nf90_put_att(ncid, ncv, 'long_name', trim(buffer)) )
824 call check( nf90_put_att(ncid, ncv, 'grid_mapping', 'crs') )
825 call check( nf90_put_att(ncid, ncv, 'coordinates', 'lat lon') )
826 
827 call check( nf90_inq_dimid(ncid, 'JMAX', nc2d(1)) )
828 call check( nf90_inq_dimid(ncid, 'IMAX', nc2d(2)) )
829 call check( nf90_def_var(ncid, 'dH_dtau', nf90_float, nc2d, ncv) )
830 buffer = 'm s-1'
831 call check( nf90_put_att(ncid, ncv, 'units', trim(buffer)) )
832 buffer = 'tendency_of_land_ice_thickness'
833 call check( nf90_put_att(ncid, ncv, 'standard_name', trim(buffer)) )
834 buffer = 'Rate of change of the ice thickness'
835 call check( nf90_put_att(ncid, ncv, 'long_name', trim(buffer)) )
836 call check( nf90_put_att(ncid, ncv, 'grid_mapping', 'crs') )
837 call check( nf90_put_att(ncid, ncv, 'coordinates', 'lat lon') )
838 
839 call check( nf90_inq_dimid(ncid, 'JMAX', nc2d(1)) )
840 call check( nf90_inq_dimid(ncid, 'IMAX', nc2d(2)) )
841 call check( nf90_def_var(ncid, 'vx_b_g', nf90_float, nc2d, ncv) )
842 buffer = 'm s-1'
843 call check( nf90_put_att(ncid, ncv, 'units', trim(buffer)) )
844 buffer = 'land_ice_base_x_velocity'
845 call check( nf90_put_att(ncid, ncv, 'standard_name', trim(buffer)) )
846 buffer = 'Horizontal velocity vx at the ice base'
847 call check( nf90_put_att(ncid, ncv, 'long_name', trim(buffer)) )
848 call check( nf90_put_att(ncid, ncv, 'grid_mapping', 'crs') )
849 call check( nf90_put_att(ncid, ncv, 'coordinates', 'lat lon') )
850 
851 call check( nf90_inq_dimid(ncid, 'JMAX', nc2d(1)) )
852 call check( nf90_inq_dimid(ncid, 'IMAX', nc2d(2)) )
853 call check( nf90_def_var(ncid, 'vy_b_g', nf90_float, nc2d, ncv) )
854 buffer = 'm s-1'
855 call check( nf90_put_att(ncid, ncv, 'units', trim(buffer)) )
856 buffer = 'land_ice_base_y_velocity'
857 call check( nf90_put_att(ncid, ncv, 'standard_name', trim(buffer)) )
858 buffer = 'Horizontal velocity vy at the ice base'
859 call check( nf90_put_att(ncid, ncv, 'long_name', trim(buffer)) )
860 call check( nf90_put_att(ncid, ncv, 'grid_mapping', 'crs') )
861 call check( nf90_put_att(ncid, ncv, 'coordinates', 'lat lon') )
862 
863 call check( nf90_inq_dimid(ncid, 'JMAX', nc2d(1)) )
864 call check( nf90_inq_dimid(ncid, 'IMAX', nc2d(2)) )
865 call check( nf90_def_var(ncid, 'vz_b', nf90_float, nc2d, ncv) )
866 buffer = 'm s-1'
867 call check( nf90_put_att(ncid, ncv, 'units', trim(buffer)) )
868 buffer = 'land_ice_base_z_velocity'
869 call check( nf90_put_att(ncid, ncv, 'standard_name', trim(buffer)) )
870 buffer = 'Vertical velocity vz at the ice base'
871 call check( nf90_put_att(ncid, ncv, 'long_name', trim(buffer)) )
872 call check( nf90_put_att(ncid, ncv, 'grid_mapping', 'crs') )
873 call check( nf90_put_att(ncid, ncv, 'coordinates', 'lat lon') )
874 
875 call check( nf90_inq_dimid(ncid, 'JMAX', nc2d(1)) )
876 call check( nf90_inq_dimid(ncid, 'IMAX', nc2d(2)) )
877 call check( nf90_def_var(ncid, 'vx_s_g', nf90_float, nc2d, ncv) )
878 buffer = 'm s-1'
879 call check( nf90_put_att(ncid, ncv, 'units', trim(buffer)) )
880 buffer = 'land_ice_surface_x_velocity'
881 call check( nf90_put_att(ncid, ncv, 'standard_name', trim(buffer)) )
882 buffer = 'Horizontal velocity vx at the ice surface'
883 call check( nf90_put_att(ncid, ncv, 'long_name', trim(buffer)) )
884 call check( nf90_put_att(ncid, ncv, 'grid_mapping', 'crs') )
885 call check( nf90_put_att(ncid, ncv, 'coordinates', 'lat lon') )
886 
887 call check( nf90_inq_dimid(ncid, 'JMAX', nc2d(1)) )
888 call check( nf90_inq_dimid(ncid, 'IMAX', nc2d(2)) )
889 call check( nf90_def_var(ncid, 'vy_s_g', nf90_float, nc2d, ncv) )
890 buffer = 'm s-1'
891 call check( nf90_put_att(ncid, ncv, 'units', trim(buffer)) )
892 buffer = 'land_ice_surface_y_velocity'
893 call check( nf90_put_att(ncid, ncv, 'standard_name', trim(buffer)) )
894 buffer = 'Horizontal velocity vy at the ice surface'
895 call check( nf90_put_att(ncid, ncv, 'long_name', trim(buffer)) )
896 call check( nf90_put_att(ncid, ncv, 'grid_mapping', 'crs') )
897 call check( nf90_put_att(ncid, ncv, 'coordinates', 'lat lon') )
898 
899 call check( nf90_inq_dimid(ncid, 'JMAX', nc2d(1)) )
900 call check( nf90_inq_dimid(ncid, 'IMAX', nc2d(2)) )
901 call check( nf90_def_var(ncid, 'vz_s', nf90_float, nc2d, ncv) )
902 buffer = 'm s-1'
903 call check( nf90_put_att(ncid, ncv, 'units', trim(buffer)) )
904 buffer = 'land_ice_surface_z_velocity'
905 call check( nf90_put_att(ncid, ncv, 'standard_name', trim(buffer)) )
906 buffer = 'Vertical velocity vz at the ice surface'
907 call check( nf90_put_att(ncid, ncv, 'long_name', trim(buffer)) )
908 call check( nf90_put_att(ncid, ncv, 'grid_mapping', 'crs') )
909 call check( nf90_put_att(ncid, ncv, 'coordinates', 'lat lon') )
910 
911 call check( nf90_inq_dimid(ncid, 'JMAX', nc2d(1)) )
912 call check( nf90_inq_dimid(ncid, 'IMAX', nc2d(2)) )
913 call check( nf90_def_var(ncid, 'temp_b', nf90_float, nc2d, ncv) )
914 buffer = 'degC'
915 call check( nf90_put_att(ncid, ncv, 'units', trim(buffer)) )
916 buffer = 'basal_temperature'
917 call check( nf90_put_att(ncid, ncv, 'standard_name', trim(buffer)) )
918 buffer = 'Temperature at the ice base'
919 call check( nf90_put_att(ncid, ncv, 'long_name', trim(buffer)) )
920 call check( nf90_put_att(ncid, ncv, 'grid_mapping', 'crs') )
921 call check( nf90_put_att(ncid, ncv, 'coordinates', 'lat lon') )
922 
923 call check( nf90_inq_dimid(ncid, 'JMAX', nc2d(1)) )
924 call check( nf90_inq_dimid(ncid, 'IMAX', nc2d(2)) )
925 call check( nf90_def_var(ncid, 'temph_b', nf90_float, nc2d, ncv) )
926 buffer = 'degC'
927 call check( nf90_put_att(ncid, ncv, 'units', trim(buffer)) )
928 buffer = 'basal_temperature_rel_to_pmp'
929 call check( nf90_put_att(ncid, ncv, 'standard_name', trim(buffer)) )
930 buffer = 'Temperature at the ice base relative to the pressure melting point'
931 call check( nf90_put_att(ncid, ncv, 'long_name', trim(buffer)) )
932 call check( nf90_put_att(ncid, ncv, 'grid_mapping', 'crs') )
933 call check( nf90_put_att(ncid, ncv, 'coordinates', 'lat lon') )
934 
935 call check( nf90_inq_dimid(ncid, 'JMAX', nc2d(1)) )
936 call check( nf90_inq_dimid(ncid, 'IMAX', nc2d(2)) )
937 call check( nf90_def_var(ncid, 'p_b_w', nf90_float, nc2d, ncv) )
938 buffer = 'Pa'
939 call check( nf90_put_att(ncid, ncv, 'units', trim(buffer)) )
940 buffer = 'basal_water_pressure'
941 call check( nf90_put_att(ncid, ncv, 'standard_name', trim(buffer)) )
942 buffer = 'Basal water pressure'
943 call check( nf90_put_att(ncid, ncv, 'long_name', trim(buffer)) )
944 call check( nf90_put_att(ncid, ncv, 'grid_mapping', 'crs') )
945 call check( nf90_put_att(ncid, ncv, 'coordinates', 'lat lon') )
946 
947 call check( nf90_inq_dimid(ncid, 'JMAX', nc2d(1)) )
948 call check( nf90_inq_dimid(ncid, 'IMAX', nc2d(2)) )
949 call check( nf90_def_var(ncid, 'H_w', nf90_float, nc2d, ncv) )
950 buffer = 'm'
951 call check( nf90_put_att(ncid, ncv, 'units', trim(buffer)) )
952 buffer = 'subglacial_water_thickness'
953 call check( nf90_put_att(ncid, ncv, 'standard_name', trim(buffer)) )
954 buffer = 'Effective thickness of subglacial water'
955 call check( nf90_put_att(ncid, ncv, 'long_name', trim(buffer)) )
956 call check( nf90_put_att(ncid, ncv, 'grid_mapping', 'crs') )
957 call check( nf90_put_att(ncid, ncv, 'coordinates', 'lat lon') )
958 
959 call check( nf90_inq_dimid(ncid, 'JMAX', nc2d(1)) )
960 call check( nf90_inq_dimid(ncid, 'IMAX', nc2d(2)) )
961 call check( nf90_def_var(ncid, 'q_gl_g', nf90_float, nc2d, ncv) )
962 buffer = 'm2 s-1'
963 call check( nf90_put_att(ncid, ncv, 'units', trim(buffer)) )
964 buffer = 'land_ice_volume_flux_across_gl'
965 call check( nf90_put_att(ncid, ncv, 'standard_name', trim(buffer)) )
966 buffer = 'Horizontal volume flux across the grounding line'
967 call check( nf90_put_att(ncid, ncv, 'long_name', trim(buffer)) )
968 call check( nf90_put_att(ncid, ncv, 'grid_mapping', 'crs') )
969 call check( nf90_put_att(ncid, ncv, 'coordinates', 'lat lon') )
970 
971 call check( nf90_def_var(ncid, 'V_tot', nf90_float, ncv) )
972 buffer = 'm3'
973 call check( nf90_put_att(ncid, ncv, 'units', trim(buffer)) )
974 buffer = 'land_ice_volume'
975 call check( nf90_put_att(ncid, ncv, 'standard_name', trim(buffer)) )
976 buffer = 'Ice volume'
977 call check( nf90_put_att(ncid, ncv, 'long_name', trim(buffer)) )
978 
979 call check( nf90_def_var(ncid, 'A_grounded', nf90_float, ncv) )
980 buffer = 'm2'
981 call check( nf90_put_att(ncid, ncv, 'units', trim(buffer)) )
982 buffer = 'land_ice_area_grounded'
983 call check( nf90_put_att(ncid, ncv, 'standard_name', trim(buffer)) )
984 buffer = 'Area covered by grounded ice'
985 call check( nf90_put_att(ncid, ncv, 'long_name', trim(buffer)) )
986 
987 call check( nf90_def_var(ncid, 'A_floating', nf90_float, ncv) )
988 buffer = 'm2'
989 call check( nf90_put_att(ncid, ncv, 'units', trim(buffer)) )
990 buffer = 'land_ice_area_floating'
991 call check( nf90_put_att(ncid, ncv, 'standard_name', trim(buffer)) )
992 buffer = 'Area covered by floating ice'
993 call check( nf90_put_att(ncid, ncv, 'long_name', trim(buffer)) )
994 
995 call check( nf90_enddef(ncid) )
996 
997 !-------- Computation of the ice volume,
998 ! the area covered by grounded ice
999 ! and the area covered by floating ice --------
1000 
1001 v_tot = 0.0_dp
1002 a_grounded = 0.0_dp
1003 a_floating = 0.0_dp
1004 
1005 do i=0, imax
1006 do j=0, jmax
1007 
1008  if ( (maske(j,i)==0).or.(maske(j,i)==3) ) & ! grounded or floating ice
1009  v_tot = v_tot + (h_c(j,i)+h_t(j,i))*area(j,i)
1010 
1011  if (maske(j,i)==0) & ! grounded ice
1012  a_grounded = a_grounded + area(j,i)
1013 
1014  if (maske(j,i)==3) & ! floating ice
1015  a_floating = a_floating + area(j,i)
1016 
1017 end do
1018 end do
1019 
1020 !-------- Convert data to real*4 --------
1021 
1022 #if ( !defined(OUT_TIMES) || OUT_TIMES==1 )
1023 time_conv = real(time,sp)
1024 #elif OUT_TIMES == 2
1025 time_conv = real((time+year_zero),sp)
1026 #else
1027 stop ' output_nc: OUT_TIMES must be either 1 or 2!'
1028 #endif
1029 
1030 delta_ts_conv = real(delta_ts,sp)
1031 glac_index_conv = real(glac_index,sp)
1032 z_sl_conv = real(z_sl,sp)
1033 v_tot_conv = real(v_tot,sp)
1034 a_grounded_conv = real(a_grounded,sp)
1035 a_floating_conv = real(a_floating,sp)
1036 h_r_conv = real(h_r,sp)
1037 
1038 do i=0, imax
1039  xi_conv(i) = real(xi(i),sp)
1040 end do
1041 
1042 do j=0, jmax
1043  eta_conv(j) = real(eta(j),sp)
1044 end do
1045 
1046 do kc=0, kcmax
1047  sigma_level_c_conv(kc) = real(eaz_c_quotient(kc),sp)
1048 end do
1049 
1050 do kt=0, ktmax
1051  sigma_level_t_conv(kt) = real(zeta_t(kt),sp)
1052 end do
1053 
1054 do kr=0, krmax
1055  sigma_level_r_conv(kr) = real(kr,sp)/real(krmax,sp)
1056 end do
1057 
1058 do i=0, imax
1059 do j=0, jmax
1060 
1061  lambda_conv(j,i) = real(lambda(j,i),sp)
1062  phi_conv(j,i) = real(phi(j,i),sp)
1063  lon_conv(j,i) = real(lambda(j,i)*pi_180_inv,sp) ! longitude in deg
1064  lon_conv(j,i) = modulo(lon_conv(j,i)+180.0_sp, 360.0_sp)-180.0_sp
1065  ! mapping to interval [-180 deg, +180 deg)
1066  lat_conv(j,i) = real(phi(j,i) *pi_180_inv,sp) ! latitute in deg
1067  if (lat_conv(j,i) > 90.0_sp) lat_conv(j,i) = 90.0_sp
1068  if (lat_conv(j,i) < -90.0_sp) lat_conv(j,i) = -90.0_sp
1069  ! constraining to interval [-90 deg, +90 deg]
1070  temp_s_conv(j,i) = real(temp_s(j,i),sp)
1071  as_perp_conv(j,i) = real(as_perp(j,i),sp)
1072  zs_conv(j,i) = real(zs(j,i),sp)
1073  zm_conv(j,i) = real(zm(j,i),sp)
1074  zb_conv(j,i) = real(zb(j,i),sp)
1075  zl_conv(j,i) = real(zl(j,i),sp)
1076  h_c_conv(j,i) = real(H_c(j,i),sp)
1077  h_t_conv(j,i) = real(H_t(j,i),sp)
1078  h_conv(j,i) = real(H_c(j,i)+H_t(j,i),sp)
1079  q_bm_conv(j,i) = real(Q_bm(j,i),sp)
1080  q_tld_conv(j,i) = real(Q_tld(j,i),sp)
1081  am_perp_conv(j,i) = real(am_perp(j,i),sp)
1082  qx_conv(j,i) = real(qx(j,i),sp)
1083  qy_conv(j,i) = real(qy(j,i),sp)
1084  dzs_dtau_conv(j,i) = real(dzs_dtau(j,i),sp)
1085  dzm_dtau_conv(j,i) = real(dzm_dtau(j,i),sp)
1086  dzb_dtau_conv(j,i) = real(dzb_dtau(j,i),sp)
1087  dzl_dtau_conv(j,i) = real(dzl_dtau(j,i),sp)
1088  dh_c_dtau_conv(j,i) = real(dH_c_dtau(j,i),sp)
1089  dh_t_dtau_conv(j,i) = real(dH_t_dtau(j,i),sp)
1090  dh_dtau_conv(j,i) = real(dH_c_dtau(j,i)+dH_t_dtau(j,i),sp)
1091  vx_b_g_conv(j,i) = real(vx_b_g(j,i),sp)
1092  vy_b_g_conv(j,i) = real(vy_b_g(j,i),sp)
1093  vz_b_conv(j,i) = real(vz_b(j,i),sp)
1094  vx_s_g_conv(j,i) = real(vx_s_g(j,i),sp)
1095  vy_s_g_conv(j,i) = real(vy_s_g(j,i),sp)
1096  vz_s_conv(j,i) = real(vz_s(j,i),sp)
1097  temp_b_conv(j,i) = real(temp_b(j,i),sp)
1098  temph_b_conv(j,i) = real(temph_b(j,i),sp)
1099  p_b_w_conv(j,i) = real(p_b_w(j,i),sp)
1100  h_w_conv(j,i) = real(H_w(j,i),sp)
1101  q_gl_g_conv(j,i) = real(q_gl_g(j,i),sp)
1102 
1103  do kr=0, krmax
1104  temp_r_conv(kr,j,i) = real(temp_r(kr,j,i),sp)
1105  end do
1106 
1107  do kt=0, ktmax
1108  vx_t_conv(kt,j,i) = real(vx_t(kt,j,i),sp)
1109  vy_t_conv(kt,j,i) = real(vy_t(kt,j,i),sp)
1110  vz_t_conv(kt,j,i) = real(vz_t(kt,j,i),sp)
1111  omega_t_conv(kt,j,i) = real(omega_t(kt,j,i),sp)
1112  age_t_conv(kt,j,i) = real(age_t(kt,j,i),sp)
1113  end do
1114 
1115  do kc=0, kcmax
1116  vx_c_conv(kc,j,i) = real(vx_c(kc,j,i),sp)
1117  vy_c_conv(kc,j,i) = real(vy_c(kc,j,i),sp)
1118  vz_c_conv(kc,j,i) = real(vz_c(kc,j,i),sp)
1119  temp_c_conv(kc,j,i) = real(temp_c(kc,j,i),sp)
1120  age_c_conv(kc,j,i) = real(age_c(kc,j,i),sp)
1121  end do
1122 
1123 end do
1124 end do
1125 
1126 !-------- Write data on file --------
1127 
1128 call check( nf90_inq_varid(ncid, 'crs', ncv) )
1129 call check( nf90_put_var(ncid, ncv, 0) )
1130 
1131 call check( nf90_inq_varid(ncid, 'time', ncv) )
1132 call check( nf90_put_var(ncid, ncv, time_conv) )
1133 
1134 if (forcing_flag == 1) then
1135 
1136  call check( nf90_inq_varid(ncid, 'delta_ts', ncv) )
1137  call check( nf90_put_var(ncid, ncv, delta_ts_conv) )
1138 
1139 else if (forcing_flag == 2) then
1140 
1141  call check( nf90_inq_varid(ncid, 'glac_index', ncv) )
1142  call check( nf90_put_var(ncid, ncv, glac_index_conv) )
1143 
1144 else if (forcing_flag == 3) then
1145 
1146  call check( nf90_inq_varid(ncid, 'glac_index', ncv) )
1147  glac_index_conv = 1.11e+11 ! dummy value
1148  call check( nf90_put_var(ncid, ncv, glac_index_conv) )
1149 
1150 end if
1151 
1152 call check( nf90_inq_varid(ncid, 'z_sl', ncv) )
1153 call check( nf90_put_var(ncid, ncv, z_sl_conv) )
1154 
1155 call check( nf90_inq_varid(ncid, 'xi', ncv) )
1156 nc1cor(1) = 1
1157 nc1cnt(1) = imax + 1
1158 call check( nf90_put_var(ncid, ncv, xi_conv, start=nc1cor, count=nc1cnt) )
1159 
1160 call check( nf90_inq_varid(ncid, 'eta', ncv) )
1161 nc1cor(1) = 1
1162 nc1cnt(1) = jmax + 1
1163 call check( nf90_put_var(ncid, ncv, eta_conv, start=nc1cor, count=nc1cnt) )
1164 
1165 call check( nf90_inq_varid(ncid, 'sigma_level_c', ncv) )
1166 nc1cor(1) = 1
1167 nc1cnt(1) = kcmax + 1
1168 call check( nf90_put_var(ncid, ncv, sigma_level_c_conv, start=nc1cor, count=nc1cnt) )
1169 
1170 call check( nf90_inq_varid(ncid, 'sigma_level_t', ncv) )
1171 nc1cor(1) = 1
1172 nc1cnt(1) = ktmax + 1
1173 call check( nf90_put_var(ncid, ncv, sigma_level_t_conv, start=nc1cor, count=nc1cnt) )
1174 
1175 call check( nf90_inq_varid(ncid, 'sigma_level_r', ncv) )
1176 nc1cor(1) = 1
1177 nc1cnt(1) = krmax + 1
1178 call check( nf90_put_var(ncid, ncv, sigma_level_r_conv, start=nc1cor, count=nc1cnt) )
1179 
1180 call check( nf90_inq_varid(ncid, 'lon', ncv) )
1181 nc2cor(1) = 1
1182 nc2cor(2) = 1
1183 nc2cnt(1) = jmax + 1
1184 nc2cnt(2) = imax + 1
1185 call check( nf90_put_var(ncid, ncv, lon_conv, start=nc2cor, count=nc2cnt) )
1186 
1187 call check( nf90_inq_varid(ncid, 'lat', ncv) )
1188 nc2cor(1) = 1
1189 nc2cor(2) = 1
1190 nc2cnt(1) = jmax + 1
1191 nc2cnt(2) = imax + 1
1192 call check( nf90_put_var(ncid, ncv, lat_conv, start=nc2cor, count=nc2cnt) )
1193 
1194 call check( nf90_inq_varid(ncid, 'lambda', ncv) )
1195 nc2cor(1) = 1
1196 nc2cor(2) = 1
1197 nc2cnt(1) = jmax + 1
1198 nc2cnt(2) = imax + 1
1199 call check( nf90_put_var(ncid, ncv, lambda_conv, start=nc2cor, count=nc2cnt) )
1200 
1201 call check( nf90_inq_varid(ncid, 'phi', ncv) )
1202 nc2cor(1) = 1
1203 nc2cor(2) = 1
1204 nc2cnt(1) = jmax + 1
1205 nc2cnt(2) = imax + 1
1206 call check( nf90_put_var(ncid, ncv, phi_conv, start=nc2cor, count=nc2cnt) )
1207 
1208 call check( nf90_inq_varid(ncid, 'temp_s', ncv) )
1209 nc2cor(1) = 1
1210 nc2cor(2) = 1
1211 nc2cnt(1) = jmax + 1
1212 nc2cnt(2) = imax + 1
1213 call check( nf90_put_var(ncid, ncv, temp_s_conv, start=nc2cor, count=nc2cnt) )
1214 
1215 call check( nf90_inq_varid(ncid, 'as_perp', ncv) )
1216 nc2cor(1) = 1
1217 nc2cor(2) = 1
1218 nc2cnt(1) = jmax + 1
1219 nc2cnt(2) = imax + 1
1220 call check( nf90_put_var(ncid, ncv, as_perp_conv, start=nc2cor, count=nc2cnt) )
1221 
1222 call check( nf90_inq_varid(ncid, 'maske', ncv) )
1223 nc2cor(1) = 1
1224 nc2cor(2) = 1
1225 nc2cnt(1) = jmax + 1
1226 nc2cnt(2) = imax + 1
1227 call check( nf90_put_var(ncid, ncv, maske, start=nc2cor, count=nc2cnt) )
1228 
1229 call check( nf90_inq_varid(ncid, 'n_cts', ncv) )
1230 nc2cor(1) = 1
1231 nc2cor(2) = 1
1232 nc2cnt(1) = jmax + 1
1233 nc2cnt(2) = imax + 1
1234 call check( nf90_put_var(ncid, ncv, n_cts, start=nc2cor, count=nc2cnt) )
1235 
1236 call check( nf90_inq_varid(ncid, 'zs', ncv) )
1237 nc2cor(1) = 1
1238 nc2cor(2) = 1
1239 nc2cnt(1) = jmax + 1
1240 nc2cnt(2) = imax + 1
1241 call check( nf90_put_var(ncid, ncv, zs_conv, start=nc2cor, count=nc2cnt) )
1242 
1243 call check( nf90_inq_varid(ncid, 'zm', ncv) )
1244 nc2cor(1) = 1
1245 nc2cor(2) = 1
1246 nc2cnt(1) = jmax + 1
1247 nc2cnt(2) = imax + 1
1248 call check( nf90_put_var(ncid, ncv, zm_conv, start=nc2cor, count=nc2cnt) )
1249 
1250 call check( nf90_inq_varid(ncid, 'zb', ncv) )
1251 nc2cor(1) = 1
1252 nc2cor(2) = 1
1253 nc2cnt(1) = jmax + 1
1254 nc2cnt(2) = imax + 1
1255 call check( nf90_put_var(ncid, ncv, zb_conv, start=nc2cor, count=nc2cnt) )
1256 
1257 call check( nf90_inq_varid(ncid, 'zl', ncv) )
1258 nc2cor(1) = 1
1259 nc2cor(2) = 1
1260 nc2cnt(1) = jmax + 1
1261 nc2cnt(2) = imax + 1
1262 call check( nf90_put_var(ncid, ncv, zl_conv, start=nc2cor, count=nc2cnt) )
1263 
1264 call check( nf90_inq_varid(ncid, 'H_c', ncv) )
1265 nc2cor(1) = 1
1266 nc2cor(2) = 1
1267 nc2cnt(1) = jmax + 1
1268 nc2cnt(2) = imax + 1
1269 call check( nf90_put_var(ncid, ncv, h_c_conv, start=nc2cor, count=nc2cnt) )
1270 
1271 call check( nf90_inq_varid(ncid, 'H_t', ncv) )
1272 nc2cor(1) = 1
1273 nc2cor(2) = 1
1274 nc2cnt(1) = jmax + 1
1275 nc2cnt(2) = imax + 1
1276 call check( nf90_put_var(ncid, ncv, h_t_conv, start=nc2cor, count=nc2cnt) )
1277 
1278 call check( nf90_inq_varid(ncid, 'H', ncv) )
1279 nc2cor(1) = 1
1280 nc2cor(2) = 1
1281 nc2cnt(1) = jmax + 1
1282 nc2cnt(2) = imax + 1
1283 call check( nf90_put_var(ncid, ncv, h_conv, start=nc2cor, count=nc2cnt) )
1284 
1285 call check( nf90_inq_varid(ncid, 'H_R', ncv) )
1286 call check( nf90_put_var(ncid, ncv, h_r_conv) )
1287 
1288 if (flag_3d_output) then
1289 
1290  call check( nf90_inq_varid(ncid, 'vx_c', ncv) )
1291  nc3cor(1) = 1
1292  nc3cor(2) = 1
1293  nc3cor(3) = 1
1294  nc3cnt(1) = kcmax + 1
1295  nc3cnt(2) = jmax + 1
1296  nc3cnt(3) = imax + 1
1297  call check( nf90_put_var(ncid, ncv, vx_c_conv, start=nc3cor, count=nc3cnt) )
1298 
1299  call check( nf90_inq_varid(ncid, 'vy_c', ncv) )
1300  nc3cor(1) = 1
1301  nc3cor(2) = 1
1302  nc3cor(3) = 1
1303  nc3cnt(1) = kcmax + 1
1304  nc3cnt(2) = jmax + 1
1305  nc3cnt(3) = imax + 1
1306  call check( nf90_put_var(ncid, ncv, vy_c_conv, start=nc3cor, count=nc3cnt) )
1307 
1308  call check( nf90_inq_varid(ncid, 'vz_c', ncv) )
1309  nc3cor(1) = 1
1310  nc3cor(2) = 1
1311  nc3cor(3) = 1
1312  nc3cnt(1) = kcmax + 1
1313  nc3cnt(2) = jmax + 1
1314  nc3cnt(3) = imax + 1
1315  call check( nf90_put_var(ncid, ncv, vz_c_conv, start=nc3cor, count=nc3cnt) )
1316 
1317  call check( nf90_inq_varid(ncid, 'vx_t', ncv) )
1318  nc3cor(1) = 1
1319  nc3cor(2) = 1
1320  nc3cor(3) = 1
1321  nc3cnt(1) = ktmax + 1
1322  nc3cnt(2) = jmax + 1
1323  nc3cnt(3) = imax + 1
1324  call check( nf90_put_var(ncid, ncv, vx_t_conv, start=nc3cor, count=nc3cnt) )
1325 
1326  call check( nf90_inq_varid(ncid, 'vy_t', ncv) )
1327  nc3cor(1) = 1
1328  nc3cor(2) = 1
1329  nc3cor(3) = 1
1330  nc3cnt(1) = ktmax + 1
1331  nc3cnt(2) = jmax + 1
1332  nc3cnt(3) = imax + 1
1333  call check( nf90_put_var(ncid, ncv, vy_t_conv, start=nc3cor, count=nc3cnt) )
1334 
1335  call check( nf90_inq_varid(ncid, 'vz_t', ncv) )
1336  nc3cor(1) = 1
1337  nc3cor(2) = 1
1338  nc3cor(3) = 1
1339  nc3cnt(1) = ktmax + 1
1340  nc3cnt(2) = jmax + 1
1341  nc3cnt(3) = imax + 1
1342  call check( nf90_put_var(ncid, ncv, vz_t_conv, start=nc3cor, count=nc3cnt) )
1343 
1344  call check( nf90_inq_varid(ncid, 'temp_c', ncv) )
1345  nc3cor(1) = 1
1346  nc3cor(2) = 1
1347  nc3cor(3) = 1
1348  nc3cnt(1) = kcmax + 1
1349  nc3cnt(2) = jmax + 1
1350  nc3cnt(3) = imax + 1
1351  call check( nf90_put_var(ncid, ncv, temp_c_conv, start=nc3cor, count=nc3cnt) )
1352 
1353  call check( nf90_inq_varid(ncid, 'omega_t', ncv) )
1354  nc3cor(1) = 1
1355  nc3cor(2) = 1
1356  nc3cor(3) = 1
1357  nc3cnt(1) = ktmax + 1
1358  nc3cnt(2) = jmax + 1
1359  nc3cnt(3) = imax + 1
1360  call check( nf90_put_var(ncid, ncv, omega_t_conv, start=nc3cor, count=nc3cnt) )
1361 
1362  call check( nf90_inq_varid(ncid, 'temp_r', ncv) )
1363  nc3cor(1) = 1
1364  nc3cor(2) = 1
1365  nc3cor(3) = 1
1366  nc3cnt(1) = krmax + 1
1367  nc3cnt(2) = jmax + 1
1368  nc3cnt(3) = imax + 1
1369  call check( nf90_put_var(ncid, ncv, temp_r_conv, start=nc3cor, count=nc3cnt) )
1370 
1371 end if
1372 
1373 call check( nf90_inq_varid(ncid, 'Q_bm', ncv) )
1374 nc2cor(1) = 1
1375 nc2cor(2) = 1
1376 nc2cnt(1) = jmax + 1
1377 nc2cnt(2) = imax + 1
1378 call check( nf90_put_var(ncid, ncv, q_bm_conv, start=nc2cor, count=nc2cnt) )
1379 
1380 call check( nf90_inq_varid(ncid, 'Q_tld', ncv) )
1381 nc2cor(1) = 1
1382 nc2cor(2) = 1
1383 nc2cnt(1) = jmax + 1
1384 nc2cnt(2) = imax + 1
1385 call check( nf90_put_var(ncid, ncv, q_tld_conv, start=nc2cor, count=nc2cnt) )
1386 
1387 call check( nf90_inq_varid(ncid, 'am_perp', ncv) )
1388 nc2cor(1) = 1
1389 nc2cor(2) = 1
1390 nc2cnt(1) = jmax + 1
1391 nc2cnt(2) = imax + 1
1392 call check( nf90_put_var(ncid, ncv, am_perp_conv, start=nc2cor, count=nc2cnt) )
1393 
1394 call check( nf90_inq_varid(ncid, 'qx', ncv) )
1395 nc2cor(1) = 1
1396 nc2cor(2) = 1
1397 nc2cnt(1) = jmax + 1
1398 nc2cnt(2) = imax + 1
1399 call check( nf90_put_var(ncid, ncv, qx_conv, start=nc2cor, count=nc2cnt) )
1400 
1401 call check( nf90_inq_varid(ncid, 'qy', ncv) )
1402 nc2cor(1) = 1
1403 nc2cor(2) = 1
1404 nc2cnt(1) = jmax + 1
1405 nc2cnt(2) = imax + 1
1406 call check( nf90_put_var(ncid, ncv, qy_conv, start=nc2cor, count=nc2cnt) )
1407 
1408 if (flag_3d_output) then
1409 
1410  call check( nf90_inq_varid(ncid, 'age_c', ncv) )
1411  nc3cor(1) = 1
1412  nc3cor(2) = 1
1413  nc3cor(3) = 1
1414  nc3cnt(1) = kcmax + 1
1415  nc3cnt(2) = jmax + 1
1416  nc3cnt(3) = imax + 1
1417  call check( nf90_put_var(ncid, ncv, age_c_conv, start=nc3cor, count=nc3cnt) )
1418 
1419  call check( nf90_inq_varid(ncid, 'age_t', ncv) )
1420  nc3cor(1) = 1
1421  nc3cor(2) = 1
1422  nc3cor(3) = 1
1423  nc3cnt(1) = ktmax + 1
1424  nc3cnt(2) = jmax + 1
1425  nc3cnt(3) = imax + 1
1426  call check( nf90_put_var(ncid, ncv, age_t_conv, start=nc3cor, count=nc3cnt) )
1427 
1428 end if
1429 
1430 call check( nf90_inq_varid(ncid, 'dzs_dtau', ncv) )
1431 nc2cor(1) = 1
1432 nc2cor(2) = 1
1433 nc2cnt(1) = jmax + 1
1434 nc2cnt(2) = imax + 1
1435 call check( nf90_put_var(ncid, ncv, dzs_dtau_conv, start=nc2cor, count=nc2cnt) )
1436 
1437 call check( nf90_inq_varid(ncid, 'dzm_dtau', ncv) )
1438 nc2cor(1) = 1
1439 nc2cor(2) = 1
1440 nc2cnt(1) = jmax + 1
1441 nc2cnt(2) = imax + 1
1442 call check( nf90_put_var(ncid, ncv, dzm_dtau_conv, start=nc2cor, count=nc2cnt) )
1443 
1444 call check( nf90_inq_varid(ncid, 'dzb_dtau', ncv) )
1445 nc2cor(1) = 1
1446 nc2cor(2) = 1
1447 nc2cnt(1) = jmax + 1
1448 nc2cnt(2) = imax + 1
1449 call check( nf90_put_var(ncid, ncv, dzb_dtau_conv, start=nc2cor, count=nc2cnt) )
1450 
1451 call check( nf90_inq_varid(ncid, 'dzl_dtau', ncv) )
1452 nc2cor(1) = 1
1453 nc2cor(2) = 1
1454 nc2cnt(1) = jmax + 1
1455 nc2cnt(2) = imax + 1
1456 call check( nf90_put_var(ncid, ncv, dzl_dtau_conv, start=nc2cor, count=nc2cnt) )
1457 
1458 call check( nf90_inq_varid(ncid, 'dH_c_dtau', ncv) )
1459 nc2cor(1) = 1
1460 nc2cor(2) = 1
1461 nc2cnt(1) = jmax + 1
1462 nc2cnt(2) = imax + 1
1463 call check( nf90_put_var(ncid, ncv, dh_c_dtau_conv, start=nc2cor, count=nc2cnt) )
1464 
1465 call check( nf90_inq_varid(ncid, 'dH_t_dtau', ncv) )
1466 nc2cor(1) = 1
1467 nc2cor(2) = 1
1468 nc2cnt(1) = jmax + 1
1469 nc2cnt(2) = imax + 1
1470 call check( nf90_put_var(ncid, ncv, dh_t_dtau_conv, start=nc2cor, count=nc2cnt) )
1471 
1472 call check( nf90_inq_varid(ncid, 'dH_dtau', ncv) )
1473 nc2cor(1) = 1
1474 nc2cor(2) = 1
1475 nc2cnt(1) = jmax + 1
1476 nc2cnt(2) = imax + 1
1477 call check( nf90_put_var(ncid, ncv, dh_dtau_conv, start=nc2cor, count=nc2cnt) )
1478 
1479 call check( nf90_inq_varid(ncid, 'vx_b_g', ncv) )
1480 nc2cor(1) = 1
1481 nc2cor(2) = 1
1482 nc2cnt(1) = jmax + 1
1483 nc2cnt(2) = imax + 1
1484 call check( nf90_put_var(ncid, ncv, vx_b_g_conv, start=nc2cor, count=nc2cnt) )
1485 
1486 call check( nf90_inq_varid(ncid, 'vy_b_g', ncv) )
1487 nc2cor(1) = 1
1488 nc2cor(2) = 1
1489 nc2cnt(1) = jmax + 1
1490 nc2cnt(2) = imax + 1
1491 call check( nf90_put_var(ncid, ncv, vy_b_g_conv, start=nc2cor, count=nc2cnt) )
1492 
1493 call check( nf90_inq_varid(ncid, 'vz_b', ncv) )
1494 nc2cor(1) = 1
1495 nc2cor(2) = 1
1496 nc2cnt(1) = jmax + 1
1497 nc2cnt(2) = imax + 1
1498 call check( nf90_put_var(ncid, ncv, vz_b_conv, start=nc2cor, count=nc2cnt) )
1499 
1500 call check( nf90_inq_varid(ncid, 'vx_s_g', ncv) )
1501 nc2cor(1) = 1
1502 nc2cor(2) = 1
1503 nc2cnt(1) = jmax + 1
1504 nc2cnt(2) = imax + 1
1505 call check( nf90_put_var(ncid, ncv, vx_s_g_conv, start=nc2cor, count=nc2cnt) )
1506 
1507 call check( nf90_inq_varid(ncid, 'vy_s_g', ncv) )
1508 nc2cor(1) = 1
1509 nc2cor(2) = 1
1510 nc2cnt(1) = jmax + 1
1511 nc2cnt(2) = imax + 1
1512 call check( nf90_put_var(ncid, ncv, vy_s_g_conv, start=nc2cor, count=nc2cnt) )
1513 
1514 call check( nf90_inq_varid(ncid, 'vz_s', ncv) )
1515 nc2cor(1) = 1
1516 nc2cor(2) = 1
1517 nc2cnt(1) = jmax + 1
1518 nc2cnt(2) = imax + 1
1519 call check( nf90_put_var(ncid, ncv, vz_s_conv, start=nc2cor, count=nc2cnt) )
1520 
1521 call check( nf90_inq_varid(ncid, 'temp_b', ncv) )
1522 nc2cor(1) = 1
1523 nc2cor(2) = 1
1524 nc2cnt(1) = jmax + 1
1525 nc2cnt(2) = imax + 1
1526 call check( nf90_put_var(ncid, ncv, temp_b_conv, start=nc2cor, count=nc2cnt) )
1527 
1528 call check( nf90_inq_varid(ncid, 'temph_b', ncv) )
1529 nc2cor(1) = 1
1530 nc2cor(2) = 1
1531 nc2cnt(1) = jmax + 1
1532 nc2cnt(2) = imax + 1
1533 call check( nf90_put_var(ncid, ncv, temph_b_conv, start=nc2cor, count=nc2cnt) )
1534 
1535 call check( nf90_inq_varid(ncid, 'p_b_w', ncv) )
1536 nc2cor(1) = 1
1537 nc2cor(2) = 1
1538 nc2cnt(1) = jmax + 1
1539 nc2cnt(2) = imax + 1
1540 call check( nf90_put_var(ncid, ncv, p_b_w_conv, start=nc2cor, count=nc2cnt) )
1541 
1542 call check( nf90_inq_varid(ncid, 'H_w', ncv) )
1543 nc2cor(1) = 1
1544 nc2cor(2) = 1
1545 nc2cnt(1) = jmax + 1
1546 nc2cnt(2) = imax + 1
1547 call check( nf90_put_var(ncid, ncv, h_w_conv, start=nc2cor, count=nc2cnt) )
1548 
1549 call check( nf90_inq_varid(ncid, 'q_gl_g', ncv) )
1550 nc2cor(1) = 1
1551 nc2cor(2) = 1
1552 nc2cnt(1) = jmax + 1
1553 nc2cnt(2) = imax + 1
1554 call check( nf90_put_var(ncid, ncv, q_gl_g_conv, start=nc2cor, count=nc2cnt) )
1555 
1556 call check( nf90_inq_varid(ncid, 'V_tot', ncv) )
1557 call check( nf90_put_var(ncid, ncv, v_tot_conv) )
1558 
1559 call check( nf90_inq_varid(ncid, 'A_grounded', ncv) )
1560 call check( nf90_put_var(ncid, ncv, a_grounded_conv) )
1561 
1562 call check( nf90_inq_varid(ncid, 'A_floating', ncv) )
1563 call check( nf90_put_var(ncid, ncv, a_floating_conv) )
1564 
1565 !-------- Close NetCDF file --------
1566 
1567 call check( nf90_sync(ncid) )
1568 
1569 call check( nf90_close(ncid) )
1570 
1571 !-------- Increase file counter --------
1572 
1573 ndat = ndat+1
1574 
1575 if (flag_3d_output) then
1576  ndat3d = ndat
1577 else
1578  ndat2d = ndat
1579 end if
1580 
1581 end subroutine output_nc
1582 !