Broken QgsPoint.project() method or wrong/outdated pyQGIS documentation

[DelazJ]

Trying the code provided in the PyQGIS docs in the Python console, I don't get the same results

p = QgsPoint( 1, 2 ) # 2D point
pr = p.project ( 1, 0 )
# pr is a 2D point: 'Point (1 3)'
pr
<QgsPoint: Point (1 3)>
pr = p.project ( 1, 0, 90 )
# pr is a 2D point: 'Point (1 3)'
pr
<QgsPoint: Point (1 3)>
pr = p.project (1, 0, 0 )
# pr is a 3D point: 'PointZ (1 2 1)'
pr
<QgsPoint: PointZ (1 2 nan)>

p = QgsPoint( QgsWkbTypes.PointZ, 1, 2, 2 ) # 3D point
pr = p.project ( 1, 0 )
# pr is a 3D point: 'PointZ (1 3 2)'
pr
<QgsPoint: PointZM (1001 2 2 2)>
pr = p.project ( 1, 0, 90 )
# pr is a 3D point: 'PointZ (1 3 2)'
pr
<QgsPoint: PointZM (1001 2 2 2)>
pr = p.project (1, 0, 0 )
# pr is a 3D point: 'PointZ (1 2 3)'
pr
<QgsPoint: PointZM (1001 1 3 2)>

QGIS version	3.16.5-Hannover	QGIS code revision	58ba7c1ed6
Compiled against Qt	5.14.2	Running against Qt	5.14.2
Compiled against GDAL/OGR	3.2.1	Running against GDAL/OGR	3.2.1
Compiled against GEOS	3.9.1-CAPI-1.14.2	Running against GEOS	3.9.1-CAPI-1.14.2
Compiled against SQLite	3.31.1	Running against SQLite	3.31.1
PostgreSQL Client Version	12.3	SpatiaLite Version	4.3.0a
QWT Version	6.1.4	QScintilla2 Version	2.11.4
Compiled against PROJ	6.3.2	Running against PROJ	Rel. 6.3.2, May 1st, 2020
OS Version	macOS 10.15
Active python plugins	QuickOSM; QPackage; qgis2web; Qgis2threejs; OSMDownloader; MemoryLayerSaver; qgis_resource_sharing; StreetView; arrayplus; processing; db_manager

[agiudiceandrea]

Hi @DelazJ, it seems to me this

p = QgsPoint( 1, 2 ) # 2D point

pr = p.project (1, 0, 0 )
# pr is a 3D point: 'PointZ (1 2 1)'
pr
<QgsPoint: PointZ (1 2 nan)>

is a bug. Anyway the behaviour seems changed since https://github.com/qgis/QGIS/pull/4712 by @m-kuhn, but I'm not sure it was intended or a mistake (see the final Note in https://github.com/qgis/QGIS/pull/4712#issue-125016805 and these two conflicting commits 7f17498caf8e274b7652efb4f541b2fb30cd5a77 and 96cf4b7026f4fca8106bd8a3a0ec8a25d37e8bf5. See also the original intended behaviour https://github.com/qgis/QGIS/pull/4058#discussion_r97871745)

This

p = QgsPoint( QgsWkbTypes.PointZ, 1, 2, 2 ) # 3D point

is wrong in the PyQGIS Docs (the same error is in QgsGeometryUtils.midpoint() examples).

It should be:

p = QgsPoint( 1, 2, 2, wkbType=QgsWkbTypes.PointZ ) # 3D point

or simply:

p = QgsPoint( 1, 2, 2 ) # 3D point

The following commands will return the correct values.

[m-kuhn]

This was intended as indicated in the final note. I still think it's reasonable, if we don't have a z we cannot calculate with it and want the user to consciusly take care of setting an initial z value.

[agiudiceandrea]

This was intended as indicated in the final note.

So, was "If the point is 2D, the Z value is assumed to be 0." introduced with 96cf4b7026f4fca8106bd8a3a0ec8a25d37e8bf5 a typo? Anyway the current behaviour seems in contrast with the original intended behaviour I agree with: "If I project a 2D point using an inclination I'd expect a 3D point result projected from an initial 0 z value" (https://github.com/qgis/QGIS/pull/4058#discussion_r97871745 by @nyalldawson) otherwise I don't understand the usefulness of this behaviour : "If a 2D point is projected a 3D point will be returned except if inclination is 90".

[m-kuhn]

For me it's a typo.

I don't agree with:

"If I project a 2D point using an inclination I'd expect a 3D point result projected from an initial 0 z value

To me this deviates from the idea of NULL values / unknown values. Like 2 + NULL which is 2 + unknown and not 2 + 0.

Analogue, if there is a point known in 2D space but it's altitude is unknown, a new point in 3D space with known X and Y and unknown Z should be returned.

Otherwise I think a new parameter should be added to specify the default 3D value explicitly (defaultZ = nan).

[agiudiceandrea]

if there is a point known in 2D space but it's altitude is unknown, a new point in 3D space with known X and Y and unknown Z should be returned

Yes, actually this seems to me a valid reason. But then, why "If a 2D point is projected a 3D point will be returned except if inclination is 90"? Why don't simply return always a 2D point, instead of return a PointZ (x+dx, y+dy, nan) if inclination != 90 and a PointZ (x+dx, y+dy) if inclination == 90?

[m-kuhn]

Very valid point. Either we always or never promote to 3D points. I'd also like to highlight that my interpretation is a purely theoretical perspective, that focuses that qgis does not try to smart out the user and rather does nothing than something unexpected. On the other hand, 0 as default will often be desired, so this needs to be easy to achieve. @lbartoletti and @nyalldawson might have some realworld usecases.

[agiudiceandrea]

Just for reference, IINM it seems QgsPoint::project with three parameters (including elevation) is only used in QgsQuadrilateral::rectangleFrom3Points