WIP

test/vitest/examples/example-01.spec.ts

@@ -11,12 +11,17 @@ import {
 /**
  * Example 1: A and B to delta
  *
+ * Given two positions A and B. Find the exact vector from A to B in meters
+ * north, east and down, and find the direction (azimuth/bearing) to B, relative
+ * to north. Use WGS-84 ellipsoid.
+ *
  * @see https://www.ffi.no/en/research/n-vector/#example_1
  */
 test("Example 1", () => {
+  // PROBLEM:
+
   // Given two positions, A and B as latitudes, longitudes and depths (relative
   // to Earth, E):
-
   const aLat = 1,
     aLon = 2,
     aDepth = 3;
@@ -38,36 +43,45 @@ test("Example 1", () => {
 
   // SOLUTION:
 
-  // Step 1: First, the given latitudes and longitudes are converted to
-  //         n-vectors:
+  // Step 1
+  //
+  // First, the given latitudes and longitudes are converted to n-vectors:
   const a = fromGeodeticCoordinates(radians(aLat), radians(aLon));
   const b = fromGeodeticCoordinates(radians(bLat), radians(bLon));
 
-  // Step 2: When the positions are given as n-vectors (and depths), it is easy
-  //         to find the delta vector decomposed in E. No ellipsoid is specified
-  //         when calling the function, thus WGS-84 (default) is used:
+  // Step 2
+  //
+  // When the positions are given as n-vectors (and depths), it is easy to find
+  // the delta vector decomposed in E. No ellipsoid is specified when calling
+  // the function, thus WGS-84 (default) is used:
   const abE = delta(a, b, aDepth, bDepth);
 
-  // Step 3: We now have the delta vector from A to B, but the three coordinates
-  //         of the vector are along the Earth coordinate frame E, while we need
-  //         the coordinates to be north, east and down. To get this, we define
-  //         a North-East-Down coordinate frame called N, and then we need the
-  //         rotation matrix (direction cosine matrix) rEN to go between E and
-  //         N. We have a simple function that calculates rEN from an n-vector,
-  //         and we use this function (using the n-vector at position A):
+  // Step 3
+  //
+  // We now have the delta vector from A to B, but the three coordinates of the
+  // vector are along the Earth coordinate frame E, while we need the
+  // coordinates to be north, east and down. To get this, we define a
+  // North-East-Down coordinate frame called N, and then we need the rotation
+  // matrix (direction cosine matrix) rEN to go between E and N. We have a
+  // simple function that calculates rEN from an n-vector, and we use this
+  // function (using the n-vector at position A):
   const rEN = toRotationMatrix(a);
 
-  // Step 4: Now the delta vector is easily decomposed in N. Since the vector is
-  //         decomposed in E, we must use rNE (rNE is the transpose of rEN):
+  // Step 4
+  //
+  // Now the delta vector is easily decomposed in N. Since the vector is
+  // decomposed in E, we must use rNE (rNE is the transpose of rEN):
   const abN = transform(transpose(rEN), abE);
 
-  // Step 5: The three components of abN are the north, east and down
-  //         displacements from A to B in meters. The azimuth is simply found
-  //         from element 1 and 2 of the vector (the north and east components):
+  // Step 5
+  //
+  // The three components of abN are the north, east and down displacements from
+  // A to B in meters. The azimuth is simply found from element 1 and 2 of the
+  // vector (the north and east components):
   const azimuth = Math.atan2(abN[1], abN[0]);
 
-  expect(abN[0]).toBeCloseTo(331730.2347808944, 8); // meters
-  expect(abN[1]).toBeCloseTo(332997.8749892695, 8); // meters
-  expect(abN[2]).toBeCloseTo(17404.27136193635, 8); // meters
+  expect(abN[0]).toBeCloseTo(331730.2347808944, 8);
+  expect(abN[1]).toBeCloseTo(332997.8749892695, 8);
+  expect(abN[2]).toBeCloseTo(17404.27136193635, 8);
   expect(azimuth).toBeCloseTo(radians(45.10926323826139), 15);
 });

test/vitest/examples/example-02.spec.ts

@@ -22,47 +22,60 @@ import {
  * @see https://www.ffi.no/en/research/n-vector/#example_2
  */
 test("Example 2", () => {
-  // delta vector from B to C, decomposed in B is given:
-  const bc: Vector = [3000, 2000, 100];
+  // PROBLEM:
 
-  // Position and orientation of B is given:
-  // unit to get unit length of vector
+  // A radar or sonar attached to a vehicle B (Body coordinate frame) measures
+  // the distance and direction to an object C. We assume that the distance and
+  // two angles measured by the sensor (typically bearing and elevation relative
+  // to B) are already converted (by converting from spherical to Cartesian
+  // coordinates) to the vector bcB (i.e. the vector from B to C, decomposed in
+  // B):
+  const bcB: Vector = [3000, 2000, 100];
+
+  // The position of B is given as an n-vector and a depth:
   const b = normalize([1, 2, 3]);
   const bDepth = -400;
-  // the three angles are yaw, pitch, and roll
-  const r = eulerZYXToRotationMatrix(radians(10), radians(20), radians(30));
 
-  // A custom reference ellipsoid is given (replacing WGS-84):
-  // (WGS-72)
+  // The orientation (attitude) of B is given as rNB, specified as yaw, pitch,
+  // roll:
+  const rNB = eulerZYXToRotationMatrix(radians(10), radians(20), radians(30));
+
+  // Use the WGS-72 ellipsoid:
   const e = WGS_72;
 
-  // Find the position of C.
+  // Find the exact position of object C as an n-vector and a depth.
 
   // SOLUTION:
 
-  // Step1: Find R_EN:
-  const rn = toRotationMatrix(b);
+  // Step 1
+  //
+  // The delta vector is given in B. It should be decomposed in E before using
+  // it, and thus we need rEB. This matrix is found from the matrices rEN and
+  // rNB, and we need to find rEN, as in Example 1:
+  const rEN = toRotationMatrix(b);
 
-  // Step2: Find R_EB, from R_EN and R_NB:
-  // Note: closest frames cancel
-  const rb = multiply(rn, r);
+  // Step 2
+  //
+  // Now, we can find rEB y using that the closest frames cancel when
+  // multiplying two rotation matrices (i.e. N is cancelled here):
+  const rEB = multiply(rEN, rNB);
 
-  // Step3: Decompose the delta vector in E:
-  // no transpose of R_EB, since the vector is in B
-  const bce = transform(rb, bc);
+  // Step 3
+  //
+  // The delta vector is now decomposed in E:
+  const bcE = transform(rEB, bcB);
 
-  // Step4: Find the position of C, using the functions that goes from one
-  // position and a delta, to a new position:
-  const [c, cDepth] = destination(b, bce, bDepth, e);
+  // Step 4
+  //
+  // It is now easy to find the position of C using destination (with custom
+  // ellipsoid overriding the default WGS-84):
+  const [c, cDepth] = destination(b, bcE, bDepth, e);
 
-  // When displaying the resulting position for humans, it is more convenient
-  // to see lat, long:
+  // Use human-friendly outputs:
   const [lat, lon] = toGeodeticCoordinates(c);
-
-  // Here we also assume that the user wants the output to be height (= -depth):
   const height = -cDepth;
 
   expect(degrees(lat)).toBeCloseTo(53.32637826433107, 13);
   expect(degrees(lon)).toBeCloseTo(63.46812343514746, 13);
-  expect(height).toBeCloseTo(406.0071960700098, 15); // meters
+  expect(height).toBeCloseTo(406.0071960700098, 15);
 });

test/vitest/examples/example-03.spec.ts

@@ -16,21 +16,29 @@ import {
  * @see https://www.ffi.no/en/research/n-vector/#example_3
  */
 test("Example 3", () => {
-  // Position B is given as p_EB_E ("ECEF-vector")
-  const pb: Vector = apply((n) => n * 6371e3, [0.71, -0.72, 0.1]); // m
+  // PROBLEM:
 
-  // Find position B as geodetic latitude, longitude and height
+  // Position B is given as an “ECEF-vector” pb (i.e. a vector from E, the
+  // center of the Earth, to B, decomposed in E):
+  const pb: Vector = apply((n) => n * 6371e3, [0.71, -0.72, 0.1]);
+
+  // Find the geodetic latitude, longitude and height, assuming WGS-84
+  // ellipsoid.
 
   // SOLUTION:
 
-  // Find n-vector from the p-vector:
+  // Step 1
+  //
+  // We have a function that converts ECEF-vectors to n-vectors:
   const [b, bDepth] = fromECEF(pb);
 
-  // Convert to lat, long and height:
+  // Step 2
+  //
+  // Find latitude, longitude and height:
   const [lat, lon] = toGeodeticCoordinates(b);
   const height = -bDepth;
 
   expect(degrees(lat)).toBeCloseTo(5.685075734513181, 14);
   expect(degrees(lon)).toBeCloseTo(-45.40066325579215, 14);
-  expect(height).toBeCloseTo(95772.10761821801, 15); // meters
+  expect(height).toBeCloseTo(95772.10761821801, 15);
 });

test/vitest/examples/example-04.spec.ts

@@ -14,22 +14,24 @@ import {
  * @see https://www.ffi.no/en/research/n-vector/#example_4
  */
 test("Example 4", () => {
-  // Position B is given with lat, long and height:
+  // PROBLEM:
+
+  // Geodetic latitude, longitude and height are given for position B:
   const bLat = 1;
   const bLon = 2;
   const bHeight = 3;
 
-  // Find the vector p_EB_E ("ECEF-vector")
+  // Find the ECEF-vector for this position.
 
   // SOLUTION:
 
-  // Step1: Convert to n-vector:
+  // Step 1: First, the given latitude and longitude are converted to n-vector:
   const b = fromGeodeticCoordinates(radians(bLat), radians(bLon));
 
-  // Step2: Find the ECEF-vector p_EB_E:
+  // Step 2: Convert to an ECEF-vector:
   const pb = toECEF(b, -bHeight);
 
-  expect(pb[0]).toBeCloseTo(6373290.277218279, 8); // meters
-  expect(pb[1]).toBeCloseTo(222560.2006747365, 8); // meters
-  expect(pb[2]).toBeCloseTo(110568.8271817859, 8); // meters
+  expect(pb[0]).toBeCloseTo(6373290.277218279, 8);
+  expect(pb[1]).toBeCloseTo(222560.2006747365, 8);
+  expect(pb[2]).toBeCloseTo(110568.8271817859, 8);
 });

test/vitest/examples/example-05.spec.ts

@@ -24,7 +24,6 @@ test.each`
   ${"or input as lat/long in deg:"} | ${fromGeodeticCoordinates(radians(88), radians(0))} | ${fromGeodeticCoordinates(radians(89), radians(-170))} | ${332.4564441053448} | ${332.4187248568097}
 `(
   "Example 5 ($label)",
-  // Position A and B are given as n_EA_E and n_EB_E:
   ({
     a,
     b,
@@ -36,16 +35,23 @@ test.each`
     expectedGCD: number;
     expectedED: number;
   }) => {
-    // m, mean Earth radius
+    // PROBLEM:
+
+    // Given two positions A and B as n-vectors, find the surface distance (i.e.
+    // great circle distance).The heights of A and B are not relevant (i.e. if
+    // they do not have zero height, we seek the distance between the points
+    // that are at the surface of the Earth, directly above/below A and B). The
+    // Euclidean distance (chord length) should also be found.
+
+    // Use Earth radius r:
     const r = 6371e3;
 
     // SOLUTION:
 
-    // The great circle distance is given by equation (16) in Gade (2010):
-    // Well conditioned for all angles:
+    // Find the great circle distance:
     const gcd = Math.atan2(norm(cross(a, b)), dot(a, b)) * r;
 
-    // The Euclidean distance is given by:
+    // Find the Euclidean distance:
     const ed = norm(apply((b, a) => b - a, b, a)) * r;
 
     expect(gcd / 1000).toBeCloseTo(expectedGCD, 10); // kilometers

test/vitest/examples/example-06.spec.ts

@@ -23,7 +23,6 @@ test.each`
   ${"or input as lat/long in deg:"} | ${fromGeodeticCoordinates(radians(89.9), radians(-150))} | ${fromGeodeticCoordinates(radians(89.9), radians(150))} | ${89.91282199988446} | ${173.4132244463705}
 `(
   "Example 6 ($label)",
-  // Position B is given at time t0 as pt0 and at time t1 as pt1:
   ({
     pt0,
     pt1,
@@ -35,16 +34,19 @@ test.each`
     expectedLat: number;
     expectedLon: number;
   }) => {
-    // The times are given as:
+    // PROBLEM:
+
+    // Given the position of B at time t0 and t1, pt0 and pt1:
     const t0 = 10;
     const t1 = 20;
-    const ti = 16; // time of interpolation
+    const ti = 16;
 
-    // Find the interpolated position at time ti, n_EB_E_ti
+    // Find an interpolated position at time ti, pti. All positions are given as
+    // n-vectors.
 
     // SOLUTION:
 
-    // Using standard interpolation:
+    // Standard interpolation can be used directly with n-vector:
     const pti = normalize(
       apply(
         (pt0, pt1) => pt0 + ((ti - t0) * (pt1 - pt0)) / (t1 - t0),
@@ -53,8 +55,7 @@ test.each`
       ),
     );
 
-    // When displaying the resulting position for humans, it is more convenient
-    // to see lat, long:
+    // Use human-friendly outputs:
     const [lat, lon] = toGeodeticCoordinates(pti);
 
     expect(degrees(lat)).toBeCloseTo(expectedLat, 12);

test/vitest/examples/example-07.spec.ts

@@ -20,7 +20,6 @@ test.each`
   ${"or input as lat/long in deg:"} | ${fromGeodeticCoordinates(radians(90), radians(0))} | ${fromGeodeticCoordinates(radians(60), radians(10))} | ${fromGeodeticCoordinates(radians(50), radians(-20))} | ${[0.3841171702926, -0.04660240548568945, 0.9221074857571395]}
 `(
   "Example 7 ($label)",
-  // Three positions A, B and C are given:
   ({
     a,
     b,
@@ -32,9 +31,16 @@ test.each`
     c: Vector;
     expectedM: Vector;
   }) => {
+    // PROBLEM:
+
+    // Three positions A, B, and C are given as n-vectors.
+    //
+    // Find the mean position, M. Note that the calculation is independent of
+    // the heights/depths of the positions.
+
     // SOLUTION:
 
-    // Find the horizontal mean position, M:
+    // The mean position is simply given by the mean n-vector:
     const m = normalize(apply((a, b, c) => a + b + c, a, b, c));
 
     expect(m[0]).toBeCloseTo(expectedM[0], 16);