diff --git a/bld/namelist_files/namelist_defaults_ctsm.xml b/bld/namelist_files/namelist_defaults_ctsm.xml
index 52d876d6e3..c9da63db74 100644
--- a/bld/namelist_files/namelist_defaults_ctsm.xml
+++ b/bld/namelist_files/namelist_defaults_ctsm.xml
@@ -677,8 +677,6 @@ lnd/clm2/surfdata_map/surfdata_1x1_mexicocityMEX_16pfts_Irrig_CMIP6_simyr2000_c1
lnd/clm2/surfdata_map/surfdata_1x1_urbanc_alpha_16pfts_Irrig_CMIP6_simyr2000_c171214.nc
-
-lnd/clm2/surfdata_map/surfdata_0.47x0.63_16pfts_Irrig_CMIP6_simyr1850_c180508.nc
lnd/clm2/surfdata_map/surfdata_360x720cru_16pfts_Irrig_CMIP6_simyr1850_c170824.nc
diff --git a/doc/ChangeLog b/doc/ChangeLog
index 54e91aba49..b2a1e1d818 100644
--- a/doc/ChangeLog
+++ b/doc/ChangeLog
@@ -1,4 +1,137 @@
===============================================================
+Tag name: ctsm1.0.dev052
+Originator(s): sacks (Bill Sacks)
+Date: Mon Jul 22 14:02:43 MDT 2019
+One-line Summary: Fix rare soil color bug in mksurfdata_map
+
+Purpose of changes
+------------------
+
+Under rare conditions, mksurfdata_map could put the default soil color
+in an output cell where there is actually more information. This tag
+fixes that issue. None of the out-of-the-box surface datasets are
+impacted by this bug, so I have not recreated any surface datasets. (I
+checked all out-of-the-box surface datasets except for
+surfdata_0.125x0.125_mp24_simyr2000_c150114.nc, because it doesn't get
+remade cleanly out of the box.)
+
+Also:
+
+- Add some unit tests for the creation of soil color in mksurfdata_map
+
+- Point to correct (existing) surface dataset for year-1850 at f05
+ resolution
+
+Bugs fixed or introduced
+------------------------
+
+Issues fixed (include CTSM Issue #):
+- Resolves ESCOMP/ctsm#4 (Minor bug in creation of soil color in
+ mksurfdata_map: points can be given the default soil color, when they
+ should have a real color)
+- Resolves ESCOMP/ctsm#765 (Year-1850 f05 surface dataset missing from
+ inputdata repository)
+
+Significant changes to scientifically-supported configurations
+--------------------------------------------------------------
+
+Does this tag change answers significantly for any of the following physics configurations?
+(Details of any changes will be given in the "Answer changes" section below.)
+
+ [Put an [X] in the box for any configuration with significant answer changes.]
+
+[ ] clm5_0
+
+[ ] ctsm5_0-nwp
+
+[ ] clm4_5
+
+Notes of particular relevance for users
+---------------------------------------
+
+Caveats for users (e.g., need to interpolate initial conditions): none
+
+Changes to CTSM's user interface (e.g., new/renamed XML or namelist variables): none
+
+Changes made to namelist defaults (e.g., changed parameter values): none
+
+Changes to the datasets (e.g., parameter, surface or initial files):
+- Year-1850 f05 surface dataset now points to a file that actually exists
+
+Substantial timing or memory changes: none
+
+Notes of particular relevance for developers: (including Code reviews and testing)
+---------------------------------------------
+NOTE: Be sure to review the steps in README.CHECKLIST.master_tags as well as the coding style in the Developers Guide
+
+Caveats for developers (e.g., code that is duplicated that requires double maintenance):
+
+Changes to tests or testing: none
+
+Code reviewed by: self
+
+
+CTSM testing:
+
+ [PASS means all tests PASS and OK means tests PASS other than expected fails.]
+
+ build-namelist tests:
+
+ cheyenne - not run
+
+ tools-tests (test/tools):
+
+ cheyenne - ok
+
+ Most tests pass, including baseline comparisons. The following tests fail, but also failed for me on master:
+
+ 019 smiS4 TSMscript_tools.sh ncl_scripts getregional_datasets.pl getregional ....................\c
+ rc=6 FAIL
+ 020 bliS4 TBLscript_tools.sh ncl_scripts getregional_datasets.pl getregional ....................\c
+ rc=4 FAIL
+ 027 smf84 TSMscript_tools.sh PTCLM PTCLMmkdata PTCLM_USUMB_clm4_5^buildtools ....................\c
+ rc=6 FAIL
+ 028 blf84 TBLscript_tools.sh PTCLM PTCLMmkdata PTCLM_USUMB_clm4_5^buildtools ....................\c
+ rc=4 FAIL
+ 029 smfc4 TSMscript_tools.sh PTCLM PTCLMmkdata PTCLM_USUMB_Cycle_clm4_5^buildtools ..............\c
+ rc=6 FAIL
+ 030 blfc4 TBLscript_tools.sh PTCLM PTCLMmkdata PTCLM_USUMB_Cycle_clm4_5^buildtools ..............\c
+ rc=4 FAIL
+ 031 smfg4 TSMscript_tools.sh PTCLM PTCLMmkdata PTCLM_USUMB_Global_clm4_5^buildtools .............\c
+ rc=6 FAIL
+ 032 blfg4 TBLscript_tools.sh PTCLM PTCLMmkdata PTCLM_USUMB_Global_clm4_5^buildtools .............\c
+ rc=4 FAIL
+
+ PTCLM testing (tools/shared/PTCLM/test):
+
+ cheyenne - not run
+
+ python testing (see instructions in python/README.md; document testing done):
+
+ (any machine) - not run
+
+ regular tests (aux_clm):
+
+ cheyenne ---- not run
+ hobart ------ not run
+
+If the tag used for baseline comparisons was NOT the previous tag, note that here:
+
+
+Answer changes
+--------------
+
+Changes answers relative to baseline: NO
+
+Detailed list of changes
+------------------------
+
+List any externals directories updated (cime, rtm, mosart, cism, fates, etc.): none
+
+Pull Requests that document the changes (include PR ids): none
+
+===============================================================
+===============================================================
Tag name: ctsm1.0.dev051
Originator(s): sacks (Bill Sacks)
Date: Fri Jul 19 13:26:25 MDT 2019
diff --git a/doc/ChangeSum b/doc/ChangeSum
index 41684e0579..701d7dda0d 100644
--- a/doc/ChangeSum
+++ b/doc/ChangeSum
@@ -1,5 +1,6 @@
Tag Who Date Summary
============================================================================================================================
+ ctsm1.0.dev052 sacks 07/22/2019 Fix rare soil color bug in mksurfdata_map
ctsm1.0.dev051 sacks 07/19/2019 Update water tracers for remainder of first stage of hydrology
ctsm1.0.dev050 slevis 07/15/2019 dz --> dz_lake bug-fix in LakeTemperatureMod.F90 line 960
ctsm1.0.dev049 erik 06/23/2019 Update mosart and intel to intel-19 on cheyenne
diff --git a/tools/mksurfdata_map/src/CMakeLists.txt b/tools/mksurfdata_map/src/CMakeLists.txt
index 4cd91e8a39..0cf2ca90b0 100644
--- a/tools/mksurfdata_map/src/CMakeLists.txt
+++ b/tools/mksurfdata_map/src/CMakeLists.txt
@@ -17,6 +17,7 @@ list(APPEND mksurfdat_sources
mkpctPftTypeMod.F90
mkutilsMod.F90
mkpftUtilsMod.F90
+ mksoilUtilsMod.F90
mkvarctl.F90
mkvarpar.F90
shr_const_mod.F90
diff --git a/tools/mksurfdata_map/src/Srcfiles b/tools/mksurfdata_map/src/Srcfiles
index 25e50d3f95..b126a5b508 100644
--- a/tools/mksurfdata_map/src/Srcfiles
+++ b/tools/mksurfdata_map/src/Srcfiles
@@ -12,6 +12,7 @@ mkpftConstantsMod.F90
mkpftUtilsMod.F90
mkpctPftTypeMod.F90
mkpftMod.F90
+mksoilUtilsMod.F90
mksoilMod.F90
mkvocefMod.F90
mksurfdat.F90
diff --git a/tools/mksurfdata_map/src/mkpftMod.F90 b/tools/mksurfdata_map/src/mkpftMod.F90
index c99d66c9ac..fb86543944 100644
--- a/tools/mksurfdata_map/src/mkpftMod.F90
+++ b/tools/mksurfdata_map/src/mkpftMod.F90
@@ -524,7 +524,8 @@ subroutine mkpft(ldomain, mapfname, fpft, ndiag, allow_no_crops, &
call gridmap_areaave(tgridmap, pctcrop_i, pctcrop_o, nodata=0._r8)
do m = 0, num_natpft
- call gridmap_areaave_scs(tgridmap, pct_nat_pft_i(:,m), pct_nat_pft_o(:,m), nodata=0._r8,src_wt=pctnatveg_i*0.01_r8,dst_wt=pctnatveg_o*0.01_r8)
+ call gridmap_areaave_scs(tgridmap, pct_nat_pft_i(:,m), pct_nat_pft_o(:,m), &
+ nodata=0._r8,src_wt=pctnatveg_i*0.01_r8,dst_wt=pctnatveg_o*0.01_r8)
do no = 1,ns_o
if (pctlnd_o(no) < 1.0e-6 .or. pctnatveg_o(no) < 1.0e-6) then
if (m == 0) then
@@ -536,7 +537,8 @@ subroutine mkpft(ldomain, mapfname, fpft, ndiag, allow_no_crops, &
enddo
end do
do m = 1, num_cft
- call gridmap_areaave_scs(tgridmap, pct_cft_i(:,m), pct_cft_o(:,m), nodata=0._r8,src_wt=pctcrop_i*0.01_r8,dst_wt=pctcrop_o*0.01_r8)
+ call gridmap_areaave_scs(tgridmap, pct_cft_i(:,m), pct_cft_o(:,m), &
+ nodata=0._r8,src_wt=pctcrop_i*0.01_r8,dst_wt=pctcrop_o*0.01_r8)
do no = 1,ns_o
if (pctlnd_o(no) < 1.0e-6 .or. pctcrop_o(no) < 1.0e-6) then
if (m == 1) then
diff --git a/tools/mksurfdata_map/src/mksoilMod.F90 b/tools/mksurfdata_map/src/mksoilMod.F90
index 1192863737..e73c611308 100644
--- a/tools/mksurfdata_map/src/mksoilMod.F90
+++ b/tools/mksurfdata_map/src/mksoilMod.F90
@@ -15,6 +15,7 @@ module mksoilMod
use shr_kind_mod, only : r8 => shr_kind_r8, r4=>shr_kind_r4
use shr_sys_mod , only : shr_sys_flush
use mkdomainMod , only : domain_checksame
+ use mksoilUtilsMod, only : mkrank, dominant_soil_color
implicit none
SAVE
@@ -43,7 +44,6 @@ module mksoilMod
! !PRIVATE DATA MEMBERS:
!
! !PRIVATE MEMBER FUNCTIONS:
- private :: mkrank
private :: mksoiltexInit ! Soil texture Initialization
private :: mksoilcolInit ! Soil color Initialization
private :: mksoilfmaxInit ! Soil fmax Initialization
@@ -597,21 +597,15 @@ subroutine mksoilcol(ldomain, mapfname, datfname, ndiag, &
!EOP
type(gridmap_type) :: tgridmap
type(domain_type) :: tdomain ! local domain
- integer, parameter :: num=2 ! set soil mapunit number
- integer :: wsti(num) ! index to 1st and 2nd largest wst
- real(r8), allocatable :: wst(:,:) ! overlap weights, by surface type
real(r8), allocatable :: gast_i(:) ! global area, by surface type
real(r8), allocatable :: gast_o(:) ! global area, by surface type
integer , allocatable :: soil_color_i(:) ! input grid: BATS soil color
- integer , allocatable :: color(:) ! 0: none; 1: some
- real(r8) :: wt ! map overlap weight
real(r8) :: sum_fldi ! global sum of dummy input fld
real(r8) :: sum_fldo ! global sum of dummy output fld
character(len=35), allocatable :: col(:) ! name of each color
- integer :: k,l,n,m,ni,no,ns_i,ns_o ! indices
+ integer :: k,l,m,ni,no,ns_i,ns_o ! indices
integer :: ncid,dimid,varid ! input netCDF id's
integer :: ier ! error status
- integer :: miss = 99999 ! missing data indicator
real(r8) :: relerr = 0.00001 ! max error: sum overlap wts ne 1
character(len=32) :: subname = 'mksoilcol'
!-----------------------------------------------------------------------
@@ -703,62 +697,14 @@ subroutine mksoilcol(ldomain, mapfname, datfname, ndiag, &
call domain_checksame( tdomain, ldomain, tgridmap )
- ! find area of overlap for each soil color for each no
-
- allocate(wst(0:nsoicol,ns_o))
- wst(0:nsoicol,:) = 0
- allocate(color(ns_o))
- color(:) = 0
-
- ! TODO: need to do a loop to determine
- ! the maximum number of over lap cells throughout the grid
- ! first get an array that is novr(ns_o) and fill this in - then set
- ! maxovr - to max(novr) - then allocate the array wst to be size of
- ! maxovr,ns_o or 0:nsoilcol,ns_o
-
- do n = 1,tgridmap%ns
- ni = tgridmap%src_indx(n)
- no = tgridmap%dst_indx(n)
- wt = tgridmap%wovr(n)
- k = soil_color_i(ni) * tdomain%mask(ni)
- wst(k,no) = wst(k,no) + wt
- if (k>0 .and. wst(k,no)>0.) then
- color(no) = 1
- wst(0,no) = 0.0
- end if
- enddo
-
- soil_color_o(:) = 0
- do no = 1,ns_o
-
- ! Rank non-zero weights by color type. wsti(1) is the most extensive
- ! color type.
-
- if (color(no) == 1) then
- call mkrank (nsoicol, wst(0:nsoicol,no), miss, wsti, num)
- soil_color_o(no) = wsti(1)
- end if
-
- ! If land but no color, set color to 15 (in older dataset generic
- ! soil color 4)
-
- if (nsoicol == 8) then
- if (soil_color_o(no)==0) soil_color_o(no) = 4
- else if (nsoicol == 20) then
- if (soil_color_o(no)==0) soil_color_o(no) = 15
- end if
-
- ! Error checks
+ ! Determine dominant soil color for each output cell
- if (soil_color_o(no) < 0 .or. soil_color_o(no) > nsoicol) then
- write (6,*) 'MKSOILCOL error: land model soil color = ', &
- soil_color_o(no),' is not valid for lon,lat = ',no
- call abort()
- end if
-
- enddo
- deallocate (wst)
- deallocate (color)
+ call dominant_soil_color( &
+ tgridmap = tgridmap, &
+ mask_i = tdomain%mask, &
+ soil_color_i = soil_color_i, &
+ nsoicol = nsoicol, &
+ soil_color_o = soil_color_o)
! Global sum of output field
@@ -983,96 +929,6 @@ subroutine mkorganic(ldomain, mapfname, datfname, ndiag, organic_o)
end subroutine mkorganic
-!-----------------------------------------------------------------------
-!BOP
-!
-! !ROUTINE: mkrank
-!
-! !INTERFACE:
-subroutine mkrank (n, a, miss, iv, num)
-!
-! !DESCRIPTION:
-! Return indices of largest [num] values in array [a]. Private method
-! only used for soil color and soil texture.
-!
-! !USES:
-!
-! !ARGUMENTS:
- implicit none
- integer , intent(in) :: n !array length
- real(r8), intent(in) :: a(0:n) !array to be ranked
- integer , intent(in) :: miss !missing data value
- integer , intent(in) :: num !number of largest values requested
- integer , intent(out):: iv(num) !index to [num] largest values in array [a]
-!
-! !CALLED FROM:
-! subroutine mksoilcol
-! subroutine mksoiltex
-!
-! !REVISION HISTORY:
-! Author: Gordon Bonan
-!
-!
-! !LOCAL VARIABLES:
-!EOP
- real(r8) a_max !maximum value in array
- integer i !array index
- real(r8) delmax !tolerance for finding if larger value
- integer m !do loop index
- integer k !do loop index
- logical exclude !true if data value has already been chosen
-!-----------------------------------------------------------------------
-
- delmax = 1.e-06
-
- ! Find index of largest non-zero number
-
- iv(1) = miss
- a_max = -9999.
-
- do i = 0, n
- if (a(i)>0. .and. (a(i)-a_max)>delmax) then
- a_max = a(i)
- iv(1) = i
- end if
- end do
-
- ! iv(1) = miss indicates no values > 0. this is an error
-
- if (iv(1) == miss) then
- write (6,*) 'MKRANK error: iv(1) = missing'
- call abort()
- end if
-
- ! Find indices of the next [num]-1 largest non-zero number.
- ! iv(m) = miss if there are no more values > 0
-
- do m = 2, num
- iv(m) = miss
- a_max = -9999.
- do i = 0, n
-
- ! exclude if data value has already been chosen
-
- exclude = .false.
- do k = 1, m-1
- if (i == iv(k)) exclude = .true.
- end do
-
- ! if not already chosen, see if it is the largest of
- ! the remaining values
-
- if (.not. exclude) then
- if (a(i)>0. .and. (a(i)-a_max)>delmax) then
- a_max = a(i)
- iv(m) = i
- end if
- end if
- end do
- end do
-
-end subroutine mkrank
-
!-----------------------------------------------------------------------
!BOP
!
diff --git a/tools/mksurfdata_map/src/mksoilUtilsMod.F90 b/tools/mksurfdata_map/src/mksoilUtilsMod.F90
new file mode 100644
index 0000000000..e8c672590d
--- /dev/null
+++ b/tools/mksurfdata_map/src/mksoilUtilsMod.F90
@@ -0,0 +1,224 @@
+module mksoilUtilsMod
+
+ !-----------------------------------------------------------------------
+ !BOP
+ !
+ ! !MODULE: mksoilUtils
+ !
+ ! !DESCRIPTION:
+ ! Lower-level utilities used in making soil data.
+ !
+ ! These are separated out from mksoilMod mainly as an aid to testing.
+ !
+ ! !REVISION HISTORY:
+ ! Author: Bill Sacks
+ !
+ !-----------------------------------------------------------------------
+ !!USES:
+ use shr_kind_mod, only : r8 => shr_kind_r8
+ use mkgridmapMod, only : gridmap_type
+
+ implicit none
+ private
+
+ !
+ ! !PUBLIC MEMBER FUNCTIONS:
+ !
+ public :: dominant_soil_color
+ public :: mkrank
+
+ !
+ ! !PRIVATE MEMBER FUNCTIONS:
+ !
+
+ !EOP
+ !===============================================================
+contains
+ !===============================================================
+
+ !-----------------------------------------------------------------------
+ subroutine dominant_soil_color(tgridmap, mask_i, soil_color_i, nsoicol, soil_color_o)
+ !
+ ! !DESCRIPTION:
+ ! Determine the dominant soil color in each output cell
+ !
+ ! !ARGUMENTS:
+ type(gridmap_type) , intent(in) :: tgridmap
+ integer , intent(in) :: mask_i(:) ! input grid: land mask (1 = land, 0 = ocean)
+ integer , intent(in) :: soil_color_i(:) ! input grid: BATS soil color
+ integer , intent(in) :: nsoicol ! number of soil colors
+ integer , intent(out) :: soil_color_o(:) ! output grid: soil color classes
+ !
+ ! !LOCAL VARIABLES:
+ integer, parameter :: num = 2 ! set soil mapunit number
+ integer :: wsti(num) ! index to 1st and 2nd largest wst
+ integer :: k, n, ni, no, ns_i, ns_o
+ real(r8) :: wt ! map overlap weight
+ real(r8), allocatable :: wst(:,:) ! overlap weights, by surface type
+ logical :: has_color ! whether this grid cell has non-zero color
+ integer, parameter :: miss = 99999 ! missing data indicator
+
+ character(len=*), parameter :: subname = 'dominant_soil_color'
+ !-----------------------------------------------------------------------
+
+ ns_i = size(mask_i)
+ if (size(soil_color_i) /= ns_i) then
+ write(6,*) subname, ' ERROR: size of soil_color_i should match size of mask_i'
+ write(6,*) 'size(mask_i), size(soil_color_i) = ', &
+ size(mask_i), size(soil_color_i)
+ call abort()
+ end if
+
+ ! find area of overlap for each soil color for each no
+
+ ns_o = size(soil_color_o)
+ allocate(wst(0:nsoicol,ns_o))
+ wst(0:nsoicol,:) = 0
+
+ ! TODO: need to do a loop to determine
+ ! the maximum number of over lap cells throughout the grid
+ ! first get an array that is novr(ns_o) and fill this in - then set
+ ! maxovr - to max(novr) - then allocate the array wst to be size of
+ ! maxovr,ns_o or 0:nsoilcol,ns_o
+
+ do n = 1,tgridmap%ns
+ ni = tgridmap%src_indx(n)
+ no = tgridmap%dst_indx(n)
+ wt = tgridmap%wovr(n)
+ k = soil_color_i(ni) * mask_i(ni)
+ wst(k,no) = wst(k,no) + wt
+ enddo
+
+ soil_color_o(:) = 0
+ do no = 1,ns_o
+
+ ! If the output cell has any non-zero-colored inputs, then set the weight of
+ ! zero-colored inputs to 0, to ensure that the zero-color is NOT dominant.
+ if (any(wst(1:nsoicol,no) > 0.)) then
+ has_color = .true.
+ wst(0,no) = 0.0
+ else
+ has_color = .false.
+ end if
+
+ ! Rank non-zero weights by color type. wsti(1) is the most extensive
+ ! color type.
+
+ if (has_color) then
+ call mkrank (nsoicol, wst(0:nsoicol,no), miss, wsti, num)
+ soil_color_o(no) = wsti(1)
+ end if
+
+ ! If land but no color, set color to 15 (in older dataset generic
+ ! soil color 4)
+
+ if (nsoicol == 8) then
+ if (soil_color_o(no)==0) then
+ soil_color_o(no) = 4
+ end if
+ else if (nsoicol == 20) then
+ if (soil_color_o(no)==0) then
+ soil_color_o(no) = 15
+ end if
+ else
+ write(6,*) 'MKSOILCOL error: unhandled nsoicol: ', nsoicol
+ call abort()
+ end if
+
+ ! Error checks
+
+ if (soil_color_o(no) < 0 .or. soil_color_o(no) > nsoicol) then
+ write (6,*) 'MKSOILCOL error: land model soil color = ', &
+ soil_color_o(no),' is not valid for lon,lat = ',no
+ call abort()
+ end if
+
+ end do
+
+ deallocate (wst)
+
+ end subroutine dominant_soil_color
+
+
+ !-----------------------------------------------------------------------
+ !BOP
+ !
+ ! !ROUTINE: mkrank
+ !
+ ! !INTERFACE:
+ subroutine mkrank (n, a, miss, iv, num)
+ !
+ ! !DESCRIPTION:
+ ! Return indices of largest [num] values in array [a].
+ !
+ ! !ARGUMENTS:
+ integer , intent(in) :: n !array length
+ real(r8), intent(in) :: a(0:n) !array to be ranked
+ integer , intent(in) :: miss !missing data value
+ integer , intent(in) :: num !number of largest values requested
+ integer , intent(out):: iv(num) !index to [num] largest values in array [a]
+ !
+ ! !REVISION HISTORY:
+ ! Author: Gordon Bonan
+ !
+ ! !LOCAL VARIABLES:
+ !EOP
+ real(r8) a_max !maximum value in array
+ integer i !array index
+ real(r8) delmax !tolerance for finding if larger value
+ integer m !do loop index
+ integer k !do loop index
+ logical exclude !true if data value has already been chosen
+ !-----------------------------------------------------------------------
+
+ delmax = 1.e-06
+
+ ! Find index of largest non-zero number
+
+ iv(1) = miss
+ a_max = -9999.
+
+ do i = 0, n
+ if (a(i)>0. .and. (a(i)-a_max)>delmax) then
+ a_max = a(i)
+ iv(1) = i
+ end if
+ end do
+
+ ! iv(1) = miss indicates no values > 0. this is an error
+
+ if (iv(1) == miss) then
+ write (6,*) 'MKRANK error: iv(1) = missing'
+ call abort()
+ end if
+
+ ! Find indices of the next [num]-1 largest non-zero number.
+ ! iv(m) = miss if there are no more values > 0
+
+ do m = 2, num
+ iv(m) = miss
+ a_max = -9999.
+ do i = 0, n
+
+ ! exclude if data value has already been chosen
+
+ exclude = .false.
+ do k = 1, m-1
+ if (i == iv(k)) exclude = .true.
+ end do
+
+ ! if not already chosen, see if it is the largest of
+ ! the remaining values
+
+ if (.not. exclude) then
+ if (a(i)>0. .and. (a(i)-a_max)>delmax) then
+ a_max = a(i)
+ iv(m) = i
+ end if
+ end if
+ end do
+ end do
+
+ end subroutine mkrank
+
+end module mksoilUtilsMod
diff --git a/tools/mksurfdata_map/src/test/CMakeLists.txt b/tools/mksurfdata_map/src/test/CMakeLists.txt
index 690bee396e..5cfb20f21b 100644
--- a/tools/mksurfdata_map/src/test/CMakeLists.txt
+++ b/tools/mksurfdata_map/src/test/CMakeLists.txt
@@ -1,4 +1,5 @@
add_subdirectory(mkpctPftType_test)
add_subdirectory(mkpftUtils_test)
add_subdirectory(mkgridmap_test)
-add_subdirectory(mkindexmap_test)
\ No newline at end of file
+add_subdirectory(mkindexmap_test)
+add_subdirectory(mksoilUtils_test)
\ No newline at end of file
diff --git a/tools/mksurfdata_map/src/test/mksoilUtils_test/CMakeLists.txt b/tools/mksurfdata_map/src/test/mksoilUtils_test/CMakeLists.txt
new file mode 100644
index 0000000000..4d94b8114b
--- /dev/null
+++ b/tools/mksurfdata_map/src/test/mksoilUtils_test/CMakeLists.txt
@@ -0,0 +1,7 @@
+set (pfunit_sources
+ test_dominant_soil_color.pf)
+
+create_pFUnit_test(mksoilUtils test_mksoilUtils_exe
+ "${pfunit_sources}" "")
+
+target_link_libraries(test_mksoilUtils_exe mksurfdat)
\ No newline at end of file
diff --git a/tools/mksurfdata_map/src/test/mksoilUtils_test/test_dominant_soil_color.pf b/tools/mksurfdata_map/src/test/mksoilUtils_test/test_dominant_soil_color.pf
new file mode 100644
index 0000000000..b506549d87
--- /dev/null
+++ b/tools/mksurfdata_map/src/test/mksoilUtils_test/test_dominant_soil_color.pf
@@ -0,0 +1,140 @@
+module test_dominant_soil_color
+
+ ! Tests of mksoilUtilsMod: dominant_soil_color
+
+ use pfunit_mod
+ use mksoilUtilsMod
+ use shr_kind_mod , only : r8 => shr_kind_r8
+ use mkgridmapMod, only : gridmap_type, gridmap_clean, for_test_create_gridmap
+
+ implicit none
+
+ @TestCase
+ type, extends(TestCase) :: tdsc
+ type(gridmap_type) :: gridmap
+ contains
+ procedure :: setUp
+ procedure :: tearDown
+ procedure :: createGridmap1dst
+ end type tdsc
+
+ real(r8), parameter :: tol = 1.e-13_r8
+
+contains
+
+ subroutine setUp(this)
+ class(tdsc), intent(inout) :: this
+ end subroutine setUp
+
+ subroutine tearDown(this)
+ class(tdsc), intent(inout) :: this
+ call gridmap_clean(this%gridmap)
+ end subroutine tearDown
+
+ subroutine createGridmap1dst(this, wovr)
+ ! Create this%gridmap with a single destination point
+ class(tdsc), intent(inout) :: this
+ real(r8), intent(in) :: wovr(:) ! overlap weights
+
+ integer :: i
+ integer :: npts
+ integer :: src_indx(size(wovr))
+ integer :: dst_indx(size(wovr))
+
+ dst_indx(:) = 1
+ npts = size(wovr)
+ src_indx(:) = [(i, i = 1, npts)]
+
+ call for_test_create_gridmap(this%gridmap, na = npts, nb = 1, ns = npts, &
+ src_indx = src_indx, dst_indx = dst_indx, wovr = wovr)
+ end subroutine createGridmap1dst
+
+ @Test
+ subroutine equalWeights(this)
+ ! Four inputs with equal weight; two of one class, one of each of two other classes
+ class(tdsc), intent(inout) :: this
+ integer :: mask_i(4)
+ integer :: soil_color_i(4)
+ integer :: soil_color_o(1)
+
+ call this%createGridmap1dst([0.25_r8, 0.25_r8, 0.25_r8, 0.25_r8])
+ mask_i(:) = 1
+ soil_color_i(:) = [1, 2, 2, 3]
+
+ call dominant_soil_color(this%gridmap, mask_i, soil_color_i, 20, soil_color_o)
+
+ @assertEqual(2, soil_color_o(1))
+ end subroutine equalWeights
+
+ @Test
+ subroutine inequalWeights(this)
+ ! Four inputs with inequal weight
+ class(tdsc), intent(inout) :: this
+ integer :: mask_i(4)
+ integer :: soil_color_i(4)
+ integer :: soil_color_o(1)
+
+ call this%createGridmap1dst([0.5_r8, 0.2_r8, 0.2_r8, 0.1_r8])
+ mask_i(:) = 1
+ soil_color_i(:) = [3, 1, 1, 2]
+
+ call dominant_soil_color(this%gridmap, mask_i, soil_color_i, 20, soil_color_o)
+
+ @assertEqual(3, soil_color_o(1))
+ end subroutine inequalWeights
+
+ @Test
+ subroutine noColor(this)
+ ! No color in input
+ class(tdsc), intent(inout) :: this
+ integer :: mask_i(4)
+ integer :: soil_color_i(4)
+ integer :: soil_color_o(1)
+
+ call this%createGridmap1dst([0.25_r8, 0.25_r8, 0.25_r8, 0.25_r8])
+ ! Some points are inside the mask with color = 0, other points are outside the mask
+ mask_i(:) = [1, 0, 0, 1]
+ soil_color_i(:) = [0, 1, 1, 0]
+
+ call dominant_soil_color(this%gridmap, mask_i, soil_color_i, 20, soil_color_o)
+
+ @assertEqual(15, soil_color_o(1))
+ end subroutine noColor
+
+ @Test
+ subroutine noColorInFirstPoints(this)
+ ! No color in the first points, but a color in the last point
+ class(tdsc), intent(inout) :: this
+ integer :: mask_i(4)
+ integer :: soil_color_i(4)
+ integer :: soil_color_o(1)
+
+ call this%createGridmap1dst([0.25_r8, 0.25_r8, 0.25_r8, 0.25_r8])
+ ! Some points are inside the mask with color = 0, other points are outside the mask
+ mask_i(:) = 1
+ soil_color_i(:) = [0, 0, 0, 1]
+
+ call dominant_soil_color(this%gridmap, mask_i, soil_color_i, 20, soil_color_o)
+
+ @assertEqual(1, soil_color_o(1))
+ end subroutine noColorInFirstPoints
+
+ @Test
+ subroutine noColorInLastPoints(this)
+ ! No color in the last points, but a color in the first point
+ class(tdsc), intent(inout) :: this
+ integer :: mask_i(4)
+ integer :: soil_color_i(4)
+ integer :: soil_color_o(1)
+
+ call this%createGridmap1dst([0.25_r8, 0.25_r8, 0.25_r8, 0.25_r8])
+ ! Some points are inside the mask with color = 0, other points are outside the mask
+ mask_i(:) = 1
+ soil_color_i(:) = [1, 0, 0, 0]
+
+ call dominant_soil_color(this%gridmap, mask_i, soil_color_i, 20, soil_color_o)
+
+ @assertEqual(1, soil_color_o(1))
+ end subroutine noColorInLastPoints
+
+end module test_dominant_soil_color