Unmanned aerial vehicle remote sensing and laser radar collaborative farmland boundary extraction method

Abstract

This invention relates to the field of farmland remote sensing monitoring technology, and discloses a method for farmland boundary extraction using a combination of UAV remote sensing and lidar. The method first acquires a set of oblique UAV photographic images of the target area, airborne lidar point cloud data, and shooting attitude parameters to extract an initial visually recognized boundary line. Then, a discretized local polar coordinate occlusion perception model is constructed to screen boundary segments to be compensated. Subsequently, candidate canopy point clouds are extracted based on the boundary segments to be compensated, and the equivalent occlusion canopy height is calculated. A multi-view occlusion shrinkage compensation vector set is constructed to determine the visually corrected boundary points. Finally, lidar reference boundary points are extracted, and the farmland physical boundary line is reconstructed after dual-track consistency verification. This invention can effectively eliminate the systematic shrinkage deviation of the boundary caused by crop canopy occlusion under oblique perspectives, ensuring the boundary positioning accuracy in scenarios such as farmland rights confirmation.

Description

Technical Field

[0001] This invention relates to the field of farmland remote sensing monitoring technology, and more specifically, to a method for farmland boundary extraction using a combination of UAV remote sensing and lidar. Background Technology

[0002] The extraction of farmland boundaries is a core technical support for farmland rights confirmation and registration, cadastral surveying, high-standard farmland construction, and refined management of agricultural resources. Its positioning accuracy directly affects the scientific nature of farmland area calculation, land ownership division, and agricultural production planning. With the development of low-altitude remote sensing technology, UAV oblique photogrammetry and airborne lidar technology have become the mainstream technical means for farmland boundary extraction due to their high-resolution and high-timeliness ground feature observation capabilities. The synergistic application of the two is considered a key path to overcome the limitations of single data sources and improve the accuracy of boundary extraction.

[0003] In the field of remote sensing monitoring of arable land, many technical solutions have been explored for the extraction of arable land information. For example, Chinese patent application CN120411786A discloses a method for cross-measurement of farmland shading by fusing UAV laser point clouds and imagery. This method first collects farmland point cloud data and farmland image data and fuses them to generate fused point cloud data. Then, it obtains second green area data by identifying first green area data and calculating the visible light vegetation index. After integration, it forms green area point cloud data and performs regionalization processing. Subsequently, it removes floating areas through vertical structure analysis to obtain farmland area data, and finally generates orthophotos of farmland areas. This method realizes the identification and removal of vegetation shading in farmland and effectively avoids misjudgment of undulating terrain or densely vegetated areas. Chinese patent application CN120913078A discloses a method and system for farmland segmentation in remote sensing images that integrates context and boundary awareness. It constructs a farmland segmentation model that includes components such as a backbone network and feature enhancement modules. Through iterative optimization of a composite loss function, it realizes the segmentation of farmland in remote sensing images, strengthens the model's ability to distinguish farmland boundaries, and is suitable for agricultural interpretation and monitoring scenarios.

[0004] However, the aforementioned existing technologies still have significant technical limitations and are unable to solve the problem of systematic boundary shrinkage under the tilted perspective of tall crop fields. Specifically, while existing technology CN120411786A can remove shading from cultivated land, its core objective is the purification of green space data. It does not establish a compensation model for the boundary geometric offset caused by crop canopy during oblique photography, and cannot correct the systematic error of the boundary "retreating" into the field due to canopy projection shading. Existing technology CN120913078A focuses on improving the semantic segmentation accuracy of remote sensing images, but does not combine the three-dimensional elevation information of lidar point clouds to construct the shading geometric inversion relationship, and does not consider the differences in shading offset under different shooting perspectives. At the same time, traditional remote sensing-lidar fusion... Many fusion methods employ simple weighted averages or confidence level fusion strategies, completely ignoring the geometric projection nature of canopy shading under tilted perspectives. This offset is a systematic error strongly correlated with crop height, shooting pitch angle, and boundary orientation azimuth angle, rather than random noise. If the correlation between the shading offset and the perspective, crop height, and boundary orientation is ignored, the fused boundary will still have an inward shrinkage deviation in tall crop areas, failing to meet the requirements for centimeter-level to decimeter-level boundary accuracy for farmland ownership confirmation. This leads to inaccurate positioning of farmland ownership boundary points and deviations in area calculation, affecting the authority of land ownership confirmation and the accuracy of agricultural resource management. Summary of the Invention

[0005] To overcome the aforementioned shortcomings of existing technologies, this invention provides a method for extracting farmland boundaries using a combination of UAV remote sensing and lidar. It achieves precise selection of boundary segments to be compensated by constructing a discretized local polar coordinate occlusion perception model. By combining lidar point cloud computing with the equivalent occlusion canopy height and constructing a multi-view occlusion shrinkage compensation vector set to determine visually corrected boundary points, data fusion and boundary reconstruction are completed through dual-track boundary consistency verification. This invention effectively eliminates the systematic shrinkage deviation of boundaries caused by crop canopy projection occlusion from tilted perspectives, achieving the extraction of farmland physical boundaries. It improves the accuracy and reliability of boundary positioning in tall crop areas, providing precise technical support for farmland rights confirmation and other related operations.

[0006] To achieve the above objectives, the present invention provides the following technical solution:

[0007] A method for extracting farmland boundaries using a combination of UAV remote sensing and lidar includes:

[0008] Acquire the UAV oblique photogrammetry image set, airborne LiDAR point cloud data and UAV shooting attitude parameters of the target area. Extract the initial visual recognition boundary line from the UAV oblique photogrammetry image set. Construct a discretized local polar coordinate occlusion perception model based on the initial visual recognition boundary line and UAV shooting attitude parameters. Use the discretized local polar coordinate occlusion perception model to screen the boundary segments to be compensated.

[0009] Candidate canopy point sets are extracted from airborne lidar point cloud data based on the boundary segment to be compensated. The equivalent occlusion canopy height is calculated based on the candidate canopy point sets. A set of multi-view occlusion shrinkage compensation vectors is constructed based on the equivalent occlusion canopy height. The visual correction boundary points are determined based on the set of visual occlusion shrinkage compensation vectors.

[0010] Extract lidar reference boundary points from airborne lidar point cloud data, perform dual-track boundary consistency verification between visual correction boundary points and lidar reference boundary points, and reconstruct the physical boundary line of cultivated land based on the verification results.

[0011] The method for extracting the initial visual recognition boundary line includes:

[0012] Orthorectification and semantic segmentation are performed on the UAV oblique photography image set to extract the boundary pixel set between cultivated land and non-cultivated land, and curve fitting is performed on the boundary pixel set to generate the initial visual recognition boundary line.

[0013] The method for constructing a discretized local polar coordinate occlusion perception model and using the discretized local polar coordinate occlusion perception model to filter the boundary segments to be compensated includes:

[0014] The initial visual recognition boundary line is subjected to equidistant divergence processing to generate a sequence of discrete boundary segments. For each discrete boundary segment in the sequence, a local polar coordinate system is constructed and occlusion condition determination is performed. Discrete boundary segments that meet the occlusion condition are marked as boundary segments to be compensated.

[0015] The method for generating the discrete boundary segment sequence includes:

[0016] A sequence of discrete boundary sampling points is obtained by performing equidistant divergent sampling along the initial visual recognition boundary line. Adjacent discrete boundary sampling points in the sequence are then connected sequentially to construct a sequence of discrete boundary segments.

[0017] The method for constructing a local polar coordinate system includes:

[0018] Select the target boundary segment from the discrete boundary segment sequence, calculate the tangential direction of the target boundary segment, and establish a local polar coordinate system with the midpoint of the target boundary segment as the origin and the tangential direction as the polar axis.

[0019] The drone shooting attitude parameters include at least the camera spatial position coordinates and the shooting pitch angle;

[0020] The method for determining the occlusion condition includes:

[0021] In the local polar coordinate system, define the inner positive normal vector pointing into the field and the outer normal vector pointing out of the field.

[0022] The camera line-of-sight vector is constructed based on the camera spatial position coordinates and the midpoint coordinates of the target boundary segment in the drone's shooting attitude parameters. The camera line-of-sight vector is then projected onto the horizontal plane of the local polar coordinate system to calculate the line-of-sight azimuth angle.

[0023] Determine whether the target boundary segment meets the occlusion condition based on the line of sight azimuth, and mark the target boundary segment that meets the occlusion condition as the boundary segment to be compensated.

[0024] The criterion for determining whether the occlusion condition is met is that the absolute value of the line-of-sight azimuth angle is less than 90°.

[0025] The method for extracting the candidate canopy point cloud includes:

[0026] Extend the boundary segment to be compensated along the direction of the inner normal vector of the field to construct a long strip boundary buffer zone. Select the point cloud data falling into the long strip boundary buffer zone from the airborne lidar point cloud data to form a candidate canopy point cloud set.

[0027] The method for calculating the equivalent obstructed canopy height includes:

[0028] Identify low-elevation points in the candidate canopy point cloud set, remove low-elevation points from the candidate canopy point cloud set, and retain the remaining point cloud as canopy points. Calculate the vertical distance from each retained canopy point to the boundary segment to be compensated, and use an inverse weighting function to calculate the equivalent shading canopy height of the boundary segment to be compensated based on the vertical distance.

[0029] The method for constructing the multi-view occlusion shrinkage compensation vector set includes:

[0030] By combining the shooting pitch angle, line of sight azimuth angle and equivalent shading canopy height in the drone shooting attitude parameters, the shading inward distance is calculated, and a single-view shading inward compensation vector with a modulus of the shading inward distance is constructed along the normal vector direction of the outer side of the field.

[0031] For the same boundary segment to be compensated, multiple single-view occlusion shrinkage compensation vectors from different images are collected to form a multi-view occlusion shrinkage compensation vector set.

[0032] The dual-track boundary consistency check includes:

[0033] The spatial Euclidean distance deviation between the visually corrected boundary points and the lidar reference boundary points is calculated. When the spatial Euclidean distance deviation is greater than the preset consistency verification threshold, a conflict arbitration mechanism based on local point cloud density is triggered to obtain dynamic fusion weights. The visually corrected boundary points and lidar reference boundary points are weighted using the dynamic fusion weights to obtain the final fused boundary points. Otherwise, equal-weight fusion is performed on the visually corrected boundary points and lidar reference boundary points to obtain the final fused boundary points.

[0034] Compared with the prior art, the beneficial effects of the present invention are as follows:

[0035] This invention accurately identifies boundary segments with canopy occlusion by constructing a discretized local polar coordinate occlusion perception model. Based on the equivalent occlusion canopy height, it calculates a multi-view occlusion shrinkage compensation vector, achieving geometric compensation for the systematic shrinkage of the boundary caused by the "projection occlusion" of the crop canopy from an inclined perspective. At the same time, by using a dual-track boundary consistency verification between visually corrected boundary points and lidar reference boundary points, and combining local point cloud density for dynamic weight allocation, it effectively eliminates the boundary positioning error caused by ignoring systematic geometric offset in traditional methods, making the extracted farmland boundary closer to the real physical boundary, and improving the accuracy and reliability of farmland boundary extraction in tall crop areas. Attached Figure Description

[0036] To more clearly illustrate the technical solutions in the embodiments of the present invention or the prior art, the drawings used in the description of the embodiments or the prior art will be briefly introduced below. Obviously, the drawings described below are only some embodiments of the present invention. For those skilled in the art, other drawings can be obtained based on these drawings without creative effort.

[0037] Figure 1 This is a flowchart of a method for extracting farmland boundaries using a combination of UAV remote sensing and lidar, provided in an embodiment of the present invention.

[0038] Figure 2 A schematic diagram illustrating the construction of a local polar coordinate system and the definition of a normal vector for a discrete boundary segment provided in an embodiment of the present invention;

[0039] Figure 3 This is a schematic diagram of dual-track boundary consistency verification and fusion provided in an embodiment of the present invention;

[0040] Figure 4 This is a functional block diagram of the farmland boundary extraction system that combines UAV remote sensing and lidar, provided in an embodiment of the present invention. Detailed Implementation

[0041] The technical solutions of the embodiments of the present invention will be clearly and completely described below with reference to the accompanying drawings. Obviously, the described embodiments are only some embodiments of the present invention, and not all embodiments. Based on the embodiments of the present invention, all other embodiments obtained by those skilled in the art without creative effort are within the scope of protection of the present invention.

[0042] Example 1

[0043] Please see Figure 1As shown, this embodiment provides a method for extracting farmland boundaries using a combination of UAV remote sensing and lidar, including:

[0044] Step S10: Obtain the UAV oblique photography image set, airborne LiDAR point cloud data and UAV shooting attitude parameters of the target area; extract the initial visual recognition boundary line from the UAV oblique photography image set; construct a discretized local polar coordinate occlusion perception model based on the initial visual recognition boundary line and UAV shooting attitude parameters; and use the discretized local polar coordinate occlusion perception model to screen the boundary segments to be compensated.

[0045] Further, step S10 includes:

[0046] Step S11: Obtain the UAV oblique photography image set, airborne LiDAR point cloud data and UAV shooting attitude parameters of the target area; perform orthorectification and semantic segmentation on the UAV oblique photography image set; extract the boundary pixel set between cultivated land and non-cultivated land; and perform curve fitting on the boundary pixel set to generate the initial visual recognition boundary line; the UAV shooting attitude parameters include at least the camera spatial position coordinates and shooting pitch angle.

[0047] Specifically, the UAV oblique photography image set refers to a series of digital image data with different shooting angles collected by a UAV equipped with a multi-angle camera array flying along a preset flight path over a target cultivated area. These images include images taken from the orthogonal direction and images taken at an oblique angle relative to the ground. The oblique angle images can capture the lateral texture information of ground objects, but when shooting areas with tall crops, the canopy will cause visual occlusion of the ground boundary. Airborne LiDAR point cloud data refers to a set of three-dimensional spatial points obtained by a laser scanning device installed on the UAV platform by emitting laser pulses and receiving reflected signals. Each point contains planar coordinates and elevation values. Laser pulses have the ability to penetrate the vegetation canopy and can obtain the elevation information of the ground below the canopy, but the penetration effect is affected by the vegetation density. UAV shooting attitude parameters refer to the flight state data recorded by the UAV when collecting each frame of image. The camera spatial position coordinates represent the three-dimensional position of the camera's optical center in the geographic coordinate system, and the shooting pitch angle represents the angle between the camera's optical axis and the horizontal plane. The shooting pitch angle is usually positive. The smaller the value, the closer the camera is to shooting horizontally, and the larger the value, the closer the camera is to shooting vertically downwards.

[0048] Orthorectification is an image processing technique that eliminates geometric distortions in images. Original oblique images, due to factors such as camera tilt, terrain undulations, and lens distortion, exhibit discrepancies between the geometric positions of features in the image and their actual geographical locations. Orthorectification uses a digital elevation model and camera intrinsic and extrinsic parameters to geometrically transform the image, ensuring that the processed image conforms to the principles of orthophoto projection, with each pixel's ground position matching its map coordinates. Using orthorectification as a preprocessing method unifies multiple images from different shooting angles into a common geographic reference frame, establishing a geometric benchmark for subsequent semantic segmentation and boundary information extraction, and eliminating positional inconsistencies caused by viewing angle differences between different images.

[0049] Semantic segmentation is a deep learning-based image analysis technique that can classify each pixel in an image, dividing the entire image into regions with different semantic meanings. In the scenario of farmland boundary extraction, after training with a large number of labeled samples, the semantic segmentation model can identify pixels belonging to farmland and non-farmland areas in the image. Non-farmland areas include land features such as roads, canals, woodlands, and buildings. The boundary pixel set refers to the pixels located at the boundary between farmland and non-farmland areas. These pixels are represented in the semantic segmentation results as locations where the categories of adjacent pixels change. The boundary pixel set is extracted using an edge detection algorithm. Curve fitting is performed on the boundary pixel set. Parametric curve equations are used to mathematically model discrete edge pixels. Fitting algorithms include polynomial fitting and B-spline fitting. The fitting process aims to minimize the sum of squares of the distances between the curve and the edge pixels, generating a continuous and smooth initial visual recognition boundary line.

[0050] The initial visual boundary line reflects the projected outline of the crop canopy on the ground, rather than the actual physical boundary of the cultivated land. When a drone photographs tall crop fields from an oblique angle, the projected position of the top of the crop canopy relative to the roots in the image shifts along the line of sight. This causes the boundary position identified in the image to shift inward relative to the actual boundary where the crop roots are located, and this shift is directly related to the crop height and the shooting tilt angle. Using oblique drone imagery as the data source for boundary recognition allows for the acquisition of rich texture details and spectral information of ground features. Semantic segmentation can accurately distinguish cultivated and non-cultivated areas at the pixel level, providing a high-resolution initial boundary outline for subsequent geometric compensation calculations. Simultaneously, acquiring airborne LiDAR point cloud data provides three-dimensional spatial data support for subsequent steps to calculate the equivalent occlusion canopy height and extract LiDAR reference boundary points. Recording the drone's shooting attitude parameters establishes a correspondence between each image frame and the camera's spatial position and attitude angle, enabling subsequent steps to independently calculate occlusion geometric parameters for each image, laying the data foundation for multi-view occlusion compensation.

[0051] Step S12: Perform equidistant divergent sampling along the initial visual recognition boundary line to obtain a sequence of discrete boundary sampling points. Connect adjacent discrete boundary sampling points in the sequence of discrete boundary sampling points in sequence to construct a sequence of discrete boundary segments.

[0052] Specifically, equidistant discretization sampling refers to the process of selecting sampling points at equal intervals along the initial visual recognition boundary line, according to a preset sampling step size. The sampling step size is determined comprehensively based on the local curvature characteristics of the farmland boundary and the computational accuracy requirements for crop canopy shading compensation: an excessively large sampling step size will cause the discretized straight line segments to fail to accurately approximate the curvature of the original curve, resulting in approximation errors in areas with large boundary curvature; an excessively small sampling step size will cause the number of sampling points to increase dramatically, increasing the complexity of subsequent calculations. For example, the sampling step size can be set to L, with a value ranging from 0.5 meters to 2 meters. The specific value is determined based on the minimum value of the boundary curvature radius. When the minimum boundary curvature radius is Rmin, the sampling step size L should satisfy the condition that L is less than the sine value of Rmin multiplied by a preset angle threshold. The preset angle threshold represents the upper limit of the tangential angle change between adjacent sampling points. For example, the preset angle threshold can be 10°. The boundary discrete sampling point sequence is a set of sampling points arranged sequentially along the initial visual recognition boundary line. Each sampling point has definite planar coordinates, and the arc distance between adjacent sampling points is equal to the sampling step size L. The sampling process starts from the beginning of the initial visual recognition boundary line and sequentially determines the position of each sampling point along the curve extension direction until the end of the boundary curve is reached. Adjacent boundary discrete sampling points in the boundary discrete sampling point sequence are connected sequentially to construct a discrete boundary segment sequence composed of multiple straight line segments connected end to end. Each straight line segment in the discrete boundary segment sequence is called a discrete boundary segment, and each discrete boundary segment is represented by a directed line segment pointing from the previous sampling point to the next sampling point. The discrete boundary segment sequence preserves the topological connection relationship and direction order of the initial visual recognition boundary line, and any two adjacent discrete boundary segments share an endpoint.

[0053] By employing equidistant discrete sampling, the continuous curve is transformed into a sequence of discrete line segments. This transforms the nonlinear occlusion calculation problem, which originally required solving on a curved surface, into a linear problem solved on planar line segments. Each discrete boundary segment can be independently analyzed for occlusion geometry. Within the scale range of a single discrete boundary segment, the curvature variation of the curve is sufficiently small, allowing the boundary of that segment to be considered a locally straight edge. This approximation has sufficient geometric accuracy when the sampling step size is sufficiently small. The construction of the discrete boundary segment sequence decomposes the complex global boundary occlusion problem into multiple independent local boundary occlusion problems. Each discrete boundary segment can establish its own local coordinate system and independently calculate the occlusion compensation. The calculation processes between the discrete boundary segments do not interfere with each other, reducing the coupling complexity of the algorithm. Step S12 transforms the initial visual recognition boundary line in the form of a continuous curve output in step S11 into the form of a sequence of discrete boundary segments, creating the conditions for segmented processing in step S13 for establishing a local polar coordinate system for each discrete boundary segment. The discrete boundary segment sequence provides a clear computational unit for subsequent steps S14 to determine the occlusion conditions and S20 to calculate the occlusion compensation vector, enabling the occlusion perception model to perform refined analysis segment by segment and avoiding the problem of inconsistent occlusion directions caused by changes in boundary orientation during overall calculation.

[0054] Step S13: Select the target boundary segment in the discrete boundary segment sequence, calculate the tangential direction of the target boundary segment, establish a local polar coordinate system with the midpoint of the target boundary segment as the origin and the tangential direction as the polar axis, and define the inner positive normal vector pointing to the inside of the field and the outer normal vector pointing to the outside of the field in the local polar coordinate system.

[0055] Specifically, see Figure 2The target boundary segment refers to the discrete boundary segment currently undergoing occlusion analysis selected from the discrete boundary segment sequence. Each discrete boundary segment in the sequence is treated as a target boundary segment during occlusion analysis, and the occlusion analysis process traverses all discrete boundary segments in the sequence. The tangential direction refers to the extension direction of the target boundary segment on the horizontal plane, represented by a unit direction vector formed by the starting point and ending point of the target boundary segment. The tangential direction is parallel to the straight line segment containing the target boundary segment. To calculate the tangential direction of the target boundary segment, the plane coordinates of the ending point are subtracted from the plane coordinates of the starting point to obtain the direction vector. This direction vector is then normalized to obtain the unit tangential vector, which has a magnitude of 1 and a direction consistent with the orientation of the target boundary segment. The local polar coordinate system is a two-dimensional coordinate system established with the target boundary segment as a reference. The origin is set at the midpoint of the target boundary segment, and the plane coordinates of the midpoint are equal to the arithmetic mean of the plane coordinates of the starting and ending points of the target boundary segment. The polar axis direction is set to the tangential direction of the target boundary segment. The polar axis extends along the tangential direction, and the positive polar axis direction is consistent with the tangential vector direction. Any point in the local polar coordinate system can be represented by two parameters: the polar radius and the polar angle. The polar radius is the distance from the point to the origin, and the polar angle is the angle between the azimuth vector of the point relative to the origin and the positive polar axis direction.

[0056] In a local polar coordinate system, define the inner positive normal vector and the outer normal vector of the field. Both vectors are perpendicular to the tangential direction of the target boundary segment. The inner positive normal vector is a unit vector pointing from the midpoint of the target boundary segment into the inner region of the cultivated field, while the outer normal vector is a unit vector pointing from the midpoint of the target boundary segment into the outer region of the cultivated field. The inner and outer normal vectors are in opposite directions, and their sum is zero. The method for determining the inner positive normal vector is as follows: Based on the tangential vector of the target boundary segment, rotate the tangential vector counterclockwise by 90° around the vertical axis to obtain a candidate normal vector. Use the position of the geometric center point of the cultivated field area enclosed by the initial visual recognition boundary line to determine whether the candidate normal vector points into the field. If the angle between the candidate normal vector and the vector pointing from the midpoint of the target boundary segment to the geometric center point is less than 90°, then the candidate normal vector is the inner positive normal vector; otherwise, reverse the candidate normal vector to obtain the inner positive normal vector. The outer normal vector of a field is the reverse of the inner normal vector. A local polar coordinate system is used instead of the global geographic coordinate system for occlusion geometry analysis. This allows occlusion determination for each target boundary segment to be performed within a local coordinate framework with itself as the reference, transforming the spatial relationship between the UAV's viewing angle and the boundary orientation into an angular relationship in the local polar coordinate system. Boundary segments with different orientations have different azimuth angles in the global geographic coordinate system. If occlusion conditions are uniformly calculated in the global coordinate system, the combination of the boundary azimuth angle and the line-of-sight azimuth angle must be considered simultaneously, making the calculation process complex. The local polar coordinate system fixes the tangential direction of the target boundary segment as the polar axis direction, normalizing the boundary orientation to a zero-degree direction in the local coordinate system. The determination of occlusion conditions only needs to consider the polar angle value of the line-of-sight vector in the local coordinate system, simplifying the mathematical expression of the occlusion conditions.

[0057] The definitions of the inner and outer normal vectors of the field provide a directional reference for subsequent steps in determining the relative position of the line of sight to the boundary segment. When the angle between the projection of the line of sight vector in the local polar coordinate system and the outer normal vector of the field is less than 90°, it indicates that the camera is observing inward from outside the field. In this case, the crop canopy is located between the line of sight and the actual boundary, resulting in an occlusion effect. When the angle between the projection of the line of sight vector in the local polar coordinate system and the inner normal vector of the field is less than 90°, it indicates that the camera is observing outward from inside the field. In this case, the actual boundary is located between the line of sight and the crop canopy, and no occlusion effect occurs. The inner normal vector of the field also provides an extension direction for constructing the elongated boundary buffer zone in step S21, and the outer normal vector of the field provides a compensation direction for constructing the single-view occlusion shrinkage compensation vector in step S23. Step S13 establishes a local polar coordinate system using the target boundary segment in the discrete boundary segment sequence output in step S12 as the object, providing a coordinate reference frame for projecting the camera line of sight vector onto the local coordinate system and calculating the line of sight azimuth in step S14. Without the establishment of the local coordinate system in step S13, step S14 cannot quantify the relative angular relationship between the viewing direction and the boundary orientation, and the determination of occlusion conditions will lack a geometric reference. The establishment of the local polar coordinate system allows boundary segments with different orientations to use a unified angle threshold to determine occlusion conditions in their respective local coordinate systems, avoiding the problem of inconsistent occlusion determination standards caused by differences in boundary orientation.

[0058] Step S14: Construct a camera line-of-sight vector based on the camera spatial position coordinates and the midpoint coordinates of the target boundary segment in the UAV shooting attitude parameters. Project the camera line-of-sight vector onto the horizontal plane of the local polar coordinate system to calculate the line-of-sight azimuth angle. Determine whether the target boundary segment meets the occlusion condition based on the line-of-sight azimuth angle. Mark the target boundary segment that meets the occlusion condition as the boundary segment to be compensated.

[0059] Specifically, the camera line-of-sight vector is a three-dimensional spatial vector representing the direction from the camera's optical center to the midpoint of the target boundary segment. When constructing the camera line-of-sight vector, the camera's spatial position coordinates are read from the UAV shooting attitude parameters obtained in step S11. These coordinates include the horizontal and vertical components of the camera's optical center in the geographic coordinate system. Subtracting the camera's spatial position coordinates from the midpoint coordinates of the target boundary segment yields a three-dimensional direction vector pointing from the camera to the midpoint of the boundary segment; this direction vector is the camera line-of-sight vector. The three components of the camera line-of-sight vector represent the projection lengths of the line of sight in the east, north, and vertical directions in the geographic coordinate system, respectively. The camera line-of-sight vector is projected onto the horizontal plane of the local polar coordinate system. The projection process retains the horizontal component of the camera line-of-sight vector while discarding the vertical component, resulting in the horizontal projection vector of the camera line-of-sight vector. The direction of the horizontal projection vector represents the camera's orientation relative to the midpoint of the target boundary segment on the horizontal plane. When calculating the line-of-sight azimuth, the angle between the horizontal projection vector and the tangential direction of the target boundary segment is obtained. The angle is calculated using the vector dot product formula: the line-of-sight azimuth equals the dot product of the horizontal projection vector and the tangential vector divided by the product of the magnitudes of the two vectors, then taking the inverse cosine function. The range of the line-of-sight azimuth is 0 to 180 degrees. To distinguish whether the line of sight is located inside or outside the field boundary segment, the sign of the dot product between the horizontal projection vector and the normal vector to the outside of the field needs to be determined: a positive dot product indicates the line of sight is incident from the outside of the field, while a negative dot product indicates the line of sight is incident from the inside of the field. The line-of-sight azimuth is then used to determine whether the target boundary segment meets the occlusion conditions. The occlusion condition is determined by whether the line-of-sight azimuth is within a preset effective occlusion range. The determination of the effective occlusion range is based on the geometric premise of occlusion: when the camera is located outside the field boundary segment and observing from the inside, the crop canopy lies between the camera's line of sight and the actual ground boundary, blocking the line of sight from reaching the ground boundary and creating a visual occlusion effect. When the camera is located inside the field boundary segment and observing from the outside, the actual ground boundary lies within the range directly observable by the camera's line of sight, and the canopy lies in the extension direction of the line of sight, not between the line of sight and the boundary, thus no visual occlusion effect occurs. The effective occlusion range is set when the absolute value of the line of sight azimuth is less than 90°. When the absolute value of the line of sight azimuth is less than 90°, the target boundary segment is considered to meet the occlusion condition; when the absolute value of the line of sight azimuth is greater than or equal to 90°, the target boundary segment is considered not to meet the occlusion condition.

[0060] Target boundary segments that meet the occlusion conditions are marked as boundary segments to be compensated. These boundary segments are those with canopy occlusion geometry under the current image capture viewpoint, and their occlusion compensation needs to be calculated in subsequent steps. Target boundary segments that do not meet the occlusion conditions do not require compensation calculation, as the positions of these boundary segments identified in the current image do not have a systematic offset from their true boundary positions due to canopy occlusion. For each discrete boundary segment in the discrete boundary segment sequence, each frame of the UAV oblique photogrammetry image set is traversed. Using the UAV capture attitude parameters corresponding to each frame, the line-of-sight azimuth angle calculation and occlusion condition determination are repeatedly performed to establish an occlusion relationship mapping table between discrete boundary segments and images, recording in which images each discrete boundary segment is marked as a boundary segment to be compensated. The line-of-sight azimuth angle is used to determine the occlusion condition, transforming the occurrence of canopy occlusion into a quantifiable angle threshold comparison problem, giving the occlusion determination process a clear mathematical expression and verifiable geometric basis. The effective shading range is set at ±90°. This threshold stems from the geometrically necessary condition for shading: the canopy can only shade the boundary when the line of sight is incident from outside the field. An absolute value of the line of sight azimuth less than 90° means that the main direction of the horizontal component of the line of sight points inwards into the field. By filtering boundary segments to be compensated based on shading conditions, compensation calculations are performed only on the necessary boundary segments in subsequent steps. This avoids invalid calculations for boundary segments without shading effects and also avoids applying incorrect compensation offsets to boundary segments without shading. Step S14 uses the local polar coordinate system established in step S13 to convert the camera line of sight vector into the line of sight azimuth. The filtered boundary segments to be compensated are then used in step S20 to extract candidate canopy point clouds and calculate shading compensation vectors. Without the shading condition determination in step S14, step S20 would be unable to distinguish which boundary segments need compensation and which do not, potentially applying incorrect offsets towards the outside of the field to boundary segments without shading, causing the boundary extraction results to deviate from the true position. The boundary segment marker to be compensated and the corresponding line-of-sight azimuth data output in step S14 provide the necessary angle parameters for calculating the occlusion indentation distance in step S23.

[0061] The discretized local polar coordinate occlusion perception model constructed in step S10 discretizes the initial visual recognition boundary line into multiple independent discrete boundary segments and establishes a local polar coordinate system centered on each discrete boundary segment. This transforms the complex occlusion analysis problem, which originally required considering both boundary orientation changes and viewpoint changes simultaneously in the global coordinate system, into a simplified angle determination problem performed independently in each local coordinate system. This problem decomposition strategy reduces the nonlinear occlusion calculation of continuous curved boundary segments to the linear occlusion calculation of discrete straight line segments. The occlusion condition determination of each discrete boundary segment depends only on its own tangential direction and the viewing azimuth angle of the current image. The calculations between each discrete boundary segment are independent of each other, making them suitable for parallel processing and significantly improving the computational efficiency in complex curved boundary scenarios. The design of the local polar coordinate system provides a unified coordinate reference for boundary segments with different orientations in their respective local coordinate systems. The tangential direction is used as the polar axis direction to normalize the boundary orientation. The determination of occlusion conditions can be applied to boundary segments with any orientation using a unified angle threshold, avoiding the complexity of setting different occlusion determination parameters for boundary segments with different orientations in the global coordinate system. The definitions of the inner and outer normal vectors of the field establish a correspondence between the spatial regions on both sides of the boundary segment and the shading effect. This allows the sign of the line-of-sight azimuth angle to directly reflect the camera's position relative to the inner and outer sides of the boundary segment, transforming the determination of shading conditions into a simple judgment of angle sign and amplitude. The shading condition filtering mechanism in step S10 accurately identifies the boundary segments requiring compensation under the current shooting angle by using a threshold judgment of the line-of-sight azimuth angle. This limits the scope of subsequent shading compensation calculations to a subset of boundary segments that genuinely require compensation. This precise filtering avoids meaningless compensation calculations for boundary segments without shading effects, reducing computational resource consumption. More importantly, it avoids applying erroneous compensation offsets to the outer side of the field segment observed from inside the field, preventing reverse deviations in boundary position due to incorrect compensation direction judgments. The boundary segments requiring compensation and their corresponding local polar coordinate system parameters output in step S10 provide a precisely defined geometric reference frame and processing unit division for the calculation of the equivalent shading canopy height and the solution of the shading compensation vector in step S20. The inner normal vector of the field defines the extension direction of the elongated boundary buffer zone constructed in step S21, ensuring that the extracted candidate canopy point clusters originate from the inner region of the field within the boundary segment. The outer normal vector of the field defines the direction of action of the single-view occlusion shrinkage compensation vector in step S23, ensuring that the compensation offset points outward to correct the systematic error of the visual boundary shrinking inward. The line-of-sight azimuth data provides a view-dependent correction factor for the calculation of the occlusion shrinkage distance in step S23, allowing the occlusion compensation amount to be adjusted according to the angle between the line of sight and the boundary segment normal, reflecting the difference in the intensity of the occlusion effect under different incident angles.The discretized local polar coordinate occlusion perception model established in step S10 decomposes the complex canopy occlusion geometry problem into a local problem that can be solved independently. While ensuring geometric accuracy, it reduces the algorithm complexity and lays a reliable geometric analysis foundation for subsequent multi-view occlusion compensation and dual-track fusion processing. This enables the entire farmland boundary extraction method to handle complex curved boundaries and multi-angle tilted images.

[0062] Step S20: Based on the boundary segment to be compensated, extract candidate canopy point cloud sets from the airborne lidar point cloud data, calculate the equivalent occlusion canopy height based on the candidate canopy point cloud sets, construct a multi-view occlusion shrinkage compensation vector set based on the equivalent occlusion canopy height, and determine the visual correction boundary points based on the visual occlusion shrinkage compensation vector set.

[0063] Further, step S20 includes:

[0064] Step S21: Extend the boundary segment to be compensated along the direction of the inner normal vector of the field to construct a long strip boundary buffer zone. Select the point cloud data falling into the long strip boundary buffer zone from the airborne lidar point cloud data to form a candidate canopy point cloud set.

[0065] Specifically, the elongated boundary buffer zone is a rectangular spatial region extending into the cultivated field from the boundary segment to be compensated as the geometric reference. It is used to filter canopy height information related to the current boundary segment to be compensated from the airborne lidar point cloud data. The process of constructing the elongated boundary buffer zone is as follows: taking the boundary segment to be compensated marked in step S14 as the starting edge, the length of the boundary segment to be compensated is the length dimension of the elongated boundary buffer zone; extending a preset buffer width W into the field along the direction of the inner normal vector of the field, the buffer width W is the width dimension of the elongated boundary buffer zone. The buffer width W is set based on the effective influence range of crop canopy occlusion on the visual boundary: if the buffer width W is too small, the number of extracted canopy point clouds will be insufficient to accurately reflect the actual height characteristics of the crop at the boundary; if the buffer width W is too large, it will introduce canopy points with large differences in height characteristics between the center area of ​​the field and the boundary, interfering with the accurate calculation of the equivalent occluded canopy height. The value of the buffer width W is related to the type of crop planted in the target cultivated area. For example, for tall crops such as corn, the buffer width W can be set to two to four times the average row spacing of the crops to ensure that the canopy information of the first two to four rows of crops at the boundary is included.

[0066] The coordinates of the four vertices of the elongated boundary buffer zone are determined through vector operations: using the starting point coordinates of the boundary segment to be compensated as the base point, the buffer width W is translated along the normal vector direction inside the field to obtain the two vertices of the elongated boundary buffer zone on the starting point side; using the ending point coordinates of the boundary segment to be compensated as the base point, the buffer width W is translated along the normal vector direction inside the field to obtain the two vertices of the elongated boundary buffer zone on the ending point side. The four vertices are connected sequentially to form a closed rectangular region boundary. The process of selecting candidate canopy point cloud sets from the airborne lidar point cloud data obtained in step S11 is as follows: traverse each point cloud data in the airborne lidar point cloud data, read the planar coordinates of the point cloud data, and determine whether the planar coordinates fall within the rectangular region of the elongated boundary buffer zone. The determination method uses the ray method or cross product method to determine the positional relationship between the point and the polygon. If the planar coordinates fall within the rectangular region, the point cloud data is included in the candidate canopy point cloud set; if the planar coordinates are outside the rectangular region, the point cloud data is discarded. A long, narrow boundary buffer zone was used for point cloud screening, limiting the calculation range of the subsequent equivalent shading canopy height to the spatial region directly related to the boundary segment to be compensated. This eliminated interference from point cloud data from the field center and other areas near the boundary segment. The design of the long, narrow boundary buffer zone extending along the normal vector direction of the inner side of the field ensures that the extracted candidate canopy point cloud sets come from the inner region of the boundary segment to be compensated. This aligns with the physical location of the canopy shading effect: the crop canopy is located inside the field of the boundary segment, and only the canopy height inside the field will obstruct the line of sight entering from outside the field. The width of the long, narrow boundary buffer zone ensures that the candidate canopy point cloud sets mainly include the height information of the first few rows of crops at the boundary. The shading effect of these front-row crops is much greater than that of the back-row crops, which is consistent with the physical characteristics of the shading effect. Step S21 uses the inner normal vector of the field defined in step S13 to determine the extension direction of the elongated boundary buffer zone, and uses the boundary segment to be compensated marked in step S14 to determine the positional reference of the elongated boundary buffer zone. Without the spatial filtering process in step S21, step S22 will not be able to obtain canopy point cloud data related to the current boundary segment to be compensated, and the calculation of the equivalent shading canopy height will lack a data source. The construction of the elongated boundary buffer zone provides step S22 with a spatially oriented candidate canopy point cloud set, and provides a unified spatial analysis area for step S31 to extract lidar reference boundary points and calculate local point cloud density.

[0067] Step S22: Remove low elevation points from the candidate canopy point cloud set, and retain the remaining point cloud as canopy points. Calculate the vertical distance from each retained canopy point to the boundary segment to be compensated, and use an inverse weighting function to calculate the equivalent shading canopy height of the boundary segment to be compensated based on the vertical distance.

[0068] Specifically, the equivalent shading canopy height is the optically equivalent height value characterizing the shading effect of the crop canopy on the line of sight at the boundary segment to be compensated, rather than the simple average or maximum height of the canopy. The process of removing low-elevation points from the candidate canopy point cloud set is as follows: all point cloud data in the candidate canopy point cloud set are sorted in ascending order by elevation value, and a preset ground point removal ratio R is set. ground Ground point removal ratio R ground Characterize the proportion of point cloud data belonging to ground points in the candidate canopy point cloud set, and sort the results by elevation (R). ground A certain proportion of the point cloud data is identified as low-elevation points and removed from the candidate canopy point cloud set; the remaining point cloud data is retained as canopy points. Ground point removal proportion R ground The setting is determined based on the probability that the lidar penetrates the crop canopy and reaches the ground. For example, in areas with dense canopies of tall crops, the ground point removal ratio R is... ground It can be set to 10% to 30%; for low-growing crop areas with sparse canopies, the ground point removal ratio R ground It can be set to 30% to 50%.

[0069] The process of calculating the perpendicular distance from each preserved canopy point to the boundary segment to be compensated is as follows: Read the planar coordinates of the preserved canopy point, and calculate the perpendicular distance from these coordinates to the line containing the boundary segment to be compensated. The perpendicular distance is calculated using the point-to-line distance formula. Let the planar coordinates of the preserved canopy point be P'. i If the starting point of the boundary segment to be compensated is A, and the ending point of the boundary segment to be compensated is B, then the canopy point P' i The vertical distance d to the boundary segment to be compensated i equal to vector with vector The magnitude of the cross product divided by the vector The modulus, the formula is: , where i is the index of the canopy point.

[0070] The process of calculating the equivalent shading canopy height based on vertical distance using an inverse weighting function is as follows: Assume there are N canopy points to be retained, and the elevation value of the i-th canopy point is h. i The vertical distance from the i-th canopy point to the boundary segment to be compensated is d. i The formula for calculating the equivalent shading canopy height Heff is:

[0071]

[0072] Among them, w i w is the weight value of the i-th canopy point. i Calculated using an inverse weighting function. ε is a preset distance smoothing constant. The setting of ε is used to avoid vertical distance d when the canopy point is exactly above the boundary segment to be compensated. i The problem of numerical overflow caused by weight values ​​approaching zero and becoming infinite can be addressed by setting the distance smoothing constant ε to one-tenth to one-fifth of the sampling step size L.

[0073] The inverse weighting function design ensures that canopy points closer to the boundary segment to be compensated contribute more to the equivalent shading canopy height, while canopy points farther from the boundary segment contribute less. This weighting mechanism aligns with the physical laws of canopy shading effects: when the line of sight is incident obliquely from outside the field, the front row of crop canopies near the boundary lies on the direct path between the line of sight and the ground boundary, resulting in the most significant shading effect; while the rear row of crop canopies further from the boundary may also obstruct the line of sight, their shading effect is partially weakened or completely blocked by the front row canopies. Traditional average height calculation methods assign the same weight to all canopy points, causing the equivalent shading canopy height to be pulled towards the overall average by the canopy height in the center of the field, failing to accurately reflect the actual shading height of the front row of crops at the boundary. The inverse weighting function ensures that the equivalent shading canopy height is primarily determined by canopy points near the boundary, more accurately characterizing the optically equivalent height that produces the shading effect. Step S22 uses the candidate canopy point cloud extracted in step S21 as input data, and the calculated equivalent shading canopy height serves as the core parameter for calculating the shading inward distance in step S23. Without the calculation of the equivalent shading canopy height in step S22, step S23 cannot determine the crop canopy shading height parameter, and the geometric calculation of the shading inward distance will lack the height dimension input. The equivalent shading canopy height transforms the three-dimensional spatial information of the LiDAR point cloud into a height parameter directly corresponding to the shading geometric model, thus realizing the data interface between LiDAR data and the visual shading compensation model.

[0074] Step S23: Combine the shooting pitch angle, line of sight azimuth angle and equivalent shading canopy height in the UAV shooting attitude parameters to calculate the shading inward distance, and construct a single-view shading inward compensation vector with a modulus of the shading inward distance along the normal vector direction of the outer side of the field.

[0075] Specifically, the occlusion inward distance is a geometric quantity representing the distance by which the visually recognized boundary shifts inward relative to the actual physical boundary of the field due to crop canopy occlusion at a specific shooting angle. The shooting pitch angle θ is read from the UAV shooting attitude parameters obtained in step S11. The shooting pitch angle θ represents the angle between the camera's optical axis and the horizontal plane. The value of the shooting pitch angle θ ranges from 0° to 90°. The smaller the value, the closer the camera is to shooting horizontally, and the more significant the occlusion effect; the larger the value, the closer the camera is to shooting vertically downward, and the weaker the occlusion effect. The line-of-sight azimuth angle φ is obtained from the calculation results in step S14. The line-of-sight azimuth angle φ represents the angle between the horizontal projection of the camera's line of sight and the tangential direction of the boundary segment to be compensated.

[0076] The calculation of the shading inward distance Δd is derived based on trigonometric geometric relationships: The crop canopy is simplified as a vertical shading body of height Heeff. When the camera observes the boundary segment obliquely at a shooting pitch angle θ, the top of the canopy will shade the ground boundary for a certain distance in the line-of-sight direction. Let Δd be the shading distance projected onto the ground by the top of the canopy in the line-of-sight direction. proj According to trigonometric relations, Δd proj =Heff×cot(θ), where cot(θ) is the cotangent of the shooting pitch angle θ. The above projection occlusion distance Δd proj This is the occlusion distance along the horizontal projection direction of the line of sight. It needs to be further projected to the normal direction perpendicular to the boundary segment to be compensated to obtain the true boundary shrinkage distance. There is an angle between the horizontal projection direction of the line of sight and the normal direction of the boundary segment to be compensated. This angle is equal to 90° minus the absolute value of the line of sight azimuth angle φ. Therefore, the formula for calculating the occlusion shrinkage distance Δd is: Δd = Heff × cot(θ) × |cos(φ)|, where |cos(φ)| is the absolute value of the cosine of the line of sight azimuth angle. When the azimuth angle φ of the line of sight is zero, the horizontal projection direction of the line of sight coincides with the normal direction of the boundary segment to be compensated, and |cos(φ)| equals one. The occlusion shrinkage distance reaches the extreme value under the current pitch angle. When the absolute value of the azimuth angle φ of the line of sight is close to 90°, the horizontal projection direction of the line of sight is nearly parallel to the tangential direction of the boundary segment to be compensated, and |cos(φ)| approaches zero. The occlusion shrinkage distance approaches zero. This is consistent with the physical fact that when the line of sight is observed along the direction of the boundary segment, the boundary segment will not have a displacement in the normal direction due to the occlusion of the canopy.

[0077] The process of constructing the single-view occlusion shrinkage compensation vector is as follows: The direction of the single-view occlusion shrinkage compensation vector is set to the field outer normal vector defined in step S13. The field outer normal vector points outward from the field, opposite to the direction in which the visual recognition boundary shrinks due to occlusion; the magnitude of the single-view occlusion shrinkage compensation vector is set to the occlusion shrinkage distance Δd. The single-view occlusion shrinkage compensation vector represents the distance and direction that needs to be offset from the midpoint of the boundary segment to be compensated on the initial visual recognition boundary line along the field outer direction, in order to correct the boundary shrinkage error caused by canopy occlusion in this frame of image. The occlusion shrinkage distance is calculated using trigonometric relationships, transforming the physical process of canopy occlusion effect into a mathematical model that can be accurately calculated, thus establishing a deterministic functional relationship between the occlusion compensation amount and the shooting angle parameters and canopy height parameters. The formula for occlusion shrinkage distance includes two angular parameters: the shooting pitch angle θ and the line-of-sight azimuth angle φ. This reflects the dual dependence of the occlusion effect on the shooting posture: the shooting pitch angle θ determines the tilt of the line of sight relative to the ground, and the line-of-sight azimuth angle φ determines the relative orientation of the line of sight to the normal of the boundary segment. The direction of the single-view occlusion shrinkage compensation vector is set to the normal vector outside the field, ensuring that the compensation offset direction is opposite to the boundary shrinkage direction, pushing the visual recognition boundary outward to approximate the real physical boundary. Step S23 comprehensively utilizes the field outside normal vector defined in step S13, the line-of-sight azimuth angle calculated in step S14, and the equivalent occlusion canopy height calculated in step S22 to construct an occlusion compensation vector for a single frame image. Without the occlusion shrinkage distance calculation in step S23, the occlusion error of the visual boundary cannot be quantified, and the multi-view consistency constraint in the subsequent step S24 will lose the input of the compensation vector. The single-view occlusion shrinkage compensation vector establishes independent occlusion compensation parameters for each frame of image, providing vector elements for step S24 to collect multi-view compensation vectors and perform consistent fusion.

[0078] Step S24: Collect multiple single-view occlusion shrinkage compensation vectors from different images for the same boundary segment to be compensated to form a multi-view occlusion shrinkage compensation vector set. Apply multi-view consistency constraints to the multi-view occlusion shrinkage compensation vector set to determine the visual correction boundary points.

[0079] Starting from the midpoint of the boundary segment to be compensated, each vector in the multi-view occlusion shrinkage compensation vector set is applied to obtain the candidate true boundary point location set. The least squares method is used to search for the spatial point with the smallest sum of squared Euclidean distances among all points in the candidate true boundary point location set, and this spatial point is determined as the visual correction boundary point.

[0080] Specifically, the multi-view occlusion shrinkage compensation vector set is a set of single-view occlusion shrinkage compensation vectors calculated from multiple images taken from different shooting angles for the same boundary segment to be compensated. Since the UAV flies along a preset flight path when acquiring oblique photogrammetric images, the same boundary segment to be compensated may appear in multiple images taken from different positions and attitudes. Each image corresponds to different camera spatial coordinates and shooting pitch angles, resulting in different line-of-sight azimuth angles and occlusion shrinkage distances. For the same boundary segment to be compensated, all images in the UAV oblique photogrammetric image set that can observe the boundary segment to be compensated are traversed, and the calculation process from steps S14 to S23 is repeated to obtain a set of single-view occlusion shrinkage compensation vectors. This set of vectors is then collected to form the multi-view occlusion shrinkage compensation vector set. Let the multi-view occlusion shrinkage compensation vector set contain M single-view occlusion shrinkage compensation vectors, denoted as V1, V2, ..., V M V M This is the Mth single-view occlusion shrinkage compensation vector.

[0081] The process of obtaining the candidate true boundary point position set by applying each vector in the multi-view occlusion shrinkage compensation vector set, starting from the midpoint of the boundary segment to be compensated, is as follows: Let the midpoint coordinate of the boundary segment to be compensated be O. The midpoint coordinate O is the center position of the boundary segment on the initial visual recognition boundary line. The j-th single-view occlusion shrinkage compensation vector V in the multi-view occlusion shrinkage compensation vector set... j Adding the coordinates of the midpoint O, we obtain the position P of the j-th candidate true boundary point. j =O+V j Perform the above operation on all M vectors in the multi-view occlusion shrinkage compensation vector set to obtain a set of candidate true boundary point positions containing M candidate true boundary point positions, denoted as P1, P2, ..., P... M P M Let P be the location of the Mth candidate true boundary point. The process of searching for the spatial point with the minimum sum of squared Euclidean distances from all points in the candidate true boundary point location set using the least squares method is as follows: Construct the objective function F(P), which is defined as the sum of squared Euclidean distances from the spatial point P to all points in the candidate true boundary point location set, i.e. , where ||PP j || represents the position P of the spatial point P to be determined and the j-th candidate true boundary point P. j The Euclidean distance between them. Find the spatial point P that minimizes the objective function F(P). Since the objective function F(P) is a quadratic convex function, the minimum point exists uniquely. Taking the partial derivative of the objective function F(P) with respect to P and setting the partial derivative to zero, we obtain the analytical solution of the minimum point as the geometric centroid of all points in the set of candidate true boundary points, i.e. The optimal spatial point P is... optimalThe visually corrected boundary point was determined. Essentially, the visually corrected boundary point is the best statistical estimate of the true physical location of the boundary segment to be compensated in three-dimensional space after the crop canopy shading effect is eliminated. It comprehensively utilizes redundant geometric information from multi-angle observations, and through spatial geometric intersection and statistical averaging, effectively corrects the systematic inward shrinkage deviation of the boundary system caused by shading in a single-view image. It has higher positioning accuracy and robustness than single-frame image recognition results.

[0082] The design of the multi-view consistency constraint is based on the following physical reasoning: the location of the real physical boundary is uniquely determined in space. If the occlusion compensation model is accurate and the shooting posture parameters of each frame are accurate, then the candidate real boundary point locations inferred from images from different perspectives should converge to the same spatial location. A single image may have posture parameter measurement errors or semantic segmentation recognition jitter, resulting in random error components in the single-view occlusion shrinkage compensation vector. These random errors are usually uncorrelated across different perspectives and conform to a normal distribution. The least squares method is used to obtain the geometric centroid of the candidate real boundary point location set. The statistical averaging effect of multi-view observations is used to reduce the impact of single-view random errors on boundary positioning accuracy, achieving a statistical improvement in boundary position accuracy. Step S24 gathers the single-view occlusion shrinkage compensation vectors calculated for different images in steps S14 to S23, and obtains the visually corrected boundary points through geometric intersection. Without the multi-view consistency constraint in step S24, boundary positioning will rely entirely on the occlusion compensation results of a single image, and the redundant information from multi-view observations cannot be used to reduce random errors.

[0083] Step S21 involves spatially oriented filtering of airborne lidar point cloud data by constructing a long, narrow boundary buffer zone. This restricts the point cloud data used for canopy height calculation to the inner region of the field directly related to the boundary segment to be compensated, eliminating interference from point cloud data in the field center and other areas near the boundary segment. This ensures that the calculated equivalent shading canopy height accurately reflects the actual height characteristics of the crop canopy at the boundary segment to be compensated. The design of the long, narrow boundary buffer zone extending along the positive normal vector of the inner side of the field aligns with the physical location of the canopy shading effect, ensuring that the extracted candidate canopy point cloud sets originate from the spatial region where shading occurs. An inverse weighting function is used to calculate the equivalent shading canopy height, ensuring that canopy points closer to the boundary segment to be compensated contribute more to the height calculation result, while canopy points farther away contribute less. This weighting mechanism conforms to the physical laws of the canopy shading effect: the shading effect of the front row crop canopy on the line of sight is much greater than that of the back row crop canopy. The inverse weighting function ensures that the equivalent occlusion canopy height is primarily determined by canopy points near the boundary, avoiding the drawback of traditional average height calculation methods that pull the equivalent height towards the overall average, thus improving the accuracy of the equivalent occlusion canopy height in representing the actual occlusion effect. The formula for calculating the occlusion retraction distance Δd is: Δd = Heff × cot(θ) × |cos(φ)|. This integrates the equivalent occlusion canopy height Heff, the shooting pitch angle θ, and the line-of-sight azimuth angle φ into the calculation of the occlusion compensation amount, allowing the compensation amount to adaptively adjust according to changes in the shooting angle. The occlusion retraction distance decreases as the shooting pitch angle θ increases and decreases as the line-of-sight azimuth angle φ deviates from the boundary normal direction, a variation consistent with the physical characteristics of the canopy occlusion effect. Multi-view consistency constraints are used to obtain visual correction boundary points. Redundant observation information from multiple images of the same boundary segment to be compensated in different viewpoints is utilized, and random errors from single-view observations are reduced through geometric centroid calculation. The multi-view fusion mechanism ensures that the boundary localization accuracy is not dependent on the observation quality of a single image, but rather determined by the comprehensive observation effect of multiple images, thus improving the robustness of the boundary localization results to single-frame image errors. The 3D spatial information of the lidar point cloud is transformed into equivalent occlusion canopy height parameters, establishing a connection between lidar data and the visual occlusion compensation model. Lidar point clouds provide crop canopy height information, while visual remote sensing provides planar position information and multi-view observation capabilities for the boundary; the two data sources complement each other in step S20. The visually corrected boundary points output in step S20 are boundary position estimates after occlusion geometric compensation and multi-view fusion, eliminating the systematic shrinkage error caused by canopy occlusion under tilted views. This provides high-quality visual boundary data for the dual-track boundary consistency verification in step S30, enabling subsequent lidar-visual fusion to be performed based on more accurate boundary estimation.

[0084] Step S30: Extract the lidar reference boundary points from the airborne lidar point cloud data and calculate the local point cloud density. Perform a dual-track boundary consistency check between the visual correction boundary points and the lidar reference boundary points. Reconstruct the physical boundary line of the cultivated land based on the check results.

[0085] Further, step S30 includes:

[0086] Step S31: Extract elevation change features from airborne lidar point cloud data. Within the elongated boundary buffer zone, extract lidar reference boundary points using elevation change features from airborne lidar point cloud data, and calculate the number of point clouds per unit area within the elongated boundary buffer zone to obtain the local point cloud density.

[0087] Specifically, the lidar reference boundary points are boundary position estimation points independently identified within the elongated boundary buffer zone using airborne lidar point cloud data. These points, along with the visually corrected boundary points obtained in step S24, form a parallel boundary localization result. The elevation abrupt change characteristic refers to the significant step change in the elevation values ​​of the lidar point cloud at the boundary of cultivated land. This is because cultivated land is used for planting crops, while non-cultivated land is occupied by roads, ditches, or bare land. The elevation values ​​of the point cloud on the cultivated land side are significantly higher than those on the non-cultivated land side. The process of extracting the lidar reference boundary points is as follows: Within the long strip-shaped boundary buffer zone constructed in step S21, the same buffer width W is extended along the opposite direction of the inner normal vector of the field to form an extended buffer area spanning both sides of the boundary segment to be compensated. This extended buffer area contains point cloud data from both the inner and outer sides of the field. The point cloud data within the extended buffer area are sorted according to the projection distance along the inner normal vector of the field, and the projection distance is divided into several equally spaced distance intervals. The width of each distance interval is set as a multiple of the average spacing of the point clouds. The average elevation of all point clouds in each distance interval is calculated to obtain the elevation distribution curve along the normal direction. The maximum value of the elevation mean change rate is searched on the elevation distribution curve. The elevation change rate is calculated by dividing the difference in the average elevation values ​​of adjacent distance intervals by the interval width. The distance interval corresponding to the maximum value of the change rate is the location where the elevation change occurs. The projection point of the elevation change location on the boundary segment to be compensated is determined as the lidar reference boundary point.

[0088] The physical basis for extracting lidar reference boundary points using elevation abrupt change features is as follows: In cultivated fields, the reflected points received by lidar mainly come from the crop canopy and a small number of pulses penetrating the canopy to reach the ground, resulting in an overall higher elevation distribution. In areas outside the fields, such as roads or bare land, the lidar reflected points come directly from the ground, resulting in a lower elevation distribution. The cultivated land boundary is precisely the dividing line between these two elevation distribution patterns; therefore, there is a correspondence between the elevation abrupt change location and the actual physical boundary location. Elevation abrupt change detection does not rely on visual recognition results from imagery; it is entirely based on the physical ranging principle of lidar, providing a boundary localization method independent of visual remote sensing.

[0089] Local point cloud density is a quantitative indicator characterizing the spatial sparseness of lidar point cloud distribution within a long, narrow boundary buffer zone, reflecting the ability of lidar pulses to penetrate the crop canopy and reach the ground. The formula for calculating local point cloud density ρ is: ρ = N ground / S, where N ground The number of point clouds identified as ground points within the long, narrow boundary buffer zone is determined using the same elevation sorting and elimination method as in step S22, selecting the top R values ​​from the elevation sorting results. ground The low-elevation points are determined as ground points; S is the planar area of ​​the elongated boundary buffer zone, which is equal to the length of the boundary segment to be compensated multiplied by the buffer width W. The physical meaning of the local point cloud density ρ is as follows: when the crop canopy is sparse or has gaps, the laser pulse can penetrate the canopy to reach the ground, resulting in a large number of ground points and a high local point cloud density ρ; when the crop canopy is dense, the laser pulse is blocked by the canopy and has difficulty reaching the ground, resulting in a sparse number of ground points and a low local point cloud density ρ. The local point cloud density ρ serves as the criterion for conflict arbitration in the subsequent step S33, indicating the reliability of the lidar reference boundary points.

[0090] Step S31 extracts the lidar reference boundary points within the elongated boundary buffer zone constructed in step S21, forming two independent boundary positioning tracks with the visually corrected boundary points obtained in step S24. The lidar track obtains elevation information based on the principle of physical ranging and is unaffected by visual occlusion, but is limited by point cloud density; the visual track identifies boundaries based on image texture, offering high resolution but subject to canopy occlusion errors. The data sources and error characteristics of the two tracks are independent, providing comparable dual boundary estimates for the consistency verification in step S32 and the conflict arbitration in step S33. Without the lidar reference boundary point extraction in step S31, the subsequent fusion process would rely solely on a single visually corrected boundary point, failing to utilize the lidar's penetration advantage to supplement boundary information in visually occluded areas, and also failing to identify and correct potential positioning errors through the dual-track verification mechanism.

[0091] Step S32, see Figure 3Calculate the spatial Euclidean distance deviation between the visually corrected boundary points and the lidar reference boundary points. When the spatial Euclidean distance deviation is greater than the preset consistency verification threshold, trigger the conflict arbitration mechanism based on local point cloud density and proceed to step S33 to obtain the dynamic fusion weights; otherwise, perform equal-weight fusion on the visually corrected boundary points and the lidar reference boundary points to obtain the final fused boundary points.

[0092] Specifically, the spatial Euclidean distance deviation is a geometric quantity characterizing the degree of spatial positional difference between the visual correction boundary point and the lidar reference boundary point. Let the planar coordinates of the visual correction boundary point obtained in step S24 be P. visual The planar coordinates of the lidar reference boundary point obtained in step S31 are P. lidar The spatial Euclidean distance deviation Δdist is calculated using the two-dimensional Euclidean distance formula: Δdist = ||P visual -P lidar ||, where ||·|| represents the magnitude operation of a vector, i.e., the arithmetic square root of the sum of the squares of the differences in the planar coordinates of two points. The spatial Euclidean distance deviation Δdist quantifies the degree of deviation between the boundary positions identified by the two techniques. The smaller the value, the closer the results of the two techniques are; the larger the value, the more significant the discrepancy between the two.

[0093] The consistency verification threshold Tcons is determined based on the combined upper limits of the positioning accuracy of the two boundary recognition methods: the residual error of visually corrected boundary points after occlusion compensation and multi-view fusion mainly comes from the measurement error of the shooting posture parameters and the edge positioning jitter of semantic segmentation. The amplitude of these error components is related to the image ground resolution and the accuracy of the posture measurement equipment; the positioning error of the lidar reference boundary points mainly comes from the discretization error of the point cloud sampling interval and the elevation change detection. The amplitude of these error components is related to the lidar point cloud density. The consistency verification threshold Tcons should be set as a multiple of the sum of the above two types of error amplitudes to ensure that the result deviation falls within the consistency threshold range when both technologies are working normally. For example, when the image ground resolution is at the centimeter level and the lidar point cloud density is tens of points per square meter, the consistency verification threshold Tcons can be set to 0.2 meters to 0.5 meters.

[0094] The execution process of dual-track boundary consistency verification is as follows: The calculated spatial Euclidean distance deviation Δdist is numerically compared with the consistency verification threshold Tcons; if Δdist is less than Tcons, it is determined that the positioning results of the visually corrected boundary point and the lidar reference boundary point are consistent, and the positioning of the two techniques in the current boundary segment is highly consistent. At this time, there is no need for conflict arbitration, and equal-weight fusion is directly performed on the two boundary points, that is, the arithmetic mean of the planar coordinates of the visually corrected boundary point and the lidar reference boundary point is taken as the final fused boundary point P.final =(P visual +P lidar If Δdist is greater than or equal to Tcons, it is determined that there is a significant discrepancy between the positioning results of the two techniques, and equal weight fusion cannot be simply adopted. It is necessary to trigger a conflict arbitration mechanism based on local point cloud density and proceed to step S33 to determine the differentiated dynamic fusion weight.

[0095] A dual-track boundary consistency verification mechanism is adopted to quickly output the fusion result when the results of the two techniques match, and trigger an arbitration process when the results differ, thus realizing conditional adaptive switching of the fusion strategy. Equal-weight fusion is reasonable when the consistency condition is met, because two independent techniques yielding similar results means that their respective systematic and random errors are small, and averaging can further reduce residual random errors. The arbitration mechanism is triggered when the spatial Euclidean distance deviation exceeds the consistency verification threshold. This criterion directly reflects the comparability of the two positioning results: forcibly averaging when the deviation is too large will introduce the positioning error of one party into the final result, reducing the fusion quality. Step S32 takes the visually corrected boundary points from step S24 and the lidar reference boundary points from step S31 as inputs, and determines whether to proceed to the conflict arbitration process in step S33 through consistency verification. Consistency verification acts as a quality control node in the fusion process, dividing the boundary segments into two categories according to the degree of matching between the two-track results: matching boundary segments can be directly fused and output, while differing boundary segments require further arbitration. This classification process avoids the quality degradation that might result from applying a uniform fusion strategy to all boundary segments, enabling the fusion process to handle boundary segments with different quality states differently. Without the consistency check in step S32, all boundary segments would enter the same fusion process, making it impossible to identify and distinguish between high-confidence boundary segments and those with discrepancies, significantly reducing the adaptability of the fusion strategy.

[0096] Step S33: Compare the local point cloud density with the effective penetration density threshold, and adaptively assign dynamic fusion weights to the visual correction boundary points and the lidar reference boundary points based on the comparison results to obtain the adaptively assigned dynamic fusion weights.

[0097] If the local point cloud density is greater than the effective penetration density threshold, the lidar reference boundary point is given a higher first weight and the visual correction boundary point is given a lower second weight. If the local point cloud density is less than or equal to the effective penetration density threshold, the visual correction boundary point is given a higher third weight and the lidar reference boundary point is given a lower fourth weight.

[0098] Specifically, the conflict arbitration mechanism is a decision-making process that, when there is a significant discrepancy between the positioning results of visually corrected boundary points and lidar reference boundary points, judges the credibility of each technical means based on local point cloud density and assigns differentiated fusion weights accordingly. The effective penetration density threshold Tpene is a critical point cloud density value that distinguishes whether lidar has successfully penetrated the crop canopy. When the local point cloud density is greater than the effective penetration density threshold, it indicates that the lidar pulse can effectively penetrate the canopy to reach the ground, and the obtained ground point cloud is sufficient to accurately characterize the boundary position, and the lidar reference boundary point has high credibility. When the local point cloud density is less than or equal to the effective penetration density threshold, it indicates that the crop canopy is too dense, and lidar has difficulty penetrating to obtain ground information, and the credibility of the lidar reference boundary point decreases. In this case, the visually corrected boundary points after occlusion geometric compensation are more reliable. The effective penetration density threshold Tpene is determined as follows: At known boundary locations under different canopy density conditions, lidar point cloud data is collected and local point cloud density is calculated. Simultaneously, the positioning error between the lidar reference boundary point and the known true boundary location is measured. A correspondence between local point cloud density and positioning error is established, and the trend of increasing positioning error as point cloud density decreases is analyzed. The point cloud density value corresponding to the point where the positioning error begins to increase significantly is determined as the effective penetration density threshold Tpene. For example, for areas with tall crops such as corn, the effective penetration density threshold Tpene can be set to 5-10 points per square meter. This range ensures that when the ground point cloud density is lower than this value, the lidar's ability to identify the boundary significantly decreases.

[0099] The execution process of conflict arbitration and dynamic weight allocation is as follows: The local point cloud density ρ calculated in step S31 is numerically compared with the effective penetration density threshold Tpene; if ρ is greater than Tpene, the lidar penetration capability is sufficient, and a higher first weight w is assigned to the lidar reference boundary point. lidar,high The second weight w with lower visual correction boundary points visual,low The weights satisfy the normalization constraint w lidar,high +w visual,low =1; If ρ is less than or equal to Tpene, the LiDAR penetration capability is insufficient, and a higher third weight w is assigned to the visual correction boundary point. visual,high The fourth weight w with lower reference boundary points of the lidar lidar,low The weights also satisfy the normalization constraint w. visual,high +w lidar,low =1. The specific value of the weight is determined based on the degree to which the local point cloud density deviates from the effective penetration density threshold. The greater the deviation, the higher the weight of the reliable data source. For example, the first weight w lidar,high and the third weight w visual,high It can be set to 0.7 to 0.9, the second weight w visual,lowand the fourth weight w lidar,low It can be set to 0.1 to 0.3.

[0100] The physical basis for conflict arbitration based on local point cloud density is that local point cloud density directly reflects the ability of the lidar to acquire ground information at the current boundary segment, and the quality of ground information acquisition determines the accuracy of elevation change detection, and thus the reliability of the lidar reference boundary point. When dense canopy leads to sparse point cloud, the number of sample points for elevation change detection decreases, statistical stability declines, and random error of the detection result increases. In contrast, although visually corrected boundary points suffer from systematic shrinkage error caused by canopy occlusion, this error has been corrected through occlusion geometric compensation in step S20. In dense canopy areas, the residual error of visually corrected boundary points after compensation may be smaller than the random error of the lidar reference boundary point. The dynamic weight allocation mechanism enables the fusion process to automatically select the data source with higher reliability at the current boundary segment as the dominant one, achieving complementary utilization of the advantages of the two technologies. Step S33 forms a data-driven decision-making collaboration relationship with steps S31 and S32. Step S33 uses the local point cloud density calculated in step S31 as the arbitration criterion to allocate weights to the divergence boundary segment triggered in step S32. Local point cloud density, as a quantitative indicator of lidar penetration capability, transforms abstract sensor performance differences into numerical parameters that can be directly used for weight calculation, giving arbitration decisions a clear physical basis rather than subjective settings. Without the conflict arbitration mechanism in step S33, the divergence boundary segments cannot obtain reasonable fusion weights, potentially leading to errors from low-reliability data sources being carried into the final fusion result, weakening the fusion quality. The dynamic weight allocation mechanism allows the fusion process to be personalized for the specific data quality status of each boundary segment, avoiding the inconsistency of fixed-weight strategies under different crop density conditions.

[0101] Step S34: Use dynamic fusion weights to perform weighted calculations on the visually corrected boundary points and the lidar reference boundary points to obtain the final fused boundary points. Perform spline interpolation on all the final fused boundary points corresponding to the discrete boundary segment sequence to generate continuous farmland physical boundary lines.

[0102] Specifically, the final fused boundary point is the weighted average position of the visually corrected boundary point and the lidar reference boundary point, representing the positioning result of the physical boundary of the cultivated land at the current boundary segment to be compensated. The weighted calculation process is as follows: Let the planar coordinates of the visually corrected boundary point be P. visual The planar coordinates of the lidar reference boundary point are P. lidar The fusion weight for visually corrected boundary points is w. visual The fusion weight of the lidar reference boundary points is w lidarThe fusion weights are derived from the equal weight allocation in step S32 or the dynamic weight allocation in step S33, satisfying the normalization constraint w. visual +w lidar =1; Planar coordinates P of the final fusion boundary point final The calculation formula is: P final =w visual ×P visual +w lidar ×P lidar The weighted average calculation biases the position of the final fusion boundary point towards the data source with higher weight. A higher weight means that the data source has a higher degree of credibility at the current boundary segment. Therefore, the weighted average result can make full use of the positioning information of the credible data source. The process of orderly aggregating all the final fusion boundary points corresponding to the discrete boundary segment sequence is as follows: The discrete boundary segment sequence constructed in step S12 contains multiple discrete boundary segments arranged sequentially according to the direction of the initial visual recognition boundary line. Each discrete boundary segment generates a final fusion boundary point after being processed in steps S20 to S33. According to the arrangement order of the discrete boundary segments in the discrete boundary segment sequence, the final fusion boundary points corresponding to each discrete boundary segment are arranged sequentially to form the final fusion boundary point sequence. The final fusion boundary point sequence retains the topological orientation of the original boundary, and the connection relationship between two adjacent final fusion boundary points in the sequence is consistent with the adjacency relationship of the original discrete boundary segments.

[0103] Spline interpolation is a mathematical method for reconstructing a discrete sequence of final fusion boundary points into a continuous, smooth curve. The implementation process of the cubic spline interpolation algorithm is as follows: Each point in the final fusion boundary point sequence is used as an interpolation node. Let the number of nodes be K, and the planar coordinates of each node be P... final,1 P final,2 , ..., P final,KA cubic polynomial curve segment is constructed between two adjacent nodes, and the entire curve is composed of K-1 cubic polynomial curve segments connected end to end. At each node, the function values, first derivatives, and second derivatives of the two adjacent curve segments are required to be continuous. These continuity conditions ensure a smooth transition of the curve at the nodes, without sharp corners or inflections. Natural boundary conditions are used, meaning the second derivative of the curve is zero at both endpoints. By solving the linear equations formed by the continuity and boundary conditions, the coefficients of each cubic polynomial segment are determined, thus obtaining the complete cubic spline interpolation curve. The implementation process of the Bézier curve fitting algorithm is similar to that of cubic spline interpolation, except that Bernstein basis functions are used as the mathematical expression of the curve, and the position of the control points is optimized by minimizing the distance between the curve and the final fused boundary point. The physical boundary line of cultivated land is a continuous and smooth curve obtained after spline interpolation, representing the true physical boundary of the target cultivated land area. The physical boundary line of cultivated land has the following geometric characteristics: the curve passes through or approaches all the final fused boundary points, ensuring consistency with the discrete positioning results; the curve has a continuous tangential direction and curvature at each node, avoiding artificial angles introduced by discretization processing, and conforming to the natural shape of the real field boundary; the curve can be interpolated at any position, providing continuous boundary position information, and meeting the accuracy requirements of subsequent area calculation and cadastral survey.

[0104] Step S34, together with steps S12, S32, and S33, forms a collaborative reconstruction relationship from discrete to continuous. Step S34 applies the fusion weights determined in steps S32 and S33 to the weighted calculation of boundary points, reintegrating the multiple boundary segments generated by the discretization in step S12 into a continuous curve. Spline interpolation eliminates the piecewise effect introduced by the discretization process on the boundary curve, restoring the continuity and smoothness of the boundary. Without the spline interpolation in step S34, the output would only be a discrete set of boundary points, unsuitable for direct use in area measurement and cadastral registration, and the visualization of the boundary would present a broken line shape rather than a natural curve shape. The physical boundary line of cultivated land, as the output of the entire method, integrates the high-resolution advantages of visual remote sensing, the systematic error correction capability of occlusion geometric compensation, the random error suppression capability of multi-view fusion, the penetration ranging advantage of lidar, and the quality control capability of density adaptive fusion, effectively eliminating the systematic inward shrinkage error of the boundary caused by canopy occlusion under tilted views.

[0105] Step S30 constructs a dual-track parallel boundary recognition architecture combining lidar and visual remote sensing. LiDAR's penetration characteristics are used to independently extract boundary location information, forming a mutually verified dual positioning result with the visually corrected boundary points obtained in step S20 after occlusion compensation. The lidar track acquires elevation information based on physical ranging principles, unaffected by visual occlusion from crop canopies, and can accurately identify boundary locations in sparse canopy areas. The visual track identifies boundaries based on image texture and spectral information, offering high resolution and continuous coverage. After occlusion geometric compensation, it can provide reliable boundary estimates in dense canopy areas. The data sources and error mechanisms of the two tracks are independent, forming a complementary redundant positioning capability. A dual-track boundary consistency verification mechanism is employed, rapidly outputting the fusion result when the visually corrected boundary points match the lidar reference boundary point positioning results, and triggering a conflict arbitration mechanism based on local point cloud density when the positioning results differ. Consistency verification categorizes boundary segments according to the degree of agreement between the two-track results, enabling the fusion strategy to differentiate between boundary segments with different quality levels: high-consistency boundary segments are further reduced by equal-weight fusion to reduce residual random error, while low-consistency boundary segments select a more reliable data source as the fusion leader through density-adaptive weight allocation. This condition-adaptive fusion strategy avoids the inconsistency problem of fixed-weight strategies under different crop density conditions, improving the adaptability of the fusion results to complex field conditions.

[0106] Step S30 uses local point cloud density as the criterion for conflict arbitration, providing a clear physical basis for weight allocation decisions. Local point cloud density directly reflects the lidar's ability to penetrate the crop canopy. High density means the lidar has successfully acquired ground information, resulting in accurate and reliable elevation change detection; low density means the lidar is blocked by the canopy, leading to missing ground information and decreased statistical stability of elevation change detection. The dynamic weight allocation mechanism automatically adjusts the contribution ratio of the two data sources based on local point cloud density. In sparse canopy areas, it tends to rely on the physical ranging accuracy of the lidar, while in dense canopy areas, it tends to rely on the visually corrected results after occlusion compensation. This adaptive weight allocation strategy achieves complementary utilization of the lidar's penetration advantage and the resolution advantage of visual remote sensing, ensuring high positioning accuracy of the fusion result under different crop density conditions. Cubic spline interpolation or Bézier curve fitting algorithms are used to reconstruct the discrete final fusion boundary points into a continuous and smooth physical boundary line of cultivated land, eliminating the piecewise effect introduced by discretization on the boundary curve and restoring the continuity and smoothness of the boundary. The spline interpolation curve has continuous tangential direction and curvature at each node, conforming to the natural shape of the actual field boundary and avoiding geometric distortion that may be introduced by the broken line shape in area calculation and cadastral registration. The physical boundary line of cultivated land integrates the high-resolution advantages of visual remote sensing, the systematic error correction capability of occlusion geometric compensation, the random error suppression capability of multi-view fusion, the penetration and ranging advantages of lidar, and the quality control capability of density adaptive fusion. This makes the final output boundary positioning accuracy significantly better than the positioning capability of a single technical means, and eliminates the systematic inward deviation caused by canopy occlusion in tall crop areas, meeting the boundary accuracy requirements of cultivated land rights confirmation measurement. The physical boundary line of cultivated land output in step S30 does not depend on the performance of a single technical means, but achieves mutual compensation of the limitations of technical means through dual-track verification and density adaptive fusion mechanism. Visual remote sensing is easily affected by occlusion in dense canopy areas, and lidar's penetration capability decreases in dense canopy areas. The limitations of the two technical means exhibit complementary characteristics in the crop density dimension. The fusion strategy in step S30 uses local point cloud density as an indicator of the reliability of the two technical means, and automatically selects the more reliable data source as the fusion leader under different crop density conditions, so that the final boundary positioning accuracy is not limited by the limitation of a single technical means, and can maintain stable and reliable positioning quality throughout the entire domain.

[0107] Example 2

[0108] This embodiment, based on Embodiment 1, provides a farmland boundary extraction system that combines UAV remote sensing and lidar, such as... Figure 4 As shown, it includes:

[0109] Occlusion perception filtering module: used to acquire UAV oblique photography image set, airborne LiDAR point cloud data and UAV shooting attitude parameters of the target area, extract initial visual recognition boundary line from UAV oblique photography image set, construct discretized local polar coordinate occlusion perception model based on initial visual recognition boundary line and UAV shooting attitude parameters, and use discretized local polar coordinate occlusion perception model to filter boundary segments to be compensated.

[0110] Visual correction solution module: Based on the boundary segment to be compensated, candidate canopy point cloud sets are extracted from the airborne lidar point cloud data. The equivalent occlusion canopy height is calculated based on the candidate canopy point cloud sets. A set of multi-view occlusion shrinkage compensation vectors is constructed based on the equivalent occlusion canopy height. The visual correction boundary points are determined based on the set of visual occlusion shrinkage compensation vectors.

[0111] Dual-track verification and reconstruction module: used to extract lidar reference boundary points from airborne lidar point cloud data, perform dual-track boundary consistency verification between visual correction boundary points and lidar reference boundary points, and reconstruct the physical boundary line of cultivated land based on the verification results.

[0112] Furthermore, in the occlusion perception and filtering module, the method of constructing a discretized local polar coordinate occlusion perception model and using the discretized local polar coordinate occlusion perception model to filter the boundary segments to be compensated includes:

[0113] The initial visual recognition boundary line is subjected to equidistant divergence processing to generate a sequence of discrete boundary segments. For each discrete boundary segment in the sequence, a local polar coordinate system is constructed and occlusion condition determination is performed. Discrete boundary segments that meet the occlusion condition are marked as boundary segments to be compensated.

[0114] Furthermore, in the visual correction solution module, the method for extracting the candidate canopy point cloud includes:

[0115] Extend the boundary segment to be compensated along the direction of the inner normal vector of the field to construct a long strip boundary buffer zone. Select the point cloud data falling into the long strip boundary buffer zone from the airborne lidar point cloud data to form a candidate canopy point cloud set.

[0116] The method for calculating the equivalent shading canopy height includes: identifying low-elevation points in the candidate canopy point cloud set; removing low-elevation points from the candidate canopy point cloud set; retaining the remaining point cloud as canopy points; calculating the vertical distance from each retained canopy point to the boundary segment to be compensated; and using an inverse weighting function to calculate the equivalent shading canopy height of the boundary segment to be compensated based on the vertical distance.

[0117] Furthermore, in the dual-track verification and reconstruction module, the dual-track boundary consistency verification includes:

[0118] The spatial Euclidean distance deviation between the visually corrected boundary points and the lidar reference boundary points is calculated. When the spatial Euclidean distance deviation is greater than the preset consistency verification threshold, a conflict arbitration mechanism based on local point cloud density is triggered to obtain dynamic fusion weights. The visually corrected boundary points and lidar reference boundary points are weighted using the dynamic fusion weights to obtain the final fused boundary points. Otherwise, equal-weight fusion is performed on the visually corrected boundary points and lidar reference boundary points to obtain the final fused boundary points.

[0119] The methods and systems of this application may be implemented in many ways. For example, they may be implemented by software, hardware, firmware, or any combination of software, hardware, and firmware. The above-described order of steps for the method is for illustrative purposes only, and the steps of the method of this application are not limited to the order specifically described above, unless otherwise specifically stated.

[0120] In addition, the parts of the technical solutions provided in the embodiments of this application that are consistent with the implementation principles of the corresponding technical solutions in the prior art have not been described in detail, so as to avoid excessive elaboration.

[0121] The specific embodiments described above further illustrate the purpose, technical solution, and beneficial effects of the present invention. It should be understood that the above descriptions are merely specific embodiments of the present invention and are not intended to limit the invention. Any modifications, equivalent substitutions, or improvements made within the spirit and principles of the present invention should be included within the scope of protection of the present invention.

Claims

1. A method for extracting farmland boundaries using a combination of UAV remote sensing and lidar, characterized in that, The method includes: Acquire a set of oblique UAV photographic images of the target area, point cloud data of airborne LiDAR, and UAV shooting attitude parameters. The UAV shooting attitude parameters include at least the camera spatial position coordinates and shooting pitch angle. Extract initial visual recognition boundary lines from the UAV oblique photographic image set, perform equidistant divergence processing on the initial visual recognition boundary lines to generate a sequence of discrete boundary segments, construct a local polar coordinate system for each discrete boundary segment in the discrete boundary segment sequence and perform occlusion condition determination, and mark the discrete boundary segments that meet the occlusion conditions as boundary segments to be compensated. Candidate canopy point sets are extracted from airborne lidar point cloud data based on the boundary segment to be compensated. The equivalent occlusion canopy height is calculated based on the candidate canopy point sets. A set of multi-view occlusion shrinkage compensation vectors is constructed based on the equivalent occlusion canopy height. The visual correction boundary points are determined based on the set of visual occlusion shrinkage compensation vectors. The method for constructing a multi-view occlusion shrinkage compensation vector set includes: calculating the line-of-sight azimuth angle and the equivalent occlusion canopy height; combining the shooting pitch angle, line-of-sight azimuth angle, and equivalent occlusion canopy height to calculate the occlusion shrinkage distance; defining a field-outer normal vector pointing to the outside of the field in a local polar coordinate system; constructing a single-view occlusion shrinkage compensation vector with a modulus equal to the occlusion shrinkage distance along the direction of the field-outer normal vector; and collecting multiple single-view occlusion shrinkage compensation vectors from different images for the same boundary segment to be compensated to form a multi-view occlusion shrinkage compensation vector set. Extract lidar reference boundary points from airborne lidar point cloud data, perform dual-track boundary consistency verification between visual correction boundary points and lidar reference boundary points, and reconstruct the physical boundary line of cultivated land based on the verification results.

2. The method for extracting farmland boundaries using a combination of UAV remote sensing and lidar as described in claim 1, characterized in that, The method for extracting the initial visual recognition boundary line includes: Orthorectification and semantic segmentation are performed on the UAV oblique photography image set to extract the boundary pixel set between cultivated land and non-cultivated land, and curve fitting is performed on the boundary pixel set to generate the initial visual recognition boundary line.

3. The method for extracting farmland boundaries using a combination of UAV remote sensing and lidar as described in claim 2, characterized in that, The method for generating the discrete boundary segment sequence includes: A sequence of discrete boundary sampling points is obtained by performing equidistant divergent sampling along the initial visual recognition boundary line. Adjacent discrete boundary sampling points in the sequence are then connected sequentially to construct a sequence of discrete boundary segments.

4. The method for extracting farmland boundaries using a combination of UAV remote sensing and lidar as described in claim 3, characterized in that, The method for constructing a local polar coordinate system includes: Select the target boundary segment from the discrete boundary segment sequence, calculate the tangential direction of the target boundary segment, and establish a local polar coordinate system with the midpoint of the target boundary segment as the origin and the tangential direction as the polar axis.

5. The method for extracting farmland boundaries using a combination of UAV remote sensing and lidar as described in claim 4, characterized in that, The method for determining the occlusion condition includes: Define a positive normal vector pointing inside the field in the local polar coordinate system; The camera line-of-sight vector is constructed based on the camera spatial position coordinates and the midpoint coordinates of the target boundary segment in the drone's shooting attitude parameters. The camera line-of-sight vector is then projected onto the horizontal plane of the local polar coordinate system to calculate the line-of-sight azimuth angle. Determine whether the target boundary segment meets the occlusion condition based on the line of sight azimuth, and mark the target boundary segment that meets the occlusion condition as the boundary segment to be compensated.

6. The method for extracting farmland boundaries using a combination of UAV remote sensing and lidar as described in claim 5, characterized in that, The criterion for determining whether the occlusion condition is met is that the absolute value of the line-of-sight azimuth angle is less than 90°.

7. The method for extracting farmland boundaries using a combination of UAV remote sensing and lidar as described in claim 6, characterized in that, The method for extracting the candidate canopy point cloud includes: Extend the boundary segment to be compensated along the direction of the inner normal vector of the field to construct a long strip boundary buffer zone. Select the point cloud data falling into the long strip boundary buffer zone from the airborne lidar point cloud data to form a candidate canopy point cloud set.

8. The method for extracting farmland boundaries using a combination of UAV remote sensing and lidar as described in claim 7, characterized in that, The method for calculating the equivalent obstructed canopy height includes: Identify low-elevation points in the candidate canopy point cloud set, remove low-elevation points from the candidate canopy point cloud set, and retain the remaining point cloud as canopy points. Calculate the vertical distance from each retained canopy point to the boundary segment to be compensated, and use an inverse weighting function to calculate the equivalent shading canopy height of the boundary segment to be compensated based on the vertical distance.

9. The method for extracting farmland boundaries using a combination of UAV remote sensing and lidar as described in claim 8, characterized in that, The dual-track boundary consistency check includes: The spatial Euclidean distance deviation between the visually corrected boundary points and the lidar reference boundary points is calculated. When the spatial Euclidean distance deviation is greater than the preset consistency verification threshold, a conflict arbitration mechanism based on local point cloud density is triggered to obtain dynamic fusion weights. The visually corrected boundary points and lidar reference boundary points are weighted using the dynamic fusion weights to obtain the final fused boundary points. Otherwise, equal-weight fusion is performed on the visually corrected boundary points and lidar reference boundary points to obtain the final fused boundary points.

Images