Unit tests for QgsDistanceArea::geodesicLine

qgis · Jan 8, 2019 · 9c7be5f · 9c7be5f
1 parent b0c3b9f
commit 9c7be5f
Show file tree

Hide file tree

Showing 2 changed files with 64 additions and 4 deletions.
diff --git a/src/core/qgsdistancearea.cpp b/src/core/qgsdistancearea.cpp
@@ -579,10 +579,14 @@ QList< QList<QgsPointXY> > QgsDistanceArea::geodesicLine( const QgsPointXY &p1,
 
       try
       {
+        QgsPointXY p;
         if ( prevLon < -120 )
-          currentPart << mCoordTransform.transform( QgsPointXY( -180, lat180 ), QgsCoordinateTransform::ReverseTransform );
+          p = mCoordTransform.transform( QgsPointXY( -180, lat180 ), QgsCoordinateTransform::ReverseTransform );
         else
-          currentPart << mCoordTransform.transform( QgsPointXY( 180, lat180 ), QgsCoordinateTransform::ReverseTransform );
+          p = mCoordTransform.transform( QgsPointXY( 180, lat180 ), QgsCoordinateTransform::ReverseTransform );
+
+        if ( std::isfinite( p.x() ) && std::isfinite( p.y() ) )
+          currentPart << p;
       }
       catch ( QgsCsException & )
       {
@@ -593,10 +597,14 @@ QList< QList<QgsPointXY> > QgsDistanceArea::geodesicLine( const QgsPointXY &p1,
       currentPart.clear();
       try
       {
+        QgsPointXY p;
         if ( lon < -120 )
-          currentPart << mCoordTransform.transform( QgsPointXY( -180, lat180 ), QgsCoordinateTransform::ReverseTransform );
+          p = mCoordTransform.transform( QgsPointXY( -180, lat180 ), QgsCoordinateTransform::ReverseTransform );
         else
-          currentPart << mCoordTransform.transform( QgsPointXY( 180, lat180 ), QgsCoordinateTransform::ReverseTransform );
+          p = mCoordTransform.transform( QgsPointXY( 180, lat180 ), QgsCoordinateTransform::ReverseTransform );
+
+        if ( std::isfinite( p.x() ) && std::isfinite( p.y() ) )
+          currentPart << p;
       }
       catch ( QgsCsException & )
       {

diff --git a/tests/src/python/test_qgsdistancearea.py b/tests/src/python/test_qgsdistancearea.py
@@ -674,5 +674,57 @@ def testGeodesicIntersectionAtDateLine(self):
                                                                   QgsPointXY(175.76717768974583578,
                                                                              8.93749416467257873)), 8.836479503936307, 5)
 
+    def testGeodesicLine(self):
+        da = QgsDistanceArea()
+        crs = QgsCoordinateReferenceSystem(4326, QgsCoordinateReferenceSystem.EpsgCrsId)
+        da.setSourceCrs(crs, QgsProject.instance().transformContext())
+        da.setEllipsoid("WGS84")
+        g = QgsGeometry.fromMultiPolylineXY(da.geodesicLine(QgsPointXY(105.4, 66.4), QgsPointXY(208.4, -77.8),
+                                                            1000000, True))
+        self.assertEqual(g.asWkt(5), 'MultiLineString ((105.4 66.4, 114.11119 58.36882, 119.52376 49.95732, 123.30625 41.35664, 126.19835 32.6479, 128.57411 23.87234, 130.647 15.05482, 132.55465 6.21309, 134.39916 -2.63822, 136.27014 -11.48549, 138.26237 -20.31469, 140.4956 -29.10966, 143.14591 -37.84912, 146.5073 -46.50015, 151.13295 -55.00229, 158.2045 -63.2234, 170.60189 -70.81368, 180 -73.89148),(-180 -73.89148, -164.6126 -76.66761, -151.6 -77.8))')
+        g = QgsGeometry.fromMultiPolylineXY(da.geodesicLine(QgsPointXY(105.4, 66.4), QgsPointXY(208.4, -77.8),
+                                                            1000000, False))
+        self.assertEqual(g.asWkt(5), 'MultiLineString ((105.4 66.4, 114.11119 58.36882, 119.52376 49.95732, 123.30625 41.35664, 126.19835 32.6479, 128.57411 23.87234, 130.647 15.05482, 132.55465 6.21309, 134.39916 -2.63822, 136.27014 -11.48549, 138.26237 -20.31469, 140.4956 -29.10966, 143.14591 -37.84912, 146.5073 -46.50015, 151.13295 -55.00229, 158.2045 -63.2234, 170.60189 -70.81368, -164.6126 -76.66761, -151.6 -77.8))')
+        g = QgsGeometry.fromMultiPolylineXY(da.geodesicLine(QgsPointXY(105.4, 66.4), QgsPointXY(208.4, -77.8),
+                                                            100000, True))
+        self.assertEqual(g.asWkt(5), 'MultiLineString ((105.4 66.4, 106.50684 65.62452, 107.54925 64.84137, 108.53251 64.05125, 109.46143 63.25477, 110.34036 62.45245, 111.17326 61.6448, 111.96369 60.83224, 112.7149 60.01516, 113.42984 59.19392, 114.11119 58.36882, 114.76141 57.54016, 115.38271 56.70818, 115.97713 55.87314, 116.54653 55.03522, 117.09261 54.19464, 117.61695 53.35156, 118.12096 52.50615, 118.60597 51.65855, 119.0732 50.8089, 119.52376 49.95732, 119.95868 49.10392, 120.37892 48.24881, 120.78537 47.39209, 121.17884 46.53385, 121.5601 45.67416, 121.92984 44.81311, 122.28874 43.95077, 122.6374 43.08721, 122.97639 42.22247, 123.30625 41.35664, 123.62747 40.48975, 123.94053 39.62186, 124.24585 38.75302, 124.54386 37.88326, 124.83494 37.01265, 125.11945 36.14121, 125.39773 35.26898, 125.67011 34.39599, 125.93688 33.52229, 126.19835 32.6479, 126.45477 31.77286, 126.7064 30.89719, 126.9535 30.02092, 127.19628 29.14407, 127.43498 28.26668, 127.66979 27.38877, 127.90093 26.51035, 128.12857 25.63146, 128.35291 24.75212, 128.57411 23.87234, 128.79234 22.99214, 129.00775 22.11156, 129.22051 21.23059, 129.43076 20.34927, 129.63863 19.46761, 129.84427 18.58563, 130.0478 17.70334, 130.24935 16.82077, 130.44904 15.93792, 130.647 15.05482, 130.84333 14.17148, 131.03815 13.28792, 131.23156 12.40414, 131.42367 11.52017, 131.61458 10.63603, 131.8044 9.75171, 131.99322 8.86725, 132.18114 7.98265, 132.36825 7.09792, 132.55465 6.21309, 132.74043 5.32816, 132.92567 4.44315, 133.11048 3.55808, 133.29493 2.67295, 133.47912 1.78778, 133.66313 0.90258, 133.84706 0.01736, 134.03098 -0.86785, 134.21498 -1.75305, 134.39916 -2.63822, 134.5836 -3.52335, 134.76839 -4.40843, 134.95362 -5.29344, 135.13937 -6.17837, 135.32575 -7.06321, 135.51283 -7.94794, 135.70071 -8.83255, 135.8895 -9.71702, 136.07928 -10.60134, 136.27014 -11.48549, 136.46221 -12.36947, 136.65557 -13.25325, 136.85033 -14.13682, 137.04659 -15.02017, 137.24448 -15.90328, 137.44411 -16.78614, 137.64558 -17.66872, 137.84904 -18.55102, 138.05459 -19.43301, 138.26237 -20.31469, 138.47252 -21.19602, 138.68518 -22.077, 138.90049 -22.9576, 139.1186 -23.83781, 139.33968 -24.71761, 139.56389 -25.59697, 139.7914 -26.47588, 140.0224 -27.35431, 140.25706 -28.23225, 140.4956 -29.10966, 140.73822 -29.98653, 140.98515 -30.86282, 141.2366 -31.73851, 141.49283 -32.61358, 141.75409 -33.488, 142.02065 -34.36173, 142.2928 -35.23474, 142.57084 -36.107, 142.8551 -36.97847, 143.14591 -37.84912, 143.44364 -38.71891, 143.74868 -39.58779, 144.06142 -40.45572, 144.38232 -41.32265, 144.71183 -42.18852, 145.05045 -43.0533, 145.39872 -43.91692, 145.7572 -44.77931, 146.12651 -45.64041, 146.5073 -46.50015, 146.90028 -47.35845, 147.3062 -48.21523, 147.72589 -49.0704, 148.16022 -49.92387, 148.61014 -50.77552, 149.0767 -51.62525, 149.56099 -52.47294, 150.06423 -53.31844, 150.58774 -54.16161, 151.13295 -55.00229, 151.7014 -55.84031, 152.29481 -56.67548, 152.91501 -57.50758, 153.56405 -58.33637, 154.24414 -59.16162, 154.95771 -59.98302, 155.70745 -60.80027, 156.49628 -61.61301, 157.32744 -62.42086, 158.2045 -63.2234, 159.13138 -64.02012, 160.11242 -64.8105, 161.15241 -65.59393, 162.25662 -66.36973, 163.43088 -67.13713, 164.68163 -67.89527, 166.01595 -68.64317, 167.44162 -69.37973, 168.96717 -70.10371, 170.60189 -70.81368, 172.35586 -71.50806, 174.2399 -72.18501, 176.26551 -72.8425, 178.4447 -73.4782, 180 -73.89148),(-180 -73.89148, -179.21026 -74.08953, -176.68721 -74.67356, -173.97467 -75.22708, -171.06257 -75.74654, -167.94325 -76.22808, -164.6126 -76.66761, -161.07142 -77.06086, -157.32679 -77.40349, -153.39333 -77.69129, -151.6 -77.8))')
+        g = QgsGeometry.fromMultiPolylineXY(da.geodesicLine(QgsPointXY(105.4, 66.4), QgsPointXY(208.4, -77.8),
+                                                            100000, False))
+        self.assertEqual(g.asWkt(5), 'MultiLineString ((105.4 66.4, 106.50684 65.62452, 107.54925 64.84137, 108.53251 64.05125, 109.46143 63.25477, 110.34036 62.45245, 111.17326 61.6448, 111.96369 60.83224, 112.7149 60.01516, 113.42984 59.19392, 114.11119 58.36882, 114.76141 57.54016, 115.38271 56.70818, 115.97713 55.87314, 116.54653 55.03522, 117.09261 54.19464, 117.61695 53.35156, 118.12096 52.50615, 118.60597 51.65855, 119.0732 50.8089, 119.52376 49.95732, 119.95868 49.10392, 120.37892 48.24881, 120.78537 47.39209, 121.17884 46.53385, 121.5601 45.67416, 121.92984 44.81311, 122.28874 43.95077, 122.6374 43.08721, 122.97639 42.22247, 123.30625 41.35664, 123.62747 40.48975, 123.94053 39.62186, 124.24585 38.75302, 124.54386 37.88326, 124.83494 37.01265, 125.11945 36.14121, 125.39773 35.26898, 125.67011 34.39599, 125.93688 33.52229, 126.19835 32.6479, 126.45477 31.77286, 126.7064 30.89719, 126.9535 30.02092, 127.19628 29.14407, 127.43498 28.26668, 127.66979 27.38877, 127.90093 26.51035, 128.12857 25.63146, 128.35291 24.75212, 128.57411 23.87234, 128.79234 22.99214, 129.00775 22.11156, 129.22051 21.23059, 129.43076 20.34927, 129.63863 19.46761, 129.84427 18.58563, 130.0478 17.70334, 130.24935 16.82077, 130.44904 15.93792, 130.647 15.05482, 130.84333 14.17148, 131.03815 13.28792, 131.23156 12.40414, 131.42367 11.52017, 131.61458 10.63603, 131.8044 9.75171, 131.99322 8.86725, 132.18114 7.98265, 132.36825 7.09792, 132.55465 6.21309, 132.74043 5.32816, 132.92567 4.44315, 133.11048 3.55808, 133.29493 2.67295, 133.47912 1.78778, 133.66313 0.90258, 133.84706 0.01736, 134.03098 -0.86785, 134.21498 -1.75305, 134.39916 -2.63822, 134.5836 -3.52335, 134.76839 -4.40843, 134.95362 -5.29344, 135.13937 -6.17837, 135.32575 -7.06321, 135.51283 -7.94794, 135.70071 -8.83255, 135.8895 -9.71702, 136.07928 -10.60134, 136.27014 -11.48549, 136.46221 -12.36947, 136.65557 -13.25325, 136.85033 -14.13682, 137.04659 -15.02017, 137.24448 -15.90328, 137.44411 -16.78614, 137.64558 -17.66872, 137.84904 -18.55102, 138.05459 -19.43301, 138.26237 -20.31469, 138.47252 -21.19602, 138.68518 -22.077, 138.90049 -22.9576, 139.1186 -23.83781, 139.33968 -24.71761, 139.56389 -25.59697, 139.7914 -26.47588, 140.0224 -27.35431, 140.25706 -28.23225, 140.4956 -29.10966, 140.73822 -29.98653, 140.98515 -30.86282, 141.2366 -31.73851, 141.49283 -32.61358, 141.75409 -33.488, 142.02065 -34.36173, 142.2928 -35.23474, 142.57084 -36.107, 142.8551 -36.97847, 143.14591 -37.84912, 143.44364 -38.71891, 143.74868 -39.58779, 144.06142 -40.45572, 144.38232 -41.32265, 144.71183 -42.18852, 145.05045 -43.0533, 145.39872 -43.91692, 145.7572 -44.77931, 146.12651 -45.64041, 146.5073 -46.50015, 146.90028 -47.35845, 147.3062 -48.21523, 147.72589 -49.0704, 148.16022 -49.92387, 148.61014 -50.77552, 149.0767 -51.62525, 149.56099 -52.47294, 150.06423 -53.31844, 150.58774 -54.16161, 151.13295 -55.00229, 151.7014 -55.84031, 152.29481 -56.67548, 152.91501 -57.50758, 153.56405 -58.33637, 154.24414 -59.16162, 154.95771 -59.98302, 155.70745 -60.80027, 156.49628 -61.61301, 157.32744 -62.42086, 158.2045 -63.2234, 159.13138 -64.02012, 160.11242 -64.8105, 161.15241 -65.59393, 162.25662 -66.36973, 163.43088 -67.13713, 164.68163 -67.89527, 166.01595 -68.64317, 167.44162 -69.37973, 168.96717 -70.10371, 170.60189 -70.81368, 172.35586 -71.50806, 174.2399 -72.18501, 176.26551 -72.8425, 178.4447 -73.4782, -179.21026 -74.08953, -176.68721 -74.67356, -173.97467 -75.22708, -171.06257 -75.74654, -167.94325 -76.22808, -164.6126 -76.66761, -161.07142 -77.06086, -157.32679 -77.40349, -153.39333 -77.69129, -151.6 -77.8))')
+
+        g = QgsGeometry.fromMultiPolylineXY(da.geodesicLine(QgsPointXY(121.4, -76.4), QgsPointXY(-121.6, 76.8),
+                                                            1000000, True))
+        self.assertEqual(g.asWkt(5), 'MultiLineString ((121.4 -76.4, 144.24671 -70.15471, 155.62067 -62.40689, 162.18341 -54.11007, 166.52394 -45.56426, 169.70696 -36.88448, 172.23571 -28.12423, 174.38053 -19.31324, 176.30526 -10.47118, 178.12301 -1.61326, 179.92499 7.24704, 180 7.61121),(-180 7.61121, -178.20071 16.0965, -176.15151 24.92076, -173.78654 33.70222, -170.88367 42.4158, -167.0472 51.01916, -161.47936 59.42708, -152.33578 67.43298, -134.83075 74.42214, -121.6 76.8))')
+        g = QgsGeometry.fromMultiPolylineXY(da.geodesicLine(QgsPointXY(121.4, -76.4), QgsPointXY(-121.6, 76.8),
+                                                            1000000, False))
+        self.assertEqual(g.asWkt(5), 'MultiLineString ((121.4 -76.4, 144.24671 -70.15471, 155.62067 -62.40689, 162.18341 -54.11007, 166.52394 -45.56426, 169.70696 -36.88448, 172.23571 -28.12423, 174.38053 -19.31324, 176.30526 -10.47118, 178.12301 -1.61326, 179.92499 7.24704, -178.20071 16.0965, -176.15151 24.92076, -173.78654 33.70222, -170.88367 42.4158, -167.0472 51.01916, -161.47936 59.42708, -152.33578 67.43298, -134.83075 74.42214, -121.6 76.8))')
+        g = QgsGeometry.fromMultiPolylineXY(da.geodesicLine(QgsPointXY(121.4, 6.4), QgsPointXY(-121.6, 2.8),
+                                                            1000000, True))
+        self.assertEqual(g.asWkt(5), 'MultiLineString ((121.4 6.4, 130.40033 7.32484, 139.43407 8.06935, 148.49637 8.61432, 157.57946 8.9455, 166.67342 9.05419, 175.76718 8.93749, 180 8.80684),(-180 8.80684, -175.15031 8.59851, -166.08891 8.04617, -157.05629 7.29488, -148.0572 6.36403, -139.09301 5.2773, -130.16168 4.06195, -121.6 2.8))')
+        g = QgsGeometry.fromMultiPolylineXY(da.geodesicLine(QgsPointXY(121.4, 6.4), QgsPointXY(-121.6, 2.8),
+                                                            1000000, False))
+        self.assertEqual(g.asWkt(5),
+                         'MultiLineString ((121.4 6.4, 130.40033 7.32484, 139.43407 8.06935, 148.49637 8.61432, 157.57946 8.9455, 166.67342 9.05419, 175.76718 8.93749, -175.15031 8.59851, -166.08891 8.04617, -157.05629 7.29488, -148.0572 6.36403, -139.09301 5.2773, -130.16168 4.06195, -121.6 2.8))')
+
+        # different ellipsoid, should be respected
+        da.setEllipsoid("Phobos2000")
+        g = QgsGeometry.fromMultiPolylineXY(da.geodesicLine(QgsPointXY(121.4, 6.4), QgsPointXY(-121.6, 2.8),
+                                                            1000000, False))
+        self.assertEqual(g.asWkt(5),
+                         'MultiLineString ((121.4 6.4, 82.62887 -29.40012, -121.6 2.8))')
+
+        da.setEllipsoid("WGS84")
+        # with reprojection
+        da.setSourceCrs(QgsCoordinateReferenceSystem('EPSG:3857'), QgsProject.instance().transformContext())
+        g = QgsGeometry.fromMultiPolylineXY(da.geodesicLine(QgsPointXY(-13536427, 14138932), QgsPointXY(13760912, -13248201),
+                                                            1000000, False))
+        self.assertEqual(g.asWkt(0),
+                         'MultiLineString ((-13536427 14138932, -16514348 11691516, -17948849 9406595, -18744235 7552985, -19255354 6014890, -19622372 4688888, -19909239 3505045, 19925702 2415579, 19712755 1385803, 19513769 388441, 19318507 -600065, 19117459 -1602293, 18899973 -2642347, 18651869 -3748726, 18351356 -4958346, 17960498 -6322823, 17404561 -7918366, 16514601 -9855937, 14851845 -12232940, 13760912 -13248201))')
+        g = QgsGeometry.fromMultiPolylineXY(da.geodesicLine(QgsPointXY(-13536427, 14138932), QgsPointXY(13760912, -13248201),
+                                                            1000000, True))
+        self.assertEqual(g.asWkt(0),
+                         'MultiLineString ((-13536427 14138932, -16514348 11691516, -17948849 9406595, -18744235 7552985, -19255354 6014890, -19622372 4688888, -19909239 3505045),(19925702 2415579, 19712755 1385803, 19513769 388441, 19318507 -600065, 19117459 -1602293, 18899973 -2642347, 18651869 -3748726, 18351356 -4958346, 17960498 -6322823, 17404561 -7918366, 16514601 -9855937, 14851845 -12232940, 13760912 -13248201))')
+
+
 if __name__ == '__main__':
     unittest.main()