1### R code from vignette source 'combine_maptools.Rnw' 2 3################################################### 4### code chunk number 1: combine_maptools.Rnw:44-47 5################################################### 6owidth <- getOption("width") 7options("width"=90) 8.PngNo <- 0 9 10 11################################################### 12### code chunk number 2: afig (eval = FALSE) 13################################################### 14## .PngNo <- .PngNo + 1; file <- paste("Fig-bitmap-", .PngNo, ".pdf", sep="") 15## pdf(file=file, width = 6.5, height = 3.5, pointsize = 12, bg = "white") 16## opar <- par(mar=c(3,3,1,1)+0.1) 17 18 19################################################### 20### code chunk number 3: afigl (eval = FALSE) 21################################################### 22## .PngNo <- .PngNo + 1; file <- paste("Fig-bitmap-", .PngNo, ".pdf", sep="") 23## pdf(file=file, width = 6.5, height = 3.5, pointsize = 12, bg = "white") 24 25 26################################################### 27### code chunk number 4: bfigl (eval = FALSE) 28################################################### 29## .PngNo <- .PngNo + 1; file <- paste("Fig-bitmap-", .PngNo, ".pdf", sep="") 30## pdf(file=file, width = 6.5, height = 5, pointsize = 12, bg = "white") 31 32 33################################################### 34### code chunk number 5: bfig (eval = FALSE) 35################################################### 36## .PngNo <- .PngNo + 1; file <- paste("Fig-bitmap-", .PngNo, ".pdf", sep="") 37## pdf(file=file, width = 6.5, height = 5, pointsize = 12, bg = "white") 38## opar <- par(mar=c(3,3,1,1)+0.1) 39 40 41################################################### 42### code chunk number 6: zfig (eval = FALSE) 43################################################### 44## par(opar) 45## dev.null <- dev.off() 46## cat("\\includegraphics[width=0.95\\textwidth]{", file, "}\n\n", sep="") 47 48 49################################################### 50### code chunk number 7: zfigl (eval = FALSE) 51################################################### 52## dev.null <- dev.off() 53## cat("\\includegraphics[width=0.95\\textwidth]{", file, "}\n\n", sep="") 54 55 56################################################### 57### code chunk number 8: combine_maptools.Rnw:101-111 58################################################### 59owd <- getwd() 60setwd(system.file("shapes", package="maptools")) 61library(maptools) 62nc90 <- readShapeSpatial("co37_d90") 63proj4string(nc90) <- CRS("+proj=longlat +datum=NAD27") 64sc90 <- readShapeSpatial("co45_d90") 65proj4string(sc90) <- CRS("+proj=longlat +datum=NAD27") 66va90 <- readShapeSpatial("co51_d90") 67proj4string(va90) <- CRS("+proj=longlat +datum=NAD27") 68setwd(owd) 69 70 71################################################### 72### code chunk number 9: combine_maptools.Rnw:117-124 73################################################### 74.PngNo <- .PngNo + 1; file <- paste("Fig-bitmap-", .PngNo, ".pdf", sep="") 75pdf(file=file, width = 6.5, height = 5, pointsize = 12, bg = "white") 76opar <- par(mar=c(3,3,1,1)+0.1) 77oopar <- par(mar=c(3,2,1,1)+0.1) 78plot(va90, xlim=c(-85,-75), ylim=c(32,40), axes=TRUE, border="grey10") 79plot(nc90, add=TRUE, border="grey40") 80plot(sc90, add=TRUE, border="grey70") 81par(oopar) 82par(opar) 83dev.null <- dev.off() 84cat("\\includegraphics[width=0.95\\textwidth]{", file, "}\n\n", sep="") 85 86 87################################################### 88### code chunk number 10: combine_maptools.Rnw:148-149 89################################################### 90library(maptools) 91 92 93################################################### 94### code chunk number 11: combine_maptools.Rnw:151-155 95################################################### 96names(sc90) 97sc90a <- spChFIDs(sc90, paste(sc90$ST, sc90$CO, sep="")) 98sc90a <- sc90a[,-(1:4)] 99names(sc90a) 100 101 102################################################### 103### code chunk number 12: combine_maptools.Rnw:157-158 104################################################### 105proj4string(sc90a) <- CRS(proj4string(sc90a)) 106 107 108################################################### 109### code chunk number 13: combine_maptools.Rnw:169-170 110################################################### 111names(nc90) 112 113 114################################################### 115### code chunk number 14: combine_maptools.Rnw:172-173 (eval = FALSE) 116################################################### 117## nc90a <- spChFIDs(nc90, paste(nc90$ST, nc90$CO, sep="")) 118 119 120################################################### 121### code chunk number 15: combine_maptools.Rnw:175-176 122################################################### 123try1 <- try(spChFIDs(nc90, paste(nc90$ST, nc90$CO, sep="")), silent=TRUE) 124 125 126################################################### 127### code chunk number 16: combine_maptools.Rnw:178-179 128################################################### 129cat(try1) 130 131 132################################################### 133### code chunk number 17: combine_maptools.Rnw:189-190 134################################################### 135table(table(paste(nc90$ST, nc90$CO, sep=""))) 136 137 138################################################### 139### code chunk number 18: combine_maptools.Rnw:216-219 140################################################### 141if (rgeosStatus()) { 142nc90a <- unionSpatialPolygons(nc90, IDs=paste(nc90$ST, nc90$CO, sep="")) 143} 144 145 146################################################### 147### code chunk number 19: combine_maptools.Rnw:229-234 148################################################### 149if (rgeosStatus()) { 150nc90_df <- as(nc90, "data.frame")[!duplicated(nc90$CO),-(1:4)] 151row.names(nc90_df) <- paste(nc90_df$ST, nc90_df$CO, sep="") 152nc90b <- SpatialPolygonsDataFrame(nc90a, nc90_df) 153} 154 155 156################################################### 157### code chunk number 20: combine_maptools.Rnw:254-263 158################################################### 159if (rgeosStatus()) { 160va90a <- spChFIDs(va90, paste(va90$ST, va90$CO, sep="")) 161va90a <- va90a[,-(1:4)] 162va90_pl <- slot(va90a, "polygons") 163va90_pla <- lapply(va90_pl, checkPolygonsHoles) 164p4sva <- CRS(proj4string(va90a)) 165vaSP <- SpatialPolygons(va90_pla, proj4string=p4sva) 166va90b <- SpatialPolygonsDataFrame(vaSP, data=as(va90a, "data.frame")) 167} 168 169 170################################################### 171### code chunk number 21: combine_maptools.Rnw:303-309 172################################################### 173if (rgeosStatus()) { 174nc_sc_va90 <- spRbind(spRbind(nc90b, sc90a), va90b) 175FIPS <- row.names(nc_sc_va90) 176str(FIPS) 177length(slot(nc_sc_va90, "polygons")) 178} 179 180 181################################################### 182### code chunk number 22: combine_maptools.Rnw:338-340 183################################################### 184t1 <- read.fwf(system.file("share/90mfips.txt", package="maptools"), skip=21, 185 widths=c(4,4,4,4,2,6,2,3,3,1,7,5,3,51), colClasses = "character") 186 187 188################################################### 189### code chunk number 23: combine_maptools.Rnw:342-348 190################################################### 191t2 <- t1[1:2004,c(1,7,8,14)] 192t3 <- t2[complete.cases(t2),] 193cnty1 <- t3[t3$V7 != " ",] 194ma1 <- t3[t3$V7 == " ",c(1,4)] 195cnty2 <- cnty1[which(!is.na(match(cnty1$V7, c("37", "45", "51")))),] 196cnty2$FIPS <- paste(cnty2$V7, cnty2$V8, sep="") 197 198 199################################################### 200### code chunk number 24: combine_maptools.Rnw:368-375 201################################################### 202if (rgeosStatus()) { 203MA_FIPS <- cnty2$V1[match(FIPS, cnty2$FIPS)] 204MA <- ma1$V14[match(MA_FIPS, ma1$V1)] 205MA_df <- data.frame(MA_FIPS=MA_FIPS, MA=MA, row.names=FIPS) 206nc_sc_va90a <- spCbind(nc_sc_va90, MA_df) 207ncscva_MA <- unionSpatialPolygons(nc_sc_va90a, nc_sc_va90a$MA_FIPS) 208} 209 210 211################################################### 212### code chunk number 25: combine_maptools.Rnw:381-392 213################################################### 214.PngNo <- .PngNo + 1; file <- paste("Fig-bitmap-", .PngNo, ".pdf", sep="") 215pdf(file=file, width = 6.5, height = 5, pointsize = 12, bg = "white") 216opar <- par(mar=c(3,3,1,1)+0.1) 217if (rgeosStatus()) { 218oopar <- par(mar=c(3,2,1,1)+0.1) 219plot(nc_sc_va90, border="grey", axes=TRUE) 220plot(ncscva_MA, lwd=2, add=TRUE) 221text(coordinates(ncscva_MA), labels=row.names(ncscva_MA), cex=0.6) 222par(oopar) 223} else { 224plot(1) 225} 226par(opar) 227dev.null <- dev.off() 228cat("\\includegraphics[width=0.95\\textwidth]{", file, "}\n\n", sep="") 229 230 231################################################### 232### code chunk number 26: combine_maptools.Rnw:410-417 233################################################### 234if (rgeosStatus()) { 235np <- sapply(slot(ncscva_MA, "polygons"), function(x) length(slot(x, "Polygons"))) 236table(np) 237MA_fips <- row.names(ncscva_MA) 238MA_name <- ma1$V14[match(MA_fips, ma1$V1)] 239data.frame(MA_fips, MA_name)[np > 1,] 240} 241 242 243################################################### 244### code chunk number 27: combine_maptools.Rnw:425-426 245################################################### 246options("width"=owidth) 247 248 249