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