[R] rounding of voronoi vertices using deldir()
Prof Brian Ripley
ripley at stats.ox.ac.uk
Thu Apr 6 08:35:32 CEST 2006
The problem is most likely your use of cat() for output. Consider
> x <- 8665540.49905558
> cat(x, "\n")
8665540
> cat(as.character(x), "\n")
8665540.49905558
> options(digits=10)
> cat(x, "\n")
8665540.499
So it would be best to do the conversions yourself, and I would
investigate using format() to do so. But just changing options(digits=)
may suffice.
The 'digits' option in deldir rounds the output, and with the size of
numbers you have 10dp is at or below representation accuracy.
There are alternative R packages for Voronoi polygons, e.g. tripack.
On Wed, 5 Apr 2006, Mike Leahy wrote:
> Hello list,
>
> I'm just getting started with using R - I have been trying over the past
> day or so to work out a method for generating voronoi polygons for
> PostGIS using SQL. I was able to put together a procedure which works
> relatively well, but is somewhat inefficient. Someone on the PostGIS
> list pointed me to the deldir() function in R, for which I can export a
> text file with x/y coordinates from a PostGIS table, and write an SQL
> script with insert statements containing PostGIS ewkt-compatible
> geometries representing the voronoi polygons generated by deldir().
>
> The script I'm using is below. I've added a very small frac parameter
> to ensure none of the points are dropped from my dataset (which was
> actually happening for me with some fairly dispersed clusters of
> points), and a rectangular window that is much larger than the dataset
> itself to ensure that the voronoi polygons extend far enough to actually
> cover all of the points in my dataset.
>
> The problem I'm having is that the vertices at the corners of these
> polygons seem to have their coordinates rounded to one decimal place.
> Based on the documentation, the 'digits' option should allow me to
> change this, but I'm not having getting anything different no matter
> what I set it to. Is there something obvious that I've missed? I
> realize that in most practical applications it's not a significant
> issue, but I'd prefer more accuracy as I may be using voronoi polygons
> to evaluate some statistical methods.
>
> Below are some sample coordinates from my points.txt file, and the
> script I'm running in R to write out the voronoi polygons. If anyone
> can see what's going, I'd be happy to hear it, as the R calculations are
> much faster than the method I'm using within PostGIS/SQL:
>
> points.txt:
> 266662.462042329 | 8665540.49905558
> 270171.836250559 | 8667802.6446983
> 268895.572741816 | 8674257.75469324
> 270054.378262961 | 8666483.37597101
> 268402.641255299 | 8664853.87941629
> 265707.056272354 | 8665434.09025432
> 269985.118229025 | 8667743.14071004
> 269282.034045422 | 8665403.39312076
> ....
>
>
> R
> library(deldir)
> points = scan(file="points.txt",what=list(x=0.0,y=0.0),sep="|")
> voro =
> deldir(points$x,points$y,digits=10,frac=0.000000000000001,rw=c(min(points$x)-abs(min(points$x)-max(points$x)),max(points$x)+abs(min(points$x)-max(points$x)),min(points$y)-abs(min(points$y)-max(points$y)),max(points$y)+abs(min(points$y)-max(points$y))))
> # generate voronoi edges
> tiles = tile.list(voro) # combine edges into polygons
> sink("voronoi.sql") # redirect output to file
> for (i in 1:length(tiles)) { # write out polygons
> tile = tiles[[i]]
> cat("insert into mytable (the_geom) values(geomfromtext('POLYGON((")
> for (j in 1:length(tile$x)) {
> cat (tile$x[[j]],' ',tile$y[[j]],",")
> }
> cat (tile$x[[1]],' ',tile$y[[1]]) #close polygon
> cat ("))',32718));\n") # add SRID and newline
> }
> sink() # output to terminal
> q()
> n
>
>
>
> Regards,
> Mike
>
> ______________________________________________
> R-help at stat.math.ethz.ch mailing list
> https://stat.ethz.ch/mailman/listinfo/r-help
> PLEASE do read the posting guide! http://www.R-project.org/posting-guide.html
>
--
Brian D. Ripley, ripley at stats.ox.ac.uk
Professor of Applied Statistics, http://www.stats.ox.ac.uk/~ripley/
University of Oxford, Tel: +44 1865 272861 (self)
1 South Parks Road, +44 1865 272866 (PA)
Oxford OX1 3TG, UK Fax: +44 1865 272595
More information about the R-help
mailing list