boost/geometry
2
0
Fork 0
mirror of https://github.com/boostorg/geometry.git synced 2026-01-30 20:02:11 +00:00
Code Issues Projects Releases Wiki Activity

[formulas] Fix direct geodesic method by performing normalization

Browse Source 
- Add minus sign for B12 evaluation
This commit is contained in:
Adeel Ahmad
2018-05-29 21:11:53 +05:00
parent 3dd6bce720
commit b8a225e1cf
1 changed files with 13 additions and 10 deletions
Download Patch File Download Diff File Expand all files Collapse all files

23
include/boost/geometry/formulas/karney_direct.hpp
View File

@@ -103,7 +103,6 @@ public:
CT const c1 = 1;
CT const c2 = 2;
CT const a = CT(get_radius<0>(spheroid));
CT const b = CT(get_radius<2>(spheroid));
CT const f = formula::flattening<CT>(spheroid);
CT const one_minus_f = c1 - f;
@@ -121,6 +120,9 @@ public:
math::sin_cos_degrees<CT>(math::round_angle<CT>(lat1), sin_beta1, cos_beta1);
sin_beta1 *= one_minus_f;
math::normalize<CT>(sin_beta1, cos_beta1);
cos_beta1 = std::max(sqrt(std::numeric_limits<CT>::min()), cos_beta1);
// Obtain alpha 0 by solving the spherical triangle.
CT const sin_alpha0
= sin_alpha1 * cos_beta1;
@@ -146,8 +148,12 @@ public:
CT const sin_tau12 = sin(tau12);
CT const cos_tau12 = cos(tau12);
CT const sin_sigma1 = sin_beta1;
CT const cos_sigma1 = cos_beta1 * cos_alpha1;
CT sin_sigma1 = sin_beta1;
CT sin_omega1 = sin_alpha0 * sin_beta1;
CT cos_sigma1, cos_omega1;
cos_sigma1 = cos_omega1 = sin_beta1 != 0 || cos_alpha1 != 0 ? cos_beta1 * cos_alpha1 : 1;
math::normalize<CT>(sin_sigma1, cos_sigma1);
CT const B11 = sin_cos_series(sin_sigma1, cos_sigma1, coeffs_C1);
CT const sin_B11 = sin(B11);
@@ -162,9 +168,9 @@ public:
CT coeffs_C1p[SeriesOrder + 1];
series_expansion::evaluate_coeffs_C1p<CT, SeriesOrder>(epsilon, coeffs_C1p);
CT const B12 = sin_cos_series(sin_tau1 * cos_tau12 + cos_tau1 * sin_tau12,
CT const B12 = - sin_cos_series(sin_tau1 * cos_tau12 + cos_tau1 * sin_tau12,
cos_tau1 * cos_tau12 - sin_tau1 * sin_tau12,
coeffs_C1p); // < 0?
coeffs_C1p);
CT const sigma12 = tau12 - (B12 - B11);
CT const sin_sigma12 = sin(sigma12);
@@ -200,9 +206,6 @@ public:
result.lat2 /= math::d2r<T>();
// Find the longitude at the second point.
CT const sin_omega1 = sin_beta1 * sin_alpha0;
CT const cos_omega1 = cos_beta1 * cos_alpha1;
CT const sin_omega2 = sin_alpha0 * sin_sigma2;
CT const cos_omega2 = cos_sigma2;
@@ -218,11 +221,11 @@ public:
// Compute the size of coefficient array.
size_t const coeffs_C3_size = (SeriesOrder * (SeriesOrder - 1)) / 2;
CT coeffs_C3x[coeffs_C3_size];
series_expansion::evaluate_coeffs_C3<double, SeriesOrder>(n, coeffs_C3x);
series_expansion::evaluate_coeffs_C3x<double, SeriesOrder>(n, coeffs_C3x);
// Evaluate C3 coefficients.
CT coeffs_C3[SeriesOrder];
math::evaluate_coeffs_var2<double, SeriesOrder>(epsilon, coeffs_C3x, coeffs_C3);
series_expansion::evaluate_coeffs_C3<double, SeriesOrder>(epsilon, coeffs_C3, coeffs_C3x);
CT const B31 = sin_cos_series(sin_sigma1, cos_sigma1, coeffs_C3);