Intersecting more than two layers at same time using QGIS

Question

Background:

At 3.20.3, there are two intersect tools available: Vector > Geoprocessing tools > Intersection and Processing Toolbox > SAGA > Features - Polygons > Intersect.

In both cases, only two layers can be intersected at a time. However, I have numerous layers that need to be intersected, and it would be much more efficient to conduct the task once, rather than two-at-a-time.

I notice that the Advanced ArcGIS Desktop 10.8.1 Intersection tool allows more than two input layers. I've searched the internet, plugins, and the documentation without success for a similar QGIS solution.

Question:

How do I intersect more than two layers at a time using QGIS?

Answer 1

Short answer: the principles

You can use QGIS expressions to generate custom overlay like intersection of as many layer as you want. Use Geometry generator or Geoemtry by expression (see here for differences).

How to do it

Let's say you have 5 polygon layers, name p1 to p5. Say you look for the areas where all five intersect. Split up this task in several steps. This is because the expression intersect() only takes two arguments (input layers), so you have to create a "nested" (combined) statement, say like (pseudocode):

intersect the (interesction of A and B) with the (interesction of C and D) 

Or more detailed:

1. Intersection of p2 with p3 (lines 4-7 in the expression below)
2. Intersection of p4 with p5 (lines 8-11)
3. Intersect the result of steps 1 and 2 (you get: intersection of layers p2, p3, p4 and p5)
4. Intersect p1 with the result of step 3.

The expression to use

The expression looks like this (see below at the bottom for detailed explanation of the functions used):

intersection (
    $geometry, 
    intersection (
        intersection (
            buffer (collect_geometries(overlay_intersects('p2', $geometry)),0), 
            buffer (collect_geometries (overlay_intersects( 'p3', $geometry)),0)
        ),
        intersection (
            buffer (collect_geometries(overlay_intersects('p4', $geometry)),0), 
            buffer (collect_geometries (overlay_intersects( 'p5', $geometry)),0)
        )
    )
)

Screenshot: the expression with Geometry generator: the diagonal X hached area is created with the expression above - the area where all 5 layers intersect:

Variant: how to adapt the expression to changing number of layers

For an even number of intersecting polygons, adapt the expression as follows (here: 6 polygons p1 to p6). As you see, I added a buffer([geometry],0) statement. Sometimes, without this, the expression returns NULL.

intersection (
    intersection (
        buffer (collect_geometries (overlay_intersects ('p5',$geometry)),0), 
        buffer (collect_geometries (overlay_intersects ('p6',$geometry)),0)
    ),
    intersection (
        intersection (
            buffer (collect_geometries (overlay_intersects ('p4',$geometry)),0), 
            buffer (collect_geometries (overlay_intersects ('p3',$geometry)),0)
        ),
        intersection(
            $geometry, 
            collect_geometries (overlay_intersects ('p2',$geometry))
        )
    )
)

Screenshot: the diagonal X hached area, representing the common space of all 6 polygons, is bordered by all 6 polgons as you can see from the color of the outer boundaries: red, blue, green, yellow, black and pink:

Explanation: How this expression works

• The intersection is done using the intersection() function, this only takes two arguments: intersection(geometry1, geometry2).
• So we have to combine intermediate steps with nested (combined) intersection functions. Geometry generator as well as Geoemtry by expression are always based on a specific layer, so to get the geometries of the other layers to combine (intersect) them, we use the overlay_intersect() function.
• This returns an array of all geometries, so to get geometries, we must convert from array to geometries using collect_geometries().
• To work properly, we must buffer these geometries by 0, using buffer() function. For the very layer we work on, it's enough to get the geometries with $geometry.

Answer 2

As your initial question focused on intersection only, without attributes, my first answer did just that. As you specified later, you're interested in the attributes as well. So this is best done using a virtual layer, based on the same principles: as st_intersection() takes just two arguments, you have to combine several st_intersection() statements.

To create a new layer with the common area of four layers p1 to p4 plus all attributes, use this query in the dialog creating the "Virtual Layer" through Layer > Add Layer > Add/Edit Virtual Layer...:

select
    st_intersection(st_intersection(a.geometry, b.geometry), st_intersection(c.geometry, d.geometry)),
    *
from
    p1 as a, p2 as b, p3 as c, p4 as d

Screenshot: the hatched area in pink is the area where all four layers intersect. It's created as a virtual layer. As you can see from the attribute tables, it does not only create the geometry, but also adds all attributes, thanks to the , * statement at the end of the first line: