Soil Water Balance (SWB2)
Loading...
Searching...
No Matches
model_initialize.F90
Go to the documentation of this file.
2
3 use iso_c_binding, only : c_int, c_float, c_double, c_bool
5 use datetime
8 use dictionary
9 use exceptions
11 use grid
12! use loop_iterate
13 use logfiles
14 use model_domain, only : model, &
15 minmaxmean, &
18
22 use parameters
24 use simulation_datetime, only : sim_dt
25 use fstring
26 use fstring_list
28 implicit none
29
30 private
31
34
36 character (len=38) :: sname
37 character (len=256) :: spathname
38 logical (c_bool) :: loptional
39 integer (c_int) :: idatatype
40 end type gridded_datasets_t
41
43 character (len=23) :: sname
44 logical (c_bool) :: loptional
45 end type methods_list_t
46
47 integer (c_int), parameter :: number_of_known_grids = 46
48 integer (c_int), parameter :: number_of_known_methods = 18
49
51
52 [ gridded_datasets_t("PRECIPITATION ", "", true, datatype_float ), &
53 gridded_datasets_t("TMIN ", "", true, datatype_float ), &
54 gridded_datasets_t("TMAX ", "", true, datatype_float ), &
55 gridded_datasets_t("AVAILABLE_WATER_CONTENT ", "", true, datatype_float ), &
56 gridded_datasets_t("REFERENCE_ET0 ", "", true, datatype_float ), &
57 gridded_datasets_t("POTENTIAL_ET ", "", true, datatype_float ), &
58 gridded_datasets_t("ACTUAL_ET ", "", true, datatype_float ), &
59 gridded_datasets_t("SOLAR_RADIATION ", "", true, datatype_float ), &
60 gridded_datasets_t("WIND_SPEED ", "", true, datatype_float ), &
61 gridded_datasets_t("RAINFALL_ZONE ", "", true, datatype_int ), &
62 gridded_datasets_t("REFERENCE_ET_ZONE ", "", true, datatype_int ), &
63 gridded_datasets_t("POTENTIAL_ET_ZONE ", "", true, datatype_int ), &
64 gridded_datasets_t("FLOW_DIRECTION ", "", true, datatype_int), &
65 gridded_datasets_t("FOG_RATIO ", "", true, datatype_float ), &
66 gridded_datasets_t("LAND_USE ", "", false, datatype_int ), &
67 gridded_datasets_t("SOILS_CODE ", "", true, datatype_int ), &
68 gridded_datasets_t("HYDROLOGIC_SOILS_GROUP ", "", false, datatype_int ), &
69 gridded_datasets_t("INITIAL_PERCENT_SOIL_MOISTURE ", "", false, datatype_float), &
70 gridded_datasets_t("INITIAL_SNOW_COVER_STORAGE ", "", true, datatype_float), &
71 gridded_datasets_t("INITIAL_CONTINUOUS_FROZEN_GROUND_INDEX", "", true, datatype_float), &
72 gridded_datasets_t("CFGI_LOWER_LIMIT ", "" ,true, datatype_float), &
73 gridded_datasets_t("CFGI_UPPER_LIMIT ", "" ,true, datatype_float), &
74 gridded_datasets_t("PERCENT_CANOPY_COVER ", "", true, datatype_float ), &
75 gridded_datasets_t("PERCENT_PERVIOUS_COVER ", "", true, datatype_float ), &
76 gridded_datasets_t("PERCENT_IMPERVIOUS_COVER ", "", true, datatype_float ), &
77 gridded_datasets_t("FRACTION_CANOPY_COVER ", "", true, datatype_float ), &
78 gridded_datasets_t("FRACTION_PERVIOUS_COVER ", "", true, datatype_float ), &
79 gridded_datasets_t("FRACTION_IMPERVIOUS_COVER ", "", true, datatype_float ), &
80 gridded_datasets_t("STEMFLOW_FRACTION ", "", true, datatype_float ), &
81 gridded_datasets_t("EVAPORATION_TO_RAINFALL_RATIO ", "", true, datatype_float ), &
82 gridded_datasets_t("RAINFALL_ADJUST_FACTOR ", "", true, datatype_float ), &
83 gridded_datasets_t("CESSPOOL_LEAKAGE ", "", true, datatype_float ), &
84 gridded_datasets_t("STORM_DRAIN_CAPTURE_FRACTION ", "", true, datatype_float ), &
85 gridded_datasets_t("WATER_BODY_LEAKAGE ", "", true, datatype_float ), &
86 gridded_datasets_t("WATER_MAIN_LEAKAGE ", "", true, datatype_float ), &
87 gridded_datasets_t("DISPOSAL_WELL_DISCHARGE ", "", true, datatype_float ), &
88 gridded_datasets_t("ANNUAL_DIRECT_NET_INFILTRATION_RATE ", "", true, datatype_float ), &
89 gridded_datasets_t("ANNUAL_SEPTIC_DISCHARGE ", "", true, datatype_float ), &
90 gridded_datasets_t("SEPTIC_DISCHARGE ", "", true, datatype_float ), &
91 gridded_datasets_t("RUNOFF_ZONE ", "", true, datatype_int ), &
92 gridded_datasets_t("POLYGON_ID ", "", true, datatype_int ), &
93 gridded_datasets_t("SOIL_STORAGE_MAX ", "", true, datatype_float ), &
94 gridded_datasets_t("PLANT_AVAILABLE_WATER ", "", true, datatype_float ), &
95 gridded_datasets_t("MAXIMUM_NET_INFILTRATION ", "", true, datatype_float ), &
96 gridded_datasets_t("IRRIGATION_MASK ", "", true, datatype_int), &
97 gridded_datasets_t("RELATIVE_HUMIDITY ", "", true, datatype_float ) ]
98
100
101 [ methods_list_t("INTERCEPTION ", false), &
102 methods_list_t("EVAPOTRANSPIRATION ", false), &
103 methods_list_t("RUNOFF ", false), &
104 methods_list_t("PRECIPITATION ", false), &
105 methods_list_t("FOG ", true), &
106 methods_list_t("AVAILABLE_WATER_CONTENT", true), &
107 methods_list_t("SOIL_STORAGE_MAX ", true), &
108 methods_list_t("SOIL_MOISTURE ", false), &
109 methods_list_t("IRRIGATION ", true), &
110 methods_list_t("CROP_COEFFICIENT ", true), &
111 methods_list_t("GROWING_DEGREE_DAY ", true), &
112 methods_list_t("DYNAMIC_LANDUSE ", true), &
113 methods_list_t("DIRECT_RECHARGE ", true), &
114 methods_list_t("DIRECT_NET_INFILTRATION", true), &
115 methods_list_t("DIRECT_SOIL_MOISTURE ", true), &
116 methods_list_t("FLOW_ROUTING ", true), &
117 methods_list_t("ROOTING_DEPTH ", true), &
118 methods_list_t("DUMP_VARIABLES ", true) ]
119
120 ! grid that will be used in the calculation of cell latitudes
121 type (general_grid_t), pointer :: pcoord_grd
122
123contains
124
125 subroutine initialize_all( output_prefix, output_dirname, data_dirname, &
126 lookup_table_dirname, weather_data_dirname )
127
128! use polygon_summarize, only : initialize_polygon_summarize
129
130 character (len=*), intent(in) :: output_prefix, output_dirname, data_dirname, &
131 lookup_table_dirname, weather_data_dirname
132
133 ! [ LOCALS ]
134 integer (c_int) :: iindex
135
136 call model%set_default_method_pointers()
137
138 ! set output directory names for NetCDF and Surfer/Arc ASCII grid output
139 call grid_set_output_directory_name( output_dirname )
140 call set_data_directory( data_dirname )
141 call set_lookup_table_directory( lookup_table_dirname )
142 call set_weather_data_directory( weather_data_dirname )
143 call set_output_directory( output_dirname )
144 call model%set_output_directory( output_dirname )
145 call set_output_prefix( output_prefix )
146
147 ! define SWB project boundary and geographic projection
149
150 ! define the start and end date for the simulation
152
153 ! specify all tables that have been defined in the control file as ***_LOOKUP_TABLE
155
156 ! scan input file entries for keywords associated with known gridded datasets
157 ! (e.g. PRECIPITATION, TMIN, TMAX, FOG_ZONE, etc.)
158
159 ! if the grid is mentioned in the control file, call the bound "initialize" method that
160 ! will wire in the type of data file. all associated methods will also be acted upon if
161 ! present in the control file (e.g. TMAX_ADD_OFFSET, TMAX_NETCDF_X_VAR, etc.)
162 do iindex = 1, ubound(known_grids, 1)
163
164 call initialize_generic_grid( skey=known_grids(iindex)%sName, &
165 spathname=known_grids(iindex)%sPathname, &
166 loptional=known_grids(iindex)%lOptional, &
167 idatatype=known_grids(iindex)%iDataType )
168
169 enddo
170
171 ! at this point the list of files that need to be incorporated into the parameter
172 ! database should be complete; munge each file in turn and add the parameter values to
173 ! the parameter dictionary
174
175 ! this should be OK to run again; duplicate entries should now be ignored
176 call params%munge_file()
177
178 ! scan the control file input for method specifications
179 ! (e.g. EVAPOTRANSPIRATION_METHOD HARGREAVES-SAMANI )
180 do iindex = 1, ubound(known_methods, 1)
181
182 call initialize_generic_method( skey=known_methods( iindex )%sName, &
183 loptional=known_methods( iindex )%lOptional )
184
185 enddo
186
187 ! check to see that there are no NULL method pointers
188 call model%preflight_check_method_pointers()
189
190 ! bring in soils (HSG, AWC), landuse, and flow direction data from native grids
191 ! and pack that data into vectors for active grid cells only
193
194
196
197 ! call each of the initialization routines associated with the chosen methods
198 call model%initialize_methods()
199
200 ! find general program options and process
202
203 ! allow user to enable or disable specific output files
205 !
206 ! ! open and prepare NetCDF files for output
207 ! call initialize_output( MODEL )
208
209! call initialize_polygon_summarize()
210
211 ! if any fatal warnings have been issued to this point, this
212 ! is the end of the simulation!
214
215 ! dump details about all gridded datasets currently in the data_catalog
216 call dat%print()
217
218 ! dump details about all dictionary keys and values to the debug logfile
219 call cf_dict%print_all( iloglevel=log_debug, lecho=false )
220
221 end subroutine initialize_all
222
223!--------------------------------------------------------------------------------------------------
224
225 !> Initialize soils, landuse, and available water content values.
226 !!
227 !! @NOTE Rather than read the raw values and initialize the SWB
228 !! variables in one step, we have two steps:
229 !! 1. read raw;
230 !! 2. initialize SWB variables.
231 !!
232 !! Two steps are needed because the call to "set_inactive_cells" requires access to the "raw"
233 !! values, not the SWB variables. Basically any cell with negative values found in the landuse, soils groups,
234 !! or awc grids results in that cell being inactivated. Once the complete collection of active cells is
235 !! found, the SWB variables can be initialized.
236 !!
237
239
241
242 call read_landuse_codes()
243
245
246 ! call to following results in the assembly of a table of maximum rooting
247 ! depths by hydrologic soil group and landuse code; max rooting depths
248 ! are also assigned to individual cells
249 call model%update_rooting_depth_table()
250
251 ! read in the 2D array of AWC values; needed so that we can use it to set inactive cells
252 call model%read_AWC_data()
253
254 ! determine which grid cells should be inactive on the basis of available water content,
255 ! hydrologic soils group, and/or soils group grids
256 call model%set_inactive_cells()
257
258 ! allocate and zero out model arrays
259 call model%initialize_arrays()
260
261 ! create 1D vectors containing row-column indices of unpacked grid
262 call model%initialize_row_column_indices()
263
264 ! read_polygon matches HWB polygon attributes with the corresponding gridcells in SWB
265! call read_polygon_id()
266
267 ! now that we know which cells are active, pack the 2D array in to 1D vector of AWC values
268 call model%init_AWC()
269
270 ! to this point, rooting depths
271 model%rooting_depth_max = pack( rooting_depth_max, model%active )
272
273 ! if the method selected reads these directly from a grid, the rooting depths previously
274 ! calculated will be overwritten and recalculated using the gridded max soil moisture and awc values
275 call model%init_soil_storage_max()
276
278
280
282
283!--------------------------------------------------------------------------------------------------
284
286
288
290
292
294
296
297 call model%init_continuous_frozen_ground_index()
298
300
301 call storm_drain_capture_initialize(model%active, model%landuse_index)
302
303 call model%initialize_growing_season()
304
305 end subroutine initialize_ancillary_values
306
307 !--------------------------------------------------------------------------------------------------
308
309 subroutine set_lookup_table_directory( lookup_table_dirname )
310
311 character(len=*), intent(in) :: lookup_table_dirname
312
313 ! setting the MODULE variable LOOKUP_TABLE_DIRECTORY_NAME, module = constants_and_conversions
314 if( len_trim(lookup_table_dirname) > 0 ) lookup_table_directory_name = lookup_table_dirname
315
316 end subroutine set_lookup_table_directory
317
318 !--------------------------------------------------------------------------------------------------
319
320 subroutine set_data_directory( data_dirname )
321
322 character(len=*), intent(in) :: data_dirname
323
324 integer (c_int) :: indx
325
326 if ( len_trim(data_dirname) > 0 ) then
327
328 do indx=1,ubound(known_grids,1)
329
330 select case (trim(known_grids(indx)%sName))
331
332 case("PRECIPITATION","TMIN","TMAX")
333
334 case default
335
336 known_grids(indx)%sPathname = trim( data_dirname )
337
338 end select
339
340 enddo
341
342 ! setting the MODULE variable DATA_DIRECTORY_NAME, module = file_operations
343 data_directory_name = data_dirname
344
345 endif
346
347 end subroutine set_data_directory
348
349 !--------------------------------------------------------------------------------------------------
350
351 subroutine set_weather_data_directory( weather_data_dirname )
352
353 character(len=*), intent(in) :: weather_data_dirname
354
355 integer (c_int) :: indx
356
357 do indx=1,ubound(known_grids,1)
358
359 select case (trim(known_grids(indx)%sName))
360
361 case("PRECIPITATION","TMIN","TMAX")
362
363 known_grids(indx)%sPathname = trim( weather_data_dirname )
364
365 case default
366
367 end select
368
369 enddo
370
371 end subroutine set_weather_data_directory
372
373!--------------------------------------------------------------------------------------------------
374
376
377 type (DATA_CATALOG_ENTRY_T), pointer :: pINITIAL_SNOW_COVER_STORAGE
378
379 ! [ LOCALS ]
380 real (c_float), allocatable :: fInitial_Snow_Cover_Storage(:)
381 integer (c_int) :: iStat
382
383 allocate ( finitial_snow_cover_storage( count( model%active ) ), stat=istat )
384
385 ! locate the data structure associated with the gridded initial_snow_cover_storage
386 pinitial_snow_cover_storage => dat%find("INITIAL_SNOW_COVER_STORAGE")
387
388 ! OK, now try the old SWB 1.0 syntax
389 if ( .not. associated( pinitial_snow_cover_storage ) ) &
390 ! locate the data structure associated with the gridded initial_snow_cover_storage
391 pinitial_snow_cover_storage => dat%find("INITIAL_SNOW_COVER")
392
393 if ( .not. associated( pinitial_snow_cover_storage ) ) then
394 call warn(smessage="An INITIAL_SNOW_COVER_STORAGE grid (or constant) was not found.", &
395 shints="Check your control file to see that a valid INITIAL_SNOW_COVER_STORAGE grid or" &
396 //" constant is specified.", lfatal=false )
397
398 model%snow_storage = 0.0_c_float
399
400 else
401
402 call pinitial_snow_cover_storage%getvalues()
403
404 ! map the 2D array of INITIAL_PERCENT_SOIL_MOISTURE values to the vector of active cells
405 finitial_snow_cover_storage = pack( pinitial_snow_cover_storage%pGrdBase%rData, model%active )
406
407 if ( minval( finitial_snow_cover_storage ) < fzero &
408 .or. maxval( finitial_snow_cover_storage ) > 300.0_c_float ) &
409 call warn(smessage="One or more initial snow cover storage values outside of " &
410 //"valid range (0 to 300)", lfatal=true )
411
412 model%snow_storage = finitial_snow_cover_storage
413
414 endif
415
416 end subroutine initialize_snow_storage
417
418!--------------------------------------------------------------------------------------------------
419
421
422 type (DATA_CATALOG_ENTRY_T), pointer :: pINITIAL_PERCENT_SOIL_MOISTURE
423
424 ! [ LOCALS ]
425 real (c_float), allocatable :: fInitial_Percent_Soil_Moisture(:)
426 integer (c_int) :: iStat
427
428 allocate ( finitial_percent_soil_moisture( count( model%active ) ), stat=istat )
429
430 ! locate the data structure associated with the gridded initial_percent_soil_moisture entries
431 pinitial_percent_soil_moisture => dat%find("INITIAL_PERCENT_SOIL_MOISTURE")
432
433 ! OK, now try the old SWB 1.0 syntax
434 if ( .not. associated( pinitial_percent_soil_moisture ) ) &
435 pinitial_percent_soil_moisture => dat%find("INITIAL_SOIL_MOISTURE")
436
437 if ( .not. associated( pinitial_percent_soil_moisture ) ) then
438 call warn(smessage="An INITIAL_PERCENT_SOIL_MOISTURE grid (or constant) was not found.", &
439 shints="Check your control file to see that a valid INITIAL_PERCENT_SOIL_MOISTURE grid or" &
440 //" constant is specified.", lfatal=true )
441 else
442
443 call pinitial_percent_soil_moisture%getvalues()
444
445 ! map the 2D array of INITIAL_PERCENT_SOIL_MOISTURE values to the vector of active cells
446 finitial_percent_soil_moisture = pack( pinitial_percent_soil_moisture%pGrdBase%rData, model%active )
447
448 if ( minval( finitial_percent_soil_moisture ) < fzero &
449 .or. maxval( finitial_percent_soil_moisture ) > 100.0_c_float ) &
450 call warn(smessage="One or more initial percent soils moisture values outside of " &
451 //"valid range (0% to 100%)", lfatal=true )
452
453 model%soil_storage = finitial_percent_soil_moisture / 100.0_c_float * model%soil_storage_max
454
455 endif
456
457 end subroutine initialize_soil_storage
458
459!--------------------------------------------------------------------------------------------------
460
462
463 ! [ LOCALS ]
464 type (DATA_CATALOG_ENTRY_T), pointer :: pPERCENT_IMPERVIOUS
465 type (DATA_CATALOG_ENTRY_T), pointer :: pPERCENT_PERVIOUS
466 type (DATA_CATALOG_ENTRY_T), pointer :: pFRACTION_IMPERVIOUS
467 type (DATA_CATALOG_ENTRY_T), pointer :: pFRACTION_PERVIOUS
468 type ( GENERAL_GRID_T ), pointer :: pTempGrd
469
470 ppercent_impervious => dat%find("PERCENT_IMPERVIOUS_COVER")
471 ppercent_pervious => dat%find("PERCENT_PERVIOUS_COVER")
472 pfraction_impervious => dat%find("FRACTION_IMPERVIOUS_COVER")
473 pfraction_pervious => dat%find("FRACTION_PERVIOUS_COVER")
474
475 if ( associated(ppercent_impervious) ) then
476
477 call ppercent_impervious%getvalues()
478
479 if (associated( ppercent_impervious%pGrdBase) ) then
480 model%pervious_fraction = pack( 1.0_c_float - ppercent_impervious%pGrdBase%rData/100.0_c_float, model%active )
481 else
482 call die("INTERNAL PROGRAMMING ERROR: attempted use of NULL pointer", __file__, __line__)
483 endif
484
485 elseif ( associated( ppercent_pervious ) ) then
486
487 call ppercent_pervious%getvalues()
488
489 if (associated( ppercent_pervious%pGrdBase) ) then
490 model%pervious_fraction = pack( (ppercent_pervious%pGrdBase%rData/100.0_c_float), model%active )
491 else
492 call die("INTERNAL PROGRAMMING ERROR: attempted use of NULL pointer", __file__, __line__)
493 endif
494
495 elseif ( associated(pfraction_impervious) ) then
496
497 call pfraction_impervious%getvalues()
498
499 if (associated( pfraction_impervious%pGrdBase) ) then
500 model%pervious_fraction = pack( 1.0_c_float - pfraction_impervious%pGrdBase%rData, model%active )
501 else
502 call die("INTERNAL PROGRAMMING ERROR: attempted use of NULL pointer", __file__, __line__)
503 endif
504
505 elseif ( associated( pfraction_pervious ) ) then
506
507 call pfraction_pervious%getvalues()
508
509 if (associated( pfraction_pervious%pGrdBase) ) then
510 model%pervious_fraction = pack( pfraction_pervious%pGrdBase%rData, model%active )
511 else
512 call die("INTERNAL PROGRAMMING ERROR: attempted use of NULL pointer", __file__, __line__)
513 endif
514
515 else
516
517 model%pervious_fraction = 1.0_c_float
518
519 endif
520
521 if ( minval( model%pervious_fraction ) < fzero &
522 .or. maxval( model%pervious_fraction ) > 1.0_c_float ) &
523 call warn(smessage="One or more percent (im)pervious cover percent/fraction values are outside of " &
524 //"valid range (0% to 100% or 0.0 to 1.0)", lfatal=true )
525
526 if ( all( model%pervious_fraction < 0.01_c_float ) ) &
527 call warn(smessage="All (im)pervious cover percent/fraction values are suspiciously low " &
528 //"(less than 1% or less than 0.01)", lfatal=true, &
529 shints="Check to see whether (im)pervious cover is expressed as a fraction (0.0-1.0)" &
530 //" or a percentage (0-100%)." )
531
532 ptempgrd => grid_create( inx=model%number_of_columns, iny=model%number_of_rows, &
533 rx0=model%X_ll, ry0=model%Y_ll, &
534 rgridcellsize=model%gridcellsize, idatatype=grid_datatype_real )
535
536 ptempgrd%rData = unpack( model%pervious_fraction, model%active, model%nodata_fill_value )
537
538 call grid_writearcgrid( sfilename="Fraction_pervious_surface__as_read_in_unitless.asc", pgrd=ptempgrd )
539
540 call grid_destroy( ptempgrd )
541
542 end subroutine initialize_percent_pervious
543
544!--------------------------------------------------------------------------------------------------
545
547
548 ! [ LOCALS ]
549 type (DATA_CATALOG_ENTRY_T), pointer :: pPERCENT_CANOPY_COVER
550 type (DATA_CATALOG_ENTRY_T), pointer :: pFRACTION_CANOPY_COVER
551 type ( GENERAL_GRID_T ), pointer :: pTempGrd
552
553 ppercent_canopy_cover => dat%find("PERCENT_CANOPY_COVER")
554 pfraction_canopy_cover => dat%find("FRACTION_CANOPY_COVER")
555
556 if ( associated(ppercent_canopy_cover) ) then
557
558 call ppercent_canopy_cover%getvalues()
559
560 if (associated( ppercent_canopy_cover%pGrdBase) ) then
561 model%canopy_cover_fraction = pack( ppercent_canopy_cover%pGrdBase%rData/100.0_c_float, model%active )
562 else
563 call die("INTERNAL PROGRAMMING ERROR: attempted use of NULL pointer", __file__, __line__)
564 endif
565
566 elseif ( associated(pfraction_canopy_cover) ) then
567
568 call pfraction_canopy_cover%getvalues()
569
570 if (associated( pfraction_canopy_cover%pGrdBase) ) then
571 model%canopy_cover_fraction = pack( pfraction_canopy_cover%pGrdBase%rData, model%active )
572 else
573 call die("INTERNAL PROGRAMMING ERROR: attempted use of NULL pointer", __file__, __line__)
574 endif
575
576 else
577
578 model%canopy_cover_fraction = 1.0_c_float
579
580 call warn("Could not find a grid or constant value for the canopy cover fraction. Using a" &
581 //" value of 1.0 for the entire model domain." )
582
583 endif
584
585 if ( minval( model%canopy_cover_fraction ) < fzero &
586 .or. maxval( model%canopy_cover_fraction ) > 1.0_c_float ) &
587 call warn(smessage="One or more percent canopy cover percent/fraction values values are outside of " &
588 //"valid range (0% to 100% or 0.0 to 1.0)", lfatal=true )
589
590 if ( all( model%canopy_cover_fraction < 0.01_c_float ) ) &
591 call warn(smessage="All canopy cover percent/fraction values are suspiciously low " &
592 //"(less than 1% or less than 0.01)", lfatal=true, &
593 shints="Check to see whether canopy cover is expressed as a fraction (0.0-1.0)" &
594 //" or a percentage (0-100%)." )
595
596 ptempgrd => grid_create( inx=model%number_of_columns, iny=model%number_of_rows, &
597 rx0=model%X_ll, ry0=model%Y_ll, &
598 rgridcellsize=model%gridcellsize, idatatype=grid_datatype_real )
599
600 ptempgrd%rData = unpack( model%canopy_cover_fraction, model%active, model%nodata_fill_value )
601
602 call grid_writearcgrid( sfilename="Fraction_canopy_cover__as_read_in_unitless.asc", pgrd=ptempgrd )
603
604 call grid_destroy( ptempgrd )
605
607
608!--------------------------------------------------------------------------------------------------
609
611
612 ! [ LOCALS ]
613 integer (c_int) :: iIndex
614 type (DATA_CATALOG_ENTRY_T), pointer :: pHSG
615
616 phsg => dat%find("HYDROLOGIC_SOILS_GROUP")
617
618 if ( associated(phsg) ) then
619
620 if (associated( phsg%pGrdBase) ) then
621 model%soil_group = pack( phsg%pGrdBase%iData, model%active )
622 else
623 call die("INTERNAL PROGRAMMING ERROR: attempted use of NULL pointer", __file__, __line__)
624 endif
625
626 else
627
628 call die("Attempted use of NULL pointer. Failed to find HYDROLOGIC_SOILS_GROUP data element.", &
629 __file__, __line__)
630
631 endif
632
633 call logs%write("Hydrologic soils groups as read into SWB data structure", ilinesbefore=1, ilinesafter=1, iloglevel=log_debug)
634
635 do iindex = 1, maxval(phsg%pGrdBase%iData)
636
637 call logs%write( ascharacter(count(model%soil_group == iindex) )//" cells belong to soils group " &
638 //ascharacter(iindex), iloglevel=log_debug )
639
640 end do
641
642 call logs%write("", ilinesbefore=1, iloglevel=log_debug)
643
645
646!--------------------------------------------------------------------------------------------------
647
649
650 ! [ LOCALS ]
651 type (DATA_CATALOG_ENTRY_T), pointer :: pHSG
652
653 phsg => dat%find("HYDROLOGIC_SOILS_GROUP")
654
655 if ( associated(phsg) ) then
656
657 call phsg%getvalues()
658 call grid_writearcgrid("Hydrologic_soil_groups__as_read_into_SWB.asc", phsg%pGrdBase )
659
660 else
661
662 call warn(smessage="HYDROLOGIC_SOILS_GROUP dataset is flawed or missing.", lfatal=true, &
663 iloglevel = log_all, shints="Check to see that a valid path and filename have" &
664 //" been ~included in the control file for the HYDROLOGIC_SOILS_GROUP dataset.", &
665 lecho = true )
666
667 endif
668
669 end subroutine read_hydrologic_soil_groups
670
671!--------------------------------------------------------------------------------------------------
672
673 subroutine read_polygon_id()
674
675
676
677 ! [ LOCALS ]
678 type (DATA_CATALOG_ENTRY_T), pointer :: pPOLYGON_ID
679 logical (c_bool) :: any_problems
680 type (GENERAL_GRID_T), pointer :: pTempGrd
681
682
683
684 ptempgrd => grid_create( inx=model%number_of_columns, iny=model%number_of_rows, &
685 rx0=model%X_ll, ry0=model%Y_ll, &
686 rgridcellsize=model%gridcellsize, idatatype=grid_datatype_real )
687
688 any_problems = true
689
690 do
691
692 ppolygon_id => dat%find("POLYGON_ID")
693
694 if ( .not. associated( ppolygon_id ) ) exit
695
696 call ppolygon_id%getvalues()
697 call grid_writearcgrid("Polygon_ID__as_read_into_SWB.asc", ppolygon_id%pGrdBase )
698
699 model%polygon_id = pack( ppolygon_id%pGrdBase%iData, model%active )
700
701 exit
702
703 enddo
704
705! call slList%append("POLY_ID")
706
707! !> Determine how many POLYGON_IDS are present
708! call PARAMS%get_parameters( slKeys=slList, iValues=polygon_id )
709! iNumberOfPolygonIDs = count( polygon_id > 0 )
710! call slList%clear()
711
712! call slList%append("rooting_depth_inches")
713
714! ! retrieve rooting depths, in inches
715! call PARAMS%get_parameters( slKeys=slList, fValues=rooting_depth_inches )
716! call slList%clear()
717
718! call slList%append("SMC")
719
720! ! retrieve soil moisture maximum content, in inches
721! call PARAMS%get_parameters( slKeys=slList, fValues=soil_moisture_storage )
722! call slList%clear()
723
724! if ( ubound( soil_moisture_storage, 1) == ubound( polygon_id, 1 ) ) any_problems = FALSE
725
726! exit
727
728! enddo
729
730! print *, any_problems, ubound( polygon_id, 1), ubound( soil_moisture_storage, 1), ubound( rooting_depth_inches, 1 )
731! call minmaxmean( soil_moisture_storage, "Max_soil_storage")
732! call minmaxmean( rooting_depth_inches, "Rooting_depth_inches")
733
734! call assert( .not. any_problems, "One or more steps failed while processing POLYGON_ID.", &
735! __FILE__, __LINE__ )
736
737! !$OMP PARALLEL DO
738
739! do index=1, ubound( polygon_id, 1)
740
741! if ( mod( index, 10 ) == 0 ) print *, "Processing polygon number "//asCharacter(index) &
742! //" of "//asCharacter( ubound( polygon_id, 1) )
743
744! where ( MODEL%polygon_id == polygon_id( index ) )
745
746! MODEL%hwb_rooting_depth = rooting_depth_inches( MODEL%landuse_index( index ) ) / 12.0_c_float
747
748! MODEL%hwb_soil_storage_max = soil_moisture_storage( index )
749
750! end where
751
752! enddo
753
754! !$OMP END PARALLEL DO
755
756! where ( MODEL%hwb_rooting_depth > 0.0_c_float )
757
758! MODEL%hwb_awc_in_per_ft = MODEL%hwb_soil_storage_max / MODEL%hwb_rooting_depth
759
760! end where
761
762! pTempGrd%rData = unpack( MODEL%hwb_rooting_depth, MODEL%active, MODEL%nodata_fill_value )
763
764! call grid_WriteArcGrid("Maximum_rooting_depth__as_read_in_from_HWB_input__feet.asc", pTempGrd )
765
766
767! pTempGrd%rData = unpack( MODEL%hwb_soil_storage_max, MODEL%active, MODEL%nodata_fill_value )
768
769! call grid_WriteArcGrid("Maximum_soil_storage__as_read_in_from_HWB_input__inches.asc", pTempGrd )
770
771
772! pTempGrd%rData = unpack( MODEL%hwb_awc_in_per_ft, MODEL%active, MODEL%nodata_fill_value )
773
774! call grid_WriteArcGrid("Available_water_content__as_read_in_from_HWB_input__inches_per_foot.asc", pTempGrd )
775
776 call grid_destroy( ptempgrd )
777
778 end subroutine read_polygon_id
779
780!--------------------------------------------------------------------------------------------------
781
782! subroutine read_landuse_codes
783!
784! ! [ LOCALS ]
785! type (DATA_CATALOG_ENTRY_T), pointer :: pLULC
786!
787! pLULC => DAT%find("LAND_USE")
788!
789! if ( associated(pLULC) ) then
790!
791! call pLULC%getvalues()
792! call grid_WriteArcGrid("Landuse_land_cover__as_read_into_SWB.asc", pLULC%pGrdBase )
793!
794! else
795!
796! call warn(sMessage="LAND_USE dataset is flawed or missing.", lFatal=TRUE, &
797! iLogLevel = LOG_ALL, sHints="Check to see that a valid path and filename have" &
798! //" been ~included in the control file for the LAND_USE dataset.", &
799! lEcho = TRUE )
800!
801! endif
802!
803! end subroutine read_landuse_codes
804!
805! !--------------------------------------------------------------------------------------------------
806
807 subroutine write_control_file( sFilename, sGridSpecification, slExtraDirectives )
808
809 character (len=*), intent(in) :: sfilename
810 character (len=*), intent(in), optional :: sgridspecification
811 type (fstring_list_t), optional :: slextradirectives
812
813 ! [ LOCALS ]
814 character (len=:), allocatable :: stext
815 integer (c_int) :: iindex
816 integer (c_int) :: icount
817 type (ascii_file_t), allocatable :: cf
818 type (dict_entry_t), pointer :: pdict
819
820 pdict => null()
821 allocate(cf)
822
823 call cf%open( sfilename = sfilename )
824
825 ! get GRID specification and move pointer past it; discard values
826 if ( present( sgridspecification ) ) then
827 call cf_dict%get_value(stext, "GRID")
828 pdict => cf_dict%get_next_entry()
829 call assert( associated( pdict ), "INTERNAL PROGRAMMING ERROR -- Attempted use of null poitner", &
830 __file__, __line__ )
831 call cf%writeLine( trim( sgridspecification ) )
832
833 else
834
835 pdict => cf_dict%get_entry(1)
836
837 endif
838
839 do while ( associated( pdict ) )
840
841 call cf_dict%get_value( stext )
842 call cf%writeLine( trim(pdict%key)//" "//stext )
843 pdict => cf_dict%get_next_entry()
844
845 enddo
846
847 if ( present( slextradirectives ) ) then
848
849 icount = slextradirectives%count
850
851 do iindex=1, icount
852
853 call cf%writeLine( trim( slextradirectives%get( iindex ) ) )
854
855 enddo
856
857 endif
858
859 call cf%close()
860
861 end subroutine write_control_file
862
863!--------------------------------------------------------------------------------------------------
864
865 subroutine read_control_file( sFilename )
866
867 character (len=*), intent(in) :: sfilename
868
869 ! [ LOCALS ]
870 character (len=256) :: srecord, skey, svalue
871 integer (c_int) :: istat
872 type (ascii_file_t), allocatable :: cf
873 integer (c_int) :: dumpfile_count
874
875 allocate(cf)
876 dumpfile_count = 0
877
878 ! open the control file and define the comment characters and delimiters to be used in
879 ! parsing the ASCII text
880 call cf%open( sfilename = sfilename, &
881 scommentchars = "#%!+=|[{(-*$", &
882 sdelimiters = "WHITESPACE", &
883 lhasheader = .false._c_bool )
884
885 ! set flag to have 'chomp' remove extra whitespace characters if found between items of interest
886 cf%remove_extra_delimiters = true
887
888 do
889
890 ! read in next line of the control file
891 srecord = cf%readLine()
892
893 if ( cf%isEOF() ) exit
894
895 ! create and allocate memory for a single dictionary entry
896 cf_entry => null()
897 allocate( cf_entry, stat=istat )
898 call assert(istat == 0, "Failed to allocate memory for dictionary object", &
899 __file__, __line__ )
900
901 ! break off key value for the current record
902 call chomp(srecord, skey, cf%sDelimiters, cf%remove_extra_delimiters )
903
904 if ( len_trim( skey ) > 0 ) then
905
906 ! there can be more than one instance of 'DUMP_VARIABLES'
907 ! need to make the key unique; tack on an instance number
908 if ( skey .strequal. "DUMP_VARIABLES" ) then
909 dumpfile_count = dumpfile_count + 1
910 skey = trim(skey)//"_"//ascharacter(dumpfile_count)
911 endif
912
913 ! first add the key value to the directory entry data structure
914 call cf_entry%add_key( skey )
915
916 ! break off first directive for the current record
917 call chomp( srecord, svalue, cf%sDelimiters, cf%remove_extra_delimiters )
918
919 do while ( len_trim( svalue ) > 0 )
920
921 ! add the next directive snippet to dictionary entry data structure
922 call cf_entry%add_value( svalue )
923 ! break off next directive for the current record
924 call chomp( srecord, svalue, cf%sDelimiters, cf%remove_extra_delimiters )
925
926 enddo
927
928 ! add the dictionary entry to the dictionary data structure
929 call cf_dict%add_entry( cf_entry )
930
931 endif
932
933 enddo
934
935 ! close the control file
936 call cf%close()
937
938 end subroutine read_control_file
939
940!--------------------------------------------------------------------------------------------------
941
942
943
944 !> Generic routine to handle intake of gridded data.
945 !!
946 !! @param[in] sKey Name for the data type being processed.
947 !! @param[in] lOptional If lOptional is TRUE, kill model run eventually IF sKey
948 !! does not match a known data type.
949 !! @param[in] iDataType Datatype as defined in @ref constants_and-conversion.F90
950 !!
951 !! This routine accepts a data grid type, for example, "PRECIPITATION", and attempts to
952 !! handle all related control file directives associated with MODEL data type. In MODEL way,
953 !! a new gridded data type may be added simply by extending the list of known data types
954 !! to the list defined in the module variable @ref KNOWN_TYPES.
955 !!
956
957 subroutine initialize_generic_grid(sKey, sPathname, lOptional, iDataType )
958
959 character (len=*), intent(in) :: sKey
960 character (len=*), intent(in) :: sPathname
961 logical (c_bool), intent(in) :: lOptional
962 integer (c_int), intent(in) :: iDataType
963
964 ! [ LOCALS ]
965 type (FSTRING_LIST_T) :: myDirectives
966 type (FSTRING_LIST_T) :: myOptions
967 integer (c_int) :: iIndex
968 character (len=512) :: sCmdText
969 character (len=512) :: sArgText
970 character (len=512) :: sArgText_1
971 character (len=512) :: sArgText_2
972 character (len=512) :: sArgText_3
973 character (len=512) :: sArgText_4
974 integer (c_int) :: iStat
975 type (DATA_CATALOG_ENTRY_T), pointer :: pENTRY
976 logical (c_bool) :: lGridPresent
977
978 pentry => null()
979 lgridpresent = false
980
981 ! obtain a string list of directives that contain the keyword
982 ! (e.g. if the key is "TMIN", MODEL might return:
983 ! "TMIN ARC_ASCII input/mygrid.asc"
984 ! "TMIN_PROJECTION_DEFINITION =Proj=latlon +datum=WGS84"
985 ! "TMIN_MINIMUM_ALLOWED_VALUE -60.0" )
986 ! the loop below then handles the specific directives in turn
987 mydirectives = cf_dict%grep_keys( skey )
988
989! call myDirectives%print
990
991 if ( mydirectives%count == 0 ) then
992
993 call logs%write("Your control file is missing gridded data relating to "//dquote(skey)//".", &
994 iloglevel=log_all, lecho=false )
995
996 if (.not. loptional) then
997 call warn("Your control file is missing gridded data relating to "//dquote(skey)//".", &
998 lfatal = true )
999
1000 endif
1001
1002 else
1003
1004 ! allocate memory for a generic data_catalog_entry
1005 allocate(pentry, stat=istat)
1006 call assert( istat == 0, "Failed to allocate memory for the "//dquote(skey)//" data structure", &
1007 __file__, __line__ )
1008
1009 scmdtext = ""
1010 sargtext = ""
1011 sargtext_1 = ""
1012 sargtext_2 = ""
1013 sargtext_3 = ""
1014 sargtext_4= ""
1015 call myoptions%clear()
1016
1017 ! process all known directives associated with key word
1018 do iindex = 1, mydirectives%count
1019
1020 ! myDirectives is a string list of all SWB directives that contain sKey
1021 ! sCmdText contains an individual directive
1022 scmdtext = mydirectives%get(iindex)
1023
1024 ! For MODEL directive, obtain the associated dictionary entries
1025 call cf_dict%get_values(scmdtext, myoptions )
1026
1027 ! most of the time, we only care about the first dictionary entry, obtained below
1028 sargtext_1 = myoptions%get(1)
1029 sargtext_2 = myoptions%get(2)
1030 sargtext_3 = myoptions%get(3)
1031 sargtext_4 = myoptions%get(4)
1032
1033 ! dictionary entries are initially space-delimited; sArgText contains
1034 ! all dictionary entries present, concatenated, with a space between entries
1035 sargtext = myoptions%get(1, myoptions%count )
1036
1037 ! echo the original directive and dictionary entries to the logfile
1038 call logs%write("> "//trim(scmdtext)//" "//trim(sargtext), ilinesbefore=1 )
1039
1040 ! first option is that the key value and directive are the same
1041 ! (.e.g. "PRECIPITATION"; no trailing underscores or modifiers )
1042 ! MODEL is a grid definition directive
1043 if ( scmdtext .strapprox. skey ) then
1044
1045 pentry%sVariableName_z = aslowercase( skey )
1046
1047 ! determine the type of grid and act appropriately
1048 if (sargtext_1 .strapprox. "CONSTANT" ) then
1049
1050 if (.not. is_numeric(sargtext_2)) &
1051 call die("Non-numeric argument supplied as a CONSTANT in your control file.", &
1052 shints="Faulty control file argument was "//squote(trim(scmdtext)//" "//trim(sargtext)))
1053
1054 select case ( idatatype )
1055
1056 case ( datatype_float )
1057
1058 call pentry%initialize( &
1059 sdescription=trim(scmdtext), &
1060 rconstant=asfloat(sargtext_2) )
1061 lgridpresent = true
1062
1063 case ( datatype_int )
1064
1065 call pentry%initialize( &
1066 sdescription=trim(scmdtext), &
1067 iconstant=asint(sargtext_2) )
1068 lgridpresent = true
1069
1070 case default
1071
1072 call die( "INTERNAL PROGRAMMING ERROR: Unhandled data type selected.", &
1073 __file__, __line__ )
1074
1075 end select
1076
1077 elseif ( (sargtext_1 .strapprox. "TABLE") &
1078 .or. (sargtext_1 .strapprox. "TABLE_LOOKUP") ) then
1079
1080 if (len_trim(sargtext_2) > 0) &
1081 call params%add_file( fix_pathname(trim(spathname)//sargtext_2))
1082
1083 select case ( idatatype )
1084
1085 case ( datatype_float )
1086
1087 call pentry%initialize( &
1088 sdescription=trim(scmdtext), &
1089 sdatecolumnname = "date", &
1090 svaluecolumnname = pentry%sVariableName_z, &
1091 stype = "float")
1092 lgridpresent = true
1093
1094 case ( datatype_int )
1095
1096 call pentry%initialize( &
1097 sdescription=trim(scmdtext), &
1098 sdatecolumnname = "date", &
1099 svaluecolumnname = pentry%sVariableName_z, &
1100 stype = "integer")
1101 lgridpresent = true
1102
1103 case default
1104
1105 call die( "INTERNAL PROGRAMMING ERROR: Unhandled data type selected.", &
1106 __file__, __line__ )
1107 end select
1108
1109 elseif ( (sargtext_1 .strapprox. "ARC_ASCII") &
1110 .or. (sargtext_1 .strapprox. "SURFER") &
1111 .or. (sargtext_1 .strapprox. "ARC_GRID") ) then
1112
1113 call pentry%initialize( &
1114 sdescription=trim(scmdtext), &
1115 sfiletype=trim(sargtext_1), &
1116 sfilename=fix_pathname(trim(spathname)//sargtext_2), &
1117 idatatype=idatatype )
1118 lgridpresent = true
1119
1120 elseif ( sargtext_1 .strapprox. "NETCDF" ) then
1121
1122 call pentry%initialize_netcdf( &
1123 sdescription=trim(scmdtext), &
1124 sfilename = fix_pathname(trim(spathname)//sargtext_2), &
1125 idatatype=idatatype )
1126 lgridpresent = true
1127
1128 ! elseif ( sArgText_1 .strapprox. "TABLE" ) then
1129
1130 ! ! add code to get the table header name and table values
1131 ! if (len_trim(sArgText_2) > 0) &
1132 ! call PARAMS%add_file( fix_pathname(trim(sPathname)//sArgText_2))
1133
1134 else
1135
1136 call warn( "Did not find a valid "//dquote(skey)//" option. Value supplied was: "//dquote(sargtext_1), &
1137 lfatal = true, shints="Valid options include "//dquote("CONSTANT")//", " &
1138 //dquote("ARC_ASCII")//", "//dquote("ARC_GRID") &
1139 //", "//dquote("SURFER")//", or "//dquote("NETCDF"), &
1140 lecho = true )
1141
1142 endif
1143
1144 elseif ( scmdtext .containssimilar. "_USE_MAJORITY_FILTER" ) then
1145
1146 call pentry%set_majority_filter_flag( true )
1147
1148 elseif ( scmdtext .containssimilar. "_MONTHNAMES_CAPITALIZED" ) then
1149
1150 pentry%iFilename_Monthname_Capitalization_Rule = file_template_capitalized_monthname
1151
1152 elseif ( scmdtext .containssimilar. "_MONTHNAMES_LOWERCASE" ) then
1153
1154 pentry%iFilename_Monthname_Capitalization_Rule = file_template_lowercase_monthname
1155
1156 elseif ( scmdtext .containssimilar. "_MONTHNAMES_UPPERCASE" ) then
1157
1158 pentry%iFilename_Monthname_Capitalization_Rule = file_template_uppercase_monthname
1159
1160 elseif ( scmdtext .containssimilar. "_CONVERSION_FACTOR" ) then
1161
1162 call pentry%set_scale(asdouble(sargtext_1))
1163
1164 elseif ( scmdtext .containssimilar. "NETCDF_X_VAR_ADD_OFFSET" ) then
1165
1166 call pentry%set_X_offset( asdouble( sargtext_1 ) )
1167
1168 elseif ( scmdtext .containssimilar. "NETCDF_Y_VAR_ADD_OFFSET" ) then
1169
1170 call pentry%set_Y_offset( asdouble( sargtext_1 ) )
1171
1172 elseif ( scmdtext .containssimilar. "_SCALE_FACTOR" ) then
1173
1174 call pentry%set_scale(asdouble(sargtext_1))
1175
1176 elseif ( scmdtext .containssimilar. "_ADD_OFFSET" ) then
1177
1178 call pentry%set_add_offset(asdouble(sargtext_1))
1179
1180 elseif ( scmdtext .containssimilar. "_SUBTRACT_OFFSET" ) then
1181
1182 call pentry%set_sub_offset(asdouble(sargtext_1))
1183
1184 elseif ( scmdtext .containssimilar. "_UNITS_KELVIN" ) then
1185
1186 call pentry%set_sub_offset(freezing_point_of_water_kelvin)
1187 call pentry%set_add_offset(freezing_point_of_water_fahrenheit)
1188 call pentry%set_scale(f_per_c)
1189
1190 elseif ( scmdtext .containssimilar. "_UNITS_CELSIUS" ) then
1191
1192 call pentry%set_add_offset(freezing_point_of_water_fahrenheit)
1193 call pentry%set_scale(f_per_c)
1194
1195 elseif ( scmdtext .containssimilar. "_UNITS_MILLIMETERS" ) then
1196
1197 call pentry%set_scale(1.0_c_double / mm_per_in)
1198
1199 elseif ( scmdtext .containssimilar. "_COORDINATE_TOLERANCE" ) then
1200
1201 call pentry%set_coordinate_tolerance( asdouble(sargtext_1))
1202
1203 elseif ( scmdtext .containssimilar. "NETCDF_X_VAR" ) then
1204
1205 pentry%sVariableName_x = trim(sargtext_1)
1206
1207 elseif ( scmdtext .containssimilar. "NETCDF_Y_VAR" ) then
1208
1209 pentry%sVariableName_y = trim(sargtext_1)
1210
1211 elseif ( (scmdtext .containssimilar. "NETCDF_Z_VAR") &
1212 .or. (scmdtext .containssimilar. "COLUMN_NAME") ) then
1213
1214 pentry%sVariableName_z = trim(sargtext_1)
1215
1216 elseif ( (scmdtext .containssimilar. "NETCDF_TIME_VAR") &
1217 .or. (scmdtext .containssimilar. "DATE_COLUMN_NAME") ) then
1218
1219 pentry%sVariableName_time = trim(sargtext_1)
1220
1221 elseif ( scmdtext .containssimilar. "NETCDF_VARIABLE_ORDER" ) then
1222
1223 call pentry%set_variable_order( aslowercase(sargtext_1) )
1224
1225 elseif ( scmdtext .containssimilar. "NETCDF_FLIP_VERTICAL") then
1226
1227 call pentry%set_grid_flip_vertical()
1228
1229 elseif ( scmdtext .containssimilar. "NETCDF_FLIP_HORIZONTAL" ) then
1230
1231 call pentry%set_grid_flip_horizontal()
1232
1233 elseif ( scmdtext .containssimilar."NETCDF_NO_AUTOMATIC_GRID_FLIPPING" ) then
1234
1235 call pentry%do_not_allow_netcdf_grid_data_flipping()
1236
1237 elseif ( scmdtext .containssimilar. "ALLOW_MISSING_FILES" ) then
1238
1239 call pentry%allow_missing_files()
1240
1241 elseif ( scmdtext .containssimilar. "NETCDF_MAKE_LOCAL_ARCHIVE" ) then
1242
1243 call pentry%set_make_local_archive(true)
1244
1245 elseif ( scmdtext .containssimilar. "_PROJECTION_DEFINITION" ) then
1246
1247 call pentry%set_source_PROJ4( trim(sargtext) )
1248
1249 elseif ( scmdtext .containssimilar. "_DATE_COLUMN" ) then
1250
1251 pentry%sDateColumnName = trim( sargtext_1 )
1252
1253 elseif ( scmdtext .containssimilar. "_VALUE_COLUMN" ) then
1254
1255 pentry%sValueColumnName = trim( sargtext_1 )
1256
1257 elseif ( scmdtext .containssimilar. "_MINIMUM_ALLOWED_VALUE" ) then
1258
1259 pentry%rMinAllowedValue = asfloat(sargtext_1)
1260
1261 elseif ( scmdtext .containssimilar. "_MAXIMUM_ALLOWED_VALUE" ) then
1262
1263 pentry%rMaxAllowedValue = asfloat(sargtext_1)
1264
1265 elseif ( scmdtext .containssimilar. "_MISSING_VALUES_CODE" ) then
1266
1267 pentry%rMissingValuesCode = asfloat(sargtext_1)
1268
1269 elseif ( scmdtext .containssimilar. "_MISSING_VALUES_OPERATOR" ) then
1270
1271 pentry%sMissingValuesOperator = trim(sargtext_1)
1272
1273 elseif ( scmdtext .containssimilar. "_MISSING_VALUES_ACTION" ) then
1274
1275 if (sargtext_1 .strapprox. "ZERO") then
1276
1277 pentry%iMissingValuesAction = missing_values_zero_out
1278
1279 elseif (sargtext_1 .strapprox. "MEAN" ) then
1280
1281 pentry%iMissingValuesAction = missing_values_replace_with_mean
1282
1283 else
1284
1285 call warn("Unknown missing value action supplied for " &
1286 //dquote(skey)//" data: "//dquote(sargtext_1) )
1287
1288 endif
1289
1290 elseif ( scmdtext .containssimilar. "_METHOD" ) then
1291
1292 ! no operation; just keep SWB quiet about METHOD and it will be included in the
1293 ! methods initialization section.
1294
1295 elseif ( scmdtext .containssimilar. "_LOOKUP_TABLE" ) then
1296
1297 ! no operation; just keep SWB quiet about LOOKUP_TABLE and it will be included in the
1298 ! lookup table initialization section.
1299
1300 else
1301
1302 call warn("Unknown directive detected in code at line "//ascharacter(__line__)//", file "//__file__ &
1303 //". ~Ignoring. Directive is: "//dquote(scmdtext), iloglevel=log_debug )
1304
1305 endif
1306
1307 enddo
1308
1309 ! if an unadorned grid specification directive was processed, then we can add the key and
1310 ! the data_catalog_entry to the data_catalog
1311 if ( lgridpresent ) call dat%add( key=skey, data=pentry )
1312
1313 pentry => null()
1314
1315 endif
1316
1317 call mydirectives%clear()
1318 call myoptions%clear()
1319
1320 end subroutine initialize_generic_grid
1321
1322!--------------------------------------------------------------------------------------------------
1323
1325
1326 ! [ LOCALS ]
1327 type (FSTRING_LIST_T) :: myOptions
1328 character (len=:), allocatable :: sArgText
1329 real (c_double) :: rX0, rX1, rY0, rY1, rGridCellSize
1330 integer (c_int) :: iNX, iNY
1331 real (c_float) :: fTempVal
1332
1333 rx1 = 0.0_c_double
1334 ry1 = 0.0_c_double
1335 rgridcellsize = 0.0_c_double
1336
1337 ! For MODEL directive, obtain the associated dictionary entries
1338 call cf_dict%get_values( "GRID", myoptions )
1339
1340 ! dictionary entries are initially space-delimited; sArgText contains
1341 ! all dictionary entries present, concatenated, with a space between entries
1342 sargtext = myoptions%get(1, myoptions%count )
1343
1344 ! echo the original directive and dictionary entries to the logfile
1345 call logs%write("> GRID "//sargtext, ilinesbefore=1 )
1346
1347 inx = asint( myoptions%get(1) )
1348 iny = asint( myoptions%get(2) )
1349 rx0 = asdouble( myoptions%get(3) )
1350 ry0 = asdouble( myoptions%get(4) )
1351
1352 if ( myoptions%count == 5 ) then
1353
1354 rgridcellsize = asdouble( myoptions%get(5) )
1355
1356 call model%initialize_grid(inx, iny, rx0, ry0, rgridcellsize)
1357
1358 rx1 = rx0 + rgridcellsize * real(inx, c_double)
1359 ry1 = ry0 + rgridcellsize * real(iny, c_double)
1360
1361 elseif ( myoptions%count == 7 ) then
1362
1363 rx1 = asdouble( myoptions%get(5) )
1364 ry1 = asdouble( myoptions%get(6) )
1365 rgridcellsize = asdouble( myoptions%get(7) )
1366
1367 ftempval = real(( rx1 - rx0 ) / real(inx, c_double), c_float)
1368
1369 call model%initialize_grid(inx, iny, rx0, ry0, rgridcellsize)
1370
1371 else
1372
1373 call warn("Grid specification is flawed or missing.", lfatal=true, iloglevel = log_all, lecho = true )
1374
1375 endif
1376
1377 call myoptions%clear()
1378
1379 ! For MODEL directive, obtain the associated dictionary entries
1380 call cf_dict%get_values( "BASE_PROJECTION_DEFINITION", myoptions )
1381
1382 if ( myoptions%get(1) .strequal. "<NA>" ) then
1383 call die(smessage="Your control file is missing a BASE_PROJECTION_DEFINITION entry.", &
1384 shints="This version of SWB requires that you add a BASE_PROJECTION_DEFINITION entry " &
1385 //"to your control file." )
1386 endif
1387
1388 ! dictionary entries are initially space-delimited; sArgText contains
1389 ! all dictionary entries present, concatenated, with a space between entries
1390 sargtext = myoptions%get(1, myoptions%count )
1391
1392 ! echo the original directive and dictionary entries to the logfile
1393 call logs%write("> BASE_PROJECTION_DEFINITION "//sargtext, ilinesbefore=1)
1394
1395 ! BNDS is a module-level data structure that will be used in other modules to
1396 ! supply bounding box information for the SWB project area
1397 bnds%iNumCols = inx
1398 bnds%iNumRows = iny
1399 bnds%fX_ll = rx0
1400 bnds%fY_ll = ry0
1401 bnds%fY_ur = ry1
1402 bnds%fX_ur = rx1
1403 bnds%fGridCellSize = rgridcellsize
1404 bnds%sPROJ4_string = trim(sargtext)
1405
1406 model%PROJ4_string = trim(sargtext)
1407
1408 end subroutine initialize_grid_options
1409
1410!--------------------------------------------------------------------------------------------------
1411
1413
1414 ! [ LOCALS ]
1415 type (FSTRING_LIST_T) :: myOptions
1416 integer (c_int) :: iIndex
1417 integer (c_int) :: jIndex
1418 character (len=:), allocatable :: sArgText
1419 character (len=:), allocatable :: sOutput
1420 logical (c_bool) :: enable_output
1421
1422 enable_output = true
1423
1424 ! For MODEL directive, obtain the associated dictionary entries
1425 call cf_dict%get_values( "OUTPUT", myoptions )
1426
1427 ! dictionary entry will look like this:
1428 ! 1) ENABLE
1429 ! 2) net_infiltration
1430 ! 3) actual_evapotranspiration
1431 ! etc.
1432
1433 ! in the logic that follows, it is assumed that the first entry in the
1434 ! set of returned dictionary values will be an action verb ('ENABLE', 'DISABLE')
1435 ! and the remaining values will be the variables to enable or disable
1436
1437 ! dictionary entries are initially space-delimited; sArgText contains
1438 ! all dictionary entries present, concatenated, with a space between entries
1439 sargtext = myoptions%get(1, myoptions%count )
1440
1441 ! echo the original directive and dictionary entries to the logfile
1442 call logs%write("> OUTPUT"//sargtext, ilinesbefore=1 )
1443
1444 ! iterate over all dictionary entries
1445 do iindex=1,myoptions%count
1446
1447 soutput = myoptions%get( iindex )
1448
1449 if ( ( soutput .strapprox. "ENABLE") .or. ( soutput .strapprox. "ACTIVE") ) then
1450
1451 enable_output = true
1452
1453 elseif ( ( soutput .strapprox. "DISABLE") .or. ( soutput .strapprox. "INACTIVE") ) then
1454
1455 enable_output = false
1456
1457 endif
1458
1459 if ( soutput .strapprox. "ALL" ) then
1460
1461 outspecs(:)%is_active = enable_output
1462
1463 else
1464
1465 do jindex=lbound( outspecs, 1), ubound( outspecs, 1)
1466
1467 if ( outspecs( jindex )%variable_name .strapprox. soutput ) then
1468 outspecs( jindex )%is_active = enable_output
1469 if ( enable_output ) call logs%write("> Enabling output for "//squote(soutput), ilinesbefore=1 )
1470 if ( .not. enable_output ) call logs%write("> Disabling output for "//squote(soutput), ilinesbefore=1 )
1471 endif
1472 enddo
1473
1474 endif
1475
1476 enddo
1477
1478 call myoptions%clear()
1479
1480 end subroutine initialize_output_options
1481
1482!--------------------------------------------------------------------------------------------------
1483
1485
1486 ! [ LOCALS ]
1487 type (FSTRING_LIST_T) :: myDirectives
1488 type (FSTRING_LIST_T) :: myOptions
1489 integer (c_int) :: iIndex
1490 character (len=:), allocatable :: sCmdText
1491 character (len=:), allocatable :: sOptionText
1492 character (len=:), allocatable :: sArgText
1493 logical (c_bool) :: lHaveStartDate
1494 logical (c_bool) :: lHaveEndDate
1495
1496 lhavestartdate = false
1497 lhaveenddate = false
1498
1499 mydirectives = cf_dict%grep_keys("DATE")
1500
1501! if ( myDirectives%count < 2 ) then
1502
1503! call warn(sMessage="Your control file seems to be missing START_DATE and/or END_DATE", &
1504! sHints="Add a START_DATE and/or END_DATE directive to your control file. Date " &
1505! //"~should be specified as mm/dd/yyyy.", lFatal = TRUE, iLogLevel = LOG_ALL, &
1506! lEcho = TRUE )
1507
1508
1509 call logs%set_loglevel( log_all )
1510 call logs%set_echo( false )
1511
1512 do iindex = 1, mydirectives%count
1513
1514 ! myDirectives is a string list of all SWB directives that contain the string "DATE"
1515 ! sCmdText contains an individual directive
1516 scmdtext = mydirectives%get(iindex)
1517
1518 ! For MODEL directive, obtain the associated dictionary entries
1519 call cf_dict%get_values(scmdtext, myoptions )
1520
1521 ! dictionary entries are initially space-delimited; sArgText contains
1522 ! all dictionary entries present, concatenated, with a space between entries
1523 sargtext = myoptions%get(1, myoptions%count )
1524
1525 ! echo the original directive and dictionary entries to the logfile
1526 call logs%write("> "//scmdtext//" "//sargtext, ilinesbefore=1 )
1527
1528 ! most of the time, we only care about the first dictionary entry, obtained below
1529 soptiontext = myoptions%get(1)
1530
1531 select case ( scmdtext )
1532
1533 case ( "START_DATE", "STARTDATE", "BEGIN_DATE" )
1534
1535 lhavestartdate = true
1536 call sim_dt%start%parseDate( soptiontext )
1537 call sim_dt%start%calcJulianDay()
1538
1539 case ( "END_DATE", "ENDDATE", "STOP_DATE" )
1540
1541 lhaveenddate = true
1542 call sim_dt%end%parseDate( soptiontext )
1543 call sim_dt%end%calcJulianDay()
1544
1545 case default
1546
1547 call warn("Unknown directive present, line "//ascharacter(__line__)//", file "//__file__ &
1548 //". Ignoring. Directive is: "//dquote(scmdtext), iloglevel=log_debug )
1549
1550 end select
1551
1552 enddo
1553
1554 if ( lhavestartdate .and. lhaveenddate ) then
1555
1556 sim_dt%curr = sim_dt%start
1557 sim_dt%iDOY = day_of_year( sim_dt%curr%getJulianDay() )
1558
1559 sim_dt%iDaysInMonth = sim_dt%curr%dayspermonth()
1560 sim_dt%iDaysInYear = sim_dt%curr%daysperyear()
1561 sim_dt%lIsLeapYear = sim_dt%curr%isLeapYear()
1562
1563 call logs%write("Model run start date set to: "//sim_dt%start%prettydate(), itab=4, lecho=true)
1564 call logs%write("Model run end date set to: "//sim_dt%end%prettydate(), itab=4, lecho=true)
1565
1566 else
1567
1568 call warn(smessage="Your control file seems to be missing START_DATE and/or END_DATE", &
1569 shints="Add a START_DATE and/or END_DATE directive to your control file. Date " &
1570 //"~should be specified as mm/dd/yyyy.", lfatal = true, iloglevel = log_all, &
1571 lecho = true )
1572
1573 endif
1574
1575 end subroutine initialize_start_and_end_dates
1576
1577!--------------------------------------------------------------------------------------------------
1578
1579 !> Find any parameter tables specified in the control file; process and store contents.
1580 !!
1581 !! Any control file entry that contains the text "LOOKUP_TABLE" is assumed to specify a parameter table
1582 !! that needs to be read in and processed.
1583 !!
1584 !! @note The entries given in the files are implicitly assumed to be in sorted order. For example, if the parameters
1585 !! pertain to landuse codes, it is assumed that all table values are given in order from lowest to highest
1586 !! landuse code.
1587
1589
1590 ! [ LOCALS ]
1591 type (FSTRING_LIST_T) :: myDirectives
1592 type (FSTRING_LIST_T) :: myOptions
1593 integer (c_int) :: iIndex
1594 character (len=:), allocatable :: sCmdText
1595 character (len=:), allocatable :: sOptionText
1596 character (len=:), allocatable :: sArgText
1597 type (PARAMETERS_T) :: PARAMS_LU_TABLE
1598 integer (c_int) :: iCount
1599
1600 icount = 0
1601
1602 mydirectives = cf_dict%grep_keys("LOOKUP_TABLE")
1603
1604 if ( mydirectives%count == 0 ) then
1605
1606 call warn("Your control file seems to be missing the required lookup table(s).", &
1607 lfatal = true )
1608
1609 else
1610
1611 call logs%set_loglevel( log_all )
1612 call logs%set_echo( false )
1613
1614 ! iterate over list of all control file statements that contain the phrase "LOOKUP_TABLE"
1615 do iindex = 1, mydirectives%count
1616
1617 ! sCmdText contains an individual directive
1618 scmdtext = mydirectives%get(iindex)
1619
1620 ! for MODEL directive, obtain the associated dictionary entries
1621 call cf_dict%get_values(scmdtext, myoptions )
1622
1623 ! dictionary entries are initially space-delimited; sArgText contains
1624 ! all dictionary entries present, concatenated, with a space between entries
1625 sargtext = myoptions%get(1, myoptions%count )
1626
1627 ! echo the original directive and dictionary entries to the logfile
1628 call logs%write("> "//scmdtext//" "//sargtext, ilinesbefore=1 )
1629
1630 ! most of the time, we only care about the first dictionary entry, obtained below
1631 soptiontext = fix_pathname( myoptions%get(1) )
1632 if (allocated(lookup_table_directory_name)) then
1633 if (len_trim(lookup_table_directory_name) > 0) &
1634 soptiontext = trim(lookup_table_directory_name)//trim(soptiontext)
1635 endif
1636
1637 if ( index(string=scmdtext, substring="LOOKUP_TABLE" ) > 0 ) then
1638
1639 call params_lu_table%add_file( fix_pathname( soptiontext ))
1640 icount = icount + 1
1641
1642 else
1643
1644 call warn("Unknown directive present, line "//ascharacter(__line__)//", file "//__file__ &
1645 //". Ignoring. Directive is: "//dquote(scmdtext), iloglevel=log_debug )
1646
1647 endif
1648
1649 enddo
1650
1651 if ( icount > 0 ) then
1652
1653 call params_lu_table%munge_file(delimiters=tab)
1654 call params_dict%print_all(sdescription="LOOKUP TABLE dictionary", &
1655 iloglevel=log_debug, lecho=false)
1656
1657 endif
1658
1659 endif
1660
1661 end subroutine initialize_lookup_tables
1662
1663!--------------------------------------------------------------------------------------------------
1664
1665 subroutine initialize_generic_method( sKey, lOptional)
1666
1667 character (len=*), intent(in) :: sKey
1668 logical (c_bool), intent(in) :: lOptional
1669
1670 ! [ LOCALS ]
1671 type (FSTRING_LIST_T) :: myDirectives
1672 type (FSTRING_LIST_T) :: myOptions
1673 integer (c_int) :: iIndex
1674 integer (c_int) :: indx
1675 character (len=:), allocatable :: sCmdText
1676! character (len=:), allocatable :: sOptionText
1677 type (FSTRING_LIST_T) :: argv_list
1678 character (len=:), allocatable :: sArgText
1679 logical (c_bool) :: lFatal
1680 integer (c_int) :: num_elements
1681
1682 ! obtain a list of control file directives whose key values contain the string sKey
1683 mydirectives = cf_dict%grep_keys( trim(skey) )
1684
1685 lfatal = .not. loptional
1686
1687 if ( mydirectives%count == 0 ) then
1688
1689 call warn("Your control file is missing any of the required directives relating to "//dquote(skey)//" method.", &
1690 lfatal = lfatal, iloglevel = log_all, lecho = true )
1691
1692 else
1693
1694 call logs%set_loglevel( log_all )
1695 call logs%set_echo( false )
1696
1697 ! repeat MODEL process for each control file directive in list
1698 do iindex = 1, mydirectives%count
1699
1700 ! sCmdText contains an individual directive (e.g. TMAX NETCDF input/tmax_%0m_%0d_%0Y.nc )
1701 scmdtext = mydirectives%get(iindex)
1702
1703 ! For MODEL directive, obtain the associated dictionary entries
1704 call cf_dict%get_values(scmdtext, myoptions )
1705
1706 ! dictionary entries are initially space-delimited; sArgText contains
1707 ! all dictionary entries present, concatenated, with a space between entries
1708 sargtext = myoptions%get(1, myoptions%count )
1709
1710 ! echo the original directive and dictionary entries to the logfile
1711 call logs%write("> "//scmdtext//" "//sargtext, ilinesbefore=1, iloglevel=log_all, lecho=false )
1712
1713 ! most of the time, we only care about the first dictionary entry, obtained below
1714! sOptionText = myOptions%get(1)
1715
1716 num_elements = myoptions%count
1717
1718 call argv_list%clear()
1719
1720 do indx=1, myoptions%count
1721 call argv_list%append( myoptions%get( indx ) )
1722 enddo
1723
1724 ! Any entry in the control file that contains the substring "METHOD" will be
1725 ! handed to the "set_method_pointers" subroutine in an attempt to wire up the correct
1726 ! process modules
1727 if ( ( scmdtext .contains. "METHOD" ) .or. ( scmdtext .contains. "DUMP" ) ) then
1728
1729 call model%set_method_pointers( trim(scmdtext), argv_list )
1730
1731 endif
1732
1733 enddo
1734
1735 endif
1736
1737 end subroutine initialize_generic_method
1738
1739!--------------------------------------------------------------------------------------------------
1740
1742
1743 ! [ LOCALS ]
1744 type (FSTRING_LIST_T) :: myDirectives
1745 type (FSTRING_LIST_T) :: myOptions
1746 integer (c_int) :: iIndex
1747 integer (c_int) :: indx
1748 character (len=:), allocatable :: sCmdText
1749! character (len=:), allocatable :: sOptionText
1750 type (FSTRING_LIST_T) :: argv_list
1751 character (len=:), allocatable :: sArgText
1752 integer (c_int) :: num_elements
1753 character (len=:), allocatable :: Option_Name
1754
1755 ! obtain a list of control file directives whose key values contain the string sKey
1756 mydirectives = cf_dict%grep_keys( "OPTION" )
1757
1758 call logs%set_loglevel( log_all )
1759 call logs%set_echo( false )
1760
1761 ! repeat MODEL process for each control file directive in list
1762 do iindex = 1, mydirectives%count
1763
1764 ! sCmdText contains an individual directive (e.g. TMAX NETCDF input/tmax_%0m_%0d_%0Y.nc )
1765 scmdtext = mydirectives%get(iindex)
1766
1767 ! For MODEL directive, obtain the associated dictionary entries
1768 call cf_dict%get_values(scmdtext, myoptions )
1769
1770 ! dictionary entries are initially space-delimited; sArgText contains
1771 ! all dictionary entries present, concatenated, with a space between entries
1772 sargtext = myoptions%get(1, myoptions%count )
1773
1774 ! echo the original directive and dictionary entries to the logfile
1775 call logs%write("> "//scmdtext//" "//sargtext, ilinesbefore=1 )
1776
1777 ! most of the time, we only care about the first dictionary entry, obtained below
1778! sOptionText = myOptions%get(1)
1779
1780 num_elements = myoptions%count
1781
1782 call argv_list%clear()
1783
1784 do indx=1, myoptions%count
1785 call argv_list%append( myoptions%get( indx ) )
1786 enddo
1787
1788 ! Any entry in the control file that contains the substring "OPTION" will be
1789 ! processed directly
1790 if ( ( scmdtext .contains. "OPTION" ) ) then
1791
1792 option_name = argv_list%get(1)
1793
1794 select case ( option_name )
1795
1796 case ( "NO_LATLON_IN_OUTPUT" )
1797
1799 call logs%write("==> LATITUDE and LONGITUDE will *not* be included in output NetCDF files.", ilinesbefore=1 )
1800
1801 end select
1802
1803 endif
1804
1805 enddo
1806
1807 end subroutine initialize_program_options
1808
1809
1810!--------------------------------------------------------------------------------------------------
1811
1812
1814
1815 ! [ LOCALS ]
1816
1817 pcoord_grd => grid_create( inx=model%number_of_columns, iny=model%number_of_rows, &
1818 rx0=model%X_ll, ry0=model%Y_ll, &
1819 rgridcellsize=model%gridcellsize, idatatype=grid_datatype_real )
1820
1821 allocate ( model%X(model%number_of_columns ) )
1822 allocate ( model%Y(model%number_of_rows ) )
1823
1824 ! call the grid routine to populate the X and Y values
1826
1827 ! populating these in order to have them available later for use in writing results to NetCDF
1828 model%X = pcoord_grd%rX( :, 1 )
1829 model%Y = pcoord_grd%rY( 1, : )
1830
1831 ! transform to unprojected (lat/lon) coordinate system
1832 call grid_transform(pgrd=pcoord_grd, sfromproj4=model%PROJ4_string, &
1833 stoproj4="+proj=lonlat +ellps=GRS80 +datum=WGS84 +no_defs" )
1834
1835 model%latitude = real(pack( pcoord_grd%rY, model%active ), c_float)
1836
1837 model%X_lon = pcoord_grd%rX
1838 model%Y_lat = pcoord_grd%rY
1839
1840 pcoord_grd%rData=real(pcoord_grd%rX, c_float)
1841 call grid_writearcgrid( sfilename="Longitude__calculated.asc", pgrd=pcoord_grd )
1842 pcoord_grd%rData=real(pcoord_grd%rY, c_float)
1843 call grid_writearcgrid( sfilename="Latitude__calculated.asc", pgrd=pcoord_grd )
1844
1845 call grid_destroy( pcoord_grd )
1846
1847 end subroutine initialize_latitude
1848
1849! !--------------------------------------------------------------------------------------------------
1850!
1851! !> Match landuse codes from table with those contained in the gridded landuse.
1852! !!
1853! !! This routine loops through all known
1854!
1855! subroutine initialize_landuse_codes()
1856!
1857! ! [ LOCALS ]
1858! integer (c_int) :: iIndex
1859! integer (c_int), allocatable :: iLandUseCodes(:)
1860! type (DATA_CATALOG_ENTRY_T), pointer :: pLULC
1861! integer (c_int) :: iIndex2
1862! integer (c_int) :: iCount
1863! integer (c_int) :: iStat
1864! logical (c_bool) :: lMatch
1865! type (FSTRING_LIST_T) :: slList
1866!
1867! call slList%append("LU_Code")
1868! call slList%append("LU_code")
1869! call slList%append("Landuse_Code")
1870! call slList%append("LULC_Code")
1871!
1872! !> Determine how many landuse codes are present
1873! call PARAMS%get_parameters( slKeys=slList, iValues=iLanduseCodes, lFatal=TRUE )
1874!
1875! ! obtain a pointer to the LAND_USE grid
1876! pLULC => DAT%find("LAND_USE")
1877!
1878! if ( associated(pLULC) ) then
1879!
1880! if (associated( pLULC%pGrdBase) ) then
1881! MODEL%landuse_code = pack( pLULC%pGrdBase%iData, MODEL%active )
1882! else
1883! call die("INTERNAL PROGRAMMING ERROR: attempted use of NULL pointer", __FILE__, __LINE__)
1884! endif
1885! else
1886! call die("Attempted use of NULL pointer. Failed to find LAND_USE data element.", &
1887! __FILE__, __LINE__)
1888! endif
1889!
1890! ! setting this to a value that is likely valid; if this is set to a negative value, a landuse
1891! ! code that is present in the grid but not in the lookup table will trigger a fatal error, however,
1892! ! processing will continue until a bounds error is triggered somewhere else in the code,
1893! MODEL%landuse_index = -9999
1894! iCount = 0
1895!
1896! ! only run through matching process if we have found a LU_Code entry in the
1897! ! parameter dictionary
1898! if ( all( iLandUseCodes >= 0 ) ) then
1899!
1900! do iIndex = 1, ubound(MODEL%landuse_code,1)
1901!
1902! lMatch = FALSE
1903!
1904! do iIndex2=1, ubound(iLandUseCodes, 1)
1905!
1906! if (MODEL%landuse_code(iIndex) == iLandUseCodes(iIndex2) ) then
1907! MODEL%landuse_index(iIndex) = iIndex2
1908! iCount = iCount + 1
1909! lMatch = TRUE
1910! exit
1911! endif
1912!
1913! enddo
1914!
1915! if ( .not. lMatch ) then
1916! call warn(sMessage="Failed to match landuse code "//asCharacter(MODEL%landuse_code(iIndex) ) &
1917! //" with a corresponding landuse code from lookup tables.", &
1918! sHints="Make sure your lookup table(s) have landuse codes corresponding to all values in " &
1919! //"the land-use grid.", lFatal=TRUE, iLogLevel=LOG_ALL, lEcho=TRUE)
1920! ! we are setting this value to a valid value. this should not cause problems with any
1921! ! calculations because we have already thrown a fatal error
1922! MODEL%landuse_index(iIndex) = 1
1923! endif
1924! enddo
1925!
1926! call LOGS%write("Matches were found between landuse grid value and table value for " &
1927! //asCharacter(iCount)//" cells out of a total of "//asCharacter(ubound(MODEL%landuse_code,1))//" active cells.", &
1928! iLinesBefore=1, iLinesAfter=1, iLogLevel=LOG_ALL)
1929!
1930! call slList%clear()
1931!
1932! endif
1933!
1934! ! if we have more than one cell for which an index value could not be assigned, trigger fatal error
1935! if ( count(MODEL%landuse_index < 0) > 0 ) then
1936! call warn(asCharacter(count(MODEL%landuse_index < 0))//" landuse codes could not be " &
1937! //" assigned a landuse index value.", lFatal=TRUE, sHints="Make sure that you have an " &
1938! //"entry in the landuse lookup table for each unique code contained in your landuse grid." )
1939! endif
1940!
1941!
1942! end subroutine initialize_landuse_codes
1943
1944!--------------------------------------------------------------------------------------------------
1945
1947
1948 integer (c_int) :: iIndex
1949 type (FSTRING_LIST_T) :: slList
1950 integer( c_int), allocatable :: iLanduseTableCodes(:)
1951 integer (c_int) :: iNumberOfLanduses
1952 real (c_float), allocatable :: SURFACE_STORAGE_MAXIMUM(:)
1953 real (c_float) :: current_surface_storage_max
1954
1955
1956 ! create list of possible table headings to look for...
1957 call sllist%append( "LU_Code" )
1958 call sllist%append( "Landuse_Lookup_Code" )
1959
1960 !> Determine how many landuse codes are present
1961 call params%get_parameters( slkeys=sllist, ivalues=ilandusetablecodes )
1962 inumberoflanduses = count( ilandusetablecodes >= 0 )
1963
1964 call sllist%clear()
1965 call sllist%append("Surface_Storage_Max")
1966 call sllist%append("Surface_Storage_Maximum")
1967
1968 model%surface_storage_max = 0.0_c_float
1969
1970 call params%get_parameters( slkeys=sllist, fvalues=surface_storage_maximum, lfatal=false )
1971
1972 if ( all( surface_storage_maximum > rtinyval ) ) then
1973
1974 do iindex=1, ubound( surface_storage_maximum, 1)
1975
1976 current_surface_storage_max = surface_storage_maximum( iindex )
1977
1978 where( model%landuse_index == iindex )
1979
1980 model%surface_storage_max = current_surface_storage_max
1981
1982 end where
1983
1984 enddo
1985
1986 endif
1987
1988 end subroutine initialize_surface_storage_max
1989
1990end module model_initialize
This module contains physical constants and convenience functions aimed at performing unit conversion...
character(len=len_trim(input_pathname)) function fix_pathname(input_pathname)
logical(c_bool), parameter, public true
real(c_double), parameter, public f_per_c
real(c_float), parameter, public rtinyval
real(c_double), parameter, public mm_per_in
character(len=:), allocatable, public lookup_table_directory_name
logical(c_bool), parameter, public false
real(c_double), parameter, public freezing_point_of_water_kelvin
real(c_double), parameter, public freezing_point_of_water_fahrenheit
integer(c_int), parameter datatype_float
impure elemental logical(c_bool) function, public is_numeric(value)
Determine if string contains numeric values.
character(len=:), allocatable, public data_directory_name
integer(c_int), parameter datatype_int
Defines the DATA_CATALOG_T data type, which contains type-bound procedures to add,...
type(data_catalog_t), public dat
DAT is a global to hold data catalog entries.
type(general_grid_t), pointer, public pgrd
integer(c_int), parameter, public missing_values_replace_with_mean
integer(c_int), parameter, public file_template_uppercase_monthname
integer(c_int), parameter, public missing_values_zero_out
integer(c_int), parameter, public file_template_capitalized_monthname
integer(c_int), parameter, public file_template_lowercase_monthname
This module contains the DATETIME_T class and associated time and date-related routines,...
Definition datetime.F90:9
integer(c_int) function, public day_of_year(ijd)
type(dict_entry_t), pointer, public cf_entry
type(dict_t), public cf_dict
subroutine, public check_for_fatal_warnings()
subroutine, public warn(smessage, smodule, iline, shints, lfatal, iloglevel, lecho)
subroutine, public die(smessage, smodule, iline, shints, scalledby, icalledbyline)
character(len=1), parameter, public tab
Definition fstring.F90:167
Provides support for input and output of gridded ASCII data, as well as for creation and destruction ...
Definition grid.F90:8
subroutine, public grid_populatexy(pgrd, rx, ry)
Definition grid.F90:2195
subroutine, public grid_destroy(pgrd)
Definition grid.F90:366
subroutine, public grid_writearcgrid(sfilename, pgrd)
Definition grid.F90:1054
integer(c_int), parameter, public grid_datatype_real
Definition grid.F90:26
subroutine, public grid_set_output_directory_name(sdirname)
Definition grid.F90:182
subroutine, public grid_transform(pgrd, sfromproj4, stoproj4, rx, ry)
Call PROJ4 to transform coordinates.
Definition grid.F90:1538
type(logfile_t), public logs
Definition logfiles.F90:62
subroutine, public initialize_landuse_codes()
Match landuse codes from table with those contained in the gridded landuse.
type(model_domain_t), public model
subroutine, public read_landuse_codes
real(c_float), dimension(:,:), allocatable, public rooting_depth_max
subroutine initialize_generic_method(skey, loptional)
type(methods_list_t), dimension(number_of_known_methods), parameter known_methods
type(general_grid_t), pointer pcoord_grd
subroutine initialize_grid_options()
subroutine initialize_ancillary_values()
subroutine, public write_control_file(sfilename, sgridspecification, slextradirectives)
subroutine, public read_control_file(sfilename)
subroutine initialize_latitude()
subroutine initialize_program_options()
subroutine, public initialize_all(output_prefix, output_dirname, data_dirname, lookup_table_dirname, weather_data_dirname)
subroutine initialize_lookup_tables()
Find any parameter tables specified in the control file; process and store contents.
subroutine initialize_snow_storage()
subroutine initialize_percent_pervious()
subroutine initialize_output_options()
integer(c_int), parameter number_of_known_methods
subroutine initialize_soil_storage()
type(gridded_datasets_t), dimension(number_of_known_grids) known_grids
subroutine set_data_directory(data_dirname)
subroutine initialize_start_and_end_dates()
subroutine set_weather_data_directory(weather_data_dirname)
subroutine read_hydrologic_soil_groups
subroutine initialize_surface_storage_max()
subroutine initialize_hydrologic_soil_groups
integer(c_int), parameter number_of_known_grids
subroutine initialize_percent_canopy_cover
subroutine initialize_soils_landuse_awc_flowdir_values()
Initialize soils, landuse, and available water content values.
subroutine initialize_generic_grid(skey, spathname, loptional, idatatype)
Generic routine to handle intake of gridded data.
subroutine set_lookup_table_directory(lookup_table_dirname)
subroutine read_polygon_id()
type(output_specs_t), dimension(ncdf_num_outputs) outspecs
Definition output.F90:30
subroutine initialize_output(cells)
Definition output.F90:117
subroutine set_output_prefix(output_prefix)
Definition output.F90:107
subroutine set_output_directory(output_dir_name)
Definition output.F90:97
subroutine set_output_latlon_option(output_includes_latlon_l)
Definition output.F90:87
type(parameters_t), public params
type(dict_t), public params_dict
Module precipitation__method_of_fragments provides support for creating synthetic daily precipitation...
type(date_range_t), public sim_dt
subroutine, public storm_drain_capture_initialize(is_cell_active, landuse_index)