ceph/src/boost/libs/geometry/test/formulas/vertex_longitude.cpp

   1 // Boost.Geometry
   2 // Unit Test
   3
   4 // Copyright (c) 2017 Oracle and/or its affiliates.
   5
   6 // Contributed and/or modified by Vissarion Fysikopoulos, on behalf of Oracle
   7
   8 // Use, modification and distribution is subject to the Boost Software License,
   9 // Version 1.0. (See accompanying file LICENSE_1_0.txt or copy at
  10 // http://www.boost.org/LICENSE_1_0.txt)
  11
  12 #include <geometry_test_common.hpp>
  13
  14 #include "test_formula.hpp"
  15 #include "vertex_longitude_cases.hpp"
  16
  17 #include <boost/geometry/formulas/vertex_latitude.hpp>
  18 #include <boost/geometry/formulas/vertex_longitude.hpp>
  19 #include <boost/geometry/formulas/vincenty_inverse.hpp>
  20 #include <boost/geometry/formulas/thomas_inverse.hpp>
  21 #include <boost/geometry/formulas/andoyer_inverse.hpp>
  22
  23 #include <boost/geometry/strategies/strategies.hpp>
  24
  25 #include <boost/geometry/util/math.hpp>
  26
  27 #define BOOST_GEOMETRY_TEST_DEBUG
  28
  29 namespace bg = boost::geometry;
  30
  31 template<typename CT>
  32 CT test_vrt_lon_sph(CT lon1r,
  33                     CT lat1r,
  34                     CT lon2r,
  35                     CT lat2r)
  36 {
  37     CT a1 = bg::formula::spherical_azimuth(lon1r, lat1r, lon2r, lat2r);
  38
  39     typedef bg::model::point<CT, 2,
  40             bg::cs::spherical_equatorial<bg::radian> > point;
  41
  42     bg::model::segment<point> segment(point(lon1r, lat1r),
  43                                       point(lon2r, lat2r));
  44     bg::model::box<point> box;
  45     bg::envelope(segment, box);
  46
  47     CT vertex_lat;
  48     CT lat_sum = lat1r + lat2r;
  49     if (lat_sum > CT(0))
  50     {
  51         vertex_lat = bg::get_as_radian<bg::max_corner, 1>(box);
  52     } else {
  53         vertex_lat = bg::get_as_radian<bg::min_corner, 1>(box);
  54     }
  55
  56     bg::strategy::azimuth::spherical<> azimuth;
  57
  58     return bg::formula::vertex_longitude
  59             <CT, bg::spherical_equatorial_tag>::
  60             apply(lon1r, lat1r,
  61                   lon2r, lat2r,
  62                   vertex_lat,
  63                   a1,
  64                   azimuth);
  65 }
  66
  67 template
  68 <
  69         template <typename, bool, bool, bool, bool, bool> class FormulaPolicy,
  70         typename CT
  71         >
  72 CT test_vrt_lon_geo(CT lon1r,
  73                     CT lat1r,
  74                     CT lon2r,
  75                     CT lat2r)
  76 {
  77     // WGS84
  78     bg::srs::spheroid<CT> spheroid(6378137.0, 6356752.3142451793);
  79
  80     typedef FormulaPolicy<CT, false, true, false, false, false> formula;
  81     CT a1 = formula::apply(lon1r, lat1r, lon2r, lat2r, spheroid).azimuth;
  82
  83     typedef bg::model::point<CT, 2, bg::cs::geographic<bg::radian> > geo_point;
  84
  85     bg::model::segment<geo_point> segment(geo_point(lon1r, lat1r),
  86                                           geo_point(lon2r, lat2r));
  87     bg::model::box<geo_point> box;
  88     bg::envelope(segment, box);
  89
  90     CT vertex_lat;
  91     CT lat_sum = lat1r + lat2r;
  92     if (lat_sum > CT(0))
  93     {
  94         vertex_lat = bg::get_as_radian<bg::max_corner, 1>(box);
  95     } else {
  96         vertex_lat = bg::get_as_radian<bg::min_corner, 1>(box);
  97     }
  98
  99     bg::strategy::azimuth::geographic<> azimuth_geographic;
 100
 101     return bg::formula::vertex_longitude
 102             <CT, bg::geographic_tag>::apply(lon1r, lat1r,
 103                                             lon2r, lat2r,
 104                                             vertex_lat,
 105                                             a1,
 106                                             azimuth_geographic);
 107
 108 }
 109
 110 void test_all(expected_results const& results)
 111 {
 112     double const d2r = bg::math::d2r<double>();
 113
 114     double lon1r = results.p1.lon * d2r;
 115     double lat1r = results.p1.lat * d2r;
 116     double lon2r = results.p2.lon * d2r;
 117     double lat2r = results.p2.lat * d2r;
 118
 119     if(lon1r > lon2r)
 120     {
 121         std::swap(lon1r, lon2r);
 122         std::swap(lat1r, lat2r);
 123     }
 124
 125     double res_an = test_vrt_lon_geo<bg::formula::andoyer_inverse>
 126             (lon1r, lat1r, lon2r, lat2r);
 127     double res_th = test_vrt_lon_geo<bg::formula::thomas_inverse>
 128             (lon1r, lat1r, lon2r, lat2r);
 129     double res_vi = test_vrt_lon_geo<bg::formula::vincenty_inverse>
 130             (lon1r, lat1r, lon2r, lat2r);
 131     double res_sh = test_vrt_lon_sph(lon1r, lat1r, lon2r, lat2r);
 132
 133     bg::math::normalize_longitude<bg::radian, double>(res_an);
 134     bg::math::normalize_longitude<bg::radian, double>(res_th);
 135     bg::math::normalize_longitude<bg::radian, double>(res_vi);
 136     bg::math::normalize_longitude<bg::radian, double>(res_sh);
 137
 138     check_one(res_an, results.andoyer * d2r, res_vi, 0.001);
 139     check_one(res_th, results.thomas * d2r, res_vi, 0.01);//in some tests thomas gives low accuracy
 140     check_one(res_vi, results.vincenty * d2r, res_vi, 0.0000001);
 141     check_one(res_sh, results.spherical * d2r, res_vi, 1);
 142 }
 143
 144 int test_main(int, char*[])
 145 {
 146
 147     for (size_t i = 0; i < expected_size; ++i)
 148     {
 149         test_all(expected[i]);
 150     }
 151
 152     return 0;
 153 }
 154