[R-sig-Geo] Shapes rearranged upon write

[Lyndon Estes]

Hi Sam,

I wrote the following function a while ago to account for this problem
of dataframe row reordering relative to the polygon IDs:

joinAttributeTable <- function(x, y, xcol, ycol) {
# Merges data frame to SpatialPolygonsDataFrame, keeping the correct
order. Code from suggestions at:
# https://stat.ethz.ch/pipermail/r-sig-geo/2008-January/003064.html
# Args:
#   x: SpatialPolygonsDataFrame
#   y: Name of data.frame to merge
#   xcol: Merge column name
#   ycol: Merge column name
# Returns: Shapefile with merged attribute table

  x$sort_id <- 1:nrow(as(x, "data.frame"))  # Column containing
original row order for later sorting
	
  x.dat <- as(x, "data.frame")  # Create new data.frame object
  x.dat2 <- merge(x.dat, y, by.x = xcol, by.y = ycol)  # Merge
  x.dat2.ord <- x.dat2[order(x.dat2$sort_id), ]  # Reorder back to original
	x2 <- x[x$sort_id %in% x.dat2$sort_id, ]  # Make new set of polygons,
dropping those which aren't in merge
  x2.dat <- as(x2, "data.frame")  # Make updated x2 into a data.frame
	row.names(x.dat2.ord) <- row.names(x2.dat)  # Reassign row.names from
original data.frame
  x2 at data <- x.dat2.ord  # Assign to shapefile the new data.frame
  return(x2)
}

I haven't used this in a while, so I hope it still works.  I think
others could suggest a more efficient, existing method for dealing
with this as well.

Cheers, Lyndon


On Wed, Nov 2, 2011 at 12:12 PM, Sam Rabin <srabin at princeton.edu> wrote:
> Hello!
>
> I have a shapefile (.shp) containing the Landsat grid coverage of South America. I also have 97 other shapefiles, each of which covers one of the Landsat grid cells (the 97 being a subset of the cells in the first shapefile). I managed to get R (via maptools) to perform a calculation on each of the 97 and put the results in a data frame. Among other things, this data frame contains a column for grid cell ID that matches up with a column in the Landsat grid shapefile. I used merge() to combine these and then saved a new shapefile, following the instructions at http://help.nceas.ucsb.edu/R:_Spatial#Append_a_second_set_of_attributes_to_a_spatial_data_file.27s_attribute_table. Here is the code I used:
>
>        require(maptools)
>        landsat1 = readShapePoly("landsat_grid")
>          # Make a backup
>          landsat_orig = landsat1
>        landsat_2_spdf = merge(landsat1 at data,newcalc,by.x="PR",by.y="pr2",all.x=TRUE,sort=FALSE)
>        landsat1 at data = landsat_2_spdf
>        writeSpatialShape(landsat1, "methodA")
>
> Looking at the new data with
>
>        View(landsat1 at data),
>
> it seemed like everything went great. The ID's matched up perfectly, and I checked some of the calculated figures and they matched, too. The polygons from landsat_grid.shp that were not included in the 97 were assigned NA for the column with the calculated figures, appropriately.
>
> Unfortunately, when I imported the new shapefile to my GIS program, all the tiles were rearranged. screenshot_expectedcoverage (http://tinyurl.com/3nbm3q7 ) shows fine-scale data for what I was expecting (Amazonia; imagine a polygon drawn around all the data there), and screenshot_methodA (http://tinyurl.com/3htvlga) shows what I got. (Only the 97 polygons are shown.)
>
> I tried changing sort to TRUE in the fourth line (and the output shapefile name to "methodB"), but the tiles are still messed up, just differently — see screenshot_methodB (http://tinyurl.com/3h95vr3).
>
> In both shape files, the attribute tables look fine, but the tiles are just totally rearranged. They're not all scattered over the world, either. In screenshot_rearrangement (which shows all tiles, not just the 97) (http://tinyurl.com/6xfabvm), you can see that they cover South America perfectly — they all line up with the polygons from landsat_grid.shp (not shown).
>
> What do I need to do to fix this?
>
> Thanks very much in advance.
>
>
>
> Sam Rabin
>     Graduate student
>         Princeton University
>         Ecology & Evolutionary Biology
>
>
>
>
>
>
>        [[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
>
>