> testLines <- mapproject(yValue, xValue, proj="lambert", param=c(50,65))
For starters, if you give the x,y values in reverse order of what the
mapproject function expects you need to label them: y=yValue,
x=xValue.
(Also, I would have expected longitudes in the Americas to be
negative, but mapproject doesn't appear to care.)
Bill Dunlap
TIBCO Software
wdunlap tibco.com
On Tue, Sep 23, 2014 at 9:17 AM, Alain Dubreuil <alain.dubre...@cae.com> wrote:
> Hi all,
>
> Based on Ray's suggestions, I have tried the following script:
>
> library(mapproj)
> library(maps)
> resuPdfFileName="C:/linesTest.pdf"
> pdf(resuPdfFileName)
> # Create a map of Canada in Lambert projection
> map(database= "worldHires","Canada", ylim=c(39,90), xlim=c(-145,-25),
> col=alpha("grey90",0.5), fill=TRUE, projection="lambert", param=c(50,65))
> # Use test position vectors to draw lines
> yValue <- c(49.0, 49.0, 60.0, 60.0)
> xValue <- c(105.0, 120.0, 120.0, 105.0)
> # Convert the test vectors into lambert and then draw the lines
> testLines <- mapproject(yValue, xValue, proj="lambert", param=c(50,65))
> lines(testLines)
> dev.off()
>
> The script draws the map of Canada, but fails to draw the lines. Please let
> me know what I'm doing wrong because I can't see it. By the way, not
> specifying the lambert projection in the call to mapproject yields different
> results than specifying it which seems contrary to the documentation (?).
>
> Thanks
>
> Alain Dubreuil
> Ottawa, Canada
>
>
> -----Original Message-----
> From: Ray Brownrigg [mailto:ray.brownr...@ecs.vuw.ac.nz]
> Sent: September-23-14 4:41 AM
> To: Alain Dubreuil; r-help@r-project.org
> Subject: Re: [R] Plotting boundary lines from shapefiles overtop a map of
> Canada
>
> On 23/09/2014 3:27 a.m., Alain Dubreuil wrote:
>> Hi. I have a requirement to plot a series of points on a map of Canada
>> along with boundaries defining search and rescue (SAR) regions. I have been
>> successful in plotting the map of Canada (Lambert projection) and the
>> points, but I have been unable thus far to plot the SAR regions on top of
>> the map. I'm at the point now where I need help to resolve the issue.
>>
>> To plot the map of Canada, I have used the following line of code:
>>
>> map(database= "worldHires","Canada", ylim=c(39,90),
>> xlim=c(-150,-25), col=alpha("grey90",0.5), fill=TRUE,
>> projection="lambert", param=c(50,65))
>>
>> Note that the ylim and xlim limits go wider that the actual coordinates of
>> Canada, but that is necessary because the SAR regions go out to sea quite a
>> distance. Also, I need the map to go all the way to the North Pole.
>>
>> To plot the points, I have used a "dummy" list of points which I will
>> eventually replace with my real data. I convert the points to the lambert
>> projection on the fly using the following lines of code:
>>
>> lon <- c(-60, -60, -80, -80.1, -90, -105, -140) #a test longitude
>> vector
>> lat <- c(42.5, 42.6, 54.6, 54.4, 75, 68.3, 60) #a test latitude
>> vector
>> coord <- mapproject(lon, lat, proj="lambert", param=c(50,65))
>> #convert points to projected lat/long
>> points(coord, add=TRUE, pch=20, cex=1.2, col=alpha("red", 0.5))
>> #plot converted points
>>
>> As stated, plotting the SAR regions has not worked thus far. The best I
>> have ever gotten is a square box around the map. I have data files that
>> list the coordinates of the SAR regions, which is a succession of up to
>> 12100 lat & long points. A colleague converted those data files into
>> shapefiles defining polygons, with the coordinates already projected to
>> Lambert. I have tried various options to plot the regions, but none have
>> worked.
>>
>> Using readOGR:
>>
>> region <- readOGR(dsn="C:/myRfolder",layer="mySARshapefile")
>> plot(region, add=TRUE, xlim=c(-150,-25),ylim=c(39,90),
>> col=alpha("lightgreen", 0.6), border=TRUE)
> You don't state which package readOGR() comes from, but it is not part of the
> maps package, which doesn't understand shapefiles, so just using
> plot() on top of a (projected) map() is unlikely to succeed.
>
> I believe what you have to do is go back to your lat/long pairs for your SAR
> regions and use mapproject() to convert them to the coordinates used by the
> plotted projection. Note that you don't need the proj="lambert"
> option when you call mapproject() after a call to map() with a projection
> because the most recent projection (and its parameters) are "remembered".
> Also I suspect (though untested) is that if you put NA pairs in between your
> lists of projected SAR regions, then you can just use lines() to draw them
> all at once.
>
> Hope this helps,
> Ray Brownrigg
>> Using read.shp and draw.shp:
>>
>> region <- read.shp("C:/myRfolder/mySARshapefile.shp")
>> draw.shape(shape=region, type="poly", col="red")
>>
>> Using readShapePoly:
>>
>> region <- readShapePoly("C:/ myRfolder/mySARshapefile.shp")
>> plot(halRegion, add=TRUE, xlim=c(-150,-25),ylim=c(39,90),
>> col=alpha("lightgreen", 0.6), border=TRUE)
>>
>> Using readShapeLines after converting the region coordinates to a Lines
>> shapefile instead of a Polygon shapefile:
>>
>> region <- readShapeLines("C:/myRfolder/mySARshapefile_lines.shp")
>> lines(region, col=alpha("black", 0.6))
>>
>> I have tried playing with spplot, but I haven't quite understood how this
>> one works yet (gives me an error message: "Error in
>> stack.SpatialPointsDataFrame(as(data, "SpatialPointsDataFrame"), : all
>> factors should have identical levels")
>>
>> I would appreciate any help or insight that you could provide to help me get
>> those boundaries drawn on-top of the country map.
>>
>> Thanks
>>
>> Alain Dubreuil
>> Ottawa, Canada
>>
>
> ______________________________________________
> 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.
______________________________________________
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.
