diff --git a/scm/etc/scripts/UFS_case_gen.py b/scm/etc/scripts/UFS_case_gen.py index 54802d3a7..c6f1dfab1 100755 --- a/scm/etc/scripts/UFS_case_gen.py +++ b/scm/etc/scripts/UFS_case_gen.py @@ -56,6 +56,8 @@ missing_variable_snow_layers = 3 missing_variable_soil_layers = 4 missing_variable_ice_layers = 2 +missing_variable_vegetation_categories = 20 +missing_variable_soil_categories = 16 # Path to the directory containing processed case input files PROCESSED_CASE_DIR = '../../data/processed_case_input' @@ -1686,7 +1688,7 @@ def get_UFS_surface_data(dir, tile, i, j, old_chgres, lam): ######################################################################################## # ######################################################################################## -def get_UFS_oro_data(dir, tile, i, j, lam): +def get_UFS_oro_data(dir, tile, i, j, lam, old_chgres): """Get the orographic data for the given tile and indices""" if lam: @@ -1721,6 +1723,9 @@ def get_UFS_oro_data(dir, tile, i, j, lam): #lake variables (optional) lake_frac_in = read_NetCDF_var(nc_file, "lake_frac", i, j) lake_depth_in = read_NetCDF_var(nc_file, "lake_depth", i, j) + + vegtype_frac_in = read_NetCDF_surface_var(nc_file, "vegetation_type_pct", i, j, old_chgres, missing_variable_vegetation_categories) + soiltype_frac_in = read_NetCDF_surface_var(nc_file, "soil_type_pct", i, j, old_chgres, missing_variable_soil_categories) # nc_file.close() @@ -1744,7 +1749,9 @@ def get_UFS_oro_data(dir, tile, i, j, lam): "oro_uf": orog_raw_in, "landfrac": land_frac_in, "lakefrac": lake_frac_in, - "lakedepth": lake_depth_in} + "lakedepth": lake_depth_in, + "vegtype_frac": vegtype_frac_in, + "soiltype_frac": soiltype_frac_in} return oro @@ -3208,6 +3215,8 @@ def write_SCM_case_file(state, surface, oro, forcing, init, case, date, forcing_ snow_dim = nc_file.createDimension('nsnow', len(surface["snicexy"])) nslsnw_dim = nc_file.createDimension('nsoil_plus_nsnow',len(surface["snicexy"]) + len(surface["stc"])) ice_dim = nc_file.createDimension('nice', len(surface["tiice"])) + vegcat_dim = nc_file.createDimension('nvegcat', len(oro["vegtype_frac"])) + soilcat_dim =nc_file.createDimension('nsoilcat', len(oro["soiltype_frac"])) # timei_var = nc_file.createVariable('t0', wp, ('t0')) @@ -3388,7 +3397,9 @@ def write_SCM_case_file(state, surface, oro, forcing, init, case, date, forcing_ {"name": "oro_uf", "type":wp, "dimd": ('t0'), "units": "m", "desc": "unfiltered orography"}, \ {"name": "landfrac", "type":wp, "dimd": ('t0'), "units": "none", "desc": "fraction of horizontal grid area occupied by land"}, \ {"name": "lakefrac", "type":wp, "dimd": ('t0'), "units": "none", "desc": "fraction of horizontal grid area occupied by lake", "default_value":0}, \ - {"name": "lakedepth", "type":wp, "dimd": ('t0'), "units": "none", "desc": "lake depth", "default_value":0}] + {"name": "lakedepth", "type":wp, "dimd": ('t0'), "units": "none", "desc": "lake depth", "default_value":0}, + {"name": "vegtype_frac", "type":wp, "dimd": ('t0', 'nvegcat'), "units": "none", "desc": "fraction of horizontal grid area occupied by given vegetation category"}, + {"name": "soiltype_frac","type":wp, "dimd": ('t0', 'nsoilcat'), "units": "none", "desc": "fraction of horizontal grid area occupied by given soil category"}] # var_nsst = [{"name": "tref", "type":wp, "dimd": ('t0'), "units": "K", "desc": "sea surface reference temperature for NSST"}, \ {"name": "z_c", "type":wp, "dimd": ('t0'), "units": "m", "desc": "sub-layer cooling thickness for NSST"}, \ @@ -3682,7 +3693,7 @@ def main(): check_IC_hist_surface_compatibility(forcing_dir, hist_i, hist_j, surface_data, lam, old_chgres, tile) #read in orographic data for the initial conditions at the IC point - oro_data = get_UFS_oro_data(in_dir, tile, IC_i, IC_j, lam) + oro_data = get_UFS_oro_data(in_dir, tile, IC_i, IC_j, lam, old_chgres) #read in the initial condition profiles diff --git a/scm/src/scm_input.F90 b/scm/src/scm_input.F90 index 2233098bf..13ca8d9be 100644 --- a/scm/src/scm_input.F90 +++ b/scm/src/scm_input.F90 @@ -13,6 +13,8 @@ module scm_input integer :: missing_snow_layers = 3 integer :: missing_soil_layers = 4 integer :: missing_ice_layers = 2 +integer :: missing_nvegcat = 20 +integer :: missing_nsoilcat = 16 contains @@ -219,6 +221,7 @@ subroutine get_case_init(scm_state, scm_input) integer :: input_nice !< number of sea ice levels in the input file integer :: input_ntimes !< number of times represented in the input file integer :: input_nsoil_plus_nsnow !< number of combined snow and soil levels in the input file + integer :: input_nvegcat, input_nsoilcat real(kind=dp) :: input_lat !< column latitude (deg) real(kind=dp) :: input_lon !< column longitude (deg) real(kind=dp) :: input_area !< surface area [m^2] @@ -458,7 +461,19 @@ subroutine get_case_init(scm_state, scm_input) input_nice = missing_ice_layers else call check(NF90_INQUIRE_DIMENSION(ncid, varID, tmpName, input_nice),"nf90_inq_dim(nice)") - end if + end if + ierr = NF90_INQ_DIMID(ncid,"nvegcat",varID) + if(ierr /= NF90_NOERR) then + input_nvegcat = missing_nvegcat + else + call check(NF90_INQUIRE_DIMENSION(ncid, varID, tmpName, input_nvegcat),"nf90_inq_dim(nvegcat)") + end if + ierr = NF90_INQ_DIMID(ncid,"nsoilcat",varID) + if(ierr /= NF90_NOERR) then + input_nsoilcat = missing_nsoilcat + else + call check(NF90_INQUIRE_DIMENSION(ncid, varID, tmpName, input_nsoilcat),"nf90_inq_dim(nsoilcat)") + end if !> - Allocate the dimension variables. allocate(input_pres(input_nlev),input_time(input_ntimes), stat=allocate_status) @@ -714,7 +729,7 @@ subroutine get_case_init(scm_state, scm_input) call check(NF90_CLOSE(NCID=ncid),"nf90_close()") - call scm_input%create(input_ntimes, input_nlev, input_nsoil, input_nsnow, input_nice) + call scm_input%create(input_ntimes, input_nlev, input_nsoil, input_nsnow, input_nice, input_nvegcat, input_nsoilcat) ! GJF already done in scm_input%create routine !scm_input%input_nlev = input_nlev @@ -1050,6 +1065,8 @@ subroutine get_case_init_DEPHY(scm_state, scm_input) real(kind=dp), allocatable :: input_landfrac(:) !< fraction of horizontal grid area occupied by land real(kind=dp), allocatable :: input_lakefrac(:) !< fraction of horizontal grid area occupied by lake real(kind=dp), allocatable :: input_lakedepth(:) !< lake depth (m) + real(kind=dp), allocatable :: input_vegtype_frac(:,:) !< fraction of horizontal grid area occupied by given vegetation category + real(kind=dp), allocatable :: input_soiltype_frac(:,:) !< fraction of horizontal grid area occupied by given soil category real(kind=dp), allocatable :: input_tvxy(:) !< vegetation temperature (K) real(kind=dp), allocatable :: input_tgxy(:) !< ground temperature for Noahmp (K) @@ -1165,7 +1182,7 @@ subroutine get_case_init_DEPHY(scm_state, scm_input) integer :: jdat(1:8), idat(1:8) !(yr, mon, day, t-zone, hr, min, sec, mil-sec) logical :: needed_for_lsm_ics, needed_for_model_ics, lev_in_altitude - integer :: input_n_init_times, input_n_forcing_times, input_n_lev, input_n_snow, input_n_ice, input_n_soil + integer :: input_n_init_times, input_n_forcing_times, input_n_lev, input_n_snow, input_n_ice, input_n_soil, input_nvegcat, input_nsoilcat missing_value_eps = missing_value + 0.01 @@ -1213,6 +1230,18 @@ subroutine get_case_init_DEPHY(scm_state, scm_input) input_n_ice = missing_ice_layers else call check(NF90_INQUIRE_DIMENSION(ncid, varID, tmpName, input_n_ice),"nf90_inq_dim(nice)") + end if + ierr = NF90_INQ_DIMID(ncid,"nvegcat",varID) + if(ierr /= NF90_NOERR) then + input_nvegcat = missing_nvegcat + else + call check(NF90_INQUIRE_DIMENSION(ncid, varID, tmpName, input_nvegcat),"nf90_inq_dim(nvegcat)") + end if + ierr = NF90_INQ_DIMID(ncid,"nsoilcat",varID) + if(ierr /= NF90_NOERR) then + input_nsoilcat = missing_nsoilcat + else + call check(NF90_INQUIRE_DIMENSION(ncid, varID, tmpName, input_nsoilcat),"nf90_inq_dim(nsoilcat)") end if !> - Allocate the dimension variables. @@ -1399,6 +1428,8 @@ subroutine get_case_init_DEPHY(scm_state, scm_input) input_landfrac ( input_n_init_times), & input_lakefrac ( input_n_init_times), & input_lakedepth ( input_n_init_times), & + input_vegtype_frac (input_nvegcat, input_n_init_times), & + input_soiltype_frac (input_nsoilcat, input_n_init_times),& stat=allocate_status) allocate(input_tvxy ( input_n_init_times), & input_tgxy ( input_n_init_times), & @@ -1543,6 +1574,8 @@ subroutine get_case_init_DEPHY(scm_state, scm_input) call NetCDF_read_var(ncid, "landfrac", needed_for_model_ics, input_landfrac) call NetCDF_read_var(ncid, "lakefrac", needed_for_model_ics, input_lakefrac) call NetCDF_read_var(ncid, "lakedepth", needed_for_model_ics, input_lakedepth) + call NetCDF_read_var(ncid, "vegtype_frac", needed_for_model_ics, input_vegtype_frac) + call NetCDF_read_var(ncid, "soiltype_frac", needed_for_model_ics, input_soiltype_frac) !NSST variables call NetCDF_read_var(ncid, "tref", needed_for_model_ics, input_tref) @@ -1805,7 +1838,7 @@ subroutine get_case_init_DEPHY(scm_state, scm_input) call check(NF90_CLOSE(NCID=ncid),"nf90_close()") - call scm_input%create(input_n_forcing_times, input_n_lev, input_n_soil, input_n_snow, input_n_ice) + call scm_input%create(input_n_forcing_times, input_n_lev, input_n_soil, input_n_snow, input_n_ice, input_nvegcat, input_nsoilcat) !fill the scm_input DDT @@ -2046,6 +2079,8 @@ subroutine get_case_init_DEPHY(scm_state, scm_input) scm_input%input_landfrac = input_landfrac(active_init_time) scm_input%input_lakefrac = input_lakefrac(active_init_time) scm_input%input_lakedepth= input_lakedepth(active_init_time) + scm_input%input_vegtype_frac = input_vegtype_frac(:,active_init_time) + scm_input%input_soiltype_frac = input_soiltype_frac(:,active_init_time) scm_input%input_tref = input_tref(active_init_time) scm_input%input_z_c = input_z_c(active_init_time) diff --git a/scm/src/scm_type_defs.F90 b/scm/src/scm_type_defs.F90 index 01629bacf..5302b7507 100644 --- a/scm/src/scm_type_defs.F90 +++ b/scm/src/scm_type_defs.F90 @@ -180,6 +180,8 @@ module scm_type_defs integer :: input_nsoil !< number of soil levels in the input file integer :: input_nsnow !< number of snow layers in the input file integer :: input_nice !< number of sea ice layers in the input file + integer :: input_nvegcat !< number of vegetation type categories + integer :: input_nsoilcat !< number of soil type categories integer :: input_ntimes !< number of times in the input file where forcing is available real(kind=dp) :: input_lat !< latitude of column center real(kind=dp) :: input_lon !< longitude of column center @@ -262,6 +264,8 @@ module scm_type_defs real(kind=dp) :: input_landfrac !< fraction of horizontal grid area occupied by land real(kind=dp) :: input_lakefrac !< fraction of horizontal grid area occupied by lake real(kind=dp) :: input_lakedepth !< lake depth (m) + real(kind=dp), allocatable :: input_vegtype_frac(:) !< fraction of horizontal grid area occupied by given vegetation category + real(kind=dp), allocatable :: input_soiltype_frac(:) !< fraction of horizontal grid area occupied by given soil category real(kind=dp) :: input_tvxy !< vegetation temperature (K) real(kind=dp) :: input_tgxy !< ground temperature for Noahmp (K) @@ -668,14 +672,16 @@ subroutine scm_state_create(scm_state, n_columns, n_levels, n_soil, n_snow, n_ti end subroutine scm_state_create - subroutine scm_input_create(scm_input, ntimes, nlev, nsoil, nsnow, nice) + subroutine scm_input_create(scm_input, ntimes, nlev, nsoil, nsnow, nice, nvegcat, nsoilcat) class(scm_input_type) :: scm_input - integer, intent(in) :: ntimes, nlev, nsoil, nsnow, nice + integer, intent(in) :: ntimes, nlev, nsoil, nsnow, nice, nvegcat, nsoilcat scm_input%input_nlev = nlev scm_input%input_nsoil = nsoil scm_input%input_nsnow = nsnow scm_input%input_nice = nice + scm_input%input_nvegcat = nvegcat + scm_input%input_nsoilcat = nsoilcat scm_input%input_ntimes = ntimes scm_input%input_lat = real_zero @@ -762,6 +768,10 @@ subroutine scm_input_create(scm_input, ntimes, nlev, nsoil, nsnow, nice) scm_input%input_landfrac = real_zero scm_input%input_lakefrac = real_zero scm_input%input_lakedepth = real_zero + allocate(scm_input%input_vegtype_frac(nvegcat)) + scm_input%input_vegtype_frac = real_zero + allocate(scm_input%input_soiltype_frac(nsoilcat)) + scm_input%input_soiltype_frac = real_zero scm_input%input_tvxy = real_zero scm_input%input_tgxy = real_zero @@ -1085,8 +1095,10 @@ subroutine physics_set(physics, scm_input, scm_state) call conditionally_set_var(scm_input%input_landfrac, physics%Sfcprop%landfrac(i), "landfrac", physics%Model%frac_grid, missing_var(17)) call conditionally_set_var(scm_input%input_lakefrac, physics%Sfcprop%lakefrac(i), "lakefrac", (physics%Model%lkm == 1), missing_var(18)) call conditionally_set_var(scm_input%input_lakedepth, physics%Sfcprop%lakedepth(i), "lakedepth", (physics%Model%lkm == 1), missing_var(19)) + call conditionally_set_var(scm_input%input_vegtype_frac(:), physics%Sfcprop%vegtype_frac(i,:), "vegtype_frac", .true., missing_var(20)) + call conditionally_set_var(scm_input%input_soiltype_frac(:), physics%Sfcprop%soiltype_frac(i,:), "soiltype_frac", .true., missing_var(21)) - n = 19 + n = 21 if ( i==1 .and. ANY( missing_var(1:n) ) ) then write(0,'(a)') "INPUT CHECK: Some missing input data was found related to (potentially non-required) orography and gravity wave drag parameters. This may lead to crashes or other strange behavior." write(0,'(a)') "Check scm_type_defs.F90/physics_set to see the names of variables that are missing, corresponding to the following indices:"