commit 98882016a1dc897a6016ebbf912c001c13dd5e28
parent 779b5fda418fa4b71c19afef9f1652af6acd4d6e
Author: eamoncaddigan <eamon.caddigan@gmail.com>
Date: Mon, 12 Sep 2016 15:45:22 -0400
More modularity between calculateSLoWC and calculateTCPA
Diffstat:
7 files changed, 84 insertions(+), 67 deletions(-)
diff --git a/NAMESPACE b/NAMESPACE
@@ -6,4 +6,5 @@ export(createTrajectory)
export(identifyLoWC)
export(identifyNMAC)
export(interpolateTrajectory)
+export(is.flattrajectory)
export(is.flighttrajectory)
diff --git a/R/calculateSLoWC.R b/R/calculateSLoWC.R
@@ -17,67 +17,52 @@
calculateSLoWC <- function(trajectory1, trajectory2) {
checkTrajectories(trajectory1, trajectory2)
- # Find the "origin" lon/lat for the encounter. Distances will be represented
- # in feet north/east from this point. Use the centroid of the trajectories.
- lon0 <- mean(c(trajectory1$longitude, trajectory2$longitude))
- lat0 <- mean(c(trajectory1$latitude, trajectory2$latitude))
-
- t1 <- trajectoryToXYZ(trajectory1, c(lon0, lat0))
- t2 <- trajectoryToXYZ(trajectory2, c(lon0, lat0))
-
- # Flat Earth approximation of aircraft position and velocity
- ac1XYZ <- t1$position
- ac2XYZ <- t2$position
+ # DAA Well Clear thresholds
+ DMOD <- 4000 # ft
+ DH_thr <- 450 # ft
+ TauMod_thr <- 35 # s
- # Convert bearing/speed to velocity vector.
- ac1Velocity <- t1$velocity
- ac2Velocity <- t2$velocity
+ if (!is.flattrajectory(trajectory1)) {
+ lon0 <- mean(c(trajectory1$longitude, trajectory2$longitude))
+ lat0 <- mean(c(trajectory1$latitude, trajectory2$latitude))
+
+ trajectory1 <- trajectoryToXYZ(trajectory1, c(lon0, lat0))
+ trajectory2 <- trajectoryToXYZ(trajectory2, c(lon0, lat0))
+ }
# Distance between aircraft
- dXYZ <- ac2XYZ - ac1XYZ
+ dXYZ <- trajectory2$position - trajectory1$position
dXYZ[, 3] <- abs(dXYZ[, 3])
- relativeVelocity <- ac2Velocity - ac1Velocity
+ relativeVelocity <- trajectory2$velocity - trajectory1$velocity
# Calculate the range
R <- sqrt(apply(dXYZ[, 1:2, drop = FALSE]^2, 1, sum))
+ # Time to the closest point of approach (s)
+ tCPA <- calculateTCPA(trajectory1, trajectory2)
+
+ # Horizontal missed distance
+ HMD <- sqrt((dXYZ[, 1] + relativeVelocity[, 1] * tCPA)^2 +
+ (dXYZ[, 2] + relativeVelocity[, 2] * tCPA)^2 )
+
# Note: the code below here is very close to a direct translation of the
# MATLAB script to R (with a few modifications to permit parallelization).
# Additional optimizations are possible.
-
- # DAA Well Clear thresholds
- DMOD <- 4000 # ft
- DH_thr <- 450 # ft
- TauMod_thr <- 35 # s
-
- dX <- dXYZ[, 1]
- dY <- dXYZ[, 2]
- dH <- dXYZ[, 3]
-
- vrX <- relativeVelocity[, 1]
- vrY <- relativeVelocity[, 2]
-
+
# Horizontal size of the Hazard Zone (TauMod_thr boundary)
- Rdot <- (dX * vrX + dY * vrY) / R;
+ Rdot <- apply(dXYZ[, 1:2, drop = FALSE] * relativeVelocity, 1, sum) / R
S <- pmax(DMOD, .5 * (sqrt( (Rdot * TauMod_thr)^2 + 4 * DMOD^2) - Rdot * TauMod_thr))
# Safeguard against x/0
S[R < 1e-4] <- DMOD
-
- # Calculate time to CPA and projected HMD
- tCPA <- -(dX * vrX + dY * vrY) / (vrX^2 + vrY^2)
- # Safegaurd against singularity
- tCPA[(vrX^2 + vrY^2) == 0 | (dX * vrX + dY * vrY) > 0] <- 0
-
- HMD <- sqrt( (dX + vrX * tCPA)^2 + (dY + vrY * tCPA)^2 )
-
+
# Three penetration components (Range, HMD, vertical separation)
- RangePen <- pmin( ( R / S) , 1)
- HMDPen <- pmin( (HMD / DMOD) , 1)
- DHPen <- pmin( (dH / DH_thr) , 1)
-
+ RangePen <- pmin((R / S), 1)
+ HMDPen <- pmin((HMD / DMOD), 1)
+ DHPen <- pmin((dXYZ[, 3] / DH_thr), 1)
+
hpen <- FGnorm(RangePen, HMDPen)
vSLoWC <- 100 * (1 - FGnorm(hpen, DHPen))
-
+
return(vSLoWC)
}
diff --git a/R/calculateTCPA.R b/R/calculateTCPA.R
@@ -15,27 +15,18 @@
calculateTCPA <- function(trajectory1, trajectory2) {
checkTrajectories(trajectory1, trajectory2)
- # Find the "origin" lon/lat for the encounter. Distances will be represented
- # in feet north/east from this point. Use the centroid of the trajectories.
- lon0 <- mean(c(trajectory1$longitude, trajectory2$longitude))
- lat0 <- mean(c(trajectory1$latitude, trajectory2$latitude))
-
- # Flat Earth approximation of aircraft position and velocity
- ac1XYZ <- cbind(lonlatToXY(trajectory1$longitude, trajectory1$latitude,
- lon0, lat0),
- trajectory1$altitude)
- ac2XYZ <- cbind(lonlatToXY(trajectory2$longitude, trajectory2$latitude,
- lon0, lat0),
- trajectory2$altitude)
-
- # Convert bearing/speed to velocity vector.
- ac1Velocity <- bearingToXY(trajectory1$bearing, trajectory1$groundspeed)
- ac2Velocity <- bearingToXY(trajectory2$bearing, trajectory2$groundspeed)
+ if (!is.flattrajectory(trajectory1)) {
+ lon0 <- mean(c(trajectory1$longitude, trajectory2$longitude))
+ lat0 <- mean(c(trajectory1$latitude, trajectory2$latitude))
+
+ trajectory1 <- trajectoryToXYZ(trajectory1, c(lon0, lat0))
+ trajectory2 <- trajectoryToXYZ(trajectory2, c(lon0, lat0))
+ }
# Distance between aircraft
- dXYZ <- ac2XYZ - ac1XYZ
+ dXYZ <- trajectory2$position - trajectory1$position
dXYZ[, 3] <- abs(dXYZ[, 3])
- relativeVelocity <- ac2Velocity - ac1Velocity
+ relativeVelocity <- trajectory2$velocity - trajectory1$velocity
# Calculate time to CPA and projected HMD
tCPADividend <- -apply(dXYZ[, 1:2, drop = FALSE] * relativeVelocity, 1, sum)
diff --git a/R/checkTrajectories.R b/R/checkTrajectories.R
@@ -1,14 +1,32 @@
-#' Check that both arguments are trajectories with the same timecourse.
+#' Check that both arguments are trajectories with the same timecourse.
#'
#' @param trajectory1 An object to test.
#' @param trajectory2 Another object to test.
#' @return TRUE if everything works. Raises an error otherwise.
+#'
+#' @details There are two options for passing this test. Either both objects can
+#' be flighttrajectories, in which case the time stamps must match, or they
+#' can both be flattrajectories, in which case they must have the same number
+#' of samples and origin.
checkTrajectories <- function(trajectory1, trajectory2) {
- if (!is.flighttrajectory(trajectory1) || !is.flighttrajectory(trajectory2)) {
- stop("Both arguments must be instances of flighttrajectory")
- }
- if (!isTRUE(all.equal(trajectory1$timestamp, trajectory2$timestamp))) {
- stop("Trajectories must have matching time stamps")
+ # Two flight trajectories
+ if (is.flighttrajectory(trajectory1)) {
+ if (!is.flighttrajectory(trajectory2)) {
+ stop("'trajectory1' is a flighttrajectory; 'trajectory2' must match this type")
+ }
+ if (!isTRUE(all.equal(trajectory1$timestamp, trajectory2$timestamp))) {
+ stop("flighttrajectories must have matching time stamps")
+ }
+ } else if (is.flattrajectory(trajectory1)) {
+ if (!is.flattrajectory(trajectory2)) {
+ stop("'trajectory1' is a flattrajectory; 'trajectory2' must match this type")
+ }
+ if (nrow(trajectory1$position) != nrow(trajectory2$position)) {
+ stop("flattrajectories must have the same number of samples")
+ }
+ if (!isTRUE(all.equal(trajectory1$origin, trajectory2$origin))) {
+ stop("flattrajectories must have the same origin")
+ }
}
return(TRUE)
}
diff --git a/R/trajectoryToXYZ.R b/R/trajectoryToXYZ.R
@@ -24,3 +24,7 @@ trajectoryToXYZ <- function(trajectory, origin) {
return(flattrajectory)
}
+
+#' Check if an object is a flighttrajectory
+#' @export
+is.flattrajectory <- function(x) inherits(x, "flattrajectory")
diff --git a/man/checkTrajectories.Rd b/man/checkTrajectories.Rd
@@ -17,4 +17,10 @@ TRUE if everything works. Raises an error otherwise.
\description{
Check that both arguments are trajectories with the same timecourse.
}
+\details{
+There are two options for passing this test. Either both objects can
+ be flighttrajectories, in which case the time stamps must match, or they
+ can both be flattrajectories, in which case they must have the same number
+ of samples and origin.
+}
diff --git a/man/is.flattrajectory.Rd b/man/is.flattrajectory.Rd
@@ -0,0 +1,12 @@
+% Generated by roxygen2: do not edit by hand
+% Please edit documentation in R/trajectoryToXYZ.R
+\name{is.flattrajectory}
+\alias{is.flattrajectory}
+\title{Check if an object is a flighttrajectory}
+\usage{
+is.flattrajectory(x)
+}
+\description{
+Check if an object is a flighttrajectory
+}
+
