Ocentric (XYZ) Lidar Data Do Not Render Properly

[jlaura]

What is the bug or the crash?

COPC formatted lidar data, with proper projection WKT projection information do not render properly. Sample data are available either from https://na-c.entwine.io/lola.copc.laz or our LOLA in XYZ https://asc-moon.s3-us-west-2.amazonaws.com/hobu/LolaRDR_0N15N_0E15E.copc.laz. These are a global and 15˚ tile of some lunar orbiter laser altimeter data. The WKT projection defined properly in the data is:

    DATUM["Moon(2015)-Sphere",
        ELLIPSOID["Moon(2015)-Sphere",1737400,0,
            LENGTHUNIT["metre",1]]],
    PRIMEM["Reference meridian",0,
        ANGLEUNIT["degree",0.0174532925199433]],
    CS[Cartesian,3],
        AXIS["(X)",geocentricX,
            ORDER[1],
            LENGTHUNIT["metre",1]],
        AXIS["(Y)",geocentricY,
            ORDER[2],
            LENGTHUNIT["metre",1]],
        AXIS["(Z)",geocentricZ,
            ORDER[3],
            LENGTHUNIT["metre",1]]]

Depending on the CRS set for the project, incorrect rendering examples are (project SRS set to the above WKT string):

and (project SRS set to IAU_2015:30100):

Steps to reproduce the issue

Add one of the above linked COPC files to QGIS as a point cloud source. The data are freely available and will stream efficiently to the installation.

Versions

3.34.x, 3.36.x

Supported QGIS version

• [X] I'm running a supported QGIS version according to the roadmap.

New profile

• [X] I tried with a new QGIS profile

Additional context

See moonlidar.com for the correct rendering of these data.

I can only surmise that the issue is some combination of body-centered body-fixed (XYZ meters) rendering and/or the fact that the lunar sphere is used with a 1737.4 km radius while the QGIS code base makes some Earth radii assumption.

[uclaros]

@jlaura Can you please describe what is wrong with the rendering? The link you provided does not have a 2D view of the data to compare. Furthermore, it uses Altitude to style the points which is not an attribute in your point cloud files.

[jlaura]

@uclaros Absolutely! Here is a link to a thumbnail that is generated from these data if they are encoded in the COPC file as lat/lon/height above the datum (which stores poorly), as opposed to XYZ encoding that they are using now. Yes, this thumbnail has been gridded from the point data and is offered only as an example of what these data should render as. I have also included that image inline.

moonlidar.com shows the correct point rendering from the first COPC file linked in the original issue. One sees the lunar topography and not the general red-to-blue, left-to-right gradient that QGIS renders.

As for styling, the points should be x,y and then stylized with the z dimension. If these were to be converted using something like pyproj, I would use:

bcbf = pyproj.Proj(proj='geocent'j, a=1737400, b=1737400)
lla = pyproj.Proj(proj='lonlat', a=1737400, b=1737400)
trans = pyproj.Transformer.from_proj(bcbf, lla)
lon, lat, elev = trans.transform(x,y,z)

where that code snippet assumes the data come in in geocentric XYZ and the desired coordinates for rendering are lat, lon, elevation.

Thanks for your interested in the issue! Please let me know if anything is unclear, I have misunderstood anything, or if I can provide more information.