Recently, while out in the field, I encountered some issues with sample plot division that conventional grid division methods could not generate. What should I do in this case? With Kimi’s help, I used QGIS Python programming to complete some advanced custom geographic analysis functions.

Problem

How to divide the lakeside into 50 grids?

Problem Breakdown

1. Extract the lake area range
2. Extract the lakeshore line
3. Generate a grid based on the lakeshore line

Extracting the Lake Area Range

I used the Zhongke Xingtou Smart Brain Engine to complete the extraction of the lake area range. The Zhongke Xingtou Smart Brain Engine provides a very convenient water body extraction function, which can effectively complete the extraction of water body ranges.

• Xingtou Earth Homepage: https://geovisearth.com/

The Xingtou Earth Smart Brain Engine offers many convenient tools, combined with high-resolution images provided by the Zhongke Xingtou Data Cloud, making it easy to complete some water body extraction tasks.

The processing speed is very fast. After extracting the water body, save the results to the cloud drive, and then download them from the cloud drive. The results can be saved directly in geojson format, allowing for direct opening in QGIS.

Extracting the Lakeshore Line

In QGIS, you can directly convert the lake area range into a line using the Polygons to Lines tool.

However, I found that the resulting line was too complex for generating a grid, so I ultimately chose to manually draw a line along the lakeshore.

First, create a gpkg database to facilitate the storage of vector data. Then right-click, New Table; a vector line is equivalent to a new table in the gpkg.

Here, I need to create a line vector.

Then drag this line vector into the layer and start editing.

Add line features, start drawing, and after finishing, click the pencil button to save and stop editing.

Generating Sample Plots Along the Lakeshore

Note that this command was generated after multiple discussions with Kimi regarding QGIS. The complete dialogue process can be referenced:

Click the link to view the conversation with Kimi: https://www.kimi.com/share/d1psv64nvj4kt3e4p0mg

Command

Source H:/Lashihai/Lashihai.gpkg|layername=WaterLine

Please help me generate 50 equidistant points based on this line, and then generate 50 bounding rectangles based on these 50 points. Please complete the Python code.

Kimi’s Response

Below is a complete script that can be run directly in the QGIS Python console. Logic:

1. Read the first line feature from the line layer (WaterLine);
2. Generate 50 equidistant points along the line length L (including the starting point, excluding the endpoint);
3. Generate square bounding rectangles with side length s = L/50 centered at each point;
4. Output the results to a new temporary layer “WaterLine_Grid50”.

# -*- coding: utf-8 -*-
"""
Generate 50 equidistant points along WaterLine and create square bounding rectangles centered at each point
"""

from qgis.core import (
    QgsProject, QgsVectorLayer, QgsFeature, QgsGeometry,
    QgsRectangle, QgsPoint, QgsField, QgsWkbTypes
)
from qgis.PyQt.QtCore import QVariant

# 1. Load the line layer
uri = r"H:/Lashihai/Lashihai.gpkg|layername=WaterLine"
line_lyr = QgsVectorLayer(uri, "WaterLine", "ogr")
if not line_lyr.isValid():
    raise RuntimeError("Failed to load line layer!")

# 2. Get the first line
line_feat = next(line_lyr.getFeatures())
line_geom = line_feat.geometry()
line_length = line_geom.length()

# 3. Calculate the square side length
side = line_length / 50.0
print(f"Square side length: {side:.2f} map units")

# 4. Create output layer
grid_lyr = QgsVectorLayer("Polygon?crs=" + line_lyr.crs().authid(),
                          "WaterLine_Grid50", "memory")
prov = grid_lyr.dataProvider()
prov.addAttributes([QgsField("id", QVariant.Int)])
grid_lyr.updateFields()

# 5. Generate 50 points → 50 squares
for i in range(50):
    # Equidistant point position (0~49)
    dist = (i * line_length) / 50.0
    pt = line_geom.interpolate(dist).asPoint()          # QgsPointXY
    # Bounding rectangle centered at this point
    xmin = pt.x() - side / 2.0
    xmax = pt.x() + side / 2.0
    ymin = pt.y() - side / 2.0
    ymax = pt.y() + side / 2.0
    rect_geom = QgsGeometry.fromRect(QgsRectangle(xmin, ymin, xmax, ymax))

    feat = QgsFeature()
    feat.setGeometry(rect_geom)
    feat.setAttributes([i])
    prov.addFeature(feat)

grid_lyr.updateExtents()
QgsProject.instance().addMapLayer(grid_lyr)
print("Generated 50 bounding rectangles and loaded them into the map.")

Usage Instructions

1. Open QGIS → Plugins → Python Console (<span>Ctrl+Alt+P</span>).
2. Paste the above code and press Enter to run.

This way, sample plots are generated along the lakeshore.

★

Reprinted from: Xu Xiaoyang Geographic Data Science

The article only represents the author’s views and is not related to this public account. Copyright belongs to the original author.

Original title:Advanced Geographic Analysis Functions Using QGIS Python Programming——Taking Lakeshore Sample Plot Generation as an Example

·END·

Text and image editing: Fan Xuesong Deng Huiling

Reviewed by: Wang Yuyao Guo Shuyang

Final review: Chu Mingruo Li Yuzhu Dai Haoyu Bi Siqi·Recommended Articles·

1. Call for Papers | Artificial Intelligence in Coastal Wetlands Remote Sensing

2. Article Appreciation | Scientific Bulletin: Ancient DNA Reveals Ancient Family Lineages and Social Structures

3. Geographic Data | 2024 Heihe River Basin Surface Process Comprehensive Observation Network Dataset

4. Journal Directory | Recommended Articles from the Journal of Earth Information Science, Issue 8, 2025