diff --git a/DESCRIPTION b/DESCRIPTION
index 0c5151978a47f47d02a0e1dfe2ea051c28b6d17a..d2444dab773613047218c7002202b443924f499b 100644
--- a/DESCRIPTION
+++ b/DESCRIPTION
@@ -1,6 +1,6 @@
 Package: ClimProjDiags
 Title: Set of Tools to Compute Various Climate Indices
-Version: 0.3.0
+Version: 0.3.1
 Authors@R: c(
     person("BSC-CNS", role = c("aut", "cph")),
     person("Nuria", "Perez-Zanon", , "nuria.perez@bsc.es", role = c("aut"), comment = c(ORCID = "0000-0001-8568-3071")),
diff --git a/NEWS.md b/NEWS.md
index dd95442f5ba4fc20e1e46608107618ec7440c725..d72e8eba87f3961e8e032cfe8eb563b16944b2b4 100644
--- a/NEWS.md
+++ b/NEWS.md
@@ -1,3 +1,6 @@
+# 0.3.1 (Release date: 2023-03-23)
+- Subset(): Prioritize the dimension names from names(dim(x)) rather than attribute 'dimensions'; If the input data doesn't have dimension names, the output doesn't have either.
+
 # 0.3.0 (Release date: 2023-02-28)
 - SelBox() and ShiftLon() to accept non-numerical data input  
 - SelBox() uses the latitude and longitude dimension name instead of index  
diff --git a/R/Subset.R b/R/Subset.R
index 6c2321519dc55a32aa8226cd136c90bc2b7b3667..d9cb384e2daf32cfb6aca0294f101db81596eff0 100644
--- a/R/Subset.R
+++ b/R/Subset.R
@@ -4,11 +4,11 @@
 #'similar way as done in the function \code{take()} in the package plyr. There
 #'are two main snprovements:\cr\cr First, the input array can have dimension 
 #'names, either in \code{names(dim(x))} or in the attribute 'dimensions'. If 
-#'both exist, the attribute 'dimensions' is prioritized. The dimensions to 
-#'subset along can be specified via the parameter \code{along} either with
-#'integer indices or either by their name.\cr\cr Second, there are additional
-#'ways to adjust which dimensions are dropped in the resulting array: either to
-#'drop all, to drop none, to drop only the ones that have been sliced or to drop
+#'both exist, \code{names(dim(x))} is prioritized. The dimensions to subset
+#'along can be specified via the parameter \code{along} either with integer
+#'indices or either by their name.\cr\cr Second, there are additional ways to
+#'adjust which dimensions are dropped in the resulting array: either to drop
+#'all, to drop none, to drop only the ones that have been sliced or to drop
 #'only the ones that have not been sliced.\cr\cr 
 #'
 #'@param x A named multidimensional array to be sliced. It can have dimension 
@@ -48,19 +48,27 @@ Subset <- function(x, along, indices, drop = FALSE) {
   }
 
   # Take the input array dimension names
-  dim_names <- attr(x, 'dimensions')
-  if (!is.character(dim_names)) {
-    dim_names <- names(dim(x))
-  } else {
-    names(dim(x)) <- dim_names
+  x_has_names_dim <- TRUE
+  dim_names <- names(dim(x))
+  if (is.null(dim_names)) {
+    dim_names <- attr(x, 'dimensions')
+    x_has_names_dim <- FALSE
   }
 
   if (!is.character(dim_names)) {
     if (any(sapply(along, is.character))) {
       stop("The input array 'x' doesn't have labels for the dimensions but the parameter 'along' contains dimension names.")
     }
+  } else {
+    if (!is.null(names(dim(x))) & !is.null(attr(x, 'dimensions'))) {
+      if (any(names(dim(x)) != attr(x, 'dimensions'))) {
+        warning("Found attribute 'dimensions' containing different dimension names from ",
+                "dim(names(x)). Use the latter one only.")
+      }
+    } else if (is.null(names(dim(x))) & !is.null(attr(x, 'dimensions'))) {
+      names(dim(x)) <- dim_names
+    }
   }
-
   # Check along
   if (any(sapply(along, function(x) !is.numeric(x) && !is.character(x))) |
       length(along) == 0) {
@@ -134,6 +142,10 @@ Subset <- function(x, along, indices, drop = FALSE) {
       }
     }
   }
+  # if names(dim(x)) doesn't exist, remove the dimension names
+  if (!x_has_names_dim) {
+    names(dim(subset)) <- NULL
+  }
 
   # Amend the final dimensions and put dimnames and attributes
   metadata <- attributes(x)
@@ -153,12 +165,14 @@ Subset <- function(x, along, indices, drop = FALSE) {
       }
     }
   } else if (is.character(dim_names) & !identical(dim_names, character(0))) {
-    names(metadata[['dim']]) <- dim_names
+    if (x_has_names_dim) {
+      names(metadata[['dim']]) <- dim_names
+    }
     if ('dimensions' %in% names(attributes(x))) {
       metadata[['dimensions']] <- dim_names
     }
   }
-
   attributes(subset) <- metadata
-  subset
+
+  return(subset)
 }
diff --git a/man/Subset.Rd b/man/Subset.Rd
index f50c5da0bb4098fb407a0335904a70147fa47c70..6b73f28bfa0c24ad41e10ad6cd943ba838392e62 100644
--- a/man/Subset.Rd
+++ b/man/Subset.Rd
@@ -31,11 +31,11 @@ This function allows to subset (i.e. slice, take a chunk of) an array, in a
 similar way as done in the function \code{take()} in the package plyr. There
 are two main snprovements:\cr\cr First, the input array can have dimension 
 names, either in \code{names(dim(x))} or in the attribute 'dimensions'. If 
-both exist, the attribute 'dimensions' is prioritized. The dimensions to 
-subset along can be specified via the parameter \code{along} either with
-integer indices or either by their name.\cr\cr Second, there are additional
-ways to adjust which dimensions are dropped in the resulting array: either to
-drop all, to drop none, to drop only the ones that have been sliced or to drop
+both exist, \code{names(dim(x))} is prioritized. The dimensions to subset
+along can be specified via the parameter \code{along} either with integer
+indices or either by their name.\cr\cr Second, there are additional ways to
+adjust which dimensions are dropped in the resulting array: either to drop
+all, to drop none, to drop only the ones that have been sliced or to drop
 only the ones that have not been sliced.\cr\cr
 }
 \examples{
diff --git a/tests/testthat/test-Subset.R b/tests/testthat/test-Subset.R
index e769453cdca7a6188caf80ebbe41e5de2bfef3dd..8144f96262ee9d76f32684565a9f1644a9979ba0 100644
--- a/tests/testthat/test-Subset.R
+++ b/tests/testthat/test-Subset.R
@@ -43,6 +43,11 @@ expect_error(
 Subset(x = dat1, 'dat', 1, drop = 'yes'),
 "Parameter 'drop' must be one of TRUE, FALSE, 'all', 'selected', 'non-selected', 'none'."
 )
+expect_warning(
+Subset(dat4, 'lat', 1),
+"Found attribute 'dimensions' containing different dimension names from dim(names(x)). Use the latter one only.",
+fixed = T
+)
 
 })
 
@@ -248,43 +253,29 @@ check.attributes = F
 )
 expect_equal(
 attributes(Subset(dat3, 'dat', 1, drop = FALSE)),
-list(dim = c(dat = 1, lat = 2, lon = 10), dimensions = c('dat', 'lat', 'lon'))
+list(dim = c(1, 2, 10), dimensions = c('dat', 'lat', 'lon'))
 )
 
-# drop: lat, lon
+# drop: dat, lat
 expect_equal(
-Subset(dat1, c('lat', 'lon'), list(1, 2), drop = TRUE),
-Subset(dat3, c('lat', 'lon'), list(1, 2), drop = TRUE),
-check.attributes = F
+attributes(Subset(dat3, c('dat', 'lat'), list(1, 1), drop = FALSE)),
+list(dim = c(1, 1, 10), dimensions = c('dat', 'lat', 'lon'))
 )
 expect_equal(
-attributes(Subset(dat3, c('lat', 'lon'), list(1, 2), drop = TRUE)),
-list(dim = 1, dimensions = c('dat', 'lat', 'lon'))
+attributes(Subset(dat3, c('dat', 'lat'), list(1, 1), drop = TRUE)),
+list(dim = c(10), dimensions = c('lon'))
 )
 
-})
-
-test_that("5. dat4", {
-
-# drop: lat
-expect_equal(
-Subset(dat4, 'latitude', 1, drop = FALSE),
-Subset(dat1, 'lat', 1, drop = FALSE),
-check.attributes = F
-)
-expect_equal(
-attributes(Subset(dat4, 'latitude', 1, drop = FALSE)),
-list(dim = c(dataset = 1, latitude = 1, longitude = 10), dimensions = c('dataset', 'latitude', 'longitude'))
-)
 # drop: lat, lon
 expect_equal(
-Subset(dat4, c('latitude', 'longitude'), list(1, 2), drop = TRUE),
 Subset(dat1, c('lat', 'lon'), list(1, 2), drop = TRUE),
+Subset(dat3, c('lat', 'lon'), list(1, 2), drop = TRUE),
 check.attributes = F
 )
 expect_equal(
-attributes(Subset(dat4, c('latitude', 'longitude'), list(1, 2), drop = TRUE)),
-list(dim = 1, dimensions = c('dataset', 'latitude', 'longitude'))
+attributes(Subset(dat3, c('lat', 'lon'), list(1, 2), drop = TRUE)),
+list(dim = 1, dimensions = c('dat', 'lat', 'lon'))
 )
 
 })
+