Skip to content

Scan and shatter are failing when reprojecting point data #113

Description

@bmcgaughey1

I am testing with NOAA COPC data that is in geographic CRS (EPSG:6318) but I suspect this problem will come up any time we want to reproject point data.

My code gets the bounding box from PDAL quickinfo, then transforms the box into UTM 8N (EPSG:26908) with extra points along edges to account for distortion in the transform (and deals with axis order). DB is created using the projected bounds for all point tiles. Shatter is called using the transformed bounding box for each file. PDAL pipeline created and fed into scan and shatter includes a filters.reprojection stage at the end of the pipeline.

Problem is that the bounds argument added by SilviMetric to the reader stage for scan and shatter uses the desired output CRS (UTM 8N), so no points are returned from the stage for further processing.

Not sure how we deal with this. I think it makes sense to do reprojection of bounds so we can use the projected extent for initialize to set the DB extent. Then we have to use the projected bounds for each tile, but SilviMetric's pipeline modification doesn't allow a bounds argument for the reader that is in a different projection from the point data.

We could pass the original geographic bounds to shatter for each point tile, but this fails since DB was created with projected bounds.

URLs for testing data (used data locally to speed up runs):
PDAL reports CRS as 4326: WGS84)
assets = [
'https://noaa-nos-coastal-lidar-pds.s3.us-east-1.amazonaws.com/laz/geoid12b/10045/20230707_TNFWI_664000_6244000.copc.laz',
'https://noaa-nos-coastal-lidar-pds.s3.us-east-1.amazonaws.com/laz/geoid12b/10045/20230707_TNFWI_664000_6243000.copc.laz',
'https://noaa-nos-coastal-lidar-pds.s3.us-east-1.amazonaws.com/laz/geoid12b/10045/20230707_TNFWI_665000_6244000.copc.laz',
'https://noaa-nos-coastal-lidar-pds.s3.us-east-1.amazonaws.com/laz/geoid12b/10045/20230707_TNFWI_665000_6243000.copc.laz'
]

pipeline given to scan and shatter (bounds is added before the filename in the reader.copc stage):
"pipeline":
[
{
"filename": "H:/NOAATestData/20230707_TNFWI_664000_6243000.copc.laz",
"override_srs": "{"type": "CompoundCRS", "name": "NAD83(2011); NAVD88 height", "components": [{"type": "GeographicCRS", "name": "NAD83(2011)", "datum": {"type": "GeodeticReferenceFrame", "name": "NAD83 (National Spatial Reference System 2011)", "ellipsoid": {"name": "GRS 1980", "semi_major_axis": 6378137, "inverse_flattening": 298.257222101}}, "coordinate_system": {"subtype": "ellipsoidal", "axis": [{"name": "Geodetic latitude", "abbreviation": "Lat", "direction": "north", "unit": "degree"}, {"name": "Geodetic longitude", "abbreviation": "Lon", "direction": "east", "unit": "degree"}]}, "id": {"authority": "EPSG", "code": 6318}}, {"type": "VerticalCRS", "name": "NAVD88 height", "datum": {"type": "VerticalReferenceFrame", "name": "North American Vertical Datum 1988"}, "coordinate_system": {"subtype": "vertical", "axis": [{"name": "Up", "abbreviation": "", "direction": "up", "unit": "metre"}]}, "id": {"authority": "EPSG", "code": 5703}}]}",
"tag": "readers_copc1",
"type": "readers.copc"
},
{
"expression": "(((Classification != 7) && (Classification != 9) && (Classification != 18)) && (Synthetic != 1) && (Withheld != 1))",
"inputs":
[
"readers_copc1"
],
"tag": "filters_expression1",
"type": "filters.expression"
},
{
"allow_extrapolation": "true",
"inputs":
[
"filters_expression1"
],
"tag": "filters_hag_delaunay1",
"type": "filters.hag_delaunay"
},
{
"expression": "HeightAboveGround >= 2.0 && HeightAboveGround <= 150.0",
"inputs":
[
"filters_hag_delaunay1"
],
"tag": "filters_expression2",
"type": "filters.expression"
},
{
"in_srs": "{"type": "CompoundCRS", "name": "NAD83(2011); NAVD88 height", "components": [{"type": "GeographicCRS", "name": "NAD83(2011)", "datum": {"type": "GeodeticReferenceFrame", "name": "NAD83 (National Spatial Reference System 2011)", "ellipsoid": {"name": "GRS 1980", "semi_major_axis": 6378137, "inverse_flattening": 298.257222101}}, "coordinate_system": {"subtype": "ellipsoidal", "axis": [{"name": "Geodetic latitude", "abbreviation": "Lat", "direction": "north", "unit": "degree"}, {"name": "Geodetic longitude", "abbreviation": "Lon", "direction": "east", "unit": "degree"}]}, "id": {"authority": "EPSG", "code": 6318}}, {"type": "VerticalCRS", "name": "NAVD88 height", "datum": {"type": "VerticalReferenceFrame", "name": "North American Vertical Datum 1988"}, "coordinate_system": {"subtype": "vertical", "axis": [{"name": "Up", "abbreviation": "", "direction": "up", "unit": "metre"}]}, "id": {"authority": "EPSG", "code": 5703}}]}",
"inputs":
[
"filters_expression2"
],
"out_srs": "EPSG:26908",
"tag": "filters_reprojection1",
"type": "filters.reprojection"
},
{
"dimensions": "HeightAboveGround=>Z",
"inputs":
[
"filters_reprojection1"
],
"tag": "filters_ferry1",
"type": "filters.ferry"
}
]
}

Metadata

Metadata

Assignees

No one assigned

    Labels

    No labels
    No labels

    Type

    No type

    Fields

    No fields configured for issues without a type.

    Projects

    No projects

    Milestone

    No milestone

    Relationships

    None yet

    Development

    No branches or pull requests

    Issue actions