diff --git a/R/CST_Subset.R b/R/CST_Subset.R
index 2a5defbf5a9a57969c80a78a045f549556cc898a..35717757c0679bec88b8cdc293b99a24425f4116 100644
--- a/R/CST_Subset.R
+++ b/R/CST_Subset.R
@@ -85,7 +85,12 @@ CST_Subset <- function(x, along, indices, drop = FALSE, var_dim = NULL,
     # Only rename coordinates that have not been dropped
     if (dim_name %in% names(x$dims)) {
       # Subset coordinate by indices
-      x$coords[[dim_name]] <- .subset_with_attrs(x$coords[[dim_name]], index)
+      if (is.null(dim(x$coords[[dim_name]])) | length(dim(x$coords[[dim_name]])) == 1) {
+        x$coords[[dim_name]] <- .subset_with_attrs(x$coords[[dim_name]], index)
+      } else {
+        x$coords[[dim_name]] <- ClimProjDiags::Subset(x$coords[[dim_name]], along = dim_name,
+                                                      indices = index)
+      }
     }
   }
   # Remove dropped coordinates
@@ -114,9 +119,14 @@ CST_Subset <- function(x, along, indices, drop = FALSE, var_dim = NULL,
     }
     if ((dim_name %in% names(x$dims)) &&
 	      (dim_name %in% names(x$attrs$Variable$metadata))) {
+      variable <- x$attrs$Variable$metadata[[dim_name]]
       # Subset coords by indices
-      x$attrs$Variable$metadata[[dim_name]] <-
- 	      .subset_with_attrs(x$attrs$Variable$metadata[[dim_name]], index)
+      if (is.null(dim(variable)) | length(dim(variable)) == 1) {
+        x$attrs$Variable$metadata[[dim_name]] <- .subset_with_attrs(variable, index)
+      } else {
+        x$attrs$Variable$metadata[[dim_name]] <- ClimProjDiags::Subset(variable, along = dim_name,
+                                                                       indices = index)
+      }
     }
   }
   # Remove metadata from variables that were dropped
@@ -136,6 +146,7 @@ CST_Subset <- function(x, along, indices, drop = FALSE, var_dim = NULL,
   return(x)
 }
 
+# Function to subset vectors with attributes
 .subset_with_attrs <- function(x, ...) {
   l <- x[...]
   x.dims <- names(dim(x))
diff --git a/tests/testthat/test-CST_Subset.R b/tests/testthat/test-CST_Subset.R
index 313673ba495427abc62e029b930dc858642763fc..cb567279f715e40de2d3b5116c2a04961fbc94c1 100644
--- a/tests/testthat/test-CST_Subset.R
+++ b/tests/testthat/test-CST_Subset.R
@@ -118,7 +118,7 @@ suppressWarnings(
                 var = c('tas', 'sfcWind'),
                 sdate = c('20170101'),
                 ensemble = indices(1),
-                time = indices(1),
+                time = indices(1:3),
                 lat = indices(1:10),
                 lon = indices(1:10),
                 synonims = list(lat = c('lat', 'latitude'),
@@ -148,7 +148,7 @@ test_that("3. Output checks with Start", {
   # Check dimensions
   expect_equal(
     dim(res8$data), 
-    c(dat = 1, var = 1, sdate = 1, ensemble = 1, time = 1, lat = 10, lon = 2)
+    c(dat = 1, var = 1, sdate = 1, ensemble = 1, time = 3, lat = 10, lon = 2)
   )
   expect_equal(
     dim(res8$data), 
@@ -156,7 +156,7 @@ test_that("3. Output checks with Start", {
   )
   expect_equal(
     dim(res10$data), 
-    c(time = 1, lat = 10, lon = 2)
+    c(time = 3, lat = 10, lon = 2)
   )
   # Check coordinates
   expect_equal(
@@ -220,3 +220,35 @@ test_that("3. Output checks with Start", {
     c(1) 
   )
 })
+
+##############################################
+
+test_that("3. Output checks with Start", {
+  res11 <- CST_Subset(exp_start, along = c("dat", "lon", 'time', 'var'),
+                      indices = list(1, 1:2, 1:2, 1), dat_dim = 'dat', 
+                      var_dim = 'var', drop = 'non-selected')
+  expect_equal(
+    dim(res11$data),
+    c(dat = 1, var = 1, time = 2, lat = 10, lon = 2)
+  )
+  expect_equal(
+    names(res11$coords),
+    names(res11$dims)
+  )
+  expect_equal(
+    dim(res11$attrs$Dates), 
+    c(time = 2)
+  )
+  expect_equal(
+    dim(res11$coords$time), 
+    c(sdate = 1, time = 2)
+  )
+  expect_equal(
+    dim(res11$attrs$source_files), 
+    c(var = 1)
+  )
+  expect_equal(
+    names(res11$attrs$Variable$metadata), 
+    c("time", "lat", "lon", "tas")
+  )
+})