[PROJ] Horner polynomial evaluation

Tue Nov 21 11:53:15 PST 2023

Dear list members, Dear Jochem,

thanks for the detailed description, it would have been very useful for 
me two years ago, when I prepared a paper, where I used Horner 
Polynomials for transforming ellipsoidal coordinates to spherical 
coordinates.
The oblique Mercator projection used in Hungary (EOV, epsg:23700), is 
defined in GIS environment as Hotine (or Swiss) Oblique Mercator.  As it 
is a two-step projection, using the pipeline operator, an alternative 
method can be defined,that coincides the original definition of the 
projection: first ellipsoidal coordinates are transformed to spherical 
coordinates on the Gauss sphere (aposphere) and in the second step the 
spherical form of the oblique Mercator projection equations are applied 
on these spherical coordinates to get the projected coordinates.
The pipeline operator enables us to apply this two step method, and we 
wanted to use the Horner method, for the first step, to transform the 
ellipsoidal coordinates to spherical coordinates. (The coefficients were 
published in the standardization document, published in 1972, in the 
calculator era).
The manual of Horner Polynomial Evaluation on proj.org states, that both 
input and output coordinates can be either geodetic, either projected, 
but we could not apply the method on geodetic coordinates in a cct 
pipeline.
As a roundabout process, in the pipeline,
1:
the input geodetic coordinates were converted to fake 'projected' 
coordinates, using the
+step +proj=eqc +a=57.295779513082323
that literally does not change(!) the figures, but proj regards the 
output as projected coordinates. The equidistant cylindrical projection 
with +a=57.29... Earth radius ('a' is actually 1 radian in degrees) 
converts 47 degree North latitude and 19 degree East longitude to 47 
meter Northings, and 19 meter Eastings...
2:
We used Horner Polynomial Evolution to modify the coordinates. Proj 
believes that it transforms projected coordinates to projected 
coordinates.
3:
We used the inverse projection, to get ellipsoidal (actually spherical) 
coordinates:
+step +proj=eqc +a=57.295779513082323 +inv
and we finally used the spherical form of the oblique Mercator 
projection.

The details are in the paper:
https://www.fig.net/resources/proceedings/fig_proceedings/fig2022/papers/ts01c/TS01C_molnar_toth_11649.pdf

Regards, Gabor

On 2023-11-21 16:23, Lesparre, Jochem via PROJ wrote:
> Dear list members,
> 
> The implementation of +proj=horner is nice. I used it to reproduce a
> 99 year old polynomial transformation of an even older Dutch CRS
> published in 1827.  However, I encountered some problems with the
> documentation on
> https://proj.org/en/9.3/operations/transformations/horner.html. The
> order of the coefficients in the formulas and in the description of
> the PROJ parameters (+fwd_u and +fwd_u) seem mutually inconsistent.
> 
> PROBLEM 1: The matrix formula (after formula 3) seems inconsistent
> with formula (1) and the actual implementation.
> 
> \begin{bmatrix}
> 
>      u_{0,1} + u_{0,2} U + ... & u_{1,0} + u_{1,1} U + u_{2,0} V + ...
> \\
> 
>      v_{1,0} + v_{1,1} V + v_{2,0} U + ... & v_{0,1} + v_{0,2} V \\
> 
> \end{bmatrix}
> 
> should be changed to:
> 
> \begin{bmatrix}
> 
>      u_{1,0} + u_{2,0} U + ... & u_{0,1} + u_{1,1} U + u_{0,2} V + ...
> \\
> 
>     v_{0,1} + v_{1,1} V + v_{0,2} U + ... & v_{1,0} + v_{2,0} V \\
> 
> \end{bmatrix}
> 
> PROBLEM 2: The description of the parameters seems inconsistent with
> formula (1) and the actual implementation.
> 
> +fwd_u=<u_11,u_12,...,u_ij,..,u_mn>
> 
> and
> 
> +fwd_v=<v_11,v_12,...,v_ij,..,v_mn>
> 
> should be changed to:
> 
> +fwd_u=<u_00,u_10,...,u_ij,..,u_nn>
> 
> and
> 
> +fwd_v=<v_00,v_01,...,v_ij,..,v_nn>
> 
> As a result it was difficult to determine the correct order of the
> coefficients for +fwd_u and +fwd_v parameters. I think it would be
> clarifying to add a table for the order of the coefficients in the
> documentation (example for +deg=3):
> 
> Order of coefficients of +fwd_u:
> 
>     | x^0 | x^1 | x^2 | x^3
> 
> y^0 | 1   | 2   | 3   | 4
> 
> y^1 | 5   | 6   | 7   | -
> 
> y^2 | 8   | 9   | -   | -
> 
> y^3 | 10  | -   | -   | -
> 
> Order of coefficients of +fwd_v:
> 
>     | x^0 | x^1 | x^2 | x^3
> 
> y^0 | 1   | 5   | 8   | 10
> 
> y^1 | 2   | 6   | 9   | -
> 
> y^2 | 3   | 7   | -   | -
> 
> y^3 | 4   | -   | -   | -
> 
> SUGGESTION FOR A DIFFERENT IMPLEMENTATION: Actually, I would have
> preferred a different order of the coefficients in the implementation,
> in a way that the order would be the same for +fwd_u and +fwd_v, and
> it would keep the first coefficients the same when increasing the
> degree (+deg) of the polynomial:
> 
> Order of coefficients of +fwd_u and +fwd_v:
> 
>     | x^0 | x^1 | x^2 | x^3
> 
> y^0 | 1   | 2   | 4   | 7
> 
> y^1 | 3   | 5   | 8   | -
> 
> y^2 | 6   | 9   | -   | -
> 
> y^3 | 10  | -   | -   | -
> 
> Regards, Jochem
> 
> Disclaimer:
> De inhoud van deze e-mail is vertrouwelijk en uitsluitend bestemd voor
> de geadresseerde(n).
> Gebruik, openbaarmaking, vermenigvuldiging, verspreiding en/of
> verstrekking van deze informatie aan derden is niet toegestaan.
> Op al onze producten en diensten zijn onze algemene
> leveringsvoorwaarden van toepassing
> [https://www.kadaster.nl/algemene-leveringsvoorwaarden].
> 
> Disclaimer:
> This email and any files transmitted with it are confidential and
> intended solely for the use of the individual or entity to whom they
> are addressed.
> If you are not the intended recipient, you are notified that
> disclosing, copying, distributing or taking any action in reliance on
> the contents of this information is strictly prohibited.
> Our general terms and conditions of delivery apply to all our products
> and services
> [https://www.kadaster.com/general-terms-and-conditions].
> _______________________________________________
> PROJ mailing list
> PROJ at lists.osgeo.org
> https://lists.osgeo.org/mailman/listinfo/proj