diff --git a/.github/workflows/R-CMD-check.yaml b/.github/workflows/R-CMD-check.yaml
index 8a9d8542..f2aec811 100644
--- a/.github/workflows/R-CMD-check.yaml
+++ b/.github/workflows/R-CMD-check.yaml
@@ -36,7 +36,7 @@ jobs:
       GITHUB_PAT: ${{ secrets.GITHUB_TOKEN }}
 
     steps:
-      - uses: actions/checkout@v3
+      - uses: actions/checkout@v4
 
       - uses: r-lib/actions/setup-r@v2
         with:
diff --git a/R/SoilProfileCollection-setters.R b/R/SoilProfileCollection-setters.R
index 82442a28..b262a98a 100644
--- a/R/SoilProfileCollection-setters.R
+++ b/R/SoilProfileCollection-setters.R
@@ -544,9 +544,10 @@ setReplaceMethod("horizons", signature(object = "SoilProfileCollection"),
   if (!inherits(value, "data.frame"))
     stop("new horizon data input value must inherit from data.frame", call. = FALSE)
 
-  # allow short-circuit
-  if (all(colnames(value) %in% horizonNames(object)) &
-      all(c(idn, hdn, hzd) %in% colnames(value)) &
+  # allow short-circuit (handles NULL and non-op without going thru merge())
+  if ((all(horizonNames(object) %in% colnames(value)) ||
+       all(colnames(value) %in% horizonNames(object))) &&
+       all(c(idn, hdn, hzd) %in% colnames(value)) &&
       nrow(value) == nrow(object)) {
     if (!all(value[[idn]] %in% profile_id(object))) {
       message("Some profile IDs in input data are not present in object and no new columns to merge. Doing nothing.")
diff --git a/R/estimatePSCS.R b/R/estimatePSCS.R
index 3ac2860f..f29f6c69 100644
--- a/R/estimatePSCS.R
+++ b/R/estimatePSCS.R
@@ -96,10 +96,13 @@ estimatePSCS = function(p, hzdesgn = "hzname", clay.attr = "clay",
   default_t <- rep(25, length(soildepth))
   default_b <- rep(100, length(soildepth))
   
+  odepth <- getMineralSoilSurfaceDepth(p, hzdesgn, simplify = FALSE)[[hz.depths[2]]]
+  ohzidx <- which(odepth > 0)
+  
   # Key part A (soils with restriction in shallow depth)
-  lt36idx <- which(soildepth <= 36)
+  lt36idx <- which(soildepth - odepth <= 36)
   if (length(lt36idx) > 0) {
-    default_t[lt36idx] <- 0
+    default_t[lt36idx] <- 0 # O horizon correction applied below
     default_b[lt36idx] <- soildepth[lt36idx]
     shallow_flag[lt36idx] <- TRUE
   }
@@ -115,8 +118,6 @@ estimatePSCS = function(p, hzdesgn = "hzname", clay.attr = "clay",
   }
 
   # Adjust PSCS range downward if organic soil material is present at surface (i.e. mineral soil surface depth > 0)
-  odepth <- getMineralSoilSurfaceDepth(p, hzdesgn, simplify = FALSE)[[hz.depths[2]]]
-  ohzidx <- which(odepth > 0)
   if (length(ohzidx) > 0) {
     
     default_t[ohzidx] <- default_t[ohzidx] + odepth[ohzidx]
diff --git a/tests/testthat/test-estimatePSCS.R b/tests/testthat/test-estimatePSCS.R
index 8eda0752..498d889d 100644
--- a/tests/testthat/test-estimatePSCS.R
+++ b/tests/testthat/test-estimatePSCS.R
@@ -6,10 +6,22 @@ data(sp1, package = 'aqp')
 depths(sp1) <- id ~ top + bottom
 site(sp1) <- ~ group
 
-p <- sp1[1]
+p <- sp1[1,]
 attr <- 'prop' # clay contents %
 
-q <- sp1[2]
+q <- sp1[2,]
+
+x <- data.frame(
+  peiid = 706300,
+  taxsubgrp = "Lithic Humicryods",
+  top = c(0, 13, 16, 18, 24, 40),
+  bottom = c(13, 16, 18, 24, 40, 65),
+  name = c("Oi", "A", "E", "Bhs", "2C", "2R"),
+  texture = c("SPM", "SIL", "SIL", "SIL", "SIL", "BR"),
+  prop = c(0, 6, 6, 6, 6, 6)
+)
+depths(x) <- peiid ~ top + bottom
+site(x) <- ~ taxsubgrp
 
 test_that("estimatePSCS()", {
 
@@ -64,6 +76,14 @@ test_that("estimatePSCS()", {
   expect_error(estimatePSCS(q2, clay.attr = 'prop', texcl.attr = "texture", hzdesgn = 'foo'))
 })
 
+test_that("estimatePSCS() thin soil profile with O horizon", {
+  expect_equal(estimatePSCS(x, clay.attr = 'prop', texcl.attr = "foo", hzdesgn = 'name'), c(13, 40))
+  expect_equal(estimatePSCS(c(q,x), clay.attr = 'prop', texcl.attr = "foo", hzdesgn = 'name'),
+               data.frame(id = c("706300", "P002"), 
+                          pscs_top = c(13, 30),
+                          pscs_bottom = c(40, 59)))
+})
+
 test_that("estimatePSCS() multiple profiles",{
   e <- estimatePSCS(sp1, clay.attr = 'prop', texcl.attr = "texture", hzdesgn = 'name')
   expect_equal(e$pscs_top, c(49, 30, 2, 32, 5, 31, 25, 27, 28))