geod.py - Add return_back_azimuth: Optional[bool] = None to `inv_in…

…termediate` and `fwd_intermediate` (`.. versionadded:: 3.5.0 return_back_azimuth`) to allow compatibility between the azimuth output of `fwd`/`fwd_intermediate` and to emit warning for the braking change; update `npts` to override the new default.

test_geod.py - update tests for the new default and add some new tests for new functionality
history.rst - document this change

docs/history.rst

@@ -4,6 +4,13 @@ Change Log
 Latest
 ------
 
+3.5.0
+-----
+- ENH: Add `return_back_azimuth: bool` to allow compatibility between the azimuth output of the following functions (issue #1163):
+    `fwd` and `fwd_intermediate`, `inv` and `inv_intermediate`,
+    Note: BREAKING CHANGE for the default value `return_back_azimuth=True` in the functions `fwd_intermediate` and `inv_intermediate`
+    to mach the default value in `fwd` and `inv`
+
 3.4.1
 -----
 - WHL: Add win32 to build_wheels matrix (pull #1169)

pyproj/geod.py

@@ -18,6 +18,7 @@
 ]
 
 import math
+import warnings
 from typing import Any, Dict, List, Optional, Tuple, Union
 
 from pyproj._geod import Geod as _Geod
@@ -505,6 +506,7 @@ def npts(
             out_lons=None,
             out_lats=None,
             out_azis=None,
+            return_back_azimuth=False,
         )
         return list(zip(res.lons, res.lats))
 
@@ -523,9 +525,11 @@ def inv_intermediate(
         out_lons: Optional[Any] = None,
         out_lats: Optional[Any] = None,
         out_azis: Optional[Any] = None,
+        return_back_azimuth: Optional[bool] = None,
     ) -> GeodIntermediateReturn:
         """
         .. versionadded:: 3.1.0
+        .. versionadded:: 3.5.0 return_back_azimuth
 
         Given a single initial point and terminus point,
         and the number of points, returns
@@ -634,12 +638,24 @@ def inv_intermediate(
             az12(s) of the intermediate point(s)
             If None then buffers would be allocated internnaly
             unless requested otherwise by the flags
+        return_back_azimuth: bool, default=True
+            if True, out_azis will store the back azimuth,
+            Otherwise, out_azis will store the forward azimuth.
 
         Returns
         -------
         GeodIntermediateReturn:
             number of points, distance and output arrays (GeodIntermediateReturn docs)
         """
+        if return_back_azimuth is None:
+            return_back_azimuth = True
+            warnings.warn(
+                "Back azimuth is being returned by default to be compatible with fwd()"
+                "This is a breaking change for pyproj 3.5+."
+                "To avoid this warning, set return_back_azimuth=True."
+                "Otherwise, to restore old behaviour, set return_back_azimuth=False."
+                "This warning will be removed in future version."
+            )
         return super()._inv_or_fwd_intermediate(
             lon1=lon1,
             lat1=lat1,
@@ -654,6 +670,7 @@ def inv_intermediate(
             out_lons=out_lons,
             out_lats=out_lats,
             out_azis=out_azis,
+            return_back_azimuth=return_back_azimuth,
         )
 
     def fwd_intermediate(
@@ -670,9 +687,11 @@ def fwd_intermediate(
         out_lons: Optional[Any] = None,
         out_lats: Optional[Any] = None,
         out_azis: Optional[Any] = None,
+        return_back_azimuth: Optional[bool] = None,
     ) -> GeodIntermediateReturn:
         """
         .. versionadded:: 3.1.0
+        .. versionadded:: 3.5.0 return_back_azimuth
 
         Given a single initial point and azimuth, number of points (npts)
         and delimiter distance between two successive points (del_s), returns
@@ -766,12 +785,24 @@ def fwd_intermediate(
             az12(s) of the intermediate point(s)
             If None then buffers would be allocated internnaly
             unless requested otherwise by the flags
+        return_back_azimuth: bool, default=True
+            if True, out_azis will store the back azimuth,
+            Otherwise, out_azis will store the forward azimuth.
 
         Returns
         -------
         GeodIntermediateReturn:
             number of points, distance and output arrays (GeodIntermediateReturn docs)
         """
+        if return_back_azimuth is None:
+            return_back_azimuth = True
+            warnings.warn(
+                "Back azimuth is being returned by default to be compatible with inv()"
+                "This is a breaking change for pyproj 3.5+."
+                "To avoid this warning, set return_back_azimuth=True."
+                "Otherwise, to restore old behaviour, set return_back_azimuth=False."
+                "This warning will be removed in future version."
+            )
         return super()._inv_or_fwd_intermediate(
             lon1=lon1,
             lat1=lat1,
@@ -786,6 +817,7 @@ def fwd_intermediate(
             out_lons=out_lons,
             out_lats=out_lats,
             out_azis=out_azis,
+            return_back_azimuth=return_back_azimuth,
         )
 
     def line_length(self, lons: Any, lats: Any, radians: bool = False) -> float:

test/test_geod.py

@@ -4,7 +4,7 @@
 import pickle
 import shutil
 import tempfile
-from contextlib import contextmanager
+from contextlib import contextmanager, nullcontext
 from itertools import permutations
 
 import numpy as np
@@ -160,32 +160,40 @@ def test_geod_inverse_transform():
             ).transpose()
 
     del_s = dist / (point_count - 1)
-    lons_a = np.empty(point_count)
-    lats_a = np.empty(point_count)
-    azis_a = np.empty(point_count)
-
-    print("test inv_intermediate (by npts) with azi output")
-    res = gg.inv_intermediate(
-        out_lons=lons_a,
-        out_lats=lats_a,
-        out_azis=azis_a,
-        lon1=lon1pt,
-        lat1=lat1pt,
-        lon2=lon2pt,
-        lat2=lat2pt,
-        npts=point_count,
-        initial_idx=0,
-        terminus_idx=0,
-    )
-    assert res.npts == point_count
-    assert_almost_equal(res.del_s, del_s)
-    assert_almost_equal(res.dist, dist)
-    assert_almost_equal(res.lons, lons)
-    assert_almost_equal(res.lats, lats)
-    assert_almost_equal(res.azis[:-1], azis12[1:])
-    assert res.lons is lons_a
-    assert res.lats is lats_a
-    assert res.azis is azis_a
+    for return_back_azimuth in [None, True, False]:
+        print(f"test inv_intermediate (by npts) with azi output {return_back_azimuth=}")
+        lons_a = np.empty(point_count)
+        lats_a = np.empty(point_count)
+        azis_a = np.empty(point_count)
+
+        with pytest.warns(
+            UserWarning
+        ) if return_back_azimuth is None else nullcontext():
+            res = gg.inv_intermediate(
+                out_lons=lons_a,
+                out_lats=lats_a,
+                out_azis=azis_a,
+                lon1=lon1pt,
+                lat1=lat1pt,
+                lon2=lon2pt,
+                lat2=lat2pt,
+                npts=point_count,
+                initial_idx=0,
+                terminus_idx=0,
+                return_back_azimuth=return_back_azimuth,
+            )
+        assert res.npts == point_count
+        assert_almost_equal(res.del_s, del_s)
+        assert_almost_equal(res.dist, dist)
+        assert_almost_equal(res.lons, lons)
+        assert_almost_equal(res.lats, lats)
+        out_azis = res.azis[:-1]
+        if return_back_azimuth in [True, None]:
+            out_azis = reverse_azimuth(out_azis)
+        assert_almost_equal(out_azis, azis12[1:])
+        assert res.lons is lons_a
+        assert res.lats is lats_a
+        assert res.azis is azis_a
 
     for flags in (GeodIntermediateFlag.AZIS_DISCARD, GeodIntermediateFlag.AZIS_KEEP):
         print("test inv_intermediate (by npts) without azi output, no buffers")
@@ -198,6 +206,7 @@ def test_geod_inverse_transform():
             initial_idx=0,
             terminus_idx=0,
             flags=flags,
+            return_back_azimuth=False,
         )
         assert res.npts == point_count
         assert_almost_equal(res.del_s, del_s)
@@ -226,41 +235,53 @@ def test_geod_inverse_transform():
             initial_idx=0,
             terminus_idx=0,
             flags=flags,
+            return_back_azimuth=False,
         )
         assert res.npts == point_count
-        assert_almost_equal(res.del_s, del_s)
-        assert_almost_equal(res.dist, dist)
-        assert_almost_equal(res.lons, lons_a)
-        assert_almost_equal(res.lats, lats_a)
         assert res.lons is lons_b
         assert res.lats is lats_b
         if flags == GeodIntermediateFlag.AZIS_DISCARD:
             assert res.azis is None
         else:
             assert_almost_equal(res.azis, azis_a)
+        assert_almost_equal(res.del_s, del_s)
+        assert_almost_equal(res.dist, dist)
+        assert_almost_equal(res.lons, lons_a)
+        assert_almost_equal(res.lats, lats_a)
 
-    print("test fwd_intermediate")
-    res = gg.fwd_intermediate(
-        out_lons=lons_b,
-        out_lats=lats_b,
-        out_azis=azis_b,
-        lon1=lon1pt,
-        lat1=lat1pt,
-        azi1=true_az12,
-        npts=point_count,
-        del_s=del_s,
-        initial_idx=0,
-        terminus_idx=0,
-    )
-    assert res.npts == point_count
-    assert_almost_equal(res.del_s, del_s)
-    assert_almost_equal(res.dist, dist)
-    assert_almost_equal(res.lons, lons_a)
-    assert_almost_equal(res.lats, lats_a)
-    assert_almost_equal(res.azis, azis_a)
-    assert res.lons is lons_b
-    assert res.lats is lats_b
-    assert res.azis is azis_b
+    for return_back_azimuth in [False, True, None]:
+        print(f"test fwd_intermediate with {return_back_azimuth=}")
+        lons_b = np.empty(point_count)
+        lats_b = np.empty(point_count)
+        azis_b = np.empty(point_count)
+
+        with pytest.warns(
+            UserWarning
+        ) if return_back_azimuth is None else nullcontext():
+            res = gg.fwd_intermediate(
+                out_lons=lons_b,
+                out_lats=lats_b,
+                out_azis=azis_b,
+                lon1=lon1pt,
+                lat1=lat1pt,
+                azi1=true_az12,
+                npts=point_count,
+                del_s=del_s,
+                initial_idx=0,
+                terminus_idx=0,
+                return_back_azimuth=return_back_azimuth,
+            )
+        assert res.npts == point_count
+        assert res.lons is lons_b
+        assert res.lats is lats_b
+        assert res.azis is azis_b
+        assert_almost_equal(res.del_s, del_s)
+        assert_almost_equal(res.dist, dist)
+        assert_almost_equal(res.lons, lons_a)
+        assert_almost_equal(res.lats, lats_a)
+        if return_back_azimuth in [True, None]:
+            azis_b = reverse_azimuth(azis_b)
+        assert_almost_equal(azis_b, azis_a)
 
     print("test inv_intermediate (by del_s)")
     for del_s_fact, flags in (
@@ -280,6 +301,7 @@ def test_geod_inverse_transform():
             initial_idx=0,
             terminus_idx=0,
             flags=flags,
+            return_back_azimuth=False,
         )
         assert res.npts == point_count
         assert_almost_equal(res.del_s, del_s)