Compute surfaces of intersection between circles in julia with LibGEOS

Specific Domains Geo
,
1

I need to compute surfaces of intersection between different circles. I am used to do it in R with the sf library and it works fine but I need to do it in Julia with the LibGEOS library and it is not as easy…

My first question would be, how to create a circle polygon? By creating the center and extend it with a buffer? For example:
circle_i = buffer(point_i, radius_i) ?

My second question is, now that I have created some circles, can I put them in a common object, like a multipolygon?

My third and last question is, if I have an ensemble of circles, can I use the intersection function to compute all the intersections of all my circles?

I hope I am clear enough, I could really use some help regarding these questions!

2

Hi, welcome! Yeah this sure is possible. Not quite as well documented as sf though.

GEOS doesn’t have circles, only the geometry types listed here. Indeed taking a buffer of a point is a good way to create a circular polygon.

Why do you want to put the circles together in a common object? It should be doable to put them together in a MultiPolygon or GeometryCollection, but in the example below a simple Vector of Polygons is probably easier.

using LibGEOS, Plots

# 3 overlapping circles in a triangle
circles = [
    buffer(Point(0.0, 0.0), 1.0),
    buffer(Point(1.0, 0.0), 1.0),
    buffer(Point(0.5, 1.0), 1.0),
]

isects = Polygon[]
for (i, g1) in enumerate(circles)
    for (j, g2) in enumerate(circles)
        # avoid counting double, and only use non empty intersections
        if i > j && LibGEOS.intersects(g1, g2)
            isect = LibGEOS.intersection(g1, g2)
            push!(isects, isect)
        end
    end
end

plot(isects)

image
image895×599 29.1 KB

3

Thank you so much for your answer, it is very helpful!

I guess for the intersection function I wanted an object like a multipolygon because I need to compute all the intersections between the circles, not only 2 by 2 intersections. With your example, I would also need to have a polygon for the intersection between the three circles. Ideally I would need as many polygons as overlapped surfaces, just like this figure below:

ex_R
ex_R747×465 44.3 KB

But since the intersection function only takes 2 polygons as arguments it might not be possible?

4

Hmm, I’m actually not sure how to do this. So in your figure there are no overlapping areas anymore right? It looks like the result of a union in QGIS. However union in GEOS also directly dissolves the result. Though probably with a combination of the other functions you can effectively do the same?

# need to go through GeoInterface, see issue #87
mp = MultiPolygon(GeoInterface.coordinates(circles))
plot(LibGEOS.unaryUnion(mp))

image
image889×590 22.1 KB

5

You should be able to write each colored surface as intersections between circles or their complements.
For example, the orange region is A \cap B \cap \neg C, where A is the top circle, B is the right circle and C is the left circle.

1 Like