add RGB values

[wiesehahn]

Hi, I want to merge RGB values from raster image to the las pointcloud.

Currently I tried:

las_rgb <- merge_spatial(las, img[[1]], attribute = "r")
las_rgb <- merge_spatial(las_rgb, img[[2]], attribute = "g")
las_rgb <- merge_spatial(las_rgb, img[[3]], attribute = "b")

las_rgb <- add_lasrgb(las_rgb, as.integer(las_rgb$r), as.integer(las_rgb$g), as.integer(las_rgb$b))

But I guess there is a better way, since the intermediate r value seems obsolete. How would I do it properly?

Also it seems strange that plot(las_rgb, color = "RGB) cant find RGB values when I do not run add_lasrgb (since I named them with lower letters), but filter_poi seems to recognize the lower letters as RGB values and renames them to capital RGB, which then results in two attributes for each R, G and B.

las_rgb <- merge_spatial(las, img[[1]], attribute = "r")
las_rgb <- merge_spatial(las_rgb, img[[2]], attribute = "g")
las_rgb <- merge_spatial(las_rgb, img[[3]], attribute = "b")
plot(las_rgb, color = "RGB")
# Error: No 'RGB' attributes found.
las_rgb <- add_lasrgb(las_rgb, as.integer(las_rgb$r), as.integer(las_rgb$g), as.integer(las_rgb$b))
plot(las_rgb, color = "RGB")
# Error in if (colmax > 255) 16 else 8 : 
#  missing value where TRUE/FALSE needed
las_rgb <- filter_poi(las_rgb, !is.na(R) & !is.na(G) & !is.na(B))
# Attribute 'r' renamed 'R' to match with default attribute names.
# Attribute 'g' renamed 'G' to match with default attribute names.
# Attribute 'b' renamed 'B' to match with default attribute names.
str(las_rgb@data)
# Classes ‘data.table’ and 'data.frame':    336191 obs. of  23 variables:
#  $ X                : num  567275 567275 567274 567274 567274 ...
#  $ Y                : num  5708885 5708885 5708885 5708885 5708885 ...
#  $ Z                : num  0 0.998 0 1.313 -0.025 ...
#  $ gpstime          : num  567074 567074 567074 567074 567074 ...
#  $ Intensity        : int  176 20 136 26 149 117 12 12 113 23 ...
#  $ ReturnNumber     : int  5 5 6 4 5 3 3 4 5 2 ...
#  $ NumberOfReturns  : int  5 6 6 5 5 3 3 5 5 4 ...
#  $ ScanDirectionFlag: int  0 0 0 0 0 0 0 0 0 0 ...
#  $ EdgeOfFlightline : int  0 0 0 0 0 0 0 0 0 0 ...
#  $ Classification   : int  2 13 2 13 20 20 13 13 20 13 ...
#  $ Synthetic_flag   : logi  FALSE FALSE FALSE FALSE FALSE FALSE ...
#  $ Keypoint_flag    : logi  FALSE FALSE FALSE FALSE FALSE FALSE ...
#  $ Withheld_flag    : logi  FALSE FALSE FALSE FALSE FALSE FALSE ...
#  $ ScanAngleRank    : int  30 30 30 30 30 30 30 30 30 30 ...
#  $ UserData         : int  50 47 51 53 53 49 122 61 55 54 ...
#  $ PointSourceID    : int  150 150 150 150 150 150 150 150 150 150 ...
#  $ Zref             : num  322 323 322 323 322 ...
#  $ R                : int  104 104 81 108 88 84 108 92 68 81 ...
#  $ G                : int  96 96 73 100 80 72 100 84 55 73 ...
#  $ B                : int  85 85 60 87 67 60 87 71 46 60 ...
#  $ R                : int  26728 26728 20817 27756 22616 21588 27756 23644 17476 20817 ...
#  $ G                : int  24672 24672 18761 25700 20560 18504 25700 21588 14135 18761 ...
#  $ B                : int  21845 21845 15420 22359 17219 15420 22359 18247 11822 15420 ...
#  - attr(*, ".internal.selfref")=<externalptr>

[Jean-Romain]

From the doc:

RasterStack, RasterBrick, multibands stars or multilayer SpatRaster must have 3 layers for RGB colors. It colorizes the point cloud with RGB values.

So this should work

las_rgb <- merge_spatial(las, img) 

Yet you found a bug with automatic naming.

[wiesehahn]

Ok, but I still have a problem, with my above mentioned aproach I get one point with NA value for whatever reason (the raster seems to be fine).

> is.na(img)
class       : SpatRaster 
dimensions  : 500, 500, 3  (nrow, ncol, nlyr)
resolution  : 0.2, 0.2  (x, y)
extent      : 567222.4, 567322.4, 5708785, 5708885  (xmin, xmax, ymin, ymax)
coord. ref. : ETRS89 / UTM zone 32N (EPSG:25832) 
source(s)   : memory
names       :     R,     G,     B 
min values  : FALSE, FALSE, FALSE 
max values  : FALSE, FALSE, FALSE 

Thats why plot() did not work and I used las_rgb <- filter_poi(las_rgb, !is.na(R) & !is.na(G) & !is.na(B)) to filter this point.

With las_rgb <- merge_spatial(las, img) I cannot get around this and get Error: Some points were associated with an NA RGB color. RGB cannot be NA in a LAS object. Colorization aborted.

[Jean-Romain]

The img either does not cover the point cloud entirely or contain NA. If not please provide a reproducible example.

[wiesehahn]

You were right I was trying to crop the image with the bbox of the las, which did not perfectly fit with default sanpping option nearest, changing to out did the trick.

> bbox <- st_bbox(las)
> 
> 
> # crop to las extent
> img1 <- crop(img, st_as_sfc(bbox))
> 
> ext(img1)
SpatExtent : 569137.4, 569147.4, 5709415.4, 5709425.2 (xmin, xmax, ymin, ymax)
> ext(las)
SpatExtent : 569137.482, 569147.434, 5709415.464, 5709425.295 (xmin, xmax, ymin, ymax)
> 
> # try to add RGB
> las_rgb <- merge_spatial(sample, img1) 
Error: Some points were associated with an NA RGB color. RGB cannot be NA in a LAS object. Colorization aborted.
> # crop to extented bbox
> img2 <- crop(img, st_as_sfc(bbox),snap="out")
> 
> ext(img2)
SpatExtent : 569137.4, 569147.6, 5709415.4, 5709425.4 (xmin, xmax, ymin, ymax)
> ext(las)
SpatExtent : 569137.482, 569147.434, 5709415.464, 5709425.295 (xmin, xmax, ymin, ymax)
> 
> # add RGB works
> las_rgb <- merge_spatial(sample, img2) 

Although I am wondering if it would be possible / useful to add RGB values where there is valid data only (not sure if NA is valid in las data though).

[Jean-Romain]

No because NA is not a writable value in a LAS file and lidR tries to ensure that LAS objects are always writable in LAS files.

[wiesehahn]

Makes sense, thanks!

[Jean-Romain]

Fixed, attributes are no longer automatically renamed when filtering a point cloud