Add invertible Horner polynomials

[georgeouzou]

This enables the user to create horner polynomial transformations
in the real number space without the need to supply "inv" parameters:
"+inv_origin", "+inv_u" and "+inv_v". We can get the inverse
transformation by only using the "fwd" parameters and solving
a system of equations iteratively. The system is usually solved in a
very small number of iterations.

A new optional parameter "+inv_tolerance" is added that can make the
iterative process stop earlier if a certain coordinate accuracy is
met. For example "+inv_tolerance=0.1" will make the iterative
process stop when the calculated coordinates differ less than
0.1 meters, if the units of the transformation destination is meters.

• Tests added
• Added clear title that can be used to generate release notes
• Fully documented, including updating docs/source/*.rst for new API

[georgeouzou]

I will also try to create docs for this.

[georgeouzou]

I also believe there are errors in the horner documentation

• Input and Output types can be Geodetic or projected coordinates
• Equation (3) seems wrong as the final coordinates are not calculated by adding the input X,Y.
• Equation (1) seems also a little off as the Σ u_i,j Uⁱ V^j needs to be Σ u_i,j U^j Vⁱ
• u and v in +fwd_u, +fwd_v, +inv_u, +inv_v, +origin_u, +origin_v can express northing/easting, easting/northing, latitude/longitude, longitude/latitude etc... Parameter documentation implies that u, v expresses northing/easting.

[kbevers]

I can't really comment on the maths but here's a few suggestions for improving the docs. I hope @busstoptaktik can find the time to give this a proper review at some point. He has the expertise in this part of the code base.

[kbevers]

Doc changes looks good to me!

[georgeouzou]

Just added complex inversion too :)

[busstoptaktik]

Certainly a useful addition to the Horner code. See, however, my comments regarding the probable un-soundness of evaluating the polynomium by brute force, rather than Horner's scheme.

Also, IOGP Publication 373-7-2 – Geomatics Guidance Note number 7, part 2 – September 2019 has a very clear section on the criteria for polynomium reversibility. It would probably be a good thing to quote and refer extensively from/to that in the documentation updates. For example, the Geomatics Guidance Note says:

Polynomial reversibility

Approximation polynomials are in a strict mathematical sense not reversible, i.e. the same polynomial coefficients cannot be used to execute the reverse transformation

However, under certain conditions, described below, a satisfactory solution for the reverse transformation may be obtained using the forward coefficient values in a single step, rather than multiple step iteration. If such a solution is possible, in the EPSG Dataset the polynomial coordinate transformation method is classified as a reversible polynomial of degree n.

A (general) polynomial transformation is reversible only when the following conditions are met.

1. The coordinates of source and target evaluation point are (numerically) the same.
2. The unit of measure of the coordinate differences in source and target coordinate reference system are the same.
3. The scaling factors applied to source and target coordinate differences are the same.
4. The spatial variation of the differences between the coordinate reference systems around any given location is sufficiently small.

Your code solves a more general problem, and hence extends on the EPSG sense of a reversible polynomium - it probably would be a good idea to explicitly state this in the documentation.

[georgeouzou]

@busstoptaktik i will make a separate commit for the docs.
As i understand we implemented the general case of the second bullet of Polynomial reversibility section (page 126)

2. By applying the polynomial transformation with the same coefficients but with their signs reversed
   and then iterate to an acceptable solution, the number of iteration steps being dependent on the
   desired accuracy. (Note that only the signs of the polynomial coefficients should be reversed and not
   the coordinates of the evaluation points or the scaling factors!) The iteration procedure is usually
   described by the information source of the polynomial transformation.

So, i will update the docs, so as to make it clear this is an iterative inverse transformation and not a "reversible polynomial" as presented in the publication.

Another thing i realized when i read the section of the polynomial transformations in the publication, is that:

• we do not support using scale factors m_S and m_T
• we use evaluation points in a different way than the usage presented in the publication. Specifically we do:
  X_T = ΔX,
  Y_T = ΔY
  instead of:
  X_T = X_S – X_SO + X_TO + ΔX,
  Y_T = Y_S – Y_SO + Y_TO + ΔY
  where ΔΧ, ΔΥ the result of the polynomium evaluation.

[georgeouzou]

@busstoptaktik i added the doc reference to the publication. Also i refactored the code as requested in the review.