Tree ring detection and ring width measurement from density profiles

Introduction

Load the library, an example density profile and then trim it and detrend it. While ring detection will work on untrimmed and undetrended density profiles, you will obtain better results by trimming and detrending the density profile.

library(densitr)
dp <- dpload(dp.file = system.file("extdata", "00010001.dpa", package = "densitr"))
dp.trimmed <- dptrim(dp)
dp.detrended  <- dpdetrend(dp.trimmed, type = "gam")

Detecting tree rings

Tree rings can identified using dprings function, which will identify local peaks and valleys inside the density profile. It will return a data frame by default, where the value is the horizontal index of the peak/valley:

rings <- dprings(dp.detrended)
head(rings)
#>    value   type amplitude
#> 1     19   peak  438.7755
#> 37   170 valley  400.5265
#> 2    248   peak  440.3054
#> 38   436 valley  407.3829
#> 3    506   peak  484.3507
#> 39   613 valley  417.8074

When dprings is called with return.plot = TRUE, the function will return a diagnostic plot. The green points are valleys, blue points are peaks and red points were either a repetition of peak or valley that was automatically excluded.

dprings(dp.detrended, return.plot = TRUE)

The dashed line on the graph is the minimum peak value limit, represented by the overall mean of the profile. You could also remove some noise by smoothing the profile first by adding the argument smooth = TRUE. This applies a LOESS regression with a given span and removes some of the noise.

dprings(dp.detrended, smooth = TRUE, return.plot = TRUE)

Different tree species will require different parameters, try to adjust either pps or threshold.sd to get better results. Best results are expected in softwoods, as they have a clearer transition between early and late wood.

Extracting ring widths

To extract ring widths from the ring detection, use get_RW on the data frame returned by dprings, which will a vector of peak-to-peak distances representing ring widths. The unit is dictated by the unit of the original density profile, which can be extracted using dp$footer$xUnit.

get_RW(rings)
#>  [1] 229 258 175 633 163 191 290 156 455 226 254 329 135 193 287 257 384 277 201
#> [20] 298 205 276 241 167 263 290 239 146 323 267 266 315 259
dp$footer$xUnit
#> [1] "1/100 mm"
length(get_RW(rings))
#> [1] 33

Overall, in this particular density profile we have detected 33 rings, the values can then be further examined as demonstrated below.

rw <- get_RW(rings)
## convert to milimetres
rw <- rw / 100
summary(rw)
#>    Min. 1st Qu.  Median    Mean 3rd Qu.    Max. 
#>   1.350   2.010   2.580   2.621   2.900   6.330