In thomasgaquiere/Maria: Simulate Different Forms of Logging on a Forest

INPUTS

data(Paracou6_2016)
data(DTMParacou)

inventory <- addtreedim(inventorycheckformat(Paracou6_2016), volumeparameters = ForestZoneVolumeParametersTable))

inventory <- treeselection(inventory, objective = 20, scenario ="manual", fuel = "2", diversification = TRUE, specieslax = FALSE,
                           objectivelax = FALSE, topography = DTMParacou, plotslope = PlotSlope, speciescriteria = SpeciesCriteria,
                           advancedloggingparameters = loggingparameters())$inventory

speciescriteria = SpeciesCriteria
scenario = "manual"
fuel = "2"
objective = 20
diversification = TRUE
directionalfelling = "2"
specieslax = FALSE
objectivelax = FALSE
topography = DTMParacou
advancedloggingparameters = loggingparameters()
# MainTrail

Create the tree polygone

# Take selected trees
SelectedTrees <- inventory %>% 
  filter(Selected == "1") %>%
  dplyr::select(idTree,DBH,TrunkHeight,TreeHeight,CrownHeight,CrownDiameter,Selected, Xutm, Yutm)

# sp::coordinates(SelectedTrees) <- ~ Xutm + Yutm
# 
# sp::proj4string(SelectedTrees) <- raster::crs(topography)
# 
# SelectedTrees <- st_as_sf(as(SelectedTrees,"SpatialPoints"))

# The Trunk
pts = list(matrix(c(SelectedTrees$Xutm[1]-(SelectedTrees$DBH[1]/100)/2, SelectedTrees$Yutm[1], 
                    SelectedTrees$Xutm[1]-(SelectedTrees$DBH[1]/100)/2, SelectedTrees$Yutm[1] + SelectedTrees$TrunkHeight[1], 
                    SelectedTrees$Xutm[1]+(SelectedTrees$DBH[1]/100)/2, SelectedTrees$Yutm[1] + SelectedTrees$TrunkHeight[1], 
                    SelectedTrees$Xutm[1]+(SelectedTrees$DBH[1]/100)/2, SelectedTrees$Yutm[1],
                    SelectedTrees$Xutm[1]-(SelectedTrees$DBH[1]/100)/2, SelectedTrees$Yutm[1]) # the return
                  ,ncol=2, byrow=TRUE)) # DBH in cm to m

(Trunk = st_polygon(pts))

# The Crown
dat <- data.frame(
  x = SelectedTrees$Xutm[1], #centroid location (x) same x than the tree
  y = SelectedTrees$Yutm[1] + SelectedTrees$TrunkHeight[1] + SelectedTrees$CrownHeight[1]/2, #centroid location (y) = trunk height + CrownHeight/2
  ex = SelectedTrees$CrownDiameter[1]/2, #Size along x-axis (CrownDiameter/2)
  ey = SelectedTrees$CrownHeight[1]/2, #Size along y-axis (CrownHeight/2)
  stringsAsFactors = FALSE
)
dat <- st_as_sf(dat, coords = c("x", "y")) # spacial dataframe

Crown <- st_ellipse(pnt = dat, ex = dat$ex, ey = dat$ey) # Create the ellipse

a <- st_difference(st_union(Crown, Trunk)) # -> multypolygon of which we keep only the points that do not overlap

ggplot() + # plot 2 polygones
  geom_sf(data = Trunk) +
  geom_sf(data = Crown) 

ggplot() +
  geom_sf(data = a) # plot the multypolygon

Treefelling orientation

(voir shéma sur google slides) 1) choper le point (Trail) de la pistele plus proche de l'emplacement de l'arbre (Foot) 2) calculer l'angle (theta) entre le vecteur que forme la hauteur de l'arbre et le vecteur opposé à celui formé entre le point "Foot" et le point "Trail", par rapport au point "Foot" 3) rotation du polygone d'un angle de 60 +/- theta °

# Create a polygon ScndTrail
pts = list(matrix(c(286520, 583134,
                    286520, 583240, 
                    286524, 583240, 
                    286524, 583134,
                    286520, 583134) # the return
                  ,ncol=2, byrow=TRUE)) # DBH en cm à passer en m

(ScndTrail <- st_polygon(pts))

ggplot() + 
  geom_sf(data = a) + # the tree
  geom_sf(data = ScndTrail) # the 2nd trail

# Find the point (Trail) on the ScndTrail closest to the location of the tree (Foot)
Foot <- st_point(c(SelectedTrees$Xutm[1],SelectedTrees$Yutm[1])) # tree foot point

NearestPoints <- st_nearest_points(Foot, ScndTrail) # from the Foot of the tree to the ScndTrail (linestring)

NearestPoint <- st_cast(NearestPoints, "POINT") # to have start (Foot) and end (Trail) points
Trail <- NearestPoint[[2]] # the point (Trail) on the ScndTrail closest to the location of the tree (Foot)

ggplot() + # plot 2 polygones
  geom_sf(data = a) + # the tree (polygone)
  geom_sf(data = Foot) + # the foot of the tree (point)
  geom_sf(data = ScndTrail) + # the 2nd trail (polygone)
  geom_sf(data = NearestPoints) + # the shortest way from the foot to the trail (linestring)
  geom_sf(data = Trail) # the point on the ScndTrail closest to the location of the tree (Foot)

# Compute theta angle
theta <- as.numeric(matlib::angle(c(Foot[1] - Foot[1], Foot[2] + SelectedTrees$TreeHeight[1]),
                                  c(Foot[1] - Trail[1], Foot[2] -Trail[2]), degree = TRUE))

Tree rotation

ref : https://github.com/jeffreyevans/spatialEco/blob/master/R/rotate.polygon.R

RotatePolygon <- function(p, angle, fixed) { # angle in the clockwise direction

  p.coords <- sf::st_coordinates(p)[,1:2] # Polygone coordinates extraction

  p.center <- suppressWarnings(sf::st_coordinates(fixed))

  rotate.coords <- function(xy, a, center) {

    co <- cos(-a * pi / 180)

    si <- sin(-a * pi / 180)

    adj <- matrix(rep(center, nrow(xy)), ncol=2, byrow=TRUE) # matrix with fixed point coordinates

    xy <- xy-adj

    cbind(co * xy[,1] - si * xy[,2],si * xy[,1] + co * xy[,2]) + adj
  }

  p.rotate <- rotate.coords(p.coords, a = angle, center = c(p.center[1], 
                                                            p.center[2]))

  Turned <- sf::st_sfc(sf::st_polygon(list(p.rotate))) # create the turned polygon

  return(Turned)

}

TurnedTrunkA <- RotatePolygon(Trunk, angle = 60 - theta, fixed = Foot) # turned trunk
TurnedCrownA <- RotatePolygon(Crown, angle = 60 - theta, fixed = Foot) # turned crown
TurnedTrunkB <- RotatePolygon(Trunk, angle = 300 - theta, fixed = Foot) # turned trunk (angle in the clockwise direction: 360-60-theta)
TurnedCrownB <- RotatePolygon(Crown, angle = 300 - theta, fixed = Foot) # turned crown (angle in the clockwise direction)

# For a random direction felling
RandomAngle <- sample(c(0:180), size = 1)
RandomTrunk <- RotatePolygon(Trunk, angle = RandomAngle, fixed = Foot) # turned trunk
RandomCrown <- RotatePolygon(Crown, angle = RandomAngle, fixed = Foot) # turned crown



ggplot() + # plot 2 polygones
  geom_sf(data = a) + # the tree (polygone)
  geom_sf(data = Foot) + # the foot of the tree (point)
  geom_sf(data = ScndTrail) + # the 2nd trail (polygone)
  geom_sf(data = NearestPoints) + # the shortest way from the foot to the trail
  geom_sf(data = Trail) + # the point on the ScndTrail closest to the location of the tree (Foot)
  geom_sf(data = TurnedTrunkA) +
  geom_sf(data = TurnedCrownA) +
  geom_sf(data = TurnedTrunkB) +
  geom_sf(data = TurnedCrownB) +
  geom_sf(data = RandomTrunk) +
  geom_sf(data = RandomCrown)