| blob | 1499cd44700de76b58fd791365bd97813b8d0bd9 |
1 !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
2 ! Module: source_data_module
3 !
4 ! Description:
5 !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
6 module source_data_module
8 use hash_module
9 use list_module
10 use module_debug
11 use misc_definitions_module
13 ! Parameters
14 integer, parameter :: RETURN_LANDMASK = 0, &
15 RETURN_DOMCAT_LM = 1, &
16 RETURN_DFDX_LM = 2, &
17 RETURN_DFDY_LM = 3, &
18 RETURN_FIELDNAME = 4, &
19 RETURN_DOMCAT = 5, &
20 RETURN_DFDX = 6, &
21 RETURN_DFDY = 7
22 integer, parameter :: MAX_LANDMASK_CATEGORIES = 100
24 ! Module variables
25 integer :: num_entries
26 integer :: next_field = 1
27 integer :: output_field_state = RETURN_LANDMASK
28 integer, pointer, dimension(:) :: source_proj, source_wordsize, source_endian, source_fieldtype, &
29 source_dest_fieldtype, source_priority, source_tile_x, source_tile_y, &
30 source_tile_z, source_tile_z_start, source_tile_z_end, source_tile_bdr, &
31 source_category_min, source_category_max, source_smooth_option, &
32 source_smooth_passes, source_output_stagger, source_row_order
33 integer :: source_iswater, source_islake, source_isice, source_isurban, source_isoilwater
34 real, pointer, dimension(:) :: source_dx, source_dy, source_known_x, source_known_y, &
35 source_known_lat, source_known_lon, source_masked, source_truelat1, source_truelat2, &
36 source_stdlon, source_scale_factor, source_missing_value, source_fill_missing
37 character (len=128), pointer, dimension(:) :: source_fieldname, source_path, source_interp_string, &
38 source_dominant_category, source_dominant_only, source_dfdx, source_dfdy, &
39 source_z_dim_name, source_units, source_descr
40 character (len=128) :: source_mminlu
41 logical, pointer, dimension(:) :: is_proj, is_wordsize, is_endian, is_fieldtype, &
42 is_dest_fieldtype, is_priority, is_tile_x, is_tile_y, is_tile_z, &
43 is_tile_z_start, is_tile_z_end, is_tile_bdr, is_category_min, &
44 is_category_max, is_masked, &
45 is_dx, is_dy, is_known_x, is_known_y, &
46 is_known_lat, is_known_lon, is_truelat1, is_truelat2, is_stdlon, &
47 is_scale_factor, is_fieldname, is_path, is_dominant_category, &
48 is_dominant_only, is_dfdx, is_dfdy, is_z_dim_name, &
49 is_smooth_option, is_smooth_passes, is_signed, is_missing_value, &
50 is_fill_missing, is_halt_missing, is_output_stagger, is_row_order, &
51 is_units, is_descr, is_subgrid
53 type (list), pointer, dimension(:) :: source_res_path, source_interp_option, source_landmask_land, &
54 source_landmask_water
55 type (hashtable) :: bad_files ! Track which files produce errors when we try to open them
56 type (hashtable) :: duplicate_fnames ! Remember which output fields we have returned
59 contains
62 !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
63 ! Name: get_datalist
64 !
65 ! Purpose:
66 !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
67 subroutine get_datalist()
69 use gridinfo_module
70 use stringutil
72 implicit none
74 ! Parameters
75 integer, parameter :: BUFSIZE = 256
77 ! Local variables
78 integer :: nparams, idx, eos, ispace, comma, i, j, funit
79 logical :: have_specification, is_used
80 character (len=128) :: res_string, path_string, interp_string, landmask_string
81 character (len=BUFSIZE) :: buffer
83 nparams = 0
84 num_entries = 0
86 do funit=10,100
87 inquire(unit=funit, opened=is_used)
88 if (.not. is_used) exit
89 end do
90 open(funit,file=trim(opt_geogrid_tbl_path)//'GEOGRID.TBL',form='formatted',status='old',err=1000)
92 !
93 ! We will first go through the file to determine how many field
94 ! specifications there are
95 !
96 10 read(funit,'(a)',end=20,err=1001) buffer
97 call despace(buffer)
99 ! Is this line a comment?
100 if (buffer(1:1) == '#') then
102 ! Are we beginning a new field specification?
103 else if (index(buffer,'=====') /= 0) then
104 if (nparams > 0) num_entries = num_entries + 1
105 nparams = 0
107 else
108 eos = index(buffer,'#')
109 if (eos /= 0) buffer(eos:BUFSIZE) = ' '
111 ! Does this line contain at least one parameter specification?
112 if (index(buffer,'=') /= 0) then
113 nparams = nparams + 1
114 end if
115 end if
116 go to 10
118 20 rewind(funit)
120 ! In case the last entry ended without a ======== line
121 if (nparams > 0) num_entries = num_entries + 1
123 call mprintf(.true.,STDOUT,'Parsed %i entries in GEOGRID.TBL', i1=num_entries)
125 !
126 ! Now that we know how many fields the user has specified, allocate
127 ! the properly sized arrays
128 !
129 allocate(source_wordsize(num_entries))
130 allocate(source_endian(num_entries))
131 allocate(source_fieldtype(num_entries))
132 allocate(source_dest_fieldtype(num_entries))
133 allocate(source_proj(num_entries))
134 allocate(source_priority(num_entries))
135 allocate(source_dx(num_entries))
136 allocate(source_dy(num_entries))
137 allocate(source_known_x(num_entries))
138 allocate(source_known_y(num_entries))
139 allocate(source_known_lat(num_entries))
140 allocate(source_known_lon(num_entries))
141 allocate(source_truelat1(num_entries))
142 allocate(source_truelat2(num_entries))
143 allocate(source_stdlon(num_entries))
144 allocate(source_fieldname(num_entries))
145 allocate(source_path(num_entries))
146 allocate(source_interp_string(num_entries))
147 allocate(source_tile_x(num_entries))
148 allocate(source_tile_y(num_entries))
149 allocate(source_tile_z(num_entries))
150 allocate(source_tile_z_start(num_entries))
151 allocate(source_tile_z_end(num_entries))
152 allocate(source_category_min(num_entries))
153 allocate(source_category_max(num_entries))
154 allocate(source_tile_bdr(num_entries))
155 allocate(source_masked(num_entries))
156 allocate(source_output_stagger(num_entries))
157 allocate(source_row_order(num_entries))
158 allocate(source_dominant_category(num_entries))
159 allocate(source_dominant_only(num_entries))
160 allocate(source_dfdx(num_entries))
161 allocate(source_dfdy(num_entries))
162 allocate(source_scale_factor(num_entries))
163 allocate(source_z_dim_name(num_entries))
164 allocate(source_units(num_entries))
165 allocate(source_descr(num_entries))
166 allocate(source_smooth_option(num_entries))
167 allocate(source_smooth_passes(num_entries))
168 allocate(source_missing_value(num_entries))
169 allocate(source_fill_missing(num_entries))
170 allocate(source_res_path(num_entries))
171 allocate(source_interp_option(num_entries))
172 allocate(source_landmask_land(num_entries))
173 allocate(source_landmask_water(num_entries))
174 do i=1,num_entries
175 call list_init(source_res_path(i))
176 call list_init(source_interp_option(i))
177 call list_init(source_landmask_land(i))
178 call list_init(source_landmask_water(i))
179 end do
181 allocate(is_wordsize(num_entries))
182 allocate(is_endian(num_entries))
183 allocate(is_fieldtype(num_entries))
184 allocate(is_dest_fieldtype(num_entries))
185 allocate(is_proj(num_entries))
186 allocate(is_priority(num_entries))
187 allocate(is_dx(num_entries))
188 allocate(is_dy(num_entries))
189 allocate(is_known_x(num_entries))
190 allocate(is_known_y(num_entries))
191 allocate(is_known_lat(num_entries))
192 allocate(is_known_lon(num_entries))
193 allocate(is_truelat1(num_entries))
194 allocate(is_truelat2(num_entries))
195 allocate(is_stdlon(num_entries))
196 allocate(is_fieldname(num_entries))
197 allocate(is_path(num_entries))
198 allocate(is_tile_x(num_entries))
199 allocate(is_tile_y(num_entries))
200 allocate(is_tile_z(num_entries))
201 allocate(is_tile_z_start(num_entries))
202 allocate(is_tile_z_end(num_entries))
203 allocate(is_category_min(num_entries))
204 allocate(is_category_max(num_entries))
205 allocate(is_tile_bdr(num_entries))
206 allocate(is_masked(num_entries))
207 allocate(is_halt_missing(num_entries))
208 allocate(is_output_stagger(num_entries))
209 allocate(is_row_order(num_entries))
210 allocate(is_dominant_category(num_entries))
211 allocate(is_dominant_only(num_entries))
212 allocate(is_dfdx(num_entries))
213 allocate(is_dfdy(num_entries))
214 allocate(is_scale_factor(num_entries))
215 allocate(is_z_dim_name(num_entries))
216 allocate(is_units(num_entries))
217 allocate(is_descr(num_entries))
218 allocate(is_smooth_option(num_entries))
219 allocate(is_smooth_passes(num_entries))
220 allocate(is_signed(num_entries))
221 allocate(is_missing_value(num_entries))
222 allocate(is_fill_missing(num_entries))
223 allocate(is_subgrid(num_entries))
225 source_wordsize=0
226 source_endian=0
227 source_fieldtype=0
228 source_dest_fieldtype=0
229 source_proj=0
230 source_priority=0
231 source_dx=0
232 source_dy=0
233 source_known_x=0
234 source_known_y=0
235 source_known_lat=0
236 source_known_lon=0
237 source_truelat1=0
238 source_truelat2=0
239 source_stdlon=0
240 source_fieldname=' '
241 source_path=' '
242 source_interp_string=' '
243 source_tile_x=0
244 source_tile_y=0
245 source_tile_z=0
246 source_tile_z_start=0
247 source_tile_z_end=0
248 source_category_min=0
249 source_category_max=0
250 source_tile_bdr=0
251 source_masked=0
252 source_output_stagger=0
253 source_row_order=0
254 source_dominant_category=' '
255 source_dominant_only=' '
256 source_dfdx=' '
257 source_dfdy=' '
258 source_scale_factor=0
259 source_z_dim_name=' '
260 source_units=' '
261 source_descr=' '
262 source_smooth_option=0
263 source_smooth_passes=0
264 source_missing_value=0
265 source_fill_missing=0
267 is_wordsize=.false.
268 is_endian=.false.
269 is_fieldtype=.false.
270 is_dest_fieldtype=.false.
271 is_proj=.false.
272 is_priority=.false.
273 is_dx=.false.
274 is_dy=.false.
275 is_known_x=.false.
276 is_known_y=.false.
277 is_known_lat=.false.
278 is_known_lon=.false.
279 is_truelat1=.false.
280 is_truelat2=.false.
281 is_stdlon=.false.
282 is_fieldname=.false.
283 is_path=.false.
284 is_tile_x=.false.
285 is_tile_y=.false.
286 is_tile_z=.false.
287 is_tile_z_start=.false.
288 is_tile_z_end=.false.
289 is_category_min=.false.
290 is_category_max=.false.
291 is_tile_bdr=.false.
292 is_masked=.false.
293 is_halt_missing=.false.
294 is_output_stagger=.false.
295 is_row_order=.false.
296 is_dominant_category=.false.
297 is_dominant_only=.false.
298 is_dfdx=.false.
299 is_dfdy=.false.
300 is_scale_factor=.false.
301 is_z_dim_name=.false.
302 is_units=.false.
303 is_descr=.false.
304 is_smooth_option=.false.
305 is_smooth_passes=.false.
306 is_signed=.false.
307 is_missing_value=.false.
308 is_fill_missing=.false.
309 is_subgrid=.false.
311 write(source_mminlu,'(a4)') 'USGS'
312 source_iswater = 16
313 source_islake = -1
314 source_isice = 24
315 source_isurban = 1
316 source_isoilwater = 14
318 !
319 ! Actually read and save the specifications
320 !
321 nparams = 0
322 i = 1
323 30 buffer = ' '
324 read(funit,'(a)',end=40,err=1001) buffer
325 call despace(buffer)
327 ! Is this line a comment?
328 if (buffer(1:1) == '#') then
329 ! Do nothing.
331 ! Are we beginning a new field specification?
332 else if (index(buffer,'=====') /= 0) then !{
333 if (nparams > 0) i = i + 1
334 nparams = 0
335 if (i <= num_entries) then
336 !BUG: Are these initializations needed now that we've added initializations for
337 ! all options after their initial allocation above?
338 is_path(i) = .false.
339 is_masked(i) = .false.
340 is_halt_missing(i) = .false.
341 is_output_stagger(i) = .false.
342 is_dominant_category(i) = .false.
343 is_dominant_only(i) = .false.
344 is_dfdx(i) = .false.
345 is_dfdy(i) = .false.
346 is_dest_fieldtype(i) = .false.
347 is_priority(i) = .false.
348 is_z_dim_name(i) = .false.
349 is_smooth_option(i) = .false.
350 is_smooth_passes(i) = .false.
351 is_fill_missing(i) = .false.
352 is_subgrid(i) = .false.
353 end if
355 else
356 ! Check whether the current line is a comment
357 if (buffer(1:1) /= '#') then
358 have_specification = .true.
359 else
360 have_specification = .false.
361 end if
363 ! If only part of the line is a comment, just turn the comment into spaces
364 eos = index(buffer,'#')
365 if (eos /= 0) buffer(eos:BUFSIZE) = ' '
367 do while (have_specification) !{
369 ! If this line has no semicolon, it may contain a single specification,
370 ! so we set have_specification = .false. to prevent the line from being
371 ! processed again and "pretend" that the last character was a semicolon
372 eos = index(buffer,';')
373 if (eos == 0) then
374 have_specification = .false.
375 eos = BUFSIZE
376 end if
378 idx = index(buffer(1:eos-1),'=')
380 if (idx /= 0) then !{
381 nparams = nparams + 1
383 if (index('name',trim(buffer(1:idx-1))) /= 0) then
384 ispace = idx+1
385 do while ((ispace < eos) .and. (buffer(ispace:ispace) /= ' '))
386 ispace = ispace + 1
387 end do
388 is_fieldname(i) = .true.
389 source_fieldname(i) = ' '
390 source_fieldname(i)(1:ispace-idx) = buffer(idx+1:ispace-1)
392 else if (index('priority',trim(buffer(1:idx-1))) /= 0) then
393 is_priority(i) = .true.
394 read(buffer(idx+1:eos-1),'(i10)') source_priority(i)
396 else if (index('dest_type',trim(buffer(1:idx-1))) /= 0) then
397 if (index('continuous',trim(buffer(idx+1:eos-1))) /= 0) then
398 is_dest_fieldtype(i) = .true.
399 source_dest_fieldtype(i) = CONTINUOUS
400 else if (index('categorical',trim(buffer(idx+1:eos-1))) /= 0) then
401 is_dest_fieldtype(i) = .true.
402 source_dest_fieldtype(i) = CATEGORICAL
403 end if
405 else if (index('interp_option',trim(buffer(1:idx-1))) /= 0) then
406 ispace = idx+1
407 do while ((ispace < eos) .and. (buffer(ispace:ispace) /= ' '))
408 ispace = ispace + 1
409 end do
410 interp_string = ' '
411 interp_string(1:ispace-idx-1) = buffer(idx+1:ispace-1)
412 ispace = index(interp_string,':')
413 if (ispace /= 0) then
414 write(res_string,'(a)') interp_string(1:ispace-1)
415 else
416 res_string = 'default'
417 end if
418 write(interp_string,'(a)') trim(interp_string(ispace+1:128))
419 if (list_search(source_interp_option(i), ckey=res_string, cvalue=interp_string)) then
420 call mprintf(.true., WARN, &
421 'In GEOGRID.TBL entry %i, multiple interpolation methods are '// &
422 'given for resolution %s. %s will be used.', &
423 i1=i, s1=trim(res_string), s2=trim(interp_string))
424 else
425 call list_insert(source_interp_option(i), ckey=res_string, cvalue=interp_string)
426 end if
428 else if (index('smooth_option',trim(buffer(1:idx-1))) /= 0) then
429 if ((index('1-2-1',trim(buffer(idx+1:eos-1))) /= 0) .and. &
430 (len_trim(buffer(idx+1:eos-1)) == len('1-2-1'))) then
431 is_smooth_option(i) = .true.
432 source_smooth_option(i) = ONETWOONE
433 else if ((index('smth-desmth',trim(buffer(idx+1:eos-1))) /= 0) .and. &
434 (len_trim(buffer(idx+1:eos-1)) == len('smth-desmth'))) then
435 is_smooth_option(i) = .true.
436 source_smooth_option(i) = SMTHDESMTH
437 else if ((index('smth-desmth_special',trim(buffer(idx+1:eos-1))) /= 0) .and. &
438 (len_trim(buffer(idx+1:eos-1)) == len('smth-desmth_special'))) then
439 is_smooth_option(i) = .true.
440 source_smooth_option(i) = SMTHDESMTH_SPECIAL
441 end if
443 else if (index('smooth_passes',trim(buffer(1:idx-1))) /= 0) then
444 is_smooth_passes(i) = .true.
445 read(buffer(idx+1:eos-1),'(i10)') source_smooth_passes(i)
447 else if (index('rel_path',trim(buffer(1:idx-1))) /= 0) then
448 ispace = idx+1
449 do while ((ispace < eos) .and. (buffer(ispace:ispace) /= ' '))
450 ispace = ispace + 1
451 end do
452 path_string = ' '
453 path_string(1:ispace-idx-1) = buffer(idx+1:ispace-1)
454 if (path_string(ispace-idx-1:ispace-idx-1) /= '/') &
455 path_string(ispace-idx:ispace-idx) = '/'
456 if (path_string(1:1) == '/') then
457 path_string(1:127) = path_string(2:128)
458 path_string(128:128) = ' '
459 end if
460 ispace = index(path_string,':')
461 if (ispace /= 0) then
462 write(res_string,'(a)') path_string(1:ispace-1)
463 else
464 res_string = 'default'
465 end if
466 write(path_string,'(a)') trim(geog_data_path)//trim(path_string(ispace+1:128))
467 if (list_search(source_res_path(i), ckey=res_string, cvalue=path_string)) then
468 call mprintf(.true., WARN, &
469 'In GEOGRID.TBL entry %i, multiple paths are given for '// &
470 'resolution %s. %s will be used.', &
471 i1=i, s1=trim(res_string), s2=trim(path_string))
472 else
473 call list_insert(source_res_path(i), ckey=res_string, cvalue=path_string)
474 end if
476 else if (index('abs_path',trim(buffer(1:idx-1))) /= 0) then
477 ispace = idx+1
478 do while ((ispace < eos) .and. (buffer(ispace:ispace) /= ' '))
479 ispace = ispace + 1
480 end do
481 path_string = ' '
482 path_string(1:ispace-idx-1) = buffer(idx+1:ispace-1)
483 if (path_string(ispace-idx-1:ispace-idx-1) /= '/') &
484 path_string(ispace-idx:ispace-idx) = '/'
485 ispace = index(path_string,':')
486 if (ispace /= 0) then
487 write(res_string,'(a)') path_string(1:ispace-1)
488 else
489 res_string = 'default'
490 end if
491 write(path_string,'(a)') trim(path_string(ispace+1:128))
492 if (list_search(source_res_path(i), ckey=res_string, cvalue=path_string)) then
493 call mprintf(.true., WARN, &
494 'In GEOGRID.TBL entry %i, multiple paths are given for '// &
495 'resolution %s. %s will be used.', &
496 i1=i, s1=trim(res_string), s2=trim(path_string))
497 else
498 call list_insert(source_res_path(i), ckey=res_string, cvalue=path_string)
499 end if
501 else if (index('output_stagger',trim(buffer(1:idx-1))) /= 0) then
502 if (index('M',trim(buffer(idx+1:eos-1))) /= 0) then
503 is_output_stagger(i) = .true.
504 source_output_stagger(i) = M
505 else if (index('U',trim(buffer(idx+1:eos-1))) /= 0) then
506 is_output_stagger(i) = .true.
507 source_output_stagger(i) = U
508 else if (index('V',trim(buffer(idx+1:eos-1))) /= 0) then
509 is_output_stagger(i) = .true.
510 source_output_stagger(i) = V
511 else if (index('HH',trim(buffer(idx+1:eos-1))) /= 0) then
512 is_output_stagger(i) = .true.
513 source_output_stagger(i) = HH
514 else if (index('VV',trim(buffer(idx+1:eos-1))) /= 0) then
515 is_output_stagger(i) = .true.
516 source_output_stagger(i) = VV
517 end if
519 else if ((index('landmask_water',trim(buffer(1:idx-1))) /= 0) .and. &
520 (len_trim(buffer(1:idx-1)) == 14)) then
521 ispace = idx+1
522 do while ((ispace < eos) .and. (buffer(ispace:ispace) /= ' '))
523 ispace = ispace + 1
524 end do
525 landmask_string = ' '
526 landmask_string(1:ispace-idx-1) = buffer(idx+1:ispace-1)
527 ispace = index(landmask_string,':')
528 if (ispace /= 0) then
529 write(res_string,'(a)') landmask_string(1:ispace-1)
530 else
531 res_string = 'default'
532 end if
533 write(landmask_string,'(a)') trim(landmask_string(ispace+1:128))
534 if (list_search(source_landmask_water(i), ckey=res_string, cvalue=landmask_string)) then
535 call mprintf(.true., WARN, &
536 'In GEOGRID.TBL entry %i, multiple landmask category specifications are given for '// &
537 'resolution %s. %s will be used.', &
538 i1=i, s1=trim(res_string), s2=trim(landmask_string))
539 else
540 call list_insert(source_landmask_water(i), ckey=res_string, cvalue=landmask_string)
541 end if
543 else if ((index('landmask_land',trim(buffer(1:idx-1))) /= 0) .and. &
544 (len_trim(buffer(1:idx-1)) == 13)) then
545 ispace = idx+1
546 do while ((ispace < eos) .and. (buffer(ispace:ispace) /= ' '))
547 ispace = ispace + 1
548 end do
549 landmask_string = ' '
550 landmask_string(1:ispace-idx-1) = buffer(idx+1:ispace-1)
551 ispace = index(landmask_string,':')
552 if (ispace /= 0) then
553 write(res_string,'(a)') landmask_string(1:ispace-1)
554 else
555 res_string = 'default'
556 end if
557 write(landmask_string,'(a)') trim(landmask_string(ispace+1:128))
558 if (list_search(source_landmask_land(i), ckey=res_string, cvalue=landmask_string)) then
559 call mprintf(.true., WARN, &
560 'In GEOGRID.TBL entry %i, multiple landmask category specifications are given for '// &
561 'resolution %s. %s will be used.', &
562 i1=i, s1=trim(res_string), s2=trim(landmask_string))
563 else
564 call list_insert(source_landmask_land(i), ckey=res_string, cvalue=landmask_string)
565 end if
567 else if ((index('masked',trim(buffer(1:idx-1))) /= 0) .and. &
568 (len_trim(buffer(1:idx-1)) == 6)) then
569 if (index('water',trim(buffer(idx+1:eos-1))) /= 0) then
570 is_masked(i) = .true.
571 source_masked(i) = 0.
572 else if (index('land',trim(buffer(idx+1:eos-1))) /= 0) then
573 is_masked(i) = .true.
574 source_masked(i) = 1.
575 end if
577 else if (index('fill_missing',trim(buffer(1:idx-1))) /= 0) then
578 is_fill_missing(i) = .true.
579 read(buffer(idx+1:eos-1),*) source_fill_missing(i)
581 else if (index('halt_on_missing',trim(buffer(1:idx-1))) /= 0) then
582 if (index('yes',trim(buffer(idx+1:eos-1))) /= 0) then
583 is_halt_missing(i) = .true.
584 else if (index('no',trim(buffer(idx+1:eos-1))) /= 0) then
585 is_halt_missing(i) = .false.
586 end if
588 else if (index('dominant_category',trim(buffer(1:idx-1))) /= 0) then
589 ispace = idx+1
590 do while ((ispace < eos) .and. (buffer(ispace:ispace) /= ' '))
591 ispace = ispace + 1
592 end do
593 is_dominant_category(i) = .true.
594 source_dominant_category(i) = ' '
595 source_dominant_category(i)(1:ispace-idx) = buffer(idx+1:ispace-1)
597 else if (index('dominant_only',trim(buffer(1:idx-1))) /= 0) then
598 ispace = idx+1
599 do while ((ispace < eos) .and. (buffer(ispace:ispace) /= ' '))
600 ispace = ispace + 1
601 end do
602 is_dominant_only(i) = .true.
603 source_dominant_only(i) = ' '
604 source_dominant_only(i)(1:ispace-idx) = buffer(idx+1:ispace-1)
606 else if (index('df_dx',trim(buffer(1:idx-1))) /= 0) then
607 ispace = idx+1
608 do while ((ispace < eos) .and. (buffer(ispace:ispace) /= ' '))
609 ispace = ispace + 1
610 end do
611 is_dfdx(i) = .true.
612 source_dfdx(i) = ' '
613 source_dfdx(i)(1:ispace-idx) = buffer(idx+1:ispace-1)
615 else if (index('df_dy',trim(buffer(1:idx-1))) /= 0) then
616 ispace = idx+1
617 do while ((ispace < eos) .and. (buffer(ispace:ispace) /= ' '))
618 ispace = ispace + 1
619 end do
620 is_dfdy(i) = .true.
621 source_dfdy(i) = ' '
622 source_dfdy(i)(1:ispace-idx) = buffer(idx+1:ispace-1)
624 else if (index('z_dim_name',trim(buffer(1:idx-1))) /= 0) then
625 ispace = idx+1
626 do while ((ispace < eos) .and. (buffer(ispace:ispace) /= ' '))
627 ispace = ispace + 1
628 end do
629 is_z_dim_name(i) = .true.
630 source_z_dim_name(i) = ' '
631 source_z_dim_name(i)(1:ispace-idx) = buffer(idx+1:ispace-1)
633 else if (index('subgrid',trim(buffer(1:idx-1))) /= 0) then
634 if (index('yes',trim(buffer(idx+1:eos-1))) /= 0) then
635 is_subgrid(i) = .true.
636 else if (index('no',trim(buffer(idx+1:eos-1))) /= 0) then
637 is_subgrid(i) = .false.
638 end if
640 else
641 call mprintf(.true., WARN, 'In GEOGRID.TBL, unrecognized option %s in '// &
642 'entry %i.',i1=idx, s1=buffer(i:idx-1))
644 end if
646 end if !} index(buffer(1:eos-1),'=') /= 0
648 buffer = buffer(eos+1:BUFSIZE)
649 end do ! while eos /= 0 }
651 end if !} index(buffer, '=====') /= 0
652 go to 30
654 40 close(funit)
656 ! Check the user specifications for gross errors
657 if ( .not. check_data_specification() ) then
658 call datalist_destroy()
659 call mprintf(.true.,ERROR,'Errors were found in either index files or GEOGRID.TBL.')
660 end if
662 call hash_init(bad_files)
664 return
666 1000 call mprintf(.true.,ERROR,'Could not open GEOGRID.TBL')
668 1001 call mprintf(.true.,ERROR,'Encountered error while reading GEOGRID.TBL')
670 end subroutine get_datalist
673 !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
674 ! Name: get_source_params
675 !
676 ! Purpose: For each field, this routine reads in the metadata in the index file
677 ! for the first available resolution of data specified by res_string. Also
678 ! based on res_string, this routine sets the interpolation sequence for the
679 ! field. This routine should be called prior to processing a field for each
680 ! domain.
681 !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
682 subroutine get_source_params(res_string)
684 use stringutil
686 implicit none
688 ! Parameters
689 integer, parameter :: BUFSIZE = 256
691 ! Arguments
692 character (len=128), intent(in) :: res_string
694 ! Local variables
695 integer :: idx, i, is, ie, ispace, eos, iquoted, funit
696 character (len=128) :: temp_data, temp_interp
697 character (len=256) :: test_fname
698 character (len=BUFSIZE) :: buffer
699 logical :: have_specification, is_used
701 ! For each entry in the GEOGRID.TBL file
702 do idx=1,num_entries
704 ! Initialize metadata
705 is_wordsize(idx) = .false.
706 is_endian(idx) = .false.
707 is_row_order(idx) = .false.
708 is_fieldtype(idx) = .false.
709 is_proj(idx) = .false.
710 is_dx(idx) = .false.
711 is_dy(idx) = .false.
712 is_known_x(idx) = .false.
713 is_known_y(idx) = .false.
714 is_known_lat(idx) = .false.
715 is_known_lon(idx) = .false.
716 is_truelat1(idx) = .false.
717 is_truelat2(idx) = .false.
718 is_stdlon(idx) = .false.
719 is_tile_x(idx) = .false.
720 is_tile_y(idx) = .false.
721 is_tile_z(idx) = .false.
722 is_tile_z_start(idx) = .false.
723 is_tile_z_end(idx) = .false.
724 is_category_min(idx) = .false.
725 is_category_max(idx) = .false.
726 is_tile_bdr(idx) = .false.
727 is_fieldname(idx) = .false.
728 is_scale_factor(idx) = .false.
729 is_units(idx) = .false.
730 is_descr(idx) = .false.
731 is_signed(idx) = .false.
732 is_missing_value(idx) = .false.
735 ! Set the interpolator sequence for the field to be the first value in res_string that matches
736 ! the resolution keyword for an interp_sequence specification
737 is = 1
738 ie = index(res_string(is:128),'+') - 1
739 if (ie <= 0) ie = 128
740 temp_interp = res_string(is:ie)
741 do while (.not. list_search(source_interp_option(idx), ckey=temp_interp, cvalue=source_interp_string(idx)) &
742 .and. is <= 128)
743 call mprintf(.true., INFORM, 'For %s, couldn''t find interpolator sequence for '// &
744 'resolution %s.', &
745 s1=trim(source_fieldname(idx)), s2=trim(temp_interp))
746 is = ie+2
747 ie = is + index(res_string(is:128),'+') - 2
748 if (ie - is <= 0) ie = 128
749 temp_interp = res_string(is:ie)
750 end do
752 if (is > 128) then
753 temp_interp = 'default'
754 if (list_search(source_interp_option(idx), ckey=temp_interp, cvalue=source_interp_string(idx))) then
755 call mprintf(.true., INFORM, 'Using default interpolator sequence for %s.', &
756 s1=trim(source_fieldname(idx)))
757 else
758 call mprintf(.true., ERROR, 'Could not find interpolator sequence for requested resolution for %s.'// &
759 ' The sources specified in namelist.wps is not listed in GEOGRID.TBL.', &
760 s1=trim(source_fieldname(idx)))
761 end if
762 else
763 call mprintf(.true., INFORM, 'Using %s interpolator sequence for %s.', &
764 s1=temp_interp, s2=trim(source_fieldname(idx)))
765 end if
767 ! Set the data source for the field to be the first value in res_string that matches
768 ! the resolution keyword for an abs_path or rel_path specification
769 is = 1
770 ie = index(res_string(is:128),'+') - 1
771 if (ie <= 0) ie = 128
772 temp_data = res_string(is:ie)
773 do while (.not. list_search(source_res_path(idx), ckey=temp_data, cvalue=source_path(idx)) &
774 .and. is <= 128)
775 call mprintf(.true., INFORM, 'For %s, couldn''t find %s data source.', &
776 s1=trim(source_fieldname(idx)), s2=trim(temp_data))
777 is = ie+2
778 ie = is + index(res_string(is:128),'+') - 2
779 if (ie - is <= 0) ie = 128
780 temp_data = res_string(is:ie)
781 end do
783 if (is > 128) then
784 temp_data = 'default'
785 if (list_search(source_res_path(idx), ckey=temp_data, cvalue=source_path(idx))) then
786 call mprintf(.true., INFORM, 'Using default data source for %s.', &
787 s1=trim(source_fieldname(idx)))
788 else
789 call mprintf(.true., ERROR, 'Could not find data resolution for requested resolution for %s. '// &
790 'The source specified in namelist.wps is not listed in GEOGRID.TBL.', &
791 s1=trim(source_fieldname(idx)))
792 end if
793 else
794 call mprintf(.true., INFORM, 'Using %s data source for %s.', &
795 s1=temp_data, s2=trim(source_fieldname(idx)))
796 end if
798 call mprintf(trim(temp_data) /= trim(temp_interp),WARN,'For %s, using %s data source with %s interpolation sequence.', &
799 s1=source_fieldname(idx), s2=temp_data, s3=temp_interp)
801 write(test_fname, '(a)') trim(source_path(idx))//'index'
803 !
804 ! Open the index file for the data source for this field, and read in metadata specs
805 !
806 do funit=10,100
807 inquire(unit=funit, opened=is_used)
808 if (.not. is_used) exit
809 end do
810 open(funit,file=trim(test_fname),form='formatted',status='old',err=1000)
812 30 buffer = ' '
813 read(funit,'(a)',end=40, err=1001) buffer
814 call despace(buffer)
816 ! Is this line a comment?
817 if (buffer(1:1) == '#') then
818 ! Do nothing.
820 else
821 have_specification = .true.
823 ! If only part of the line is a comment, just turn the comment into spaces
824 eos = index(buffer,'#')
825 if (eos /= 0) buffer(eos:BUFSIZE) = ' '
827 do while (have_specification) !{
829 ! If this line has no semicolon, it may contain a single specification,
830 ! so we set have_specification = .false. to prevent the line from being
831 ! processed again and pretend that the last character was a semicolon
832 eos = index(buffer,';')
833 if (eos == 0) then
834 have_specification = .false.
835 eos = BUFSIZE
836 end if
838 i = index(buffer(1:eos-1),'=')
840 if (i /= 0) then !{
842 if (index('projection',trim(buffer(1:i-1))) /= 0) then
843 if (index('lambert',trim(buffer(i+1:eos-1))) /= 0) then
844 is_proj(idx) = .true.
845 source_proj(idx) = PROJ_LC
846 else if (index('polar_wgs84',trim(buffer(i+1:eos-1))) /= 0 .and. &
847 len_trim('polar_wgs84') == len_trim(buffer(i+1:eos-1))) then
848 is_proj(idx) = .true.
849 source_proj(idx) = PROJ_PS_WGS84
850 else if (index('albers_nad83',trim(buffer(i+1:eos-1))) /= 0 .and. &
851 len_trim('albers_nad83') == len_trim(buffer(i+1:eos-1))) then
852 is_proj(idx) = .true.
853 source_proj(idx) = PROJ_ALBERS_NAD83
854 else if (index('polar',trim(buffer(i+1:eos-1))) /= 0 .and. &
855 len_trim('polar') == len_trim(buffer(i+1:eos-1))) then
856 is_proj(idx) = .true.
857 source_proj(idx) = PROJ_PS
858 else if (index('mercator',trim(buffer(i+1:eos-1))) /= 0) then
859 is_proj(idx) = .true.
860 source_proj(idx) = PROJ_MERC
861 else if (index('regular_ll',trim(buffer(i+1:eos-1))) /= 0) then
862 is_proj(idx) = .true.
863 source_proj(idx) = PROJ_LATLON
864 end if
866 else if (index('type',trim(buffer(1:i-1))) /= 0) then
867 if (index('continuous',trim(buffer(i+1:eos-1))) /= 0) then
868 is_fieldtype(idx) = .true.
869 source_fieldtype(idx) = CONTINUOUS
870 else if (index('categorical',trim(buffer(i+1:eos-1))) /= 0) then
871 is_fieldtype(idx) = .true.
872 source_fieldtype(idx) = CATEGORICAL
873 end if
875 else if (index('signed',trim(buffer(1:i-1))) /= 0) then
876 if (index('yes',trim(buffer(i+1:eos-1))) /= 0) then
877 is_signed(idx) = .true.
878 else if (index('no',trim(buffer(i+1:eos-1))) /= 0) then
879 is_signed(idx) = .false.
880 end if
882 else if (index('units',trim(buffer(1:i-1))) /= 0) then
883 ispace = i+1
884 iquoted = 0
885 do while (((ispace < eos) .and. (buffer(ispace:ispace) /= ' ')) .or. (iquoted == 1))
886 if (buffer(ispace:ispace) == '"' .or. buffer(ispace:ispace) == '''') iquoted = mod(iquoted+1,2)
887 ispace = ispace + 1
888 end do
889 is_units(idx) = .true.
890 source_units(idx) = ' '
891 if (buffer(i+1:i+1) == '"' .or. buffer(i+1:i+1) == '''') i = i + 1
892 if (buffer(ispace-1:ispace-1) == '"' .or. buffer(ispace-1:ispace-1) == '''') ispace = ispace - 1
893 source_units(idx)(1:ispace-i) = buffer(i+1:ispace-1)
895 else if (index('description',trim(buffer(1:i-1))) /= 0) then
896 ispace = i+1
897 iquoted = 0
898 do while (((ispace < eos) .and. (buffer(ispace:ispace) /= ' ')) .or. (iquoted == 1))
899 if (buffer(ispace:ispace) == '"' .or. buffer(ispace:ispace) == '''') iquoted = mod(iquoted+1,2)
900 ispace = ispace + 1
901 end do
902 is_descr(idx) = .true.
903 source_descr(idx) = ' '
904 if (buffer(i+1:i+1) == '"' .or. buffer(i+1:i+1) == '''') i = i + 1
905 if (buffer(ispace-1:ispace-1) == '"' .or. buffer(ispace-1:ispace-1) == '''') ispace = ispace - 1
906 source_descr(idx)(1:ispace-i) = buffer(i+1:ispace-1)
908 else if (index('mminlu',trim(buffer(1:i-1))) /= 0) then
909 ispace = i+1
910 iquoted = 0
911 do while (((ispace < eos) .and. (buffer(ispace:ispace) /= ' ')) .or. (iquoted == 1))
912 if (buffer(ispace:ispace) == '"' .or. buffer(ispace:ispace) == '''') iquoted = mod(iquoted+1,2)
913 ispace = ispace + 1
914 end do
915 if (buffer(i+1:i+1) == '"' .or. buffer(i+1:i+1) == '''') i = i + 1
916 if (buffer(ispace-1:ispace-1) == '"' .or. buffer(ispace-1:ispace-1) == '''') ispace = ispace - 1
917 source_mminlu(1:ispace-i) = buffer(i+1:ispace-1)
919 else if (index('iswater',trim(buffer(1:i-1))) /= 0) then
920 read(buffer(i+1:eos-1),*) source_iswater
922 else if (index('islake',trim(buffer(1:i-1))) /= 0) then
923 read(buffer(i+1:eos-1),*) source_islake
925 else if (index('isice',trim(buffer(1:i-1))) /= 0) then
926 read(buffer(i+1:eos-1),*) source_isice
928 else if (index('isurban',trim(buffer(1:i-1))) /= 0) then
929 read(buffer(i+1:eos-1),*) source_isurban
931 else if (index('isoilwater',trim(buffer(1:i-1))) /= 0) then
932 read(buffer(i+1:eos-1),*) source_isoilwater
934 else if (index('dx',trim(buffer(1:i-1))) /= 0) then
935 is_dx(idx) = .true.
936 read(buffer(i+1:eos-1),*) source_dx(idx)
938 else if (index('dy',trim(buffer(1:i-1))) /= 0) then
939 is_dy(idx) = .true.
940 read(buffer(i+1:eos-1),*) source_dy(idx)
942 else if (index('known_x',trim(buffer(1:i-1))) /= 0) then
943 is_known_x(idx) = .true.
944 read(buffer(i+1:eos-1),*) source_known_x(idx)
946 else if (index('known_y',trim(buffer(1:i-1))) /= 0) then
947 is_known_y(idx) = .true.
948 read(buffer(i+1:eos-1),*) source_known_y(idx)
950 else if (index('known_lat',trim(buffer(1:i-1))) /= 0) then
951 is_known_lat(idx) = .true.
952 read(buffer(i+1:eos-1),*) source_known_lat(idx)
954 else if (index('known_lon',trim(buffer(1:i-1))) /= 0) then
955 is_known_lon(idx) = .true.
956 read(buffer(i+1:eos-1),*) source_known_lon(idx)
958 else if (index('stdlon',trim(buffer(1:i-1))) /= 0) then
959 is_stdlon(idx) = .true.
960 read(buffer(i+1:eos-1),*) source_stdlon(idx)
962 else if (index('truelat1',trim(buffer(1:i-1))) /= 0) then
963 is_truelat1(idx) = .true.
964 read(buffer(i+1:eos-1),*) source_truelat1(idx)
966 else if (index('truelat2',trim(buffer(1:i-1))) /= 0) then
967 is_truelat2(idx) = .true.
968 read(buffer(i+1:eos-1),*) source_truelat2(idx)
970 else if (index('wordsize',trim(buffer(1:i-1))) /= 0) then
971 is_wordsize(idx) = .true.
972 read(buffer(i+1:eos-1),'(i10)') source_wordsize(idx)
974 else if (index('endian',trim(buffer(1:i-1))) /= 0) then
975 if (index('big',trim(buffer(i+1:eos-1))) /= 0) then
976 is_endian(idx) = .true.
977 source_endian(idx) = BIG_ENDIAN
978 else if (index('little',trim(buffer(i+1:eos-1))) /= 0) then
979 is_endian(idx) = .true.
980 source_endian(idx) = LITTLE_ENDIAN
981 else
982 call mprintf(.true.,WARN,'Invalid value for keyword ''endian'' '// &
983 'specified in index file. BIG_ENDIAN will be used.')
984 end if
986 else if (index('row_order',trim(buffer(1:i-1))) /= 0) then
987 if (index('bottom_top',trim(buffer(i+1:eos-1))) /= 0) then
988 is_row_order(idx) = .true.
989 source_row_order(idx) = BOTTOM_TOP
990 else if (index('top_bottom',trim(buffer(i+1:eos-1))) /= 0) then
991 is_row_order(idx) = .true.
992 source_row_order(idx) = TOP_BOTTOM
993 end if
995 else if (index('tile_x',trim(buffer(1:i-1))) /= 0) then
996 is_tile_x(idx) = .true.
997 read(buffer(i+1:eos-1),'(i10)') source_tile_x(idx)
999 else if (index('tile_y',trim(buffer(1:i-1))) /= 0) then
1000 is_tile_y(idx) = .true.
1001 read(buffer(i+1:eos-1),'(i10)') source_tile_y(idx)
1003 else if (index('tile_z',trim(buffer(1:i-1))) /= 0) then
1004 is_tile_z(idx) = .true.
1005 read(buffer(i+1:eos-1),'(i10)') source_tile_z(idx)
1007 else if (index('tile_z_start',trim(buffer(1:i-1))) /= 0) then
1008 is_tile_z_start(idx) = .true.
1009 read(buffer(i+1:eos-1),'(i10)') source_tile_z_start(idx)
1011 else if (index('tile_z_end',trim(buffer(1:i-1))) /= 0) then
1012 is_tile_z_end(idx) = .true.
1013 read(buffer(i+1:eos-1),'(i10)') source_tile_z_end(idx)
1015 else if (index('category_min',trim(buffer(1:i-1))) /= 0) then
1016 is_category_min(idx) = .true.
1017 read(buffer(i+1:eos-1),'(i10)') source_category_min(idx)
1019 else if (index('category_max',trim(buffer(1:i-1))) /= 0) then
1020 is_category_max(idx) = .true.
1021 read(buffer(i+1:eos-1),'(i10)') source_category_max(idx)
1023 else if (index('tile_bdr',trim(buffer(1:i-1))) /= 0) then
1024 is_tile_bdr(idx) = .true.
1025 read(buffer(i+1:eos-1),'(i10)') source_tile_bdr(idx)
1027 else if (index('missing_value',trim(buffer(1:i-1))) /= 0) then
1028 is_missing_value(idx) = .true.
1029 read(buffer(i+1:eos-1),*) source_missing_value(idx)
1031 else if (index('scale_factor',trim(buffer(1:i-1))) /= 0) then
1032 is_scale_factor(idx) = .true.
1033 read(buffer(i+1:eos-1),*) source_scale_factor(idx)
1035 else
1036 call mprintf(.true., WARN, 'In %s, unrecognized option %s in entry %i.', &
1037 s1=trim(test_fname), s2=buffer(i:i-1), i1=i)
1039 end if
1041 end if !} index(buffer(1:eos-1),'=') /= 0
1043 buffer = buffer(eos+1:BUFSIZE)
1044 end do ! while eos /= 0 }
1046 end if !} index(buffer, '=====') /= 0
1048 go to 30
1050 40 close(funit)
1052 end do
1054 return
1056 1000 call mprintf(.true.,ERROR,'Could not open %s', s1=trim(test_fname))
1057 1001 call mprintf(.true.,ERROR,'Encountered error while reading %s', s1=trim(test_fname))
1059 end subroutine get_source_params
1062 !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
1063 ! Name: datalist_destroy()
1064 !
1065 ! Purpose: This routine deallocates any memory that was allocated by the
1066 ! get_datalist() subroutine.
1067 !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
1068 subroutine datalist_destroy()
1070 implicit none
1072 ! Local variables
1073 integer :: i
1075 if (associated(source_wordsize)) then
1076 deallocate(source_wordsize)
1077 deallocate(source_endian)
1078 deallocate(source_fieldtype)
1079 deallocate(source_dest_fieldtype)
1080 deallocate(source_proj)
1081 deallocate(source_priority)
1082 deallocate(source_dx)
1083 deallocate(source_dy)
1084 deallocate(source_known_x)
1085 deallocate(source_known_y)
1086 deallocate(source_known_lat)
1087 deallocate(source_known_lon)
1088 deallocate(source_truelat1)
1089 deallocate(source_truelat2)
1090 deallocate(source_stdlon)
1091 deallocate(source_fieldname)
1092 deallocate(source_path)
1093 deallocate(source_interp_string)
1094 deallocate(source_tile_x)
1095 deallocate(source_tile_y)
1096 deallocate(source_tile_z)
1097 deallocate(source_tile_z_start)
1098 deallocate(source_tile_z_end)
1099 deallocate(source_tile_bdr)
1100 deallocate(source_category_min)
1101 deallocate(source_category_max)
1102 deallocate(source_masked)
1103 deallocate(source_output_stagger)
1104 deallocate(source_row_order)
1105 deallocate(source_dominant_category)
1106 deallocate(source_dominant_only)
1107 deallocate(source_dfdx)
1108 deallocate(source_dfdy)
1109 deallocate(source_scale_factor)
1110 deallocate(source_z_dim_name)
1111 deallocate(source_smooth_option)
1112 deallocate(source_smooth_passes)
1113 deallocate(source_units)
1114 deallocate(source_descr)
1115 deallocate(source_missing_value)
1116 deallocate(source_fill_missing)
1117 do i=1,num_entries
1118 call list_destroy(source_res_path(i))
1119 call list_destroy(source_interp_option(i))
1120 call list_destroy(source_landmask_land(i))
1121 call list_destroy(source_landmask_water(i))
1122 end do
1123 deallocate(source_res_path)
1124 deallocate(source_interp_option)
1125 deallocate(source_landmask_land)
1126 deallocate(source_landmask_water)
1128 deallocate(is_wordsize)
1129 deallocate(is_endian)
1130 deallocate(is_fieldtype)
1131 deallocate(is_dest_fieldtype)
1132 deallocate(is_proj)
1133 deallocate(is_priority)
1134 deallocate(is_dx)
1135 deallocate(is_dy)
1136 deallocate(is_known_x)
1137 deallocate(is_known_y)
1138 deallocate(is_known_lat)
1139 deallocate(is_known_lon)
1140 deallocate(is_truelat1)
1141 deallocate(is_truelat2)
1142 deallocate(is_stdlon)
1143 deallocate(is_fieldname)
1144 deallocate(is_path)
1145 deallocate(is_tile_x)
1146 deallocate(is_tile_y)
1147 deallocate(is_tile_z)
1148 deallocate(is_tile_z_start)
1149 deallocate(is_tile_z_end)
1150 deallocate(is_tile_bdr)
1151 deallocate(is_category_min)
1152 deallocate(is_category_max)
1153 deallocate(is_masked)
1154 deallocate(is_halt_missing)
1155 deallocate(is_output_stagger)
1156 deallocate(is_row_order)
1157 deallocate(is_dominant_category)
1158 deallocate(is_dominant_only)
1159 deallocate(is_dfdx)
1160 deallocate(is_dfdy)
1161 deallocate(is_scale_factor)
1162 deallocate(is_z_dim_name)
1163 deallocate(is_smooth_option)
1164 deallocate(is_smooth_passes)
1165 deallocate(is_signed)
1166 deallocate(is_units)
1167 deallocate(is_descr)
1168 deallocate(is_missing_value)
1169 deallocate(is_fill_missing)
1170 end if
1172 call hash_destroy(bad_files)
1174 end subroutine datalist_destroy
1177 !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
1178 ! Name: reset_next_field
1179 !
1180 ! Purpose: To reset the pointer to the next field in the list of fields
1181 ! specified by the user.
1182 !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
1183 subroutine reset_next_field()
1185 implicit none
1187 next_field = 1
1189 end subroutine reset_next_field
1192 !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
1193 ! Name: get_next_fieldname
1194 !
1195 ! Purpose: Calling this routine results in field_name being set to the name of
1196 ! the field currently pointed to. If istatus /= 0 upon return, an error
1197 ! occurred, and the value of field_name is undefined.
1198 !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
1199 subroutine get_next_fieldname(field_name, istatus)
1201 implicit none
1203 ! Arguments
1204 integer, intent(out) :: istatus
1205 character (len=128), intent(out) :: field_name
1207 istatus = 1
1209 if (next_field <= num_entries) then
1211 field_name = source_fieldname(next_field)
1212 next_field = next_field + 1
1213 istatus = 0
1215 end if
1217 end subroutine get_next_fieldname
1220 !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
1221 ! Name: get_next_output_fieldname
1222 !
1223 ! Purpose:
1224 !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
1225 recursive subroutine get_next_output_fieldname(nest_num, field_name, ndims, &
1226 min_cat, max_cat, &
1227 istagger, memorder, dimnames, units, &
1228 description, sr_x, sr_y, istatus)
1230 use gridinfo_module
1232 implicit none
1234 #include "wrf_io_flags.h"
1236 ! Arguments
1237 integer, intent(in) :: nest_num
1238 integer, intent(out) :: istatus, ndims, istagger, min_cat, max_cat
1239 integer, intent(out) :: sr_x, sr_y
1240 character (len=128), intent(out) :: memorder, field_name, units, description
1241 character (len=128), dimension(3), intent(out) :: dimnames
1243 ! Local variables
1244 integer :: temp_fieldtype
1245 integer, dimension(MAX_LANDMASK_CATEGORIES) :: landmask
1246 logical :: is_water_mask, is_dom_only
1247 character (len=128) :: domcat_name, dfdx_name, dfdy_name
1248 character (len=256) :: temphash
1250 istatus = 1
1252 if (output_field_state == RETURN_LANDMASK) then
1253 call hash_init(duplicate_fnames)
1254 call get_landmask_field(geog_data_res(nest_num), field_name, is_water_mask, landmask, istatus)
1255 if (istatus == 0) then
1256 temphash(129:256) = ' '
1257 temphash(1:128) = field_name(1:128)
1258 call hash_insert(duplicate_fnames, temphash)
1259 call get_domcategory_name(field_name, domcat_name, is_dom_only, istatus)
1260 ! We will only save the dominant category
1261 if (is_dom_only .and. (istatus == 0)) then
1262 output_field_state = RETURN_DOMCAT_LM
1263 call get_next_output_fieldname(nest_num, field_name, ndims, &
1264 min_cat, max_cat, istagger, &
1265 memorder, dimnames, units, description, &
1266 sr_x, sr_y, istatus)
1267 return
1268 else
1269 ndims = 2
1270 min_cat = 1
1271 max_cat = 1
1272 temp_fieldtype = iget_fieldtype(field_name, istatus)
1273 if (istatus == 0) then
1274 if (temp_fieldtype == CONTINUOUS) then
1275 call get_max_levels(field_name, min_cat, max_cat, istatus)
1276 else if (temp_fieldtype == CATEGORICAL) then
1277 call get_max_categories(field_name, min_cat, max_cat, istatus)
1278 end if
1279 if (max_cat - min_cat > 0) ndims = 3
1280 end if
1281 call get_output_stagger(field_name, istagger, istatus)
1282 if (istagger == M) then
1283 dimnames(1) = 'west_east'
1284 dimnames(2) = 'south_north'
1285 else if (istagger == U) then
1286 dimnames(1) = 'west_east_stag'
1287 dimnames(2) = 'south_north'
1288 else if (istagger == V) then
1289 dimnames(1) = 'west_east'
1290 dimnames(2) = 'south_north_stag'
1291 else if (istagger == HH) then
1292 dimnames(1) = 'west_east'
1293 dimnames(2) = 'south_north'
1294 else if (istagger == VV) then
1295 dimnames(1) = 'west_east'
1296 dimnames(2) = 'south_north'
1297 end if
1298 if (ndims == 2) then
1299 memorder = 'XY '
1300 dimnames(3) = ' '
1301 else if (ndims == 3) then
1302 memorder = 'XYZ'
1303 call get_z_dim_name(field_name, dimnames(3), istatus)
1304 istatus = 0
1305 else
1306 memorder = ' '
1307 dimnames(3) = ' '
1308 end if
1309 call get_subgrid_dim_name(nest_num, field_name, dimnames(1:2), sr_x, sr_y, istatus)
1310 call get_source_units(field_name, 1, units, istatus)
1311 if (istatus /= 0) units = '-'
1312 call get_source_descr(field_name, 1, description, istatus)
1313 if (istatus /= 0) description = '-'
1314 istatus = 0
1315 output_field_state = RETURN_DOMCAT_LM
1316 end if
1317 else
1318 output_field_state = RETURN_FIELDNAME
1319 call get_next_output_fieldname(nest_num, field_name, ndims, &
1320 min_cat, max_cat, istagger, &
1321 memorder, dimnames, units, description, &
1322 sr_x, sr_y, istatus)
1323 return
1324 end if
1326 else if (output_field_state == RETURN_FIELDNAME) then
1327 call get_next_fieldname(field_name, istatus)
1328 temphash(129:256) = ' '
1329 temphash(1:128) = field_name(1:128)
1330 if (istatus == 0 .and. (.not. hash_search(duplicate_fnames, temphash))) then
1331 call hash_insert(duplicate_fnames, temphash)
1332 call get_domcategory_name(field_name, domcat_name, is_dom_only, istatus)
1333 ! We will only save the dominant category
1334 if (is_dom_only .and. (istatus == 0)) then
1335 output_field_state = RETURN_DOMCAT
1336 call get_next_output_fieldname(nest_num, field_name, ndims, &
1337 min_cat, max_cat, istagger, &
1338 memorder, dimnames, units, description, &
1339 sr_x, sr_y, istatus)
1340 return
1342 ! Return the fractional field
1343 else
1344 ndims = 2
1345 min_cat = 1
1346 max_cat = 1
1347 temp_fieldtype = iget_fieldtype(field_name, istatus)
1348 if (istatus == 0) then
1349 if (temp_fieldtype == CONTINUOUS) then
1350 call get_max_levels(field_name, min_cat, max_cat, istatus)
1351 else if (temp_fieldtype == CATEGORICAL) then
1352 call get_max_categories(field_name, min_cat, max_cat, istatus)
1353 end if
1354 if (max_cat - min_cat > 0) ndims = 3
1355 end if
1356 call get_output_stagger(field_name, istagger, istatus)
1357 if (istagger == M) then
1358 dimnames(1) = 'west_east'
1359 dimnames(2) = 'south_north'
1360 else if (istagger == U) then
1361 dimnames(1) = 'west_east_stag'
1362 dimnames(2) = 'south_north'
1363 else if (istagger == V) then
1364 dimnames(1) = 'west_east'
1365 dimnames(2) = 'south_north_stag'
1366 else if (istagger == HH) then
1367 dimnames(1) = 'west_east'
1368 dimnames(2) = 'south_north'
1369 else if (istagger == VV) then
1370 dimnames(1) = 'west_east'
1371 dimnames(2) = 'south_north'
1372 end if
1373 if (ndims == 2) then
1374 memorder = 'XY '
1375 dimnames(3) = ' '
1376 else if (ndims == 3) then
1377 memorder = 'XYZ'
1378 call get_z_dim_name(field_name, dimnames(3), istatus)
1379 istatus = 0
1380 else
1381 memorder = ' '
1382 dimnames(3) = ' '
1383 end if
1384 call get_subgrid_dim_name(nest_num, field_name, dimnames(1:2), sr_x, sr_y, istatus)
1385 call get_source_units(field_name, 1, units, istatus)
1386 if (istatus /= 0) units = '-'
1387 call get_source_descr(field_name, 1, description, istatus)
1388 if (istatus /= 0) description = '-'
1389 istatus = 0
1390 output_field_state = RETURN_DOMCAT
1391 end if
1392 else if (istatus /= 0) then
1393 output_field_state = RETURN_LANDMASK
1394 call hash_destroy(duplicate_fnames)
1395 return
1396 else if (hash_search(duplicate_fnames, temphash)) then
1397 call get_next_output_fieldname(nest_num, field_name, ndims, &
1398 min_cat, max_cat, istagger, &
1399 memorder, dimnames, units, description, &
1400 sr_x, sr_y, istatus)
1401 return
1402 end if
1404 else if (output_field_state == RETURN_DOMCAT .or. &
1405 output_field_state == RETURN_DOMCAT_LM ) then
1406 if (output_field_state == RETURN_DOMCAT) then
1407 next_field = next_field - 1
1408 call get_next_fieldname(field_name, istatus)
1409 else
1410 call get_landmask_field(geog_data_res(nest_num), field_name, is_water_mask, landmask, istatus)
1411 end if
1412 if (istatus == 0) then
1413 call get_domcategory_name(field_name, domcat_name, is_dom_only, istatus)
1414 if (istatus == 0) then
1415 ndims = 2
1416 min_cat = 1
1417 max_cat = 1
1418 call get_output_stagger(field_name, istagger, istatus)
1419 if (istagger == M) then
1420 dimnames(1) = 'west_east'
1421 dimnames(2) = 'south_north'
1422 else if (istagger == U) then
1423 dimnames(1) = 'west_east_stag'
1424 dimnames(2) = 'south_north'
1425 else if (istagger == V) then
1426 dimnames(1) = 'west_east'
1427 dimnames(2) = 'south_north_stag'
1428 else if (istagger == HH) then
1429 dimnames(1) = 'west_east'
1430 dimnames(2) = 'south_north'
1431 else if (istagger == VV) then
1432 dimnames(1) = 'west_east'
1433 dimnames(2) = 'south_north'
1434 end if
1435 dimnames(3) = ' '
1436 memorder = 'XY '
1438 call get_subgrid_dim_name(nest_num, field_name, dimnames(1:2), sr_x, sr_y, istatus)
1439 field_name = domcat_name
1440 units = 'category'
1441 description = 'Dominant category'
1442 if (output_field_state == RETURN_DOMCAT) then
1443 output_field_state = RETURN_DFDX
1444 else
1445 output_field_state = RETURN_DFDX_LM
1446 end if
1447 else
1448 if (output_field_state == RETURN_DOMCAT) then
1449 output_field_state = RETURN_DFDX
1450 else
1451 output_field_state = RETURN_DFDX_LM
1452 end if
1453 call get_next_output_fieldname(nest_num, field_name, ndims, &
1454 min_cat, max_cat, istagger, &
1455 memorder, dimnames, units, description, &
1456 sr_x, sr_y, istatus)
1457 return
1458 end if
1459 else
1460 call mprintf(.true., ERROR, 'get_next_output_fieldname(): In state DOMCAT, '// &
1461 'but no field name is found.')
1462 end if
1464 else if (output_field_state == RETURN_DFDX .or. &
1465 output_field_state == RETURN_DFDX_LM) then
1466 if (output_field_state == RETURN_DFDX) then
1467 next_field = next_field - 1
1468 call get_next_fieldname(field_name, istatus)
1469 else
1470 call get_landmask_field(geog_data_res(nest_num), field_name, is_water_mask, landmask, istatus)
1471 end if
1472 if (istatus == 0) then
1473 call get_dfdx_name(field_name, dfdx_name, istatus)
1474 if (istatus == 0) then
1475 ndims = 2
1476 min_cat = 1
1477 max_cat = 1
1478 temp_fieldtype = iget_fieldtype(field_name, istatus)
1479 if (istatus == 0) then
1480 if (temp_fieldtype == CONTINUOUS) then
1481 call get_max_levels(field_name, min_cat, max_cat, istatus)
1482 else if (temp_fieldtype == CATEGORICAL) then
1483 call get_max_categories(field_name, min_cat, max_cat, istatus)
1484 end if
1485 if (max_cat - min_cat > 0) ndims = 3
1486 end if
1487 call get_output_stagger(field_name, istagger, istatus)
1488 if (istagger == M) then
1489 dimnames(1) = 'west_east'
1490 dimnames(2) = 'south_north'
1491 else if (istagger == U) then
1492 dimnames(1) = 'west_east_stag'
1493 dimnames(2) = 'south_north'
1494 else if (istagger == V) then
1495 dimnames(1) = 'west_east'
1496 dimnames(2) = 'south_north_stag'
1497 else if (istagger == HH) then
1498 dimnames(1) = 'west_east'
1499 dimnames(2) = 'south_north'
1500 else if (istagger == VV) then
1501 dimnames(1) = 'west_east'
1502 dimnames(2) = 'south_north'
1503 end if
1504 if (ndims == 2) then
1505 memorder = 'XY '
1506 dimnames(3) = ' '
1507 else if (ndims == 3) then
1508 memorder = 'XYZ'
1509 call get_z_dim_name(field_name, dimnames(3), istatus)
1510 istatus = 0
1511 else
1512 memorder = ' '
1513 dimnames(3) = ' '
1514 end if
1515 field_name = dfdx_name
1516 units = '-'
1518 call get_subgrid_dim_name(nest_num, field_name, dimnames(1:2), sr_x, sr_y, istatus)
1519 description = 'df/dx'
1520 if (output_field_state == RETURN_DFDX) then
1521 output_field_state = RETURN_DFDY
1522 else
1523 output_field_state = RETURN_DFDY_LM
1524 end if
1525 else
1526 if (output_field_state == RETURN_DFDX) then
1527 output_field_state = RETURN_DFDY
1528 else
1529 output_field_state = RETURN_DFDY_LM
1530 end if
1531 call get_next_output_fieldname(nest_num, field_name, ndims, &
1532 min_cat, max_cat, istagger, &
1533 memorder, dimnames, units, description, &
1534 sr_x, sr_y, istatus)
1535 return
1536 end if
1537 else
1538 call mprintf(.true., ERROR, 'get_next_output_fieldname(): In state DFDX, '// &
1539 'but no field name is found.')
1540 end if
1542 else if (output_field_state == RETURN_DFDY .or. &
1543 output_field_state == RETURN_DFDY_LM) then
1544 if (output_field_state == RETURN_DFDY) then
1545 next_field = next_field - 1
1546 call get_next_fieldname(field_name, istatus)
1547 else
1548 call get_landmask_field(geog_data_res(nest_num), field_name, is_water_mask, landmask, istatus)
1549 end if
1550 if (istatus == 0) then
1551 call get_dfdy_name(field_name, dfdy_name, istatus)
1552 if (istatus == 0) then
1553 ndims = 2
1554 min_cat = 1
1555 max_cat = 1
1556 temp_fieldtype = iget_fieldtype(field_name, istatus)
1557 if (istatus == 0) then
1558 if (temp_fieldtype == CONTINUOUS) then
1559 call get_max_levels(field_name, min_cat, max_cat, istatus)
1560 else if (temp_fieldtype == CATEGORICAL) then
1561 call get_max_categories(field_name, min_cat, max_cat, istatus)
1562 end if
1563 if (max_cat - min_cat > 0) ndims = 3
1564 end if
1565 call get_output_stagger(field_name, istagger, istatus)
1566 if (istagger == M) then
1567 dimnames(1) = 'west_east'
1568 dimnames(2) = 'south_north'
1569 else if (istagger == U) then
1570 dimnames(1) = 'west_east_stag'
1571 dimnames(2) = 'south_north'
1572 else if (istagger == V) then
1573 dimnames(1) = 'west_east'
1574 dimnames(2) = 'south_north_stag'
1575 else if (istagger == HH) then
1576 dimnames(1) = 'west_east'
1577 dimnames(2) = 'south_north'
1578 else if (istagger == VV) then
1579 dimnames(1) = 'west_east'
1580 dimnames(2) = 'south_north'
1581 end if
1582 if (ndims == 2) then
1583 memorder = 'XY '
1584 dimnames(3) = ' '
1585 else if (ndims == 3) then
1586 memorder = 'XYZ'
1587 call get_z_dim_name(field_name, dimnames(3), istatus)
1588 istatus = 0
1589 else
1590 memorder = ' '
1591 dimnames(3) = ' '
1592 end if
1594 call get_subgrid_dim_name(nest_num, field_name, dimnames(1:2), sr_x, sr_y, istatus)
1595 field_name = dfdy_name
1596 units = '-'
1597 description = 'df/dy'
1598 output_field_state = RETURN_FIELDNAME
1599 else
1600 output_field_state = RETURN_FIELDNAME
1601 call get_next_output_fieldname(nest_num, field_name, ndims, &
1602 min_cat, max_cat, istagger, &
1603 memorder, dimnames, units, description, &
1604 sr_x, sr_y, istatus)
1605 return
1606 end if
1607 else
1608 call mprintf(.true., ERROR, 'get_next_output_fieldname(): In state DFDY, but no '// &
1609 'field name is found.')
1610 end if
1612 end if
1614 end subroutine get_next_output_fieldname
1617 !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
1618 ! Name: get_landmask_field
1619 !
1620 ! Purpose: To return the name of the field from which the landmask is to be
1621 ! computed. If no error occurs, is_water_mask is .true. if the landmask
1622 ! value specifies the value of water, and .false. if the landmask value
1623 ! specifies the value of land.
1624 !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
1625 subroutine get_landmask_field(res_string, landmask_name, is_water_mask, landmask, istatus)
1627 implicit none
1629 ! Arguments
1630 character (len=128), intent(in) :: res_string
1631 integer, dimension(:), intent(out) :: landmask
1632 integer, intent(out) :: istatus
1633 logical, intent(out) :: is_water_mask
1634 character (len=128), intent(out) :: landmask_name
1636 ! Local variables
1637 integer :: j
1638 integer :: ilen
1639 integer :: idx
1640 integer :: is, ie, sos, eos, comma
1641 character (len=128) :: temp_res, mask_cat_string
1643 istatus = 1
1645 do idx=1,num_entries
1647 if (list_length(source_landmask_land(idx)) > 0) then
1648 is = 1
1649 ie = index(res_string(is:128),'+') - 1
1650 if (ie <= 0) ie = 128
1651 temp_res = res_string(is:ie)
1652 do while (.not. list_search(source_landmask_land(idx), ckey=temp_res, cvalue=mask_cat_string) &
1653 .and. is <= 128)
1654 is = ie+2
1655 ie = is + index(res_string(is:128),'+') - 2
1656 if (ie - is <= 0) ie = 128
1657 temp_res = res_string(is:ie)
1658 end do
1660 if (is > 128) then
1661 temp_res = 'default'
1662 if (list_search(source_landmask_land(idx), ckey=temp_res, cvalue=mask_cat_string)) then
1663 is_water_mask = .false.
1664 landmask_name = source_fieldname(idx)
1665 istatus = 0
1666 end if
1667 else
1668 is_water_mask = .false.
1669 landmask_name = source_fieldname(idx)
1670 istatus = 0
1671 end if
1673 end if
1675 ! Note: The following cannot be an else-if, since different resolutions of data may
1676 ! specify, alternately, a land or a water mask, and in general we need to search
1677 ! both lists
1679 if (list_length(source_landmask_water(idx)) > 0) then
1680 is = 1
1681 ie = index(res_string(is:128),'+') - 1
1682 if (ie <= 0) ie = 128
1683 temp_res = res_string(is:ie)
1684 do while (.not. list_search(source_landmask_water(idx), ckey=temp_res, cvalue=mask_cat_string) &
1685 .and. is <= 128)
1686 is = ie+2
1687 ie = is + index(res_string(is:128),'+') - 2
1688 if (ie - is <= 0) ie = 128
1689 temp_res = res_string(is:ie)
1690 end do
1692 if (is > 128) then
1693 temp_res = 'default'
1694 if (list_search(source_landmask_water(idx), ckey=temp_res, cvalue=mask_cat_string)) then
1695 is_water_mask = .true.
1696 landmask_name = source_fieldname(idx)
1697 istatus = 0
1698 end if
1699 else
1700 is_water_mask = .true.
1701 landmask_name = source_fieldname(idx)
1702 istatus = 0
1703 end if
1704 end if
1706 if (istatus == 0) then
1707 j = 1
1708 sos = 0
1709 eos = 128
1710 comma = index(mask_cat_string(sos+1:eos-1),',')
1711 do while (comma > 0 .and. j < MAX_LANDMASK_CATEGORIES)
1712 read(mask_cat_string(sos+1:sos+comma-1),'(i10)') landmask(j)
1713 sos = sos + comma
1714 comma = index(mask_cat_string(sos+1:eos-1),',')
1715 j = j + 1
1716 end do
1717 read(mask_cat_string(sos+1:eos-1),'(i10)') landmask(j)
1718 j = j + 1
1719 if (j <= MAX_LANDMASK_CATEGORIES) then ! Terminate list with a flag value
1720 landmask(j) = INVALID
1721 end if
1722 exit
1723 end if
1725 end do
1727 end subroutine get_landmask_field
1730 !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
1731 ! Name: get_missing_value
1732 !
1733 ! Pupose: Return the value used in the source data to indicate missing data.
1734 !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
1735 subroutine get_missing_value(fieldnm, ilevel, rmissing, istatus)
1737 implicit none
1739 ! Arguments
1740 integer, intent(in) :: ilevel
1741 integer, intent(out) :: istatus
1742 real, intent(out) :: rmissing
1743 character (len=128), intent(in) :: fieldnm
1745 ! Local variables
1746 integer :: idx
1748 istatus = 1
1750 do idx=1,num_entries
1751 if ((index(source_fieldname(idx),trim(fieldnm)) /= 0) .and. &
1752 (len_trim(source_fieldname(idx)) == len_trim(fieldnm)) .and. &
1753 (source_priority(idx) == ilevel)) then
1755 if (is_missing_value(idx)) then
1756 rmissing = source_missing_value(idx)
1757 istatus = 0
1758 exit
1759 end if
1761 end if
1762 end do
1764 end subroutine get_missing_value
1767 !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
1768 ! Name: get_source_units
1769 !
1770 ! Pupose: Return a string giving the units of the specified source data.
1771 !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
1772 subroutine get_source_units(fieldnm, ilevel, cunits, istatus)
1774 implicit none
1776 ! Arguments
1777 integer, intent(in) :: ilevel
1778 integer, intent(out) :: istatus
1779 character (len=128), intent(in) :: fieldnm
1780 character (len=128), intent(out) :: cunits
1782 ! Local variables
1783 integer :: idx
1785 istatus = 1
1787 do idx=1,num_entries
1788 if ((index(source_fieldname(idx),trim(fieldnm)) /= 0) .and. &
1789 (len_trim(source_fieldname(idx)) == len_trim(fieldnm)) .and. &
1790 (source_priority(idx) == ilevel)) then
1792 if (is_units(idx)) then
1793 cunits = source_units(idx)
1794 istatus = 0
1795 exit
1796 end if
1798 end if
1799 end do
1801 end subroutine get_source_units
1804 !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
1805 ! Name: get_source_descr
1806 !
1807 ! Pupose: Return a string giving a description of the specified source data.
1808 !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
1809 subroutine get_source_descr(fieldnm, ilevel, descr, istatus)
1811 implicit none
1813 ! Arguments
1814 integer, intent(in) :: ilevel
1815 integer, intent(out) :: istatus
1816 character (len=128), intent(in) :: fieldnm
1817 character (len=128), intent(out) :: descr
1819 ! Local variables
1820 integer :: idx
1822 istatus = 1
1824 do idx=1,num_entries
1825 if ((index(source_fieldname(idx),trim(fieldnm)) /= 0) .and. &
1826 (len_trim(source_fieldname(idx)) == len_trim(fieldnm)) .and. &
1827 (source_priority(idx) == ilevel)) then
1829 if (is_units(idx)) then
1830 descr = source_descr(idx)
1831 istatus = 0
1832 exit
1833 end if
1835 end if
1836 end do
1838 end subroutine get_source_descr
1841 !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
1842 ! Name: get_missing_fill_value
1843 !
1844 ! Pupose: Return the value to fill missing points with.
1845 !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
1846 subroutine get_missing_fill_value(fieldnm, rmissing, istatus)
1848 implicit none
1850 ! Arguments
1851 integer, intent(out) :: istatus
1852 real, intent(out) :: rmissing
1853 character (len=128), intent(in) :: fieldnm
1855 ! Local variables
1856 integer :: idx
1858 istatus = 1
1860 do idx=1,num_entries
1861 if ((index(source_fieldname(idx),trim(fieldnm)) /= 0) .and. &
1862 (len_trim(source_fieldname(idx)) == len_trim(fieldnm)) ) then
1864 if (is_fill_missing(idx)) then
1865 rmissing = source_fill_missing(idx)
1866 istatus = 0
1867 exit
1868 end if
1870 end if
1871 end do
1873 end subroutine get_missing_fill_value
1876 !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
1877 ! Name: get_halt_on_missing
1878 !
1879 ! Pupose: Returns 1 if the program should halt upon encountering a missing
1880 ! value in the final output field, and 0 otherwise.
1881 !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
1882 subroutine get_halt_on_missing(fieldnm, halt, istatus)
1884 implicit none
1886 ! Arguments
1887 integer, intent(out) :: istatus
1888 logical, intent(out) :: halt
1889 character (len=128), intent(in) :: fieldnm
1891 ! Local variables
1892 integer :: idx
1894 istatus = 0
1895 halt = .false.
1897 do idx=1,num_entries
1898 if ((index(source_fieldname(idx),trim(fieldnm)) /= 0) .and. &
1899 (len_trim(source_fieldname(idx)) == len_trim(fieldnm)) ) then
1901 if (is_halt_missing(idx)) halt = .true.
1903 end if
1904 end do
1906 end subroutine get_halt_on_missing
1909 !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
1910 ! Name: get_masked_value
1911 !
1912 ! Pupose: If the field is to be masked by the landmask, returns 0 if the field
1913 ! is masked over water and 1 if the field is masked over land. If no mask is
1914 ! to be applied, -1 is returned.
1915 !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
1916 subroutine get_masked_value(fieldnm, ilevel, masked, istatus)
1918 implicit none
1920 ! Arguments
1921 integer, intent(in) :: ilevel
1922 integer, intent(out) :: istatus
1923 real, intent(out) :: masked
1924 character (len=128), intent(in) :: fieldnm
1926 ! Local variables
1927 integer :: idx
1929 istatus = 0
1930 masked = -1.
1932 do idx=1,num_entries
1933 if ((index(source_fieldname(idx),trim(fieldnm)) /= 0) .and. &
1934 (len_trim(source_fieldname(idx)) == len_trim(fieldnm)) .and. &
1935 (source_priority(idx) == ilevel)) then
1937 if (is_masked(idx)) then
1938 masked = source_masked(idx)
1939 exit
1940 end if
1942 end if
1943 end do
1945 end subroutine get_masked_value
1948 !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
1949 ! Name: get_max_levels
1950 !
1951 ! Purpose: Returns the number of levels for the field given by fieldnm.
1952 ! The number of levels will generally be specified for continuous fields,
1953 ! whereas min/max category will generally be specified for categorical ones.
1954 !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
1955 subroutine get_max_levels(fieldnm, min_level, max_level, istatus)
1957 implicit none
1959 ! Arguments
1960 integer, intent(out) :: min_level, max_level, istatus
1961 character (len=128), intent(in) :: fieldnm
1963 ! Local variables
1964 integer :: idx
1965 logical :: have_found_field
1967 have_found_field = .false.
1968 istatus = 1
1970 do idx=1,num_entries
1971 if ((index(source_fieldname(idx),trim(fieldnm)) /= 0) .and. &
1972 (len_trim(source_fieldname(idx)) == len_trim(fieldnm))) then
1974 if (is_dest_fieldtype(idx) .and. (source_dest_fieldtype(idx) /= CONTINUOUS)) then
1975 call mprintf(.true., WARN, 'In GEOGRID.TBL, destination field type for %s is '// &
1976 'not continuous and min/max levels specified.', s1=trim(fieldnm))
1977 end if
1978 if (.not. have_found_field) then
1979 if (is_tile_z_start(idx) .and. is_tile_z_end(idx)) then
1980 have_found_field = .true.
1981 istatus = 0
1982 min_level = source_tile_z_start(idx)
1983 max_level = source_tile_z_end(idx)
1984 else if (is_tile_z(idx)) then
1985 have_found_field = .true.
1986 istatus = 0
1987 min_level = 1
1988 max_level = source_tile_z(idx)
1989 end if
1991 if (.not. (is_tile_z_start(idx) .and. is_tile_z_end(idx))) then
1992 if (is_tile_z_start(idx) .or. is_tile_z_end(idx)) then
1993 call mprintf(.true., ERROR, 'In GEOGRID.TBL, only one of tile_z_start '// &
1994 'and tile_z_end specified for entry %i.',i1=idx)
1995 end if
1996 end if
1997 else
1998 if (is_tile_z_start(idx) .and. is_tile_z_end(idx)) then
1999 if (source_tile_z_start(idx) < min_level) min_level = source_tile_z_start(idx)
2000 if (source_tile_z_end(idx) > max_level) max_level = source_tile_z_end(idx)
2001 else if (is_tile_z(idx)) then
2002 if (min_level > 1) min_level = 1
2003 if (source_tile_z(idx) > max_level) max_level = source_tile_z(idx)
2004 end if
2006 if (.not. (is_tile_z_start(idx) .and. is_tile_z_end(idx))) then
2007 if (is_tile_z_start(idx) .or. is_tile_z_end(idx)) then
2008 call mprintf(.true., ERROR, 'In GEOGRID.TBL, only one of tile_z_start '// &
2009 'and tile_z_end specified for entry %i.',i1=idx)
2010 end if
2011 end if
2012 end if
2014 end if
2015 end do
2017 end subroutine get_max_levels
2020 !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
2021 ! Name: get_source_levels
2022 !
2023 ! Purpose: Return the min and max z-index for the source data for fieldname
2024 ! at a specified priority level (compared with the min/max level over
2025 ! all priority levels, as given by get_max_levels).
2026 !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
2027 subroutine get_source_levels(fieldnm, ilevel, min_level, max_level, istatus)
2029 implicit none
2031 ! Arguments
2032 integer, intent(in) :: ilevel
2033 integer, intent(out) :: min_level, max_level, istatus
2034 character (len=128), intent(in) :: fieldnm
2036 ! Local variables
2037 integer :: idx
2039 istatus = 1
2041 do idx=1,num_entries
2042 if ((index(source_fieldname(idx),trim(fieldnm)) /= 0) .and. &
2043 (len_trim(source_fieldname(idx)) == len_trim(fieldnm))) then
2044 if (ilevel == source_priority(idx)) then
2046 if (is_tile_z_start(idx) .and. is_tile_z_end(idx)) then
2047 istatus = 0
2048 min_level = source_tile_z_start(idx)
2049 max_level = source_tile_z_end(idx)
2050 else if (is_tile_z(idx)) then
2051 istatus = 0
2052 min_level = 1
2053 max_level = source_tile_z(idx)
2054 end if
2056 if (.not. (is_tile_z_start(idx) .and. is_tile_z_end(idx))) then
2057 if (is_tile_z_start(idx) .or. is_tile_z_end(idx)) then
2058 call mprintf(.true., ERROR, 'In GEOGRID.TBL, only one of tile_z_start '// &
2059 'and tile_z_end specified for entry %i.',i1=idx)
2060 end if
2061 end if
2063 end if
2064 end if
2065 end do
2067 end subroutine get_source_levels
2070 !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
2071 ! Name: get_max_categories
2072 !
2073 ! Purpose: Returns the minimum category and the maximum category for the field
2074 ! given by fieldnm.
2075 ! Min/max category will generally be specified for categorical fields,
2076 ! whereas the number of levels will generally be specified for continuous
2077 ! fields.
2078 !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
2079 subroutine get_max_categories(fieldnm, min_category, max_category, istatus)
2081 implicit none
2083 ! Arguments
2084 integer, intent(out) :: min_category, max_category, istatus
2085 character (len=128), intent(in) :: fieldnm
2087 ! Local variables
2088 integer :: idx
2089 logical :: have_found_field
2091 have_found_field = .false.
2092 istatus = 1
2094 do idx=1,num_entries
2095 if ((index(source_fieldname(idx),trim(fieldnm)) /= 0) .and. &
2096 (len_trim(source_fieldname(idx)) == len_trim(fieldnm))) then
2098 if (is_dest_fieldtype(idx) .and. (source_dest_fieldtype(idx) /= CATEGORICAL)) then
2099 call mprintf(.true., WARN, &
2100 'In GEOGRID.TBL, cannot get min/max categories for continuous '// &
2101 'field %s at entry %i. Perhaps the user has requested to '// &
2102 'perform a strange operation on the field.', s1=trim(fieldnm), i1=idx)
2103 end if
2104 if (.not. have_found_field) then
2105 if (is_category_min(idx) .and. is_category_max(idx)) then
2106 have_found_field = .true.
2107 istatus = 0
2108 min_category = source_category_min(idx)
2109 max_category = source_category_max(idx)
2110 else if (is_category_min(idx) .or. is_category_max(idx)) then
2111 call mprintf(.true., ERROR, 'In GEOGRID.TBL, only one of min_category and '// &
2112 'max_category specified for entry %i.',i1=idx)
2113 end if
2114 else
2115 if (is_category_min(idx) .and. is_category_max(idx)) then
2116 if (source_category_min(idx) < min_category) min_category = source_category_min(idx)
2117 if (source_category_max(idx) > max_category) max_category = source_category_max(idx)
2118 else if (is_category_min(idx) .or. is_category_max(idx)) then
2119 call mprintf(.true., ERROR, 'In GEOGRID.TBL, only one of min_category and '// &
2120 'max_category specified for entry %i.',i1=idx)
2121 end if
2122 end if
2124 end if
2125 end do
2127 end subroutine get_max_categories
2130 !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
2131 ! Name: get_source_categories
2132 !
2133 ! Purpose: Return the min and max category for the source data for fieldname
2134 ! at a specified priority level (compared with the min/max category over
2135 ! all priority levels, as given by get_max_categories).
2136 !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
2137 subroutine get_source_categories(fieldnm, ilevel, min_category, max_category, istatus)
2139 implicit none
2141 ! Arguments
2142 integer, intent(in) :: ilevel
2143 integer, intent(out) :: min_category, max_category, istatus
2144 character (len=128), intent(in) :: fieldnm
2146 ! Local variables
2147 integer :: idx
2149 istatus = 1
2151 do idx=1,num_entries
2152 if ((index(source_fieldname(idx),trim(fieldnm)) /= 0) .and. &
2153 (len_trim(source_fieldname(idx)) == len_trim(fieldnm))) then
2154 if (ilevel == source_priority(idx)) then
2156 if (is_category_min(idx) .and. is_category_max(idx)) then
2157 istatus = 0
2158 min_category = source_category_min(idx)
2159 max_category = source_category_max(idx)
2160 else if (is_category_min(idx) .or. is_category_max(idx)) then
2161 call mprintf(.true., ERROR, 'In GEOGRID.TBL, only one of min_category '// &
2162 'and max_category specified for entry %i.',i1=idx)
2163 end if
2165 end if
2166 end if
2167 end do
2169 end subroutine get_source_categories
2172 !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
2173 ! Name: get_domcategory_name
2174 !
2175 ! Purpose:
2176 !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
2177 subroutine get_domcategory_name(fieldnm, domcat_name, ldominant_only, istatus)
2179 implicit none
2181 ! Arguments
2182 integer, intent(out) :: istatus
2183 logical, intent(out) :: ldominant_only
2184 character (len=128), intent(in) :: fieldnm
2185 character (len=128), intent(out) :: domcat_name
2187 ! Local variables
2188 integer :: idx
2190 istatus = 1
2191 ldominant_only = .false.
2193 do idx=1,num_entries
2194 if ((index(source_fieldname(idx),trim(fieldnm)) /= 0) .and. &
2195 (len_trim(source_fieldname(idx)) == len_trim(fieldnm))) then
2197 if (is_dominant_category(idx)) then
2198 domcat_name = source_dominant_category(idx)
2199 istatus = 0
2200 if (is_dominant_only(idx)) then
2201 call mprintf(.true., WARN, 'In GEOGRID.TBL, both dominant_category and '// &
2202 'dominant_only are specified in entry %i. Using specification '// &
2203 'for dominant_category.',i1=idx)
2204 is_dominant_only(idx) = .false.
2205 end if
2206 exit
2208 else if (is_dominant_only(idx)) then
2209 domcat_name = source_dominant_only(idx)
2210 ldominant_only = .true.
2211 istatus = 0
2212 exit
2213 end if
2215 end if
2216 end do
2218 end subroutine get_domcategory_name
2221 !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
2222 ! Name: get_dfdx_name
2223 !
2224 ! Purpose:
2225 !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
2226 subroutine get_dfdx_name(fieldnm, dfdx_name, istatus)
2228 implicit none
2230 ! Arguments
2231 integer, intent(out) :: istatus
2232 character (len=128), intent(in) :: fieldnm
2233 character (len=128), intent(out) :: dfdx_name
2235 ! Local variables
2236 integer :: idx
2238 istatus = 1
2240 do idx=1,num_entries
2241 if ((index(source_fieldname(idx),trim(fieldnm)) /= 0) .and. &
2242 (len_trim(source_fieldname(idx)) == len_trim(fieldnm))) then
2244 if (is_dfdx(idx)) then
2245 dfdx_name = source_dfdx(idx)
2246 istatus = 0
2247 exit
2248 end if
2250 end if
2251 end do
2253 end subroutine get_dfdx_name
2256 !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
2257 ! Name: get_dfdy_name
2258 !
2259 ! Purpose:
2260 !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
2261 subroutine get_dfdy_name(fieldnm, dfdy_name, istatus)
2263 implicit none
2265 ! Arguments
2266 integer, intent(out) :: istatus
2267 character (len=128), intent(in) :: fieldnm
2268 character (len=128), intent(out) :: dfdy_name
2270 ! Local variables
2271 integer :: idx
2273 istatus = 1
2275 do idx=1,num_entries
2276 if ((index(source_fieldname(idx),trim(fieldnm)) /= 0) .and. &
2277 (len_trim(source_fieldname(idx)) == len_trim(fieldnm))) then
2279 if (is_dfdy(idx)) then
2280 dfdy_name = source_dfdy(idx)
2281 istatus = 0
2282 exit
2283 end if
2285 end if
2286 end do
2288 end subroutine get_dfdy_name
2291 !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
2292 ! Name: get_z_dim_name
2293 !
2294 ! Purpose:
2295 !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
2296 subroutine get_z_dim_name(fieldnm, z_dim, istatus)
2298 implicit none
2300 ! Arguments
2301 integer, intent(out) :: istatus
2302 character (len=128), intent(in) :: fieldnm
2303 character (len=128), intent(out) :: z_dim
2305 ! Local variables
2306 integer :: idx
2308 istatus = 1
2310 do idx=1,num_entries
2311 if ((index(source_fieldname(idx),trim(fieldnm)) /= 0) .and. &
2312 (len_trim(source_fieldname(idx)) == len_trim(fieldnm))) then
2313 if (is_z_dim_name(idx)) then
2314 z_dim = source_z_dim_name(idx)
2315 istatus = 0
2316 exit
2317 end if
2318 end if
2319 end do
2321 end subroutine get_z_dim_name
2324 !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
2325 ! Name: get_field_scale_factor
2326 !
2327 ! Purpose:
2328 !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
2329 subroutine get_field_scale_factor(fieldnm, ilevel, scale_factor, istatus)
2331 implicit none
2333 ! Arguments
2334 integer, intent(in) :: ilevel
2335 integer, intent(out) :: istatus
2336 real, intent(out) :: scale_factor
2337 character (len=128), intent(in) :: fieldnm
2339 ! Local variables
2340 integer :: idx
2342 istatus = 1
2344 do idx=1,num_entries
2345 if ((index(source_fieldname(idx),trim(fieldnm)) /= 0) .and. &
2346 (len_trim(source_fieldname(idx)) == len_trim(fieldnm)) .and. &
2347 (ilevel == source_priority(idx))) then
2349 if (is_scale_factor(idx)) then
2350 scale_factor = source_scale_factor(idx)
2351 istatus = 0
2352 end if
2354 end if
2355 end do
2357 end subroutine get_field_scale_factor
2360 !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
2361 ! Name: get_output_stagger
2362 !
2363 ! Pupose:
2364 !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
2365 subroutine get_output_stagger(fieldnm, istagger, istatus)
2367 use gridinfo_module
2369 implicit none
2371 ! Arguments
2372 integer, intent(out) :: istatus, istagger
2373 character (len=128), intent(in) :: fieldnm
2375 ! Local variables
2376 integer :: idx
2378 istatus = 1
2379 do idx=1,num_entries
2380 if ((index(source_fieldname(idx),trim(fieldnm)) /= 0) .and. &
2381 (len_trim(source_fieldname(idx)) == len_trim(fieldnm))) then
2383 if (is_output_stagger(idx)) then
2384 istatus = 0
2385 istagger = source_output_stagger(idx)
2386 exit
2387 else
2388 if (gridtype == 'C') then
2389 istatus = 0
2390 istagger = M
2391 exit
2392 else if (gridtype == 'E') then
2393 istatus = 0
2394 istagger = HH
2395 exit
2396 end if
2397 end if
2399 end if
2400 end do
2402 end subroutine get_output_stagger
2405 !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
2406 ! Name: get_subgrid_dim_name
2407 !
2408 ! Pupose:
2409 !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
2410 subroutine get_subgrid_dim_name(nest_num, field_name, dimnames, &
2411 sub_x, sub_y, istatus)
2413 use gridinfo_module
2415 implicit none
2416 integer, intent(in) :: nest_num
2417 integer, intent(out) :: sub_x, sub_y, istatus
2418 character(len=128), intent(in) :: field_name
2419 character(len=128), dimension(2), intent(inout) :: dimnames
2420 integer :: idx, nlen
2422 sub_x = 1
2423 sub_y = 1
2425 istatus = 0
2426 do idx=1,num_entries
2427 if ((index(source_fieldname(idx),trim(field_name)) /= 0) .and. &
2428 (len_trim(source_fieldname(idx)) == len_trim(field_name))) then
2429 if (is_subgrid(idx)) then
2430 istatus = 0
2431 if (is_output_stagger(idx)) then
2432 call mprintf(.true.,ERROR,'Cannot use subgrids on variables with staggered grids')
2433 end if
2434 dimnames(1) = trim(dimnames(1))//"_subgrid"
2435 dimnames(2) = trim(dimnames(2))//"_subgrid"
2436 sub_x = subgrid_ratio_x(nest_num)
2437 sub_y = subgrid_ratio_y(nest_num)
2438 end if
2439 end if
2440 end do
2442 end subroutine get_subgrid_dim_name
2445 !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
2446 ! Name: get_interp_option
2447 !
2448 ! Pupose:
2449 !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
2450 subroutine get_interp_option(fieldnm, ilevel, interp_opt, istatus)
2452 implicit none
2454 ! Arguments
2455 integer, intent(in) :: ilevel
2456 integer, intent(out) :: istatus
2457 character (len=128), intent(in) :: fieldnm
2458 character (len=128), intent(out) :: interp_opt
2460 ! Local variables
2461 integer :: idx
2463 istatus = 1
2464 do idx=1,num_entries
2465 if ((index(source_fieldname(idx),trim(fieldnm)) /= 0) .and. &
2466 (len_trim(source_fieldname(idx)) == len_trim(fieldnm))) then
2467 if (ilevel == source_priority(idx)) then
2469 interp_opt = source_interp_string(idx)
2470 istatus = 0
2471 exit
2473 end if
2474 end if
2475 end do
2477 end subroutine get_interp_option
2480 !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
2481 ! Name: get_gcel_threshold
2482 !
2483 ! Pupose:
2484 !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
2485 subroutine get_gcell_threshold(interp_opt, threshold, istatus)
2487 implicit none
2489 ! Arguments
2490 integer, intent(out) :: istatus
2491 real, intent(out) :: threshold
2492 character (len=128), intent(in) :: interp_opt
2494 ! Local variables
2495 integer :: i, p1, p2
2497 istatus = 1
2498 threshold = 1.0
2500 i = index(interp_opt,'average_gcell')
2501 if (i /= 0) then
2502 ! Move the "average_gcell" option to the beginning
2503 ! if (i /= 1) then
2504 ! p1 =
2505 ! end if
2507 ! Check for a threshold
2508 p1 = index(interp_opt(i:128),'(')
2509 p2 = index(interp_opt(i:128),')')
2510 if (p1 /= 0 .and. p2 /= 0) then
2511 read(interp_opt(p1+1:p2-1),*,err=1000) threshold
2512 else
2513 call mprintf(.true., WARN, 'Problem with specified threshold '// &
2514 'for average_gcell interp option. Setting threshold to 0.0.')
2515 threshold = 0.0
2516 end if
2517 end if
2518 istatus = 0
2520 return
2522 1000 call mprintf(.true., ERROR, 'Threshold option to average_gcell interpolator '// &
2523 'must be a real number, enclosed in parentheses immediately '// &
2524 'after keyword "average_gcell"')
2526 end subroutine get_gcell_threshold
2529 !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
2530 ! Name: get_smooth_option
2531 !
2532 ! Pupose:
2533 !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
2534 subroutine get_smooth_option(fieldnm, smooth_opt, smooth_passes, istatus)
2536 implicit none
2538 ! Arguments
2539 integer, intent(out) :: istatus, smooth_opt, smooth_passes
2540 character (len=128), intent(in) :: fieldnm
2542 ! Local variables
2543 integer :: idx
2545 istatus = 1
2547 do idx=1,num_entries
2548 if ((index(source_fieldname(idx),trim(fieldnm)) /= 0) .and. &
2549 (len_trim(source_fieldname(idx)) == len_trim(fieldnm))) then
2551 if (is_smooth_option(idx)) then
2552 istatus = 0
2553 smooth_opt = source_smooth_option(idx)
2555 if (is_smooth_passes(idx)) then
2556 smooth_passes = source_smooth_passes(idx)
2557 else
2558 smooth_passes = 1
2559 end if
2561 exit
2562 end if
2564 end if
2565 end do
2567 end subroutine get_smooth_option
2570 !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
2571 ! Name: iget_fieldtype
2572 !
2573 ! Purpose:
2574 !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
2575 function iget_fieldtype(fieldnm, istatus)
2577 implicit none
2579 ! Arguments
2580 integer, intent(out) :: istatus
2581 character (len=128), intent(in) :: fieldnm
2583 ! Local variables
2584 integer :: idx
2586 ! Return value
2587 integer :: iget_fieldtype
2589 istatus = 1
2591 do idx=1,num_entries
2592 if ((index(source_fieldname(idx),trim(fieldnm)) /= 0) .and. &
2593 (len_trim(source_fieldname(idx)) == len_trim(fieldnm))) then
2595 if (is_dest_fieldtype(idx)) then
2596 iget_fieldtype = source_dest_fieldtype(idx)
2597 istatus = 0
2598 exit
2599 end if
2601 end if
2602 end do
2604 end function iget_fieldtype
2607 !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
2608 ! Name: iget_source_fieldtype
2609 !
2610 ! Purpose: Given a resolution, in degrees, and the name of a field, this
2611 ! function returns the type (categorical, continuous, etc.) of the source
2612 ! data that will be used. This may, in general, depend on the field name
2613 ! and the resolution; for example, near 30 second resolution, land use data
2614 ! may come from a categorical field, whereas for lower resolutions, it may
2615 ! come from arrays of land use fractions for each category.
2616 !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
2617 function iget_source_fieldtype(fieldnm, ilevel, istatus)
2619 implicit none
2621 ! Arguments
2622 integer, intent(in) :: ilevel
2623 integer, intent(out) :: istatus
2624 character (len=128), intent(in) :: fieldnm
2626 ! Return value
2627 integer :: iget_source_fieldtype
2629 ! Local variables
2630 integer :: idx
2632 ! Find information about the source tiles for the specified fieldname
2633 istatus = 1
2634 do idx=1,num_entries
2635 if ((index(source_fieldname(idx),trim(fieldnm)) /= 0) .and. &
2636 (len_trim(source_fieldname(idx)) == len_trim(fieldnm))) then
2637 if (ilevel == source_priority(idx)) then
2638 istatus = 0
2639 iget_source_fieldtype = source_fieldtype(idx)
2640 exit
2641 end if
2642 end if
2643 end do
2645 end function iget_source_fieldtype
2648 !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
2649 ! Name: get_data_tile
2650 !
2651 ! Purpose:
2652 !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
2653 subroutine get_data_tile(xlat, xlon, ilevel, field_name, &
2654 file_name, array, start_x_dim, end_x_dim, start_y_dim, &
2655 end_y_dim, start_z_dim, end_z_dim, npts_bdr, &
2656 istatus)
2658 implicit none
2660 ! Arguments
2661 integer, intent(in) :: ilevel
2662 integer, intent(out) :: istatus
2663 integer, intent(out) :: start_x_dim, end_x_dim, &
2664 start_y_dim, end_y_dim, &
2665 start_z_dim, end_z_dim, &
2666 npts_bdr
2667 real, intent(in) :: xlat, xlon ! Location that tile should contain
2668 real, pointer, dimension(:,:,:) :: array ! The array to be allocated by this routine
2669 character (len=128), intent(in) :: field_name
2670 character (len=256), intent(out) :: file_name
2672 ! Local variables
2673 integer :: j, k
2674 integer :: local_wordsize, local_endian, sign_convention, irow_order, strlen
2675 integer :: xdim,ydim,zdim
2676 real :: scalefac
2677 real, allocatable, dimension(:) :: temprow
2679 call get_tile_fname(file_name, xlat, xlon, ilevel, field_name, istatus)
2681 if (index(file_name, 'OUTSIDE') /= 0) then
2682 istatus = 1
2683 return
2684 else if (hash_search(bad_files, file_name)) then
2685 istatus = 1
2686 return
2687 end if
2689 call get_tile_dimensions(xlat, xlon, start_x_dim, end_x_dim, start_y_dim, end_y_dim, &
2690 start_z_dim, end_z_dim, npts_bdr, local_wordsize, local_endian, &
2691 sign_convention, ilevel, field_name, istatus)
2693 xdim = (end_x_dim-start_x_dim+1)
2694 ydim = (end_y_dim-start_y_dim+1)
2695 zdim = (end_z_dim-start_z_dim+1)
2697 if (associated(array)) deallocate(array)
2698 allocate(array(xdim,ydim,zdim))
2700 call get_row_order(field_name, ilevel, irow_order, istatus)
2701 if (istatus /= 0) irow_order = BOTTOM_TOP
2703 call s_len(file_name,strlen)
2705 scalefac = 1.0
2707 call read_geogrid(file_name, strlen, array, xdim, ydim, zdim, &
2708 sign_convention, local_endian, scalefac, local_wordsize, istatus)
2710 if (irow_order == TOP_BOTTOM) then
2711 allocate(temprow(xdim))
2712 do k=1,zdim
2713 do j=1,ydim
2714 if (ydim-j+1 <= j) exit
2715 temprow(1:xdim) = array(1:xdim,j,k)
2716 array(1:xdim,j,k) = array(1:xdim,ydim-j+1,k)
2717 array(1:xdim,ydim-j+1,k) = temprow(1:xdim)
2718 end do
2719 end do
2720 deallocate(temprow)
2721 end if
2723 if (istatus /= 0) then
2724 start_x_dim = INVALID
2725 start_y_dim = INVALID
2726 end_x_dim = INVALID
2727 end_y_dim = INVALID
2729 call hash_insert(bad_files, file_name)
2730 end if
2732 end subroutine get_data_tile
2735 !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
2736 ! Name: get_row_order
2737 !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
2738 subroutine get_row_order(fieldnm, ilevel, irow_order, istatus)
2740 implicit none
2742 ! Arguments
2743 integer, intent(in) :: ilevel
2744 character (len=128), intent(in) :: fieldnm
2745 integer, intent(out) :: irow_order, istatus
2747 ! Local variables
2748 integer :: idx
2750 istatus = 1
2751 do idx=1,num_entries
2752 if ((index(source_fieldname(idx),trim(fieldnm)) /= 0) .and. &
2753 (len_trim(source_fieldname(idx)) == len_trim(fieldnm))) then
2754 if (ilevel == source_priority(idx)) then
2755 if (is_row_order(idx)) then
2756 irow_order = source_row_order(idx)
2757 istatus = 0
2758 exit
2759 end if
2760 end if
2761 end if
2762 end do
2764 end subroutine get_row_order
2767 !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
2768 ! Name: get_tile_dimensions
2769 !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
2770 subroutine get_tile_dimensions(xlat, xlon, start_x_dim, end_x_dim, start_y_dim, end_y_dim, &
2771 start_z_dim, end_z_dim, npts_bdr, bytes_per_datum, endianness, &
2772 sign_convention, ilevel, fieldnm, istatus)
2774 use llxy_module
2776 implicit none
2778 ! Arguments
2779 integer, intent(in) :: ilevel
2780 integer, intent(out) :: start_x_dim, end_x_dim, &
2781 start_y_dim, end_y_dim, &
2782 start_z_dim, end_z_dim, &
2783 npts_bdr, &
2784 bytes_per_datum, endianness, &
2785 sign_convention, istatus
2786 real, intent(in) :: xlat, xlon
2787 character (len=128), intent(in) :: fieldnm
2789 ! Local variables
2790 integer :: idx, current_domain
2791 real :: rx, ry
2793 istatus = 1
2794 do idx=1,num_entries
2795 if ((index(source_fieldname(idx),trim(fieldnm)) /= 0) .and. &
2796 (len_trim(source_fieldname(idx)) == len_trim(fieldnm))) then
2797 if (ilevel == source_priority(idx)) then
2798 istatus = 0
2799 exit
2800 end if
2801 end if
2802 end do
2804 if (istatus /= 0) then
2805 start_x_dim = 1
2806 start_y_dim = 1
2807 start_z_dim = 1
2808 end_x_dim = 1
2809 end_y_dim = 1
2810 end_z_dim = 1
2811 bytes_per_datum = 0
2812 return
2813 end if
2815 current_domain = iget_selected_domain()
2816 call select_domain(SOURCE_PROJ)
2817 call lltoxy(xlat, xlon, rx, ry, M)
2818 call select_domain(current_domain)
2820 start_x_dim = source_tile_x(idx) * nint(real(floor((rx-0.5) / real(source_tile_x(idx))))) + 1
2821 end_x_dim = start_x_dim + source_tile_x(idx) - 1
2823 start_y_dim = source_tile_y(idx) * nint(real(floor((ry-0.5) / real(source_tile_y(idx))))) + 1
2824 end_y_dim = start_y_dim + source_tile_y(idx) - 1
2826 if (is_tile_z_start(idx) .and. is_tile_z_end(idx)) then
2827 start_z_dim = source_tile_z_start(idx)
2828 end_z_dim = source_tile_z_end(idx)
2829 else if (is_tile_z(idx)) then
2830 start_z_dim = 1
2831 end_z_dim = source_tile_z(idx)
2832 end if
2834 if (.not. (is_tile_z_start(idx) .and. is_tile_z_end(idx))) then
2835 if (is_tile_z_start(idx) .or. is_tile_z_end(idx)) then
2836 call mprintf(.true., ERROR, 'In GEOGRID.TBL, only one of tile_z_start and '// &
2837 'tile_z_end specified for entry %i.',i1=idx)
2838 end if
2839 end if
2841 if (is_tile_bdr(idx)) then
2842 npts_bdr = source_tile_bdr(idx)
2843 else
2844 npts_bdr = 0
2845 end if
2847 start_x_dim = start_x_dim - npts_bdr
2848 end_x_dim = end_x_dim + npts_bdr
2849 start_y_dim = start_y_dim - npts_bdr
2850 end_y_dim = end_y_dim + npts_bdr
2852 if (is_wordsize(idx)) then
2853 bytes_per_datum = source_wordsize(idx)
2854 else
2855 bytes_per_datum = 1
2856 call mprintf(.true.,ERROR,'In GEOGRID.TBL, no wordsize specified for data in entry %i.',i1=idx)
2857 end if
2859 if (is_endian(idx)) then
2860 endianness = source_endian(idx)
2861 else
2862 endianness = BIG_ENDIAN
2863 end if
2865 if (is_signed(idx)) then
2866 sign_convention = 1
2867 else
2868 sign_convention = 0
2869 end if
2871 end subroutine get_tile_dimensions
2874 !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
2875 ! Name: get_tile_fname
2876 !
2877 ! Purpose:
2878 !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
2879 subroutine get_tile_fname(test_fname, xlat, xlon, ilevel, fieldname, istatus)
2881 use llxy_module
2882 use gridinfo_module
2884 implicit none
2886 ! Arguments
2887 integer, intent(in) :: ilevel
2888 integer, intent(out) :: istatus
2889 real, intent(in) :: xlat, xlon
2890 character (len=*), intent(in) :: fieldname
2891 character (len=256), intent(out) :: test_fname
2893 ! Local variables
2894 integer :: current_domain, idx
2895 real :: rx, ry
2897 istatus = 1
2898 write(test_fname, '(a)') 'TILE.OUTSIDE.DOMAIN'
2900 do idx=1,num_entries
2901 if ((index(source_fieldname(idx),trim(fieldname)) /= 0) .and. &
2902 (len_trim(source_fieldname(idx)) == len_trim(fieldname))) then
2903 if (ilevel == source_priority(idx)) then
2904 istatus = 0
2905 exit
2906 end if
2907 end if
2908 end do
2910 if (istatus /= 0) return
2912 current_domain = iget_selected_domain()
2913 call select_domain(SOURCE_PROJ)
2914 call lltoxy(xlat, xlon, rx, ry, M)
2915 call select_domain(current_domain)
2917 ! rx = real(source_tile_x(idx)) * real(floor((rx-0.5*source_dx(idx))/ real(source_tile_x(idx)))) + 1.0
2918 ! ry = real(source_tile_y(idx)) * real(floor((ry-0.5*source_dy(idx))/ real(source_tile_y(idx)))) + 1.0
2919 rx = real(source_tile_x(idx)) * real(floor((rx-0.5) / real(source_tile_x(idx)))) + 1.0
2920 ry = real(source_tile_y(idx)) * real(floor((ry-0.5) / real(source_tile_y(idx)))) + 1.0
2922 if (rx > 0. .and. ry > 0.) then
2923 write(test_fname, '(a,i5.5,a1,i5.5,a1,i5.5,a1,i5.5)') trim(source_path(idx)), &
2924 nint(rx),'-',nint(rx)+source_tile_x(idx)-1,'.',nint(ry),'-',nint(ry)+source_tile_y(idx)-1
2925 end if
2927 end subroutine get_tile_fname
2930 !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
2931 ! Name: get_source_resolution
2932 !
2933 ! Purpose:
2934 !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
2935 subroutine get_source_resolution(fieldnm, ilevel, src_dx, src_dy, istatus)
2937 implicit none
2939 ! Arguments
2940 integer, intent(in) :: ilevel
2941 integer, intent(out) :: istatus
2942 real, intent(out) :: src_dx, src_dy
2943 character (len=*), intent(in) :: fieldnm
2945 ! Local variables
2946 integer :: idx
2948 istatus = 1
2949 do idx=1,num_entries
2950 if ((index(source_fieldname(idx),trim(fieldnm)) /= 0) .and. &
2951 (len_trim(source_fieldname(idx)) == len_trim(fieldnm))) then
2952 if (ilevel == source_priority(idx)) then
2953 if (is_dx(idx) .and. is_dy(idx)) then
2954 src_dx = source_dx(idx)
2955 src_dy = source_dy(idx)
2956 if (source_proj(idx) /= PROJ_LATLON) then
2957 src_dx = src_dx / 111000.
2958 src_dy = src_dy / 111000.
2959 end if
2960 istatus = 0
2961 exit
2962 end if
2963 end if
2964 end if
2965 end do
2967 end subroutine get_source_resolution
2970 !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
2971 ! Name: get_data_projection
2972 !
2973 ! Purpose: To acquire the parameters necessary in defining the grid on which
2974 ! the user-specified data for field 'fieldnm' are given.
2975 !
2976 ! NOTES: If the routine successfully acquires values for all necessary
2977 ! parameters, istatus is set to 0. In case of an error,
2978 ! OR IF THE USER HAS NOT SPECIFIED A TILE OF DATA FOR FIELDNM,
2979 ! istatus is set to 1.
2980 !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
2981 subroutine get_data_projection(fieldnm, iproj, stand_lon, truelat1, truelat2, &
2982 dxkm, dykm, known_x, known_y, known_lat, known_lon, ilevel, istatus)
2984 implicit none
2986 ! Arguments
2987 integer, intent(in) :: ilevel
2988 integer, intent(out) :: iproj, istatus
2989 real, intent(out) :: stand_lon, truelat1, truelat2, dxkm, dykm, &
2990 known_x, known_y, known_lat, known_lon
2991 character (len=*), intent(in) :: fieldnm
2993 ! Local variables
2994 integer :: idx
2996 istatus = 1
2997 do idx=1,num_entries
2998 if ((index(source_fieldname(idx),trim(fieldnm)) /= 0) .and. &
2999 (len_trim(source_fieldname(idx)) == len_trim(fieldnm))) then
3000 if (ilevel == source_priority(idx)) then
3001 istatus = 0
3002 if (is_proj(idx)) then
3003 iproj = source_proj(idx)
3004 else
3005 iproj = 1
3006 call mprintf(.true., ERROR, &
3007 'In GEOGRID.TBL, no specification for projection in entry %i.',i1=idx)
3008 end if
3009 if (is_known_x(idx)) then
3010 known_x = source_known_x(idx)
3011 else
3012 known_x = 1.
3013 call mprintf(.true., ERROR, &
3014 'In GEOGRID.TBL, no specification for known_x in entry %i.',i1=idx)
3015 end if
3016 if (is_known_y(idx)) then
3017 known_y = source_known_y(idx)
3018 else
3019 known_y = 1.
3020 call mprintf(.true., ERROR, &
3021 'In GEOGRID.TBL, no specification for known_y in entry %i.',i1=idx)
3022 end if
3023 if (is_known_lat(idx)) then
3024 known_lat = source_known_lat(idx)
3025 else
3026 known_lat = 1.
3027 call mprintf(.true., ERROR, &
3028 'In GEOGRID.TBL, no specification for known_lat in entry %i.',i1=idx)
3029 end if
3030 if (is_known_lon(idx)) then
3031 known_lon = source_known_lon(idx)
3032 else
3033 known_lon = 1.
3034 call mprintf(.true., ERROR, &
3035 'In GEOGRID.TBL, no specification for known_lon in entry %i.',i1=idx)
3036 end if
3037 if (is_truelat1(idx)) then
3038 truelat1 = source_truelat1(idx)
3039 else if (is_proj(idx) .and. source_proj(idx) /= PROJ_LATLON) then
3040 truelat1 = 1.
3041 call mprintf(.true., WARN, &
3042 'In GEOGRID.TBL, no specification for truelat1 in entry %i.',i1=idx)
3043 end if
3044 if (is_truelat2(idx)) then
3045 truelat2 = source_truelat2(idx)
3046 else if (is_proj(idx) .and. source_proj(idx) /= PROJ_LATLON) then
3047 truelat2 = 1.
3048 call mprintf(.true., WARN, &
3049 'In GEOGRID.TBL, no specification for truelat2 in entry %i.',i1=idx)
3050 end if
3051 if (is_stdlon(idx)) then
3052 stand_lon = source_stdlon(idx)
3053 else if (is_proj(idx) .and. source_proj(idx) /= PROJ_LATLON) then
3054 stand_lon = 1.
3055 call mprintf(.true., WARN, &
3056 'In GEOGRID.TBL, no specification for stdlon in entry %i.',i1=idx)
3057 end if
3058 if (is_dx(idx)) then
3059 dxkm = source_dx(idx)
3060 else
3061 dxkm = 1.
3062 call mprintf(.true., ERROR, &
3063 'In GEOGRID.TBL, no specification for dx in entry %i.',i1=idx)
3064 end if
3065 if (is_dy(idx)) then
3066 dykm = source_dy(idx)
3067 else
3068 dykm = 1.
3069 call mprintf(.true., ERROR, &
3070 'In GEOGRID.TBL, no specification for dy in entry %i.',i1=idx)
3071 end if
3072 exit
3073 end if
3074 end if
3075 end do
3077 end subroutine get_data_projection
3080 !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
3081 ! Name: check_data_specification
3082 !
3083 ! Purpose: To check for obvious errors in the user source data specifications.
3084 ! Returns .true. if specification passes all checks, and .false. otherwise.
3085 ! For failing checks, diagnostic messages are printed.
3086 !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
3087 function check_data_specification( )
3089 implicit none
3091 ! Return value
3092 logical :: check_data_specification
3094 ! Local variables
3095 integer :: i, j, istatus
3096 integer, pointer, dimension(:) :: priorities
3097 real :: rmissing
3098 logical :: begin_priority, halt
3099 character (len=128) :: cur_name
3101 check_data_specification = .false.
3103 ! Check that each specification has a name, priority level, and path
3104 do i=1,num_entries
3105 if (.not. is_fieldname(i)) then
3106 call mprintf(.true., ERROR, &
3107 'In GEOGRID.TBL, specification %i does not have a name.',i1=i)
3108 end if
3109 if (.not. is_priority(i)) then
3110 call mprintf(.true., ERROR, &
3111 'In GEOGRID.TBL, specification %i does not have a priority.',i1=i)
3112 end if
3113 if (list_length(source_res_path(i)) == 0) then
3114 call mprintf(.true., ERROR, &
3115 'In GEOGRID.TBL, no path (relative or absolute) is specified '// &
3116 'for entry %i.',i1=i)
3117 end if
3118 end do
3120 ! The fill_missing and halt_on_missing options are mutually exclusive
3121 do i=1,num_entries
3122 call get_halt_on_missing(source_fieldname(i), halt, istatus)
3123 call get_missing_fill_value(source_fieldname(i), rmissing, istatus)
3124 if (halt .and. (istatus == 0)) then
3125 call mprintf(.true., ERROR, 'In GEOGRID.TBL, the halt_on_missing and fill_missing '// &
3126 'options are mutually exclusive, but both are given for field %s', &
3127 s1=trim(source_fieldname(i)))
3128 end if
3129 end do
3131 ! Check that the field from which landmask is calculated is not output on a staggering
3132 do i=1,num_entries
3133 if (list_length(source_landmask_land(i)) > 0 .or. list_length(source_landmask_water(i)) > 0) then
3134 if (is_output_stagger(i)) then
3135 if (source_output_stagger(i) /= M) then
3136 call mprintf(.true., ERROR, 'In GEOGRID.TBL, landmask cannot be derived from '// &
3137 'a field that is computed on a staggered grid at entry %i.',i1=i)
3138 end if
3139 end if
3140 end if
3141 end do
3143 ! Also check that any field that is to be masked by the landmask is not output on a staggering
3144 do i=1,num_entries
3145 if (is_masked(i) .and. is_output_stagger(i)) then
3146 if (source_output_stagger(i) /= M) then
3147 call mprintf(.true., ERROR, 'In GEOGRID.TBL, landmask cannot be used with '// &
3148 'a field that is computed on a staggered grid at entry %i.',i1=i)
3149 end if
3150 end if
3151 end do
3153 allocate(priorities(num_entries))
3155 ! Now check that priorities for each source are unique and in the interval [1,n], n <= num_entries
3156 do i=1,num_entries
3157 priorities = 0
3158 cur_name = source_fieldname(i)
3159 do j=1,num_entries
3160 if (source_fieldname(j) == cur_name) then
3162 if (source_priority(j) > num_entries .or. source_priority(j) < 1) then
3163 call mprintf(.true., ERROR, 'In GEOGRID.TBL, priorities for %s do not '// &
3164 'form a sequence 1,2,...,n.', s1=trim(cur_name))
3166 else
3167 if (priorities(source_priority(j)) == 1) then
3168 call mprintf(.true., ERROR, 'In GEOGRID.TBL, more than one entry for %s '// &
3169 'has priority %i.', s1=trim(cur_name), i1=source_priority(j))
3171 else
3172 priorities(source_priority(j)) = 1
3173 end if
3174 end if
3176 end if
3177 end do
3179 begin_priority = .false.
3180 do j=num_entries,1,-1
3181 if (.not.begin_priority .and. priorities(j) == 1) then
3182 begin_priority = .true.
3183 else if (begin_priority .and. priorities(j) == 0) then
3184 call mprintf(.true., ERROR, 'In GEOGRID.TBL, no entry for %s has '// &
3185 'priority %i, but an entry has priority %i.', &
3186 s1=trim(cur_name), i1=j, i2=j+1)
3187 end if
3188 end do
3189 end do
3191 deallocate(priorities)
3193 ! Units must match for all priority levels of a field
3194 do i=1,num_entries
3195 if (source_priority(i) == 1) then
3196 do j=1,num_entries
3197 if ((source_fieldname(i) == source_fieldname(j)) .and. &
3198 (source_units(i) /= source_units(j))) then
3199 call mprintf(.true., ERROR, 'In GEOGRID.TBL, units for %s at entry %i '// &
3200 'do not match units at entry %i (%s)', &
3201 s1=trim(source_fieldname(i)), i1=j, i2=i, s2=trim(source_units(i)))
3202 end if
3203 end do
3204 end if
3205 end do
3207 ! Make sure that user has not asked to calculate landmask from a continuous field
3208 do i=1,num_entries
3209 if (is_dest_fieldtype(i)) then
3210 if (source_dest_fieldtype(i) == CONTINUOUS) then
3211 if (list_length(source_landmask_water(i)) > 0 .or. list_length(source_landmask_land(i)) > 0) then
3212 call mprintf(.true., ERROR, 'In GEOGRID.TBL, landmask cannot be '// &
3213 'calculated from a continuous destination field at entry %i.',i1=i)
3214 end if
3215 end if
3216 end if
3217 end do
3219 ! If either min_category or max_category is specified, then both must be specified
3220 do i=1,num_entries
3221 if (is_category_min(i) .or. is_category_max(i)) then
3222 if (.not. is_category_min(i)) then
3223 call mprintf(.true., ERROR, 'In GEOGRID.TBL, for index file of data at '// &
3224 'entry %i, category_max is specified, but category_min is '// &
3225 'not. Both must be specified.',i1=i)
3226 else if (.not. is_category_max(i)) then
3227 call mprintf(.true., ERROR, 'In GEOGRID.TBL, for index file of data at '// &
3228 'entry %i, category_min is specified, but category_max is '// &
3229 'not. Both must be specified.',i1=i)
3230 end if
3231 end if
3232 end do
3234 ! For continuous data, (category_max - category_min + 1) should equal tile_z
3235 do i=1,num_entries
3236 if (is_fieldtype(i)) then
3237 if (source_fieldtype(i) == CONTINUOUS) then
3238 if (is_category_max(i) .and. is_category_min(i) .and. is_tile_z(i)) then
3239 if (source_tile_z(i) /= (source_category_max(i) - source_category_min(i) + 1)) then
3240 call mprintf(.true., ERROR, 'In GEOGRID.TBL, tile_z must equal '// &
3241 '(category_max - category_min + 1) at entry %i.',i1=i)
3242 end if
3243 else if (is_category_max(i) .and. is_category_min(i) .and. &
3244 is_tile_z_start(i) .and. is_tile_z_end(i)) then
3245 if (source_tile_z_end(i) /= source_category_max(i) .or. &
3246 source_tile_z_start(i) /= source_category_min(i)) then
3247 call mprintf(.true., ERROR, 'In GEOGRID.TBL, tile_z_end must equal '// &
3248 'category_max, and tile_z_start must equal category_min '// &
3249 'at entry %i.',i1=i)
3250 end if
3251 end if
3252 end if
3253 end if
3254 end do
3256 ! Make sure that user has not named a dominant category or computed slope field
3257 ! the same as a fractional field
3258 do i=1,num_entries
3259 if (source_dominant_category(i) == source_fieldname(i)) then
3260 call mprintf(.true., ERROR, 'In GEOGRID.TBL, dominant category cannot have '// &
3261 'the same name as the field at entry %i.',i1=i)
3262 end if
3264 do j=1,num_entries
3265 if (.not. is_dominant_only(i)) then
3266 if (is_dfdx(j)) then
3267 if (source_dfdx(j) == source_fieldname(i)) then
3268 call mprintf(.true., ERROR, 'In GEOGRID.TBL, field name at entry %i '// &
3269 'cannot have the same name as the slope field df_dx at entry %i.', &
3270 i1=i, i2=j)
3271 end if
3272 end if
3273 if (is_dfdy(j)) then
3274 if (source_dfdy(j) == source_fieldname(i)) then
3275 call mprintf(.true., ERROR, 'In GEOGRID.TBL, field name at entry %i '// &
3276 'cannot have the same name as the slope field df_dy at entry %i.', &
3277 i1=i, i2=j)
3278 end if
3279 end if
3280 if (is_dfdx(j) .and. is_dominant_category(i)) then
3281 if (source_dfdx(j) == source_dominant_category(i)) then
3282 call mprintf(.true., ERROR, 'In GEOGRID.TBL, dominant field name at '// &
3283 'entry %i cannot have the same name as the slope field df_dx '// &
3284 'at entry %i.',i1=i, i2=j)
3285 end if
3286 end if
3287 if (is_dfdy(j) .and. is_dominant_category(i)) then
3288 if (source_dfdy(j) == source_dominant_category(i)) then
3289 call mprintf(.true., ERROR, 'In GEOGRID.TBL, dominant field name at '// &
3290 'entry %i cannot have the same name as the slope field '// &
3291 'df_dy at entry %i.',i1=i, i2=j)
3292 end if
3293 end if
3294 else
3295 if (is_dfdx(j)) then
3296 if (source_dfdx(j) == source_dominant_only(i)) then
3297 call mprintf(.true., ERROR, 'In GEOGRID.TBL, dominant field name at '// &
3298 'entry %i cannot have the same name as the slope field '// &
3299 'df_dx at entry %i.',i1=i, i2=j)
3300 end if
3301 end if
3302 if (is_dfdy(j)) then
3303 if (source_dfdy(j) == source_dominant_only(i)) then
3304 call mprintf(.true., ERROR, 'In GEOGRID.TBL, dominant field name at '// &
3305 'entry %i cannot have the same name as the slope field '// &
3306 'df_dy at entry %i.',i1=i, i2=j)
3307 end if
3308 end if
3309 end if
3310 if (i /= j) then
3311 if (is_dfdx(i)) then
3312 if (is_dfdx(j)) then
3313 if (source_dfdx(j) == source_dfdx(i)) then
3314 call mprintf(.true., ERROR, 'In GEOGRID.TBL, slope field df_dx at '// &
3315 'entry %i cannot have the same name as the slope '// &
3316 'field df_dx at entry %i.',i1=i, i2=j)
3317 end if
3318 end if
3319 if (is_dfdy(j)) then
3320 if (source_dfdy(j) == source_dfdx(i)) then
3321 call mprintf(.true., ERROR, 'In GEOGRID.TBL, slope field df_dx at '// &
3322 'entry %i cannot have the same name as the slope field '// &
3323 'df_dy at entry %i.',i1=i, i2=j)
3324 end if
3325 end if
3326 end if
3327 if (is_dfdy(i)) then
3328 if (is_dfdx(j)) then
3329 if (source_dfdx(j) == source_dfdy(i)) then
3330 call mprintf(.true., ERROR, 'In GEOGRID.TBL, slope field df_dy at '// &
3331 'entry %i cannot have the same name as the slope field '// &
3332 'df_dx at entry %i.',i1=i, i2=j)
3333 end if
3334 end if
3335 if (is_dfdy(j)) then
3336 if (source_dfdy(j) == source_dfdy(i)) then
3337 call mprintf(.true., ERROR, 'In GEOGRID.TBL, slope field df_dy at '// &
3338 'entry %i cannot have the same name as the slope field '// &
3339 'df_dy at entry %i.',i1=i, i2=j)
3340 end if
3341 end if
3342 end if
3343 if (is_dominant_category(i)) then
3344 if (source_dominant_category(i) == source_fieldname(j)) then ! Possible exception
3345 if (.not. (is_dominant_only(j) .and. (source_dominant_category(i) /= source_dominant_only(j)))) then
3346 call mprintf(.true., ERROR, 'In GEOGRID.TBL, dominant category at '// &
3347 'entry %i cannot have the same name as the field at '// &
3348 'entry %i.',i1=i, i2=j)
3349 end if
3350 else if (is_dominant_category(j) .and. &
3351 (source_dominant_category(i) == source_dominant_category(j))) then
3352 call mprintf(.true., ERROR, 'In GEOGRID.TBL, dominant category at entry '// &
3353 '%i cannot have the same name as dominant category at '// &
3354 'entry %i.',i1=i, i2=j)
3355 else if (is_dominant_only(j) .and. &
3356 (source_dominant_category(i) == source_dominant_only(j))) then
3357 call mprintf(.true., ERROR, 'In GEOGRID.TBL, dominant category at '// &
3358 'entry %i cannot have the same name as dominant_only '// &
3359 'category at entry %i.',i1=i, i2=j)
3360 end if
3361 else if (is_dominant_only(i)) then
3362 if (source_dominant_only(i) == source_fieldname(j)) then ! Possible exception
3363 if (.not. (is_dominant_only(j) .and. (source_dominant_only(i) /= source_dominant_only(j)))) then
3364 call mprintf(.true., ERROR, 'In GEOGRID.TBL, dominant_only category '// &
3365 'at entry %i cannot have the same name as the field at '// &
3366 'entry %i.',i1=i, i2=j)
3367 end if
3368 else if (is_dominant_category(j) .and. &
3369 (source_dominant_only(i) == source_dominant_category(j))) then
3370 call mprintf(.true., ERROR, 'In GEOGRID.TBL, dominant_only category '// &
3371 'at entry %i cannot have the same name as dominant '// &
3372 'category at entry %i.',i1=i, i2=j)
3373 else if (is_dominant_only(j) .and. &
3374 (source_dominant_only(i) == source_dominant_only(j))) then
3375 call mprintf(.true., ERROR, 'In GEOGRID.TBL, dominant_only category '// &
3376 'at entry %i cannot have the same name as dominant_only '// &
3377 'category at entry %i.',i1=i, i2=j)
3378 end if
3379 end if
3380 end if
3381 end do
3382 end do
3384 check_data_specification = .true.
3385 end function check_data_specification
3388 !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
3389 ! Name: s_len
3390 !
3391 ! Purpose: This routine receives a fortran string, and returns the number of
3392 ! characters in the string before the first "space" is encountered. It
3393 ! considers ascii characters 33 to 126 to be valid characters, and ascii
3394 ! 0 to 32, and 127 to be "space" characters.
3395 !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
3396 subroutine s_len(string, s_length)
3398 implicit none
3400 ! Arguments
3401 character (len=*), intent(in) :: string
3402 integer, intent(out) :: s_length
3404 ! Local variables
3405 integer :: i, len_str, aval
3406 logical :: space
3408 space = .false.
3409 i = 1
3410 len_str = len(string)
3411 s_length = len_str
3412 do while ((i .le. len_str) .and. (.not. space))
3413 aval = ichar(string(i:i))
3414 if ((aval .lt. 33) .or. (aval .gt. 126)) then
3415 s_length = i - 1
3416 space = .true.
3417 endif
3418 i = i + 1
3419 enddo
3421 end subroutine s_len
3423 end module source_data_module