coordinate transformation correction

R/coordinates.R

@@ -20,20 +20,42 @@ lin2vec0 <- function(azi, inc) {
 
 vec2lin0 <- function(x, y, z) {
   n <- vnorm(cbind(x, y, z)) # normalized vector
-  nz <- sapply(n[, 3], function(x) ifelse(x < 0, -x, x))
-  cbind(
-    azimuth = atan2d(n[, 2], n[, 1]) %% 360,
-    plunge = asind(nz)
-  )
+  # nz <- sapply(n[, 3], function(x) ifelse(x < 0, -x, x))
+  nz <- n[, 3]
+  # cbind(
+  azimuth = atan2d(n[, 2], n[, 1])
+  plunge = asind(nz)
+  # )
+
+  res <- mapply(correct_inc, azi = azimuth, inc = plunge) |> t()
+  rownames(res) <- rownames(x)
+  colnames(res) <- c('azimuth', 'plunge')
+  res
 }
 
 vec2fol0 <- function(x, y, z) {
   n <- vnorm(cbind(x, y, z)) # normalized vector
-  nz <- sapply(n[, 3], function(x) ifelse(x < 0, -x, x))
-  cbind(
-    dip_direction = (atan2d(n[, 2], n[, 1]) + 180) %% 360,
-    dip = 90 - asind(nz)
-  )
+  # nz <- sapply(n[, 3], function(x) ifelse(x < 0, -x, x))
+  nz <- n[, 3]
+
+  dip_direction <- (atan2d(n[, 2], n[, 1]) + 180)
+  dip <- 90 - asind(nz)
+
+  res <- mapply(correct_inc, azi = dip_direction, inc = dip) |> t()
+  rownames(res) <- rownames(x)
+  colnames(res) <- c('dip_direction', 'dip')
+  res
+}
+
+correct_inc <- function(azi, inc) {
+  if (inc > 90) {
+    inc <- 180 - 90
+    azi <- (azi + 180) %% 360
+  } else if (inc < 0) {
+    inc <- abs(inc)
+    azi <- (azi + 180) %% 360
+  }
+  c(azi, inc)
 }
 
 

R/gg_stereonet.R

@@ -25,10 +25,17 @@
 #'   ggplot2::geom_point(data = gg(x), ggplot2::aes(x, y), color = "red") +
 #'   ggplot2::geom_path(data = ggl(x), ggplot2::aes(x, y), color = "red")
 #'
-#' x <- Line(120, 5)
+#' x2 <- Line(120, 5)
 #' ggstereo() +
-#'   ggplot2::geom_point(data = gg(x), ggplot2::aes(x, y), color = "darkgreen") +
-#'   ggplot2::geom_path(data = ggl(x, d = 8), ggplot2::aes(x, y, group = group), color = "darkgreen")
+#'   ggplot2::geom_point(data = gg(x2), ggplot2::aes(x, y), color = "darkgreen") +
+#'   ggplot2::geom_path(data = ggl(x2, d = 8), ggplot2::aes(x, y, group = group), color = "darkgreen")
+#'   
+#' x3 <- Plane(137, 71)
+#' ggstereo() +
+#'   ggplot2::geom_point(data = gg(x3), ggplot2::aes(x, y), color = "darkgreen") +
+#'   ggplot2::geom_path(data = ggl(x3, d = 90), ggplot2::aes(x, y, group = group), color = "darkgreen", lwd = 1) +
+#'   ggplot2::geom_path(data = ggl(x3, d = 90 + 11), ggplot2::aes(x, y, group = group, color = 'sde <90')) +
+#'   ggplot2::geom_path(data = ggl(x3, d = 90 - 11), ggplot2::aes(x, y, group = group, color = 'sde >90'))
 #'   }
 NULL
 
@@ -60,7 +67,7 @@ ggl <- function(x, ..., d = 90, n = 1e3) {
   stopifnot(is.spherical(x))
   if (n %% 2 > 0) n <- n + 1
   if (is.plane(x) | is.fault(x)) {
-    x[, 1] <- 180 + x[, 1]
+    #x[, 1] <- 180 + x[, 1]
     x[, 2] <- 90 - x[, 2]
   }
 
@@ -73,9 +80,9 @@ ggl <- function(x, ..., d = 90, n = 1e3) {
 
   zaxis <- c(0, 0, 1)
 
-  res <- matrix(ncol = 4, nrow = n * nx) |>
+  res <- matrix(ncol = 3, nrow = n * nx) |>
     as.data.frame()
-  colnames(res) <- c("x", "y", "cond", "id")
+  colnames(res) <- c("x", "y",  "id")
 
   for (i in seq_along(x[, 1])) {
     D <- cbind(azimuth = seq(0, 360, l = n), plunge = rep(90 - d[i], n)) |>
@@ -88,6 +95,13 @@ ggl <- function(x, ..., d = 90, n = 1e3) {
     D1 <- vrotate(D, zaxis, deg2rad(strike))
     rotangle <- vangle(zaxis, as.line(x[i, ]))
 
+    if(d[i] < 90 & is.plane(x)){
+      d[i] <- 180 - d[i]
+      D1 <- -D1
+      x[i, 1] <- x[i, 1] + 180
+      x[i, 2] <- 90- x[i, 2]
+    }
+
     if (d[i] < 90) {
       k <- -1
     } else {
@@ -98,35 +112,36 @@ ggl <- function(x, ..., d = 90, n = 1e3) {
     D_rot <- vrotate(D1, rotaxis, k * rotangle) |>
       vec2line()
 
-    D_fixed <- list(az = D_rot[, 1], inc = D_rot[, 2]) |>
-      fix_inc()
+    D_fixed <- list(az = D_rot[, 1], inc = D_rot[, 2]) |> fix_inc()
 
-    if (d[i] < 90 & d[i] > as.line(x[i, ])[2]) {
+    if (d[i] != 90 & d[i] > as.line(x[i, ])[2]) {
       # upper hemisphere
-      flag <- TRUE
-      D_rotrot <- vrotate(D_rot, zaxis, deg2rad(180))
+      #flag <- TRUE
+      #D_rotrot <- vrotate(D_rot, zaxis, deg2rad(180))
 
-      D_rotrot_fixed <- list(az = D_rotrot[, 1], inc = D_rotrot[, 2]) |>
-        fix_inc()
+      # D_rotrot_fixed <- list(az = D_rotrot[, 1], inc = D_rotrot[, 2]) |>
+      #   fix_inc()
+      D_rotrot_fixed <- D_fixed
 
-      prec <- sqrt(.Machine$double.eps)
+      #prec <- sqrt(.Machine$double.eps)
       dangle <- vangle(Line(D_fixed$az, D_fixed$inc), as.line(x[i, ]))
-      cond <- d[i] >= (dangle + prec) | d[i] <= (dangle - prec)
+      #cond <- d[i] >= (dangle + prec) | d[i] <= (dangle - prec)
+      cond <- dplyr::near(d[i], dangle)
 
       D_fixed$az[cond] <- D_rotrot_fixed$az[cond]
       D_fixed$inc[cond] <- D_rotrot_fixed$inc[cond]
-    } else {
-      cond <- rep(FALSE, n)
-    }
+    } #else {
+      #cond <- rep(FALSE, n)
+    #}
 
-    D_fixed2 <- data.frame(x = 180 - D_fixed$az, y = D_fixed$inc)
+    D_fixed2 <- data.frame(x = (180 - D_fixed$az), y = D_fixed$inc)
 
     if (d[i] < 90) {
       D_fixed3 <- D_fixed2
-      D_fixed3$cond <- cond
+      #D_fixed3$cond <- cond
     } else {
       D_fixed3 <- utils::tail(D_fixed2, n = n / 2)
-      D_fixed3$cond <- FALSE
+      #D_fixed3$cond <- FALSE
     }
 
     D_fixed3$id <- i
@@ -137,7 +152,7 @@ ggl <- function(x, ..., d = 90, n = 1e3) {
   }
   res |>
     dplyr::as_tibble() |>
-    dplyr::mutate(group = paste0(as.character(id), as.character(cond))) |>
+    dplyr::mutate(group = as.character(id)) |>
     dplyr::left_join(
       xdf, dplyr::join_by(id)
     )
@@ -239,7 +254,7 @@ ggstereo <- function(data = NULL, mapping = aes(), earea = TRUE, centercross = T
   ggplot(data = data, mapping = mapping) +
     #theme_void() +
     theme(
-      title = element_text(element_text(face = "bold")),
+      title = element_text(face = "bold"),
       panel.background = element_blank(), 
       panel.border = element_blank(),
       panel.grid = element_blank(),
@@ -335,7 +350,7 @@ vmf_kerncontour <- function(u, hw = NULL, kernel_method = c("cross", "rot"), ngr
 #'
 #' @references Garcia Portugues, E. (2013). Exact risk improvement of
 #' bandwidth selectors for kernel density estimation with directional data.
-#' Electronic Journal of Statistics, 7, 1655<U+2013>1685.
+#' Electronic Journal of Statistics, 7, 1655-1685.
 #'
 #' @import ggplot2
 #' @importFrom Directional vmf.kerncontour euclid vmfkde.tune

R/math.R

@@ -770,22 +770,22 @@ ortensor <- function(x, norm = TRUE) {
 
 #' Eigenvalues and Eigenvectors of a Set of Vectors
 #'
-#' Decomposition of Orientation tensor
+#' Decomposition of Orientation Tensor Eigenvectors and Eigenvalues
 #'
 #' @param x numeric. Can be three element vector, three column array, or an
 #' object of class `"line"` or `"plane"`
-#' @param scaled logical. Whether the eigenvectors should be scaled by the
-#' eigenvalues (only effective if `x` is in Cartesian coordinates)
+#' @param scaled logical. Whether the Eigenvectors should be scaled by the
+#' Eigenvalues (only effective if `x` is in Cartesian coordinates).
 #' @returns list containing
 #' \describe{
 #' \item{`values`}{Eigenvalues}
-#' \item{`vectors`}{Eigenvectors in coordinates system of `x`}
+#' \item{`vectors`}{Eigenvectors in coordinate system of `x`}
 #' }
-#' @seealso [ortensor()]
+#' @seealso [ortensor()], [eigen()]
 #' @export
 #' @examples
 #' set.seed(1)
-#' mu <- Line(120, 50)
+#' mu <- rvmf(n = 1) |>  vec2line()
 #' x <- rfb(100, mu = mu, k = 1, A = diag(c(10, 0, 0)))
 #' x_eigen <- or_eigen(x)
 #' x_eigen
@@ -794,18 +794,10 @@ ortensor <- function(x, norm = TRUE) {
 #' stereo_point(mu, lab = "mu", col = 4)
 #' stereo_point(x_eigen$vectors, col = c(1, 2, 3), lab = c("E1", "E2", "E3"))
 or_eigen <- function(x, scaled = FALSE) {
-  x_or <- ortensor(x, norm = TRUE)
-  x_eigen <- eigen(x_or)
+  x_or <- ortensor(x, norm = FALSE)
+  x_eigen <- eigen(x_or, symmetric = TRUE)
   x_eigen$vectors <- t(x_eigen$vectors)
 
-  x_eigen$vectors[1, ] <- -(x_eigen$vectors[1, ])
-
-  # x_svd <- svd(x_or) # Singular Value Decomposition of a Matrix
-  # x_eigen <- list(
-  #   values = x_svd$d,
-  #   vectors = t(x_svd$u)
-  # )
-
   if (scaled) {
     x_eigen$vectors[, 1] * x_eigen$values[1]
     x_eigen$vectors[, 2] * x_eigen$values[2]
@@ -851,28 +843,28 @@ or_eigen <- function(x, scaled = FALSE) {
 #' @returns list
 #' 
 #' @references 
-#' Flinn, Derek. "On the statistical analysis of fabric diagrams." Geological Journal 3.2 (1963): 247-253.
+#' Flinn, Derek.(1963): "On the statistical analysis of fabric diagrams." Geological Journal 3.2: 247-253.
 #' 
-#' Kirschvink, J. (1980). The least-squares line and plane and the analysis of palaeomagnetic data. Geophysical Journal International, 62(3), 699-718.
+#' Kirschvink, J. (1980): The least-squares line and plane and the analysis of palaeomagnetic data. Geophysical Journal International, 62(3), 699-718.
 #'  
-#' Lisle, Richard J. "The use of the orientation tensor for the description and statistical testing of fabrics." Journal of Structural Geology 7.1 (1985): 115-117.
+#' Lisle, Richard J.  (1985): "The use of the orientation tensor for the description and statistical testing of fabrics." Journal of Structural Geology 7.1: 115-117.
 #' 
-#' Lode, Walter (1926) "Versuche über den Einfluß der mittleren Hauptspannung auf das Fließen der Metalle Eisen, Kupfer und Nickel“ 
+#' Lode, Walter (1926): "Versuche über den Einfluß der mittleren Hauptspannung auf das Fließen der Metalle Eisen, Kupfer und Nickel“ 
 #'  (*"Experiments on the influence of the mean principal stress on the flow of the metals iron, copper and nickel"*], Zeitschrift für Physik, vol. 36 (November), pp. 913–939, DOI: 10.1007/BF01400222
 #' 
-#' Mardia, Kantilal Varichand. "Statistics of directional data." Journal of the Royal Statistical Society Series B: Statistical Methodology 37.3 (1975): 349-371.
+#' Mardia, Kantilal Varichand. (1975): "Statistics of directional data." Journal of the Royal Statistical Society Series B: Statistical Methodology 37.3: 349-371.
 #' 
-#' Nadai, A., and Hodge, P. G., Jr. (December 1, 1963). "Theory of Flow and Fracture of Solids, vol. II." ASME. J. Appl. Mech. December 1963; 30(4): 640. https://doi.org/10.1115/1.3636654
+#' Nadai, A., and Hodge, P. G., Jr. (1963): "Theory of Flow and Fracture of Solids, vol. II." ASME. J. Appl. Mech. December 1963; 30(4): 640. https://doi.org/10.1115/1.3636654
 #' 
-#' Ramsay, John G. "Folding and fracturing of rocks." Mc Graw Hill Book Company 568 (1967). 
+#' Ramsay, John G. (1967): "Folding and fracturing of rocks." Mc Graw Hill Book Company 568. 
 #' 
-#' Vollmer, Frederick W. "An application of eigenvalue methods to structural domain analysis." Geological Society of America Bulletin 102.6 (1990): 786-791.
+#' Vollmer, Frederick W. (1990): "An application of eigenvalue methods to structural domain analysis." Geological Society of America Bulletin 102.6: 786-791.
 #' 
-#' Vollmer, Frederick W. "Representing Progressive Fabric Paths on a Triangular Plot Using a Fabric Density Index and Crystal Axes Eigenvector Barycenters." Geological Society of America Abstracts. Vol. 52. 2020.
+#' Vollmer, Frederick W. (2020): "Representing Progressive Fabric Paths on a Triangular Plot Using a Fabric Density Index and Crystal Axes Eigenvector Barycenters." Geological Society of America Abstracts. Vol. 52.
 #' 
-#' Watterson, Juan. Homogeneous deformation of the gneisses of Vesterland, south-west Greenland. No. 78. CA Reitzel, 1968.
+#' Watterson, Juan. (1968): "Homogeneous deformation of the gneisses of Vesterland, south-west Greenland". No. 78. CA Reitzel.
 #' 
-#' Woodcock, N. H. "Specification of fabric shapes using an eigenvalue method." Geological Society of America Bulletin 88.9 (1977): 1231-1236.
+#' Woodcock, N. H.  (1977): "Specification of fabric shapes using an eigenvalue method." Geological Society of America Bulletin 88.9: 1231-1236.
 #' 
 #' @examples
 #' set.seed(1)