This is awesome help! Thank you very much Tom! I will give it a go. Mike On Jul 17, 2006, at 10:04 PM, Mulholland, Tom wrote:
> Yes you are correct in your assumption that if you are not using > the "maps" package you need to have your own boundary files. There > is the "National Atlas of the United States" produced by the USGS > which I recall has downloads of a variety of shapefiles. http://www- > atlas.usgs.gov/ is the mainpage but administrative boundaries are > on http://www-atlas.usgs.gov/atlasftp.html? > openChapters=chpbound#chpbound in compressed files. I'm sure I have > previously found them uncompressed, but it's not too hard to > decompress tham as you only have to do it once. > > > So once you have the shapefiles you can read them using the > maptools package with the readShapePoly set of commands. > > such as > > US <- readShapePoly("name of shape file") > #make sure you have the correct transformation > # set the CRS for the shapefile you have downloaded. > # You need to check the metadata associated with your shapefile > > # I typically use sommething like > # proj4string(x) <- CRS("+proj=longlat +datum=NAD83") > # Then reproject it into the format you wnat to use > # x <- transform(x,CRS("+proj=utm +zone=50")) > # } > > # Having done that > > plot(US,...) > > #Then you add your data to the map > plot(S, pch = 20, add = TRUE) > > Yes you can extract the coordinates from maps > >> test <- map("world","UK",xlim = c(-8.08932,3.4), ylim = c >> (49.11741,61.54963)) >> str(test) > List of 4 > $ x : num [1:296] -6.24 -6.30 -6.23 -6.38 -6.46 ... > $ y : num [1:296] 58.5 58.3 58.2 58.2 58.1 ... > $ range: num [1:4] -8.09 1.73 49.25 60.84 > $ names: chr [1:25] "UK:Scotland:Isle of Lewis" "UK:Guernsey" > "UK:Great Britain" "UK:Scotland:Shetland Islands:Unst" ... > - attr(*, "class")= chr "map" > > but it looks like you would them have to process them into an sp > format. However the world map is not particularly detailed. It's > possible the US county data is of more precision. > > If you just use maps then you should be able to use "points" to > plot your data > > map('usa') > points(S$x, S$y, pch = 20) > > Not having delved into the map package I don't know exactly how it > handles projections but it looks like > > map('usa', projection = "yourCRS") > > would be sufficient. > > I hope this helps. > > Tom Mulholland > Senior Demographer > Applied Research and Modelling > State and Regional Policy > Department for Planning and Infrastructure > Perth, Western Australia > +61 (08) 9264 7936 > > >> >> -----Original Message----- >> From: [EMAIL PROTECTED] >> [mailto:[EMAIL PROTECTED] Behalf Of >> Michael Fuller >> Sent: Tuesday, 18 July 2006 5:56 AM >> To: r-sig-geo@stat.math.ethz.ch >> Subject: [R-sig-Geo] geo-newbie Q: overlay lat-long data onto >> projected >> map? >> >> >> I've spent the last 5 days trying to figure out how to overlay the >> positions of geographic data on top of a map. I've scoured the man >> pages for the maps, mapproj, maptools, mapdata, and the huge sp >> package. I've experimented with countless permutations of commands. >> So far I've learned a lot about map projections, but I still have no >> ideas how to display my data the way I want! (sigh) >> >> The data are the latitude-longitude coordinates of 1,599 >> survey sites >> located across the US. I want to plot these sites on a map of >> the US. >> Seems like a simple task. I've tried converting them to a >> SpatialPoints object and worked out how to project the sites using a >> given map projection (i.e. Bonne). But how do I correctly overlay >> (i.e. display) the data points onto a similarly projected map of the >> US? It might be possible to do it using overlay {sp} but for that it >> seems I would need the coordinates of the geopolitical boundary for >> the map of the US. Can these values be extracted somehow from the US >> map? >> >> Here's how I generated the SpatialPoints object: >> >> #x = longitude values in decimal degrees >> #y = latitude values in decimal degrees >> #create SpatialPoints object using CRS-class Bonne projection >> centered on longitude 100 >> S <- SpatialPoints(cbind(x,y),CRS("+proj=bonne +lon_0=100")) >> >> #this next part plots just the data points (cloud of points) >> plot(S,pch=20) >> >> How do I overlay the same points on top of a map of the US such that >> they are correctly positioned wrt coordinates? >> >> Also, the R-sig-geo info page does show an URL for a searchable >> archive of R-sig-geo. Does one exist? >> >> Many thanks for your time. >> Mike >> ______________________________ >> Michael M Fuller, Ph.D. >> The Institute for Environmental Modeling >> University of Tennessee >> 1416 Circle Drive >> Knoxville, TN 37996-1610 >> >> >> >> >> [[alternative HTML version deleted]] >> >> _______________________________________________ >> R-sig-Geo mailing list >> R-sig-Geo@stat.math.ethz.ch >> https://stat.ethz.ch/mailman/listinfo/r-sig-geo >> _______________________________________________ R-sig-Geo mailing list R-sig-Geo@stat.math.ethz.ch https://stat.ethz.ch/mailman/listinfo/r-sig-geo