| blob | dc2ae902977227195e7ffee414d4e9e196e4ec76 |
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, source_geotiff_file
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, is_geotiff
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 allocate(source_geotiff_file(num_entries))
175 do i=1,num_entries
176 call list_init(source_res_path(i))
177 call list_init(source_interp_option(i))
178 call list_init(source_landmask_land(i))
179 call list_init(source_landmask_water(i))
180 end do
182 allocate(is_wordsize(num_entries))
183 allocate(is_endian(num_entries))
184 allocate(is_fieldtype(num_entries))
185 allocate(is_dest_fieldtype(num_entries))
186 allocate(is_proj(num_entries))
187 allocate(is_priority(num_entries))
188 allocate(is_dx(num_entries))
189 allocate(is_dy(num_entries))
190 allocate(is_known_x(num_entries))
191 allocate(is_known_y(num_entries))
192 allocate(is_known_lat(num_entries))
193 allocate(is_known_lon(num_entries))
194 allocate(is_truelat1(num_entries))
195 allocate(is_truelat2(num_entries))
196 allocate(is_stdlon(num_entries))
197 allocate(is_fieldname(num_entries))
198 allocate(is_path(num_entries))
199 allocate(is_tile_x(num_entries))
200 allocate(is_tile_y(num_entries))
201 allocate(is_tile_z(num_entries))
202 allocate(is_tile_z_start(num_entries))
203 allocate(is_tile_z_end(num_entries))
204 allocate(is_category_min(num_entries))
205 allocate(is_category_max(num_entries))
206 allocate(is_tile_bdr(num_entries))
207 allocate(is_masked(num_entries))
208 allocate(is_halt_missing(num_entries))
209 allocate(is_output_stagger(num_entries))
210 allocate(is_row_order(num_entries))
211 allocate(is_dominant_category(num_entries))
212 allocate(is_dominant_only(num_entries))
213 allocate(is_dfdx(num_entries))
214 allocate(is_dfdy(num_entries))
215 allocate(is_scale_factor(num_entries))
216 allocate(is_z_dim_name(num_entries))
217 allocate(is_units(num_entries))
218 allocate(is_descr(num_entries))
219 allocate(is_smooth_option(num_entries))
220 allocate(is_smooth_passes(num_entries))
221 allocate(is_signed(num_entries))
222 allocate(is_missing_value(num_entries))
223 allocate(is_fill_missing(num_entries))
224 allocate(is_subgrid(num_entries))
225 allocate(is_geotiff(num_entries))
227 source_wordsize=0
228 source_endian=0
229 source_fieldtype=0
230 source_dest_fieldtype=0
231 source_proj=0
232 source_priority=0
233 source_dx=0
234 source_dy=0
235 source_known_x=0
236 source_known_y=0
237 source_known_lat=0
238 source_known_lon=0
239 source_truelat1=0
240 source_truelat2=0
241 source_stdlon=0
242 source_fieldname=' '
243 source_path=' '
244 source_interp_string=' '
245 source_tile_x=0
246 source_tile_y=0
247 source_tile_z=0
248 source_tile_z_start=0
249 source_tile_z_end=0
250 source_category_min=0
251 source_category_max=0
252 source_tile_bdr=0
253 source_masked=0
254 source_output_stagger=0
255 source_row_order=0
256 source_dominant_category=' '
257 source_dominant_only=' '
258 source_dfdx=' '
259 source_dfdy=' '
260 source_scale_factor=0
261 source_z_dim_name=' '
262 source_units=' '
263 source_descr=' '
264 source_smooth_option=0
265 source_smooth_passes=0
266 source_missing_value=0
267 source_fill_missing=0
269 is_wordsize=.false.
270 is_endian=.false.
271 is_fieldtype=.false.
272 is_dest_fieldtype=.false.
273 is_proj=.false.
274 is_priority=.false.
275 is_dx=.false.
276 is_dy=.false.
277 is_known_x=.false.
278 is_known_y=.false.
279 is_known_lat=.false.
280 is_known_lon=.false.
281 is_truelat1=.false.
282 is_truelat2=.false.
283 is_stdlon=.false.
284 is_fieldname=.false.
285 is_path=.false.
286 is_tile_x=.false.
287 is_tile_y=.false.
288 is_tile_z=.false.
289 is_tile_z_start=.false.
290 is_tile_z_end=.false.
291 is_category_min=.false.
292 is_category_max=.false.
293 is_tile_bdr=.false.
294 is_masked=.false.
295 is_halt_missing=.false.
296 is_output_stagger=.false.
297 is_row_order=.false.
298 is_dominant_category=.false.
299 is_dominant_only=.false.
300 is_dfdx=.false.
301 is_dfdy=.false.
302 is_scale_factor=.false.
303 is_z_dim_name=.false.
304 is_units=.false.
305 is_descr=.false.
306 is_smooth_option=.false.
307 is_smooth_passes=.false.
308 is_signed=.false.
309 is_missing_value=.false.
310 is_fill_missing=.false.
311 is_subgrid=.false.
313 write(source_mminlu,'(a4)') 'USGS'
314 source_iswater = 16
315 source_islake = -1
316 source_isice = 24
317 source_isurban = 1
318 source_isoilwater = 14
320 !
321 ! Actually read and save the specifications
322 !
323 nparams = 0
324 i = 1
325 30 buffer = ' '
326 read(funit,'(a)',end=40,err=1001) buffer
327 call despace(buffer)
329 ! Is this line a comment?
330 if (buffer(1:1) == '#') then
331 ! Do nothing.
333 ! Are we beginning a new field specification?
334 else if (index(buffer,'=====') /= 0) then !{
335 if (nparams > 0) i = i + 1
336 nparams = 0
337 if (i <= num_entries) then
338 !BUG: Are these initializations needed now that we've added initializations for
339 ! all options after their initial allocation above?
340 is_path(i) = .false.
341 is_masked(i) = .false.
342 is_halt_missing(i) = .false.
343 is_output_stagger(i) = .false.
344 is_dominant_category(i) = .false.
345 is_dominant_only(i) = .false.
346 is_dfdx(i) = .false.
347 is_dfdy(i) = .false.
348 is_dest_fieldtype(i) = .false.
349 is_priority(i) = .false.
350 is_z_dim_name(i) = .false.
351 is_smooth_option(i) = .false.
352 is_smooth_passes(i) = .false.
353 is_fill_missing(i) = .false.
354 is_subgrid(i) = .false.
355 end if
357 else
358 ! Check whether the current line is a comment
359 if (buffer(1:1) /= '#') then
360 have_specification = .true.
361 else
362 have_specification = .false.
363 end if
365 ! If only part of the line is a comment, just turn the comment into spaces
366 eos = index(buffer,'#')
367 if (eos /= 0) buffer(eos:BUFSIZE) = ' '
369 do while (have_specification) !{
371 ! If this line has no semicolon, it may contain a single specification,
372 ! so we set have_specification = .false. to prevent the line from being
373 ! processed again and "pretend" that the last character was a semicolon
374 eos = index(buffer,';')
375 if (eos == 0) then
376 have_specification = .false.
377 eos = BUFSIZE
378 end if
380 idx = index(buffer(1:eos-1),'=')
382 if (idx /= 0) then !{
383 nparams = nparams + 1
385 if (index('name',trim(buffer(1:idx-1))) /= 0) then
386 ispace = idx+1
387 do while ((ispace < eos) .and. (buffer(ispace:ispace) /= ' '))
388 ispace = ispace + 1
389 end do
390 is_fieldname(i) = .true.
391 source_fieldname(i) = ' '
392 source_fieldname(i)(1:ispace-idx) = buffer(idx+1:ispace-1)
394 else if (index('priority',trim(buffer(1:idx-1))) /= 0) then
395 is_priority(i) = .true.
396 read(buffer(idx+1:eos-1),'(i10)') source_priority(i)
398 else if (index('dest_type',trim(buffer(1:idx-1))) /= 0) then
399 if (index('continuous',trim(buffer(idx+1:eos-1))) /= 0) then
400 is_dest_fieldtype(i) = .true.
401 source_dest_fieldtype(i) = CONTINUOUS
402 else if (index('categorical',trim(buffer(idx+1:eos-1))) /= 0) then
403 is_dest_fieldtype(i) = .true.
404 source_dest_fieldtype(i) = CATEGORICAL
405 end if
407 else if (index('interp_option',trim(buffer(1:idx-1))) /= 0) then
408 ispace = idx+1
409 do while ((ispace < eos) .and. (buffer(ispace:ispace) /= ' '))
410 ispace = ispace + 1
411 end do
412 interp_string = ' '
413 interp_string(1:ispace-idx-1) = buffer(idx+1:ispace-1)
414 ispace = index(interp_string,':')
415 if (ispace /= 0) then
416 write(res_string,'(a)') interp_string(1:ispace-1)
417 else
418 res_string = 'default'
419 end if
420 write(interp_string,'(a)') trim(interp_string(ispace+1:128))
421 if (list_search(source_interp_option(i), ckey=res_string, cvalue=interp_string)) then
422 call mprintf(.true., WARN, &
423 'In GEOGRID.TBL entry %i, multiple interpolation methods are '// &
424 'given for resolution %s. %s will be used.', &
425 i1=i, s1=trim(res_string), s2=trim(interp_string))
426 else
427 call list_insert(source_interp_option(i), ckey=res_string, cvalue=interp_string)
428 end if
430 else if (index('smooth_option',trim(buffer(1:idx-1))) /= 0) then
431 if ((index('1-2-1',trim(buffer(idx+1:eos-1))) /= 0) .and. &
432 (len_trim(buffer(idx+1:eos-1)) == len('1-2-1'))) then
433 is_smooth_option(i) = .true.
434 source_smooth_option(i) = ONETWOONE
435 else if ((index('smth-desmth',trim(buffer(idx+1:eos-1))) /= 0) .and. &
436 (len_trim(buffer(idx+1:eos-1)) == len('smth-desmth'))) then
437 is_smooth_option(i) = .true.
438 source_smooth_option(i) = SMTHDESMTH
439 else if ((index('smth-desmth_special',trim(buffer(idx+1:eos-1))) /= 0) .and. &
440 (len_trim(buffer(idx+1:eos-1)) == len('smth-desmth_special'))) then
441 is_smooth_option(i) = .true.
442 source_smooth_option(i) = SMTHDESMTH_SPECIAL
443 end if
445 else if (index('smooth_passes',trim(buffer(1:idx-1))) /= 0) then
446 is_smooth_passes(i) = .true.
447 read(buffer(idx+1:eos-1),'(i10)') source_smooth_passes(i)
449 else if (index('rel_path',trim(buffer(1:idx-1))) /= 0) then
450 ispace = idx+1
451 do while ((ispace < eos) .and. (buffer(ispace:ispace) /= ' '))
452 ispace = ispace + 1
453 end do
454 path_string = ' '
455 path_string(1:ispace-idx-1) = buffer(idx+1:ispace-1)
456 if (path_string(ispace-idx-1:ispace-idx-1) /= '/') &
457 path_string(ispace-idx:ispace-idx) = '/'
458 if (path_string(1:1) == '/') then
459 path_string(1:127) = path_string(2:128)
460 path_string(128:128) = ' '
461 end if
462 ispace = index(path_string,':')
463 if (ispace /= 0) then
464 write(res_string,'(a)') path_string(1:ispace-1)
465 else
466 res_string = 'default'
467 end if
468 write(path_string,'(a)') trim(geog_data_path)//trim(path_string(ispace+1:128))
469 if (list_search(source_res_path(i), ckey=res_string, cvalue=path_string)) then
470 call mprintf(.true., WARN, &
471 'In GEOGRID.TBL entry %i, multiple paths are given for '// &
472 'resolution %s. %s will be used.', &
473 i1=i, s1=trim(res_string), s2=trim(path_string))
474 else
475 call list_insert(source_res_path(i), ckey=res_string, cvalue=path_string)
476 end if
478 else if (index('abs_path',trim(buffer(1:idx-1))) /= 0) then
479 ispace = idx+1
480 do while ((ispace < eos) .and. (buffer(ispace:ispace) /= ' '))
481 ispace = ispace + 1
482 end do
483 path_string = ' '
484 path_string(1:ispace-idx-1) = buffer(idx+1:ispace-1)
485 if (path_string(ispace-idx-1:ispace-idx-1) /= '/') &
486 path_string(ispace-idx:ispace-idx) = '/'
487 ispace = index(path_string,':')
488 if (ispace /= 0) then
489 write(res_string,'(a)') path_string(1:ispace-1)
490 else
491 res_string = 'default'
492 end if
493 write(path_string,'(a)') trim(path_string(ispace+1:128))
494 if (list_search(source_res_path(i), ckey=res_string, cvalue=path_string)) then
495 call mprintf(.true., WARN, &
496 'In GEOGRID.TBL entry %i, multiple paths are given for '// &
497 'resolution %s. %s will be used.', &
498 i1=i, s1=trim(res_string), s2=trim(path_string))
499 else
500 call list_insert(source_res_path(i), ckey=res_string, cvalue=path_string)
501 end if
503 else if (index('output_stagger',trim(buffer(1:idx-1))) /= 0) then
504 if (index('M',trim(buffer(idx+1:eos-1))) /= 0) then
505 is_output_stagger(i) = .true.
506 source_output_stagger(i) = M
507 else if (index('U',trim(buffer(idx+1:eos-1))) /= 0) then
508 is_output_stagger(i) = .true.
509 source_output_stagger(i) = U
510 else if (index('V',trim(buffer(idx+1:eos-1))) /= 0) then
511 is_output_stagger(i) = .true.
512 source_output_stagger(i) = V
513 else if (index('HH',trim(buffer(idx+1:eos-1))) /= 0) then
514 is_output_stagger(i) = .true.
515 source_output_stagger(i) = HH
516 else if (index('VV',trim(buffer(idx+1:eos-1))) /= 0) then
517 is_output_stagger(i) = .true.
518 source_output_stagger(i) = VV
519 end if
521 else if ((index('landmask_water',trim(buffer(1:idx-1))) /= 0) .and. &
522 (len_trim(buffer(1:idx-1)) == 14)) then
523 ispace = idx+1
524 do while ((ispace < eos) .and. (buffer(ispace:ispace) /= ' '))
525 ispace = ispace + 1
526 end do
527 landmask_string = ' '
528 landmask_string(1:ispace-idx-1) = buffer(idx+1:ispace-1)
529 ispace = index(landmask_string,':')
530 if (ispace /= 0) then
531 write(res_string,'(a)') landmask_string(1:ispace-1)
532 else
533 res_string = 'default'
534 end if
535 write(landmask_string,'(a)') trim(landmask_string(ispace+1:128))
536 if (list_search(source_landmask_water(i), ckey=res_string, cvalue=landmask_string)) then
537 call mprintf(.true., WARN, &
538 'In GEOGRID.TBL entry %i, multiple landmask category specifications are given for '// &
539 'resolution %s. %s will be used.', &
540 i1=i, s1=trim(res_string), s2=trim(landmask_string))
541 else
542 call list_insert(source_landmask_water(i), ckey=res_string, cvalue=landmask_string)
543 end if
545 else if ((index('landmask_land',trim(buffer(1:idx-1))) /= 0) .and. &
546 (len_trim(buffer(1:idx-1)) == 13)) then
547 ispace = idx+1
548 do while ((ispace < eos) .and. (buffer(ispace:ispace) /= ' '))
549 ispace = ispace + 1
550 end do
551 landmask_string = ' '
552 landmask_string(1:ispace-idx-1) = buffer(idx+1:ispace-1)
553 ispace = index(landmask_string,':')
554 if (ispace /= 0) then
555 write(res_string,'(a)') landmask_string(1:ispace-1)
556 else
557 res_string = 'default'
558 end if
559 write(landmask_string,'(a)') trim(landmask_string(ispace+1:128))
560 if (list_search(source_landmask_land(i), ckey=res_string, cvalue=landmask_string)) then
561 call mprintf(.true., WARN, &
562 'In GEOGRID.TBL entry %i, multiple landmask category specifications are given for '// &
563 'resolution %s. %s will be used.', &
564 i1=i, s1=trim(res_string), s2=trim(landmask_string))
565 else
566 call list_insert(source_landmask_land(i), ckey=res_string, cvalue=landmask_string)
567 end if
569 else if ((index('masked',trim(buffer(1:idx-1))) /= 0) .and. &
570 (len_trim(buffer(1:idx-1)) == 6)) then
571 if (index('water',trim(buffer(idx+1:eos-1))) /= 0) then
572 is_masked(i) = .true.
573 source_masked(i) = 0.
574 else if (index('land',trim(buffer(idx+1:eos-1))) /= 0) then
575 is_masked(i) = .true.
576 source_masked(i) = 1.
577 end if
579 else if (index('fill_missing',trim(buffer(1:idx-1))) /= 0) then
580 is_fill_missing(i) = .true.
581 read(buffer(idx+1:eos-1),*) source_fill_missing(i)
583 else if (index('halt_on_missing',trim(buffer(1:idx-1))) /= 0) then
584 if (index('yes',trim(buffer(idx+1:eos-1))) /= 0) then
585 is_halt_missing(i) = .true.
586 else if (index('no',trim(buffer(idx+1:eos-1))) /= 0) then
587 is_halt_missing(i) = .false.
588 end if
590 else if (index('dominant_category',trim(buffer(1:idx-1))) /= 0) then
591 ispace = idx+1
592 do while ((ispace < eos) .and. (buffer(ispace:ispace) /= ' '))
593 ispace = ispace + 1
594 end do
595 is_dominant_category(i) = .true.
596 source_dominant_category(i) = ' '
597 source_dominant_category(i)(1:ispace-idx) = buffer(idx+1:ispace-1)
599 else if (index('dominant_only',trim(buffer(1:idx-1))) /= 0) then
600 ispace = idx+1
601 do while ((ispace < eos) .and. (buffer(ispace:ispace) /= ' '))
602 ispace = ispace + 1
603 end do
604 is_dominant_only(i) = .true.
605 source_dominant_only(i) = ' '
606 source_dominant_only(i)(1:ispace-idx) = buffer(idx+1:ispace-1)
608 else if (index('df_dx',trim(buffer(1:idx-1))) /= 0) then
609 ispace = idx+1
610 do while ((ispace < eos) .and. (buffer(ispace:ispace) /= ' '))
611 ispace = ispace + 1
612 end do
613 is_dfdx(i) = .true.
614 source_dfdx(i) = ' '
615 source_dfdx(i)(1:ispace-idx) = buffer(idx+1:ispace-1)
617 else if (index('df_dy',trim(buffer(1:idx-1))) /= 0) then
618 ispace = idx+1
619 do while ((ispace < eos) .and. (buffer(ispace:ispace) /= ' '))
620 ispace = ispace + 1
621 end do
622 is_dfdy(i) = .true.
623 source_dfdy(i) = ' '
624 source_dfdy(i)(1:ispace-idx) = buffer(idx+1:ispace-1)
626 else if (index('z_dim_name',trim(buffer(1:idx-1))) /= 0) then
627 ispace = idx+1
628 do while ((ispace < eos) .and. (buffer(ispace:ispace) /= ' '))
629 ispace = ispace + 1
630 end do
631 is_z_dim_name(i) = .true.
632 source_z_dim_name(i) = ' '
633 source_z_dim_name(i)(1:ispace-idx) = buffer(idx+1:ispace-1)
635 else if (index('subgrid',trim(buffer(1:idx-1))) /= 0) then
636 ispace = idx+1
637 do while ((ispace < eos) .and. (buffer(ispace:ispace) /= ' '))
638 ispace = ispace + 1
639 end do
640 if (index('yes',trim(buffer(idx+1:eos-1))) /= 0) then
641 is_subgrid(i) = .true.
642 else if (index('no',trim(buffer(idx+1:eos-1))) /= 0) then
643 is_subgrid(i) = .false.
644 end if
646 else
647 call mprintf(.true., WARN, 'In GEOGRID.TBL, unrecognized option %s in '// &
648 'entry %i.',i1=idx, s1=buffer(i:idx-1))
650 end if
652 end if !} index(buffer(1:eos-1),'=') /= 0
654 buffer = buffer(eos+1:BUFSIZE)
655 end do ! while eos /= 0 }
657 end if !} index(buffer, '=====') /= 0
658 go to 30
660 40 close(funit)
662 ! Check the user specifications for gross errors
663 if ( .not. check_data_specification() ) then
664 call datalist_destroy()
665 call mprintf(.true.,ERROR,'Errors were found in either index files or GEOGRID.TBL.')
666 end if
668 call hash_init(bad_files)
670 return
672 1000 call mprintf(.true.,ERROR,'Could not open GEOGRID.TBL')
674 1001 call mprintf(.true.,ERROR,'Encountered error while reading GEOGRID.TBL')
676 end subroutine get_datalist
679 !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
680 ! Name: get_source_params
681 !
682 ! Purpose: For each field, this routine reads in the metadata in the index file
683 ! for the first available resolution of data specified by res_string. Also
684 ! based on res_string, this routine sets the interpolation sequence for the
685 ! field. This routine should be called prior to processing a field for each
686 ! domain.
687 !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
688 subroutine get_source_params(res_string)
690 use stringutil
691 #ifdef _HAS_GEOTIFF
692 use geotiff_module, only : open_geotiff,merge_geotiff_header
693 #endif
695 implicit none
697 ! Parameters
698 integer, parameter :: BUFSIZE = 256
700 ! Arguments
701 character (len=128), intent(in) :: res_string
703 ! Local variables
704 integer :: idx, i, is, ie, ispace, eos, iquoted, funit
705 character (len=128) :: temp_data, temp_interp
706 character (len=256) :: test_fname
707 character (len=BUFSIZE) :: buffer
708 logical :: have_specification, is_used
710 ! For each entry in the GEOGRID.TBL file
711 do idx=1,num_entries
713 ! Initialize metadata
714 is_wordsize(idx) = .false.
715 is_endian(idx) = .false.
716 is_row_order(idx) = .false.
717 is_fieldtype(idx) = .false.
718 is_proj(idx) = .false.
719 is_dx(idx) = .false.
720 is_dy(idx) = .false.
721 is_known_x(idx) = .false.
722 is_known_y(idx) = .false.
723 is_known_lat(idx) = .false.
724 is_known_lon(idx) = .false.
725 is_truelat1(idx) = .false.
726 is_truelat2(idx) = .false.
727 is_stdlon(idx) = .false.
728 is_tile_x(idx) = .false.
729 is_tile_y(idx) = .false.
730 is_tile_z(idx) = .false.
731 is_tile_z_start(idx) = .false.
732 is_tile_z_end(idx) = .false.
733 is_category_min(idx) = .false.
734 is_category_max(idx) = .false.
735 is_tile_bdr(idx) = .false.
736 is_fieldname(idx) = .false.
737 is_scale_factor(idx) = .false.
738 is_units(idx) = .false.
739 is_descr(idx) = .false.
740 is_signed(idx) = .false.
741 is_missing_value(idx) = .false.
742 is_geotiff(idx) = .false.
745 ! Set the interpolator sequence for the field to be the first value in res_string that matches
746 ! the resolution keyword for an interp_sequence specification
747 is = 1
748 ie = index(res_string(is:128),'+') - 1
749 if (ie <= 0) ie = 128
750 temp_interp = res_string(is:ie)
751 do while (.not. list_search(source_interp_option(idx), ckey=temp_interp, cvalue=source_interp_string(idx)) &
752 .and. is <= 128)
753 call mprintf(.true., INFORM, 'For %s, couldn''t find interpolator sequence for '// &
754 'resolution %s.', &
755 s1=trim(source_fieldname(idx)), s2=trim(temp_interp))
756 is = ie+2
757 ie = is + index(res_string(is:128),'+') - 2
758 if (ie - is <= 0) ie = 128
759 temp_interp = res_string(is:ie)
760 end do
762 if (is > 128) then
763 temp_interp = 'default'
764 if (list_search(source_interp_option(idx), ckey=temp_interp, cvalue=source_interp_string(idx))) then
765 call mprintf(.true., INFORM, 'Using default interpolator sequence for %s.', &
766 s1=trim(source_fieldname(idx)))
767 else
768 call mprintf(.true., ERROR, 'Could not find interpolator sequence for requested resolution for %s.'// &
769 ' The sources specified in namelist.wps is not listed in GEOGRID.TBL.', &
770 s1=trim(source_fieldname(idx)))
771 end if
772 else
773 call mprintf(.true., INFORM, 'Using %s interpolator sequence for %s.', &
774 s1=temp_interp, s2=trim(source_fieldname(idx)))
775 end if
777 ! Set the data source for the field to be the first value in res_string that matches
778 ! the resolution keyword for an abs_path or rel_path specification
779 is = 1
780 ie = index(res_string(is:128),'+') - 1
781 if (ie <= 0) ie = 128
782 temp_data = res_string(is:ie)
783 do while (.not. list_search(source_res_path(idx), ckey=temp_data, cvalue=source_path(idx)) &
784 .and. is <= 128)
785 call mprintf(.true., INFORM, 'For %s, couldn''t find %s data source.', &
786 s1=trim(source_fieldname(idx)), s2=trim(temp_data))
787 is = ie+2
788 ie = is + index(res_string(is:128),'+') - 2
789 if (ie - is <= 0) ie = 128
790 temp_data = res_string(is:ie)
791 end do
793 if (is > 128) then
794 temp_data = 'default'
795 if (list_search(source_res_path(idx), ckey=temp_data, cvalue=source_path(idx))) then
796 call mprintf(.true., INFORM, 'Using default data source for %s.', &
797 s1=trim(source_fieldname(idx)))
798 else
799 call mprintf(.true., ERROR, 'Could not find data resolution for requested resolution for %s. '// &
800 'The source specified in namelist.wps is not listed in GEOGRID.TBL.', &
801 s1=trim(source_fieldname(idx)))
802 end if
803 else
804 call mprintf(.true., INFORM, 'Using %s data source for %s.', &
805 s1=temp_data, s2=trim(source_fieldname(idx)))
806 end if
808 call mprintf(trim(temp_data) /= trim(temp_interp),WARN,'For %s, using %s data source with %s interpolation sequence.', &
809 s1=source_fieldname(idx), s2=temp_data, s3=temp_interp)
811 write(test_fname, '(a)') trim(source_path(idx))//'index'
813 !
814 ! Open the index file for the data source for this field, and read in metadata specs
815 !
816 do funit=10,100
817 inquire(unit=funit, opened=is_used)
818 if (.not. is_used) exit
819 end do
820 open(funit,file=trim(test_fname),form='formatted',status='old',err=1000)
822 30 buffer = ' '
823 read(funit,'(a)',end=40, err=1001) buffer
824 call despace(buffer)
826 ! Is this line a comment?
827 if (buffer(1:1) == '#') then
828 ! Do nothing.
830 else
831 have_specification = .true.
833 ! If only part of the line is a comment, just turn the comment into spaces
834 eos = index(buffer,'#')
835 if (eos /= 0) buffer(eos:BUFSIZE) = ' '
837 do while (have_specification) !{
839 ! If this line has no semicolon, it may contain a single specification,
840 ! so we set have_specification = .false. to prevent the line from being
841 ! processed again and pretend that the last character was a semicolon
842 eos = index(buffer,';')
843 if (eos == 0) then
844 have_specification = .false.
845 eos = BUFSIZE
846 end if
848 i = index(buffer(1:eos-1),'=')
850 if (i /= 0) then !{
852 if (index('projection',trim(buffer(1:i-1))) /= 0) then
853 if (index('lambert',trim(buffer(i+1:eos-1))) /= 0) then
854 is_proj(idx) = .true.
855 source_proj(idx) = PROJ_LC
856 else if (index('polar_wgs84',trim(buffer(i+1:eos-1))) /= 0 .and. &
857 len_trim('polar_wgs84') == len_trim(buffer(i+1:eos-1))) then
858 is_proj(idx) = .true.
859 source_proj(idx) = PROJ_PS_WGS84
860 else if (index('albers_nad83',trim(buffer(i+1:eos-1))) /= 0 .and. &
861 len_trim('albers_nad83') == len_trim(buffer(i+1:eos-1))) then
862 is_proj(idx) = .true.
863 source_proj(idx) = PROJ_ALBERS_NAD83
864 else if (index('polar',trim(buffer(i+1:eos-1))) /= 0 .and. &
865 len_trim('polar') == len_trim(buffer(i+1:eos-1))) then
866 is_proj(idx) = .true.
867 source_proj(idx) = PROJ_PS
868 else if (index('mercator',trim(buffer(i+1:eos-1))) /= 0) then
869 is_proj(idx) = .true.
870 source_proj(idx) = PROJ_MERC
871 else if (index('regular_ll',trim(buffer(i+1:eos-1))) /= 0) then
872 is_proj(idx) = .true.
873 source_proj(idx) = PROJ_LATLON
874 end if
876 else if (index('type',trim(buffer(1:i-1))) /= 0) then
877 if (index('continuous',trim(buffer(i+1:eos-1))) /= 0) then
878 is_fieldtype(idx) = .true.
879 source_fieldtype(idx) = CONTINUOUS
880 else if (index('categorical',trim(buffer(i+1:eos-1))) /= 0) then
881 is_fieldtype(idx) = .true.
882 source_fieldtype(idx) = CATEGORICAL
883 end if
885 else if (index('signed',trim(buffer(1:i-1))) /= 0) then
886 if (index('yes',trim(buffer(i+1:eos-1))) /= 0) then
887 is_signed(idx) = .true.
888 else if (index('no',trim(buffer(i+1:eos-1))) /= 0) then
889 is_signed(idx) = .false.
890 end if
892 else if (index('units',trim(buffer(1:i-1))) /= 0) then
893 ispace = i+1
894 iquoted = 0
895 do while (((ispace < eos) .and. (buffer(ispace:ispace) /= ' ')) .or. (iquoted == 1))
896 if (buffer(ispace:ispace) == '"' .or. buffer(ispace:ispace) == '''') iquoted = mod(iquoted+1,2)
897 ispace = ispace + 1
898 end do
899 is_units(idx) = .true.
900 source_units(idx) = ' '
901 if (buffer(i+1:i+1) == '"' .or. buffer(i+1:i+1) == '''') i = i + 1
902 if (buffer(ispace-1:ispace-1) == '"' .or. buffer(ispace-1:ispace-1) == '''') ispace = ispace - 1
903 source_units(idx)(1:ispace-i) = buffer(i+1:ispace-1)
905 else if (index('description',trim(buffer(1:i-1))) /= 0) then
906 ispace = i+1
907 iquoted = 0
908 do while (((ispace < eos) .and. (buffer(ispace:ispace) /= ' ')) .or. (iquoted == 1))
909 if (buffer(ispace:ispace) == '"' .or. buffer(ispace:ispace) == '''') iquoted = mod(iquoted+1,2)
910 ispace = ispace + 1
911 end do
912 is_descr(idx) = .true.
913 source_descr(idx) = ' '
914 if (buffer(i+1:i+1) == '"' .or. buffer(i+1:i+1) == '''') i = i + 1
915 if (buffer(ispace-1:ispace-1) == '"' .or. buffer(ispace-1:ispace-1) == '''') ispace = ispace - 1
916 source_descr(idx)(1:ispace-i) = buffer(i+1:ispace-1)
918 else if (index('mminlu',trim(buffer(1:i-1))) /= 0) then
919 ispace = i+1
920 iquoted = 0
921 do while (((ispace < eos) .and. (buffer(ispace:ispace) /= ' ')) .or. (iquoted == 1))
922 if (buffer(ispace:ispace) == '"' .or. buffer(ispace:ispace) == '''') iquoted = mod(iquoted+1,2)
923 ispace = ispace + 1
924 end do
925 if (buffer(i+1:i+1) == '"' .or. buffer(i+1:i+1) == '''') i = i + 1
926 if (buffer(ispace-1:ispace-1) == '"' .or. buffer(ispace-1:ispace-1) == '''') ispace = ispace - 1
927 source_mminlu(1:ispace-i) = buffer(i+1:ispace-1)
929 else if (index('iswater',trim(buffer(1:i-1))) /= 0) then
930 read(buffer(i+1:eos-1),*) source_iswater
932 else if (index('islake',trim(buffer(1:i-1))) /= 0) then
933 read(buffer(i+1:eos-1),*) source_islake
935 else if (index('isice',trim(buffer(1:i-1))) /= 0) then
936 read(buffer(i+1:eos-1),*) source_isice
938 else if (index('isurban',trim(buffer(1:i-1))) /= 0) then
939 read(buffer(i+1:eos-1),*) source_isurban
941 else if (index('isoilwater',trim(buffer(1:i-1))) /= 0) then
942 read(buffer(i+1:eos-1),*) source_isoilwater
944 else if (index('dx',trim(buffer(1:i-1))) /= 0) then
945 is_dx(idx) = .true.
946 read(buffer(i+1:eos-1),*) source_dx(idx)
948 else if (index('dy',trim(buffer(1:i-1))) /= 0) then
949 is_dy(idx) = .true.
950 read(buffer(i+1:eos-1),*) source_dy(idx)
952 else if (index('known_x',trim(buffer(1:i-1))) /= 0) then
953 is_known_x(idx) = .true.
954 read(buffer(i+1:eos-1),*) source_known_x(idx)
956 else if (index('known_y',trim(buffer(1:i-1))) /= 0) then
957 is_known_y(idx) = .true.
958 read(buffer(i+1:eos-1),*) source_known_y(idx)
960 else if (index('known_lat',trim(buffer(1:i-1))) /= 0) then
961 is_known_lat(idx) = .true.
962 read(buffer(i+1:eos-1),*) source_known_lat(idx)
964 else if (index('known_lon',trim(buffer(1:i-1))) /= 0) then
965 is_known_lon(idx) = .true.
966 read(buffer(i+1:eos-1),*) source_known_lon(idx)
968 else if (index('stdlon',trim(buffer(1:i-1))) /= 0) then
969 is_stdlon(idx) = .true.
970 read(buffer(i+1:eos-1),*) source_stdlon(idx)
972 else if (index('truelat1',trim(buffer(1:i-1))) /= 0) then
973 is_truelat1(idx) = .true.
974 read(buffer(i+1:eos-1),*) source_truelat1(idx)
976 else if (index('truelat2',trim(buffer(1:i-1))) /= 0) then
977 is_truelat2(idx) = .true.
978 read(buffer(i+1:eos-1),*) source_truelat2(idx)
980 else if (index('wordsize',trim(buffer(1:i-1))) /= 0) then
981 is_wordsize(idx) = .true.
982 read(buffer(i+1:eos-1),'(i10)') source_wordsize(idx)
984 else if (index('endian',trim(buffer(1:i-1))) /= 0) then
985 if (index('big',trim(buffer(i+1:eos-1))) /= 0) then
986 is_endian(idx) = .true.
987 source_endian(idx) = BIG_ENDIAN
988 else if (index('little',trim(buffer(i+1:eos-1))) /= 0) then
989 is_endian(idx) = .true.
990 source_endian(idx) = LITTLE_ENDIAN
991 else
992 call mprintf(.true.,WARN,'Invalid value for keyword ''endian'' '// &
993 'specified in index file. BIG_ENDIAN will be used.')
994 end if
996 else if (index('row_order',trim(buffer(1:i-1))) /= 0) then
997 if (index('bottom_top',trim(buffer(i+1:eos-1))) /= 0) then
998 is_row_order(idx) = .true.
999 source_row_order(idx) = BOTTOM_TOP
1000 else if (index('top_bottom',trim(buffer(i+1:eos-1))) /= 0) then
1001 is_row_order(idx) = .true.
1002 source_row_order(idx) = TOP_BOTTOM
1003 end if
1005 else if (index('tile_x',trim(buffer(1:i-1))) /= 0) then
1006 is_tile_x(idx) = .true.
1007 read(buffer(i+1:eos-1),'(i10)') source_tile_x(idx)
1009 else if (index('tile_y',trim(buffer(1:i-1))) /= 0) then
1010 is_tile_y(idx) = .true.
1011 read(buffer(i+1:eos-1),'(i10)') source_tile_y(idx)
1013 else if (index('tile_z',trim(buffer(1:i-1))) /= 0) then
1014 is_tile_z(idx) = .true.
1015 read(buffer(i+1:eos-1),'(i10)') source_tile_z(idx)
1017 else if (index('tile_z_start',trim(buffer(1:i-1))) /= 0) then
1018 is_tile_z_start(idx) = .true.
1019 read(buffer(i+1:eos-1),'(i10)') source_tile_z_start(idx)
1021 else if (index('tile_z_end',trim(buffer(1:i-1))) /= 0) then
1022 is_tile_z_end(idx) = .true.
1023 read(buffer(i+1:eos-1),'(i10)') source_tile_z_end(idx)
1025 else if (index('category_min',trim(buffer(1:i-1))) /= 0) then
1026 is_category_min(idx) = .true.
1027 read(buffer(i+1:eos-1),'(i10)') source_category_min(idx)
1029 else if (index('category_max',trim(buffer(1:i-1))) /= 0) then
1030 is_category_max(idx) = .true.
1031 read(buffer(i+1:eos-1),'(i10)') source_category_max(idx)
1033 else if (index('tile_bdr',trim(buffer(1:i-1))) /= 0) then
1034 is_tile_bdr(idx) = .true.
1035 read(buffer(i+1:eos-1),'(i10)') source_tile_bdr(idx)
1037 else if (index('missing_value',trim(buffer(1:i-1))) /= 0) then
1038 is_missing_value(idx) = .true.
1039 read(buffer(i+1:eos-1),*) source_missing_value(idx)
1041 else if (index('scale_factor',trim(buffer(1:i-1))) /= 0) then
1042 is_scale_factor(idx) = .true.
1043 read(buffer(i+1:eos-1),*) source_scale_factor(idx)
1045 else if (index('geotiff',trim(buffer(1:i-1))) /= 0) then
1046 is_geotiff(idx) = .true.
1047 read(buffer(i+1:eos-1),*) source_geotiff_file(idx)
1049 else
1050 call mprintf(.true., WARN, 'In %s, unrecognized option %s in entry %i.', &
1051 s1=trim(test_fname), s2=buffer(1:i-1), i1=i)
1053 end if
1055 end if !} index(buffer(1:eos-1),'=') /= 0
1057 buffer = buffer(eos+1:BUFSIZE)
1058 end do ! while eos /= 0 }
1060 end if !} index(buffer, '=====') /= 0
1062 go to 30
1064 40 close(funit)
1067 #ifdef _HAS_GEOTIFF
1068 if(is_geotiff(idx)) then
1069 if(source_geotiff_file(idx)(1:1) .ne. '/')then
1070 source_geotiff_file(idx)=TRIM(source_path(idx))//TRIM(source_geotiff_file(idx))
1071 endif
1072 call open_geotiff(source_geotiff_file(idx))
1073 call merge_geotiff_header(source_geotiff_file(idx),is_proj(idx),source_proj(idx),&
1074 is_fieldtype(idx),source_fieldtype(idx), &
1075 is_units(idx),source_units(idx), &
1076 is_descr(idx),source_descr(idx),is_dx(idx), &
1077 source_dx(idx),is_dy(idx),source_dy(idx), &
1078 is_known_x(idx),source_known_x(idx), &
1079 is_known_y(idx),source_known_y(idx), &
1080 is_known_lat(idx),source_known_lat(idx), &
1081 is_known_lon(idx),source_known_lon(idx), &
1082 is_stdlon(idx),source_stdlon(idx), &
1083 is_truelat1(idx),source_truelat1(idx), &
1084 is_truelat2(idx),source_truelat2(idx), &
1085 is_row_order(idx),source_row_order(idx), &
1086 is_tile_x(idx),source_tile_x(idx), &
1087 is_tile_y(idx),source_tile_y(idx), &
1088 is_tile_z(idx),source_tile_z(idx), &
1089 is_tile_z_start(idx),source_tile_z_start(idx), &
1090 is_tile_z_end(idx), source_tile_z_end(idx), &
1091 is_category_min(idx),source_category_min(idx), &
1092 is_category_max(idx),source_category_max(idx), &
1093 is_tile_bdr(idx),source_tile_bdr(idx), &
1094 is_missing_value(idx),source_missing_value(idx))
1095 is_wordsize(idx)=.true.
1096 source_wordsize(idx)=4
1097 is_scale_factor(idx)=.true.
1098 source_scale_factor(idx)=1.
1099 call mprintf(.true.,INFORM,'For geotiff file, %s, using the following parameters:', &
1100 s1=TRIM(source_geotiff_file(idx)))
1101 if(is_descr(idx)) &
1102 call mprintf(.true.,INFORM,'description=%s',s1=TRIM(source_descr(idx)))
1103 if(is_units(idx)) &
1104 call mprintf(.true.,INFORM,'units=%s',s1=TRIM(source_units(idx)))
1105 if(is_proj(idx)) &
1106 call mprintf(.true.,INFORM,'proj=%i',i1=source_proj(idx))
1107 if(is_stdlon(idx)) &
1108 call mprintf(.true.,INFORM,'stdlon=%f',f1=source_stdlon(idx))
1109 if(is_truelat1(idx)) &
1110 call mprintf(.true.,INFORM,'truelat1=%f',f1=source_truelat1(idx))
1111 if(is_truelat2(idx)) &
1112 call mprintf(.true.,INFORM,'truelat2=%f',f1=source_truelat2(idx))
1113 if(is_dx(idx)) &
1114 call mprintf(.true.,INFORM,'dx=%f',f1=source_dx(idx))
1115 if(is_dy(idx)) &
1116 call mprintf(.true.,INFORM,'dy=%f',f1=source_dy(idx))
1117 if(is_known_x(idx)) &
1118 call mprintf(.true.,INFORM,'known_x=%f',f1=source_known_x(idx))
1119 if(is_known_y(idx)) &
1120 call mprintf(.true.,INFORM,'known_y=%f',f1=source_known_y(idx))
1121 if(is_known_lon(idx)) &
1122 call mprintf(.true.,INFORM,'known_lon=%f',f1=source_known_lon(idx))
1123 if(is_known_lat(idx)) &
1124 call mprintf(.true.,INFORM,'known_lat=%f',f1=source_known_lat(idx))
1125 if(is_fieldtype(idx)) &
1126 call mprintf(.true.,INFORM,'fieldtype=%i',i1=source_fieldtype(idx))
1127 if(is_category_min(idx)) &
1128 call mprintf(.true.,INFORM,'category_min=%i',i1=source_category_min(idx))
1129 if(is_category_max(idx)) &
1130 call mprintf(.true.,INFORM,'category_max=%i',i1=source_category_max(idx))
1131 if(is_row_order(idx)) &
1132 call mprintf(.true.,INFORM,'row_order=%i',i1=source_row_order(idx))
1133 if(is_tile_x(idx)) &
1134 call mprintf(.true.,INFORM,'tile_x=%i',i1=source_tile_x(idx))
1135 if(is_tile_y(idx)) &
1136 call mprintf(.true.,INFORM,'tile_y=%i',i1=source_tile_y(idx))
1137 if(is_tile_z(idx)) &
1138 call mprintf(.true.,INFORM,'tile_z=%i',i1=source_tile_z(idx))
1139 if(is_tile_bdr(idx)) &
1140 call mprintf(.true.,INFORM,'tile_bdr=%i',i1=source_tile_bdr(idx))
1141 if(is_missing_value(idx)) &
1142 call mprintf(.true.,INFORM,'missing_value=%f',f1=source_missing_value(idx))
1143 endif
1144 #else
1145 if(is_geotiff(idx)) then
1146 call mprintf(.true.,ERROR,'A geotiff file has been specified in %s, '// &
1147 'but geogrid has not been compiled with geotiff '// &
1148 'support.',s1=trim(test_fname))
1149 endif
1150 #endif
1152 end do
1154 return
1156 1000 call mprintf(.true.,ERROR,'Could not open %s', s1=trim(test_fname))
1157 1001 call mprintf(.true.,ERROR,'Encountered error while reading %s', s1=trim(test_fname))
1159 end subroutine get_source_params
1162 !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
1163 ! Name: datalist_destroy()
1164 !
1165 ! Purpose: This routine deallocates any memory that was allocated by the
1166 ! get_datalist() subroutine.
1167 !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
1168 subroutine datalist_destroy()
1169 #ifdef _HAVE_GEOTIFF
1170 use geotiff_module, only : destroy_all
1171 #endif
1172 implicit none
1174 ! Local variables
1175 integer :: i
1177 #ifdef _HAVE_GEOTIFF
1178 call destroy_all()
1179 #endif
1181 if (associated(source_wordsize)) then
1182 deallocate(source_wordsize)
1183 deallocate(source_endian)
1184 deallocate(source_fieldtype)
1185 deallocate(source_dest_fieldtype)
1186 deallocate(source_proj)
1187 deallocate(source_priority)
1188 deallocate(source_dx)
1189 deallocate(source_dy)
1190 deallocate(source_known_x)
1191 deallocate(source_known_y)
1192 deallocate(source_known_lat)
1193 deallocate(source_known_lon)
1194 deallocate(source_truelat1)
1195 deallocate(source_truelat2)
1196 deallocate(source_stdlon)
1197 deallocate(source_fieldname)
1198 deallocate(source_path)
1199 deallocate(source_interp_string)
1200 deallocate(source_tile_x)
1201 deallocate(source_tile_y)
1202 deallocate(source_tile_z)
1203 deallocate(source_tile_z_start)
1204 deallocate(source_tile_z_end)
1205 deallocate(source_tile_bdr)
1206 deallocate(source_category_min)
1207 deallocate(source_category_max)
1208 deallocate(source_masked)
1209 deallocate(source_output_stagger)
1210 deallocate(source_row_order)
1211 deallocate(source_dominant_category)
1212 deallocate(source_dominant_only)
1213 deallocate(source_dfdx)
1214 deallocate(source_dfdy)
1215 deallocate(source_scale_factor)
1216 deallocate(source_z_dim_name)
1217 deallocate(source_smooth_option)
1218 deallocate(source_smooth_passes)
1219 deallocate(source_units)
1220 deallocate(source_descr)
1221 deallocate(source_missing_value)
1222 deallocate(source_fill_missing)
1223 do i=1,num_entries
1224 call list_destroy(source_res_path(i))
1225 call list_destroy(source_interp_option(i))
1226 call list_destroy(source_landmask_land(i))
1227 call list_destroy(source_landmask_water(i))
1228 end do
1229 deallocate(source_res_path)
1230 deallocate(source_interp_option)
1231 deallocate(source_landmask_land)
1232 deallocate(source_landmask_water)
1234 deallocate(is_wordsize)
1235 deallocate(is_endian)
1236 deallocate(is_fieldtype)
1237 deallocate(is_dest_fieldtype)
1238 deallocate(is_proj)
1239 deallocate(is_priority)
1240 deallocate(is_dx)
1241 deallocate(is_dy)
1242 deallocate(is_known_x)
1243 deallocate(is_known_y)
1244 deallocate(is_known_lat)
1245 deallocate(is_known_lon)
1246 deallocate(is_truelat1)
1247 deallocate(is_truelat2)
1248 deallocate(is_stdlon)
1249 deallocate(is_fieldname)
1250 deallocate(is_path)
1251 deallocate(is_tile_x)
1252 deallocate(is_tile_y)
1253 deallocate(is_tile_z)
1254 deallocate(is_tile_z_start)
1255 deallocate(is_tile_z_end)
1256 deallocate(is_tile_bdr)
1257 deallocate(is_category_min)
1258 deallocate(is_category_max)
1259 deallocate(is_masked)
1260 deallocate(is_halt_missing)
1261 deallocate(is_output_stagger)
1262 deallocate(is_row_order)
1263 deallocate(is_dominant_category)
1264 deallocate(is_dominant_only)
1265 deallocate(is_dfdx)
1266 deallocate(is_dfdy)
1267 deallocate(is_scale_factor)
1268 deallocate(is_z_dim_name)
1269 deallocate(is_smooth_option)
1270 deallocate(is_smooth_passes)
1271 deallocate(is_signed)
1272 deallocate(is_units)
1273 deallocate(is_descr)
1274 deallocate(is_missing_value)
1275 deallocate(is_fill_missing)
1276 deallocate(is_subgrid)
1277 end if
1279 call hash_destroy(bad_files)
1281 end subroutine datalist_destroy
1284 !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
1285 ! Name: reset_next_field
1286 !
1287 ! Purpose: To reset the pointer to the next field in the list of fields
1288 ! specified by the user.
1289 !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
1290 subroutine reset_next_field()
1292 implicit none
1294 next_field = 1
1296 end subroutine reset_next_field
1299 !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
1300 ! Name: get_next_fieldname
1301 !
1302 ! Purpose: Calling this routine results in field_name being set to the name of
1303 ! the field currently pointed to. If istatus /= 0 upon return, an error
1304 ! occurred, and the value of field_name is undefined.
1305 !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
1306 subroutine get_next_fieldname(field_name, istatus)
1308 implicit none
1310 ! Arguments
1311 integer, intent(out) :: istatus
1312 character (len=128), intent(out) :: field_name
1314 istatus = 1
1316 if (next_field <= num_entries) then
1318 field_name = source_fieldname(next_field)
1319 next_field = next_field + 1
1320 istatus = 0
1322 end if
1324 end subroutine get_next_fieldname
1327 !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
1328 ! Name: get_next_output_fieldname
1329 !
1330 ! Purpose:
1331 !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
1332 recursive subroutine get_next_output_fieldname(nest_num, field_name, ndims, &
1333 min_cat, max_cat, &
1334 istagger, memorder, dimnames, units, &
1335 description, sr_x, sr_y, istatus)
1337 use gridinfo_module
1339 implicit none
1341 #include "wrf_io_flags.h"
1343 ! Arguments
1344 integer, intent(in) :: nest_num
1345 integer, intent(out) :: istatus, ndims, istagger, min_cat, max_cat
1346 integer, intent(out) :: sr_x, sr_y
1347 character (len=128), intent(out) :: memorder, field_name, units, description
1348 character (len=128), dimension(3), intent(out) :: dimnames
1350 ! Local variables
1351 integer :: temp_fieldtype
1352 integer, dimension(MAX_LANDMASK_CATEGORIES) :: landmask
1353 logical :: is_water_mask, is_dom_only
1354 character (len=128) :: domcat_name, dfdx_name, dfdy_name
1355 character (len=256) :: temphash
1357 istatus = 1
1359 if (output_field_state == RETURN_LANDMASK) then
1360 call hash_init(duplicate_fnames)
1361 call get_landmask_field(geog_data_res(nest_num), field_name, is_water_mask, landmask, istatus)
1362 if (istatus == 0) then
1363 temphash(129:256) = ' '
1364 temphash(1:128) = field_name(1:128)
1365 call hash_insert(duplicate_fnames, temphash)
1366 call get_domcategory_name(field_name, domcat_name, is_dom_only, istatus)
1367 ! We will only save the dominant category
1368 if (is_dom_only .and. (istatus == 0)) then
1369 output_field_state = RETURN_DOMCAT_LM
1370 call get_next_output_fieldname(nest_num, field_name, ndims, &
1371 min_cat, max_cat, istagger, &
1372 memorder, dimnames, units, description, &
1373 sr_x, sr_y, istatus)
1374 return
1375 else
1376 ndims = 2
1377 min_cat = 1
1378 max_cat = 1
1379 temp_fieldtype = iget_fieldtype(field_name, istatus)
1380 if (istatus == 0) then
1381 if (temp_fieldtype == CONTINUOUS) then
1382 call get_max_levels(field_name, min_cat, max_cat, istatus)
1383 else if (temp_fieldtype == CATEGORICAL) then
1384 call get_max_categories(field_name, min_cat, max_cat, istatus)
1385 end if
1386 if (max_cat - min_cat > 0) ndims = 3
1387 end if
1388 call get_output_stagger(field_name, istagger, istatus)
1389 if (istagger == M) then
1390 dimnames(1) = 'west_east'
1391 dimnames(2) = 'south_north'
1392 else if (istagger == U) then
1393 dimnames(1) = 'west_east_stag'
1394 dimnames(2) = 'south_north'
1395 else if (istagger == V) then
1396 dimnames(1) = 'west_east'
1397 dimnames(2) = 'south_north_stag'
1398 else if (istagger == HH) then
1399 dimnames(1) = 'west_east'
1400 dimnames(2) = 'south_north'
1401 else if (istagger == VV) then
1402 dimnames(1) = 'west_east'
1403 dimnames(2) = 'south_north'
1404 end if
1405 if (ndims == 2) then
1406 memorder = 'XY '
1407 dimnames(3) = ' '
1408 else if (ndims == 3) then
1409 memorder = 'XYZ'
1410 call get_z_dim_name(field_name, dimnames(3), istatus)
1411 istatus = 0
1412 else
1413 memorder = ' '
1414 dimnames(3) = ' '
1415 end if
1416 call get_subgrid_dim_name(nest_num, field_name, dimnames(1:2), sr_x, sr_y, istatus)
1417 call get_source_units(field_name, 1, units, istatus)
1418 if (istatus /= 0) units = '-'
1419 call get_source_descr(field_name, 1, description, istatus)
1420 if (istatus /= 0) description = '-'
1421 istatus = 0
1422 output_field_state = RETURN_DOMCAT_LM
1423 end if
1424 else
1425 output_field_state = RETURN_FIELDNAME
1426 call get_next_output_fieldname(nest_num, field_name, ndims, &
1427 min_cat, max_cat, istagger, &
1428 memorder, dimnames, units, description, &
1429 sr_x, sr_y, istatus)
1430 return
1431 end if
1433 else if (output_field_state == RETURN_FIELDNAME) then
1434 call get_next_fieldname(field_name, istatus)
1435 temphash(129:256) = ' '
1436 temphash(1:128) = field_name(1:128)
1437 if (istatus == 0 .and. (.not. hash_search(duplicate_fnames, temphash))) then
1438 call hash_insert(duplicate_fnames, temphash)
1439 call get_domcategory_name(field_name, domcat_name, is_dom_only, istatus)
1440 ! We will only save the dominant category
1441 if (is_dom_only .and. (istatus == 0)) then
1442 output_field_state = RETURN_DOMCAT
1443 call get_next_output_fieldname(nest_num, field_name, ndims, &
1444 min_cat, max_cat, istagger, &
1445 memorder, dimnames, units, description, &
1446 sr_x, sr_y, istatus)
1447 return
1449 ! Return the fractional field
1450 else
1451 ndims = 2
1452 min_cat = 1
1453 max_cat = 1
1454 temp_fieldtype = iget_fieldtype(field_name, istatus)
1455 if (istatus == 0) then
1456 if (temp_fieldtype == CONTINUOUS) then
1457 call get_max_levels(field_name, min_cat, max_cat, istatus)
1458 else if (temp_fieldtype == CATEGORICAL) then
1459 call get_max_categories(field_name, min_cat, max_cat, istatus)
1460 end if
1461 if (max_cat - min_cat > 0) ndims = 3
1462 end if
1463 call get_output_stagger(field_name, istagger, istatus)
1464 if (istagger == M) then
1465 dimnames(1) = 'west_east'
1466 dimnames(2) = 'south_north'
1467 else if (istagger == U) then
1468 dimnames(1) = 'west_east_stag'
1469 dimnames(2) = 'south_north'
1470 else if (istagger == V) then
1471 dimnames(1) = 'west_east'
1472 dimnames(2) = 'south_north_stag'
1473 else if (istagger == HH) then
1474 dimnames(1) = 'west_east'
1475 dimnames(2) = 'south_north'
1476 else if (istagger == VV) then
1477 dimnames(1) = 'west_east'
1478 dimnames(2) = 'south_north'
1479 end if
1480 if (ndims == 2) then
1481 memorder = 'XY '
1482 dimnames(3) = ' '
1483 else if (ndims == 3) then
1484 memorder = 'XYZ'
1485 call get_z_dim_name(field_name, dimnames(3), istatus)
1486 istatus = 0
1487 else
1488 memorder = ' '
1489 dimnames(3) = ' '
1490 end if
1491 call get_subgrid_dim_name(nest_num, field_name, dimnames(1:2), sr_x, sr_y, istatus)
1492 call get_source_units(field_name, 1, units, istatus)
1493 if (istatus /= 0) units = '-'
1494 call get_source_descr(field_name, 1, description, istatus)
1495 if (istatus /= 0) description = '-'
1496 istatus = 0
1497 output_field_state = RETURN_DOMCAT
1498 end if
1499 else if (istatus /= 0) then
1500 output_field_state = RETURN_LANDMASK
1501 call hash_destroy(duplicate_fnames)
1502 return
1503 else if (hash_search(duplicate_fnames, temphash)) then
1504 call get_next_output_fieldname(nest_num, field_name, ndims, &
1505 min_cat, max_cat, istagger, &
1506 memorder, dimnames, units, description, &
1507 sr_x, sr_y, istatus)
1508 return
1509 end if
1511 else if (output_field_state == RETURN_DOMCAT .or. &
1512 output_field_state == RETURN_DOMCAT_LM ) then
1513 if (output_field_state == RETURN_DOMCAT) then
1514 next_field = next_field - 1
1515 call get_next_fieldname(field_name, istatus)
1516 else
1517 call get_landmask_field(geog_data_res(nest_num), field_name, is_water_mask, landmask, istatus)
1518 end if
1519 if (istatus == 0) then
1520 call get_domcategory_name(field_name, domcat_name, is_dom_only, istatus)
1521 if (istatus == 0) then
1522 ndims = 2
1523 min_cat = 1
1524 max_cat = 1
1525 call get_output_stagger(field_name, istagger, istatus)
1526 if (istagger == M) then
1527 dimnames(1) = 'west_east'
1528 dimnames(2) = 'south_north'
1529 else if (istagger == U) then
1530 dimnames(1) = 'west_east_stag'
1531 dimnames(2) = 'south_north'
1532 else if (istagger == V) then
1533 dimnames(1) = 'west_east'
1534 dimnames(2) = 'south_north_stag'
1535 else if (istagger == HH) then
1536 dimnames(1) = 'west_east'
1537 dimnames(2) = 'south_north'
1538 else if (istagger == VV) then
1539 dimnames(1) = 'west_east'
1540 dimnames(2) = 'south_north'
1541 end if
1542 dimnames(3) = ' '
1543 memorder = 'XY '
1545 call get_subgrid_dim_name(nest_num, field_name, dimnames(1:2), sr_x, sr_y, istatus)
1546 field_name = domcat_name
1547 units = 'category'
1548 description = 'Dominant category'
1549 if (output_field_state == RETURN_DOMCAT) then
1550 output_field_state = RETURN_DFDX
1551 else
1552 output_field_state = RETURN_DFDX_LM
1553 end if
1554 else
1555 if (output_field_state == RETURN_DOMCAT) then
1556 output_field_state = RETURN_DFDX
1557 else
1558 output_field_state = RETURN_DFDX_LM
1559 end if
1560 call get_next_output_fieldname(nest_num, field_name, ndims, &
1561 min_cat, max_cat, istagger, &
1562 memorder, dimnames, units, description, &
1563 sr_x, sr_y, istatus)
1564 return
1565 end if
1566 else
1567 call mprintf(.true., ERROR, 'get_next_output_fieldname(): In state DOMCAT, '// &
1568 'but no field name is found.')
1569 end if
1571 else if (output_field_state == RETURN_DFDX .or. &
1572 output_field_state == RETURN_DFDX_LM) then
1573 if (output_field_state == RETURN_DFDX) then
1574 next_field = next_field - 1
1575 call get_next_fieldname(field_name, istatus)
1576 else
1577 call get_landmask_field(geog_data_res(nest_num), field_name, is_water_mask, landmask, istatus)
1578 end if
1579 if (istatus == 0) then
1580 call get_dfdx_name(field_name, dfdx_name, istatus)
1581 if (istatus == 0) then
1582 ndims = 2
1583 min_cat = 1
1584 max_cat = 1
1585 temp_fieldtype = iget_fieldtype(field_name, istatus)
1586 if (istatus == 0) then
1587 if (temp_fieldtype == CONTINUOUS) then
1588 call get_max_levels(field_name, min_cat, max_cat, istatus)
1589 else if (temp_fieldtype == CATEGORICAL) then
1590 call get_max_categories(field_name, min_cat, max_cat, istatus)
1591 end if
1592 if (max_cat - min_cat > 0) ndims = 3
1593 end if
1594 call get_output_stagger(field_name, istagger, istatus)
1595 if (istagger == M) then
1596 dimnames(1) = 'west_east'
1597 dimnames(2) = 'south_north'
1598 else if (istagger == U) then
1599 dimnames(1) = 'west_east_stag'
1600 dimnames(2) = 'south_north'
1601 else if (istagger == V) then
1602 dimnames(1) = 'west_east'
1603 dimnames(2) = 'south_north_stag'
1604 else if (istagger == HH) then
1605 dimnames(1) = 'west_east'
1606 dimnames(2) = 'south_north'
1607 else if (istagger == VV) then
1608 dimnames(1) = 'west_east'
1609 dimnames(2) = 'south_north'
1610 end if
1611 if (ndims == 2) then
1612 memorder = 'XY '
1613 dimnames(3) = ' '
1614 else if (ndims == 3) then
1615 memorder = 'XYZ'
1616 call get_z_dim_name(field_name, dimnames(3), istatus)
1617 istatus = 0
1618 else
1619 memorder = ' '
1620 dimnames(3) = ' '
1621 end if
1622 field_name = dfdx_name
1623 units = '-'
1625 call get_subgrid_dim_name(nest_num, field_name, dimnames(1:2), sr_x, sr_y, istatus)
1626 description = 'df/dx'
1627 if (output_field_state == RETURN_DFDX) then
1628 output_field_state = RETURN_DFDY
1629 else
1630 output_field_state = RETURN_DFDY_LM
1631 end if
1632 else
1633 if (output_field_state == RETURN_DFDX) then
1634 output_field_state = RETURN_DFDY
1635 else
1636 output_field_state = RETURN_DFDY_LM
1637 end if
1638 call get_next_output_fieldname(nest_num, field_name, ndims, &
1639 min_cat, max_cat, istagger, &
1640 memorder, dimnames, units, description, &
1641 sr_x, sr_y, istatus)
1642 return
1643 end if
1644 else
1645 call mprintf(.true., ERROR, 'get_next_output_fieldname(): In state DFDX, '// &
1646 'but no field name is found.')
1647 end if
1649 else if (output_field_state == RETURN_DFDY .or. &
1650 output_field_state == RETURN_DFDY_LM) then
1651 if (output_field_state == RETURN_DFDY) then
1652 next_field = next_field - 1
1653 call get_next_fieldname(field_name, istatus)
1654 else
1655 call get_landmask_field(geog_data_res(nest_num), field_name, is_water_mask, landmask, istatus)
1656 end if
1657 if (istatus == 0) then
1658 call get_dfdy_name(field_name, dfdy_name, istatus)
1659 if (istatus == 0) then
1660 ndims = 2
1661 min_cat = 1
1662 max_cat = 1
1663 temp_fieldtype = iget_fieldtype(field_name, istatus)
1664 if (istatus == 0) then
1665 if (temp_fieldtype == CONTINUOUS) then
1666 call get_max_levels(field_name, min_cat, max_cat, istatus)
1667 else if (temp_fieldtype == CATEGORICAL) then
1668 call get_max_categories(field_name, min_cat, max_cat, istatus)
1669 end if
1670 if (max_cat - min_cat > 0) ndims = 3
1671 end if
1672 call get_output_stagger(field_name, istagger, istatus)
1673 if (istagger == M) then
1674 dimnames(1) = 'west_east'
1675 dimnames(2) = 'south_north'
1676 else if (istagger == U) then
1677 dimnames(1) = 'west_east_stag'
1678 dimnames(2) = 'south_north'
1679 else if (istagger == V) then
1680 dimnames(1) = 'west_east'
1681 dimnames(2) = 'south_north_stag'
1682 else if (istagger == HH) then
1683 dimnames(1) = 'west_east'
1684 dimnames(2) = 'south_north'
1685 else if (istagger == VV) then
1686 dimnames(1) = 'west_east'
1687 dimnames(2) = 'south_north'
1688 end if
1689 if (ndims == 2) then
1690 memorder = 'XY '
1691 dimnames(3) = ' '
1692 else if (ndims == 3) then
1693 memorder = 'XYZ'
1694 call get_z_dim_name(field_name, dimnames(3), istatus)
1695 istatus = 0
1696 else
1697 memorder = ' '
1698 dimnames(3) = ' '
1699 end if
1701 call get_subgrid_dim_name(nest_num, field_name, dimnames(1:2), sr_x, sr_y, istatus)
1702 field_name = dfdy_name
1703 units = '-'
1704 description = 'df/dy'
1705 output_field_state = RETURN_FIELDNAME
1706 else
1707 output_field_state = RETURN_FIELDNAME
1708 call get_next_output_fieldname(nest_num, field_name, ndims, &
1709 min_cat, max_cat, istagger, &
1710 memorder, dimnames, units, description, &
1711 sr_x, sr_y, istatus)
1712 return
1713 end if
1714 else
1715 call mprintf(.true., ERROR, 'get_next_output_fieldname(): In state DFDY, but no '// &
1716 'field name is found.')
1717 end if
1719 end if
1721 end subroutine get_next_output_fieldname
1724 !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
1725 ! Name: get_landmask_field
1726 !
1727 ! Purpose: To return the name of the field from which the landmask is to be
1728 ! computed. If no error occurs, is_water_mask is .true. if the landmask
1729 ! value specifies the value of water, and .false. if the landmask value
1730 ! specifies the value of land.
1731 !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
1732 subroutine get_landmask_field(res_string, landmask_name, is_water_mask, landmask, istatus)
1734 implicit none
1736 ! Arguments
1737 character (len=128), intent(in) :: res_string
1738 integer, dimension(:), intent(out) :: landmask
1739 integer, intent(out) :: istatus
1740 logical, intent(out) :: is_water_mask
1741 character (len=128), intent(out) :: landmask_name
1743 ! Local variables
1744 integer :: j
1745 integer :: ilen
1746 integer :: idx
1747 integer :: is, ie, sos, eos, comma
1748 character (len=128) :: temp_res, mask_cat_string
1750 istatus = 1
1752 do idx=1,num_entries
1754 if (list_length(source_landmask_land(idx)) > 0) then
1755 is = 1
1756 ie = index(res_string(is:128),'+') - 1
1757 if (ie <= 0) ie = 128
1758 temp_res = res_string(is:ie)
1759 do while (.not. list_search(source_landmask_land(idx), ckey=temp_res, cvalue=mask_cat_string) &
1760 .and. is <= 128)
1761 is = ie+2
1762 ie = is + index(res_string(is:128),'+') - 2
1763 if (ie - is <= 0) ie = 128
1764 temp_res = res_string(is:ie)
1765 end do
1767 if (is > 128) then
1768 temp_res = 'default'
1769 if (list_search(source_landmask_land(idx), ckey=temp_res, cvalue=mask_cat_string)) then
1770 is_water_mask = .false.
1771 landmask_name = source_fieldname(idx)
1772 istatus = 0
1773 end if
1774 else
1775 is_water_mask = .false.
1776 landmask_name = source_fieldname(idx)
1777 istatus = 0
1778 end if
1780 end if
1782 ! Note: The following cannot be an else-if, since different resolutions of data may
1783 ! specify, alternately, a land or a water mask, and in general we need to search
1784 ! both lists
1786 if (list_length(source_landmask_water(idx)) > 0) then
1787 is = 1
1788 ie = index(res_string(is:128),'+') - 1
1789 if (ie <= 0) ie = 128
1790 temp_res = res_string(is:ie)
1791 do while (.not. list_search(source_landmask_water(idx), ckey=temp_res, cvalue=mask_cat_string) &
1792 .and. is <= 128)
1793 is = ie+2
1794 ie = is + index(res_string(is:128),'+') - 2
1795 if (ie - is <= 0) ie = 128
1796 temp_res = res_string(is:ie)
1797 end do
1799 if (is > 128) then
1800 temp_res = 'default'
1801 if (list_search(source_landmask_water(idx), ckey=temp_res, cvalue=mask_cat_string)) then
1802 is_water_mask = .true.
1803 landmask_name = source_fieldname(idx)
1804 istatus = 0
1805 end if
1806 else
1807 is_water_mask = .true.
1808 landmask_name = source_fieldname(idx)
1809 istatus = 0
1810 end if
1811 end if
1813 if (istatus == 0) then
1814 j = 1
1815 sos = 0
1816 eos = 128
1817 comma = index(mask_cat_string(sos+1:eos-1),',')
1818 do while (comma > 0 .and. j < MAX_LANDMASK_CATEGORIES)
1819 read(mask_cat_string(sos+1:sos+comma-1),'(i10)') landmask(j)
1820 sos = sos + comma
1821 comma = index(mask_cat_string(sos+1:eos-1),',')
1822 j = j + 1
1823 end do
1824 read(mask_cat_string(sos+1:eos-1),'(i10)') landmask(j)
1825 j = j + 1
1826 if (j <= MAX_LANDMASK_CATEGORIES) then ! Terminate list with a flag value
1827 landmask(j) = INVALID
1828 end if
1829 exit
1830 end if
1832 end do
1834 end subroutine get_landmask_field
1837 !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
1838 ! Name: get_missing_value
1839 !
1840 ! Pupose: Return the value used in the source data to indicate missing data.
1841 !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
1842 subroutine get_missing_value(fieldnm, ilevel, rmissing, istatus)
1844 implicit none
1846 ! Arguments
1847 integer, intent(in) :: ilevel
1848 integer, intent(out) :: istatus
1849 real, intent(out) :: rmissing
1850 character (len=128), intent(in) :: fieldnm
1852 ! Local variables
1853 integer :: idx
1855 istatus = 1
1857 do idx=1,num_entries
1858 if ((index(source_fieldname(idx),trim(fieldnm)) /= 0) .and. &
1859 (len_trim(source_fieldname(idx)) == len_trim(fieldnm)) .and. &
1860 (source_priority(idx) == ilevel)) then
1862 if (is_missing_value(idx)) then
1863 rmissing = source_missing_value(idx)
1864 istatus = 0
1865 exit
1866 end if
1868 end if
1869 end do
1871 end subroutine get_missing_value
1874 !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
1875 ! Name: get_source_units
1876 !
1877 ! Pupose: Return a string giving the units of the specified source data.
1878 !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
1879 subroutine get_source_units(fieldnm, ilevel, cunits, istatus)
1881 implicit none
1883 ! Arguments
1884 integer, intent(in) :: ilevel
1885 integer, intent(out) :: istatus
1886 character (len=128), intent(in) :: fieldnm
1887 character (len=128), intent(out) :: cunits
1889 ! Local variables
1890 integer :: idx
1892 istatus = 1
1894 do idx=1,num_entries
1895 if ((index(source_fieldname(idx),trim(fieldnm)) /= 0) .and. &
1896 (len_trim(source_fieldname(idx)) == len_trim(fieldnm)) .and. &
1897 (source_priority(idx) == ilevel)) then
1899 if (is_units(idx)) then
1900 cunits = source_units(idx)
1901 istatus = 0
1902 exit
1903 end if
1905 end if
1906 end do
1908 end subroutine get_source_units
1911 !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
1912 ! Name: get_source_descr
1913 !
1914 ! Pupose: Return a string giving a description of the specified source data.
1915 !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
1916 subroutine get_source_descr(fieldnm, ilevel, descr, istatus)
1918 implicit none
1920 ! Arguments
1921 integer, intent(in) :: ilevel
1922 integer, intent(out) :: istatus
1923 character (len=128), intent(in) :: fieldnm
1924 character (len=128), intent(out) :: descr
1926 ! Local variables
1927 integer :: idx
1929 istatus = 1
1931 do idx=1,num_entries
1932 if ((index(source_fieldname(idx),trim(fieldnm)) /= 0) .and. &
1933 (len_trim(source_fieldname(idx)) == len_trim(fieldnm)) .and. &
1934 (source_priority(idx) == ilevel)) then
1936 if (is_units(idx)) then
1937 descr = source_descr(idx)
1938 istatus = 0
1939 exit
1940 end if
1942 end if
1943 end do
1945 end subroutine get_source_descr
1948 !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
1949 ! Name: get_missing_fill_value
1950 !
1951 ! Pupose: Return the value to fill missing points with.
1952 !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
1953 subroutine get_missing_fill_value(fieldnm, rmissing, istatus)
1955 implicit none
1957 ! Arguments
1958 integer, intent(out) :: istatus
1959 real, intent(out) :: rmissing
1960 character (len=128), intent(in) :: fieldnm
1962 ! Local variables
1963 integer :: idx
1965 istatus = 1
1967 do idx=1,num_entries
1968 if ((index(source_fieldname(idx),trim(fieldnm)) /= 0) .and. &
1969 (len_trim(source_fieldname(idx)) == len_trim(fieldnm)) ) then
1971 if (is_fill_missing(idx)) then
1972 rmissing = source_fill_missing(idx)
1973 istatus = 0
1974 exit
1975 end if
1977 end if
1978 end do
1980 end subroutine get_missing_fill_value
1983 !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
1984 ! Name: get_halt_on_missing
1985 !
1986 ! Pupose: Returns 1 if the program should halt upon encountering a missing
1987 ! value in the final output field, and 0 otherwise.
1988 !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
1989 subroutine get_halt_on_missing(fieldnm, halt, istatus)
1991 implicit none
1993 ! Arguments
1994 integer, intent(out) :: istatus
1995 logical, intent(out) :: halt
1996 character (len=128), intent(in) :: fieldnm
1998 ! Local variables
1999 integer :: idx
2001 istatus = 0
2002 halt = .false.
2004 do idx=1,num_entries
2005 if ((index(source_fieldname(idx),trim(fieldnm)) /= 0) .and. &
2006 (len_trim(source_fieldname(idx)) == len_trim(fieldnm)) ) then
2008 if (is_halt_missing(idx)) halt = .true.
2010 end if
2011 end do
2013 end subroutine get_halt_on_missing
2016 !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
2017 ! Name: get_masked_value
2018 !
2019 ! Pupose: If the field is to be masked by the landmask, returns 0 if the field
2020 ! is masked over water and 1 if the field is masked over land. If no mask is
2021 ! to be applied, -1 is returned.
2022 !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
2023 subroutine get_masked_value(fieldnm, ilevel, masked, istatus)
2025 implicit none
2027 ! Arguments
2028 integer, intent(in) :: ilevel
2029 integer, intent(out) :: istatus
2030 real, intent(out) :: masked
2031 character (len=128), intent(in) :: fieldnm
2033 ! Local variables
2034 integer :: idx
2036 istatus = 0
2037 masked = -1.
2039 do idx=1,num_entries
2040 if ((index(source_fieldname(idx),trim(fieldnm)) /= 0) .and. &
2041 (len_trim(source_fieldname(idx)) == len_trim(fieldnm)) .and. &
2042 (source_priority(idx) == ilevel)) then
2044 if (is_masked(idx)) then
2045 masked = source_masked(idx)
2046 exit
2047 end if
2049 end if
2050 end do
2052 end subroutine get_masked_value
2055 !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
2056 ! Name: get_max_levels
2057 !
2058 ! Purpose: Returns the number of levels for the field given by fieldnm.
2059 ! The number of levels will generally be specified for continuous fields,
2060 ! whereas min/max category will generally be specified for categorical ones.
2061 !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
2062 subroutine get_max_levels(fieldnm, min_level, max_level, istatus)
2064 implicit none
2066 ! Arguments
2067 integer, intent(out) :: min_level, max_level, istatus
2068 character (len=128), intent(in) :: fieldnm
2070 ! Local variables
2071 integer :: idx
2072 logical :: have_found_field
2074 have_found_field = .false.
2075 istatus = 1
2077 do idx=1,num_entries
2078 if ((index(source_fieldname(idx),trim(fieldnm)) /= 0) .and. &
2079 (len_trim(source_fieldname(idx)) == len_trim(fieldnm))) then
2081 if (is_dest_fieldtype(idx) .and. (source_dest_fieldtype(idx) /= CONTINUOUS)) then
2082 call mprintf(.true., WARN, 'In GEOGRID.TBL, destination field type for %s is '// &
2083 'not continuous and min/max levels specified.', s1=trim(fieldnm))
2084 end if
2085 if (.not. have_found_field) then
2086 if (is_tile_z_start(idx) .and. is_tile_z_end(idx)) then
2087 have_found_field = .true.
2088 istatus = 0
2089 min_level = source_tile_z_start(idx)
2090 max_level = source_tile_z_end(idx)
2091 else if (is_tile_z(idx)) then
2092 have_found_field = .true.
2093 istatus = 0
2094 min_level = 1
2095 max_level = source_tile_z(idx)
2096 end if
2098 if (.not. (is_tile_z_start(idx) .and. is_tile_z_end(idx))) then
2099 if (is_tile_z_start(idx) .or. is_tile_z_end(idx)) then
2100 call mprintf(.true., ERROR, 'In GEOGRID.TBL, only one of tile_z_start '// &
2101 'and tile_z_end specified for entry %i.',i1=idx)
2102 end if
2103 end if
2104 else
2105 if (is_tile_z_start(idx) .and. is_tile_z_end(idx)) then
2106 if (source_tile_z_start(idx) < min_level) min_level = source_tile_z_start(idx)
2107 if (source_tile_z_end(idx) > max_level) max_level = source_tile_z_end(idx)
2108 else if (is_tile_z(idx)) then
2109 if (min_level > 1) min_level = 1
2110 if (source_tile_z(idx) > max_level) max_level = source_tile_z(idx)
2111 end if
2113 if (.not. (is_tile_z_start(idx) .and. is_tile_z_end(idx))) then
2114 if (is_tile_z_start(idx) .or. is_tile_z_end(idx)) then
2115 call mprintf(.true., ERROR, 'In GEOGRID.TBL, only one of tile_z_start '// &
2116 'and tile_z_end specified for entry %i.',i1=idx)
2117 end if
2118 end if
2119 end if
2121 end if
2122 end do
2124 end subroutine get_max_levels
2127 !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
2128 ! Name: get_source_levels
2129 !
2130 ! Purpose: Return the min and max z-index for the source data for fieldname
2131 ! at a specified priority level (compared with the min/max level over
2132 ! all priority levels, as given by get_max_levels).
2133 !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
2134 subroutine get_source_levels(fieldnm, ilevel, min_level, max_level, istatus)
2136 implicit none
2138 ! Arguments
2139 integer, intent(in) :: ilevel
2140 integer, intent(out) :: min_level, max_level, istatus
2141 character (len=128), intent(in) :: fieldnm
2143 ! Local variables
2144 integer :: idx
2146 istatus = 1
2148 do idx=1,num_entries
2149 if ((index(source_fieldname(idx),trim(fieldnm)) /= 0) .and. &
2150 (len_trim(source_fieldname(idx)) == len_trim(fieldnm))) then
2151 if (ilevel == source_priority(idx)) then
2153 if (is_tile_z_start(idx) .and. is_tile_z_end(idx)) then
2154 istatus = 0
2155 min_level = source_tile_z_start(idx)
2156 max_level = source_tile_z_end(idx)
2157 else if (is_tile_z(idx)) then
2158 istatus = 0
2159 min_level = 1
2160 max_level = source_tile_z(idx)
2161 end if
2163 if (.not. (is_tile_z_start(idx) .and. is_tile_z_end(idx))) then
2164 if (is_tile_z_start(idx) .or. is_tile_z_end(idx)) then
2165 call mprintf(.true., ERROR, 'In GEOGRID.TBL, only one of tile_z_start '// &
2166 'and tile_z_end specified for entry %i.',i1=idx)
2167 end if
2168 end if
2170 end if
2171 end if
2172 end do
2174 end subroutine get_source_levels
2177 !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
2178 ! Name: get_max_categories
2179 !
2180 ! Purpose: Returns the minimum category and the maximum category for the field
2181 ! given by fieldnm.
2182 ! Min/max category will generally be specified for categorical fields,
2183 ! whereas the number of levels will generally be specified for continuous
2184 ! fields.
2185 !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
2186 subroutine get_max_categories(fieldnm, min_category, max_category, istatus)
2188 implicit none
2190 ! Arguments
2191 integer, intent(out) :: min_category, max_category, istatus
2192 character (len=128), intent(in) :: fieldnm
2194 ! Local variables
2195 integer :: idx
2196 logical :: have_found_field
2198 have_found_field = .false.
2199 istatus = 1
2201 do idx=1,num_entries
2202 if ((index(source_fieldname(idx),trim(fieldnm)) /= 0) .and. &
2203 (len_trim(source_fieldname(idx)) == len_trim(fieldnm))) then
2205 if (is_dest_fieldtype(idx) .and. (source_dest_fieldtype(idx) /= CATEGORICAL)) then
2206 call mprintf(.true., WARN, &
2207 'In GEOGRID.TBL, cannot get min/max categories for continuous '// &
2208 'field %s at entry %i. Perhaps the user has requested to '// &
2209 'perform a strange operation on the field.', s1=trim(fieldnm), i1=idx)
2210 end if
2211 if (.not. have_found_field) then
2212 if (is_category_min(idx) .and. is_category_max(idx)) then
2213 have_found_field = .true.
2214 istatus = 0
2215 min_category = source_category_min(idx)
2216 max_category = source_category_max(idx)
2217 else if (is_category_min(idx) .or. is_category_max(idx)) then
2218 call mprintf(.true., ERROR, 'In GEOGRID.TBL, only one of min_category and '// &
2219 'max_category specified for entry %i.',i1=idx)
2220 end if
2221 else
2222 if (is_category_min(idx) .and. is_category_max(idx)) then
2223 if (source_category_min(idx) < min_category) min_category = source_category_min(idx)
2224 if (source_category_max(idx) > max_category) max_category = source_category_max(idx)
2225 else if (is_category_min(idx) .or. is_category_max(idx)) then
2226 call mprintf(.true., ERROR, 'In GEOGRID.TBL, only one of min_category and '// &
2227 'max_category specified for entry %i.',i1=idx)
2228 end if
2229 end if
2231 end if
2232 end do
2234 end subroutine get_max_categories
2237 !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
2238 ! Name: get_source_categories
2239 !
2240 ! Purpose: Return the min and max category for the source data for fieldname
2241 ! at a specified priority level (compared with the min/max category over
2242 ! all priority levels, as given by get_max_categories).
2243 !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
2244 subroutine get_source_categories(fieldnm, ilevel, min_category, max_category, istatus)
2246 implicit none
2248 ! Arguments
2249 integer, intent(in) :: ilevel
2250 integer, intent(out) :: min_category, max_category, istatus
2251 character (len=128), intent(in) :: fieldnm
2253 ! Local variables
2254 integer :: idx
2256 istatus = 1
2258 do idx=1,num_entries
2259 if ((index(source_fieldname(idx),trim(fieldnm)) /= 0) .and. &
2260 (len_trim(source_fieldname(idx)) == len_trim(fieldnm))) then
2261 if (ilevel == source_priority(idx)) then
2263 if (is_category_min(idx) .and. is_category_max(idx)) then
2264 istatus = 0
2265 min_category = source_category_min(idx)
2266 max_category = source_category_max(idx)
2267 else if (is_category_min(idx) .or. is_category_max(idx)) then
2268 call mprintf(.true., ERROR, 'In GEOGRID.TBL, only one of min_category '// &
2269 'and max_category specified for entry %i.',i1=idx)
2270 end if
2272 end if
2273 end if
2274 end do
2276 end subroutine get_source_categories
2279 !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
2280 ! Name: get_domcategory_name
2281 !
2282 ! Purpose:
2283 !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
2284 subroutine get_domcategory_name(fieldnm, domcat_name, ldominant_only, istatus)
2286 implicit none
2288 ! Arguments
2289 integer, intent(out) :: istatus
2290 logical, intent(out) :: ldominant_only
2291 character (len=128), intent(in) :: fieldnm
2292 character (len=128), intent(out) :: domcat_name
2294 ! Local variables
2295 integer :: idx
2297 istatus = 1
2298 ldominant_only = .false.
2300 do idx=1,num_entries
2301 if ((index(source_fieldname(idx),trim(fieldnm)) /= 0) .and. &
2302 (len_trim(source_fieldname(idx)) == len_trim(fieldnm))) then
2304 if (is_dominant_category(idx)) then
2305 domcat_name = source_dominant_category(idx)
2306 istatus = 0
2307 if (is_dominant_only(idx)) then
2308 call mprintf(.true., WARN, 'In GEOGRID.TBL, both dominant_category and '// &
2309 'dominant_only are specified in entry %i. Using specification '// &
2310 'for dominant_category.',i1=idx)
2311 is_dominant_only(idx) = .false.
2312 end if
2313 exit
2315 else if (is_dominant_only(idx)) then
2316 domcat_name = source_dominant_only(idx)
2317 ldominant_only = .true.
2318 istatus = 0
2319 exit
2320 end if
2322 end if
2323 end do
2325 end subroutine get_domcategory_name
2328 !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
2329 ! Name: get_dfdx_name
2330 !
2331 ! Purpose:
2332 !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
2333 subroutine get_dfdx_name(fieldnm, dfdx_name, istatus)
2335 implicit none
2337 ! Arguments
2338 integer, intent(out) :: istatus
2339 character (len=128), intent(in) :: fieldnm
2340 character (len=128), intent(out) :: dfdx_name
2342 ! Local variables
2343 integer :: idx
2345 istatus = 1
2347 do idx=1,num_entries
2348 if ((index(source_fieldname(idx),trim(fieldnm)) /= 0) .and. &
2349 (len_trim(source_fieldname(idx)) == len_trim(fieldnm))) then
2351 if (is_dfdx(idx)) then
2352 dfdx_name = source_dfdx(idx)
2353 istatus = 0
2354 exit
2355 end if
2357 end if
2358 end do
2360 end subroutine get_dfdx_name
2363 !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
2364 ! Name: get_dfdy_name
2365 !
2366 ! Purpose:
2367 !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
2368 subroutine get_dfdy_name(fieldnm, dfdy_name, istatus)
2370 implicit none
2372 ! Arguments
2373 integer, intent(out) :: istatus
2374 character (len=128), intent(in) :: fieldnm
2375 character (len=128), intent(out) :: dfdy_name
2377 ! Local variables
2378 integer :: idx
2380 istatus = 1
2382 do idx=1,num_entries
2383 if ((index(source_fieldname(idx),trim(fieldnm)) /= 0) .and. &
2384 (len_trim(source_fieldname(idx)) == len_trim(fieldnm))) then
2386 if (is_dfdy(idx)) then
2387 dfdy_name = source_dfdy(idx)
2388 istatus = 0
2389 exit
2390 end if
2392 end if
2393 end do
2395 end subroutine get_dfdy_name
2398 !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
2399 ! Name: get_z_dim_name
2400 !
2401 ! Purpose:
2402 !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
2403 subroutine get_z_dim_name(fieldnm, z_dim, istatus)
2405 implicit none
2407 ! Arguments
2408 integer, intent(out) :: istatus
2409 character (len=128), intent(in) :: fieldnm
2410 character (len=128), intent(out) :: z_dim
2412 ! Local variables
2413 integer :: idx
2415 istatus = 1
2417 do idx=1,num_entries
2418 if ((index(source_fieldname(idx),trim(fieldnm)) /= 0) .and. &
2419 (len_trim(source_fieldname(idx)) == len_trim(fieldnm))) then
2420 if (is_z_dim_name(idx)) then
2421 z_dim = source_z_dim_name(idx)
2422 istatus = 0
2423 exit
2424 end if
2425 end if
2426 end do
2428 end subroutine get_z_dim_name
2431 !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
2432 ! Name: get_field_scale_factor
2433 !
2434 ! Purpose:
2435 !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
2436 subroutine get_field_scale_factor(fieldnm, ilevel, scale_factor, istatus)
2438 implicit none
2440 ! Arguments
2441 integer, intent(in) :: ilevel
2442 integer, intent(out) :: istatus
2443 real, intent(out) :: scale_factor
2444 character (len=128), intent(in) :: fieldnm
2446 ! Local variables
2447 integer :: idx
2449 istatus = 1
2451 do idx=1,num_entries
2452 if ((index(source_fieldname(idx),trim(fieldnm)) /= 0) .and. &
2453 (len_trim(source_fieldname(idx)) == len_trim(fieldnm)) .and. &
2454 (ilevel == source_priority(idx))) then
2456 if (is_scale_factor(idx)) then
2457 scale_factor = source_scale_factor(idx)
2458 istatus = 0
2459 end if
2461 end if
2462 end do
2464 end subroutine get_field_scale_factor
2467 !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
2468 ! Name: get_output_stagger
2469 !
2470 ! Pupose:
2471 !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
2472 subroutine get_output_stagger(fieldnm, istagger, istatus)
2474 use gridinfo_module
2476 implicit none
2478 ! Arguments
2479 integer, intent(out) :: istatus, istagger
2480 character (len=128), intent(in) :: fieldnm
2482 ! Local variables
2483 integer :: idx
2485 istatus = 1
2486 do idx=1,num_entries
2487 if ((index(source_fieldname(idx),trim(fieldnm)) /= 0) .and. &
2488 (len_trim(source_fieldname(idx)) == len_trim(fieldnm))) then
2490 if (is_output_stagger(idx)) then
2491 istatus = 0
2492 istagger = source_output_stagger(idx)
2493 exit
2494 else
2495 if (gridtype == 'C') then
2496 istatus = 0
2497 istagger = M
2498 exit
2499 else if (gridtype == 'E') then
2500 istatus = 0
2501 istagger = HH
2502 exit
2503 end if
2504 end if
2506 end if
2507 end do
2509 end subroutine get_output_stagger
2512 !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
2513 ! Name: get_subgrid_dim_name
2514 !
2515 ! Pupose:
2516 !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
2517 subroutine get_subgrid_dim_name(nest_num, field_name, dimnames, &
2518 sub_x, sub_y, istatus)
2520 use gridinfo_module
2522 implicit none
2523 integer, intent(in) :: nest_num
2524 integer, intent(out) :: sub_x, sub_y, istatus
2525 character(len=128), intent(in) :: field_name
2526 character(len=128), dimension(2), intent(inout) :: dimnames
2527 integer :: idx, nlen
2529 sub_x = 1
2530 sub_y = 1
2532 istatus = 0
2533 do idx=1,num_entries
2534 if ((index(source_fieldname(idx),trim(field_name)) /= 0) .and. &
2535 (len_trim(source_fieldname(idx)) == len_trim(field_name))) then
2536 if (is_subgrid(idx)) then
2537 istatus = 0
2538 if (is_output_stagger(idx)) then
2539 call mprintf(.true.,ERROR,'Cannot use subgrids on variables with staggered grids')
2540 end if
2541 dimnames(1) = trim(dimnames(1))//"_subgrid"
2542 dimnames(2) = trim(dimnames(2))//"_subgrid"
2543 sub_x = subgrid_ratio_x(nest_num)
2544 sub_y = subgrid_ratio_y(nest_num)
2545 end if
2546 end if
2547 end do
2549 end subroutine get_subgrid_dim_name
2552 !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
2553 ! Name: get_interp_option
2554 !
2555 ! Pupose:
2556 !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
2557 subroutine get_interp_option(fieldnm, ilevel, interp_opt, istatus)
2559 implicit none
2561 ! Arguments
2562 integer, intent(in) :: ilevel
2563 integer, intent(out) :: istatus
2564 character (len=128), intent(in) :: fieldnm
2565 character (len=128), intent(out) :: interp_opt
2567 ! Local variables
2568 integer :: idx
2570 istatus = 1
2571 do idx=1,num_entries
2572 if ((index(source_fieldname(idx),trim(fieldnm)) /= 0) .and. &
2573 (len_trim(source_fieldname(idx)) == len_trim(fieldnm))) then
2574 if (ilevel == source_priority(idx)) then
2576 interp_opt = source_interp_string(idx)
2577 istatus = 0
2578 exit
2580 end if
2581 end if
2582 end do
2584 end subroutine get_interp_option
2587 !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
2588 ! Name: get_gcel_threshold
2589 !
2590 ! Pupose:
2591 !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
2592 subroutine get_gcell_threshold(interp_opt, threshold, istatus)
2594 implicit none
2596 ! Arguments
2597 integer, intent(out) :: istatus
2598 real, intent(out) :: threshold
2599 character (len=128), intent(in) :: interp_opt
2601 ! Local variables
2602 integer :: i, p1, p2
2604 istatus = 1
2605 threshold = 1.0
2607 i = index(interp_opt,'average_gcell')
2608 if (i /= 0) then
2609 ! Move the "average_gcell" option to the beginning
2610 ! if (i /= 1) then
2611 ! p1 =
2612 ! end if
2614 ! Check for a threshold
2615 p1 = index(interp_opt(i:128),'(')
2616 p2 = index(interp_opt(i:128),')')
2617 if (p1 /= 0 .and. p2 /= 0) then
2618 read(interp_opt(p1+1:p2-1),*,err=1000) threshold
2619 else
2620 call mprintf(.true., WARN, 'Problem with specified threshold '// &
2621 'for average_gcell interp option. Setting threshold to 0.0.')
2622 threshold = 0.0
2623 end if
2624 end if
2625 istatus = 0
2627 return
2629 1000 call mprintf(.true., ERROR, 'Threshold option to average_gcell interpolator '// &
2630 'must be a real number, enclosed in parentheses immediately '// &
2631 'after keyword "average_gcell"')
2633 end subroutine get_gcell_threshold
2636 !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
2637 ! Name: get_smooth_option
2638 !
2639 ! Pupose:
2640 !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
2641 subroutine get_smooth_option(fieldnm, smooth_opt, smooth_passes, istatus)
2643 implicit none
2645 ! Arguments
2646 integer, intent(out) :: istatus, smooth_opt, smooth_passes
2647 character (len=128), intent(in) :: fieldnm
2649 ! Local variables
2650 integer :: idx
2652 istatus = 1
2654 do idx=1,num_entries
2655 if ((index(source_fieldname(idx),trim(fieldnm)) /= 0) .and. &
2656 (len_trim(source_fieldname(idx)) == len_trim(fieldnm))) then
2658 if (is_smooth_option(idx)) then
2659 istatus = 0
2660 smooth_opt = source_smooth_option(idx)
2662 if (is_smooth_passes(idx)) then
2663 smooth_passes = source_smooth_passes(idx)
2664 else
2665 smooth_passes = 1
2666 end if
2668 exit
2669 end if
2671 end if
2672 end do
2674 end subroutine get_smooth_option
2677 !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
2678 ! Name: iget_fieldtype
2679 !
2680 ! Purpose:
2681 !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
2682 function iget_fieldtype(fieldnm, istatus)
2684 implicit none
2686 ! Arguments
2687 integer, intent(out) :: istatus
2688 character (len=128), intent(in) :: fieldnm
2690 ! Local variables
2691 integer :: idx
2693 ! Return value
2694 integer :: iget_fieldtype
2696 istatus = 1
2698 do idx=1,num_entries
2699 if ((index(source_fieldname(idx),trim(fieldnm)) /= 0) .and. &
2700 (len_trim(source_fieldname(idx)) == len_trim(fieldnm))) then
2702 if (is_dest_fieldtype(idx)) then
2703 iget_fieldtype = source_dest_fieldtype(idx)
2704 istatus = 0
2705 exit
2706 end if
2708 end if
2709 end do
2711 end function iget_fieldtype
2714 !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
2715 ! Name: iget_source_fieldtype
2716 !
2717 ! Purpose: Given a resolution, in degrees, and the name of a field, this
2718 ! function returns the type (categorical, continuous, etc.) of the source
2719 ! data that will be used. This may, in general, depend on the field name
2720 ! and the resolution; for example, near 30 second resolution, land use data
2721 ! may come from a categorical field, whereas for lower resolutions, it may
2722 ! come from arrays of land use fractions for each category.
2723 !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
2724 function iget_source_fieldtype(fieldnm, ilevel, istatus)
2726 implicit none
2728 ! Arguments
2729 integer, intent(in) :: ilevel
2730 integer, intent(out) :: istatus
2731 character (len=128), intent(in) :: fieldnm
2733 ! Return value
2734 integer :: iget_source_fieldtype
2736 ! Local variables
2737 integer :: idx
2739 ! Find information about the source tiles for the specified fieldname
2740 istatus = 1
2741 do idx=1,num_entries
2742 if ((index(source_fieldname(idx),trim(fieldnm)) /= 0) .and. &
2743 (len_trim(source_fieldname(idx)) == len_trim(fieldnm))) then
2744 if (ilevel == source_priority(idx)) then
2745 istatus = 0
2746 iget_source_fieldtype = source_fieldtype(idx)
2747 exit
2748 end if
2749 end if
2750 end do
2752 end function iget_source_fieldtype
2755 !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
2756 ! Name: get_data_tile
2757 !
2758 ! Purpose:
2759 !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
2760 subroutine get_data_tile(xlat, xlon, ilevel, field_name, &
2761 file_name, array, start_x_dim, end_x_dim, start_y_dim, &
2762 end_y_dim, start_z_dim, end_z_dim, npts_bdr, &
2763 istatus)
2764 #ifdef _HAS_GEOTIFF
2765 use geotiff_module, only : read_geogrid_tile
2766 #endif
2767 implicit none
2769 ! Arguments
2770 integer, intent(in) :: ilevel
2771 integer, intent(out) :: istatus
2772 integer, intent(out) :: start_x_dim, end_x_dim, &
2773 start_y_dim, end_y_dim, &
2774 start_z_dim, end_z_dim, &
2775 npts_bdr
2776 real, intent(in) :: xlat, xlon ! Location that tile should contain
2777 real, pointer, dimension(:,:,:) :: array ! The array to be allocated by this routine
2778 character (len=128), intent(in) :: field_name
2779 character (len=256), intent(out) :: file_name
2781 ! Local variables
2782 integer :: j, k
2783 integer :: local_wordsize, local_endian, sign_convention, irow_order, strlen
2784 integer :: xdim,ydim,zdim
2785 real :: scalefac
2786 real, allocatable, dimension(:) :: temprow
2787 logical :: geotiff
2788 character(len=256) :: gtiff_file
2790 geotiff=is_geotiff_file(field_name,ilevel)
2791 call get_tile_fname(file_name, xlat, xlon, ilevel, field_name, istatus)
2792 if (.not.geotiff .and. index(file_name, 'OUTSIDE') /= 0) then
2793 istatus = 1
2794 return
2795 else if (hash_search(bad_files, file_name)) then
2796 istatus = 1
2797 return
2798 end if
2800 call get_tile_dimensions(xlat, xlon, start_x_dim, end_x_dim, start_y_dim, end_y_dim, &
2801 start_z_dim, end_z_dim, npts_bdr, local_wordsize, local_endian, &
2802 sign_convention, ilevel, field_name, istatus)
2804 xdim = (end_x_dim-start_x_dim+1)
2805 ydim = (end_y_dim-start_y_dim+1)
2806 zdim = (end_z_dim-start_z_dim+1)
2808 if (associated(array)) deallocate(array)
2809 allocate(array(xdim,ydim,zdim))
2811 call get_row_order(field_name, ilevel, irow_order, istatus)
2812 if (istatus /= 0) irow_order = BOTTOM_TOP
2814 call s_len(file_name,strlen)
2816 scalefac = 1.0
2818 if(geotiff) then
2819 #ifdef _HAS_GEOTIFF
2820 call get_geotiff_filename(field_name,ilevel,gtiff_file,istatus)
2821 if(istatus.ne.0)then
2822 call mprintf(.true.,ERROR,"Could not find geotiff file for "// &
2823 TRIM(field_name))
2824 endif
2825 call read_geogrid_tile(gtiff_file,start_x_dim,end_x_dim, &
2826 start_y_dim,end_y_dim, &
2827 start_z_dim,end_z_dim,array,istatus)
2828 #else
2829 call mprintf(.true.,ERROR,"GEOTIFF file specified for %s, but geogrid not "// &
2830 "compiled with geotiff support.",s1=trim(field_name))
2831 #endif
2832 else
2833 call read_geogrid(file_name, strlen, array, xdim, ydim, zdim, &
2834 sign_convention, local_endian, scalefac, local_wordsize, istatus)
2835 endif
2837 if (irow_order == TOP_BOTTOM) then
2838 allocate(temprow(xdim))
2839 do k=1,zdim
2840 do j=1,ydim
2841 if (ydim-j+1 <= j) exit
2842 temprow(1:xdim) = array(1:xdim,j,k)
2843 array(1:xdim,j,k) = array(1:xdim,ydim-j+1,k)
2844 array(1:xdim,ydim-j+1,k) = temprow(1:xdim)
2845 end do
2846 end do
2847 deallocate(temprow)
2848 end if
2850 if (istatus /= 0) then
2851 start_x_dim = INVALID
2852 start_y_dim = INVALID
2853 end_x_dim = INVALID
2854 end_y_dim = INVALID
2856 call hash_insert(bad_files, file_name)
2857 end if
2859 end subroutine get_data_tile
2862 !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
2863 ! Name: get_row_order
2864 !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
2865 subroutine get_row_order(fieldnm, ilevel, irow_order, istatus)
2867 implicit none
2869 ! Arguments
2870 integer, intent(in) :: ilevel
2871 character (len=128), intent(in) :: fieldnm
2872 integer, intent(out) :: irow_order, istatus
2874 ! Local variables
2875 integer :: idx
2877 istatus = 1
2878 do idx=1,num_entries
2879 if ((index(source_fieldname(idx),trim(fieldnm)) /= 0) .and. &
2880 (len_trim(source_fieldname(idx)) == len_trim(fieldnm))) then
2881 if (ilevel == source_priority(idx)) then
2882 if (is_row_order(idx)) then
2883 irow_order = source_row_order(idx)
2884 istatus = 0
2885 exit
2886 end if
2887 end if
2888 end if
2889 end do
2891 end subroutine get_row_order
2893 logical function is_geotiff_file(fieldname,ilevel)
2894 implicit none
2895 character(len=*),intent(in)::fieldname
2896 integer,intent(in)::ilevel
2897 integer::idx
2898 is_geotiff_file=.false.
2899 do idx=1,num_entries
2900 if ((index(source_fieldname(idx),trim(fieldname)) /= 0) .and. &
2901 (len_trim(source_fieldname(idx)) == len_trim(fieldname)) .and. &
2902 (ilevel == source_priority(idx))) then
2903 is_geotiff_file = is_geotiff(idx)
2904 exit
2905 endif
2906 enddo
2907 end function is_geotiff_file
2910 !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
2911 ! Name: get_tile_dimensions
2912 !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
2913 subroutine get_tile_dimensions(xlat, xlon, start_x_dim, end_x_dim, start_y_dim, end_y_dim, &
2914 start_z_dim, end_z_dim, npts_bdr, bytes_per_datum, endianness, &
2915 sign_convention, ilevel, fieldnm, istatus)
2917 use llxy_module
2919 implicit none
2921 ! Arguments
2922 integer, intent(in) :: ilevel
2923 integer, intent(out) :: start_x_dim, end_x_dim, &
2924 start_y_dim, end_y_dim, &
2925 start_z_dim, end_z_dim, &
2926 npts_bdr, &
2927 bytes_per_datum, endianness, &
2928 sign_convention, istatus
2929 real, intent(in) :: xlat, xlon
2930 character (len=128), intent(in) :: fieldnm
2932 ! Local variables
2933 integer :: idx, current_domain
2934 real :: rx, ry
2936 istatus = 1
2937 do idx=1,num_entries
2938 if ((index(source_fieldname(idx),trim(fieldnm)) /= 0) .and. &
2939 (len_trim(source_fieldname(idx)) == len_trim(fieldnm))) then
2940 if (ilevel == source_priority(idx)) then
2941 istatus = 0
2942 exit
2943 end if
2944 end if
2945 end do
2947 if (istatus /= 0) then
2948 start_x_dim = 1
2949 start_y_dim = 1
2950 start_z_dim = 1
2951 end_x_dim = 1
2952 end_y_dim = 1
2953 end_z_dim = 1
2954 bytes_per_datum = 0
2955 return
2956 end if
2958 current_domain = iget_selected_domain()
2959 call select_domain(SOURCE_PROJ)
2960 call lltoxy(xlat, xlon, rx, ry, M)
2961 call select_domain(current_domain)
2963 start_x_dim = source_tile_x(idx) * nint(real(floor((rx-0.5) / real(source_tile_x(idx))))) + 1
2964 end_x_dim = start_x_dim + source_tile_x(idx) - 1
2966 start_y_dim = source_tile_y(idx) * nint(real(floor((ry-0.5) / real(source_tile_y(idx))))) + 1
2967 end_y_dim = start_y_dim + source_tile_y(idx) - 1
2969 if (is_tile_z_start(idx) .and. is_tile_z_end(idx)) then
2970 start_z_dim = source_tile_z_start(idx)
2971 end_z_dim = source_tile_z_end(idx)
2972 else if (is_tile_z(idx)) then
2973 start_z_dim = 1
2974 end_z_dim = source_tile_z(idx)
2975 end if
2977 if (.not. (is_tile_z_start(idx) .and. is_tile_z_end(idx))) then
2978 if (is_tile_z_start(idx) .or. is_tile_z_end(idx)) then
2979 call mprintf(.true., ERROR, 'In GEOGRID.TBL, only one of tile_z_start and '// &
2980 'tile_z_end specified for entry %i.',i1=idx)
2981 end if
2982 end if
2984 if (is_tile_bdr(idx)) then
2985 npts_bdr = source_tile_bdr(idx)
2986 else
2987 npts_bdr = 0
2988 end if
2990 start_x_dim = start_x_dim - npts_bdr
2991 end_x_dim = end_x_dim + npts_bdr
2992 start_y_dim = start_y_dim - npts_bdr
2993 end_y_dim = end_y_dim + npts_bdr
2995 if (is_wordsize(idx)) then
2996 bytes_per_datum = source_wordsize(idx)
2997 else
2998 bytes_per_datum = 1
2999 call mprintf(.true.,ERROR,'In GEOGRID.TBL, no wordsize specified for data in entry %i.',i1=idx)
3000 end if
3002 if (is_endian(idx)) then
3003 endianness = source_endian(idx)
3004 else
3005 endianness = BIG_ENDIAN
3006 end if
3008 if (is_signed(idx)) then
3009 sign_convention = 1
3010 else
3011 sign_convention = 0
3012 end if
3014 end subroutine get_tile_dimensions
3017 !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
3018 ! Name: get_tile_fname
3019 !
3020 ! Purpose:
3021 !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
3022 subroutine get_tile_fname(test_fname, xlat, xlon, ilevel, fieldname, istatus)
3024 use llxy_module
3025 use gridinfo_module
3027 implicit none
3029 ! Arguments
3030 integer, intent(in) :: ilevel
3031 integer, intent(out) :: istatus
3032 real, intent(in) :: xlat, xlon
3033 character (len=*), intent(in) :: fieldname
3034 character (len=256), intent(out) :: test_fname
3036 ! Local variables
3037 integer :: current_domain, idx
3038 real :: rx, ry
3040 istatus = 1
3041 write(test_fname, '(a)') 'TILE.OUTSIDE.DOMAIN'
3043 do idx=1,num_entries
3044 if ((index(source_fieldname(idx),trim(fieldname)) /= 0) .and. &
3045 (len_trim(source_fieldname(idx)) == len_trim(fieldname))) then
3046 if (ilevel == source_priority(idx)) then
3047 istatus = 0
3048 exit
3049 end if
3050 end if
3051 end do
3053 if (istatus /= 0) return
3055 current_domain = iget_selected_domain()
3056 call select_domain(SOURCE_PROJ)
3057 call lltoxy(xlat, xlon, rx, ry, M)
3058 call select_domain(current_domain)
3060 ! rx = real(source_tile_x(idx)) * real(floor((rx-0.5*source_dx(idx))/ real(source_tile_x(idx)))) + 1.0
3061 ! ry = real(source_tile_y(idx)) * real(floor((ry-0.5*source_dy(idx))/ real(source_tile_y(idx)))) + 1.0
3062 rx = real(source_tile_x(idx)) * real(floor((rx-0.5) / real(source_tile_x(idx)))) + 1.0
3063 ry = real(source_tile_y(idx)) * real(floor((ry-0.5) / real(source_tile_y(idx)))) + 1.0
3065 if (rx > 0. .and. ry > 0.) then
3066 write(test_fname, '(a,i5.5,a1,i5.5,a1,i5.5,a1,i5.5)') trim(source_path(idx)), &
3067 nint(rx),'-',nint(rx)+source_tile_x(idx)-1,'.',nint(ry),'-',nint(ry)+source_tile_y(idx)-1
3068 end if
3070 end subroutine get_tile_fname
3073 !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
3074 ! Name: get_source_resolution
3075 !
3076 ! Purpose:
3077 !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
3078 subroutine get_source_resolution(fieldnm, ilevel, src_dx, src_dy, istatus)
3080 implicit none
3082 ! Arguments
3083 integer, intent(in) :: ilevel
3084 integer, intent(out) :: istatus
3085 real, intent(out) :: src_dx, src_dy
3086 character (len=*), intent(in) :: fieldnm
3088 ! Local variables
3089 integer :: idx
3091 istatus = 1
3092 do idx=1,num_entries
3093 if ((index(source_fieldname(idx),trim(fieldnm)) /= 0) .and. &
3094 (len_trim(source_fieldname(idx)) == len_trim(fieldnm))) then
3095 if (ilevel == source_priority(idx)) then
3096 if (is_dx(idx) .and. is_dy(idx)) then
3097 src_dx = source_dx(idx)
3098 src_dy = source_dy(idx)
3099 if (source_proj(idx) /= PROJ_LATLON) then
3100 src_dx = src_dx / 111000.
3101 src_dy = src_dy / 111000.
3102 end if
3103 istatus = 0
3104 exit
3105 end if
3106 end if
3107 end if
3108 end do
3110 end subroutine get_source_resolution
3113 !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
3114 ! Name: get_data_projection
3115 !
3116 ! Purpose: To acquire the parameters necessary in defining the grid on which
3117 ! the user-specified data for field 'fieldnm' are given.
3118 !
3119 ! NOTES: If the routine successfully acquires values for all necessary
3120 ! parameters, istatus is set to 0. In case of an error,
3121 ! OR IF THE USER HAS NOT SPECIFIED A TILE OF DATA FOR FIELDNM,
3122 ! istatus is set to 1.
3123 !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
3124 subroutine get_data_projection(fieldnm, iproj, stand_lon, truelat1, truelat2, &
3125 dxkm, dykm, known_x, known_y, known_lat, known_lon, ilevel, istatus)
3127 implicit none
3129 ! Arguments
3130 integer, intent(in) :: ilevel
3131 integer, intent(out) :: iproj, istatus
3132 real, intent(out) :: stand_lon, truelat1, truelat2, dxkm, dykm, &
3133 known_x, known_y, known_lat, known_lon
3134 character (len=*), intent(in) :: fieldnm
3136 ! Local variables
3137 integer :: idx
3139 istatus = 1
3140 do idx=1,num_entries
3141 if ((index(source_fieldname(idx),trim(fieldnm)) /= 0) .and. &
3142 (len_trim(source_fieldname(idx)) == len_trim(fieldnm))) then
3143 if (ilevel == source_priority(idx)) then
3144 istatus = 0
3145 if (is_proj(idx)) then
3146 iproj = source_proj(idx)
3147 else
3148 iproj = 1
3149 call mprintf(.true., ERROR, &
3150 'In GEOGRID.TBL, no specification for projection in entry %i.',i1=idx)
3151 end if
3152 if (is_known_x(idx)) then
3153 known_x = source_known_x(idx)
3154 else
3155 known_x = 1.
3156 call mprintf(.true., ERROR, &
3157 'In GEOGRID.TBL, no specification for known_x in entry %i.',i1=idx)
3158 end if
3159 if (is_known_y(idx)) then
3160 known_y = source_known_y(idx)
3161 else
3162 known_y = 1.
3163 call mprintf(.true., ERROR, &
3164 'In GEOGRID.TBL, no specification for known_y in entry %i.',i1=idx)
3165 end if
3166 if (is_known_lat(idx)) then
3167 known_lat = source_known_lat(idx)
3168 else
3169 known_lat = 1.
3170 call mprintf(.true., ERROR, &
3171 'In GEOGRID.TBL, no specification for known_lat in entry %i.',i1=idx)
3172 end if
3173 if (is_known_lon(idx)) then
3174 known_lon = source_known_lon(idx)
3175 else
3176 known_lon = 1.
3177 call mprintf(.true., ERROR, &
3178 'In GEOGRID.TBL, no specification for known_lon in entry %i.',i1=idx)
3179 end if
3180 if (is_truelat1(idx)) then
3181 truelat1 = source_truelat1(idx)
3182 else if (is_proj(idx) .and. source_proj(idx) /= PROJ_LATLON) then
3183 truelat1 = 1.
3184 call mprintf(.true., WARN, &
3185 'In GEOGRID.TBL, no specification for truelat1 in entry %i.',i1=idx)
3186 end if
3187 if (is_truelat2(idx)) then
3188 truelat2 = source_truelat2(idx)
3189 else if (is_proj(idx) .and. source_proj(idx) /= PROJ_LATLON) then
3190 truelat2 = 1.
3191 call mprintf(.true., WARN, &
3192 'In GEOGRID.TBL, no specification for truelat2 in entry %i.',i1=idx)
3193 end if
3194 if (is_stdlon(idx)) then
3195 stand_lon = source_stdlon(idx)
3196 else if (is_proj(idx) .and. source_proj(idx) /= PROJ_LATLON) then
3197 stand_lon = 1.
3198 call mprintf(.true., WARN, &
3199 'In GEOGRID.TBL, no specification for stdlon in entry %i.',i1=idx)
3200 end if
3201 if (is_dx(idx)) then
3202 dxkm = source_dx(idx)
3203 else
3204 dxkm = 1.
3205 call mprintf(.true., ERROR, &
3206 'In GEOGRID.TBL, no specification for dx in entry %i.',i1=idx)
3207 end if
3208 if (is_dy(idx)) then
3209 dykm = source_dy(idx)
3210 else
3211 dykm = 1.
3212 call mprintf(.true., ERROR, &
3213 'In GEOGRID.TBL, no specification for dy in entry %i.',i1=idx)
3214 end if
3215 exit
3216 end if
3217 end if
3218 end do
3220 end subroutine get_data_projection
3222 subroutine get_geotiff_filename(field_name,ilevel,file_name,istatus)
3223 implicit none
3224 character(len=*),intent(in)::field_name
3225 integer,intent(in)::ilevel
3226 character(len=*),intent(out)::file_name
3227 integer,intent(out)::istatus
3228 integer::idx
3229 istatus=1
3230 do idx=1,num_entries
3231 if ((index(source_fieldname(idx),trim(field_name)) /= 0) .and. &
3232 (len_trim(source_fieldname(idx)) == len_trim(field_name))) then
3233 if (ilevel == source_priority(idx)) then
3234 if (is_geotiff(idx)) then
3235 file_name=trim(source_geotiff_file(idx))
3236 istatus = 0
3237 exit
3238 end if
3239 end if
3240 end if
3241 end do
3242 end subroutine get_geotiff_filename
3244 !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
3245 ! Name: check_data_specification
3246 !
3247 ! Purpose: To check for obvious errors in the user source data specifications.
3248 ! Returns .true. if specification passes all checks, and .false. otherwise.
3249 ! For failing checks, diagnostic messages are printed.
3250 !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
3251 function check_data_specification( )
3253 implicit none
3255 ! Return value
3256 logical :: check_data_specification
3258 ! Local variables
3259 integer :: i, j, istatus
3260 integer, pointer, dimension(:) :: priorities
3261 real :: rmissing
3262 logical :: begin_priority, halt
3263 character (len=128) :: cur_name
3265 check_data_specification = .false.
3267 ! Check that each specification has a name, priority level, and path
3268 do i=1,num_entries
3269 if (.not. is_fieldname(i)) then
3270 call mprintf(.true., ERROR, &
3271 'In GEOGRID.TBL, specification %i does not have a name.',i1=i)
3272 end if
3273 if (.not. is_priority(i)) then
3274 call mprintf(.true., ERROR, &
3275 'In GEOGRID.TBL, specification %i does not have a priority.',i1=i)
3276 end if
3277 if (list_length(source_res_path(i)) == 0) then
3278 call mprintf(.true., ERROR, &
3279 'In GEOGRID.TBL, no path (relative or absolute) is specified '// &
3280 'for entry %i.',i1=i)
3281 end if
3282 end do
3284 ! The fill_missing and halt_on_missing options are mutually exclusive
3285 do i=1,num_entries
3286 call get_halt_on_missing(source_fieldname(i), halt, istatus)
3287 call get_missing_fill_value(source_fieldname(i), rmissing, istatus)
3288 if (halt .and. (istatus == 0)) then
3289 call mprintf(.true., ERROR, 'In GEOGRID.TBL, the halt_on_missing and fill_missing '// &
3290 'options are mutually exclusive, but both are given for field %s', &
3291 s1=trim(source_fieldname(i)))
3292 end if
3293 end do
3295 ! Check that the field from which landmask is calculated is not output on a staggering
3296 do i=1,num_entries
3297 if (list_length(source_landmask_land(i)) > 0 .or. list_length(source_landmask_water(i)) > 0) then
3298 if (is_output_stagger(i)) then
3299 if (source_output_stagger(i) /= M) then
3300 call mprintf(.true., ERROR, 'In GEOGRID.TBL, landmask cannot be derived from '// &
3301 'a field that is computed on a staggered grid at entry %i.',i1=i)
3302 end if
3303 end if
3304 end if
3305 end do
3307 ! Also check that any field that is to be masked by the landmask is not output on a staggering
3308 do i=1,num_entries
3309 if (is_masked(i) .and. is_output_stagger(i)) then
3310 if (source_output_stagger(i) /= M) then
3311 call mprintf(.true., ERROR, 'In GEOGRID.TBL, landmask cannot be used with '// &
3312 'a field that is computed on a staggered grid at entry %i.',i1=i)
3313 end if
3314 end if
3315 end do
3317 allocate(priorities(num_entries))
3319 ! Now check that priorities for each source are unique and in the interval [1,n], n <= num_entries
3320 do i=1,num_entries
3321 priorities = 0
3322 cur_name = source_fieldname(i)
3323 do j=1,num_entries
3324 if (source_fieldname(j) == cur_name) then
3326 if (source_priority(j) > num_entries .or. source_priority(j) < 1) then
3327 call mprintf(.true., ERROR, 'In GEOGRID.TBL, priorities for %s do not '// &
3328 'form a sequence 1,2,...,n.', s1=trim(cur_name))
3330 else
3331 if (priorities(source_priority(j)) == 1) then
3332 call mprintf(.true., ERROR, 'In GEOGRID.TBL, more than one entry for %s '// &
3333 'has priority %i.', s1=trim(cur_name), i1=source_priority(j))
3335 else
3336 priorities(source_priority(j)) = 1
3337 end if
3338 end if
3340 end if
3341 end do
3343 begin_priority = .false.
3344 do j=num_entries,1,-1
3345 if (.not.begin_priority .and. priorities(j) == 1) then
3346 begin_priority = .true.
3347 else if (begin_priority .and. priorities(j) == 0) then
3348 call mprintf(.true., ERROR, 'In GEOGRID.TBL, no entry for %s has '// &
3349 'priority %i, but an entry has priority %i.', &
3350 s1=trim(cur_name), i1=j, i2=j+1)
3351 end if
3352 end do
3353 end do
3355 deallocate(priorities)
3357 ! Units must match for all priority levels of a field
3358 do i=1,num_entries
3359 if (source_priority(i) == 1) then
3360 do j=1,num_entries
3361 if ((source_fieldname(i) == source_fieldname(j)) .and. &
3362 (source_units(i) /= source_units(j))) then
3363 call mprintf(.true., ERROR, 'In GEOGRID.TBL, units for %s at entry %i '// &
3364 'do not match units at entry %i (%s)', &
3365 s1=trim(source_fieldname(i)), i1=j, i2=i, s2=trim(source_units(i)))
3366 end if
3367 end do
3368 end if
3369 end do
3371 ! Make sure that user has not asked to calculate landmask from a continuous field
3372 do i=1,num_entries
3373 if (is_dest_fieldtype(i)) then
3374 if (source_dest_fieldtype(i) == CONTINUOUS) then
3375 if (list_length(source_landmask_water(i)) > 0 .or. list_length(source_landmask_land(i)) > 0) then
3376 call mprintf(.true., ERROR, 'In GEOGRID.TBL, landmask cannot be '// &
3377 'calculated from a continuous destination field at entry %i.',i1=i)
3378 end if
3379 end if
3380 end if
3381 end do
3383 ! If either min_category or max_category is specified, then both must be specified
3384 do i=1,num_entries
3385 if (is_category_min(i) .or. is_category_max(i)) then
3386 if (.not. is_category_min(i)) then
3387 call mprintf(.true., ERROR, 'In GEOGRID.TBL, for index file of data at '// &
3388 'entry %i, category_max is specified, but category_min is '// &
3389 'not. Both must be specified.',i1=i)
3390 else if (.not. is_category_max(i)) then
3391 call mprintf(.true., ERROR, 'In GEOGRID.TBL, for index file of data at '// &
3392 'entry %i, category_min is specified, but category_max is '// &
3393 'not. Both must be specified.',i1=i)
3394 end if
3395 end if
3396 end do
3398 ! For continuous data, (category_max - category_min + 1) should equal tile_z
3399 do i=1,num_entries
3400 if (is_fieldtype(i)) then
3401 if (source_fieldtype(i) == CONTINUOUS) then
3402 if (is_category_max(i) .and. is_category_min(i) .and. is_tile_z(i)) then
3403 if (source_tile_z(i) /= (source_category_max(i) - source_category_min(i) + 1)) then
3404 call mprintf(.true., ERROR, 'In GEOGRID.TBL, tile_z must equal '// &
3405 '(category_max - category_min + 1) at entry %i.',i1=i)
3406 end if
3407 else if (is_category_max(i) .and. is_category_min(i) .and. &
3408 is_tile_z_start(i) .and. is_tile_z_end(i)) then
3409 if (source_tile_z_end(i) /= source_category_max(i) .or. &
3410 source_tile_z_start(i) /= source_category_min(i)) then
3411 call mprintf(.true., ERROR, 'In GEOGRID.TBL, tile_z_end must equal '// &
3412 'category_max, and tile_z_start must equal category_min '// &
3413 'at entry %i.',i1=i)
3414 end if
3415 end if
3416 end if
3417 end if
3418 end do
3420 ! Make sure that user has not named a dominant category or computed slope field
3421 ! the same as a fractional field
3422 do i=1,num_entries
3423 if (source_dominant_category(i) == source_fieldname(i)) then
3424 call mprintf(.true., ERROR, 'In GEOGRID.TBL, dominant category cannot have '// &
3425 'the same name as the field at entry %i.',i1=i)
3426 end if
3428 do j=1,num_entries
3429 if (.not. is_dominant_only(i)) then
3430 if (is_dfdx(j)) then
3431 if (source_dfdx(j) == source_fieldname(i)) then
3432 call mprintf(.true., ERROR, 'In GEOGRID.TBL, field name at entry %i '// &
3433 'cannot have the same name as the slope field df_dx at entry %i.', &
3434 i1=i, i2=j)
3435 end if
3436 end if
3437 if (is_dfdy(j)) then
3438 if (source_dfdy(j) == source_fieldname(i)) then
3439 call mprintf(.true., ERROR, 'In GEOGRID.TBL, field name at entry %i '// &
3440 'cannot have the same name as the slope field df_dy at entry %i.', &
3441 i1=i, i2=j)
3442 end if
3443 end if
3444 if (is_dfdx(j) .and. is_dominant_category(i)) then
3445 if (source_dfdx(j) == source_dominant_category(i)) then
3446 call mprintf(.true., ERROR, 'In GEOGRID.TBL, dominant field name at '// &
3447 'entry %i cannot have the same name as the slope field df_dx '// &
3448 'at entry %i.',i1=i, i2=j)
3449 end if
3450 end if
3451 if (is_dfdy(j) .and. is_dominant_category(i)) then
3452 if (source_dfdy(j) == source_dominant_category(i)) then
3453 call mprintf(.true., ERROR, 'In GEOGRID.TBL, dominant field name at '// &
3454 'entry %i cannot have the same name as the slope field '// &
3455 'df_dy at entry %i.',i1=i, i2=j)
3456 end if
3457 end if
3458 else
3459 if (is_dfdx(j)) then
3460 if (source_dfdx(j) == source_dominant_only(i)) then
3461 call mprintf(.true., ERROR, 'In GEOGRID.TBL, dominant field name at '// &
3462 'entry %i cannot have the same name as the slope field '// &
3463 'df_dx at entry %i.',i1=i, i2=j)
3464 end if
3465 end if
3466 if (is_dfdy(j)) then
3467 if (source_dfdy(j) == source_dominant_only(i)) then
3468 call mprintf(.true., ERROR, 'In GEOGRID.TBL, dominant field name at '// &
3469 'entry %i cannot have the same name as the slope field '// &
3470 'df_dy at entry %i.',i1=i, i2=j)
3471 end if
3472 end if
3473 end if
3474 if (i /= j) then
3475 if (is_dfdx(i)) then
3476 if (is_dfdx(j)) then
3477 if (source_dfdx(j) == source_dfdx(i)) then
3478 call mprintf(.true., ERROR, 'In GEOGRID.TBL, slope field df_dx at '// &
3479 'entry %i cannot have the same name as the slope '// &
3480 'field df_dx at entry %i.',i1=i, i2=j)
3481 end if
3482 end if
3483 if (is_dfdy(j)) then
3484 if (source_dfdy(j) == source_dfdx(i)) then
3485 call mprintf(.true., ERROR, 'In GEOGRID.TBL, slope field df_dx at '// &
3486 'entry %i cannot have the same name as the slope field '// &
3487 'df_dy at entry %i.',i1=i, i2=j)
3488 end if
3489 end if
3490 end if
3491 if (is_dfdy(i)) then
3492 if (is_dfdx(j)) then
3493 if (source_dfdx(j) == source_dfdy(i)) then
3494 call mprintf(.true., ERROR, 'In GEOGRID.TBL, slope field df_dy at '// &
3495 'entry %i cannot have the same name as the slope field '// &
3496 'df_dx at entry %i.',i1=i, i2=j)
3497 end if
3498 end if
3499 if (is_dfdy(j)) then
3500 if (source_dfdy(j) == source_dfdy(i)) then
3501 call mprintf(.true., ERROR, 'In GEOGRID.TBL, slope field df_dy at '// &
3502 'entry %i cannot have the same name as the slope field '// &
3503 'df_dy at entry %i.',i1=i, i2=j)
3504 end if
3505 end if
3506 end if
3507 if (is_dominant_category(i)) then
3508 if (source_dominant_category(i) == source_fieldname(j)) then ! Possible exception
3509 if (.not. (is_dominant_only(j) .and. (source_dominant_category(i) /= source_dominant_only(j)))) then
3510 call mprintf(.true., ERROR, 'In GEOGRID.TBL, dominant category at '// &
3511 'entry %i cannot have the same name as the field at '// &
3512 'entry %i.',i1=i, i2=j)
3513 end if
3514 else if (is_dominant_category(j) .and. &
3515 (source_dominant_category(i) == source_dominant_category(j))) then
3516 call mprintf(.true., ERROR, 'In GEOGRID.TBL, dominant category at entry '// &
3517 '%i cannot have the same name as dominant category at '// &
3518 'entry %i.',i1=i, i2=j)
3519 else if (is_dominant_only(j) .and. &
3520 (source_dominant_category(i) == source_dominant_only(j))) then
3521 call mprintf(.true., ERROR, 'In GEOGRID.TBL, dominant category at '// &
3522 'entry %i cannot have the same name as dominant_only '// &
3523 'category at entry %i.',i1=i, i2=j)
3524 end if
3525 else if (is_dominant_only(i)) then
3526 if (source_dominant_only(i) == source_fieldname(j)) then ! Possible exception
3527 if (.not. (is_dominant_only(j) .and. (source_dominant_only(i) /= source_dominant_only(j)))) then
3528 call mprintf(.true., ERROR, 'In GEOGRID.TBL, dominant_only category '// &
3529 'at entry %i cannot have the same name as the field at '// &
3530 'entry %i.',i1=i, i2=j)
3531 end if
3532 else if (is_dominant_category(j) .and. &
3533 (source_dominant_only(i) == source_dominant_category(j))) then
3534 call mprintf(.true., ERROR, 'In GEOGRID.TBL, dominant_only category '// &
3535 'at entry %i cannot have the same name as dominant '// &
3536 'category at entry %i.',i1=i, i2=j)
3537 else if (is_dominant_only(j) .and. &
3538 (source_dominant_only(i) == source_dominant_only(j))) then
3539 call mprintf(.true., ERROR, 'In GEOGRID.TBL, dominant_only category '// &
3540 'at entry %i cannot have the same name as dominant_only '// &
3541 'category at entry %i.',i1=i, i2=j)
3542 end if
3543 end if
3544 end if
3545 end do
3546 end do
3548 check_data_specification = .true.
3549 end function check_data_specification
3552 !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
3553 ! Name: s_len
3554 !
3555 ! Purpose: This routine receives a fortran string, and returns the number of
3556 ! characters in the string before the first "space" is encountered. It
3557 ! considers ascii characters 33 to 126 to be valid characters, and ascii
3558 ! 0 to 32, and 127 to be "space" characters.
3559 !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
3560 subroutine s_len(string, s_length)
3562 implicit none
3564 ! Arguments
3565 character (len=*), intent(in) :: string
3566 integer, intent(out) :: s_length
3568 ! Local variables
3569 integer :: i, len_str, aval
3570 logical :: space
3572 space = .false.
3573 i = 1
3574 len_str = len(string)
3575 s_length = len_str
3576 do while ((i .le. len_str) .and. (.not. space))
3577 aval = ichar(string(i:i))
3578 if ((aval .lt. 33) .or. (aval .gt. 126)) then
3579 s_length = i - 1
3580 space = .true.
3581 endif
3582 i = i + 1
3583 enddo
3585 end subroutine s_len
3587 end module source_data_module