Skip to content

Commit f9ea18b

Browse files
authored
Merge pull request #19 from dastrobu/readme-updates-for-release
add example to README on azimuths computation
2 parents fcbe6eb + 808c07b commit f9ea18b

File tree

3 files changed

+48
-33
lines changed

3 files changed

+48
-33
lines changed

README.md

Lines changed: 6 additions & 2 deletions
Original file line numberDiff line numberDiff line change
@@ -19,7 +19,11 @@ Here is an example to compute the distance between two points (the poles in this
1919
import vincenty
2020
let d = try distance((lat: Double.pi / 2,lon: 0), (lat: -Double.pi / 2, lon: 0))
2121

22-
and that's it.
22+
To compute azimuths (also known as initial and final bearings)
23+
24+
let (d, (a, b)) = try solveInverse((lat: Double.pi / 2,lon: 0), (lat: -Double.pi / 2, lon: 0))
25+
26+
where `(a, b)` are the azimuths.
2327

2428
<!-- START doctoc generated TOC please keep comment here to allow auto update -->
2529
<!-- DON'T EDIT THIS SECTION, INSTEAD RE-RUN doctoc TO UPDATE -->
@@ -47,7 +51,7 @@ There are no dependencies on macOS.
4751
```swift
4852
let package = Package(
4953
dependencies: [
50-
.package(url: "https://github.com/dastrobu/vincenty.git", from: "1.0.4"),
54+
.package(url: "https://github.com/dastrobu/vincenty.git", from: "1.1.0"),
5155
]
5256
)
5357
```

Sources/vincenty/vincenty.swift

Lines changed: 19 additions & 20 deletions
Original file line numberDiff line numberDiff line change
@@ -62,12 +62,12 @@ public func distance(_ x: (lat: Double, lon: Double),
6262
/// - Throws:
6363
/// - `VincentyError.notConverged` if the distance computation does not converge within `maxIter` iterations.
6464
public func solveInverse(_ x: (lat: Double, lon: Double),
65-
_ y: (lat: Double, lon: Double),
66-
tol: Double = 1e-12,
67-
maxIter: UInt = 200,
68-
ellipsoid: (a: Double, f: Double) = wgs84
69-
) throws -> (distance:Double,azimuths:(Double,Double)) {
70-
65+
_ y: (lat: Double, lon: Double),
66+
tol: Double = 1e-12,
67+
maxIter: UInt = 200,
68+
ellipsoid: (a: Double, f: Double) = wgs84
69+
) throws -> (distance: Double, azimuths: (Double, Double)) {
70+
7171
assert(tol > 0, "tol '\(tol)' ≤ 0")
7272

7373
// validate lat and lon values
@@ -151,43 +151,42 @@ public func solveInverse(_ x: (lat: Double, lon: Double),
151151
- 1.0 / 6.0 * b * cos_2sigma
152152
* (-3 + 4 * sin_sigma * sin_sigma)
153153
* (-3 + 4 * cos_2sigma * cos_2sigma))))
154-
154+
155155
let distance = B * a * (sigma - delta_sigma)
156-
156+
157157
//Azimuth calculations:
158158
let sinSq_sigma = q * q + p * p
159159
// note special handling of exactly antipodal points where sin²σ = 0 (due to discontinuity
160160
// atan2(0, 0) = 0 but atan2(ε, 0) = π/2 / 90°) - in which case bearing is always meridional,
161161
// due north (or due south!)
162162
// α = azimuths of the geodesic; α2 the direction P₁ P₂ produced
163-
let a1 = abs(sinSq_sigma) < Double.leastNonzeroMagnitude ? 0 : atan2(cos_u_y*sin(lambda), cos_u_x*sin_u_y-sin_u_x*cos_u_y*cos(lambda))
164-
let a2 = abs(sinSq_sigma) < Double.leastNonzeroMagnitude ? Double.pi : atan2(cos_u_x*sin(lambda), -sin_u_x*cos_u_y+cos_u_x*sin_u_y*cos(lambda))
165-
163+
let a1 = abs(sinSq_sigma) < Double.leastNonzeroMagnitude ? 0 : atan2(cos_u_y * sin(lambda), cos_u_x * sin_u_y - sin_u_x * cos_u_y * cos(lambda))
164+
let a2 = abs(sinSq_sigma) < Double.leastNonzeroMagnitude ? Double.pi : atan2(cos_u_x * sin(lambda), -sin_u_x * cos_u_y + cos_u_x * sin_u_y * cos(lambda))
165+
166166
let initialTrueTrack = abs(distance) < Double.leastNonzeroMagnitude ? Double.nan : wrap2pi(a1)
167167
let finalTrueTrack = abs(distance) < Double.leastNonzeroMagnitude ? Double.nan : wrap2pi(a2)
168-
168+
169169
return (distance: distance, azimuths: (initialTrueTrack, finalTrueTrack))
170-
170+
171171
}
172172

173173
/* Source: https://www.movable-type.co.uk/scripts/geodesy/docs/dms.js.html */
174174

175-
private func wrap2pi(_ radians:Double) -> Double
176-
{
175+
private func wrap2pi(_ radians: Double) -> Double {
177176
// avoid rounding due to arithmetic ops if within range
178-
guard radians < 0 || radians >= 2*Double.pi else {
177+
guard radians < 0 || radians >= 2 * Double.pi else {
179178
return radians
180179
}
181-
180+
182181
// bearing wrapping requires a sawtooth wave function with a vertical offset equal to the
183182
// amplitude and a corresponding phase shift; this changes the general sawtooth wave function from
184183
// f(x) = (2ax/p - p/2) % p - a
185184
// to
186185
// f(x) = (2ax/p) % p
187186
// where a = amplitude, p = period, % = modulo; however, Swift '%' is a remainder operator
188187
// not a modulo operator - for modulo, replace 'x%n' with '((x%n)+n)%n'
189-
let x = radians, a = Double.pi, p = 2*Double.pi
190-
191-
return ((2*a*x/p).truncatingRemainder(dividingBy: p)+p).truncatingRemainder(dividingBy: p)
188+
let x = radians, a = Double.pi, p = 2 * Double.pi
189+
190+
return ((2 * a * x / p).truncatingRemainder(dividingBy: p) + p).truncatingRemainder(dividingBy: p)
192191
}
193192

Tests/vincentyTests/vincentyTests.swift

Lines changed: 23 additions & 11 deletions
Original file line numberDiff line numberDiff line change
@@ -13,10 +13,12 @@ private let pi: Double = Double.pi
1313
private func rad<T: BinaryFloatingPoint>(fromDegree d: T) -> T {
1414
return d * T.pi / 180
1515
}
16+
1617
/// - Returns: radians converted to degrees
1718
private func deg<T: BinaryFloatingPoint>(fromRadian r: T) -> T {
1819
return r * 180 / T.pi
1920
}
21+
2022
/// - Returns: meters converted to nautical miles
2123
private func nm<T: BinaryFloatingPoint>(fromMeters m: T) -> T {
2224
return m / 1852.0
@@ -39,7 +41,17 @@ private extension BinaryFloatingPoint {
3941
}
4042

4143
final class VincentyTests: XCTestCase {
42-
44+
45+
/// example code as shown in the Readme.md
46+
func testRreadmeExamples() {
47+
do {
48+
let d = try! vincenty.distance((lat: Double.pi / 2, lon: 0), (lat: -Double.pi / 2, lon: 0))
49+
}
50+
do {
51+
let (d, (a, b)) = try! solveInverse((lat: Double.pi / 2, lon: 0), (lat: -Double.pi / 2, lon: 0))
52+
}
53+
54+
}
4355

4456
func testShortcutForEqualPoints() {
4557
// make sure, points are equal and not identical, to check if the shortcut works correctly
@@ -98,7 +110,7 @@ final class VincentyTests: XCTestCase {
98110
x = (lat: 0.asRad, lon: -0.5.asRad)
99111
y = (lat: 0.asRad, lon: 0.5.asRad)
100112
XCTAssertEqual(try! vincenty.solveInverse(x, y).distance, 111319.491, accuracy: delta)
101-
113+
102114
//Test Cardinals
103115
x = (lat: 0.0, lon: 0.0)
104116
y = (lat: pi/2,lon: 0.0) //north pole
@@ -110,23 +122,23 @@ final class VincentyTests: XCTestCase {
110122
(_, azimuths: (initialTrueTrack, finalTrueTrack)) = try! vincenty.solveInverse(x, y)
111123
XCTAssertEqual(initialTrueTrack, Double.pi/2, accuracy: delta)
112124
XCTAssertEqual(finalTrueTrack, Double.pi/2, accuracy: delta)
113-
125+
114126
y = (lat: -pi/2,lon: 0.0) //south pole
115127
(_, azimuths: (initialTrueTrack, finalTrueTrack)) = try! vincenty.solveInverse(x, y)
116128
XCTAssertEqual(initialTrueTrack, Double.pi, accuracy: delta)
117129
XCTAssertEqual(finalTrueTrack, Double.pi, accuracy: delta)
118-
130+
119131
y = (lat: 0.0,lon: -pi/2) //west
120132
(_, azimuths: (initialTrueTrack, finalTrueTrack)) = try! vincenty.solveInverse(x, y)
121133
XCTAssertEqual(initialTrueTrack, 3*Double.pi/2, accuracy: delta)
122134
XCTAssertEqual(finalTrueTrack, 3*Double.pi/2, accuracy: delta)
123-
135+
124136
}
125-
137+
126138
/// Test against A330 FMS
127139
let fmsAcc = 0.49 //within half nm or degree
128140
func testNavigationAccurracy() {
129-
141+
130142
var x: (lat: Double, lon: Double), y: (lat: Double, lon: Double)
131143

132144
//Urabi to Bumat
@@ -135,8 +147,8 @@ final class VincentyTests: XCTestCase {
135147
var (distance, azimuths: (initialTrueTrack, _)) = try! vincenty.solveInverse(x, y)
136148
XCTAssertEqual(distance.inNm, 197, accuracy: fmsAcc)
137149
XCTAssertEqual(initialTrueTrack.asDegrees, 058, accuracy: fmsAcc)
138-
139-
150+
151+
140152
//Dacey is N5933.6 / W12604.5
141153
x = (lat: (59+33.6/60).asRad, lon: -(126+04.5/60).asRad)
142154
//MCT is N5321.4 / W00215.7
@@ -149,9 +161,9 @@ final class VincentyTests: XCTestCase {
149161
print("vincenty: \(distance), geodesic: \(gDist), delta: \(fabs(distance - gDist))")
150162
XCTAssertEqual(distance, gDist, accuracy: 1e-3)
151163
XCTAssertEqual(initialTrueTrack.asDegrees, 036, accuracy: fmsAcc)
152-
164+
153165
}
154-
166+
155167

156168
/// use geodesic as reference and test vincenty distance.
157169
func testAgainstGeodesic() {

0 commit comments

Comments
 (0)