randortho.R
##
## r a n d o r t h o . R Generate Random Orthogonal Matrix
##
randortho <- function(n, type = c("orthonormal", "unitary")) {
stopifnot(is.numeric(n), length(n) == 1,
floor(n) == ceiling(n), n >= 1)
if (n == 1)
return(matrix(1, 1, 1))
type <- match.arg(type)
if (type == "orthonormal") {
z <- randn(n, n) / sqrt(2.0)
} else {
z <- (randn(n, n) + 1i * randn(n, n)) / sqrt(2.0)
}
# QR decomposition for real or complex matrices
Z <- qr(z)
q <- qr.Q(Z); r <- qr.R(Z)
d <- diag(r)
ph <- d/abs(d)
q %*% diag(ph) %*% q
}
# rortho <- function(n) {
# warning("Function 'rortho' is deprecated; use 'randortho' instead.",
# call. = FALSE, immediate. = TRUE)
# randortho(n, type = "orthonormal")
# }