OK, the underlying problem is that the projection code (i.e. the mapproj package plus some
parts of the maps package) does not clip its output to the specified limits.
Let me explain:
1) the maps databases are made up of line segments which are combined to form
polygons
2) when you call map() with limits specified, all line segements containing *at least one
point* within the limits is returned
3) when there is no projection, the plot() function quite naturally clips the output to
the specified rectangular limits
4) when there is a projection, the limits are no longer rectangular, and so this clipping
no longer occurs, and all points are projected
Now with your code:
1) when you use the (default) "world" database, the line segments returned include a large
part of the US-Canada border, plus all of Lake Winnipeg
2) if you use the "usa" database, you only get the US-Canada border but much more than
that part specified by your limits
3) if you use the "state" database, you get much closer to what you "want" (as defined by
the limits), but you still get some parts of line segments which are outside the limits.
In order to get what you seem to want (although if you are looking at "states", why do you
go such a long way north of the 49th parallel?), you can do something like this:
library(maps)
require("mapproj")
longlatLimit <- c(-107, -93, 40, 52)
par(plt=c(0, 1, 0, 1), cex=1, cex.main=1) # Set plotting parameters
mypoly <- list(x=c(-107, -93, -93, -107, -107), y=c(40, 40, 52, 52, 40))
mymap <- map("state", plot=F, xlim=longlatLimit[1:2], ylim=longlatLimit[3:4])
whichxy <- in.polygon(mypoly, list(x=mymap$x, y=mymap$y))
whichxy[is.na(mymap$x)] <- NA
mymap$x <- mymap$x[whichxy]
mymap$y <- mymap$y[whichxy]
map(mymap, projection="azequalarea", type="n")
bound<-c(floor(longlatLimit[1]), ceiling(longlatLimit[2]),
floor(longlatLimit[3]), ceiling(longlatLimit[4]))
map.grid(lim=bound, col="light grey")
map(mymap, projection="azequalarea", add=T)
HTH
Ray Brownrigg
Dr. Alireza Zolfaghari wrote:
Ray,
You are right, it is not having any data for states, but in fact the
grids are related to states. If you run the code you sent, you could see
there is a blank area on top and left. All I want is to locate the grid
map in the middle of window. Once I do this I can add the boundary data.
Thanks,
Alireza
On Fri, Mar 13, 2009 at 9:00 AM, Ray Brownrigg
<ray.brownr...@ecs.vuw.ac.nz <mailto:ray.brownr...@ecs.vuw.ac.nz>> wrote:
I think what you are missing is that the default map database is
"world".
If you use:
library(maps)
require("mapproj")
longlatLimit<-c(-107,-93,40,52)
par(plt=c(0,1,0,1),cex=1,cex.main=1) #Set plotting parameters
map(regions="USA", projection="azequalarea",
type="n",xlim=longlatLimit[1:2],ylim=longlatLimit[3:4])
bound<-c(floor(longlatLimit[1]), ceiling(longlatLimit[2]),
floor(longlatLimit[3]), ceiling(longlatLimit[4]))
map.grid(lim=bound,col="light grey")
I think you will get what you want (but what you want is not exactly
clear).
[I am not sure why this happens, but it seems to be an artifact of
the projection code.]
[All your examples have been very confusing, since you talk about
plotting the "states", but in fact your example code plots only a grid.]
HTH
Ray Brownrigg
Dr. Alireza Zolfaghari wrote:
Hi list,
I have a real problem with plotting US state map. When I try to
plot the
northern state, there will be some blank space in the top of
graph (see case
1 example), and when I plot southern states, there will be a
blank space in
the bottom of plot (see case 2). I spent almost 2 days to figure
out a
solution, but could not. Would you help me if you know what the
problem is?
Regards,
Alireza
#################################################
#case 1
library(maps)
require("mapproj")
longlatLimit<-c(-107,-93,40,52)
par(plt=c(0,1,0,1),cex=1,cex.main=1) #Set plotting parameters
#map(projection="azequalarea",
type="n",xlim=longlatLimit[1:2],ylim=longlatLimit[3:4])
map(projection="azequalarea",
type="n",xlim=longlatLimit[1:2],ylim=longlatLimit[3:4])
bound<-c(floor(longlatLimit[1]), ceiling(longlatLimit[2]),
floor(longlatLimit[3]), ceiling(longlatLimit[4]))
map.grid(lim=bound,col="light grey")
#################################################
#case 2
library(maps)
require("mapproj")
longlatLimit<-c(-107,-93,25,37)
par(plt=c(0,1,0,1),cex=1,cex.main=1) #Set plotting parameters
#map(projection="azequalarea",
type="n",xlim=longlatLimit[1:2],ylim=longlatLimit[3:4])
map(projection="azequalarea",
type="n",xlim=longlatLimit[1:2],ylim=longlatLimit[3:4])
bound<-c(floor(longlatLimit[1]), ceiling(longlatLimit[2]),
floor(longlatLimit[3]), ceiling(longlatLimit[4]))
map.grid(lim=bound,col="light grey")
[[alternative HTML version deleted]]
______________________________________________
R-help@r-project.org <mailto:R-help@r-project.org> mailing list
https://stat.ethz.ch/mailman/listinfo/r-help
PLEASE do read the posting guide
http://www.R-project.org/posting-guide.html
<http://www.r-project.org/posting-guide.html>
and provide commented, minimal, self-contained, reproducible code.
______________________________________________
R-help@r-project.org mailing list
https://stat.ethz.ch/mailman/listinfo/r-help
PLEASE do read the posting guide http://www.R-project.org/posting-guide.html
and provide commented, minimal, self-contained, reproducible code.
