Create an inscribed square inside a polygon

Question

I have polygons of the administrative boundaries and want to create a square polygon inside of them, which will have the maximum area and will be fully covered by these boundaries. In other words: how to create an inscribed square inside a polygon.

The example is below

Answer 1

This solution is based on the 4 white dots in the following screenshots. They are the closest points to the or pole of inaccessability for each section of the outer boundary of the polygon. Section here means: divide the polygon in four parts (north/east/south/west) in a diagonal-like manner. From this rectangle, get the shorter side and reduce the longer side to this length.

From the four white points, create a bounding box to get the blue rectangle. From this, it's ease to get the red square, based on its diagonal: starting at the or pole of inaccessability, going in an angle of 45 degrees until it intersects the boundary of the rectangle.

Find a step-by-step description below.

Remarks:

Better than using centroid () (as propagated in the first version of this solution) is pole_of_inaccessibility( ) (see here for the concept of pole of inaccessability). As you can see on the following screenshot, some country shapes (like the contiguous territories of Italy or especially Croatia) are irregularily shaped in a way that the centroid does not fall in the place with most space around it inside the country borders: but we need as much space as possible to get the biggest square possible inside the territory. Thus the pole of inaccessability is better. As the first version of this solution was besed on centroid, the screenshots inserted in the steps are still based on this - with the shape used for this, centroid or pole of inaccessability are almost the same and thus does not affect much how the solution looks - for Italy or Croatia, however, it does!

Blue dots: pole of inaccessability, white dots: centroids:

Centroid of Norway even falls outside the territory of Norway: it is inside of Sweden, whereas the pole of inaccessability by definition is alwasy inside the territory:

Step by step description:

For each step, use the output of the previous step, if not stated otherwise:

1. Use Geometry by expression to create a bounding box, subdivided in 4 sections centered on the pole of inaccessibility with this expression:

collect_geometries (
    array_foreach (
        array (1,2,3,4),
        make_triangle (
            pole_of_inaccessibility ($geometry,1),
            point_n( bounds( $geometry),@element),
            point_n( bounds( $geometry),@element+1)
        )
    )
)

Orange (covered): original polygon; blue: polygons created with the expression:

2. Menu Vector / Geometry Tools / Multipart to singleparts

3. Menu Vector / Geoprocessing Tools / Intersection, set the polygon as Input layer and the singleparts from step 2 as Overlay layer

4. Menu Vector / Geometry tools / Polygons to lines

5. Add a unique id to the layer create in step 4, using field calculator with $id or @row_number to create a field id

6. Menu procesing / Toolbox / Explode lines

7. Select by expression with this expression - replace 'poly' with the name of your original polygon layer: disjoint ($geometry, boundary ( geometry ( get_feature_by_id ('poly',1))))

   Result: the diagonals are selected:

8. Toggle editing, delete selected (diagonals)

9. Menu Processing / Toolbox / Aggregate, set the id in the field Group by expression (NULL to group all features) and under Source Expression select first value as Aggregate Function for the id.

   

10. Use Geometry by expression to create the pole of inaccessibility of the original polygon with this expression: pole_of_inaccessibility ($geometry,1)

11. Use Geometry by expression on the pole of inaccessibility layer with the following expression to create the 4 closest point to the pole of inaccessibility : each for every 4 sections:

collect_geometries (
    array_foreach (
        array (1,2,3,4),
        closest_point (
            geometry (
                get_feature_by_id (
                    'Aggregiert',
                    @element
                )
            ),
            $geometry
        )
    )
)

Centroid (or pole of inaccessability, blue) and the 4 closest point, one per each section:

12. Use Geometry by expression with the following expression: bounds($geometry): you get the rectangle from the very first screenshot above.

13. To convert the rectangle to a square, use Geometry by expression with the following expression - adapt the value in line 11 from 1000 in my example to your needs - it can't be too long, so if unsure, just add several zeros:

with_variable (
    'diagonal',
    with_variable (
        'diagonal_half',
        make_line (
            intersection (
                make_line (
                    centroid ($geometry),
                    project (
                        centroid ($geometry),
                        1000,
                        radians (45)
                    )
                ),
                boundary ($geometry)
            ),
            centroid($geometry)
        ),
        extend (
            @diagonal_half,
            0,
            length (@diagonal_half)
        )
    ),
    make_square (
        start_point (@diagonal),
        end_point (@diagonal)
    )
)