[R-sig-Geo] Seeking example of making polygons from non-uniform grid

Waichler, Scott R Scott.Waichler at pnl.gov
Tue Jun 27 02:45:13 CEST 2006


> Does this help? 
> 
> ## bottom left corner of bottom left cell llX <- 0 llY <- 0.5
> 
> ## cell width and height
> cellX <- 1
> cellY <- 0.3

Michael, thanks for the example, but does this not assume that cellX is
constant across the grid, and likewise for cellY?

After my post I found a Roger Bivand example of how to process a few
irregular polygons given knowledge of the vertices.  I adapted this for
a rectangular grid situation, where the cells are varying sizes.  My
code to do this is below.  If anyone desires to make this more
efficient, please have at it.

  library(sp)
  library(spgpc)

  dx <- c(1, 2, 2, 1, 3, 3, 1, 1)
  dz <- c(2, 1, 1, 3, 3, 1, 2, 2)
  x.origin <- 0
  xnb <- c(x.origin, x.origin + cumsum(dx))  # x cell boundaries
  z.origin <- 0
  znb <- c(z.origin, z.origin + cumsum(dz))  # z nodal boundaries
  num.x <- length(dx)
  num.z <- length(dz)
  num.cells <- num.x * num.z

  k <- 0
  for(j in 1:num.z) {
    for(i in 1:num.x) {
      k <- k + 1
      ll.x <- xnb[i]  # lower left corner
      ll.y <- znb[j]
      ul.x <- ll.x     # upper left corner
      ul.y <- znb[(j+1)]
      ur.x <- xnb[(i+1)]  # upper right corner
      ur.y <- ul.y  
      lr.x <- ur.x  # lower right corner
      lr.y <- znb[j]
      this.polygon <-
        Polygon(matrix(c(ll.x, ll.y, ul.x, ul.y, ur.x, ur.y, lr.x, lr.y,
ll.x, ll.y), byrow=T, ncol=2),
                hole=F)

      cmd <- paste(sep="", 'Srs', k, ' <- Polygons(list(this.polygon),
"s', k, '")')
      eval(parse(text=cmd))
    }
  }
  
  # Assign values to the cells
  cell.values <- c(1, 1, 2, 1, 3, 3, 3, 3, 2, 2, 2, 1, 1, 2, 2, 2, 2, 2,
2, 3, 3,
                   2, 2, 3, 3, 3, 1, 1, 1, 3, 3, 3, 3, 3, 1, 1, 1, 3, 2,
3, 2, 2,
                   3, 3, 3, 3, 1, 2, 3, 3, 2, 2, 2, 2, 1, 3, 1, 1, 2, 2,
2, 2, 1, 3)

  cmd.SR <- 'SR <- SpatialPolygons(list(Srs1'
  for(k in 2:num.cells) cmd.SR <- paste(sep="", cmd.SR, ', Srs', k)
  cmd.SR <- paste(sep="", cmd.SR, ')', ', cell.values', ')')
  eval(parse(text=cmd.SR))
  SR.new <- lapply(getSpPpolygonsSlot(SR), checkPolygonsHoles)

  # Plot the individual cells/polygons
  plot(as.SpatialPolygons.PolygonsList(SR.new), col = cell.values, pbg =
"white")

  # Now union contiguous like-valued cells into minimum number of
polygons.
  SR.union <-
unionSpatialPolygons(as.SpatialPolygons.PolygonsList(SR.new),
cell.values)
  plot(SR.union, col=1:3)

Thanks,
Scott Waichler
Pacific Northwest National Laboratory
scott.waichler _at_ pnl.gov




More information about the R-sig-Geo mailing list