add soil data and update dependencies

pyproject.toml

@@ -8,31 +8,31 @@ license="BSD-3-Clause"
 packages = [{include = "terrapyn"}]
 
 [tool.poetry.dependencies]
-python = ">=3.11,<3.12"
+python = ">=3.10,<3.12"
 xarray = "^2024"
-pandas = "^2.0.0"
-geopandas = "^0.14"
-numpy = "^1.26.0"
-scipy = "^1.12.0"
+pandas = "^2"
+geopandas = "*"
+numpy = "*"
+scipy = "*"
 dask = "^2024"
-matplotlib = "^3.8.0"
-bottleneck = "^1.3.0"
-shapely = "^2.0.0"
-netcdf4 = "^1.6.0"
-ipykernel = "^6.29.0"
+matplotlib = "^3"
+bottleneck = "^1"
+shapely = "^2"
+netcdf4 = "*"
+ipykernel = "^6"
 pyarrow = "^15"
-fiona = "^1.9.0"
-tqdm = "^4.66.0"
+fiona = "*"
+tqdm = "*"
 jupyterlab = "^4.1"
-rioxarray = "^0.15"
-earthengine-api = "^0.1.394"
-polygon-geohasher = "^0.0.1"
-geemap = "^0.32.0"
+rioxarray = "*"
+earthengine-api = "^0.1"
+polygon-geohasher = "*"
+geemap = "*"
 
 [tool.poetry.group.dev.dependencies]
 pytest = "^7.3"
 pytest-cov = "^4.0"
-isort = "^5.12"
+isort = "^5"
 black = "^23.3"
 pre-commit = "^3.2"
 freezegun = "^1.2"

terrapyn/ee/data.py

@@ -109,8 +109,8 @@ def soilgrids(
         # All images have the same horizons
         horizon_top = [0, 5, 15, 30, 60, 100]
         horizon_bottom = [5, 15, 30, 60, 100, 200]
-        horizon_weights = tp.ee.stats._horizon_weights(
-            horizon_top=horizon_top, horizon_bottom=horizon_bottom, topsoil_depth=depth
+        horizon_weights = tp.ee.stats._layer_weights(
+            layer_tops=horizon_top, layer_bottoms=horizon_bottom, total_depth=depth
         )
 
         mean_images = []
@@ -122,3 +122,98 @@ def soilgrids(
         return ee.Image(mean_images)
     else:
         return ee.Image(images)
+
+
+def usda_soil_class(sand: ee.Image = None, silt: ee.Image = None, clay: ee.Image = None) -> ee.Image:
+    """
+    Takes images of sand/silt/clay percentage and returns USDA soil class (texture triangle)
+    https://www.nrcs.usda.gov/sites/default/files/2022-09/The-Soil-Survey-Manual.pdf page 122
+
+    Output Image has the same mask as `sand` image.
+
+    1  "Sa": "sand"
+    2  "LoSa": "loamy sand"
+    3  "SaLo": "sandy loam"
+    4  "Lo": "loam"
+    5  "SiLo": "silty loam"
+    6  "Si": "silt"
+    7  "SaClLo": "sandy clay loam"
+    8  "ClLo": "clay loam"
+    9  "SiClLo": "silty clay loam"
+    10  "SaCl": "sandy clay"
+    11  "SiCl": "silty clay"
+    12  "Cl": "clay"
+
+
+    Args:
+        sand: Image of sand percentage
+        silt: Image of silt percentage
+        clay: Image of clay percentage
+
+    Returns:
+        Image of soil class
+    """
+    # Initialize soil class image with crs and crsTransform from `sand` image,
+    # and set to be the same as that applied to `sand` image
+    soil_class = (
+        ee.Image(0)
+        .setDefaultProjection(crs=sand.projection().wkt(), crsTransform=sand.getInfo()["bands"][0]["crs_transform"])
+        .updateMask(sand)
+    )
+
+    # Sand - Material has more than 85 percent sand, and the percentage of silt plus
+    # 1.5 times the percentage of clay is less than 15.
+    soil_class = soil_class.where(sand.gt(85).And(silt.add(clay.multiply(1.5)).lt(15)), 1)
+    # Loamy sands - Material has between 70 and 90 percent sand, the
+    # percentage of silt plus 1.5 times the percentage of clay is 15 or more, and
+    # the percentage of silt plus twice the percentage of clay is less than 30.
+    soil_class = soil_class.where(
+        sand.gte(70).And(sand.lte(90)).And(silt.add(clay.multiply(1.5)).gte(15)).And(silt.add(clay.multiply(2)).lt(30)),
+        2,
+    )
+    # Sandy loams - Material has 7 to less than 20 percent clay and more
+    # than 52 percent sand, and the percentage of silt plus twice the percentage
+    # of clay is 30 or more; OR material has less than 7 percent clay and less
+    # than 50 percent silt, and the percentage of silt plus twice the percentage
+    # of clay is 30 or more.
+    soil_class = soil_class.where(
+        (clay.gte(7).And(clay.lt(20)).And(sand.gt(52)).And(silt.add(clay.multiply(2)).gte(30))).Or(
+            (clay.lt(7).And(silt.lt(50)).And(silt.add(clay.multiply(2)).gte(30)))
+        ),
+        3,
+    )
+    # Loam - Material has 7 to less than 27 percent clay, 28 to less than
+    # 50 percent silt, and 52 percent or less sand.
+    soil_class = soil_class.where(clay.gte(7).And(clay.lt(27)).And(silt.gte(28)).And(silt.lt(50)).And(sand.lte(52)), 4)
+    # Silt loam - Material has 50 percent or more silt and 12 to less than
+    # 27 percent clay; OR material has 50 to less than 80 percent silt and less
+    # than 12 percent clay.
+    soil_class = soil_class.where(
+        (silt.gte(50).And(clay.gte(12)).And(clay.lt(27))).Or((silt.gte(50).And(silt.lt(80)).And(clay.lt(12)))), 5
+    )
+    # Silt - Material has 80 percent or more silt and less than 12 percent
+    # clay.
+    soil_class = soil_class.where(silt.gte(80).And(clay.lt(12)), 6)
+    # Sandy clay loam - Material has 20 to less than 35 percent clay, less
+    # than 28 percent silt, and more than 45 percent sand.
+    soil_class = soil_class.where(clay.gte(20).And(clay.lt(35)).And(silt.lt(28)).And(sand.gt(45)), 7)
+    # Clay loam - Material has 27 to less than 40 percent clay and more
+    # than 20 to 45 percent sand.
+    soil_class = soil_class.where((clay.gte(27).And(clay.lt(40)).And(sand.gt(20)).And(sand.lte(45))), 8)
+    # Silty clay loam - Material has 27 to less than 40 percent clay and 20
+    # percent or less sand.
+    soil_class = soil_class.where((clay.gte(27).And(clay.lt(40)).And(sand.lte(20))), 9)
+    # Sandy clay - Material has 35 percent or more clay and more than
+    # 45 percent sand.
+    soil_class = soil_class.where(clay.gte(35).And(sand.gt(45)), 10)
+    # Silty clay - Material has 40 percent or more clay and 40 percent or
+    # more silt.
+    soil_class = soil_class.where(clay.gte(40).And(silt.gte(40)), 11)
+    # Clay - Material has 40 percent or more clay, 45 percent or less sand,
+    # and less than 40 percent silt.
+    soil_class = soil_class.where(clay.gte(40).And(sand.lte(45)).And(silt.lt(40)), 12)
+
+    # Mask out null values
+    soil_class = soil_class.selfMask()
+
+    return soil_class.rename("texture")

terrapyn/ee/stats.py

@@ -236,11 +236,13 @@ def weighted_mean(img, weight_dict, output_bandname="mean"):
     return tp.ee.utils.scale_image_bands(img, weight_dict=weights_dict).reduce("sum").rename(output_bandname)
 
 
-def _horizon_weights(horizon_top, horizon_bottom, topsoil_depth=30):
-    """calculate horizon weight based on depth range overlap"""
+def _layer_weights(layer_tops: T.Iterable, layer_bottoms: T.Iterable, total_depth: float = 30) -> T.List:
+    """
+    Calculate a weighting for each layer based on thickness and range overlap, typically used for soil horizons
+    """
     weights = []
-    for top, bottom in zip(horizon_top, horizon_bottom):
-        overlap = min(topsoil_depth, bottom) - max(0, top)
+    for top, bottom in zip(layer_tops, layer_bottoms):
+        overlap = min(total_depth, bottom) - max(0, top)
         # if depth range doesn't overlap, weight = 0
         weights.append(overlap if overlap > 0 else 0)
     return weights

terrapyn/params/__init__.py

@@ -1,3 +1,3 @@
-from . import etp, solar
+from . import etp, soil, solar
 
-__all__ = ["solar", "etp"]
+__all__ = ["solar", "etp", "soil"]

terrapyn/params/soil.py

@@ -0,0 +1,94 @@
+# USDA soil textures: https://www.nrcs.usda.gov/sites/default/files/2022-09/The-Soil-Survey-Manual.pdf page 122
+USDA_SOIL_CLASS = {
+    1: "sand",
+    2: "loamy sand",
+    3: "sandy loam",
+    4: "loam",
+    5: "silt loam",
+    6: "silt",
+    7: "sandy clay loam",
+    8: "clay loam",
+    9: "silty clay loam",
+    10: "sandy clay",
+    11: "silty clay",
+    12: "clay",
+}
+
+
+def usda_soil_class(sand=None, silt=None, clay=None):
+    """
+    Takes floats of sand/silt/clay percentage and returns USDA soil class (texture triangle)
+    https://www.nrcs.usda.gov/sites/default/files/2022-09/The-Soil-Survey-Manual.pdf page 122
+
+    1  "Sa": "sand"
+    2  "LoSa": "loamy sand"
+    3  "SaLo": "sandy loam"
+    4  "Lo": "loam"
+    5  "SiLo": "silty loam"
+    6  "Si": "silt"
+    7  "SaClLo": "sandy clay loam"
+    8  "ClLo": "clay loam"
+    9  "SiClLo": "silty clay loam"
+    10  "SaCl": "sandy clay"
+    11  "SiCl": "silty clay"
+    12  "Cl": "clay"
+
+    Returns:
+        integer soil class
+    """
+    # Sand - Material has more than 85 percent sand, and the percentage of silt plus
+    # 1.5 times the percentage of clay is less than 15.
+    if (sand > 85) and (silt + 1.5 * clay < 15):
+        soil_class = 1
+    # Loamy sands - Material has between 70 and 90 percent sand, the
+    # percentage of silt plus 1.5 times the percentage of clay is 15 or more, and
+    # the percentage of silt plus twice the percentage of clay is less than 30.
+    elif (sand >= 70) and (sand <= 90) and (silt + 1.5 * clay >= 15) and (silt + 2 * clay < 30):
+        soil_class = 2
+    # Sandy loams - Material has 7 to less than 20 percent clay and more
+    # than 52 percent sand, and the percentage of silt plus twice the percentage
+    # of clay is 30 or more; OR material has less than 7 percent clay and less
+    # than 50 percent silt, and the percentage of silt plus twice the percentage
+    # of clay is 30 or more.
+    elif ((clay >= 7) and (clay < 20) and (sand > 52) and (silt + 2 * clay >= 30)) or (
+        (clay < 7) and (silt < 50) and (silt + 2 * clay >= 30)
+    ):
+        soil_class = 3
+    # Loam - Material has 7 to less than 27 percent clay, 28 to less than
+    # 50 percent silt, and 52 percent or less sand.
+    elif (clay >= 7) and (clay < 27) and (silt >= 28) and (silt < 50) and (sand <= 52):
+        soil_class = 4
+    # Silt loam - Material has 50 percent or more silt and 12 to less than
+    # 27 percent clay; OR material has 50 to less than 80 percent silt and less
+    # than 12 percent clay.
+    elif ((silt >= 50) and (clay >= 12) and (clay < 27)) or ((silt >= 50) and (silt < 80) and (clay < 12)):
+        soil_class = 5
+    # Silt - Material has 80 percent or more silt and less than 12 percent
+    # clay.
+    elif (silt >= 80) and (clay < 12):
+        soil_class = 6
+    # Sandy clay loam - Material has 20 to less than 35 percent clay, less
+    # than 28 percent silt, and more than 45 percent sand.
+    elif (clay >= 20) and (clay < 35) and (silt < 28) and (sand > 45):
+        soil_class = 7
+    # Clay loam - Material has 27 to less than 40 percent clay and more
+    # than 20 to 45 percent sand.
+    elif (clay >= 27) and (clay < 40) and (sand > 20) and (sand <= 45):
+        soil_class = 8
+    # Silty clay loam - Material has 27 to less than 40 percent clay and 20
+    # percent or less sand.
+    elif (clay >= 27) and (clay < 40) and (sand <= 20):
+        soil_class = 9
+    # Sandy clay - Material has 35 percent or more clay and more than
+    # 45 percent sand.
+    elif (clay >= 35) and (sand > 45):
+        soil_class = 10
+    # Silty clay - Material has 40 percent or more clay and 40 percent or
+    # more silt.
+    elif (clay >= 40) and (silt >= 40):
+        soil_class = 11
+    # Clay - Material has 40 percent or more clay, 45 percent or less sand,
+    # and less than 40 percent silt.
+    elif (clay >= 40) and (sand <= 45) and (silt < 40):
+        soil_class = 12
+    return soil_class