Greetings:
I am using the PBSMapping package calcVoronoi() routine to calculate
ROIs surrounding elements in
an imported ESRI point shape file. Point coords are in decimal lat/long.
here is the gist of the code.....
sites <- read.csv("ThreeCountyPoints.csv")
coordinates(sites) = c("LonDD","LatDD")
xy=coordinates(sites)
events <-
as.EventData(data.frame(EID=1:len,X=xy[,1],Y=xy[,2],projection="LL"))
Vpolys <-calcVoronoi(events)
plot.new()
plotPolys(Vpolys) # the voronois look good.
polyData <- calcArea(Vpolys)
#
# polyData's Z attribute contains numbers in the magnitude.01 - .0001
#
Two questions:
1) What are the units of polyData.Z , and how can I convert them to km**2 ?
If I try to convert the units to UTM, as suggested in the PBSMapping
user guide, as below,
attr(Vpolys,"zone") <- 9
VpolysUTM <- convUL(Vpolys)
I get the error:
--> Error in convUL(Vpolys) : Missing or invalid projection attribute.
2) Can you suggest a way in which I could clip the Voronois in the
polySet using
as a clipping polygon the county boundary polygons contained in a
polyset
(TcPoly), generated from an incoming polygon shapefile (also
lat/lon DD):
TriCounty <-read.shape("TriCountyPolyDD.shp")
TcPolys = Map2poly(TriCounty)
# here are the contents of the TriCountyPolyDD shapefile:
# GEOGCS["GCS_North_American_1983",DATUM["D_North_American_1983",
#
SPHEROID["GRS_1980",6378137.0,298.257222101]],PRIMEM["Greenwich",0.0],UNIT["Degree",0.0174532925199433]]
Please see attached jpg for the 'big picture'
Many thanks for any suggestions!
Regards,
Rick Reeves
UCSB / NCEAS