[pdal] pdal tindex --merge problem with polygon clippping
Stefan Ziegler
stefan.ziegler.de at gmail.com
Fri Aug 14 04:01:14 PDT 2015
Hi
I created a small test case here:
http://www.catais.org/tmp/pdal_tindex.zip
It contains a small las file and the tileindex file. This is what I've done:
1) Create tileindex:
find . -iname "*.las" | pdal tindex tileindex.shp --a_srs EPSG:21781
--t_srs EPSG:21781 --fast_boundary --stdin
2) Commented out line 161 to allow "--t_srs". Pdal now does *not* reproject
data anymore to EPSG:4326. But when trying to write the las file it stops
with the following error:
PDAL: Unable to convert scaled value (2.04e+10) to int32 for dimension 'X'
when writing LAS/LAZ file tmp.las.
3) Commented out writerOptions.add() methods for scales and offsets on line
445 - 450 and added setCommonOptions(writerOptions).
Now I get the desired result by executing:
pdal tindex --scale "1 1 1" --merge tileindex.shp --lyr_name tileindex
--a_srs EPSG:21781 --t_srs EPSG:21781 tmp.las
So, this works for me but definitely introduces some evil stuff :) See also
attached patch to exactly see what I've done.
regards
Stefan
On Fri, Aug 14, 2015 at 6:16 AM, Stefan Ziegler <stefan.ziegler.de at gmail.com
> wrote:
> find . -iname "*.laz" | pdal tindex tileindex.shp --a_srs EPSG:21781
> --t_srs EPSG:21781 --fast_boundary --stdin
>
> ogrinfo -al tileindex.shp:
>
> Layer name: tileindex
> Geometry: Polygon
> Feature Count: 6
> Extent: (592000.000000, 226000.000000) - (594000.000000, 230000.000000)
> Layer SRS WKT:
> PROJCS["CH1903_LV03",
> GEOGCS["GCS_CH1903",
> DATUM["CH1903",
> SPHEROID["Bessel_1841",6377397.155,299.1528128]],
> PRIMEM["Greenwich",0],
> UNIT["Degree",0.017453292519943295]],
> PROJECTION["Hotine_Oblique_Mercator_Azimuth_Center"],
> PARAMETER["latitude_of_center",46.95240555555556],
> PARAMETER["longitude_of_center",7.439583333333333],
> PARAMETER["azimuth",90],
> PARAMETER["scale_factor",1],
> PARAMETER["false_easting",600000],
> PARAMETER["false_northing",200000],
> PARAMETER["rectified_grid_angle",90],
> UNIT["Meter",1]]
> location: String (254.0)
> srs: String (254.0)
> modified: Date (10.0)
> created: Date (10.0)
> OGRFeature(tileindex):0
> location (String) = /home/stefan/tmp/lidar/srs/./LAS_593227.laz
> srs (String) = EPSG:21781
> modified (Date) = 2015/08/13
> created (Date) = 2015/08/13
> POLYGON ((593000 227000,593000 228000,594000 228000,594000 227000,593000
> 227000))
>
> pdal info LAS_593227.laz --metadata:
>
> {
> "filename": "LAS_593227.laz",
> "metadata":
> {
> "comp_spatialreference": "PROJCS[\"CH1903 /
> LV03\",GEOGCS[\"CH1903\",DATUM[\"CH1903\",SPHEROID[\"Bessel
> 1841\",6377397.155,299.1528128,AUTHORITY[\"EPSG\",\"7004\"]],TOWGS84[674.4,15.1,405.3,0,0,0,0],AUTHORITY[\"EPSG\",\"6149\"]],PRIMEM[\"Greenwich\",0,AUTHORITY[\"EPSG\",\"8901\"]],UNIT[\"degree\",0.0174532925199433,AUTHORITY[\"EPSG\",\"9122\"]],AUTHORITY[\"EPSG\",\"4149\"]],PROJECTION[\"Hotine_Oblique_Mercator_Azimuth_Center\"],PARAMETER[\"latitude_of_center\",46.95240555555556],PARAMETER[\"longitude_of_center\",7.439583333333333],PARAMETER[\"azimuth\",90],PARAMETER[\"rectified_grid_angle\",90],PARAMETER[\"scale_factor\",1],PARAMETER[\"false_easting\",600000],PARAMETER[\"false_northing\",200000],UNIT[\"metre\",1,AUTHORITY[\"EPSG\",\"9001\"]],AXIS[\"Y\",EAST],AXIS[\"X\",NORTH],AUTHORITY[\"EPSG\",\"21781\"]]",
> "compressed": true,
> "count": 14724002,
> "creation_doy": 224,
> "creation_year": 2015,
> "dataformat_id": 3,
> "dataoffset": 2181,
> "filesource_id": 0,
> "global_encoding": "AAA=",
> "header_size": 227,
> "major_version": 1,
> "maxx": 593999.98999999999,
> "maxy": 227999.98999999999,
> "maxz": 1667.6600000000001,
> "minor_version": 2,
> "minx": 593000,
> "miny": 227000,
> "minz": 622.80000000000007,
> "offset_x": 0,
> "offset_y": 0,
> "offset_z": 0,
> "project_id": "00000000-0000-0000-0000-000000000000",
> "scale_x": 0.01,
> "scale_y": 0.01,
> "scale_z": 0.01,
> "software_id": "PDAL 1.0.0.b1 (e412bd)",
> "spatialreference": "PROJCS[\"CH1903 /
> LV03\",GEOGCS[\"CH1903\",DATUM[\"CH1903\",SPHEROID[\"Bessel
> 1841\",6377397.155,299.1528128,AUTHORITY[\"EPSG\",\"7004\"]],TOWGS84[674.4,15.1,405.3,0,0,0,0],AUTHORITY[\"EPSG\",\"6149\"]],PRIMEM[\"Greenwich\",0,AUTHORITY[\"EPSG\",\"8901\"]],UNIT[\"degree\",0.0174532925199433,AUTHORITY[\"EPSG\",\"9122\"]],AUTHORITY[\"EPSG\",\"4149\"]],PROJECTION[\"Hotine_Oblique_Mercator_Azimuth_Center\"],PARAMETER[\"latitude_of_center\",46.95240555555556],PARAMETER[\"longitude_of_center\",7.439583333333333],PARAMETER[\"azimuth\",90],PARAMETER[\"rectified_grid_angle\",90],PARAMETER[\"scale_factor\",1],PARAMETER[\"false_easting\",600000],PARAMETER[\"false_northing\",200000],UNIT[\"metre\",1,AUTHORITY[\"EPSG\",\"9001\"]],AXIS[\"Y\",EAST],AXIS[\"X\",NORTH],AUTHORITY[\"EPSG\",\"21781\"]]",
> "system_id": "PDAL"
> },
> "pdal_version": "PDAL 1.0.0.b1 (e412bd)"
> }
>
> In the mergeFile() method on line 413 you assign the output srs for the
> reprojection stage. But this cannot be set by the user and will therefore
> always be EPSG:4326. As far as I understand.
>
> regards
> Stefan
>
>
>
>
> On Fri, Aug 14, 2015 at 12:10 AM, Andrew Bell <andrew.bell.ia at gmail.com>
> wrote:
>
>> On Thu, Aug 13, 2015 at 2:46 PM, Stefan Ziegler <
>> stefan.ziegler.de at gmail.com> wrote:
>>
>>> Hi
>>>
>>> I'm trying to merge some las files with the following command:
>>>
>>> pdal tindex --merge /home/stefan/tmp/lidar/srs/tileindex.gpkg --lyr_name
>>> tileindex tmp.las
>>>
>>> The input las files have an srs assigned (EPSG:21781). But the merged
>>> one is EPSG:4326.
>>>
>>
>> Can you use ogrinfo to verify the SRS of the layer? Can you provide the
>> command used to create the index file?
>>
>> Thanks,
>>
>> --
>> Andrew Bell
>> andrew.bell.ia at gmail.com
>>
>
>
-------------- next part --------------
An HTML attachment was scrubbed...
URL: <http://lists.osgeo.org/pipermail/pdal/attachments/20150814/3dc57517/attachment-0001.html>
-------------- next part --------------
A non-text attachment was scrubbed...
Name: TIndexKernel.patch
Type: text/x-patch
Size: 1365 bytes
Desc: not available
URL: <http://lists.osgeo.org/pipermail/pdal/attachments/20150814/3dc57517/attachment-0001.bin>
More information about the pdal
mailing list