1 | ! $Id: regr_lat_time_climoz_m.F90 3065 2017-11-10 13:25:09Z fhourdin $ |
---|
2 | module regr_lat_time_climoz_m |
---|
3 | |
---|
4 | ! Author: Lionel GUEZ |
---|
5 | |
---|
6 | implicit none |
---|
7 | |
---|
8 | private |
---|
9 | public regr_lat_time_climoz |
---|
10 | |
---|
11 | contains |
---|
12 | |
---|
13 | subroutine regr_lat_time_climoz(read_climoz, regr_lat) |
---|
14 | |
---|
15 | ! "regr_lat_time_climoz" stands for "regrid latitude time |
---|
16 | ! climatology ozone". |
---|
17 | |
---|
18 | ! This procedure reads a climatology of ozone from a NetCDF file, |
---|
19 | ! regrids it in latitude and time, and writes the regridded field |
---|
20 | ! to a new NetCDF file. |
---|
21 | |
---|
22 | ! The input field depends on time, pressure level and latitude. |
---|
23 | |
---|
24 | ! If the input field has missing values, they must be signaled by |
---|
25 | ! the "missing_value" attribute. |
---|
26 | |
---|
27 | ! We assume that the input field is a step function of latitude |
---|
28 | ! and that the input latitude coordinate gives the centers of steps. |
---|
29 | ! Regridding in latitude is made by averaging, with a cosine of |
---|
30 | ! latitude factor. |
---|
31 | ! The target LMDZ latitude grid is the "scalar" grid: "rlatu". |
---|
32 | ! The values of "rlatu" are taken to be the centers of intervals. |
---|
33 | |
---|
34 | ! We assume that in the input file: |
---|
35 | |
---|
36 | ! -- Latitude is in degrees. |
---|
37 | |
---|
38 | ! -- Latitude and pressure are strictly monotonic (as all NetCDF |
---|
39 | ! coordinate variables should be). |
---|
40 | |
---|
41 | ! -- The time coordinate is in ascending order (even though we do |
---|
42 | ! not use its values). |
---|
43 | ! The input file may contain either values for 12 months or values |
---|
44 | ! for 14 months. |
---|
45 | ! If there are 14 months then we assume that we have (in that order): |
---|
46 | ! December, January, February, ..., November, December, January |
---|
47 | |
---|
48 | ! -- Missing values are contiguous, at the bottom of |
---|
49 | ! the vertical domain and at the latitudinal boundaries. |
---|
50 | |
---|
51 | ! If values are all missing at a given latitude and date, then we |
---|
52 | ! replace those missing values by values at the closest latitude, |
---|
53 | ! equatorward, with valid values. |
---|
54 | ! Then, at each latitude and each date, the missing values are replaced |
---|
55 | ! by the lowest valid value above missing values. |
---|
56 | |
---|
57 | ! Regridding in time is by linear interpolation. |
---|
58 | ! Monthly values are processed to get daily values, on the basis |
---|
59 | ! of a 360-day calendar. |
---|
60 | ! If there are 14 months, we use the first December value to |
---|
61 | ! interpolate values between January 1st and mid-January. |
---|
62 | ! We use the last January value to interpolate values between |
---|
63 | ! mid-December and end of December. |
---|
64 | ! If there are only 12 months in the input file then we assume |
---|
65 | ! periodicity for interpolation at the beginning and at the end of the |
---|
66 | ! year. |
---|
67 | use mod_grid_phy_lmdz, ONLY : nbp_lat, grid_type, unstructured |
---|
68 | use regr1_step_av_m, only: regr1_step_av |
---|
69 | use regr3_lint_m, only: regr3_lint |
---|
70 | use netcdf95, only: handle_err, nf95_close, nf95_get_att, nf95_gw_var, & |
---|
71 | nf95_inq_dimid, nf95_inq_varid, nf95_inquire_dimension, nf95_open, & |
---|
72 | nf95_put_var |
---|
73 | use netcdf, only: nf90_get_att, nf90_get_var, nf90_noerr, nf90_nowrite |
---|
74 | use assert_m, only: assert |
---|
75 | use regular_lonlat_mod, only : boundslat_reg, south |
---|
76 | use nrtype, only: pi |
---|
77 | use regular_lonlat_mod, only : lat_reg |
---|
78 | implicit none |
---|
79 | integer, intent(in):: read_climoz ! read ozone climatology |
---|
80 | ! Allowed values are 1 and 2 |
---|
81 | ! 1: read a single ozone climatology that will be used day and night |
---|
82 | ! 2: read two ozone climatologies, the average day and night |
---|
83 | ! climatology and the daylight climatology |
---|
84 | |
---|
85 | ! Variables local to the procedure: |
---|
86 | LOGICAL, intent(in), OPTIONAL :: regr_lat |
---|
87 | |
---|
88 | integer n_plev ! number of pressure levels in the input data |
---|
89 | integer n_lat ! number of latitudes in the input data |
---|
90 | integer n_month ! number of months in the input data |
---|
91 | integer n_lat_out ! number of latitudes in the output data |
---|
92 | |
---|
93 | real, pointer:: latitude(:) |
---|
94 | ! (of input data, converted to rad, sorted in strictly ascending order) |
---|
95 | |
---|
96 | real, allocatable:: lat_in_edg(:) |
---|
97 | ! (edges of latitude intervals for input data, in rad, in strictly |
---|
98 | ! ascending order) |
---|
99 | |
---|
100 | real, pointer:: plev(:) |
---|
101 | ! pressure levels of input data, sorted in strictly ascending |
---|
102 | ! order, converted to hPa |
---|
103 | |
---|
104 | logical desc_lat ! latitude in descending order in the input file |
---|
105 | logical desc_plev ! pressure levels in descending order in the input file |
---|
106 | |
---|
107 | real, allocatable:: o3_in(:, :, :, :) |
---|
108 | ! (n_lat, n_plev, n_month, read_climoz) |
---|
109 | ! ozone climatologies from the input file |
---|
110 | ! "o3_in(j, k, :, :)" is at latitude "latitude(j)" and pressure |
---|
111 | ! level "plev(k)". |
---|
112 | ! Third dimension is month index, first value may be December or January. |
---|
113 | ! "o3_in(:, :, :, 1)" is for the day- night average, "o3_in(:, :, :, 2)" |
---|
114 | ! is for daylight. |
---|
115 | |
---|
116 | real missing_value |
---|
117 | |
---|
118 | real, allocatable:: o3_regr_lat(:, :, :, :) |
---|
119 | ! (nbp_lat, n_plev, 0:13, read_climoz) |
---|
120 | ! mean of "o3_in" over a latitude interval of LMDZ |
---|
121 | ! First dimension is latitude interval. |
---|
122 | ! The latitude interval for "o3_regr_lat(j,:, :, :)" contains "rlatu(j)". |
---|
123 | ! If "j" is between 2 and "nbp_lat-1" then the interval is: |
---|
124 | ! [rlatv(j), rlatv(j-1)] |
---|
125 | ! If "j" is 1 or "nbp_lat" then the interval is: |
---|
126 | ! [rlatv(1), pi / 2] |
---|
127 | ! or: |
---|
128 | ! [- pi / 2, rlatv(nbp_lat-1)] |
---|
129 | ! respectively. |
---|
130 | ! "o3_regr_lat(:, k, :, :)" is for pressure level "plev(k)". |
---|
131 | ! Third dimension is month number, 1 for January. |
---|
132 | ! "o3_regr_lat(:, :, :, 1)" is average day and night, |
---|
133 | ! "o3_regr_lat(:, :, :, 2)" is for daylight. |
---|
134 | |
---|
135 | real, allocatable:: o3_out(:, :, :, :) |
---|
136 | ! (nbp_lat, n_plev, 360, read_climoz) |
---|
137 | ! regridded ozone climatology |
---|
138 | ! "o3_out(j, k, l, :)" is at latitude "rlatu(j)", pressure |
---|
139 | ! level "plev(k)" and date "January 1st 0h" + "tmidday(l)", in a |
---|
140 | ! 360-day calendar. |
---|
141 | ! "o3_out(:, :, :, 1)" is average day and night, |
---|
142 | ! "o3_out(:, :, :, 2)" is for daylight. |
---|
143 | |
---|
144 | integer j, k, l,m |
---|
145 | |
---|
146 | ! For NetCDF: |
---|
147 | integer ncid_in, ncid_out ! IDs for input and output files |
---|
148 | integer varid_plev, varid_time, varid, ncerr, dimid |
---|
149 | character(len=80) press_unit ! pressure unit |
---|
150 | |
---|
151 | integer varid_in(read_climoz), varid_out(read_climoz) |
---|
152 | ! index 1 is for average ozone day and night, index 2 is for |
---|
153 | ! daylight ozone. |
---|
154 | |
---|
155 | real, parameter:: tmidmonth(0:13) = (/(-15. + 30. * l, l = 0, 13)/) |
---|
156 | ! (time to middle of month, in days since January 1st 0h, in a |
---|
157 | ! 360-day calendar) |
---|
158 | ! (We add values -15 and 375 so that, for example, day 3 of the year is |
---|
159 | ! interpolated between the December and the January value.) |
---|
160 | |
---|
161 | real, parameter:: tmidday(360) = (/(l + 0.5, l = 0, 359)/) |
---|
162 | ! (time to middle of day, in days since January 1st 0h, in a |
---|
163 | ! 360-day calendar) |
---|
164 | logical :: regr_lat_ |
---|
165 | REAL,ALLOCATABLE :: lat_out(:) |
---|
166 | |
---|
167 | !--------------------------------- |
---|
168 | IF (PRESENT(regr_lat)) THEN |
---|
169 | regr_lat_=regr_lat |
---|
170 | ELSE |
---|
171 | regr_lat_=.TRUE. |
---|
172 | ENDIF |
---|
173 | |
---|
174 | print *, "Call sequence information: regr_lat_time_climoz" |
---|
175 | call assert(read_climoz == 1 .or. read_climoz == 2, "regr_lat_time_climoz") |
---|
176 | |
---|
177 | call nf95_open("climoz.nc", nf90_nowrite, ncid_in) |
---|
178 | |
---|
179 | ! Get coordinates from the input file: |
---|
180 | |
---|
181 | call nf95_inq_varid(ncid_in, "latitude", varid) |
---|
182 | call nf95_gw_var(ncid_in, varid, latitude) |
---|
183 | ! Convert from degrees to rad, because we will take the sine of latitude: |
---|
184 | latitude = latitude / 180. * pi |
---|
185 | n_lat = size(latitude) |
---|
186 | ! We need to supply the latitudes to "regr1_step_av" in |
---|
187 | ! ascending order, so invert order if necessary: |
---|
188 | desc_lat = latitude(1) > latitude(n_lat) |
---|
189 | if (desc_lat) latitude = latitude(n_lat:1:-1) |
---|
190 | |
---|
191 | ! Compute edges of latitude intervals: |
---|
192 | allocate(lat_in_edg(n_lat + 1)) |
---|
193 | lat_in_edg(1) = - pi / 2 |
---|
194 | forall (j = 2:n_lat) lat_in_edg(j) = (latitude(j - 1) + latitude(j)) / 2 |
---|
195 | lat_in_edg(n_lat + 1) = pi / 2 |
---|
196 | |
---|
197 | call nf95_inq_varid(ncid_in, "plev", varid) |
---|
198 | call nf95_gw_var(ncid_in, varid, plev) |
---|
199 | n_plev = size(plev) |
---|
200 | ! We only need the pressure coordinate to copy it to the output file. |
---|
201 | ! The program "gcm" will assume that pressure levels are in |
---|
202 | ! ascending order in the regridded climatology so invert order if |
---|
203 | ! necessary: |
---|
204 | desc_plev = plev(1) > plev(n_plev) |
---|
205 | if (desc_plev) plev = plev(n_plev:1:-1) |
---|
206 | call nf95_get_att(ncid_in, varid, "units", press_unit) |
---|
207 | if (press_unit == "Pa") then |
---|
208 | ! Convert to hPa: |
---|
209 | plev = plev / 100. |
---|
210 | elseif (press_unit /= "hPa") then |
---|
211 | print *, "regr_lat_time_climoz: the only recognized units are Pa " & |
---|
212 | // "and hPa." |
---|
213 | stop 1 |
---|
214 | end if |
---|
215 | |
---|
216 | IF (regr_lat_) THEN |
---|
217 | n_lat_out=nbp_lat |
---|
218 | ALLOCATE(lat_out(n_lat_out)) |
---|
219 | lat_out=lat_reg |
---|
220 | ELSE |
---|
221 | n_lat_out=n_lat |
---|
222 | ALLOCATE(lat_out(n_lat_out)) |
---|
223 | lat_out=latitude |
---|
224 | ENDIF |
---|
225 | |
---|
226 | ! Create the output file and get the variable IDs: |
---|
227 | call prepare_out(ncid_in, n_lat_out, lat_out, n_plev, ncid_out, varid_out, varid_plev, & |
---|
228 | varid_time) |
---|
229 | |
---|
230 | deallocate(latitude) ! pointer |
---|
231 | |
---|
232 | ! Write remaining coordinate variables: |
---|
233 | call nf95_put_var(ncid_out, varid_plev, plev) |
---|
234 | call nf95_put_var(ncid_out, varid_time, tmidday) |
---|
235 | |
---|
236 | deallocate(plev) ! pointer |
---|
237 | |
---|
238 | ! Get the number of months: |
---|
239 | call nf95_inq_dimid(ncid_in, "time", dimid) |
---|
240 | call nf95_inquire_dimension(ncid_in, dimid, nclen=n_month) |
---|
241 | |
---|
242 | allocate(o3_in(n_lat, n_plev, n_month, read_climoz)) |
---|
243 | |
---|
244 | call nf95_inq_varid(ncid_in, "tro3", varid_in(1)) |
---|
245 | ncerr = nf90_get_var(ncid_in, varid_in(1), o3_in(:, :, :, 1)) |
---|
246 | call handle_err("regr_lat_time_climoz nf90_get_var tro3", ncerr, ncid_in) |
---|
247 | |
---|
248 | if (read_climoz == 2) then |
---|
249 | call nf95_inq_varid(ncid_in, "tro3_daylight", varid_in(2)) |
---|
250 | ncerr = nf90_get_var(ncid_in, varid_in(2), o3_in(:, :, :, 2)) |
---|
251 | call handle_err("regr_lat_time_climoz nf90_get_var tro3_daylight", & |
---|
252 | ncerr, ncid_in, varid_in(2)) |
---|
253 | end if |
---|
254 | |
---|
255 | if (desc_lat) o3_in = o3_in(n_lat:1:-1, :, :, :) |
---|
256 | if (desc_plev) o3_in = o3_in(:, n_plev:1:-1, :, :) |
---|
257 | |
---|
258 | do m = 1, read_climoz |
---|
259 | ncerr = nf90_get_att(ncid_in, varid_in(m), "missing_value", & |
---|
260 | missing_value) |
---|
261 | if (ncerr == nf90_noerr) then |
---|
262 | do l = 1, n_month |
---|
263 | ! Take care of latitudes where values are all missing: |
---|
264 | |
---|
265 | ! Next to the south pole: |
---|
266 | j = 1 |
---|
267 | do while (o3_in(j, 1, l, m) == missing_value) |
---|
268 | j = j + 1 |
---|
269 | end do |
---|
270 | if (j > 1) o3_in(:j-1, :, l, m) = & |
---|
271 | spread(o3_in(j, :, l, m), dim=1, ncopies=j-1) |
---|
272 | |
---|
273 | ! Next to the north pole: |
---|
274 | j = n_lat |
---|
275 | do while (o3_in(j, 1, l, m) == missing_value) |
---|
276 | j = j - 1 |
---|
277 | end do |
---|
278 | if (j < n_lat) o3_in(j+1:, :, l, m) = & |
---|
279 | spread(o3_in(j, :, l, m), dim=1, ncopies=n_lat-j) |
---|
280 | |
---|
281 | ! Take care of missing values at high pressure: |
---|
282 | do j = 1, n_lat |
---|
283 | ! Find missing values, starting from top of atmosphere |
---|
284 | ! and going down. |
---|
285 | ! We have already taken care of latitudes full of |
---|
286 | ! missing values so the highest level has a valid value. |
---|
287 | k = 2 |
---|
288 | do while (o3_in(j, k, l, m) /= missing_value .and. k < n_plev) |
---|
289 | k = k + 1 |
---|
290 | end do |
---|
291 | ! Replace missing values with the valid value at the |
---|
292 | ! lowest level above missing values: |
---|
293 | if (o3_in(j, k, l, m) == missing_value) & |
---|
294 | o3_in(j, k:n_plev, l, m) = o3_in(j, k-1, l, m) |
---|
295 | end do |
---|
296 | end do |
---|
297 | else |
---|
298 | print *, "regr_lat_time_climoz: field ", m, & |
---|
299 | ", no missing value attribute" |
---|
300 | end if |
---|
301 | end do |
---|
302 | |
---|
303 | call nf95_close(ncid_in) |
---|
304 | |
---|
305 | allocate(o3_out(n_lat_out, n_plev, 360, read_climoz)) |
---|
306 | IF (regr_lat_) THEN |
---|
307 | allocate(o3_regr_lat(nbp_lat, n_plev, 0:13, read_climoz)) |
---|
308 | |
---|
309 | ! Regrid in latitude: |
---|
310 | ! We average with respect to sine of latitude, which is |
---|
311 | ! equivalent to weighting by cosine of latitude: |
---|
312 | if (n_month == 12) then |
---|
313 | print *, & |
---|
314 | "Found 12 months in ozone climatologies, assuming periodicity..." |
---|
315 | o3_regr_lat(nbp_lat:1:-1, :, 1:12, :) = regr1_step_av(o3_in, & |
---|
316 | xs=sin(lat_in_edg), xt=sin((/- pi / 2, boundslat_reg(nbp_lat-1:1:-1,south), pi / 2/))) |
---|
317 | ! (invert order of indices in "o3_regr_lat" because "rlatu" is |
---|
318 | ! in descending order) |
---|
319 | |
---|
320 | ! Duplicate January and December values, in preparation of time |
---|
321 | ! interpolation: |
---|
322 | o3_regr_lat(:, :, 0, :) = o3_regr_lat(:, :, 12, :) |
---|
323 | o3_regr_lat(:, :, 13, :) = o3_regr_lat(:, :, 1, :) |
---|
324 | else |
---|
325 | print *, "Using 14 months in ozone climatologies..." |
---|
326 | o3_regr_lat(nbp_lat:1:-1, :, :, :) = regr1_step_av(o3_in, & |
---|
327 | xs=sin(lat_in_edg), xt=sin((/- pi / 2, boundslat_reg(nbp_lat-1:1:-1,south), pi / 2/))) |
---|
328 | ! (invert order of indices in "o3_regr_lat" because "rlatu" is |
---|
329 | ! in descending order) |
---|
330 | end if |
---|
331 | ! Regrid in time by linear interpolation: |
---|
332 | o3_out = regr3_lint(o3_in, tmidmonth, tmidday) |
---|
333 | ELSE |
---|
334 | ! Regrid in time by linear interpolation: |
---|
335 | o3_out = regr3_lint(o3_in, tmidmonth, tmidday) |
---|
336 | ENDIF |
---|
337 | |
---|
338 | ! Write to file: |
---|
339 | |
---|
340 | do m = 1, read_climoz |
---|
341 | IF (grid_type==unstructured) THEN |
---|
342 | ! Doing spatial interpolation from XIOS need to have some point in longitude |
---|
343 | ! waiting zonal mean operation from XIOS |
---|
344 | call nf95_put_var(ncid_out, varid_out(m), SPREAD(o3_out(n_lat_out:1:-1, :, :, m),1,4)) |
---|
345 | ELSE |
---|
346 | call nf95_put_var(ncid_out, varid_out(m), o3_out(n_lat_out:1:-1, :, :, m)) |
---|
347 | ENDIF |
---|
348 | |
---|
349 | ! (The order of "rlatu" is inverted in the output file) |
---|
350 | end do |
---|
351 | |
---|
352 | call nf95_close(ncid_out) |
---|
353 | |
---|
354 | end subroutine regr_lat_time_climoz |
---|
355 | |
---|
356 | !******************************************** |
---|
357 | |
---|
358 | subroutine prepare_out(ncid_in, n_lat, lat, n_plev, ncid_out, varid_out, varid_plev, & |
---|
359 | varid_time) |
---|
360 | |
---|
361 | ! This subroutine creates the NetCDF output file, defines |
---|
362 | ! dimensions and variables, and writes one of the coordinate variables. |
---|
363 | |
---|
364 | use netcdf95, only: nf95_create, nf95_def_dim, nf95_def_var, & |
---|
365 | nf95_put_att, nf95_enddef, nf95_copy_att, nf95_put_var |
---|
366 | use netcdf, only: nf90_clobber,nf90_netcdf4, nf90_float, nf90_global, nf90_unlimited |
---|
367 | use nrtype, only: pi |
---|
368 | use mod_grid_phy_lmdz, ONLY : grid_type, unstructured |
---|
369 | |
---|
370 | integer, intent(in):: ncid_in, n_lat, n_plev |
---|
371 | real, intent(in):: lat(:) |
---|
372 | integer, intent(out):: ncid_out, varid_plev, varid_time |
---|
373 | |
---|
374 | integer, intent(out):: varid_out(:) ! dim(1 or 2) |
---|
375 | ! "varid_out(1)" is for average ozone day and night, |
---|
376 | ! "varid_out(2)" is for daylight ozone. |
---|
377 | |
---|
378 | ! Variables local to the procedure: |
---|
379 | |
---|
380 | integer ncerr |
---|
381 | integer dimid_rlatu, dimid_rlonv, dimid_plev, dimid_time |
---|
382 | integer varid_rlonv, varid_rlatu |
---|
383 | |
---|
384 | !--------------------------- |
---|
385 | |
---|
386 | print *, "Call sequence information: prepare_out" |
---|
387 | |
---|
388 | call nf95_create("climoz_LMDZ.nc", NF90_CLOBBER + NF90_NETCDF4, ncid_out) |
---|
389 | |
---|
390 | ! Dimensions: |
---|
391 | call nf95_def_dim(ncid_out, "time", nf90_unlimited, dimid_time) |
---|
392 | call nf95_def_dim(ncid_out, "plev", n_plev, dimid_plev) |
---|
393 | call nf95_def_dim(ncid_out, "rlatu", n_lat, dimid_rlatu) |
---|
394 | if (grid_type==unstructured) call nf95_def_dim(ncid_out, "rlonv", 4, dimid_rlonv) |
---|
395 | |
---|
396 | ! Define coordinate variables: |
---|
397 | |
---|
398 | call nf95_def_var(ncid_out, "time", nf90_float, dimid_time, varid_time) |
---|
399 | call nf95_put_att(ncid_out, varid_time, "units", "days since 2000-1-1") |
---|
400 | call nf95_put_att(ncid_out, varid_time, "calendar", "360_day") |
---|
401 | call nf95_put_att(ncid_out, varid_time, "standard_name", "time") |
---|
402 | |
---|
403 | call nf95_def_var(ncid_out, "plev", nf90_float, dimid_plev, varid_plev) |
---|
404 | call nf95_put_att(ncid_out, varid_plev, "units", "millibar") |
---|
405 | call nf95_put_att(ncid_out, varid_plev, "standard_name", "air_pressure") |
---|
406 | call nf95_put_att(ncid_out, varid_plev, "long_name", "air pressure") |
---|
407 | |
---|
408 | call nf95_def_var(ncid_out, "rlatu", nf90_float, dimid_rlatu, varid_rlatu) |
---|
409 | call nf95_put_att(ncid_out, varid_rlatu, "units", "degrees_north") |
---|
410 | call nf95_put_att(ncid_out, varid_rlatu, "standard_name", "latitude") |
---|
411 | |
---|
412 | if (grid_type==unstructured) then |
---|
413 | call nf95_def_var(ncid_out, "rlonv", nf90_float, dimid_rlonv, varid_rlonv) |
---|
414 | call nf95_put_att(ncid_out, varid_rlonv, "units", "degrees_east") |
---|
415 | call nf95_put_att(ncid_out, varid_rlonv, "standard_name", "longitude") |
---|
416 | endif |
---|
417 | |
---|
418 | ! Define the primary variables: |
---|
419 | |
---|
420 | if (grid_type==unstructured) then |
---|
421 | call nf95_def_var(ncid_out, "tro3", nf90_float, & |
---|
422 | (/dimid_rlonv, dimid_rlatu, dimid_plev, dimid_time/), varid_out(1)) |
---|
423 | else |
---|
424 | call nf95_def_var(ncid_out, "tro3", nf90_float, & |
---|
425 | (/dimid_rlatu, dimid_plev, dimid_time/), varid_out(1)) |
---|
426 | endif |
---|
427 | |
---|
428 | call nf95_put_att(ncid_out, varid_out(1), "long_name", & |
---|
429 | "ozone mole fraction") |
---|
430 | call nf95_put_att(ncid_out, varid_out(1), "standard_name", & |
---|
431 | "mole_fraction_of_ozone_in_air") |
---|
432 | |
---|
433 | if (size(varid_out) == 2) then |
---|
434 | if (grid_type==unstructured) then |
---|
435 | call nf95_def_var(ncid_out, "tro3_daylight", nf90_float, & |
---|
436 | (/dimid_rlonv, dimid_rlatu, dimid_plev, dimid_time/), varid_out(2)) |
---|
437 | else |
---|
438 | call nf95_def_var(ncid_out, "tro3_daylight", nf90_float, & |
---|
439 | (/ dimid_rlatu, dimid_plev, dimid_time/), varid_out(2)) |
---|
440 | endif |
---|
441 | |
---|
442 | call nf95_put_att(ncid_out, varid_out(2), "long_name", & |
---|
443 | "ozone mole fraction in daylight") |
---|
444 | end if |
---|
445 | |
---|
446 | ! Global attributes: |
---|
447 | |
---|
448 | ! The following commands, copying attributes, may fail. |
---|
449 | ! That is OK. |
---|
450 | ! It should just mean that the attribute is not defined in the input file. |
---|
451 | |
---|
452 | call nf95_copy_att(ncid_in, nf90_global, "Conventions", ncid_out, & |
---|
453 | nf90_global, ncerr) |
---|
454 | call handle_err_copy_att("Conventions") |
---|
455 | |
---|
456 | call nf95_copy_att(ncid_in, nf90_global, "title", ncid_out, nf90_global, & |
---|
457 | ncerr) |
---|
458 | call handle_err_copy_att("title") |
---|
459 | |
---|
460 | call nf95_copy_att(ncid_in, nf90_global, "institution", ncid_out, & |
---|
461 | nf90_global, ncerr) |
---|
462 | call handle_err_copy_att("institution") |
---|
463 | |
---|
464 | call nf95_copy_att(ncid_in, nf90_global, "source", ncid_out, nf90_global, & |
---|
465 | ncerr) |
---|
466 | call handle_err_copy_att("source") |
---|
467 | |
---|
468 | call nf95_put_att(ncid_out, nf90_global, "comment", "Regridded for LMDZ") |
---|
469 | |
---|
470 | call nf95_enddef(ncid_out) |
---|
471 | |
---|
472 | ! Write one of the coordinate variables: |
---|
473 | call nf95_put_var(ncid_out, varid_rlatu, lat(n_lat:1:-1) / pi * 180.) |
---|
474 | ! (convert from rad to degrees and sort in ascending order) |
---|
475 | if (grid_type==unstructured) call nf95_put_var(ncid_out, varid_rlonv, (/ 0., 90., 180.,270. /)) |
---|
476 | |
---|
477 | contains |
---|
478 | |
---|
479 | subroutine handle_err_copy_att(att_name) |
---|
480 | |
---|
481 | use netcdf, only: nf90_noerr, nf90_strerror |
---|
482 | |
---|
483 | character(len=*), intent(in):: att_name |
---|
484 | |
---|
485 | !---------------------------------------- |
---|
486 | |
---|
487 | if (ncerr /= nf90_noerr) then |
---|
488 | print *, "regr_lat_time_climoz_m prepare_out nf95_copy_att " & |
---|
489 | // att_name // " -- " // trim(nf90_strerror(ncerr)) |
---|
490 | end if |
---|
491 | |
---|
492 | end subroutine handle_err_copy_att |
---|
493 | |
---|
494 | end subroutine prepare_out |
---|
495 | |
---|
496 | end module regr_lat_time_climoz_m |
---|