1 | PROGRAM bdy_reorder |
---|
2 | |
---|
3 | !=============================================================================== |
---|
4 | ! A routine to reorder old BDY data files to make them compatible with NEMO 3.4. |
---|
5 | ! |
---|
6 | ! This routine has 2 modes: |
---|
7 | ! 1. If no template file is given then it will re-order the data in the input |
---|
8 | ! file so that it is in order of increasing nbr. |
---|
9 | ! 2. If a template file is given it will re-order the data in the input file |
---|
10 | ! be consistent with the order of the data in the template file. This is |
---|
11 | ! useful for making old barotropic and baroclinic data files consistent. |
---|
12 | ! For older versions of NEMO the barotropic boundary data did not have to |
---|
13 | ! be in the same order as the baroclinic boundary data. The rimwidth |
---|
14 | ! value can be different in the input file and template file; in this case |
---|
15 | ! the routine will just re-order the data that the two files have in common. |
---|
16 | ! |
---|
17 | ! The routine is mainly for re-ordering BDY data files, but can also be used to |
---|
18 | ! re-order BDY coordinate files if ln_coordinates is set to .true. |
---|
19 | ! |
---|
20 | ! Author: Dave Storkey Aug 2011 |
---|
21 | ! Bug notifications etc to: dave.storkey@metoffice.gov.uk |
---|
22 | ! |
---|
23 | !=============================================================================== |
---|
24 | |
---|
25 | USE netcdf |
---|
26 | |
---|
27 | IMPLICIT NONE |
---|
28 | |
---|
29 | INTEGER,PARAMETER :: numnam=11 |
---|
30 | INTEGER,PARAMETER :: sp=SELECTED_REAL_KIND(6,37) |
---|
31 | INTEGER,PARAMETER :: dp=SELECTED_REAL_KIND(12,307) |
---|
32 | INTEGER :: chunksize = 32000000 |
---|
33 | |
---|
34 | INTEGER :: jpbgrd, iostat, ncid_in, ncid_out, ncid_template, nbrid(4), nbrid_template(4) |
---|
35 | INTEGER :: dimids(10), unlimitedDimId |
---|
36 | INTEGER :: xtype, dimid, varid, ndims, ndims_var, nvars, natts, nblen(4), nblen_template(4), dimlen(10) |
---|
37 | INTEGER :: lenvar(10), igrid, ib, ib1, ir, idim, jgrid, jv, icount, attid, len1, idepth, itime |
---|
38 | INTEGER :: idim_time, idim_xb, idim_yb, idim_depth |
---|
39 | INTEGER :: nbr_min, nbr_max, nbr_min_template, nbr_max_template, strlen, nbr_match |
---|
40 | |
---|
41 | INTEGER, ALLOCATABLE, DIMENSION(:,:) :: nbr, imap |
---|
42 | INTEGER, ALLOCATABLE, DIMENSION(:,:) :: nbi, nbj |
---|
43 | INTEGER, ALLOCATABLE, DIMENSION(:,:) :: nbi_template |
---|
44 | INTEGER, ALLOCATABLE, DIMENSION(:,:) :: nbj_template |
---|
45 | INTEGER, ALLOCATABLE, DIMENSION(:,:) :: nbr_template |
---|
46 | |
---|
47 | INTEGER, ALLOCATABLE, DIMENSION(:) :: nbi_extract, nbj_extract |
---|
48 | |
---|
49 | INTEGER, ALLOCATABLE, DIMENSION(:,:,:,:) :: varin_int, varout_int |
---|
50 | REAL(sp), ALLOCATABLE, DIMENSION(:,:,:,:) :: varin_float, varout_float |
---|
51 | REAL(dp), ALLOCATABLE, DIMENSION(:,:,:,:) :: varin_dble, varout_dble |
---|
52 | |
---|
53 | LOGICAL :: ln_coordinates |
---|
54 | LOGICAL :: ln_template |
---|
55 | |
---|
56 | CHARACTER(LEN=nf90_max_name) :: attname, dimname, varname, time, date, zone, timestamp |
---|
57 | CHARACTER(LEN=nf90_max_name) :: file_in, file_out, file_template, nbi_name, nbj_name, nbr_name |
---|
58 | CHARACTER(LEN=1),DIMENSION(4) :: cgrid |
---|
59 | CHARACTER(LEN=1) :: end_letter |
---|
60 | |
---|
61 | NAMELIST/nam_bdy_reorder/ file_in, file_out, file_template, ln_coordinates, nbr_match |
---|
62 | |
---|
63 | !============================================================================= |
---|
64 | |
---|
65 | !----------------------------------------------------------------------------- |
---|
66 | ! 0. Read namelist and initialise parameters |
---|
67 | !----------------------------------------------------------------------------- |
---|
68 | |
---|
69 | file_in = '' |
---|
70 | file_out = '' |
---|
71 | file_template = '' |
---|
72 | ln_coordinates = .false. |
---|
73 | nbr_match = 0 |
---|
74 | |
---|
75 | OPEN( UNIT=numnam, FILE='nam_bdy_reorder', FORM='FORMATTED', STATUS='OLD' ) |
---|
76 | READ( numnam, nam_bdy_reorder ) |
---|
77 | CLOSE( UNIT=numnam ) |
---|
78 | |
---|
79 | IF( ln_coordinates ) THEN |
---|
80 | WRITE(*,*) 'Input file is a coordinates.bdy.nc file.' |
---|
81 | jpbgrd = 4 |
---|
82 | ELSE |
---|
83 | WRITE(*,*) 'Input file is a boundary data file.' |
---|
84 | jpbgrd = 1 |
---|
85 | ENDIF |
---|
86 | |
---|
87 | IF( LEN_TRIM(file_template) > 0 ) THEN |
---|
88 | WRITE(*,*) 'Reordering data to match template file.' |
---|
89 | ln_template = .true. |
---|
90 | ELSE |
---|
91 | WRITE(*,*) 'Reordering data in order of increasing nbr values.' |
---|
92 | ln_template = .false. |
---|
93 | ENDIF |
---|
94 | |
---|
95 | cgrid= (/'t','u','v','f'/) |
---|
96 | |
---|
97 | !----------------------------------------------------------------------------- |
---|
98 | ! 1. Read in nbr variables from file and create mapping. |
---|
99 | !----------------------------------------------------------------------------- |
---|
100 | |
---|
101 | ! 1.1 Open input files |
---|
102 | |
---|
103 | iostat = nf90_open( TRIM(file_in), nf90_nowrite, ncid_in ) |
---|
104 | IF (iostat /= nf90_noerr) THEN |
---|
105 | WRITE(6,*) TRIM(nf90_strerror(iostat)) |
---|
106 | STOP |
---|
107 | ENDIF |
---|
108 | iostat = nf90_inquire( ncid_in, ndims, nvars, natts ) |
---|
109 | |
---|
110 | IF( ln_template ) THEN |
---|
111 | iostat = nf90_open( TRIM(file_template), nf90_nowrite, ncid_template ) |
---|
112 | IF (iostat /= nf90_noerr) THEN |
---|
113 | WRITE(6,*) TRIM(nf90_strerror(iostat)) |
---|
114 | STOP |
---|
115 | ENDIF |
---|
116 | ENDIF |
---|
117 | |
---|
118 | DO igrid=1,jpbgrd |
---|
119 | |
---|
120 | ! 1.2 Find dimensions of data in input files and allocate arrays |
---|
121 | |
---|
122 | IF( ln_coordinates ) THEN |
---|
123 | nbr_name = 'nbr'//cgrid(igrid) |
---|
124 | ELSE |
---|
125 | nbr_name = 'nbrdta' |
---|
126 | ENDIF |
---|
127 | |
---|
128 | iostat = nf90_inq_varid( ncid_in, nbr_name , nbrid(igrid) ) |
---|
129 | IF (iostat /= nf90_noerr) THEN |
---|
130 | WRITE(6,*) TRIM(nf90_strerror(iostat)) |
---|
131 | STOP |
---|
132 | ENDIF |
---|
133 | iostat = nf90_inquire_variable( ncid_in, nbrid(igrid), ndims=ndims, dimids=dimids ) |
---|
134 | IF( ndims .ne. 2 ) THEN |
---|
135 | WRITE(*,*) 'ERROR : ',TRIM(nbr_name),' does not have exactly 2 dimensions.' |
---|
136 | WRITE(*,*) ' in file ',TRIM(file_in) |
---|
137 | STOP |
---|
138 | ENDIF |
---|
139 | iostat = nf90_inquire_dimension( ncid_in, dimids(1), len=nblen(igrid)) |
---|
140 | iostat = nf90_inquire_dimension( ncid_in, dimids(2), len=len1) |
---|
141 | WRITE(*,*) 'nblen(igrid), len1 : ',nblen(igrid), len1 |
---|
142 | IF( len1 .ne. 1 ) THEN |
---|
143 | WRITE(*,*) 'ERROR : second dimension of ',TRIM(nbr_name),' does not have length 1.' |
---|
144 | WRITE(*,*) ' in file ',TRIM(file_in) |
---|
145 | WRITE(*,*) ' len1 = ',len1 |
---|
146 | STOP |
---|
147 | ENDIF |
---|
148 | |
---|
149 | IF( ln_template ) THEN |
---|
150 | |
---|
151 | iostat = nf90_inq_varid( ncid_template, nbr_name , nbrid_template(igrid) ) |
---|
152 | IF (iostat /= nf90_noerr) THEN |
---|
153 | WRITE(6,*) TRIM(nf90_strerror(iostat)) |
---|
154 | STOP |
---|
155 | ENDIF |
---|
156 | iostat = nf90_inquire_variable( ncid_template, nbrid_template(igrid), ndims=ndims, dimids=dimids ) |
---|
157 | IF( ndims .ne. 2 ) THEN |
---|
158 | WRITE(*,*) 'ERROR : ',TRIM(nbr_name),' does not have exactly 2 dimensions.' |
---|
159 | WRITE(*,*) ' in file ',TRIM(file_template) |
---|
160 | STOP |
---|
161 | ENDIF |
---|
162 | iostat = nf90_inquire_dimension( ncid_template, dimids(1), len=nblen_template(igrid)) |
---|
163 | iostat = nf90_inquire_dimension( ncid_template, dimids(2), len=len1) |
---|
164 | WRITE(*,*) 'nblen_template(igrid), len1 : ',nblen_template(igrid), len1 |
---|
165 | IF( len1 .ne. 1 ) THEN |
---|
166 | WRITE(*,*) 'ERROR : second dimension of ',TRIM(nbr_name),' does not have length 1.' |
---|
167 | WRITE(*,*) ' in file ',TRIM(file_template) |
---|
168 | WRITE(*,*) ' len1 = ',len1 |
---|
169 | STOP |
---|
170 | ENDIF |
---|
171 | |
---|
172 | ENDIF ! ln_template |
---|
173 | |
---|
174 | ENDDO ! jpbgrd |
---|
175 | |
---|
176 | ALLOCATE( nbr(MAXVAL(nblen),1), imap(MAXVAL(nblen),jpbgrd) ) |
---|
177 | IF( ln_template ) THEN |
---|
178 | ALLOCATE( nbi(MAXVAL(nblen),1), nbj(MAXVAL(nblen),1) ) |
---|
179 | ALLOCATE( nbi_template(MAXVAL(nblen_template),1) ) |
---|
180 | ALLOCATE( nbj_template(MAXVAL(nblen_template),1) ) |
---|
181 | ALLOCATE( nbr_template(MAXVAL(nblen_template),1) ) |
---|
182 | ENDIF |
---|
183 | |
---|
184 | ! 1.3 Read in nbr variables and generate mapping. |
---|
185 | |
---|
186 | nbr(:,:) = -1 |
---|
187 | imap(:,:) = -1 |
---|
188 | IF( ln_template ) THEN |
---|
189 | nbi(:,:) = -1 |
---|
190 | nbj(:,:) = -1 |
---|
191 | nbi_template(:,:) = -1 |
---|
192 | nbj_template(:,:) = -1 |
---|
193 | nbr_template(:,:) = -1 |
---|
194 | ENDIF |
---|
195 | |
---|
196 | WRITE(*,*) '>>> Generating map.' |
---|
197 | DO igrid=1,jpbgrd |
---|
198 | |
---|
199 | iostat = nf90_get_var( ncid_in, nbrid(igrid), nbr(1:nblen(igrid),1:1) ) |
---|
200 | IF (iostat /= nf90_noerr) THEN |
---|
201 | WRITE(6,*) TRIM(nf90_strerror(iostat)) |
---|
202 | STOP |
---|
203 | ENDIF |
---|
204 | |
---|
205 | nbr_min = MINVAL( nbr(1:nblen(igrid),1) ) |
---|
206 | nbr_max = MAXVAL( nbr(1:nblen(igrid),1) ) |
---|
207 | IF( nbr_min .ne. 1 ) THEN |
---|
208 | WRITE(*,*) 'Something wrong with input file ',file_in |
---|
209 | WRITE(*,*) 'MIN(nbr) /= 1' |
---|
210 | STOP |
---|
211 | ENDIF |
---|
212 | |
---|
213 | IF( ln_template ) THEN |
---|
214 | |
---|
215 | IF( ln_coordinates ) THEN |
---|
216 | nbi_name = 'nbi'//cgrid(igrid) |
---|
217 | nbj_name = 'nbj'//cgrid(igrid) |
---|
218 | ELSE |
---|
219 | nbi_name = 'nbidta' |
---|
220 | nbj_name = 'nbjdta' |
---|
221 | ENDIF |
---|
222 | |
---|
223 | ! read in nbi, nbj values from input file |
---|
224 | iostat = nf90_inq_varid( ncid_in, nbi_name, varid ) |
---|
225 | IF (iostat == nf90_noerr) iostat = nf90_get_var( ncid_in, varid, nbi(1:nblen(igrid),1:1) ) |
---|
226 | IF (iostat /= nf90_noerr) THEN |
---|
227 | WRITE(6,*) TRIM(nf90_strerror(iostat)) |
---|
228 | STOP |
---|
229 | ENDIF |
---|
230 | iostat = nf90_inq_varid( ncid_in, nbj_name, varid ) |
---|
231 | IF (iostat == nf90_noerr) iostat = nf90_get_var( ncid_in, varid, nbj(1:nblen(igrid),1:1) ) |
---|
232 | IF (iostat /= nf90_noerr) THEN |
---|
233 | WRITE(6,*) TRIM(nf90_strerror(iostat)) |
---|
234 | STOP |
---|
235 | ENDIF |
---|
236 | |
---|
237 | ! read in nbi, nbj, nbr values from template file |
---|
238 | iostat = nf90_inq_varid( ncid_template, nbi_name, varid ) |
---|
239 | IF (iostat == nf90_noerr) iostat = nf90_get_var( ncid_template, varid, nbi_template(1:nblen_template(igrid),1:1) ) |
---|
240 | IF (iostat /= nf90_noerr) THEN |
---|
241 | WRITE(6,*) TRIM(nf90_strerror(iostat)) |
---|
242 | STOP |
---|
243 | ENDIF |
---|
244 | iostat = nf90_inq_varid( ncid_template, nbj_name, varid ) |
---|
245 | IF (iostat == nf90_noerr) iostat = nf90_get_var( ncid_template, varid, nbj_template(1:nblen_template(igrid),1:1) ) |
---|
246 | IF (iostat /= nf90_noerr) THEN |
---|
247 | WRITE(6,*) TRIM(nf90_strerror(iostat)) |
---|
248 | STOP |
---|
249 | ENDIF |
---|
250 | iostat = nf90_get_var( ncid_template, nbrid_template(igrid), nbr_template(1:nblen_template(igrid),1:1) ) |
---|
251 | IF (iostat /= nf90_noerr) THEN |
---|
252 | WRITE(6,*) TRIM(nf90_strerror(iostat)) |
---|
253 | STOP |
---|
254 | ENDIF |
---|
255 | |
---|
256 | nbr_min_template = MINVAL( nbr_template(1:nblen_template(igrid),1) ) |
---|
257 | nbr_max_template = MAXVAL( nbr_template(1:nblen_template(igrid),1) ) |
---|
258 | IF( nbr_min_template .ne. 1 ) THEN |
---|
259 | WRITE(*,*) 'Something wrong with template file ',file_template |
---|
260 | WRITE(*,*) 'MIN(nbr) /= 1' |
---|
261 | STOP |
---|
262 | ENDIF |
---|
263 | |
---|
264 | IF( nbr_match .lt. 1 ) nbr_match = MIN(nbr_max, nbr_max_template) |
---|
265 | |
---|
266 | ! initialise imap to the identity mapping |
---|
267 | DO ib = 1, nblen(igrid) |
---|
268 | imap(ib, igrid) = ib |
---|
269 | ENDDO |
---|
270 | |
---|
271 | ! allocate "extract" arrays |
---|
272 | icount = 0 |
---|
273 | DO ib = 1, nblen_template(igrid) |
---|
274 | IF( nbr_template(ib,1) .eq. 1 ) THEN |
---|
275 | icount = icount+1 |
---|
276 | ENDIF |
---|
277 | ENDDO |
---|
278 | ALLOCATE( nbi_extract(icount), nbj_extract(icount) ) |
---|
279 | |
---|
280 | DO ir = 1, nbr_match |
---|
281 | |
---|
282 | ! extract values from template array for this ir value |
---|
283 | icount = 0 |
---|
284 | DO ib = 1, nblen_template(igrid) |
---|
285 | IF( nbr_template(ib,1) .eq. ir ) THEN |
---|
286 | icount = icount+1 |
---|
287 | nbi_extract(icount) = nbi_template(ib,1) |
---|
288 | nbj_extract(icount) = nbj_template(ib,1) |
---|
289 | ENDIF |
---|
290 | ENDDO |
---|
291 | |
---|
292 | ! work out the mapping for this ir value |
---|
293 | icount = 1 |
---|
294 | DO ib = 1, nblen(igrid) |
---|
295 | IF( nbr(ib,1) .eq. ir ) THEN |
---|
296 | DO ib1 = 1, nblen(igrid) |
---|
297 | IF( nbi(ib1,1) .eq. nbi_extract(icount) .and. & |
---|
298 | & nbj(ib1,1) .eq. nbj_extract(icount) ) THEN |
---|
299 | imap(ib,igrid) = ib1 |
---|
300 | icount = icount + 1 |
---|
301 | EXIT |
---|
302 | ENDIF |
---|
303 | ENDDO ! ib1 |
---|
304 | ENDIF |
---|
305 | ENDDO ! ib |
---|
306 | |
---|
307 | ENDDO ! ir |
---|
308 | |
---|
309 | DEALLOCATE( nbi_extract, nbj_extract ) |
---|
310 | |
---|
311 | ELSE |
---|
312 | |
---|
313 | icount = 0 |
---|
314 | DO ir = nbr_min, nbr_max |
---|
315 | DO ib = 1, nblen(igrid) |
---|
316 | IF( nbr(ib,1) .eq. ir ) THEN |
---|
317 | icount = icount + 1 |
---|
318 | imap(icount,igrid) = ib |
---|
319 | ENDIF |
---|
320 | ENDDO |
---|
321 | ENDDO |
---|
322 | |
---|
323 | ENDIF |
---|
324 | |
---|
325 | ENDDO ! jpbgrd |
---|
326 | |
---|
327 | !----------------------------------------------------------------------------- |
---|
328 | ! 2. Open output file and copy dimensions and attributes across |
---|
329 | !----------------------------------------------------------------------------- |
---|
330 | |
---|
331 | iostat = nf90_inquire( ncid_in, ndims, nvars, natts ) |
---|
332 | |
---|
333 | ! 2.1 Create the output file |
---|
334 | |
---|
335 | WRITE(*,*) '>>> Initialising output file.' |
---|
336 | iostat = nf90_create( TRIM(file_out), nf90_64bit_offset, ncid_out, chunksize=chunksize) |
---|
337 | |
---|
338 | ! 2.2 Copy the dimensions into the output file. |
---|
339 | |
---|
340 | iostat = nf90_inquire( ncid_in, unlimitedDimId = unlimitedDimId ) |
---|
341 | dimlen(:) = 1 |
---|
342 | DO idim = 1, ndims |
---|
343 | iostat = nf90_inquire_dimension(ncid_in, idim, dimname, dimlen(idim)) |
---|
344 | IF (idim == unlimitedDimId) THEN |
---|
345 | iostat = nf90_def_dim( ncid_out, dimname, nf90_unlimited, dimid) |
---|
346 | idim_time = idim |
---|
347 | ELSE |
---|
348 | iostat = nf90_def_dim( ncid_out, dimname, dimlen(idim), dimid) |
---|
349 | IF( INDEX(dimname,'x') .gt. 0 ) THEN |
---|
350 | idim_xb = idim |
---|
351 | ELSE IF( INDEX(dimname,'y') .gt. 0 ) THEN |
---|
352 | idim_yb = idim |
---|
353 | ELSE IF( INDEX(dimname,'depth') .gt. 0 .or. INDEX(dimname,'z') .gt. 0 ) THEN |
---|
354 | idim_depth = idim |
---|
355 | ELSE |
---|
356 | WRITE(*,*) 'ERROR: Unrecognised dimension : ',dimname |
---|
357 | STOP |
---|
358 | ENDIF |
---|
359 | ENDIF |
---|
360 | END DO |
---|
361 | |
---|
362 | |
---|
363 | ! 2.2 Copy the global attributes into the output file. |
---|
364 | ! Also need to change the file_name attribute and the TimeStamp attribute. |
---|
365 | |
---|
366 | DO attid = 1, natts |
---|
367 | iostat = nf90_inq_attname( ncid_in, nf90_global, attid, attname ) |
---|
368 | WRITE(6,*)'>>> Copying attribute '//TRIM(attname)//' into destination file...' |
---|
369 | iostat = nf90_copy_att( ncid_in, nf90_global, attname, ncid_out, nf90_global ) |
---|
370 | END DO |
---|
371 | iostat = nf90_put_att( ncid_out, nf90_global, "file_name", TRIM(file_out)) |
---|
372 | CALL DATE_AND_TIME ( date=date, time=time, zone=zone ) |
---|
373 | timestamp = date(7:8) // "/" // date(5:6) // "/" // date(1:4) // " " // & |
---|
374 | time(1:2) // ":" // time(3:4) // ":" // time(5:6) // " " // & |
---|
375 | zone |
---|
376 | iostat = nf90_put_att( ncid_out, nf90_global, "TimeStamp", timestamp) |
---|
377 | |
---|
378 | ! 2.3 Copy the variable definitions and attributes into the output file. |
---|
379 | |
---|
380 | DO jv = 1, nvars |
---|
381 | iostat = nf90_inquire_variable( ncid_in, jv, varname, xtype, ndims_var, dimids, natts ) |
---|
382 | iostat = nf90_def_var( ncid_out, varname, xtype, dimids(1:ndims_var), varid ) |
---|
383 | IF (natts > 0) THEN |
---|
384 | DO attid = 1, natts |
---|
385 | iostat = nf90_inq_attname(ncid_in, varid, attid, attname) |
---|
386 | iostat = nf90_copy_att( ncid_in, varid, attname, ncid_out, varid ) |
---|
387 | END DO |
---|
388 | ENDIF |
---|
389 | END DO |
---|
390 | |
---|
391 | ! 2.4 End definitions in output file |
---|
392 | |
---|
393 | iostat = nf90_enddef( ncid_out ) |
---|
394 | |
---|
395 | !----------------------------------------------------------------------------- |
---|
396 | ! 3. Read in variables from input file, re-order and write to output file |
---|
397 | !----------------------------------------------------------------------------- |
---|
398 | |
---|
399 | IF( ln_coordinates ) THEN |
---|
400 | ALLOCATE( varin_int(MAXVAL(nblen),1,1,1), varout_int(MAXVAL(nblen),1,1,1) ) |
---|
401 | ALLOCATE( varin_float(MAXVAL(nblen),1,1,1), varout_float(MAXVAL(nblen),1,1,1) ) |
---|
402 | ALLOCATE( varin_dble(MAXVAL(nblen),1,1,1), varout_dble(MAXVAL(nblen),1,1,1) ) |
---|
403 | ELSE |
---|
404 | SELECT CASE( ndims ) |
---|
405 | CASE( 2 ) |
---|
406 | ALLOCATE( varin_int(dimlen(idim_xb),dimlen(idim_yb),1,1) ) |
---|
407 | ALLOCATE( varout_int(dimlen(idim_xb),dimlen(idim_yb),1,1) ) |
---|
408 | ALLOCATE( varin_float(dimlen(idim_xb),dimlen(idim_yb),1,1) ) |
---|
409 | ALLOCATE( varout_float(dimlen(idim_xb),dimlen(idim_yb),1,1) ) |
---|
410 | ALLOCATE( varin_dble(dimlen(idim_xb),dimlen(idim_yb),1,1) ) |
---|
411 | ALLOCATE( varout_dble(dimlen(idim_xb),dimlen(idim_yb),1,1) ) |
---|
412 | CASE ( 3 ) |
---|
413 | ALLOCATE( varin_int(dimlen(idim_xb),dimlen(idim_yb),1,dimlen(idim_time)) ) |
---|
414 | ALLOCATE( varout_int(dimlen(idim_xb),dimlen(idim_yb),1,dimlen(idim_time)) ) |
---|
415 | ALLOCATE( varin_float(dimlen(idim_xb),dimlen(idim_yb),1,dimlen(idim_time)) ) |
---|
416 | ALLOCATE( varout_float(dimlen(idim_xb),dimlen(idim_yb),1,dimlen(idim_time)) ) |
---|
417 | ALLOCATE( varin_dble(dimlen(idim_xb),dimlen(idim_yb),1,dimlen(idim_time)) ) |
---|
418 | ALLOCATE( varout_dble(dimlen(idim_xb),dimlen(idim_yb),1,dimlen(idim_time)) ) |
---|
419 | CASE ( 4 ) |
---|
420 | ALLOCATE( varin_int(dimlen(idim_xb),dimlen(idim_yb),dimlen(idim_depth),dimlen(idim_time)) ) |
---|
421 | ALLOCATE( varout_int(dimlen(idim_xb),dimlen(idim_yb),dimlen(idim_depth),dimlen(idim_time)) ) |
---|
422 | ALLOCATE( varin_float(dimlen(idim_xb),dimlen(idim_yb),dimlen(idim_depth),dimlen(idim_time)) ) |
---|
423 | ALLOCATE( varout_float(dimlen(idim_xb),dimlen(idim_yb),dimlen(idim_depth),dimlen(idim_time)) ) |
---|
424 | ALLOCATE( varin_dble(dimlen(idim_xb),dimlen(idim_yb),dimlen(idim_depth),dimlen(idim_time)) ) |
---|
425 | ALLOCATE( varout_dble(dimlen(idim_xb),dimlen(idim_yb),dimlen(idim_depth),dimlen(idim_time)) ) |
---|
426 | CASE DEFAULT |
---|
427 | WRITE(*,*) 'ERROR : Can only cope with 2, 3 or 4 dimensions in the boundary data files.' |
---|
428 | WRITE(*,*) ' This file appears to have ',ndims,' dimensions.' |
---|
429 | STOP |
---|
430 | END SELECT |
---|
431 | ENDIF |
---|
432 | |
---|
433 | DO jv = 1, nvars |
---|
434 | iostat = nf90_inquire_variable( ncid_in, jv, varname, xtype, ndims_var, dimids ) |
---|
435 | DO idim = 1, ndims_var |
---|
436 | lenvar(idim) = dimlen(dimids(idim)) |
---|
437 | ENDDO |
---|
438 | |
---|
439 | IF( ndims_var .eq. 1 ) THEN |
---|
440 | WRITE(*,*) '>>> Copying coordinate variable ',TRIM(varname) |
---|
441 | ELSE |
---|
442 | WRITE(*,*) '>>> Reordering variable ',TRIM(varname) |
---|
443 | ENDIF |
---|
444 | ! Error check here |
---|
445 | |
---|
446 | IF( ln_coordinates ) THEN |
---|
447 | strlen = len_trim(varname) |
---|
448 | end_letter = varname(strlen:strlen) |
---|
449 | jgrid = -1 |
---|
450 | DO igrid = 1,4 |
---|
451 | IF( end_letter .eq. cgrid(igrid) ) jgrid = igrid |
---|
452 | END DO |
---|
453 | IF( jgrid .lt. 0 ) THEN |
---|
454 | WRITE(*,*) 'ERROR : Could not identify grid for variable ',TRIM(varname) |
---|
455 | WRITE(*,*) ' varname : ',TRIM(varname),'!' |
---|
456 | WRITE(*,*) ' strlen : ',strlen |
---|
457 | WRITE(*,*) ' end_letter : ',end_letter |
---|
458 | WRITE(*,*) ' cgrid(1) : ',cgrid(1) |
---|
459 | WRITE(*,*) ' cgrid(2) : ',cgrid(2) |
---|
460 | WRITE(*,*) ' cgrid(3) : ',cgrid(3) |
---|
461 | WRITE(*,*) ' cgrid(4) : ',cgrid(4) |
---|
462 | STOP |
---|
463 | ENDIF |
---|
464 | ELSE |
---|
465 | jgrid=1 |
---|
466 | ENDIF |
---|
467 | |
---|
468 | SELECT CASE(xtype) |
---|
469 | |
---|
470 | CASE( NF90_INT ) |
---|
471 | SELECT CASE(ndims_var) |
---|
472 | CASE( 1 ) |
---|
473 | WRITE(*,*) 'This is a 1D integer.' |
---|
474 | ! Assume this is a depth or time coordinate variable and copy across unchanged. |
---|
475 | iostat = nf90_get_var( ncid_in, jv, varin_int(1:lenvar(1),1,1,1) ) |
---|
476 | iostat = nf90_put_var( ncid_out, jv, varout_int(1:lenvar(1),1,1,1) ) |
---|
477 | CASE( 2 ) |
---|
478 | WRITE(*,*) 'This is a 2D integer.' |
---|
479 | iostat = nf90_get_var( ncid_in, jv, varin_int(1:lenvar(1),1:1,1,1) ) |
---|
480 | DO ib = 1,lenvar(1) |
---|
481 | varout_int(ib,1,1,1) = varin_int(imap(ib,jgrid),1,1,1) |
---|
482 | END DO |
---|
483 | iostat = nf90_put_var( ncid_out, jv, varout_int(1:lenvar(1),1:1,1,1) ) |
---|
484 | CASE( 3 ) |
---|
485 | WRITE(*,*) 'This is a 3D integer.' |
---|
486 | ! Assume third dimension is time. |
---|
487 | iostat = nf90_get_var( ncid_in, jv, varin_int(1:lenvar(1),1:1,1,1:lenvar(3)) ) |
---|
488 | DO itime = 1,lenvar(3) |
---|
489 | DO ib = 1,lenvar(1) |
---|
490 | varout_int(ib,1,1,itime) = varin_int(imap(ib,jgrid),1,1,itime) |
---|
491 | END DO |
---|
492 | END DO |
---|
493 | iostat = nf90_put_var( ncid_out, jv, varout_int(1:lenvar(1),1:1,1,1:lenvar(3)) ) |
---|
494 | CASE( 4 ) |
---|
495 | WRITE(*,*) 'This is a 4D integer.' |
---|
496 | ! Assume third and fourth dimensions are depth and time respectively. |
---|
497 | iostat = nf90_get_var( ncid_in, jv, varin_int(1:lenvar(1),1:1,1:lenvar(3),1:lenvar(4)) ) |
---|
498 | DO itime = 1,lenvar(4) |
---|
499 | DO idepth = 1,lenvar(3) |
---|
500 | DO ib = 1,lenvar(1) |
---|
501 | varout_int(ib,1,idepth,itime) = varin_int(imap(ib,jgrid),1,idepth,itime) |
---|
502 | END DO |
---|
503 | END DO |
---|
504 | END DO |
---|
505 | iostat = nf90_put_var( ncid_out, jv, varout_int(1:lenvar(1),1:1,1:lenvar(3),1:lenvar(4)) ) |
---|
506 | CASE DEFAULT |
---|
507 | WRITE(*,*) 'ERROR : Variable with ',ndims_var,' dimensions. Case not coded.' |
---|
508 | STOP |
---|
509 | END SELECT |
---|
510 | |
---|
511 | CASE( NF90_FLOAT ) |
---|
512 | SELECT CASE(ndims_var) |
---|
513 | CASE( 1 ) |
---|
514 | WRITE(*,*) 'This is a 1D float.' |
---|
515 | ! Assume this is a depth or time coordinate variable and copy across unchanged. |
---|
516 | iostat = nf90_get_var( ncid_in, jv, varin_float(1:lenvar(1),1,1,1) ) |
---|
517 | iostat = nf90_put_var( ncid_out, jv, varout_float(1:lenvar(1),1,1,1) ) |
---|
518 | CASE( 2 ) |
---|
519 | WRITE(*,*) 'This is a 2D float.' |
---|
520 | iostat = nf90_get_var( ncid_in, jv, varin_float(1:lenvar(1),1:1,1,1) ) |
---|
521 | DO ib = 1,lenvar(1) |
---|
522 | varout_float(ib,1,1,1) = varin_float(imap(ib,jgrid),1,1,1) |
---|
523 | END DO |
---|
524 | iostat = nf90_put_var( ncid_out, jv, varout_float(1:lenvar(1),1:1,1,1) ) |
---|
525 | CASE( 3 ) |
---|
526 | WRITE(*,*) 'This is a 3D float.' |
---|
527 | ! Assume third dimension is time. |
---|
528 | iostat = nf90_get_var( ncid_in, jv, varin_float(1:lenvar(1),1:1,1,1:lenvar(3)) ) |
---|
529 | DO itime = 1,lenvar(3) |
---|
530 | DO ib = 1,lenvar(1) |
---|
531 | varout_float(ib,1,1,itime) = varin_float(imap(ib,jgrid),1,1,itime) |
---|
532 | END DO |
---|
533 | END DO |
---|
534 | iostat = nf90_put_var( ncid_out, jv, varout_float(1:lenvar(1),1:1,1,1:lenvar(3)) ) |
---|
535 | CASE( 4 ) |
---|
536 | WRITE(*,*) 'This is a 4D float.' |
---|
537 | ! Assume third and fourth dimensions are depth and time respectively. |
---|
538 | iostat = nf90_get_var( ncid_in, jv, varin_float(1:lenvar(1),1:1,1:lenvar(3),1:lenvar(4)) ) |
---|
539 | DO itime = 1,lenvar(4) |
---|
540 | DO idepth = 1,lenvar(3) |
---|
541 | DO ib = 1,lenvar(1) |
---|
542 | varout_float(ib,1,idepth,itime) = varin_float(imap(ib,jgrid),1,idepth,itime) |
---|
543 | END DO |
---|
544 | END DO |
---|
545 | END DO |
---|
546 | iostat = nf90_put_var( ncid_out, jv, varout_float(1:lenvar(1),1:1,1:lenvar(3),1:lenvar(4)) ) |
---|
547 | CASE DEFAULT |
---|
548 | WRITE(*,*) 'ERROR : Variable with ',ndims_var,' dimensions. Case not coded.' |
---|
549 | STOP |
---|
550 | END SELECT |
---|
551 | |
---|
552 | CASE( NF90_DOUBLE ) |
---|
553 | SELECT CASE(ndims_var) |
---|
554 | CASE( 1 ) |
---|
555 | WRITE(*,*) 'This is a 1D double.' |
---|
556 | ! Assume this is a depth or time coordinate variable and copy across unchanged. |
---|
557 | iostat = nf90_get_var( ncid_in, jv, varin_dble(1:lenvar(1),1,1,1) ) |
---|
558 | iostat = nf90_put_var( ncid_out, jv, varout_dble(1:lenvar(1),1,1,1) ) |
---|
559 | CASE( 2 ) |
---|
560 | WRITE(*,*) 'This is a 2D double.' |
---|
561 | iostat = nf90_get_var( ncid_in, jv, varin_dble(1:lenvar(1),1:1,1,1) ) |
---|
562 | DO ib = 1,lenvar(1) |
---|
563 | varout_dble(ib,1,1,1) = varin_dble(imap(ib,jgrid),1,1,1) |
---|
564 | END DO |
---|
565 | iostat = nf90_put_var( ncid_out, jv, varout_dble(1:lenvar(1),1:1,1,1) ) |
---|
566 | CASE( 3 ) |
---|
567 | WRITE(*,*) 'This is a 3D double.' |
---|
568 | ! Assume third dimension is time. |
---|
569 | iostat = nf90_get_var( ncid_in, jv, varin_dble(1:lenvar(1),1:1,1,1:lenvar(3)) ) |
---|
570 | DO itime = 1,lenvar(3) |
---|
571 | DO ib = 1,lenvar(1) |
---|
572 | varout_dble(ib,1,1,itime) = varin_dble(imap(ib,jgrid),1,1,itime) |
---|
573 | END DO |
---|
574 | END DO |
---|
575 | iostat = nf90_put_var( ncid_out, jv, varout_dble(1:lenvar(1),1:1,1,1:lenvar(3)) ) |
---|
576 | CASE( 4 ) |
---|
577 | WRITE(*,*) 'This is a 4D double.' |
---|
578 | ! Assume third and fourth dimensions are depth and time respectively. |
---|
579 | iostat = nf90_get_var( ncid_in, jv, varin_dble(1:lenvar(1),1:1,1:lenvar(3),1:lenvar(4)) ) |
---|
580 | DO itime = 1,lenvar(4) |
---|
581 | DO idepth = 1,lenvar(3) |
---|
582 | DO ib = 1,lenvar(1) |
---|
583 | varout_dble(ib,1,idepth,itime) = varin_dble(imap(ib,jgrid),1,idepth,itime) |
---|
584 | END DO |
---|
585 | END DO |
---|
586 | END DO |
---|
587 | iostat = nf90_put_var( ncid_out, jv, varout_dble(1:lenvar(1),1:1,1:lenvar(3),1:lenvar(4)) ) |
---|
588 | CASE DEFAULT |
---|
589 | WRITE(*,*) 'ERROR : Variable with ',ndims_var,' dimensions. Case not coded.' |
---|
590 | STOP |
---|
591 | END SELECT |
---|
592 | |
---|
593 | CASE DEFAULT |
---|
594 | WRITE(*,*) 'ERROR : Unrecognised data type.' |
---|
595 | STOP |
---|
596 | |
---|
597 | END SELECT |
---|
598 | |
---|
599 | END DO ! jv |
---|
600 | |
---|
601 | !----------------------------------------------------------------------------- |
---|
602 | ! 4. Close input and output files |
---|
603 | !----------------------------------------------------------------------------- |
---|
604 | |
---|
605 | iostat = nf90_close( ncid_out ) |
---|
606 | iostat = nf90_close( ncid_in ) |
---|
607 | |
---|
608 | |
---|
609 | END PROGRAM bdy_reorder |
---|