4

I am aware of how to create centroids and how to divide a polygon by line, but not how to create those lines from each centroid.

I was wondering if anybody knows a technique/algorithm/plugin that helps with splitting a polygon by its centroid, like in pizza slices.

I have some territories that I would like to split by their centroid. Once I have these "pizza slices" of each polygon, I will intersect the slices with a border polygon (created through a Buffer), which aims to perform a spatial RDD (Regression Discontinuity Design) to see whether the points are located in one or the other side of the border. I attach here an example of what I mean (and that I have created manually on Paint to show you).

The red points are the point units whose location I aim to analyze. The yellow dots I have painted myself, represent the centroids of each polygon.

enter image description here

The blue lines are what I would be aiming at. I know that the Voronoi polygons will divide an area of points in equal-size polygons, but that is not what I need.

Could anybody give me a hand on this?

3
  • 1
    So, what keeps you from creating those centroids, moving them in those fixed directions, creating lines from the moved/porjected centroids and splitting your polygons by lines?
    – Erik
    Commented Dec 5, 2023 at 14:57
  • Hi Erik, Because I do not know that intermediate step that you call "moving the centroids in fixed directions". I had in mind to use that tool of splitting the polygon by lines, but not how to create first the lines from the centroid. Thanks for answering so fast!
    – Alex_P
    Commented Dec 5, 2023 at 15:02
  • 2
    e.g. gis.stackexchange.com/questions/458121/…
    – Erik
    Commented Dec 5, 2023 at 15:50

3 Answers 3

6

I can suggest another workflow, perhaps a bit similar to @BERA's answer, however, without conversion from polygons to lines.

Let's assume there is a polygon layer called 'poly' with eight features in it, see the image below.

input

The CRS of this data set is the EPSG:5348.

Note: Centroids were only shown for visual purposes.

Step 1. Use "Geometry by expression" with this expression:

line_merge(
    collect_geometries(
        array_foreach(
            generate_series(0, 135, 45),
            extend(
                intersection(
                    rotate(
                        make_line(
                            project(
                                centroid($geometry),
                                distance(
                                    centroid($geometry),
                                    point_on_surface(
                                        boundary(
                                            minimal_circle($geometry)
                                            )
                                        )
                                    ),
                                radians(0)
                                ),
                            project(
                                centroid($geometry),
                                distance(
                                    centroid($geometry),
                                    point_on_surface(
                                        boundary(
                                            minimal_circle($geometry)
                                            )
                                        )
                                    ),
                                radians(180)
                                )
                            ),
                        @element,
                        center:=centroid($geometry)
                        ),
                    make_polygon(
                        array_first(
                            geometries_to_array(
                                order_parts(
                                    boundary($geometry),
                                    'area(make_polygon(boundary($geometry)))',
                                    False
                                    )
                                )
                            )
                        )
                    ),
                    0.01,
                    0.01
                    )
                )
            )
        )

step1

Step 2. Apply the "Multipart to singleparts" tool to previously created lines.

Step 3. Proceed with the "Split with lines" geoalgorithm to split original polygons with single-part lines and get the output like this:

result

After Step 3, application of the "Multipart to singleparts" tool can be feasible.

6

You can:

  1. Geometry By Expression with the expression: collect_geometries( array_foreach(array:=generate_series(start:=30, stop:=360, step:=60),expression:=intersection($geometry, wedge_buffer( center:=centroid($geometry), azimuth:=@element, width:=60, outer_radius:=3000)))) It creates six wedge buffers with the centroid as center and collects them into a multipart polygon.
  2. Multipart To Singleparts to split the multipart polygon to single polygons.
  3. Polygons to lines
  4. Split the polygon layer with the lines

enter image description here

3

Thanks @BERA for your answer.

It was indeed of help.

I just suggest here a piece of code, inspired by https://www.luisalucchese.com/post/qgis-the-wedge-tool/ . It needs a bunch of centroids created as shapefiles in the "centroid_folder", and you can specify of how many meters (or degrees, depending on your CRS) in OUTER_RADIUS.

import os
numsectors = 16
width = 360.0 / numsectors
centroid_folder = 'C:/Users/.../'
output_folder = 'C:/Users/.../'

# Create the output folder if it doesn't exist
os.makedirs(output_folder, exist_ok=True)



# Loop through each shapefile in the centroid folder
for shp_file in os.listdir(centroid_folder):

    if shp_file.endswith('.shp'):

        input_shp_path = os.path.join(centroid_folder, shp_file)


        # Create a unique output filename for each input shapefile
        output_filename = os.path.splitext(shp_file)[0] + '_segmented.shp'
            output_shp_path = os.path.join(output_folder, output_filename)

        for k in range(0, numsectors):

            azim = width / 2 + k * width
            param_wedge = {
                'AZIMUTH': azim,
                'INNER_RADIUS': 0,
                'INPUT': input_shp_path,
                'OUTER_RADIUS': 95000,
                'OUTPUT': output_shp_path.replace(".shp", f"_sector_{k}.shp"),
                'WIDTH': width

            }

            layer_wedge = processing.run("native:wedgebuffers", param_wedge)

            # Load the created shapefile as a vector layer
            wedge_layer = QgsVectorLayer(layer_wedge['OUTPUT'], f"Sector {k}", 'ogr')

            # Add the vector layer to the project
            QgsProject.instance().addMapLayer(wedge_layer)

print("Wedge buffers created and saved in:", output_folder)
1
  • 1
    Thank you @L.Lucchese for a great tutorial/blog!
    – Taras
    Commented Dec 11, 2023 at 16:00

Not the answer you're looking for? Browse other questions tagged or ask your own question.