Update 2 - Calculating Orangic Matter Inventory

Calculated a new difference layer as I had realized my 2014 DSM was using meters instead of feet, as my placeholder crs I had set for the 2014 LAS was using meters. Then, remodeled the y = 0.15x +0.12 for organic matter as a relation to accretion rate. into x = 6.6667y - 0.8 to find organic matter inventory in gram per cm^2. Initially had to convert my difference layer to centimeters by multiplying by 30.48 then divided by 3 to get my annual accretion rate in cm/yr. My output results were different from the ones in the article by a large margin. This could be due to a combination of factors including the way each year’s DSM was interpolated. For one, 2017 had more missing ground points in the point cloud than 2014, which may be why 2017’s DSM had more spaces of higher Z value than 2014.

t2_tin_2014 <- raster("t2_2014_tin_cropped.tif")
t2_tin_2017 <-raster("t2_2017_tin.tif")
t2_tin_diff <- raster("t2_tin_depo.tif")
#converting 2014 DSM to feet
t2_tin_2014_ft <- t2_tin_2014 * 3.28
writeRaster(t2_tin_2014_ft, 't2_tin_2014_ft.tif', overwrite=TRUE)

t2_tin_diff2 <- t2_tin_2017 - t2_tin_2014_ft
plot(t2_tin_diff2, main='accretion difference between 2017-2014(ft)')

writeRaster(t2_tin_diff2,'t2_tin_diff2.tif', overwrite=TRUE)
# change to cm per year to match article metric
t2_tin_diff2_yr <- (t2_tin_diff2 / 3) * 30.48
plot(t2_tin_diff2_yr, main='Annual change in cm')

writeRaster(t2_tin_diff2_yr,'t2_tin_diff2_yr.tif',overwrite=TRUE)

hist(t2_tin_diff2_yr,
     breaks = 6,
     main = 'annual accretion change 2014-2017',
     col='blue', 
     xlab='elevation change(cm)')

OMI2 <- 6.6667 * t2_tin_diff2_yr - 0.8

writeRaster(OMI2,'OMI.tif',overwrite=TRUE)

mapviewOptions(basemaps = c('Esri.WorldImagery', 'OpenStreetMap'))
mapview(OMI2)
hist(OMI2,
     breaks = 6,
     main = 'Calculated Organic Matter Inventory(g cm^2)',
     col='green', 
     xlab='OMI(g/cm^2')