Repository: fatiando/verde

Issue:
Projections getting more than two coordinates
**Description of the problem**

When using any subclass of `BaseGridder` to predict values through the `grid`, `scatter` or `profile` methods, we can optionally pass a `projection` in order to get, for example, a regular grid on geographic coordinates even if the gridder works in Cartesian.
These `projections` are expected to be functions that takes two coordinates (e.g. `longitude` and `latitude`) and return the projected coordinates (`easting` and `northing`, following the same example). But if we also pass some `extra_coords`, the `projection` will get the entire set of coordinates rather only the horizontal ones.


Raised on fatiando/harmonica#170
<!--
Please be as detailed as you can when describing an issue. The more information
we have, the easier it will be for us to track this down.
-->



**Full code that generated the error**

<!--
Include any data files or inputs required to run the code. It really helps if
we can run the code on our own machines.
-->

```python
import pyproj
import numpy as np
import verde as vd


# Fetch data
data = vd.datasets.fetch_texas_wind()

# Project coordinates
projection = pyproj.Proj("merc")
easting, northing = projection(data.longitude.values, data.latitude.values)

# Create a new random heights coordinates
height = np.random.random(easting.size)
coordinates = (easting, northing, height)

# Fit the spline gridder
spline = vd.Spline(damping=1e-10, mindist=100e3)
spline.fit(coordinates, data.air_temperature_c)

# Interpolate into a regular grid defined in geographic coordinates by also passing
# a extra coordinate
region = (
    data.longitude.min(),
    data.longitude.max(),
    data.latitude.min(),
    data.latitude.max(),
)
spline.grid(spacing=0.1, region=region, projection=projection, extra_coords=1)

```


**Full error message**

```
Traceback (most recent call last):
  File "projections.py", line 29, in <module>
    spline.grid(spacing=0.1, region=region, projection=projection, extra_coords=1)
  File "/home/santi/git/verde/verde/base/base_classes.py", line 303, in grid
    data = check_data(self.predict(projection(*coordinates)))
  File "/home/santi/.miniconda3/envs/verde/lib/python3.7/site-packages/pyproj/proj.py", line 193, in __call__
    lon, lat = args
ValueError: too many values to unpack (expected 2)
```



**System information**

* Operating system: Manjaro
* Python installation (Anaconda, system, ETS): Anaconda
* Version of Python: 3.7.6
* Version of this package: v1.4.0+7.g644037e
* If using conda, paste the output of `conda list` below:

<details>
<summary>output of conda list</summary>
<pre>

# Name                    Version                   Build  Channel
_libgcc_mutex             0.1                        main    conda-forge
alabaster                 0.7.12                     py_0    conda-forge
appdirs                   1.4.3                      py_1    conda-forge
astroid                   2.3.3                    py37_0    conda-forge
attrs                     19.3.0                     py_0    conda-forge
babel                     2.7.0                      py_0    conda-forge
black                     19.10b0                  py37_0    conda-forge
bokeh                     1.4.0                    py37_0    conda-forge
bzip2                     1.0.8                h516909a_2    conda-forge
ca-certificates           2020.4.5.1           hecc5488_0    conda-forge
cartopy                   0.17.0          py37h423102d_1009    conda-forge
certifi                   2020.4.5.1       py37hc8dfbb8_0    conda-forge
cffi                      1.13.2           py37h8022711_0    conda-forge
chardet                   3.0.4                 py37_1003    conda-forge
click                     7.0                        py_0    conda-forge
cloudpickle               1.2.2                      py_1    conda-forge
coverage                  5.1              py37h8f50634_0    conda-forge
cryptography              2.8              py37h72c5cf5_1    conda-forge
cycler                    0.10.0                     py_2    conda-forge
cytoolz                   0.10.1           py37h516909a_0    conda-forge
dask                      2.9.0                      py_0    conda-forge
dask-core                 2.9.0                      py_0    conda-forge
dbus                      1.13.6               he372182_0    conda-forge
distributed               2.9.0                      py_0    conda-forge
docutils                  0.15.2                   py37_0    conda-forge
entrypoints               0.3                   py37_1000    conda-forge
expat                     2.2.5             he1b5a44_1004    conda-forge
flake8                    3.7.9            py37hc8dfbb8_1    conda-forge
fontconfig                2.13.1            h86ecdb6_1001    conda-forge
freetype                  2.10.0               he983fc9_1    conda-forge
fsspec                    0.6.1                      py_0    conda-forge
geos                      3.8.0                he1b5a44_0    conda-forge
gettext                   0.19.8.1          hc5be6a0_1002    conda-forge
glib                      2.58.3          py37h6f030ca_1002    conda-forge
gst-plugins-base          1.14.5               h0935bb2_0    conda-forge
gstreamer                 1.14.5               h36ae1b5_0    conda-forge
heapdict                  1.0.1                      py_0    conda-forge
icu                       64.2                 he1b5a44_1    conda-forge
idna                      2.8                   py37_1000    conda-forge
imagesize                 1.1.0                      py_0    conda-forge
importlib_metadata        1.2.0                    py37_0    conda-forge
isort                     4.3.21                   py37_0    conda-forge
jinja2                    2.10.3                     py_0    conda-forge
joblib                    0.14.0                     py_0    conda-forge
jpeg                      9c                h14c3975_1001    conda-forge
kiwisolver                1.1.0            py37hc9558a2_0    conda-forge
lazy-object-proxy         1.4.3            py37h516909a_0    conda-forge
ld_impl_linux-64          2.33.1               h53a641e_7    conda-forge
libblas                   3.8.0               14_openblas    conda-forge
libcblas                  3.8.0               14_openblas    conda-forge
libclang                  9.0.0           default_hde54327_4    conda-forge
libffi                    3.2.1             he1b5a44_1006    conda-forge
libgcc-ng                 9.2.0                hdf63c60_0    conda-forge
libgfortran-ng            7.3.0                hdf63c60_2    conda-forge
libiconv                  1.15              h516909a_1005    conda-forge
liblapack                 3.8.0               14_openblas    conda-forge
libllvm8                  8.0.1                hc9558a2_0    conda-forge
libllvm9                  9.0.0                hc9558a2_3    conda-forge
libopenblas               0.3.7                h5ec1e0e_5    conda-forge
libpng                    1.6.37               hed695b0_0    conda-forge
libstdcxx-ng              9.2.0                hdf63c60_0    conda-forge
libtiff                   4.1.0                hc3755c2_1    conda-forge
libuuid                   2.32.1            h14c3975_1000    conda-forge
libxcb                    1.13              h14c3975_1002    conda-forge
libxkbcommon              0.9.1                hebb1f50_0    conda-forge
libxml2                   2.9.10               hee79883_0    conda-forge
llvmlite                  0.30.0           py37h8b12597_1    conda-forge
locket                    0.2.0                      py_2    conda-forge
lz4-c                     1.8.3             he1b5a44_1001    conda-forge
markupsafe                1.1.1            py37h516909a_0    conda-forge
matplotlib                3.1.2                    py37_1    conda-forge
matplotlib-base           3.1.2            py37h250f245_1    conda-forge
mccabe                    0.6.1                      py_1    conda-forge
more-itertools            8.0.2                      py_0    conda-forge
msgpack-python            0.6.2            py37hc9558a2_0    conda-forge
mypy_extensions           0.4.3                    py37_0    conda-forge
ncurses                   6.1               hf484d3e_1002    conda-forge
nose                      1.3.7                 py37_1003    conda-forge
nspr                      4.24                 he1b5a44_0    conda-forge
nss                       3.47                 he751ad9_0    conda-forge
numba                     0.46.0           py37hb3f55d8_1    conda-forge
numpy                     1.18.4           py37h8960a57_0    conda-forge
numpydoc                  1.0.0              pyh9f0ad1d_0    conda-forge
olefile                   0.46                       py_0    conda-forge
openssl                   1.1.1g               h516909a_0    conda-forge
owslib                    0.19.0                     py_1    conda-forge
packaging                 19.2                       py_0    conda-forge
pandas                    1.0.3            py37h0da4684_1    conda-forge
partd                     1.0.0                      py_0    conda-forge
pathspec                  0.6.0                      py_0    conda-forge
pcre                      8.43                 he1b5a44_0    conda-forge
pillow                    6.2.1            py37hd70f55b_1    conda-forge
pip                       20.1.1                     py_1    conda-forge
pluggy                    0.13.0                   py37_0    conda-forge
pooch                     1.1.1                      py_0    conda-forge
proj                      6.2.1                hc80f0dc_0    conda-forge
psutil                    5.6.7            py37h516909a_0    conda-forge
pthread-stubs             0.4               h14c3975_1001    conda-forge
py                        1.8.0                      py_0    conda-forge
pycodestyle               2.5.0                      py_0    conda-forge
pycparser                 2.19                     py37_1    conda-forge
pyepsg                    0.4.0                      py_0    conda-forge
pyflakes                  2.1.1                      py_0    conda-forge
pygments                  2.5.2                      py_0    conda-forge
pykdtree                  1.3.1           py37h03ebfcd_1003    conda-forge
pylint                    2.4.4                    py37_0    conda-forge
pyopenssl                 19.1.0                   py37_0    conda-forge
pyparsing                 2.4.5                      py_0    conda-forge
pyproj                    2.4.2.post1      py37h12732c1_0    conda-forge
pyqt                      5.12.3           py37hcca6a23_1    conda-forge
pyqt5-sip                 4.19.18                  pypi_0    pypi
pyqtwebengine             5.12.1                   pypi_0    pypi
pyshp                     2.1.0                      py_0    conda-forge
pysocks                   1.7.1                    py37_0    conda-forge
pytest                    5.4.2            py37hc8dfbb8_0    conda-forge
pytest-cov                2.9.0              pyh9f0ad1d_0    conda-forge
pytest-mpl                0.11                       py_1    conda-forge
python                    3.7.6           h8356626_5_cpython    conda-forge
python-dateutil           2.8.1                      py_0    conda-forge
python_abi                3.7                     1_cp37m    conda-forge
pytz                      2019.3                     py_0    conda-forge
pyyaml                    5.2              py37h516909a_0    conda-forge
qt                        5.12.5               hd8c4c69_1    conda-forge
readline                  8.0                  hf8c457e_0    conda-forge
regex                     2019.12.9        py37h516909a_0    conda-forge
requests                  2.22.0                   py37_1    conda-forge
scikit-learn              0.23.1           py37h8a51577_0    conda-forge
scipy                     1.4.1            py37ha3d9a3c_3    conda-forge
setuptools                42.0.2                   py37_0    conda-forge
shapely                   1.6.4           py37h5d51c17_1007    conda-forge
six                       1.13.0                   py37_0    conda-forge
snowballstemmer           2.0.0                      py_0    conda-forge
sortedcontainers          2.1.0                      py_0    conda-forge
sphinx                    2.2.1                      py_0    conda-forge
sphinx-gallery            0.7.0                      py_0    conda-forge
sphinx_rtd_theme          0.4.3                      py_0    conda-forge
sphinxcontrib-applehelp   1.0.1                      py_0    conda-forge
sphinxcontrib-devhelp     1.0.1                      py_0    conda-forge
sphinxcontrib-htmlhelp    1.0.2                      py_0    conda-forge
sphinxcontrib-jsmath      1.0.1                      py_0    conda-forge
sphinxcontrib-qthelp      1.0.2                      py_0    conda-forge
sphinxcontrib-serializinghtml 1.1.3                      py_0    conda-forge
sqlite                    3.30.1               hcee41ef_0    conda-forge
tblib                     1.4.0                      py_0    conda-forge
threadpoolctl             2.0.0              pyh5ca1d4c_0    conda-forge
tk                        8.6.10               hed695b0_0    conda-forge
toml                      0.10.0                     py_0    conda-forge
toolz                     0.10.0                     py_0    conda-forge
tornado                   6.0.3            py37h516909a_0    conda-forge
typed-ast                 1.4.0            py37h516909a_0    conda-forge
typing_extensions         3.7.4.1                  py37_0    conda-forge
urllib3                   1.25.7                   py37_0    conda-forge
verde                     1.2.0+4.g3f2c362           dev_0    <develop>
wcwidth                   0.1.7                      py_1    conda-forge
wheel                     0.33.6                   py37_0    conda-forge
wrapt                     1.11.2           py37h516909a_0    conda-forge
xarray                    0.15.1                     py_0    conda-forge
xorg-libxau               1.0.9                h14c3975_0    conda-forge
xorg-libxdmcp             1.1.3                h516909a_0    conda-forge
xz                        5.2.4             h14c3975_1001    conda-forge
yaml                      0.2.2                h516909a_1    conda-forge
zict                      1.0.0                      py_0    conda-forge
zipp                      0.6.0                      py_0    conda-forge
zlib                      1.2.11            h516909a_1006    conda-forge
zstd                      1.4.4                h3b9ef0a_1    conda-forge


</pre>
</details>


Discussion:

Patch:
diff --git a/verde/base/base_classes.py b/verde/base/base_classes.py
index c0e944f..40e15fb 100644
--- a/verde/base/base_classes.py
+++ b/verde/base/base_classes.py
@@ -300,7 +300,9 @@ class BaseGridder(BaseEstimator):
         if projection is None:
             data = check_data(self.predict(coordinates))
         else:
-            data = check_data(self.predict(projection(*coordinates)))
+            data = check_data(
+                self.predict(project_coordinates(coordinates, projection))
+            )
         data_names = get_data_names(data, data_names)
         coords = {dims[1]: coordinates[0][0, :], dims[0]: coordinates[1][:, 0]}
         attrs = {"metadata": "Generated by {}".format(repr(self))}
@@ -378,7 +380,9 @@ class BaseGridder(BaseEstimator):
         if projection is None:
             data = check_data(self.predict(coordinates))
         else:
-            data = check_data(self.predict(projection(*coordinates)))
+            data = check_data(
+                self.predict(project_coordinates(coordinates, projection))
+            )
         data_names = get_data_names(data, data_names)
         columns = [(dims[0], coordinates[1]), (dims[1], coordinates[0])]
         columns.extend(zip(data_names, data))
@@ -480,14 +484,14 @@ class BaseGridder(BaseEstimator):
         # geographic coordinates since we don't do actual distances on a
         # sphere).
         if projection is not None:
-            point1 = projection(*point1)
-            point2 = projection(*point2)
+            point1 = project_coordinates(point1, projection)
+            point2 = project_coordinates(point2, projection)
         coordinates, distances = profile_coordinates(point1, point2, size, **kwargs)
         data = check_data(self.predict(coordinates))
         # Project the coordinates back to have geographic coordinates for the
         # profile but Cartesian distances.
         if projection is not None:
-            coordinates = projection(*coordinates, inverse=True)
+            coordinates = project_coordinates(coordinates, projection, inverse=True)
         data_names = get_data_names(data, data_names)
         columns = [
             (dims[0], coordinates[1]),
@@ -506,6 +510,42 @@ class BaseGridder(BaseEstimator):
         return self.dims
 
 
+def project_coordinates(coordinates, projection, **kwargs):
+    """
+    Apply projection to given coordiantes
+
+    Allows to apply projections to any number of coordinates, assuming
+    that the first ones are ``longitude`` and ``latitude``.
+
+    Examples
+    --------
+
+    >>> # Define a custom projection function
+    >>> def projection(lon, lat, inverse=False):
+    ...     "Simple projection of geographic coordinates"
+    ...     radius = 1000
+    ...     if inverse:
+    ...         return (lon / radius, lat / radius)
+    ...     return (lon * radius, lat * radius)
+
+    >>> # Apply the projection to a set of coordinates containing:
+    >>> # longitude, latitude and height
+    >>> coordinates = (1., -2., 3.)
+    >>> project_coordinates(coordinates, projection)
+    (1000.0, -2000.0, 3.0)
+
+    >>> # Apply the inverse projection
+    >>> coordinates = (-500.0, 1500.0, -19.0)
+    >>> project_coordinates(coordinates, projection, inverse=True)
+    (-0.5, 1.5, -19.0)
+
+    """
+    proj_coordinates = projection(*coordinates[:2], **kwargs)
+    if len(coordinates) > 2:
+        proj_coordinates += tuple(coordinates[2:])
+    return proj_coordinates
+
+
 def get_data_names(data, data_names):
     """
     Get default names for data fields if none are given based on the data.


Test patch:
diff --git a/verde/tests/test_base.py b/verde/tests/test_base.py
index d637bd4..ad8ea59 100644
--- a/verde/tests/test_base.py
+++ b/verde/tests/test_base.py
@@ -201,6 +201,70 @@ def test_basegridder_projection():
     npt.assert_allclose(prof.distance, distance_true)
 
 
+def test_basegridder_projection_multiple_coordinates():
+    "Test BaseGridder when passing in a projection with multiple coordinates"
+
+    # Lets say the projection is doubling the coordinates
+    def proj(lon, lat, inverse=False):
+        "Project from the new coordinates to the original"
+        if inverse:
+            return (lon / 2, lat / 2)
+        return (lon * 2, lat * 2)
+
+    # Values in "geographic" coordinates
+    region = (0, 10, -10, -5)
+    shape = (51, 31)
+    angular, linear = 2, 100
+    coordinates = scatter_points(region, 1000, random_state=0, extra_coords=(1, 2))
+    data = angular * coordinates[0] + linear
+    # Project before passing to our Cartesian gridder
+    proj_coordinates = proj(coordinates[0], coordinates[1]) + coordinates[2:]
+    grd = PolyGridder().fit(proj_coordinates, data)
+
+    # Check the estimated coefficients
+    # The grid is estimated in projected coordinates (which are twice as large)
+    # so the rate of change (angular) should be half to get the same values.
+    npt.assert_allclose(grd.coefs_, [linear, angular / 2])
+
+    # The actual values for a grid
+    coordinates_true = grid_coordinates(region, shape, extra_coords=(13, 17))
+    data_true = angular * coordinates_true[0] + linear
+
+    # Check the scatter
+    scat = grd.scatter(
+        region, 1000, random_state=0, projection=proj, extra_coords=(13, 17)
+    )
+    npt.assert_allclose(scat.scalars, data)
+    npt.assert_allclose(scat.easting, coordinates[0])
+    npt.assert_allclose(scat.northing, coordinates[1])
+
+    # Check the grid
+    grid = grd.grid(region, shape, projection=proj, extra_coords=(13, 17))
+    npt.assert_allclose(grid.scalars.values, data_true)
+    npt.assert_allclose(grid.easting.values, coordinates_true[0][0, :])
+    npt.assert_allclose(grid.northing.values, coordinates_true[1][:, 0])
+
+    # Check the profile
+    prof = grd.profile(
+        (region[0], region[-1]),
+        (region[1], region[-1]),
+        shape[1],
+        projection=proj,
+        extra_coords=(13, 17),
+    )
+    npt.assert_allclose(prof.scalars, data_true[-1, :])
+    # Coordinates should still be evenly spaced since the projection is a
+    # multiplication.
+    npt.assert_allclose(prof.easting, coordinates_true[0][0, :])
+    npt.assert_allclose(prof.northing, coordinates_true[1][-1, :])
+    # Distance should still be in the projected coordinates. If the projection
+    # is from geographic, we shouldn't be returning distances in degrees but in
+    # projected meters. The distances will be evenly spaced in unprojected
+    # coordinates.
+    distance_true = np.linspace(region[0] * 2, region[1] * 2, shape[1])
+    npt.assert_allclose(prof.distance, distance_true)
+
+
 def test_check_fit_input():
     "Make sure no exceptions are raised for standard cases"
     size = 20