[R-sig-Geo] From spatial points to "spatial pictures"

[boB Rudis]

Maurizio provided his data set offline. The code + data + sample plots
are in this gist:
https://gist.github.com/hrbrmstr/e29cf9138e480db9e67d

The data is just as he described it, so it's really just a matter of
creating some polygons from the the points spline interpolation:

    library(ggplot2)
    library(ggthemes)
    library(sp)

    #' Make smooth polys from points
    #'
    #' From: http://gis.stackexchange.com/a/24929/29544
    #'
    #' @param xy n x 2 matrix with n >= k
    #' @param k # of vertices to expand points to
    #' @param ... other args passed to \code{spline()}
    spline.poly <- function(xy, vertices, k=3, ...) {

      # Wrap k vertices around each end.
      n <- dim(xy)[1]
      if (k >= 1) {
        data <- rbind(xy[(n-k+1):n,], xy, xy[1:k, ])
      } else {
        data <- xy
      }

      # Spline the x and y coordinates.
      data.spline <- spline(1:(n+2*k), data[,1], n=vertices, ...)
      x <- data.spline$x
      x1 <- data.spline$y
      x2 <- spline(1:(n+2*k), data[,2], n=vertices, ...)$y

      # Retain only the middle part.
      cbind(x1, x2)[k < x & x <= n+k, ]
    }

    # read in the data
    trees <- read.csv("CrownsExample.csv")

    # sort data by tree & point
    trees <- trees[with(trees, order(TREE, POINT)),]

    # plot just points colored by tree to make sure points look ok
    ggplot(trees, aes(x=X, y=Y, color=factor(TREE))) + coord_equal() +
geom_point()

    # small function to avoid anonymous fnction in `by` make a polygon
from our points
    make_tree_poly <- function(tree) {
      tree_pts <- spline.poly(as.matrix(tree[, c("X", "Y")]), 100)
      Polygons(list(Polygon(tree_pts)), ID=unique(tree$TREE))
    }

    # for each tree, make a set of smoothed polygons then make the whole thing
    # a SpatialPolygonsDataFrame so more spatial ops can be done on it (or more
    # data stored with it)
    tree_polys <- SpatialPolygonsDataFrame(
      SpatialPolygons(unclass(by(trees, trees$TREE, make_tree_poly)),
                      proj4string=CRS("+proj=utm +zone=32 +ellps=GRS80
+units=m +no_defs ")),
      data.frame(id=unique(trees$TREE), row.names=unique(trees$TREE)))

    # plot the new polygon tree crowns colored by tree
    tree_map <- fortify(tree_polys)
    ggplot() +
      geom_map(data=tree_map, map=tree_map,
               aes(x=long, y=lat, fill=id, map_id=id, color=id)) +
      geom_point(data=trees, aes(x=X, y=Y)) + # original points just
to show it's accurate
      coord_equal() +
      theme_map()

On Wed, Aug 5, 2015 at 9:56 AM, Maurizio Marchi
<mauriziomarchi85 at gmail.com> wrote:
> Hi everybody,
> I'm maurizio Marchi and I'm working as a post-doc at the Forestry Research
> Centre of Arezzo (Italy) and I have a question concerning spatial points
> management.
> During my research activities I have to collect data about the crown
> projection on the ground of trees. The aim is to "measure" the crown area
> and to draw a sort of picture on GIS software (QGIS in my case) to analyze
> the area and the perimeter of the crowns compared to the "free area" of a
> plot (generally a circular plot) to have an idea of how much photosynthetic
> radiation is able to reach the soil. My working steps are:
> 1) During the survey on the field each crown of each tree is initially
> "drawn" collecting 8 points that have to be connected to become a spatial
> polygon:
>
> ###---example of table to be processed
> TREE   POINT   X   Y
>     1            1       x    y
>     1            2       x    y
>     1            3       x    y
>     1            4       x    y
>     1            5       x    y
>     1            6       x    y
>     1            7       x    y
>     1            8       x    y
>     2            1       x    y
>     2            2       x    y
>     2            3       x    y
>     2            4       x    y
>     2            5       x    y
>     2            6       x    y
>     2            7       x    y
>     2            8       x    y
> [......................................]
>     N            8       x    y
>
> 2) the points,sorted by Tree become the vertices of the polygons
>
> 3) The polygons have to be processed to obtain a smoothed polygons,
> circumscribing the polygonal crowns (e.g.
> http://grasswiki.osgeo.org/wiki/V.generalize_tutorial)
>
> I generally do it on QGIS and GRASS but is there something more appropriate
> on R to do this? for example with rgdal or sp package?
>
> Many thanks
>
>
> --
> Maurizio Marchi
> Calenzano (FI) - Italy
> ID Skype: maurizioxyz
> Ubuntu 14.04 LTS
> linux user 552742
>
>         [[alternative HTML version deleted]]
>
> _______________________________________________
> R-sig-Geo mailing list
> R-sig-Geo at r-project.org
> https://stat.ethz.ch/mailman/listinfo/r-sig-geo