1=head1 NAME
2
3RhumbSolve -- perform rhumb line calculations
4
5=head1 SYNOPSIS
6
7B<RhumbSolve> [ B<-i> | B<-L> I<lat1> I<lon1> I<azi12> ]
8[ B<-e> I<a> I<f> ]
9[ B<-d> | B<-:> ] [ B<-w> ] [ B<-p> I<prec> ] [ B<-s> ]
10[ B<--comment-delimiter> I<commentdelim> ]
11[ B<--version> | B<-h> | B<--help> ]
12[ B<--input-file> I<infile> | B<--input-string> I<instring> ]
13[ B<--line-separator> I<linesep> ]
14[ B<--output-file> I<outfile> ]
15
16=head1 DESCRIPTION
17
18The path with constant heading between two points on the ellipsoid at
19(I<lat1>, I<lon1>) and (I<lat2>, I<lon2>) is called the rhumb line or
20loxodrome.  Its length is I<s12> and the rhumb line has a forward
21azimuth I<azi12> along its length.  Also computed is I<S12> is the area
22between the rhumb line from point 1 to point 2 and the equator; i.e., it
23is the area, measured counter-clockwise, of the geodesic quadrilateral
24with corners (I<lat1>,I<lon1>), (0,I<lon1>), (0,I<lon2>), and
25(I<lat2>,I<lon2>).  A point at a pole is treated as a point a tiny
26distance away from the pole on the given line of longitude.  The
27longitude becomes indeterminate when a rhumb line passes through a pole,
28and B<RhumbSolve> reports NaNs for the longitude and the area in this
29case.
30
31B<NOTE:> the rhumb line is B<not> the shortest path between two points;
32that is the geodesic and it is calculated by GeodSolve(1).
33
34B<RhumbSolve> operates in one of three modes:
35
36=over
37
38=item 1.
39
40By default, B<RhumbSolve> accepts lines on the standard input containing
41I<lat1> I<lon1> I<azi12> I<s12> and prints I<lat2> I<lon2> I<S12> on
42standard output.  This is the direct calculation.
43
44=item 2.
45
46With the B<-i> command line argument, B<RhumbSolve> performs the inverse
47calculation.  It reads lines containing I<lat1> I<lon1> I<lat2> I<lon2>
48and prints the values of I<azi12> I<s12> I<S12> for the corresponding
49shortest rhumb lines.  If the end points are on opposite meridians,
50there are two shortest rhumb lines and the east-going one is chosen.
51
52=item 3.
53
54Command line arguments B<-L> I<lat1> I<lon1> I<azi12> specify a rhumb
55line.  B<RhumbSolve> then accepts a sequence of I<s12> values (one per
56line) on standard input and prints I<lat2> I<lon2> I<S12> for each.
57This generates a sequence of points on a rhumb line.
58
59=back
60
61=head1 OPTIONS
62
63=over
64
65=item B<-i>
66
67perform an inverse calculation (see 2 above).
68
69=item B<-L> I<lat1> I<lon1> I<azi12>
70
71line mode (see 3 above); generate a sequence of points along the rhumb
72line specified by I<lat1> I<lon1> I<azi12>.  The B<-w> flag can be used
73to swap the default order of the 2 geographic coordinates, provided that
74it appears before B<-L>.  (B<-l> is an alternative, deprecated, spelling
75of this flag.)
76
77=item B<-e> I<a> I<f>
78
79specify the ellipsoid via the equatorial radius, I<a> and
80the flattening, I<f>.  Setting I<f> = 0 results in a sphere.  Specify
81I<f> E<lt> 0 for a prolate ellipsoid.  A simple fraction, e.g., 1/297,
82is allowed for I<f>.  By default, the WGS84 ellipsoid is used, I<a> =
836378137 m, I<f> = 1/298.257223563.
84
85=item B<-d>
86
87output angles as degrees, minutes, seconds instead of decimal degrees.
88
89=item B<-:>
90
91like B<-d>, except use : as a separator instead of the d, ', and "
92delimiters.
93
94=item B<-w>
95
96on input and output, longitude precedes latitude (except that on input
97this can be overridden by a hemisphere designator, I<N>, I<S>, I<E>,
98I<W>).
99
100=item B<-p> I<prec>
101
102set the output precision to I<prec> (default 3); I<prec> is the
103precision relative to 1 m.  See L</PRECISION>.
104
105=item B<-s>
106
107By default, the rhumb line calculations are carried out exactly in terms
108of elliptic integrals.  This includes the use of the addition theorem
109for elliptic integrals to compute the divided difference of the
110isometric and rectifying latitudes.  If B<-s> is supplied this divided
111difference is computed using Krueger series for the transverse Mercator
112projection which is only accurate for |I<f>| E<lt> 0.01.  See
113L</ACCURACY>.
114
115=item B<--comment-delimiter> I<commentdelim>
116
117set the comment delimiter to I<commentdelim> (e.g., "#" or "//").  If
118set, the input lines will be scanned for this delimiter and, if found,
119the delimiter and the rest of the line will be removed prior to
120processing and subsequently appended to the output line (separated by a
121space).
122
123=item B<--version>
124
125print version and exit.
126
127=item B<-h>
128
129print usage and exit.
130
131=item B<--help>
132
133print full documentation and exit.
134
135=item B<--input-file> I<infile>
136
137read input from the file I<infile> instead of from standard input; a file
138name of "-" stands for standard input.
139
140=item B<--input-string> I<instring>
141
142read input from the string I<instring> instead of from standard input.
143All occurrences of the line separator character (default is a semicolon)
144in I<instring> are converted to newlines before the reading begins.
145
146=item B<--line-separator> I<linesep>
147
148set the line separator character to I<linesep>.  By default this is a
149semicolon.
150
151=item B<--output-file> I<outfile>
152
153write output to the file I<outfile> instead of to standard output; a
154file name of "-" stands for standard output.
155
156=back
157
158=head1 INPUT
159
160B<RhumbSolve> measures all angles in degrees, all lengths (I<s12>) in
161meters, and all areas (I<S12>) in meters^2.  On input angles (latitude,
162longitude, azimuth, arc length) can be as decimal degrees or degrees,
163minutes, seconds.  For example, C<40d30>, C<40d30'>, C<40:30>, C<40.5d>,
164and C<40.5> are all equivalent.  By default, latitude precedes longitude
165for each point (the B<-w> flag switches this convention); however on
166input either may be given first by appending (or prepending) I<N> or
167I<S> to the latitude and I<E> or I<W> to the longitude.  Azimuths are
168measured clockwise from north; however this may be overridden with I<E>
169or I<W>.
170
171For details on the allowed formats for angles, see the C<GEOGRAPHIC
172COORDINATES> section of GeoConvert(1).
173
174=head1 PRECISION
175
176I<prec> gives precision of the output with I<prec> = 0 giving 1 m
177precision, I<prec> = 3 giving 1 mm precision, etc.  I<prec> is the
178number of digits after the decimal point for lengths.  For decimal
179degrees, the number of digits after the decimal point is I<prec> + 5.
180For DMS (degree, minute, seconds) output, the number of digits after the
181decimal point in the seconds component is I<prec> + 1.  The minimum
182value of I<prec> is 0 and the maximum is 10.
183
184=head1 ERRORS
185
186An illegal line of input will print an error message to standard output
187beginning with C<ERROR:> and causes B<RhumbSolve> to return an exit code
188of 1.  However, an error does not cause B<RhumbSolve> to terminate;
189following lines will be converted.
190
191=head1 ACCURACY
192
193The algorithm used by B<RhumbSolve> uses exact formulas for converting
194between the latitude, rectifying latitude (I<mu>), and isometric
195latitude (I<psi>).  These formulas are accurate for any value of the
196flattening.  The computation of rhumb lines involves the ratio (I<psi1>
197- I<psi2>) / (I<mu1> - I<mu2>) and this is subject to large round-off
198errors if I<lat1> is close to I<lat2>.  So this ratio is computed using
199divided differences using one of two methods: by default, this uses the
200addition theorem for elliptic integrals (accurate for all values of
201I<f>); however, with the B<-s> options, it is computed using the series
202expansions used by TransverseMercatorProj(1) for the conversions between
203rectifying and conformal latitudes (accurate for |I<f>| E<lt> 0.01).
204For the WGS84 ellipsoid, the error is about 10 nanometers using either
205method.
206
207=head1 EXAMPLES
208
209Route from JFK Airport to Singapore Changi Airport:
210
211   echo 40:38:23N 073:46:44W 01:21:33N 103:59:22E |
212   RhumbSolve -i -: -p 0
213
214   103:34:58.2 18523563
215
216N.B. This is B<not> the route typically taken by aircraft because it's
217considerably longer than the geodesic given by GeodSolve(1).
218
219Waypoints on the route at intervals of 2000km:
220
221   for ((i = 0; i <= 20; i += 2)); do echo ${i}000000;done |
222   RhumbSolve -L 40:38:23N 073:46:44W 103:34:58.2 -: -p 0
223
224   40:38:23.0N 073:46:44.0W 0
225   36:24:30.3N 051:28:26.4W 9817078307821
226   32:10:26.8N 030:20:57.3W 18224745682005
227   27:56:13.2N 010:10:54.2W 25358020327741
228   23:41:50.1N 009:12:45.5E 31321269267102
229   19:27:18.7N 027:59:22.1E 36195163180159
230   15:12:40.2N 046:17:01.1E 40041499143669
231   10:57:55.9N 064:12:52.8E 42906570007050
232   06:43:07.3N 081:53:28.8E 44823504180200
233   02:28:16.2N 099:24:54.5E 45813843358737
234   01:46:36.0S 116:52:59.7E 45888525219677
235
236=head1 SEE ALSO
237
238GeoConvert(1), GeodSolve(1), TransverseMercatorProj(1).
239
240An online version of this utility is availbable at
241L<https://geographiclib.sourceforge.io/cgi-bin/RhumbSolve>.
242
243The Wikipedia page, Rhumb line,
244L<https://en.wikipedia.org/wiki/Rhumb_line>.
245
246=head1 AUTHOR
247
248B<RhumbSolve> was written by Charles Karney.
249
250=head1 HISTORY
251
252B<RhumbSolve> was added to GeographicLib,
253L<https://geographiclib.sourceforge.io>, in version 1.37.
254