DERIVEDPROJCRS with 3D Helmert deriving conversion: transform to base CRS works, but to EPSG:4326 fails with "Mismatched units"

Start with the unit-class compatibility checks in pipeline.cpp, then trace how derived CRS operations are created and how the deriving conversion relates to the base and derived CRS. Reproduce with the supplied WKT and compare transformations to EPSG:7856 and EPSG:4326. Done when incompatible deriving-conversion units are consistently rejected or warned about, including in the single-step case.

AI usage

• AI (Copilot or something similar) supported the creation of this ticket. See our policy about AI tool use. Use of AI tools must be indicated.
• This ticket has been (mostly) written by myself (human), and is aimed at being read by humans. We do not accept AI generated wall-of-text

(The reproducer and all numbers were produced and checked by hand; an AI assistant helped draft and edit the text.)

Example of problem

site.wkt is a DERIVEDPROJCRS over EPSG:7856 (GDA2020 / MGA zone 56) whose DERIVINGCONVERSION is a PROJ-based seven-parameter +proj=helmert. This is a real WKT definition we received from a user - the version here has been modified with synthetic parameters for confidentiality.

The WKT definition is wrong: the 3D Helmert works on geocentric cartesian coordinates, not on projected easting/northing. The issue is not that PROJ rejects it, but that it only partly does:

DERIVEDPROJCRS["Example Site Grid",
    BASEPROJCRS["GDA2020 / MGA zone 56",
        BASEGEOGCRS["GDA2020",
            DATUM["Geocentric Datum of Australia 2020",
                ELLIPSOID["GRS 1980",6378137,298.257222101,
                    LENGTHUNIT["metre",1]]],
            PRIMEM["Greenwich",0,
                ANGLEUNIT["degree",0.0174532925199433]],
            ID["EPSG",7844]],
        CONVERSION["Map Grid of Australia zone 56",
            METHOD["Transverse Mercator",
                ID["EPSG",9807]],
            PARAMETER["Latitude of natural origin",0,
                ANGLEUNIT["degree",0.0174532925199433],
                ID["EPSG",8801]],
            PARAMETER["Longitude of natural origin",153,
                ANGLEUNIT["degree",0.0174532925199433],
                ID["EPSG",8802]],
            PARAMETER["Scale factor at natural origin",0.9996,
                SCALEUNIT["unity",1],
                ID["EPSG",8805]],
            PARAMETER["False easting",500000,
                LENGTHUNIT["metre",1],
                ID["EPSG",8806]],
            PARAMETER["False northing",10000000,
                LENGTHUNIT["metre",1],
                ID["EPSG",8807]],
            ID["EPSG",17356]],
        ID["EPSG",7856]],
    DERIVINGCONVERSION["Site calibration as PROJ-based",
        METHOD["PROJ-based operation method: +proj=pipeline +step +inv +proj=helmert +x=100 +y=200 +z=-30 +rx=1 +ry=1 +rz=-1 +s=-5 +convention=position_vector"]],
    CS[Cartesian,2],
        AXIS["site east (x)",east,
            ORDER[1],
            LENGTHUNIT["metre",1]],
        AXIS["site north (y)",north,
            ORDER[2],
            LENGTHUNIT["metre",1]]]

The CRS parses without any errors:

$ projinfo "$(cat site.wkt)" -o WKT2:2019 -q | head -1
DERIVEDPROJCRS["Example Site Grid",

Transformation to the base CRS is a lone helmert step. It instantiates and runs, and the 2D input comes out with a non-zero Z: easting/northing were treated as geocentric X/Y with Z = 0.

$ projinfo -s "$(cat site.wkt)" -t EPSG:7856 -o PROJ -q
+proj=helmert +x=100 +y=200 +z=-30 +rx=1 +ry=1 +rz=-1 +s=-5
      +convention=position_vector

$ echo "500100 7284600" | cs2cs -d 4 "$(cat site.wkt)" EPSG:7856
500232.8161    7284761.1525 2.8920

Transformation to EPSG:4326 or EPSG:28356 (GDA94 / MGA zone 56) is a pipeline whose first two steps are the helmert and the inverse utm. projinfo -o PROJ exports it, but it cannot be instantiated:

$ projinfo -s "$(cat site.wkt)" -t EPSG:4326 -o PROJ -q
+proj=pipeline
  +step +proj=helmert +x=100 +y=200 +z=-30 +rx=1 +ry=1 +rz=-1 +s=-5
        +convention=position_vector
  +step +inv +proj=utm +zone=56 +south +ellps=GRS80
  +step +proj=unitconvert +xy_in=rad +xy_out=deg
  +step +proj=axisswap +order=2,1

$ echo "500100 7284600" | cs2cs -d 6 "$(cat site.wkt)" EPSG:4326
pipeline: Pipeline: Mismatched units between step 1 and 2
[abbreviated stderr: four further unit-mismatch diagnostics omitted, one of them for steps 2 and 3]
*    * inf

$ echo "500100 7284600" | cs2cs -d 4 "$(cat site.wkt)" EPSG:28356
pipeline: Pipeline: Mismatched units between step 1 and 2
*    * inf

pyproj behaves the same: Transformer.from_crs(crs, "EPSG:7856") works, "EPSG:4326" and "EPSG:28356" raise
Input is not a transformation.: (Internal Proj Error: pipeline: Pipeline: Mismatched units between step 1 and 2).

Problem description

The "Mismatched units" check in pipeline.cpp only compares adjacent steps: the 7-parameter helmert declares
PJ_IO_UNITS_CARTESIAN, the inverse utm step expects PJ_IO_UNITS_PROJECTED, so a pipeline chaining them is
rejected. Nothing compares the deriving conversion with the CRS it is attached to. As a result:

• The single-step operation to the base CRS has no step boundary to check, so it runs and produces wrong numbers.
• A multi-step operation is only caught when a neighbour declares an incompatible class. For example, a pipeline from this CRS to a second derived CRS using +proj=affine also runs, because affine declares PJ_IO_UNITS_WHATEVER and inherits the Helmert's cartesian class during pipeline construction.
• The problem only shows up when transforming to a third CRS. This is a problem for a system that stores the accepted CRS along with the data but cannot transform it later.

I am not asking for the 3D helmert to accept projected coordinates - A site calibration on E/N should use a
projected-space operation (+proj=affine, or the 2D helmert with +theta, which declares projected units). Instead, this request is for the incompatibility to be detected consistently, and ideally early:

1. Check the deriving conversion's input unit class against the base CRS (and its output against the derived CRS) when the operation is created, so the single-step case fails the same way the multi-step case does.
2. Ideally report it at parse time (proj_create, projinfo <wkt>), so authors of such definitions find out immediately.

Expected Output

Either proj_create() rejects (or at least warns about) the CRS, or the transformation to EPSG:7856 fails with the
same "Mismatched units" error as the transformation to EPSG:4326, instead of returning 500232.8161 7284761.1525 2.8920.

For comparison, the same CRS with a projected-space +proj=affine deriving conversion (using illustrative parameters, not an equivalent of the Helmert above) transforms to EPSG:7856, EPSG:4326 and EPSG:28356 without error, and Z stays 0 for 2D input. So the step check itself is fine; what is missing is a check between the deriving conversion and the CRS it belongs to.

Environment Information

• PROJ 9.4.0 (Ubuntu 24.04 proj-bin/proj-data 9.4.0-1build2), PROJ 9.8.1 (pyproj 3.8.0 wheel) and PROJ 9.9.0 (built from the release tarball, ENABLE_TIFF=OFF ENABLE_CURL=OFF). Same results on all three; the output above
  is from 9.9.0. PROJ 9.2.0 fails the same transformations with less diagnostic output: cs2cs prints only * * inf for EPSG:4326 and aborts with "No inverse operation" for EPSG:28356, while cct on the exported pipeline reports the same "Mismatched units between step 1 and 2".
• Results are identical with PROJ_NETWORK on or off and with or without the Australian grid files installed: for the EPSG:7856 to EPSG:28356 leg PROJ selects the 7-parameter Helmert EPSG:8048 (0.01 m) ahead of the grid-based
  operations (0.05 m).
• Ubuntu 24.04 under WSL2, x86_64.

Installation method

apt-get (9.4.0), pip wheel (9.8.1), from source (9.9.0).