transform_vector is wrong / not linear for some projections

[csubich]

Description

I'm using (attempting to use) cartopy.crs.transform_vectors to transform map-coordinate wind vectors from one projection to another, currently from a rotated lat-lon projection to an ordinary lat-lon projection. In so doing, I could never quite get the wind speeds to match up after converting from (presumably) degrees/sec back to meters/second.

In drilling down to simpler and simpler reproductions, I've found that for at least cartopy.crs.RotatedPole, transform_vector isn't even linear. I could understand some scale-based difference between the output of transform_vector and a finite difference approximation via transform_points, but I can't think of any reason that the transformation should behave nonlinearly.

Code to reproduce

import cartopy
import numpy as np
print(cartopy.__version__)
latlon_crs = cartopy.crs.PlateCarree()
rot_crs = cartopy.crs.RotatedPole(pole_longitude=180,pole_latitude=45.0)
azero = np.array([0],dtype=np.float32)
aone = np.array([1],dtype=np.float32)

print('[1,0] transformed', np.concatenate(latlon_crs.transform_vectors(rot_crs,azero,azero,aone,azero)))
print('[0,1] transformed', np.concatenate(latlon_crs.transform_vectors(rot_crs,azero,azero,azero,aone)))
print('[1,1] transformed', np.concatenate(latlon_crs.transform_vectors(rot_crs,azero,azero,aone,aone)))

# Via finite difference
def do_fd(x,y):
    fd_eps = 1e-2
    fd_origin = rot_crs.transform_points(latlon_crs,azero,azero)[0,:2]    
    fd_delta = (rot_crs.transform_points(latlon_crs,azero+fd_eps*x,azero+fd_eps*y)[0,:2] - fd_origin)/fd_eps
    return fd_delta

print('---')
    
print('[1,0] via finite difference', do_fd(1,0))
print('[0,1] via finite difference', do_fd(0,1))
print('[1,1] via finite difference', do_fd(1,1))

Output

0.22.0
[1,0] transformed [ 1.00000000e+00 -6.17066907e-06]
[0,1] transformed [-6.1880869e-08  1.0000000e+00]
[1,1] transformed [1.15470764 0.81648649]
---
[1,0] via finite difference [1.41421352e+00 8.72664586e-05]
[0,1] via finite difference [0.         0.99999998]
[1,1] via finite difference [1.41396673 1.00008723]

Full environment definition

Operating system

Ubuntu 22.04.2

Cartopy version

0.22

conda list

# packages in environment at /local/drive1/subichc/miniconda/envs/cartopy_break:
#
# Name                    Version                   Build  Channel
_libgcc_mutex             0.1                 conda_forge    conda-forge
_openmp_mutex             4.5                       2_gnu    conda-forge
brotli                    1.1.0                hd590300_1    conda-forge
brotli-bin                1.1.0                hd590300_1    conda-forge
bzip2                     1.0.8                h7f98852_4    conda-forge
c-ares                    1.20.1               hd590300_0    conda-forge
ca-certificates           2023.7.22            hbcca054_0    conda-forge
cartopy                   0.22.0          py311h320fe9a_0    conda-forge
certifi                   2023.7.22          pyhd8ed1ab_0    conda-forge
contourpy                 1.1.1           py311h9547e67_1    conda-forge
cycler                    0.12.1             pyhd8ed1ab_0    conda-forge
fonttools                 4.43.1          py311h459d7ec_0    conda-forge
freetype                  2.12.1               h267a509_2    conda-forge
geos                      3.12.0               h59595ed_0    conda-forge
keyutils                  1.6.1                h166bdaf_0    conda-forge
kiwisolver                1.4.5           py311h9547e67_1    conda-forge
krb5                      1.21.2               h659d440_0    conda-forge
lcms2                     2.15                 hb7c19ff_3    conda-forge
ld_impl_linux-64          2.40                 h41732ed_0    conda-forge
lerc                      4.0.0                h27087fc_0    conda-forge
libblas                   3.9.0           19_linux64_openblas    conda-forge
libbrotlicommon           1.1.0                hd590300_1    conda-forge
libbrotlidec              1.1.0                hd590300_1    conda-forge
libbrotlienc              1.1.0                hd590300_1    conda-forge
libcblas                  3.9.0           19_linux64_openblas    conda-forge
libcurl                   8.4.0                hca28451_0    conda-forge
libdeflate                1.19                 hd590300_0    conda-forge
libedit                   3.1.20191231         he28a2e2_2    conda-forge
libev                     4.33                 h516909a_1    conda-forge
libexpat                  2.5.0                hcb278e6_1    conda-forge
libffi                    3.4.2                h7f98852_5    conda-forge
libgcc-ng                 13.2.0               h807b86a_2    conda-forge
libgfortran-ng            13.2.0               h69a702a_2    conda-forge
libgfortran5              13.2.0               ha4646dd_2    conda-forge
libgomp                   13.2.0               h807b86a_2    conda-forge
libjpeg-turbo             3.0.0                hd590300_1    conda-forge
liblapack                 3.9.0           19_linux64_openblas    conda-forge
libnghttp2                1.52.0               h61bc06f_0    conda-forge
libnsl                    2.0.1                hd590300_0    conda-forge
libopenblas               0.3.24          pthreads_h413a1c8_0    conda-forge
libpng                    1.6.39               h753d276_0    conda-forge
libsqlite                 3.43.2               h2797004_0    conda-forge
libssh2                   1.11.0               h0841786_0    conda-forge
libstdcxx-ng              13.2.0               h7e041cc_2    conda-forge
libtiff                   4.6.0                ha9c0a0a_2    conda-forge
libuuid                   2.38.1               h0b41bf4_0    conda-forge
libwebp-base              1.3.2                hd590300_0    conda-forge
libxcb                    1.15                 h0b41bf4_0    conda-forge
libzlib                   1.2.13               hd590300_5    conda-forge
matplotlib-base           3.8.0           py311h54ef318_2    conda-forge
munkres                   1.1.4              pyh9f0ad1d_0    conda-forge
ncurses                   6.4                  hcb278e6_0    conda-forge
numpy                     1.26.0          py311h64a7726_0    conda-forge
openjpeg                  2.5.0                h488ebb8_3    conda-forge
openssl                   3.1.3                hd590300_0    conda-forge
packaging                 23.2               pyhd8ed1ab_0    conda-forge
pillow                    10.1.0          py311ha6c5da5_0    conda-forge
pip                       23.3               pyhd8ed1ab_0    conda-forge
proj                      9.3.0                h1d62c97_1    conda-forge
pthread-stubs             0.4               h36c2ea0_1001    conda-forge
pyparsing                 3.1.1              pyhd8ed1ab_0    conda-forge
pyproj                    3.6.1           py311h1facc83_2    conda-forge
pyshp                     2.3.1              pyhd8ed1ab_0    conda-forge
python                    3.11.6          hab00c5b_0_cpython    conda-forge
python-dateutil           2.8.2              pyhd8ed1ab_0    conda-forge
python_abi                3.11                    4_cp311    conda-forge
readline                  8.2                  h8228510_1    conda-forge
setuptools                68.2.2             pyhd8ed1ab_0    conda-forge
shapely                   2.0.2           py311he06c224_0    conda-forge
six                       1.16.0             pyh6c4a22f_0    conda-forge
sqlite                    3.43.2               h2c6b66d_0    conda-forge
tk                        8.6.13               h2797004_0    conda-forge
tzdata                    2023c                h71feb2d_0    conda-forge
wheel                     0.41.2             pyhd8ed1ab_0    conda-forge
xorg-libxau               1.0.11               hd590300_0    conda-forge
xorg-libxdmcp             1.1.3                h7f98852_0    conda-forge
xz                        5.2.6                h166bdaf_0    conda-forge
zstd                      1.5.5                hfc55251_0    conda-forge

pip list

Package         Version
--------------- ---------
Cartopy         0.22.0
certifi         2023.7.22
contourpy       1.1.1
cycler          0.12.1
fonttools       4.43.1
kiwisolver      1.4.5
matplotlib      3.8.0
munkres         1.1.4
numpy           1.26.0
packaging       23.2
Pillow          10.1.0
pip             23.3
pyparsing       3.1.1
pyproj          3.6.1
pyshp           2.3.1
python-dateutil 2.8.2
setuptools      68.2.2
shapely         2.0.2
six             1.16.0
wheel           0.41.2

[csubich]

Reading /cartopy/blob/main/lib/cartopy/crs.py#L438, , transform_vectors makes the assumption that the transformation is a rotation in 2D Cartesian space, where the numerical vector magnitude must be preserved.

This assumption is not physical. Projecting a vector from one basis (the tangent space of one projection) to another is a general linear transformation, and there's no guarantee that vector magnitudes should be preserved. Hell, if I define 'PlateCareeHalfDegree' which is just a lat-lon projection denominated in 1/2-degree increments, then the vector angles should remain unchanged while the vector magnitudes should double/halve as appropriate.

[greglucas]

The vectors are assumed to be in the "projection" coordinates, not necessarily the geographic coordinates.

See this PR for an attempt to rectify some of this: #1926 which I think is related to your issues here?

[csubich]

This is more than just a question of projection coordinates. transform_vectors , as a transformation between vector spaces (projection or otherwise), should be linear in its components, but it's not.

PR #1926 is more related to the unintuitive behaviour of defining vectors as Δdegree on a lat/lon plot rather than the more commonly-used and intuitive Δdistance (m/s, for example).

The reproducing issue here treats the vectors as Δdeg throughout: Δ1-deg at the equator of the rotated lat-lon projection is Δ√2-deg on the PlateCaree projection. The finite difference calculation gets this correct, but transform_vectors doesn't because it normalizes output vectors to have the same magnitude as the input vectors.