3

I have a point cloud with extra attributes:

class        : LAS (LASFˇ v1.4)
point format : 6
memory       : 2 Gb 
extent       :7498400, 7499200, 5779500, 5780000 (xmin, xmax, ymin, ymax)
coord. ref.  : NA 
area         : 399973.7 units² (convex hull)
points       : 19517984 points
density      : 48.8 points/units²
names        : X Y Z gpstime Intensity ReturnNumber NumberOfReturns ScanDirectionFlag EdgeOfFlightline Classification ScannerChannel Synthetic_flag Keypoint_flag Withheld_flag Overlap_flag ScanAngle UserData PointSourceID Amplitude Reflectance Deviation 
File signature:           LASFˇ 
File source ID:           4001 
Project ID - GUID:        df4425f5-d89d-437c-baa3-cd71726ed8c6 
Version:                  1.4
System identifier:        LAStools (c) by rapidlasso GmbH 
Generating software:      lasheight (190321) commercial 
File creation d/y:        267/2018
header size:              375 
Offset to point data:     2238 
Num. var. length record:  5 
Point data format:        6 
Point data record length: 36 
Num. of point records:    19517984 
Num. of points by return: 
Scale factor X Y Z:       0.01 0.01 0.01 
Offset X Y Z:             0 0 0 
min X Y Z:                7498400 5779500 -1.17 
max X Y Z:                7499200 5780000 25.94 
Variable length records: 
   Variable length record 1 of 5 
       Description:  GeoKeyDirectoryTag (mandatory) 
       Tags:
          Key 1024 tiff_tag_location 0 count 1 value offset 1 
          Key 1025 tiff_tag_location 0 count 1 value offset 2 
          Key 1026 tiff_tag_location 34737 count 13 value offset 0 
          Key 2048 tiff_tag_location 0 count 1 value offset 32767 
          Key 2049 tiff_tag_location 34737 count 7 value offset 26 
          Key 2050 tiff_tag_location 0 count 1 value offset 32767 
          Key 2051 tiff_tag_location 0 count 1 value offset 32767 
          Key 2052 tiff_tag_location 0 count 1 value offset 9001 
          Key 2054 tiff_tag_location 0 count 1 value offset 9102 
          Key 2056 tiff_tag_location 0 count 1 value offset 32767 
          Key 2057 tiff_tag_location 34736 count 1 value offset 0 
          Key 2059 tiff_tag_location 34736 count 1 value offset 1 
          Key 2060 tiff_tag_location 0 count 1 value offset 9102 
          Key 2061 tiff_tag_location 34736 count 1 value offset 2 
          Key 3072 tiff_tag_location 0 count 1 value offset 32767 
          Key 3073 tiff_tag_location 34737 count 13 value offset 13 
          Key 3074 tiff_tag_location 0 count 1 value offset 32767 
          Key 3075 tiff_tag_location 0 count 1 value offset 1 
          Key 3076 tiff_tag_location 0 count 1 value offset 9001 
          Key 3080 tiff_tag_location 34736 count 1 value offset 4 
          Key 3081 tiff_tag_location 34736 count 1 value offset 3 
          Key 3082 tiff_tag_location 34736 count 1 value offset 5 
          Key 3083 tiff_tag_location 34736 count 1 value offset 6 
          Key 3092 tiff_tag_location 34736 count 1 value offset 7 
          Key 4099 tiff_tag_location 0 count 1 value offset 9001 
   Variable length record 2 of 5 
       Description:  GeoDoubleParamsTag (optional) 
       data:                 6378137 298.2572 0 0 21 7500000 0 0.999923 
   Variable length record 3 of 5 
       Description:  GeoASCIIParamsTag (optional) 
       data:                 Poland2000s7|Poland2000s7|ETRS89| 
   Variable length record 4 of 5 
       Description:  OGC COORDINATE SYSTEM WKT 
   Variable length record 5 of 5 
       Description:  RIEGL Extra Bytes 
       Extra Bytes Description:
          Amplitude: Echo signal amplitude [dB]
          Reflectance: Echo signal reflectance [dB]
          Deviation: Pulse shape deviation

When I check Classification attribute: 

> max(las@data$Classification)
    [1] 5

I have also a TIF raster with treeIDs. I am using following algorithm to assign treeIDs to point cloud:

las_tree = lasmergespatial(las, tif_crop, "treeID")
treeID = las_tree@data[["treeID"]]
treeID = as.numeric(treeID)
las_tree = lasaddextrabytes(las_tree, x=treeID, name="treeID", desc="ID of a tree")

Now, las_tree has summary:

class        : LAS (LASFˇ v1.4)
point format : 6
memory       : 2.2 Gb 
extent       :7498400, 7499200, 5779500, 5780000 (xmin, xmax, ymin, ymax)
coord. ref.  : NA 
area         : 399973.7 units² (convex hull)
points       : 19517984 points
density      : 48.8 points/units²
names        : X Y Z gpstime Intensity ReturnNumber NumberOfReturns ScanDirectionFlag EdgeOfFlightline Classification ScannerChannel Synthetic_flag Keypoint_flag Withheld_flag Overlap_flag ScanAngle UserData PointSourceID Amplitude Reflectance Deviation treeID 
File signature:           LASFˇ 
File source ID:           4001 
Project ID - GUID:        df4425f5-d89d-437c-baa3-cd71726ed8c6 
Version:                  1.4
System identifier:        LAStools (c) by rapidlasso GmbH 
Generating software:      lasheight (190321) commercial 
File creation d/y:        267/2018
header size:              375 
Offset to point data:     2238 
Num. var. length record:  5 
Point data format:        6 
Point data record length: 36 
Num. of point records:    19517984 
Num. of points by return: 
Scale factor X Y Z:       0.01 0.01 0.01 
Offset X Y Z:             0 0 0 
min X Y Z:                7498400 5779500 -1.17 
max X Y Z:                7499200 5780000 25.94 
Variable length records: 
   Variable length record 1 of 5 
       Description:  GeoKeyDirectoryTag (mandatory) 
       Tags:
          Key 1024 tiff_tag_location 0 count 1 value offset 1 
          Key 1025 tiff_tag_location 0 count 1 value offset 2 
          Key 1026 tiff_tag_location 34737 count 13 value offset 0 
          Key 2048 tiff_tag_location 0 count 1 value offset 32767 
          Key 2049 tiff_tag_location 34737 count 7 value offset 26 
          Key 2050 tiff_tag_location 0 count 1 value offset 32767 
          Key 2051 tiff_tag_location 0 count 1 value offset 32767 
          Key 2052 tiff_tag_location 0 count 1 value offset 9001 
          Key 2054 tiff_tag_location 0 count 1 value offset 9102 
          Key 2056 tiff_tag_location 0 count 1 value offset 32767 
          Key 2057 tiff_tag_location 34736 count 1 value offset 0 
          Key 2059 tiff_tag_location 34736 count 1 value offset 1 
          Key 2060 tiff_tag_location 0 count 1 value offset 9102 
          Key 2061 tiff_tag_location 34736 count 1 value offset 2 
          Key 3072 tiff_tag_location 0 count 1 value offset 32767 
          Key 3073 tiff_tag_location 34737 count 13 value offset 13 
          Key 3074 tiff_tag_location 0 count 1 value offset 32767 
          Key 3075 tiff_tag_location 0 count 1 value offset 1 
          Key 3076 tiff_tag_location 0 count 1 value offset 9001 
          Key 3080 tiff_tag_location 34736 count 1 value offset 4 
          Key 3081 tiff_tag_location 34736 count 1 value offset 3 
          Key 3082 tiff_tag_location 34736 count 1 value offset 5 
          Key 3083 tiff_tag_location 34736 count 1 value offset 6 
          Key 3092 tiff_tag_location 34736 count 1 value offset 7 
          Key 4099 tiff_tag_location 0 count 1 value offset 9001 
   Variable length record 2 of 5 
       Description:  GeoDoubleParamsTag (optional) 
       data:                 6378137 298.2572 0 0 21 7500000 0 0.999923 
   Variable length record 3 of 5 
       Description:  GeoASCIIParamsTag (optional) 
       data:                 Poland2000s7|Poland2000s7|ETRS89| 
   Variable length record 4 of 5 
       Description:  OGC COORDINATE SYSTEM WKT 
   Variable length record 5 of 5 
       Description:  rlas extra bytes 
       Extra Bytes Description:
          Amplitude: Echo signal amplitude [dB]
          Reflectance: Echo signal reflectance [dB]
          Deviation: Pulse shape deviation
          treeID: ID of a tree

Checking Classification attribute again:

> max(las_tree@data$Classification)
[1] 5

When I write las_tree with writeLAS:

writeLAS(las_tree, savedFilename)

and read it again with: 

las_read <- readLAS(paste(path, savedFilename, sep="\\"))

seems like I lost Classification attribute values:

> max(las_read@data$Classification)
[1] 0

Summary of las_read:

class        : LAS (LASFˇ v1.4)
point format : 6
memory       : 2.2 Gb 
extent       :7498400, 7499200, 5779500, 5780000 (xmin, xmax, ymin, ymax)
coord. ref.  : NA 
area         : 399973.7 units² (convex hull)
points       : 19517984 points
density      : 48.8 points/units²
names        : X Y Z gpstime Intensity ReturnNumber NumberOfReturns ScanDirectionFlag EdgeOfFlightline Classification ScannerChannel Synthetic_flag Keypoint_flag Withheld_flag Overlap_flag ScanAngle UserData PointSourceID Amplitude Reflectance Deviation treeID 
File signature:           LASFˇ 
File source ID:           4001 
Project ID - GUID:        df4425f5-d89d-437c-baa3-cd71726ed8c6 
Version:                  1.4
System identifier:         
Generating software:      rlas R package 
File creation d/y:        0/2018
header size:              375 
Offset to point data:     2034 
Num. var. length record:  3 
Point data format:        6 
Point data record length: 44 
Num. of point records:    19517984 
Num. of points by return: 
Scale factor X Y Z:       0.01 0.01 0.01 
Offset X Y Z:             0 0 0 
min X Y Z:                7498400 5779500 -1.17 
max X Y Z:                7499200 5780000 25.94 
Variable length records: 
   Variable length record 1 of 3 
       Description:  by LAStools of rapidlasso GmbH 
       Tags:
          Key 3072 tiff_tag_location 0 count 1 value offset 32767 
   Variable length record 2 of 3 
       Description:  by LAStools of rapidlasso GmbH 
   Variable length record 3 of 3 
       Description:  by LAStools of rapidlasso GmbH 
       Extra Bytes Description:
          Amplitude: Echo signal amplitude [dB]
          Reflectance: Echo signal reflectance [dB]
          Deviation: Pulse shape deviation
          treeID: ID of a tree

Here's the full processing:

las = readLAS(filename_path, filter = "-keep_class 2 3 4 5")
tif_crop = raster::crop(tif, raster::extent(las)*1.1)
las_tree = lasmergespatial(las, tif_crop, "treeID")
treeID = las_tree@data[["treeID"]]
treeID = as.numeric(treeID)
las_tree = lasaddextrabytes(las_tree, x=treeID, name="treeID", desc="ID of a tree")
writeLAS(las_tree, savedFilename)
las_read <- readLAS(savedFilename)

Is that a lidR (rlas or liblas) error or I am doing something wrong?

Improve this question
6
  • Are you sure that the second summary of las_tree is the good one? treeID is missing into the header but it is present into las_read. I tried to reproduce your error without success. Are you sure you have an up-to-date version of rlas? LAS 1.4 is supported for not a long time
    – JRR
    Apr 15, 2019 at 16:04
  • rlas supports more properly LAS 1.4 PRF 6 since v1.3.0 (Release date: 2019-02-03)
    – JRR
    Apr 15, 2019 at 17:03
  • Thank you for your answer @JRR. I am using rlas 1.3.1 and lidR 2.0.2 I have edited the summary of las_tree. Nonetheless the classification is reset to 0 after saving. As an extra info I can say that CloudCompare doesnt see classification attribute at all on saved cloud. According to LasInfo (from Lastools) after saving min-max classification attribute values is 0. It used to be 2-5.
    – James Ch
    Apr 16, 2019 at 7:12
  • I reproduced. You found a bug. Their is something weird that happened with LAS 1.4 prf 6. Your code (or similar) works well with a LAS 1.3 file for example. It is not related to extabytes. Simply read/write a LAS 1.4 prf 6 file is enough to get Classification lost
    – JRR
    Apr 16, 2019 at 10:12
  • 1
    See issue github.com/Jean-Romain/rlas/issues/36
    – JRR
    Apr 16, 2019 at 10:24

1 Answer 1

Reset to default
4

The lidR package relies on the rlas package to read and write las file. The rlas package has a recent support of LAS 1.4 files (v1.3.0 release date: 2019-02-03). Moreover the point record formats >6 are a bit different than former point formats. Your code is correct and you actually found a bug in function write.las from rlas that occurs with point format 6 (at least). The bug you found has been fixed in rlas 1.3.2 (release date: 2019-04-26).

Otherwise your code is correct but it can be improved. The following is working and enought:

las_tree = merge_spatial(las, tif_crop, "treeID")
las_tree = add_lasattribute(las_tree, name="treeID", desc="ID of a tree")
Improve this answer
0

Your Answer

By clicking “Post Your Answer”, you agree to our terms of service, privacy policy and cookie policy

Not the answer you're looking for? Browse other questions tagged or ask your own question.