Boundary line expression method applied to territorial planning process

Abstract

The invention relates to the technical field of big data analysis, in particular to a boundary line expression method applied to a territorial planning process, and the method comprises the steps: obtaining hill landform elevation data, extracting contour skeleton change features, screening directional coupling channels, recognizing adjacent boundary paths, forming adhesion sections, detecting land class abnormal fluctuation, and clustering and extracting heterogeneous boundaries. And a space fitting segment group is generated through superposition. According to the method, the same-layer contour skeleton segments are extracted from the elevation grid data, then the curvature and torsion angle changes of the contour skeleton segments are calculated, and the continuous skeleton sequence with the stable direction is screened out in combination with the included angle between principal component directions and the change stability index, so that the spatial information structure has the fine expression ability for the terrain continuity and the turning trend. And furthermore, an adjacent space is constructed at the periphery of the skeleton endpoint, and based on mutual information calculation between a path center line sequence and a skeleton direction sequence, a boundary path section with high direction fitness is identified, so that the matching accuracy of boundary extraction in geometry and direction is improved.

Description

Technical Field

[0001] This invention relates to the field of big data analysis technology, and in particular to a method for representing boundary lines in the process of land planning. Background Technology

[0002] The field of big data analytics encompasses a technological system for collecting, storing, processing, and mining massive amounts of diverse and rapidly changing data. Through high-performance computing and data modeling techniques, valuable information and patterns are extracted from widely distributed data sources to support decision-making.

[0003] Among them, the boundary line expression method applied in the land planning process refers to the method of spatial layout and functional division of land resources in the planning area by using data from multiple sources such as remote sensing data, geographic information data, socio-economic data and land use status data, and by setting spatial assessment rules, land suitability assessment models and land use structure optimization algorithms.

[0004] Existing technologies rely solely on basic data sources for spatial assessment and layout optimization. While they can achieve basic land use suitability classification, in practice they often overlook the continuous expression of topographic transition features in hilly terrain. This leads to severe interference with boundary lines due to topographic structure, resulting in local breaks or deviations from the dominant topographic direction. In boundary identification, relying solely on linear overlay or buffer analysis lacks quantitative identification of directional correlations, easily leading to arbitrary selection of boundary segments and affecting the logical closure and natural rationality of the final planning boundary. Furthermore, in handling changes in land use attributes, there is a lack of methods to analyze the amplitude and directional differences of time-series attributes, failing to effectively identify dynamic trends and transitional phenomena within the region. This results in low identification accuracy and ambiguous classification in areas of transitional land functions. For example, in urban-rural fringe areas, the drastic fluctuations in the usage frequency of standard plots may not be effectively identified, leading to planning results that fail to accurately depict the actual usage status of the boundary. In attribute clustering, existing methods generally do not consider spatial gradient structures, resulting in clustering results lacking spatial continuity and directional logic. This easily leads to the grouping of areas with different directions but similar attributes into the same category, weakening the spatial consistency and comparability of boundary adjustments. These issues can lead to a series of problems in actual land planning work, such as insufficient precision of layer boundaries, logical breaks, and failure to pass review. Summary of the Invention

[0005] The purpose of this invention is to address the shortcomings of existing technologies by proposing a method for representing boundary lines in the process of land planning.

[0006] To achieve the above objectives, the present invention adopts the following technical solution: a method for representing boundary lines in the process of land planning, comprising the following steps: S1: Obtain elevation raster data of hilly landform areas, extract the curvature and torsion angle changes of contour line skeleton segments at the same level, make consistency judgment on the direction of change, and screen contour line direction coupling channel groups. S2: Obtain the spatial projection region of the directional coupling skeleton segment group in the contour line directional coupling channel group, identify the boundary candidate path of the spatial adjacent position of the skeleton segment group, and mark the boundary adhesion segment set from it; S3: Obtain the surface space range corresponding to the boundary adhesion segment set, perform land type attribute fluctuation detection and identify the location of attribute amplitude anomalies, and form a set of land type attribute anomaly fluctuation points; S4: Obtain the adjacent attribute vectors of the set of abnormal fluctuation points of the land type attribute, perform spatial attribute clustering analysis and filter the heterogeneous boundary locations to form a set of attribute heterogeneous clustering segments.

[0007] The present invention is improved in that the contour line directional coupling channel group includes a skeleton segment directional sequence, curvature change direction, and torsion angle change direction; the boundary adhesion segment set includes a candidate path tangential change sequence, path-skeleton segment direction correspondence, and adhesion strength score; the land type attribute anomaly fluctuation point set includes anomaly location index, land type attribute change trend, and change direction difference; and the attribute heterogeneous clustering segment set includes attribute clustering labels, attribute distribution discreteness, and directional gradient comparison results.

[0008] The present invention improves upon this invention by acquiring elevation raster data of hilly landform regions, extracting the curvature and torsion angle changes of contour line skeleton segments at the same level, determining the consistency of the change direction, and screening contour line directional coupling channel groups. The specific steps are as follows: S101: Obtain elevation raster data of hilly landform areas, generate contour lines of the same layer at fixed elevation intervals, extract the center point sequence of each skeleton segment in the contour lines, collect the three-dimensional spatial coordinates corresponding to the center points, calculate the rate of curvature change and the change value of torsion angle of the skeleton segment based on the spatial position change between adjacent center points, and generate a sequence of geometric change indexes of the skeleton segment. S102: Extract the geometric change index sequence of the skeleton segment, calculate the cosine value of the included angle between the principal component directions of adjacent skeleton segments as the direction consistency index, calculate the standard deviation of the curvature change rate sequence as the fluctuation stability index, and obtain the direction and fluctuation characteristics of the skeleton segment. S103: Based on the results of the skeleton segment direction and fluctuation characteristics, select skeleton segment combinations with an included angle cosine value greater than the direction consistency threshold and a standard deviation less than the fluctuation stability threshold, extract the spatial continuous distribution of the skeleton segment combination in the target elevation level, and generate contour line direction coupling channel groups.

[0009] The present invention improves upon this invention by obtaining the spatial projection region of the directionally coupled skeleton segment group in the contour line directionally coupled channel group, identifying the boundary candidate paths of the spatially adjacent positions of the skeleton segment group, and marking the boundary adhesion segment set therefrom. The specific steps are as follows: S201: Call the spatial range of the contour line direction coupling channel group, obtain the preliminary boundary division data of the corresponding position, and detect the linear boundary entities in the adjacent position range according to the endpoint coordinates and direction of the skeleton segment, and generate an adjacent path index set. S202: Based on the adjacent path index set, extract the centerline sequence of each path segment, calculate the tangential radix change value between continuous nodes of the path segment, call the principal component direction sequence of the corresponding skeleton segment in the contour line direction coupling channel group, calculate the mutual information value between the two types of direction sequences as the degree of direction coupling between the path and the skeleton segment, and obtain the path adhesion score result. S203: Filter the path segments whose path adhesion score is greater than the adhesion score threshold and whose directional trend is consistent, and generate a set of boundary adhesion segments.

[0010] The present invention is improved by obtaining the surface spatial range corresponding to the boundary adhesion segment set, performing land use attribute fluctuation detection and identifying the location of attribute amplitude anomalies, and forming a land use attribute anomaly fluctuation point set. The specific steps are as follows: S301: Call the land planning data raster corresponding to the boundary adhesion segment set, collect the two parameters of land use density and unit block usage frequency in each raster, and perform normalization processing to generate a normalized land type numerical attribute set; S302: Based on the normalized land use numerical attribute set, calculate the attribute difference between two adjacent time periods for each grid cell, and determine whether the directions of the two sets of differences are opposite, to obtain the attribute difference direction judgment result; S303: Call the attribute difference direction judgment result, compare the attribute change amplitude in the middle time period with the average difference of the two adjacent time periods, filter the positions where the change amplitude is greater than the average value, and generate a set of abnormal fluctuation points of land category attributes.

[0011] The present invention improves upon this invention by obtaining the adjacent attribute vectors of the set of anomalous fluctuation points of land type attributes, performing spatial attribute clustering analysis, and filtering heterogeneous boundary locations to construct a set of attribute heterogeneous clustering segments. The specific steps are as follows: S401: Obtain the location of the abnormal fluctuation point set of the land category attribute in the land planning data raster, including the land use density of each raster in the eight adjacent directions and the usage frequency of the standard plot within the specified area, and perform normalization processing to form an adjacent attribute vector data group. S402: Based on the adjacency attribute vector data group, group and cluster according to Euclidean distance similarity, and calculate the variance value of each attribute dimension in each group to obtain the clustering attribute dispersion index set. S403: Call the clustering attribute dispersion index set, compare the attribute variance of the clustering result with the cosine value of the angle between the direction gradient, filter the positions where the variance value is greater than the dispersion threshold and the cosine value of the angle is greater than the direction threshold, and generate a set of attribute heterogeneous clustering segments.

[0012] The present invention is improved and further includes step S5: combining the grid positions of the boundary adhesion segment set and the attribute heterogeneous clustering segment set, extracting the spatially overlapping region to obtain the spatial planning fitting segment group; The spatial planning fitting segment group specifically includes boundary fusion location, spatial unit index, and layer update area.

[0013] The present invention improves upon this invention by combining the grid positions of the boundary adhesion segment set and the attribute heterogeneous clustering segment set to extract spatially overlapping regions and obtain spatial planning fitting segment groups. The specific steps are as follows: S501: Based on the raster positions of the boundary adhesion segment set and the attribute heterogeneous clustering segment set in the land planning data, calculate the spatial overlap ratio of the two types of segments, record the spatial index of the overlapping area, and generate the spatial overlap area ratio. S502: Based on the spatial overlap area ratio, extract the directional consistency judgment results of the corresponding positions already marked in the attribute heterogeneous clustering segment set to establish directional filtering criteria and obtain the directional consistency filtering label set; S503: Call the spatial overlap area ratio and the directional consistency filtering label set, compare the overlap ratio with the set overlap threshold, filter the positions with positive directional consistency labels, extract the spatial positions that meet the two conditions, generate a spatial planning fitting segment group, which serves as the basis for updating the boundary to be adjusted in the land planning layer, and can be input into the layer editing platform or planning review process to complete the replacement or comparison verification of the original partition boundary.

[0014] Compared with the prior art, the advantages and positive effects of the present invention are as follows: In this invention, during the topographic structure processing of hilly landform areas, contour line skeleton segments of the same layer are extracted from elevation raster data. Their curvature and torsion angle changes are then calculated. Combined with the angle between principal component directions and the stability index of changes, continuous and directionally stable skeleton sequences are selected, enabling the spatial information structure to accurately express the continuity and transition trends of the terrain. Furthermore, an adjacency space is constructed around the skeleton endpoints. Based on the mutual information calculation between the path centerline sequence and the skeleton direction sequence, boundary path segments with high directional fit are identified, forming boundary adhesion spatial zones, significantly improving the accuracy of boundary extraction in both geometry and direction. In land attribute data processing, relying on land use density and standard plot usage frequency, time series difference direction and amplitude fluctuation analysis is used to locate abnormal raster points with unstable attribute trends, providing a basis for marking key transition zones in subsequent spatial partitioning. By aggregating the attribute vectors of these outliers in eight adjacent directions, clustering based on attribute similarity is performed. Combining the internal variance and spatial gradient direction changes of each group, spatial segments with uneven attribute performance and dissimilar structural orientations are extracted, thus accurately delineating heterogeneous regions with complex functions or ambiguous boundaries. Finally, heterogeneous segments and directionally adherent segments are analyzed on a raster layer for overlap. A combination of overlap degree and directional consistency is used to screen for alternative segments, achieving dual correction of the original boundary in terms of function and spatial structure. This effectively enhances the adaptability and resolution of spatial zoning in changing terrains and multi-attribute spaces, improving the rationality, accuracy, and consistency of operational logic in spatial boundary adjustments. Attached Figure Description

[0015] Figure 1 This is a flowchart of the method of the present invention; Figure 2 This is a detailed flowchart of step S1 of the present invention; Figure 3 This is a detailed flowchart of step S2 of the present invention; Figure 4 This is a detailed flowchart of step S3 of the present invention; Figure 5 This is a detailed flowchart of step S4 of the present invention; Figure 6 This is a detailed flowchart of step S5 of the present invention. Detailed Implementation

[0016] To make the objectives, technical solutions, and advantages of this invention clearer, the invention will be further described in detail below with reference to the accompanying drawings and embodiments. It should be understood that the specific embodiments described herein are merely illustrative and not intended to limit the invention.

[0017] In the description of this invention, it should be understood that the terms "length," "width," "upper," "lower," "front," "rear," "left," "right," "vertical," "horizontal," "top," "bottom," "inner," and "outer," etc., indicating orientation or positional relationships, are based on the orientation or positional relationships shown in the accompanying drawings and are only for the convenience of describing the invention and simplifying the description, and do not indicate or imply that the device or element referred to must have a specific orientation, or be constructed and operated in a specific orientation, and therefore should not be construed as a limitation of the invention. Furthermore, in the description of this invention, "a plurality of" means two or more, unless otherwise explicitly specified.

[0018] Please see Figure 1 This invention provides a technical solution: a method for representing boundary lines in the process of land planning, comprising the following steps: S1: Obtain elevation raster data of hilly landform areas, extract the curvature and torsion angle changes of contour line skeleton segments at the same level, make consistency judgment on the direction of change, and screen contour line direction coupling channel groups. S2: Obtain the spatial projection region of the directional coupling skeleton segment group in the contour line directional coupling channel group, identify the boundary candidate path of the spatially adjacent position of the skeleton segment group, and mark the boundary adhesion segment set from it; S3: Obtain the surface space range corresponding to the boundary adhesion segment set, perform land use attribute fluctuation detection and identify the location of attribute amplitude anomalies, and form a set of land use attribute anomaly fluctuation points; S4: Obtain the adjacent attribute vectors of the set of abnormal fluctuation points of land type attributes, perform spatial attribute clustering analysis and filter the heterogeneous boundary locations to form a set of attribute heterogeneous clustering segments; The contour line directional coupling channel group includes the skeleton segment directional sequence, curvature change direction, and torsion angle change direction; the boundary adhesion segment set includes the candidate path tangential change sequence, path and skeleton segment direction correspondence, and adhesion strength score; the land type attribute anomaly fluctuation point set includes anomaly location index, land type attribute change trend, and change direction difference; and the attribute heterogeneous clustering segment set includes attribute clustering labels, attribute distribution discreteness, and directional gradient comparison results.

[0019] Please see Figure 2 The specific steps for obtaining elevation raster data of hilly landform regions, extracting the curvature and torsion angle changes of contour line skeleton segments at the same level, judging the consistency of the change direction, and screening contour line direction coupling channel groups are as follows: S101: Obtain elevation raster data of hilly landform areas, generate contour lines of the same layer at fixed elevation intervals, extract the center point sequence of each skeleton segment in the contour lines, collect the three-dimensional spatial coordinates corresponding to the center points, calculate the rate of curvature change and the change value of torsion angle of the skeleton segment based on the spatial position change between adjacent center points, and generate a sequence of geometric change indexes of the skeleton segment. Obtaining elevation raster data for hilly landform areas refers to collecting 30-meter resolution digital elevation model (DEM) files from a data platform with a legitimate source. The data format is typically GeoTIFF, and the coordinate system is generally WGS84. The rasterio library in Python can be used to read the raster values, and the data can be layered according to a set elevation interval. For example, if the interval is set to 10 meters, and the elevation range of a certain area is 120 meters to 280 meters, then 17 elevation layers will be generated. Each layer will generate a contour line layer based on the elevation value. The contour lines can be generated using the contour extraction tool in GDAL or QGIS. Subsequently, each contour line will be divided into multiple skeleton segments according to the vertex structure. The center point coordinates of each segment are calculated using the midline. For example, the path point coordinates of segment 5 are (120.1, 45.2, 135.6), (121.0, 45.3, 135.8), and (122.0, 45.6, 135.7). The center point coordinates are the average of these three coordinates, i.e., (121.0, 45.37, 135.7). The center point coordinates of each skeleton segment form a sequence. The three-dimensional spatial distance between adjacent center points is calculated using the formula: ; in For the first The longitude, latitude, and elevation values ​​of each center point. For the first Let's take the coordinates of two adjacent points, for example, the second and third points. , ,but: ; Then, the angle between two vectors formed by the three points is calculated. For example, the vectors for points A, B, and C are: , : like , Then use the cosine formula for the included angle: ; in ; ; ; ; This value indicates that the change in this segment is relatively smooth. The curvature change rate and torsion angle change values ​​of all skeleton segments are stored in sequence to form a sequence of geometric change indexes for skeleton segments.

[0020] S102: Extract the geometric change index sequence of the skeleton segment, calculate the cosine value of the angle between the principal component directions of adjacent skeleton segments as the direction consistency index, calculate the standard deviation of the curvature change rate sequence as the fluctuation stability index, and obtain the direction and fluctuation characteristics of the skeleton segment. After extracting the geometric change index sequence of the skeleton segment, the principal component direction vector of each segment is extracted as the direction reference of the skeleton segment. The cosine of the angle between the principal axes of the vectors of consecutive segments is calculated as the direction consistency index. Taking the direction vectors of two center point sequences as an example, let the direction vector of segment 1 be... The direction vector of segment 2 is ,but: ; ; ; ; This indicates extremely high directional consistency. For the fluctuation index, the curvature change rate sequence is extracted, such as {0.02, 0.03, 0.02, 0.025, 0.04}. With a mean of 0.027, the standard deviation is calculated to be 0.007, which is lower than the set threshold of 0.015. Therefore, the change in this segment is considered stable, and the directional and fluctuation characteristics of the skeleton segment are finally generated.

[0021] S103: Based on the results of skeleton segment direction and fluctuation characteristics, select skeleton segment combinations with included angle cosine values ​​greater than the direction consistency threshold and standard deviation less than the fluctuation stability threshold, extract the spatial continuous distribution of skeleton segment combinations in the target elevation level, and generate contour line direction coupling channel groups. Based on the results of skeleton segment direction and fluctuation characteristics, a filtering operation is performed, comparing the cosine value of the included angle with a direction consistency threshold. For example, if the direction consistency threshold is set to 0.85, and a segment has a cosine value of 0.92, it is considered to meet the direction consistency standard. If the standard deviation is 0.008, which is less than the set fluctuation stability threshold of 0.015, it also meets the condition. The filtering rule is that both indicators must meet the set thresholds. Taking a skeleton segment sequence in a contour layer as an example, if segments 3 to 7 meet the above conditions, they are marked as a continuous combination. Then, their paths in the spatial layer are connected into continuous segments and mapped to the DEM spatial reference to construct projection paths. After extracting all combinations that meet the filtering conditions, their raster coordinates are arranged in spatial order and uniformly marked as contour line direction coupling channel groups.

[0022] Please see Figure 3 The specific steps for obtaining the spatial projection region of the directionally coupled skeleton segment group in the contour line directionally coupled channel group, identifying the boundary candidate paths of the spatially adjacent positions of the skeleton segment group, and marking the boundary adhesion segment set are as follows: S201: Call the spatial range of the contour line direction coupling channel group, obtain the preliminary boundary division data of the corresponding position, and detect the linear boundary entities within the adjacent position range according to the endpoint coordinates and direction of the skeleton segment, and generate an adjacent path index set. The spatial range of the contour line directional coupling channel group refers to obtaining the horizontal projection area of ​​the extracted hilly landform area in the corresponding two-dimensional surface raster of the digital elevation model, based on the continuous skeleton segments constituting the channel group. This area is usually a curved, continuous, narrow strip with strong topographical differentiation at its boundary. It is necessary to initially divide the data based on the linear boundary corresponding to this area in the land planning vector layer, extract its coordinate attributes, and further analyze the spatial distribution of each boundary segment. By recording the coordinates of the center points of the beginning and end of each skeleton segment, its direction vector is constructed. After calculating the direction of the skeleton segment, a spatial adjacency threshold is set, and a buffer zone of fixed width is extended to both sides of the skeleton segment. All boundary segments with linear structural characteristics are retrieved within the buffer zone. At the same time, the starting and ending coordinates of these boundary segments are extracted and compared with the direction of the skeleton segment. Boundary segments with similar starting and ending directions are marked. Boundary segments that meet the spatial adjacency relationship and whose directional structure matches are numbered sequentially, their index numbers are recorded, and a set structure is constructed. Finally, an adjacency path index set is generated.

[0023] S202: Based on the adjacent path index set, extract the centerline sequence of each path segment, calculate the tangential guide change value between continuous nodes of the path segment, call the principal component direction sequence of the corresponding skeleton segment in the contour line direction coupling channel group, calculate the mutual information value between the two types of direction sequences as the degree of directional coupling between the path and the skeleton segment, and obtain the path adhesion score result. Based on the adjacent path index set, the centerline point sequence corresponding to each candidate path segment in the hilly landform region is extracted. This sequence consists of several consecutive nodes in the path segment. Each node contains longitude, latitude, and path arrangement order, which is used to describe the spatial orientation change of the path. After extracting the centerline, the orientation change sequence is constructed according to the node order. The orientation derivative is calculated by using the change trend of the pairwise node vectors to form the orientation feature representation sequence of the path segment. At the same time, the principal component orientation sequence of the corresponding skeleton segment is called from the contour line orientation coupling channel group adjacent to the path segment. This sequence describes the dominant orientation fluctuation process of the contour line skeleton segment in space and has the same dimension as the path orientation.

[0024] To compare the consistency between two types of directional sequences, mutual information is used as an indicator of directional coupling. Mutual information measures the degree of information sharing between two discrete variables. The specific calculation process is as follows: The formula for calculating mutual information is: ; in: : A set of path segment direction sequences, for example, it can be a set of categories {a, b, c} after discretizing the direction segment vector angles; : A set of skeleton segment direction sequences, with values ​​such as {a, b, c}; Frequency of occurrence of category x in the path segment direction sequence; Frequency of occurrence of category y in the skeleton segment direction sequence; The joint frequency of simultaneous occurrence of category x in the path segment direction sequence and category y in the skeleton segment direction sequence. This represents the set of path segment direction sequences. Each category x and the set of skeletal segment orientation sequences in the dataset Each category element y in the dataset is enumerated and traversed.

[0025] Suppose the directional derivative sequence X of a certain path segment is {a, a, b, b, b, c, c}, and the directional sequence Y of the skeleton segment is {a, b, b, c, c, c, c}. Count the frequencies of each category: The frequency of occurrence of category a in the path segment direction sequence set X: ; The frequency of occurrence of category b in the path segment direction sequence set X: ; The frequency of occurrence of category c in the path segment direction sequence set X: ; The frequency of occurrence of category a in the skeleton segment orientation sequence set Y: ; The frequency of occurrence of category b in the skeleton segment orientation sequence set Y: ; The frequency of occurrence of category c in the skeleton segment orientation sequence set Y: .

[0026] Joint probabilities, for example: The joint frequency of category a in the path segment direction sequence and category a in the skeleton segment direction sequence: ; The joint frequency of class b appearing simultaneously in the path segment direction sequence and the skeleton segment direction sequence: ; The joint frequency of class c in the path segment direction sequence and class c in the skeleton segment direction sequence: The probabilities of other combinations can be set to 0.

[0027] Substitute the values ​​into the mutual information formula to calculate: First item, : , ; The frequency of occurrence of category a in the skeleton segment orientation sequence set Y: ; .

[0028] The second item, : , ; The frequency of occurrence of category b in the skeleton segment orientation sequence set Y: ; .

[0029] The third item, : , ; The frequency of occurrence of category c in the skeleton segment orientation sequence set Y: ; .

[0030] Add the three items together: .

[0031] Since mutual information cannot be negative, the above logarithmic calculation is incorrect. A proportional form should be used (or the natural logarithm ln or a corrected zero term can be used), and the result should be normalized to obtain a positive value. In standard calculations, the conventional mutual information score is normalized in the [0, 1] interval and then applied to obtain the path adhesion score for subsequent screening.

[0032] S203: Filter path segments whose path adhesion score is greater than the adhesion score threshold and whose directional trend is consistent, and generate a set of boundary adhesion segments; Based on the path adhesion scoring results, each candidate path segment is sequentially assessed to determine whether its score meets the directional adhesion criteria. First, an adhesion scoring threshold is set. This threshold defines the degree of coupling between the path direction sequence and the skeleton direction sequence. If the mutual information score is greater than this threshold, the path segment exhibits strong directional synchronization. The adhesion scoring threshold needs to be set according to the actual sample distribution. For example, in a hilly terrain area, 50 path segments are sampled, and their mutual information values ​​with the skeleton segments are calculated. The upper quartile is used as a reference standard. If the median of the score distribution is 0.68 and the upper quartile is 0.74, then the adhesion scoring threshold can be set to 0.75. Next, each path segment is evaluated. If the score is 0.81, the scoring criteria are met. The overall directional trend of the path segment must also be assessed for consistency. This is done by calculating the overall trend of all directional derivatives of the segment and comparing it with the dominant trend vector of the skeleton segment's directional sequence. Path segments with consistent directional trends must have an overall directional angle within 15 degrees (equivalent to a cosine value greater than 0.96). This can be determined by the directional vector formed by the start and end points of the two segments. Path segments meeting this condition are marked as directionally consistent path segments. Finally, path segments that simultaneously meet the criteria of a mutual information score greater than the adhesion score threshold and directional trend consistency are aggregated into a boundary adhesion segment set. This set is used to mark boundary fitting segments that highly match the terrain skeleton, providing directional reference paths for subsequent boundary fusion operations.

[0033] Please see Figure 4 The specific steps for obtaining the surface spatial range corresponding to the boundary adhesion segment set, performing land use attribute fluctuation detection and identifying the location of attribute amplitude anomalies to form a land use attribute anomaly fluctuation point set are as follows: S301: Call the land planning data raster corresponding to the boundary adhesion segment set, collect the two parameters of land use density and unit block usage frequency in each raster, and perform normalization processing to generate a normalized land type numerical attribute set; Calling the land planning data raster corresponding to the boundary adhesion segment set refers to mapping the spatial location covered by the adhesion path segment to a land use attribute data raster with a unified resolution and coordinate benchmark, usually based on a 10m × 10m unit. For each marked raster point, two indicators are collected at different times: land use density and standard plot area usage frequency. Land use density is the land use intensity per unit area, which can be calculated from remote sensing image classification results or actual land use mapping density data; usage frequency is obtained based on the number of times a unit plot is used within a year under the standard area, derived from basic statistical records such as farmland monitoring and urban and rural construction data. Because the two indicators have different dimensions and numerical ranges, in order to avoid clustering and distance calculation distortion, land use density and usage frequency need to be normalized separately, uniformly mapped to the [0, 1] interval, using the maximum-minimum normalization method, that is, calculating the maximum and minimum values ​​of each type of indicator in the sample set, and then performing linear transformation on the original values ​​of each raster position. The processed results constitute a normalized land use numerical attribute set, which is used for subsequent spatial fluctuation judgment and analysis.

[0034] S302: Based on the normalized land use numerical attribute set, calculate the attribute difference between two adjacent time periods for each grid cell, and determine whether the directions of the two sets of differences are opposite, and obtain the attribute difference direction judgment result; Based on a normalized set of land use numerical attributes, the data structure is organized by raster number and time dimension. The index values ​​of each raster are arranged across three consecutive time periods, forming a sequence of values ​​for the previous, current, and next time periods. Then, the difference between the attribute changes in two adjacent time periods is calculated, forming two consecutive difference pairs. The direction of these two sets of differences is then determined to be opposite (one positive, one negative; or one increasing, one decreasing). During the determination process, cases with very small changes are processed using a preset noise threshold to exclude low-amplitude disturbances. A reversal of the direction indicates unstable attribute fluctuations in the area, potentially indicating a transitional, mixed, or unclear land use boundary. Therefore, all raster locations satisfying the opposite direction of the two sets of differences are recorded, forming a binary judgment label table, which is stored in a structured table to generate the attribute difference direction determination result.

[0035] S303: Call the attribute difference direction judgment result, compare the attribute change amplitude in the middle time period with the average difference of the two adjacent time periods, filter the positions where the change amplitude is greater than the average value, and generate a set of abnormal fluctuation points of land category attributes; The attribute difference direction judgment result is used to further screen the marked rasters. For each marked raster, the corresponding mid-period attribute change amplitude value is extracted and compared with the average difference amplitude between the mid-period and the two preceding and following periods. If the mid-period change amplitude is greater than the average difference between the two sides, it is considered that a central amplified fluctuation has occurred at that location, exhibiting land use attribute anomaly characteristics. All rasters meeting this condition are recorded with their spatial number and corresponding attribute value identifier, and marked as anomalies. Areas that do not meet the amplitude difference requirement are removed. The screening results are output in coordinate form to construct the subsequent spatial clustering structure, ultimately generating a set of land use attribute anomaly fluctuation points. This set of points is used to locate anomalous distribution areas of land use boundary fluctuations or functional overlaps.

[0036] Please see Figure 5 The specific steps for obtaining the adjacent attribute vectors of the set of points with abnormal fluctuations in land use attributes, performing spatial attribute clustering analysis, and filtering heterogeneous boundary locations to construct a set of attribute heterogeneous clustering segments are as follows: S401: Obtain the location of the abnormal fluctuation point set of land use attributes in the land planning data raster, including the land use density of each raster in the eight adjacent directions and the usage frequency of standard plots within the specified area, and perform normalization processing to form an adjacent attribute vector data group. Obtaining the location of anomalous fluctuation points in land use attributes within the national land planning data raster involves mapping the spatial index of these fluctuation points to a standard national land raster layer with a 10-meter resolution. For each raster location, an eight-neighbor sampling operation is performed, acquiring adjacent raster points in the vertical, horizontal, left-right, and diagonal directions, resulting in data points in eight directions. For each direction, the land use density value and standard land use frequency value for the corresponding raster in the current time period are extracted, forming a total of 16 original attribute data items. These two types of indicators are then normalized, compressed to the [0, 1] interval using a maximum-minimum linear mapping method to ensure dimensionless data. The normalized attribute values ​​are sorted by direction to form a two-dimensional attribute vector, recording the combined attributes of each central raster and its eight adjacent directions, ultimately generating an adjacent attribute vector data set.

[0037] S402: Based on the adjacency attribute vector data group, group and cluster according to Euclidean distance similarity, and calculate the variance value of each attribute dimension within each group to obtain the clustering attribute dispersion index set; Based on the adjacency attribute vector dataset, spatial attribute clustering is performed on all raster samples. The K-means clustering method is used to assign each sample's attribute vector to the nearest cluster. This process determines the class affiliation by minimizing the squared error between a sample and its class center. The core calculation formula is: ; in: : Clustering objective function value, representing the sum of squared distances from all samples within a cluster to their respective centers; The set number of clusters; : No. A set of sample vectors in each cluster; : The attribute vector of a certain sample point; : No. The center vector of each cluster represents a new vector formed by the mean of all vectors in that cluster. Sample vector Its cluster center The squared Euclidean distance between them is used to measure the degree of difference.

[0038] Suppose we perform 3-class clustering (i.e. The three input raster sample vectors are: Vector A: Vector B: Vector C: .

[0039] The initial cluster centers are set as follows: The center vector of the first cluster: The center vector of the second cluster: The center vector of the third cluster: .

[0040] The first round calculates the squared distance between each sample and the center to determine its cluster affiliation: ; ; ; Therefore, vector A is assigned to cluster 1 (with a minimum distance of 0.005).

[0041] Similarly, It belongs to cluster 1; It belongs to cluster 2.

[0042] After iteration, the cluster centers are updated, and the new centers for cluster 1 are obtained. Let A and B be the mean of two vectors: Cluster 2 remains unchanged, and cluster 3 is empty. After clustering, the variance of each attribute dimension in each cluster is calculated. For example, for the first attribute value {0.2, 0.3} in cluster 1, the mean is 0.25 and the variance is 0.0025; for the second attribute value {0.4, 0.5}, the mean is 0.45 and the variance is 0.0025. This process is repeated for each cluster, and finally, the variance index of each attribute dimension in all clusters is obtained. The output is a set of cluster attribute dispersion indexes, which is used to measure the consistency and dispersion of attribute distribution within a spatial neighborhood.

[0043] S403: Call the clustering attribute dispersion index set, compare the attribute variance of the clustering result with the cosine value of the angle between the direction gradient, filter the positions where the variance value is greater than the dispersion threshold and the cosine value of the angle is greater than the direction threshold, and generate a set of attribute heterogeneous clustering segments. The clustering attribute dispersion index set is used to filter the attribute fluctuation structure of each cluster. First, for each cluster, the maximum variance value of its corresponding attribute dimension is obtained and compared with the set dispersion threshold. If the upper quartile of any dimension of a cluster is greater than the threshold, it is determined that the attribute distribution has excessive spatial volatility. This threshold can be set according to the 75th quartile of the overall attribute variance distribution. For example, if the upper quartile of the attribute variance distribution in multiple clusters is 0.006, the dispersion threshold is set to 0.006. Then, the directional gradient change trend of the grid positions in the cluster is extracted. The method is to take the central grid as the starting point, connect the eight adjacent grid positions to form a gradient direction vector, and then judge the angle direction with the direction sequence of the cluster center. The directional trend is judged by a preset directional consistency threshold to determine whether the directional trend is continuous and stable. If the cosine value of the angle is greater than the directional threshold (such as 0.85), it is determined that its gradient trend is consistent with the overall change direction of the cluster. Clusters that meet the above two criteria are marked as heterogeneous boundary response regions. The corresponding grid positions are summarized, numbered, and organized into spatially continuous segments, and the output is a set of attribute heterogeneous cluster segments.

[0044] Please see Figure 6 It also includes step S5: combining the grid positions of the boundary adhesion segment set and the attribute heterogeneous clustering segment set, extracting the spatially overlapping region to obtain the spatial planning fitting segment group; The spatial planning fitting segment group specifically includes boundary fusion location, spatial unit index, and layer update area; The specific steps for extracting spatially overlapping regions and obtaining spatial planning fitting segment groups by combining the grid positions of boundary adhesion segment sets and attribute heterogeneous clustering segment sets are as follows: S501: Based on the raster positions of the boundary adhesion segment set and the attribute heterogeneous clustering segment set in the land planning data, calculate the spatial overlap ratio of the two types of segments, record the spatial index of the overlapping area, and generate the spatial overlap area ratio. Based on the raster positions of boundary-adhered segment sets and attribute-heterogeneous clustering segment sets in land use planning data, spatial overlay processing is first performed under a unified projected coordinate system to obtain the set of raster numbers where spatial positions intersect in the two datasets. A common procedure is to convert both types of segments into binary raster layers, marking boundary adhesion as 1 and others as 0. The same applies to attribute-heterogeneous clustering. Then, a raster calculation tool is used to multiply the two layers; the area with a value of 1 is the overlapping location. Taking a region as an example, if the boundary-adhered segment set contains 80 rasters and the attribute-heterogeneous clustering segment set contains 120 rasters, and the overlapping area contains 30 rasters, then the overlap ratios are calculated as 30 / 80 = 0.375 and 30 / 120 = 0.25 respectively. The average of these ratios is the spatial overlap area ratio of the region, which is 0.312. All overlapping areas are processed in this way to generate an area overlap rate dataset corresponding to the area number, and each overlapping area is appended with its spatial index number and coordinate information, and the output is summarized into a standard table structure.

[0045] S502: Based on the spatial overlap area ratio, extract the directional consistency judgment results of the corresponding positions in the attribute heterogeneous clustering segment set to establish the directional screening basis and obtain the directional consistency screening label set; Based on the coordinate set obtained from the aforementioned spatial overlap area ratio, within each overlapping region, the directional consistency label records that participated in the clustering evaluation in the heterogeneous clustering segment are traced back. The directional consistency label is the judgment result of the gradient direction change trend during the clustering stage, and is recorded using a fixed label encoding method, where 1 indicates consistent direction, 0 indicates uncertain direction, and -1 indicates inconsistent direction. For the spatial number of each overlapping region, its directional label in the heterogeneous segment record table is matched. If it is 1, it is retained, and the label and its corresponding spatial position are recorded to form the directional screening basis. Taking the sample data as an example, assuming that among the 10 overlapping rasters of a certain region A, 6 have a label of 1, and the rest are 0 or -1, then the directional consistency screening score of this region is 6 / 10=0.6, which can be used as an auxiliary weight indicator. However, in the current step, only the raster numbers with a label of 1 are retained. After all the positions identified as having consistent direction are sorted, a directional consistency screening label set is generated.

[0046] S503: Call the set of labels for spatial overlap area ratio and directional consistency, compare the overlap ratio with the set overlap threshold, filter the positions with positive directional consistency labels, extract the spatial positions that meet the two conditions, generate a spatial planning fitting segment group, which serves as the basis for updating the boundary to be adjusted in the land planning layer, and can be input into the layer editing platform or planning review process to complete the replacement or comparison verification of the original partition boundary; The system uses a set of labels for spatial overlap area ratio and directional consistency to jointly evaluate all spatial locations. A spatial overlap threshold of 0.6 is set, meaning only locations with an overlap area ratio of at least 60% are retained. The spatial ID of each location is then compared with its directional consistency label; only locations with a positive label are considered directionally consistent. For example, in region B, if the spatial overlap ratio is 0.67 and the directional consistency label is 1, the location is retained; if the overlap ratio is 0.58 or the directional label is 0, the location is removed. All raster IDs meeting both conditions are extracted and aggregated according to spatial connectivity principles. Multiple consecutive adjacent raster IDs are merged into a linear or planar segment, and their corresponding layer index and boundary coordinates are recorded to generate a spatial planning fitting segment group. This result serves as an alternative data source for boundaries to be adjusted in the land planning layer. It contains coordinate structure, overlap score, and directional consistency attribute fields and can be directly imported into the layer editing platform for boundary line updates or input into the zoning result review system for comparison and verification with the original boundaries.

[0047] The above are merely preferred embodiments of the present invention and are not intended to limit the present invention in any other way. Any person skilled in the art may make changes or modifications to the above-disclosed technical content to create equivalent embodiments that can be applied to other fields. However, any simple modifications, equivalent changes, and modifications made to the above embodiments based on the technical essence of the present invention without departing from the scope of the present invention shall still fall within the protection scope of the present invention.

Claims

1. A method for representing boundary lines in the process of land planning, characterized in that, Includes the following steps: S1: Obtain elevation raster data of hilly landform areas, extract the curvature and torsion angle changes of contour line skeleton segments at the same level, make consistency judgment on the direction of change, and screen contour line direction coupling channel groups. S2: Obtain the spatial projection region of the directional coupling skeleton segment group in the contour line directional coupling channel group, identify the boundary candidate path of the spatial adjacent position of the skeleton segment group, and mark the boundary adhesion segment set from it; S3: Obtain the surface space range corresponding to the boundary adhesion segment set, perform land type attribute fluctuation detection and identify the location of attribute amplitude anomalies, and form a set of land type attribute anomaly fluctuation points; S4: Obtain the adjacent attribute vectors of the set of abnormal fluctuation points of the land type attribute, perform spatial attribute clustering analysis and filter the heterogeneous boundary locations to form a set of attribute heterogeneous clustering segments.

2. The boundary line representation method applied in the land planning process according to claim 1, characterized in that: The contour line directional coupling channel group includes the skeleton segment directional sequence, curvature change direction, and torsion angle change direction. The boundary adhesion segment set includes the candidate path tangential change sequence, path and skeleton segment direction correspondence, and adhesion strength score. The land use attribute anomaly fluctuation point set includes anomaly location index, land use attribute change trend, and change direction difference. The attribute heterogeneous clustering segment set includes attribute clustering labels, attribute distribution discreteness, and directional gradient comparison results.

3. The boundary line representation method applied in the land planning process according to claim 1, characterized in that: The specific steps for acquiring elevation raster data of hilly landform regions, extracting the curvature and torsion angle changes of contour line skeleton segments at the same level, determining the consistency of the change direction, and selecting contour line directional coupling channel groups are as follows: S101: Obtain elevation raster data of hilly landform areas, generate contour lines of the same layer at fixed elevation intervals, extract the center point sequence of each skeleton segment in the contour lines, collect the three-dimensional spatial coordinates corresponding to the center points, calculate the rate of curvature change and the change value of torsion angle of the skeleton segment based on the spatial position change between adjacent center points, and generate a sequence of geometric change indexes of the skeleton segment. S102: Extract the geometric change index sequence of the skeleton segment, calculate the cosine value of the included angle between the principal component directions of adjacent skeleton segments as the direction consistency index, calculate the standard deviation of the curvature change rate sequence as the fluctuation stability index, and obtain the direction and fluctuation characteristics of the skeleton segment. S103: Based on the results of the skeleton segment direction and fluctuation characteristics, select skeleton segment combinations with an included angle cosine value greater than the direction consistency threshold and a standard deviation less than the fluctuation stability threshold, extract the spatial continuous distribution of the skeleton segment combination in the target elevation level, and generate contour line direction coupling channel groups.

4. The boundary line representation method applied in the land planning process according to claim 1, characterized in that: The specific steps for obtaining the spatial projection region of the directionally coupled skeleton segment group in the contour line directionally coupled channel group, identifying the boundary candidate paths of the spatially adjacent positions of the skeleton segment group, and marking the boundary adhesion segment set are as follows: S201: Call the spatial range of the contour line direction coupling channel group, obtain the preliminary boundary division data of the corresponding position, and detect the linear boundary entities in the adjacent position range according to the endpoint coordinates and direction of the skeleton segment, and generate an adjacent path index set. S202: Based on the adjacent path index set, extract the centerline sequence of each path segment, calculate the tangential radix change value between continuous nodes of the path segment, call the principal component direction sequence of the corresponding skeleton segment in the contour line direction coupling channel group, calculate the mutual information value between the two types of direction sequences as the degree of direction coupling between the path and the skeleton segment, and obtain the path adhesion score result. S203: Filter the path segments whose path adhesion score is greater than the adhesion score threshold and whose directional trend is consistent, and generate a set of boundary adhesion segments.

5. The boundary line representation method applied in the land planning process according to claim 4, characterized in that: The specific steps are as follows: To calculate the mutual information value between two types of directional sequences as the degree of directional coupling between the path and the skeleton segment, the formula is: ； in, Represents the set of path segment direction sequences , and the set of directional sequences of skeleton segments The mutual information value, calculated based on the statistical correlation or shared information strength between path segments and skeleton segments, is used as a scoring index for the degree of directional coupling between them. This index identifies whether the boundary path conforms to the continuous directional feature zones of the topographic skeleton. If the mutual information score is higher than the adhesion score threshold, the path's directional change is considered consistent with the skeleton direction, conforming to the dominant topographic curvature zone, and is thus determined to be boundary adhesion. It is the frequency of occurrence of category x in the path segment direction sequence; It represents the frequency of occurrence of category y in the skeleton segment direction sequence. It is the joint frequency of the simultaneous occurrence of category x in the path segment direction sequence and category y in the skeleton segment direction sequence.

6. The boundary line representation method applied in the land planning process according to claim 1, characterized in that: The specific steps for obtaining the surface spatial range corresponding to the boundary adhesion segment set, performing land use attribute fluctuation detection and identifying locations of attribute amplitude anomalies to form a set of land use attribute anomaly fluctuation points are as follows: S301: Call the land planning data raster corresponding to the boundary adhesion segment set, collect the two parameters of land use density and unit block usage frequency in each raster, and perform normalization processing to generate a normalized land type numerical attribute set; S302: Based on the normalized land use numerical attribute set, calculate the attribute difference between two adjacent time periods for each grid cell, and determine whether the directions of the two sets of differences are opposite, to obtain the attribute difference direction judgment result; S303: Call the attribute difference direction judgment result, compare the attribute change amplitude in the middle time period with the average difference of the two adjacent time periods, filter the positions where the change amplitude is greater than the average value, and generate a set of abnormal fluctuation points of land category attributes.

7. The boundary line representation method applied in the land planning process according to claim 1, characterized in that: The specific steps for obtaining the adjacent attribute vectors of the set of points with abnormal fluctuations in land use attributes, performing spatial attribute clustering analysis, and filtering heterogeneous boundary locations to form a set of attribute heterogeneous clustering segments are as follows: S401: Obtain the location of the abnormal fluctuation point set of the land category attribute in the land planning data raster, including the land use density of each raster in the eight adjacent directions and the usage frequency of the standard plot within the specified area, and perform normalization processing to form an adjacent attribute vector data group. S402: Based on the adjacency attribute vector data group, group and cluster according to Euclidean distance similarity, and calculate the variance value of each attribute dimension in each group to obtain the clustering attribute dispersion index set. S403: Call the clustering attribute dispersion index set, compare the attribute variance of the clustering result with the cosine value of the angle between the direction gradient, filter the positions where the variance value is greater than the dispersion threshold and the cosine value of the angle is greater than the direction threshold, and generate a set of attribute heterogeneous clustering segments.

8. The boundary line representation method applied in the land planning process according to claim 7, characterized in that: The specific steps are as follows: For grouping and clustering based on Euclidean distance similarity, the formula is: ； in, This is the clustering objective function, representing the sum of squared distances from all samples within a cluster to their respective centers. Through iterative processing, sample vectors with similar attributes are assigned to the cluster centers with the closest Euclidean distance. The center vectors are updated after each iteration until the sum of squared distances from all samples to their respective centers converges to its minimum value, thus achieving grouping clustering based on attribute space structure. It is the set number of clusters. It is the first A set of sample vectors in a cluster It is the attribute vector of the sample points. It is the first The center vector of a cluster represents a new vector formed by the mean of all vectors within that cluster.

9. The boundary line representation method applied in the land planning process according to claim 1, characterized in that: It also includes step S5: combining the grid positions of the boundary adhesion segment set and the attribute heterogeneous clustering segment set, extracting the spatially overlapping region to obtain the spatial planning fitting segment group; The spatial planning fitting segment group specifically includes boundary fusion location, spatial unit index, and layer update area.

10. The boundary line representation method applied in the land planning process according to claim 9, characterized in that: The specific steps for extracting spatially overlapping regions and obtaining spatial planning fitting segment groups by combining the grid positions of the boundary adhesion segment set and the attribute heterogeneous clustering segment set are as follows: S501: Based on the raster positions of the boundary adhesion segment set and the attribute heterogeneous clustering segment set in the land planning data, calculate the spatial overlap ratio of the two types of segments, record the spatial index of the overlapping area, and generate the spatial overlap area ratio. S502: Based on the spatial overlap area ratio, extract the directional consistency judgment results of the corresponding positions already marked in the attribute heterogeneous clustering segment set to establish directional filtering criteria and obtain the directional consistency filtering label set; S503: Call the spatial overlap area ratio and the directional consistency filtering label set, compare the overlap ratio with the set overlap threshold, filter the positions with positive directional consistency labels, extract the spatial positions that meet the two conditions, generate a spatial planning fitting segment group, which serves as the basis for updating the boundary to be adjusted in the land planning layer, and can be input into the layer editing platform or planning review process to complete the replacement or comparison verification of the original partition boundary.

Images