Hi there, I want to integrate a function over an irregular polygon. Is there any function which can implement this easily? Otherwise, I am thinking of divide the polygon into very small rectangles and use "adapt" to approximate it. Do you have any suggestions to get the fine division? Any advice is appreciated. Haiyong
Haiyong, There may be better ways, but this what I'd do. (And I'm not an expert on this.) (a) surround the polygon with a rectangle, (b) define, via an indicator function, a new function that is equal to your desired function within the polygon, and zero outside it, (c) use adapt() to integrate the new function over the whole rectangle. The tricky part is (b). How difficult this is depends on how complicated the polygon is. If it's convex then it can be represented by a set of inequalities Ax >= 0 and Bx <= 0. Ted. Haiyong Xu wrote on 02/15/2007 01:06 PM:> Hi there, > > I want to integrate a function over an irregular polygon. Is there > any function which can implement this easily? Otherwise, I am > thinking of divide the polygon into very small rectangles and use > "adapt" to approximate it. Do you have any suggestions to get the > fine division? Any advice is appreciated. > > Haiyong > > ______________________________________________ > 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 > and provide commented, minimal, self-contained, reproducible code. >-- Dr E.A. Catchpole Visiting Fellow Univ of New South Wales at ADFA, Canberra, Australia _ and University of Kent, Canterbury, England 'v' - www.pems.adfa.edu.au/~ecatchpole / \ - fax: +61 2 6268 8786 m m - ph: +61 2 6268 8895
I'm still pretty ignorant about R, but I think it might be possible to
work out an algorithm using cross products. First you would want to
subdivide the polygon into convex polygons. I haven't tried to do
that before, but it looks like it might be possible by looking at the
sign of cross products of vectors between vertices. (In other words,
pick a vertex, and then start working your way around the polygon and
pay attention to the sign of cross products of vectors from the
starting vertex to successive vertices.) Once you have convex
polygons, you can calculate the area using cross products of vectors
from some point (e.g. the origin) to adjacent vertices of the polygon.
I think that probably most computer graphics texts would have such an
algorithm.
How to implement that in R is not something I can answer, but it
doesn't sound hard.
--Paul
On 2/14/07, Haiyong Xu <xuhy at ucla.edu> wrote:> Hi there,
>
> I want to integrate a function over an irregular polygon. Is there
> any function which can implement this easily? Otherwise, I am
> thinking of divide the polygon into very small rectangles and use
> "adapt" to approximate it. Do you have any suggestions to get the
> fine division? Any advice is appreciated.
>
> Haiyong
>
> ______________________________________________
> 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
> and provide commented, minimal, self-contained, reproducible code.
>
If you can integrate over a trapezoid (with horizontal base and
vertical sides) then it's fairly easy to integrate over any (non-self
intersecting) polygon. ``Draw'' a line which is below the polygon.
For each segment in the polygon, integrate over trapezoid whose base
is on the aforesaid line and whose top edge is the aforesaid edge of
the polygon. If the edge in question goes from left to right,
multiply by -1. (Leave alone if the edge goes from right to left;
vertical edges contribute 0.) Add up the results.
This amounts to applying Green's Theorem in the plane.
I don't know of any code that does this, but it shouldn't be
hard to write. (*Given* that you can integrate over the
trapezoids.)
Another approach:
library(spatstat)
W <- owin(poly=X) # X is your polygon, represented as a list
# with components ``x'' and ``y''
giving
# the vertices of the polygon in *anti*
# clockwise order. The first vertex should
# NOT be repeated at the end.
IM <- as.im(foo,W)
int <- summary(IM)$integral # Based on a 100 x 100 pixellation
# of the bounding box of the polygon.
oop <- spatstat.options(npixel=500)
IM <- as.im(foo,W)
int <- summary(IM)$integral # Based on a 500 x 500 pixellation
# of the bounding box of the polygon.
spatstat.options(oop) # Set npixel back to the default.
HTH.
cheers,
Rolf Turner
rolf at math.unb.ca
On 2/14/2007 9:06 PM, Haiyong Xu wrote:> Hi there, > > I want to integrate a function over an irregular polygon. Is there > any function which can implement this easily? Otherwise, I am > thinking of divide the polygon into very small rectangles and use > "adapt" to approximate it. Do you have any suggestions to get the > fine division? Any advice is appreciated.If you can integrate over a triangle, you will soon be able to use gpclib to triangulate the polygon. (The code was already there in C; I wrote an interface to it and Roger Peng has packaged it up to send to CRAN.) Wait for version 1.4 or greater to appear before downloading. Duncan Murdoch