module initInterpMindist ! ------------------------------------------------------------------------ ! This module contains most of the "interesting" logic of initInterp, in terms of ! finding the input column (or landunit, patch, etc.) to use as a template for each ! output column (etc.). ! ! This is in a separate module to facilitate unit testing, since the full initInterp ! involves some awkward dependencies. ! ------------------------------------------------------------------------ #include "shr_assert.h" use shr_kind_mod , only: r8 => shr_kind_r8 use shr_log_mod , only : errMsg => shr_log_errMsg use clm_varctl , only: iulog use abortutils , only: endrun use spmdMod , only: masterproc use clm_varcon , only: spval, re use glcBehaviorMod , only: glc_behavior_type implicit none private save ! Public methods public :: set_mindist public :: set_single_match ! Public types type, public :: subgrid_special_indices_type integer :: ipft_not_vegetated integer :: icol_vegetated_or_bare_soil integer :: ilun_vegetated_or_bare_soil integer :: ilun_crop integer :: ilun_landice_multiple_elevation_classes contains procedure :: is_vegetated_landunit ! returns true if the given landunit type is natural veg or crop end type subgrid_special_indices_type type, public :: subgrid_type character(len=16) :: name ! pft, column, landunit, gridcell integer , pointer :: ptype(:) => null() ! used for patch type integer , pointer :: ctype(:) => null() ! used for patch or col type integer , pointer :: ltype(:) => null() ! used for pft, col or lun type real(r8), pointer :: topoglc(:) => null() real(r8), pointer :: lat(:) real(r8), pointer :: lon(:) real(r8), pointer :: coslat(:) contains procedure :: print_point ! print info about one point end type subgrid_type ! Private methods private :: set_glcmec_must_be_same_type private :: set_icemec_adjustable_type private :: do_fill_missing_with_natveg private :: is_sametype private :: is_baresoil character(len=*), parameter, private :: sourcefile = & __FILE__ contains !----------------------------------------------------------------------- subroutine print_point(this, index, unit) ! ! !DESCRIPTION: ! Print info about one point in a subgrid_type object ! ! !USES: ! ! !ARGUMENTS: class(subgrid_type), intent(in) :: this integer , intent(in) :: index integer , intent(in) :: unit ! unit to which we should write the info ! ! !LOCAL VARIABLES: character(len=*), parameter :: subname = 'print_point' !----------------------------------------------------------------------- write(unit,*) 'subgrid level, index = ',& this%name, index write(unit,*) 'lat, lon = ', this%lat(index), ', ', this%lon(index) if (associated(this%ltype)) then write(unit,*) 'ltype: ', this%ltype(index) end if if (associated(this%ctype)) then write(unit,*) 'ctype: ', this%ctype(index) end if if (associated(this%ptype)) then write(unit,*) 'ptype: ', this%ptype(index) end if end subroutine print_point !======================================================================= subroutine set_mindist(begi, endi, bego, endo, activei, activeo, subgridi, subgrido, & subgrid_special_indices, glc_behavior, glc_elevclasses_same, & fill_missing_with_natveg, mindist_index) ! -------------------------------------------------------------------- ! arguments integer , intent(in) :: begi, endi integer , intent(in) :: bego, endo logical , intent(in) :: activei(begi:endi) logical , intent(in) :: activeo(bego:endo) type(subgrid_type) , intent(in) :: subgridi type(subgrid_type) , intent(in) :: subgrido type(subgrid_special_indices_type), intent(in) :: subgrid_special_indices type(glc_behavior_type), intent(in) :: glc_behavior ! True if the number and bounds of glacier elevation classes are the same in the ! input and output; false otherwise. logical , intent(in) :: glc_elevclasses_same ! If false: if an output type cannot be found in the input, code aborts ! If true: if an output type cannot be found in the input, fill with closest natural ! veg column (using bare soil for patch-level variables) ! ! NOTE: always treated as true for natural veg and crop landunits/columns/patches in ! the output - e.g., if we can't find the right column type to fill crop, we always ! use the closest natural veg column, regardless of the value of this flag. logical , intent(in) :: fill_missing_with_natveg integer , intent(out) :: mindist_index(bego:endo) ! ! local variables real(r8) :: dx,dy real(r8) :: distmin,dist,hgtdiffmin,hgtdiff integer :: nsizei, nsizeo integer :: ni,no,nmin,ier,n,noloc logical :: topoglc_present logical :: closest ! Whether two glcmec columns/patches must be the same column/patch type to be ! considered the same type. This is only valid for glcmec points, and is only valid ! for subgrid name = 'pft' or 'column'. logical :: glcmec_must_be_same_type_o(bego:endo) ! -------------------------------------------------------------------- if (associated(subgridi%topoglc) .and. associated(subgrido%topoglc)) then topoglc_present = .true. else topoglc_present = .false. end if mindist_index(bego:endo) = 0 distmin = spval call set_glcmec_must_be_same_type(bego=bego, endo=endo, dimname=subgrido%name, & glc_elevclasses_same = glc_elevclasses_same, glc_behavior=glc_behavior, & glcmec_must_be_same_type_o=glcmec_must_be_same_type_o) !$OMP PARALLEL DO PRIVATE (ni,no,n,nmin,distmin,dx,dy,dist,closest,hgtdiffmin,hgtdiff) do no = bego,endo ! Only interpolate onto active points. Otherwise, the mere act of running ! init_interp (e.g., of a file onto itself) would lead to changes in a bunch of ! inactive points - e.g., going from their cold start initial conditions to some ! spunup initial conditions (from the closest active point of that type). This ! could potentially lead to different behavior in a transient run, if those points ! later became active; that's undesirable. if (activeo(no)) then nmin = 0 distmin = spval hgtdiffmin = spval do ni = begi,endi if (activei(ni)) then if (is_sametype(ni = ni, no = no, & subgridi = subgridi, subgrido = subgrido, & subgrid_special_indices = subgrid_special_indices, & glcmec_must_be_same_type = glcmec_must_be_same_type_o(no), & veg_patch_just_considers_ptype = .true.)) then dy = abs(subgrido%lat(no)-subgridi%lat(ni))*re dx = abs(subgrido%lon(no)-subgridi%lon(ni))*re * & 0.5_r8*(subgrido%coslat(no)+subgridi%coslat(ni)) dist = dx*dx + dy*dy if (topoglc_present) then hgtdiff = abs(subgridi%topoglc(ni) - subgrido%topoglc(no)) end if closest = .false. if ( dist < distmin ) then closest = .true. distmin = dist nmin = ni if (topoglc_present) then hgtdiffmin = hgtdiff end if end if if (.not. closest) then ! In *some* cases: For glc_mec points, we first find the closest ! point in lat-lon space, without consideration for column or patch ! type (above). Then, within that closest point, we find the closest ! column in topographic space; this second piece is done here. Note ! that this ordering means that we could choose a column with a very ! different topographic height from the target column, if it is ! closer in lat-lon space than any similar-height columns. if (topoglc_present) then hgtdiff = abs(subgridi%topoglc(ni) - subgrido%topoglc(no)) if ((dist == distmin) .and. (hgtdiff < hgtdiffmin)) then closest = .true. hgtdiffmin = hgtdiff distmin = dist nmin = ni end if end if end if end if end if end do ! If output type is not contained in input dataset, then use closest bare soil, ! if this point is one for which we fill missing with natveg. if ( distmin == spval .and. & do_fill_missing_with_natveg( & fill_missing_with_natveg, no, subgrido, subgrid_special_indices)) then do ni = begi, endi if (activei(ni)) then if ( is_baresoil(ni, subgridi, subgrid_special_indices)) then dy = abs(subgrido%lat(no)-subgridi%lat(ni))*re dx = abs(subgrido%lon(no)-subgridi%lon(ni))*re * & 0.5_r8*(subgrido%coslat(no)+subgridi%coslat(ni)) dist = dx*dx + dy*dy if ( dist < distmin )then distmin = dist nmin = ni end if end if end if end do end if ! Error conditions if ( distmin == spval )then write(iulog,*) 'ERROR initInterp set_mindist: & &Cannot find any input points matching output point:' call subgrido%print_point(no, iulog) write(iulog,*) ' ' write(iulog,*) 'Consider rerunning with the following in user_nl_clm:' write(iulog,*) 'init_interp_fill_missing_with_natveg = .true.' write(iulog,*) 'However, note that this will fill all missing types in the output' write(iulog,*) 'with the closest natural veg column in the input' write(iulog,*) '(using bare soil for patch-level variables).' write(iulog,*) 'So, you should consider whether that is what you want.' call endrun(msg=errMsg(sourcefile, __LINE__)) end if mindist_index(no) = nmin end if ! end if activeo block end do !$OMP END PARALLEL DO end subroutine set_mindist !----------------------------------------------------------------------- subroutine set_single_match(begi, endi, bego, endo, activeo, subgridi, subgrido, & subgrid_special_indices, glc_behavior, glc_elevclasses_same, & mindist_index) ! ! !DESCRIPTION: ! This is an alternative to set_mindist that can work when there is a single matching ! input point for each output point (or possibly 0 matching inputs for an inactive ! output point). In this method, we only look for points that share the same lat/lon ! as the output point. ! ! This method is intended to be used for the case where we want to copy areas from the ! input file to the output file (rather than maintaining areas from the surface ! dataset); some of the logic here is set up the way it is in order to support this ! copying of areas. ! ! !ARGUMENTS: integer , intent(in) :: begi, endi integer , intent(in) :: bego, endo logical , intent(in) :: activeo(bego:endo) type(subgrid_type) , intent(in) :: subgridi type(subgrid_type) , intent(in) :: subgrido type(subgrid_special_indices_type), intent(in) :: subgrid_special_indices type(glc_behavior_type), intent(in) :: glc_behavior ! True if the number and bounds of glacier elevation classes are the same in the ! input and output; false otherwise. ! ! Must be true for this method! (Otherwise we may not be able to find a single unique ! input point for each output point.) logical , intent(in) :: glc_elevclasses_same integer , intent(out) :: mindist_index(bego:endo) ! ! !LOCAL VARIABLES: integer :: ni, no, ni_match logical :: found real(r8) :: dx, dy logical :: ni_sametype ! Whether two glcmec columns/patches must be the same column/patch type to be ! considered the same type. This is only valid for glcmec points, and is only valid ! for subgrid name = 'pft' or 'column'. logical :: glcmec_must_be_same_type_o(bego:endo) ! Tolerance in lat/lon for considering a point to be at the same location real(r8) :: same_point_tol = 1.e-14_r8 character(len=*), parameter :: subname = 'set_single_match' !----------------------------------------------------------------------- if (.not. glc_elevclasses_same) then write(iulog,*) subname// & ' ERROR: For this init_interp method, the number and bounds of' write(iulog,*) 'glacier elevation classes must be the same in input and output' call endrun(msg=subname//' ERROR: glc_elevclasses_same must be true for this method') end if call set_glcmec_must_be_same_type(bego=bego, endo=endo, dimname=subgrido%name, & glc_elevclasses_same = glc_elevclasses_same, glc_behavior=glc_behavior, & glcmec_must_be_same_type_o=glcmec_must_be_same_type_o) !$OMP PARALLEL DO PRIVATE (ni,no,ni_match,found,dx,dy,ni_sametype) do no = bego, endo ni_match = 0 found = .false. do ni = begi, endi dx = abs(subgrido%lon(no)-subgridi%lon(ni)) dy = abs(subgrido%lat(no)-subgridi%lat(ni)) if (dx < same_point_tol .and. dy < same_point_tol) then ni_sametype = is_sametype(ni = ni, no = no, & subgridi = subgridi, subgrido = subgrido, & subgrid_special_indices = subgrid_special_indices, & glcmec_must_be_same_type = glcmec_must_be_same_type_o(no), & veg_patch_just_considers_ptype = .false.) if (ni_sametype) then if (found) then write(iulog,*) subname// & ' ERROR: found multiple input points matching output point' call subgrido%print_point(no, iulog) write(iulog,*) ' ' write(iulog,*) 'For this init_interp_method: for a given output point,' write(iulog,*) 'we expect to find exactly one input point at the same' write(iulog,*) 'location with the same type.' write(iulog,*) '(This most likely indicates a problem in CTSM, such as' write(iulog,*) 'having multiple columns with the same column type on a' write(iulog,*) 'given gridcell.)' call endrun(msg=subname// & ' ERROR: found multiple input points matching output point') else found = .true. ni_match = ni end if end if end if end do if (found) then mindist_index(no) = ni_match else mindist_index(no) = 0 if (activeo(no)) then write(iulog,*) subname// & ' ERROR: cannot find any input points matching output point' call subgrido%print_point(no, iulog) write(iulog,*) ' ' write(iulog,*) 'For this init_interp_method: for a given active output point,' write(iulog,*) 'we expect to find exactly one input point at the same' write(iulog,*) 'location with the same type.' write(iulog,*) 'Note that this requires the input and output grids to be identical.' write(iulog,*) '(If you need to interpolate to a different resolution, then' write(iulog,*) 'you need to use a different init_interp_method.)' call endrun(msg=subname//' ERROR: cannot find any input points matching output point') end if end if end do !$OMP END PARALLEL DO end subroutine set_single_match !----------------------------------------------------------------------- subroutine set_glcmec_must_be_same_type(bego, endo, dimname, & glc_elevclasses_same, glc_behavior, & glcmec_must_be_same_type_o) ! ! !DESCRIPTION: ! Sets the glcmec_must_be_same_type_o array for each output icemec point ! ! This array will be set to true for icemec output columns/patches for which the ! column/patch type must match the input column/patch type to be considered the same ! type. This is only valid for icemec points, and is only valid for dimname = 'pft' or ! 'column' - for others, the value is undefined. ! ! This assumes that bego and endo match the bounds that are used elsewhere in the ! model - i.e., for the decomposition within init_interp to match the model's ! decomposition. ! ! !ARGUMENTS: integer , intent(in) :: bego ! beginning index for output points integer , intent(in) :: endo ! ending index for output points character(len=*) , intent(in) :: dimname ! 'pft', 'column', etc. ! True if the number and bounds of glacier elevation classes are the same in the ! input and output; false otherwise. logical, intent(in) :: glc_elevclasses_same type(glc_behavior_type), intent(in) :: glc_behavior logical, intent(out) :: glcmec_must_be_same_type_o( bego: ) ! see description above ! ! !LOCAL VARIABLES: integer :: no ! Whether each output icemec point has adjustable column type; only valid for icemec ! points, and only valid for subgrid name = pft or column logical :: icemec_adjustable_type_o(bego:endo) character(len=*), parameter :: subname = 'set_glcmec_must_be_same_type' !----------------------------------------------------------------------- SHR_ASSERT_ALL((ubound(glcmec_must_be_same_type_o) == (/endo/)), errMsg(sourcefile, __LINE__)) if (.not. glc_elevclasses_same) then ! If the number or bounds of the elevation classes differ between input and ! output, then we ignore the column and patch types of glcmec points when ! looking for the same type - instead using closest topographic height as a ! tie-breaker between equidistant columns/patches. glcmec_must_be_same_type_o(bego:endo) = .false. else call set_icemec_adjustable_type(bego=bego, endo=endo, dimname=dimname, & glc_behavior=glc_behavior, icemec_adjustable_type_o=icemec_adjustable_type_o) do no = bego, endo if (icemec_adjustable_type_o(no)) then ! If glcmec points in this output cell have adjustable type, then we ignore ! the column and patch types of glcmec points when looking for the same ! type: we want to find the closest glcmec point in lat-lon space without ! regards for column/patch type, because the column/patch types may change ! at runtime. glcmec_must_be_same_type_o(no) = .false. else ! Otherwise, we require the column and patch types to be the same between ! input and output for this glcmec output point, as is the case for most ! other landunits. This is important for a case with interpolation to give ! bit-for-bit answers with a case without interpolation (since glcmec ! topographic heights can change after initialization, so we can't always ! rely on the point with closest topographic height to be the "right" one to ! pick as the source for interpolation). glcmec_must_be_same_type_o(no) = .true. end if end do end if end subroutine set_glcmec_must_be_same_type !----------------------------------------------------------------------- subroutine set_icemec_adjustable_type(bego, endo, dimname, glc_behavior, & icemec_adjustable_type_o) ! ! !DESCRIPTION: ! Sets the icemec_adjustable_type_o array for each output icemec point ! ! This array will be set to true for icemec points that have adjustable column type ! and false for icemec points that do not. The value will be undefined for non-icemec ! points. ! ! This can only be called for the output, not the input! ! ! The output array is only valid for dimname = pft or column; for others, the output ! values are undefined. ! ! This assumes that bego and endo match the bounds that are used elsewhere in the ! model - i.e., for the decomposition within init_interp to match the model's ! decomposition. ! ! !ARGUMENTS: integer , intent(in) :: bego ! beginning index for output points integer , intent(in) :: endo ! ending index for output points character(len=*) , intent(in) :: dimname ! 'pft', 'column', etc. type(glc_behavior_type) , intent(in) :: glc_behavior logical , intent(out) :: icemec_adjustable_type_o( bego: ) ! see documentation above ! ! !LOCAL VARIABLES: character(len=*), parameter :: subname = 'set_icemec_adjustable_type' !----------------------------------------------------------------------- SHR_ASSERT_ALL((ubound(icemec_adjustable_type_o) == (/endo/)), errMsg(sourcefile, __LINE__)) select case (dimname) case ('pft') call glc_behavior%patches_have_dynamic_type_array(bego, endo, & icemec_adjustable_type_o(bego:endo)) case ('column') call glc_behavior%cols_have_dynamic_type_array(bego, endo, & icemec_adjustable_type_o(bego:endo)) case ('landunit', 'gridcell') ! Do nothing: icemec_adjustable_type_o will be left undefined case default call endrun(subname//' ERROR: unexpected dimname: '//trim(dimname)//& errMsg(sourcefile, __LINE__)) end select end subroutine set_icemec_adjustable_type !----------------------------------------------------------------------- function do_fill_missing_with_natveg(fill_missing_with_natveg, & no, subgrido, subgrid_special_indices) ! ! !DESCRIPTION: ! Returns true if the given output point, if missing, should be filled with the ! closest natural veg point. ! ! !ARGUMENTS: logical :: do_fill_missing_with_natveg ! function result ! whether we should fill ALL missing points with natveg logical, intent(in) :: fill_missing_with_natveg integer , intent(in) :: no type(subgrid_type), intent(in) :: subgrido type(subgrid_special_indices_type), intent(in) :: subgrid_special_indices ! ! !LOCAL VARIABLES: character(len=*), parameter :: subname = 'do_fill_missing_with_natveg' !----------------------------------------------------------------------- if (subgrido%name == 'gridcell') then ! It makes no sense to try to fill missing with natveg for gridcell-level values do_fill_missing_with_natveg = .false. else if (fill_missing_with_natveg) then ! User has asked for all missing points to be filled with natveg do_fill_missing_with_natveg = .true. else if (subgrid_special_indices%is_vegetated_landunit(subgrido%ltype(no))) then ! Even if user hasn't asked for it, we fill missing vegetated points (natural veg ! and crop) with the closest natveg point. This is mainly to support the common ! use case of interpolating non-crop to crop, but also supports adding a new PFT ! type. do_fill_missing_with_natveg = .true. else do_fill_missing_with_natveg = .false. end if end function do_fill_missing_with_natveg !======================================================================= logical function is_sametype (ni, no, subgridi, subgrido, subgrid_special_indices, & glcmec_must_be_same_type, veg_patch_just_considers_ptype) ! -------------------------------------------------------------------- ! arguments integer , intent(in) :: ni integer , intent(in) :: no type(subgrid_type), intent(in) :: subgridi type(subgrid_type), intent(in) :: subgrido type(subgrid_special_indices_type), intent(in) :: subgrid_special_indices ! Whether two glcmec columns/patches must be the same column/patch type to be ! considered the same type. This is only valid for glcmec points, and is only valid ! for subgrid name = 'pft' or 'column'. logical, intent(in) :: glcmec_must_be_same_type ! For vegetated patches (natural veg or crop): if veg_patch_just_considers_ptype is ! true, then we consider two vegetated patches to be the same type if they have the ! same ptype, without regard to their column and landunit types (as long as both the ! input and the output are on either the natural veg or crop landunit). This is ! needed to handle the generic crop properly when interpolating from non-crop to ! crop, or vice versa. ! ! If false, then they need to have the same column and landunit types, too (as is the ! general case). logical, intent(in) :: veg_patch_just_considers_ptype ! -------------------------------------------------------------------- is_sametype = .false. if (trim(subgridi%name) == 'pft' .and. trim(subgrido%name) == 'pft') then if ( .not. glcmec_must_be_same_type .and. & subgridi%ltype(ni) == subgrid_special_indices%ilun_landice_multiple_elevation_classes .and. & subgrido%ltype(no) == subgrid_special_indices%ilun_landice_multiple_elevation_classes) then is_sametype = .true. else if (veg_patch_just_considers_ptype .and. & subgrid_special_indices%is_vegetated_landunit(subgrido%ltype(no))) then ! See comment attached to the declaration of veg_patch_just_considers_ptype for ! rationale for this logic. ! ! TODO(wjs, 2015-09-15) If we ever allow the same PFT to appear on multiple ! columns within a given grid cell, then this logic will need to be made ! somewhat more complex: e.g., preferably take something from the same column ! type, but if we can't find anything from the same column type, then ignore ! column type. if (subgrid_special_indices%is_vegetated_landunit(subgridi%ltype(ni)) .and. & subgridi%ptype(ni) == subgrido%ptype(no)) then is_sametype = .true. end if else if (subgridi%ptype(ni) == subgrido%ptype(no) .and. & subgridi%ctype(ni) == subgrido%ctype(no) .and. & subgridi%ltype(ni) == subgrido%ltype(no)) then is_sametype = .true. end if else if (trim(subgridi%name) == 'column' .and. trim(subgrido%name) == 'column') then if ( .not. glcmec_must_be_same_type .and. & subgridi%ltype(ni) == subgrid_special_indices%ilun_landice_multiple_elevation_classes .and. & subgrido%ltype(no) == subgrid_special_indices%ilun_landice_multiple_elevation_classes ) then is_sametype = .true. else if (subgridi%ctype(ni) == subgrido%ctype(no) .and. & subgridi%ltype(ni) == subgrido%ltype(no)) then is_sametype = .true. end if else if (trim(subgridi%name) == 'landunit' .and. trim(subgrido%name) == 'landunit') then if (subgridi%ltype(ni) == subgrido%ltype(no)) then is_sametype = .true. end if else if (trim(subgridi%name) == 'gridcell' .and. trim(subgrido%name) == 'gridcell') then is_sametype = .true. else if (masterproc) then write(iulog,*)'ERROR interpinic: is_sametype check on input and output type not supported' write(iulog,*)'typei = ',trim(subgridi%name) write(iulog,*)'typeo = ',trim(subgrido%name) end if call endrun(msg=errMsg(sourcefile, __LINE__)) end if end function is_sametype !======================================================================= logical function is_baresoil (n, subgrid, subgrid_special_indices) ! -------------------------------------------------------------------- ! arguments integer , intent(in) :: n type(subgrid_type), intent(in) :: subgrid type(subgrid_special_indices_type), intent(in) :: subgrid_special_indices ! -------------------------------------------------------------------- is_baresoil = .false. if (subgrid%name == 'pft') then if (subgrid%ptype(n) == subgrid_special_indices%ipft_not_vegetated .and. & subgrid%ctype(n) == subgrid_special_indices%icol_vegetated_or_bare_soil .and. & subgrid%ltype(n) == subgrid_special_indices%ilun_vegetated_or_bare_soil) then is_baresoil = .true. end if else if (subgrid%name == 'column') then if (subgrid%ctype(n) == subgrid_special_indices%icol_vegetated_or_bare_soil .and. & subgrid%ltype(n) == subgrid_special_indices%ilun_vegetated_or_bare_soil) then is_baresoil = .true. end if else if (subgrid%name == 'landunit') then if (subgrid%ltype(n) == subgrid_special_indices%ilun_vegetated_or_bare_soil) then is_baresoil = .true. end if else if (masterproc) then write(iulog,*)'ERROR interpinic: is_baresoil subgrid type ',subgrid%name,' not supported' end if call endrun(msg=errMsg(sourcefile, __LINE__)) end if end function is_baresoil !----------------------------------------------------------------------- function is_vegetated_landunit(this, ltype) ! ! !DESCRIPTION: ! Returns true if the given landunit type is vegetated: either natural veg or crop ! ! !USES: ! ! !ARGUMENTS: logical :: is_vegetated_landunit ! function result class(subgrid_special_indices_type), intent(in) :: this integer, intent(in) :: ltype ! landunit type of interest ! ! !LOCAL VARIABLES: character(len=*), parameter :: subname = 'is_vegetated_landunit' !----------------------------------------------------------------------- if (ltype == this%ilun_vegetated_or_bare_soil .or. & ltype == this%ilun_crop) then is_vegetated_landunit = .true. else is_vegetated_landunit = .false. end if end function is_vegetated_landunit end module initInterpMindist