library(“lidR”,“ggplot2”,“dplyr”,“terra”)

data2014<-readLAS(“2014.las”) #Read in 2014 tile header(data2014)

data2017<-readLAS(“2017.las”) #Read in 2017 tile header(data2017) #Tells me 2017 CRS is NAD83 / New York Long Island (EPSG: 2263)

crs(data2014)<-CRS(“EPSG:6347”) #SP function to define CRS in 2014 tile.

last_returns_2017<-filter_poi(data2017, ReturnNumber == NumberOfReturns)

dem_2017<-rasterize_terrain(last_returns_2017, res = 1, algorithm = tin())

crs(dem_2017)<-crs(data2017)

dem_2017_df <- as.data.frame(dem_2017, xy = TRUE) colnames(dem_2017_df)[3] <- “elevation” plot_dem_2017<-ggplot(dem_2017_df) + geom_raster(aes(x = x, y = y, fill = elevation)) + scale_fill_viridis_c(name = “Elevation (m)”) + coord_sf(crs = crs(dem_2017)) + labs(title = “Big Egg Marsh Digital Elevation Model, 2017”) + theme_bw() print(plot_dem_2017)

last_returns_2014<-filter_poi(data2014, ReturnNumber == NumberOfReturns)

dem_2014<-rasterize_terrain(last_returns_2014, res = 1, algorithm = tin())

crs(dem_2014)<-crs(data2014) ##I think here is where I will have to reproject?? Because below this produces an error of “invalid crs”)

dem_2014_df <- as.data.frame(dem_2014, xy = TRUE) colnames(dem_2014_df)[3] <- “elevation” plot_dem_2014<-ggplot(dem_2014_df) + geom_raster(aes(x = x, y = y, fill = elevation)) + scale_fill_viridis_c(name = “Elevation (m)”) + coord_sf(crs = crs(dem_2014)) + labs(title = “Big Egg Marsh Digital Elevation Model, 2014”) + theme_bw() print(plot_dem_2014)