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