CN110796741B - Airborne laser sounding point cloud filtering method based on bidirectional cloth simulation - Google Patents

Abstract

The invention discloses an airborne laser sounding point cloud filtering method based on bidirectional cloth simulation, which belongs to the technical field of ocean mapping, and comprises the steps of firstly acquiring airborne laser sounding data and calculating the coordinates of a submarine laser point; then eliminating the negative abnormal data of the seabed by constructing a transmission type iteration trend surface and expressing the continuity of the seabed topography; and finally, constructing a filtering surface according to the bidirectional cloth simulation correction model to finish airborne laser sounding point cloud filtering. The method realizes high-precision and strong-robustness filtering of the airborne laser sounding point cloud, and effectively solves the problems that the existing method can not identify negative abnormal data of the airborne laser sounding and is easy to generate excessive filtering.

Description

Airborne laser sounding point cloud filtering method based on bidirectional cloth simulation

Technical Field

The invention belongs to the technical field of ocean mapping, and particularly relates to an airborne laser sounding point cloud filtering method based on bidirectional cloth simulation.

Background

The airborne laser sounding system has the characteristics of high measurement accuracy, high measurement point density, high working efficiency, strong maneuverability, measurement continuity and the like, is particularly suitable for rapid detection of complex terrains such as shallow water areas, island reef nearby areas and the like, and can realize seamless splicing of coastline water and underwater terrains. Has important significance for meeting urgent demands of ocean, mapping, water conservancy, traffic, outcrossing, navy and other aspects. The airborne laser sounding system adopts 532nm blue-green laser with stronger water body penetrating ability as a submarine detection means, when the blue-green laser reaches the sea surface, one part of laser beams are reflected along the original path, and the other part of laser beams penetrate the sea surface to enter the water body. However, when a large amount of floating algae, fish shoals or submarine vegetation exist in the shallow sea water body, laser beams entering the water body cannot reach the sea bottom and are reflected in advance, so that a large amount of irregular noise points are generated, and the authenticity expression of the submarine topography is affected. The submarine topography is a direct basis for understanding the external shape of the earth, the movement of submarine structures and the submarine evolution, and is also important basic data in the aspects of ocean economic development, ocean scientific research, ocean military application and the like. Therefore, in order to extract the real submarine topography, it is necessary to study an effective submarine laser point cloud filtering method.

The current mainstream filtering methods are mainly divided into multi-beam water depth data filtering and land airborne laser scanning data filtering, but the applicability of the methods to airborne laser sounding data is not verified, and the methods related to airborne laser sounding point cloud data filtering are relatively few, so that an airborne laser sounding point cloud filtering method based on bidirectional cloth simulation is necessary to be provided, and high-precision and strong-robustness filtering of airborne laser sounding point clouds is realized.

Disclosure of Invention

Aiming at the technical problems in the prior art, the invention provides the airborne laser sounding point cloud filtering method based on bidirectional distribution simulation, which is reasonable in design, solves the problems that the airborne laser sounding negative abnormal data cannot be identified and excessive filtering is easy to generate in the prior art, fully reserves the characteristics of seafloor micro-topography, ensures the integrity of seafloor topography data, and has better applicability and robustness for the airborne laser sounding data.

In order to achieve the above purpose, the present invention adopts the following technical scheme:

an airborne laser sounding point cloud filtering method based on bidirectional cloth simulation comprises the following steps:

step 1: acquiring airborne laser sounding point cloud data, and calculating the coordinates of the submarine laser points;

step 2: fitting by adopting a quadric LM (Levenberg-Marquardt) algorithm, constructing a transfer iteration trend surface, removing submarine negative abnormal data, and expressing submarine continuity;

step 3: and constructing a bidirectional distribution simulation correction model, and performing bidirectional filtering of the submarine laser point cloud according to the correction model to obtain a final filtering surface.

Preferably, in step 2, the construction method of the transitive iterative trend surface is as follows:

step 2.1: quadric LM algorithm fitting: assuming that a point P is taken as a center, r is taken as a radius, searching all submarine laser points in a given radius area, and a quadric surface model formed by fitting all laser points in the area is shown as a formula (1):

wherein, (x) _j ,y _j ,z _j ) For the point set N _j Geographic coordinates of the laser point within (j=1, 2, …, k) in the local coordinate system; a. b, c, d, e and f are quadric fitting parameters;

step 2.2: at point P ₁ Taking r as a radius as a center, searching all laser points in the area and performing quadric surface LM fitting; calculating the water depth average value of laser points in curved surfaces

Standard deviation sigma of water depth _H Find and +.>

The difference is greater than 2 sigma _H The laser points are marked as negative abnormal points and stored into a suspicious point set, and the water depth values of all suspicious points are calculated into a fitted quadric surface to ensure a discrete point set N ₁ The method does not contain coarse differences and ensures the integrity of sounding data;

step 2.3: find and point P ₁ The nearest point, designated as point P ₂ Find (N) ₂ -N ₁ ∩N ₂ ) Point concentration

The difference is greater than 2 sigma _H The laser points are marked as negative abnormal points and stored into a suspicious point set, and the water depth values of all suspicious points are calculated into a fitted quadric surface to ensure a discrete point set N ₂ The method does not contain coarse differences;

step 2.4: find and point P ₂ The nearest point, excluding the completed point, is denoted as point P ₃ And (2.3) repeating the step until all submarine laser points are traversed, storing no laser points into a suspicious point set, and obtaining the final submarine topography trend surface when the fitted topography trend surface reaches a stable state.

Preferably, in step 3, the specific construction method of the bidirectional cloth simulation correction model is as follows:

step 3.1: setting filtering parameters: the grid resolution GR, namely the distance between two adjacent nodes; the stiffness parameter IMRI controls the softness of the simulated cloth; the maximum iteration number MI is used for judging whether the simulation process is terminated; judging classification thresholds of the sea bottom points and the non-sea bottom points according to the self-adaptive distance threshold H;

step 3.2: respectively performing forward filtering and backward filtering on the submarine laser points, namely respectively inverting the submarine laser points and not inverting the submarine laser points, and placing the simulated cloth right above the laser points and keeping a certain distance from the laser points; projecting the distribution node and the laser point to the same horizontal plane, performing grid-networking, determining the corresponding relation between the distribution node and the laser point, and recording the elevation values of the distribution node and the laser point;

step 3.3: the simulated cloth is enabled to do free falling motion, and the position and the stress effect of the cloth node meet the formula (2):

wherein m represents node Z _l Is the mass of (3); x represents Z _l Position at time t; g represents gravitational acceleration; c _e Representing the elastic deformation coefficient of the node spring; e represents Z _l A set of neighbors;

representing node Z _l And node Z _n Distance at time t>

Representing node Z _l And node Z _n Original distance between>

Representing Z _l Pointing to Z _n Is a unit vector of (a);

step 3.4: in the process of simulating the descending of cloth, performing collision detection on cloth nodes;

assuming that the collision state is defined as C, the collision constraint is determined by the position and velocity constraints of the cloth node, and there are two time steps:

/>

where dt is the time step; XC (t), XC (t+dt) and XC (t+2dt) are the positions of the collision at times t, t+dt and t+2dt, respectively; XC ' (t), XC ' (t+dt) and XC ' (t+2dt) are the movement velocities of the collisions at times t, t+dt and t+2dt, respectively; XC "(t) and XC" (t+dt) are the motion accelerations of the collision at times t and t+dt, respectively; Δxc "(t) and Δxc" (t+dt) are corrected accelerations of the collision at times t and t+dt, respectively, and thus equation (4) can be obtained:

wherein Δxc "(t) can be converted into additional impact forces acting on the distribution node;

step 3.5: judging the position relation between the distribution node and the submarine laser points on the basis of collision detection;

if the elevation value of the cloth node is H _cp Less than the elevation value H of the submarine laser spot _sp Moving the node to the corresponding laser point location and marking the node as an immovable point; otherwise calculate the node inside the nodeDisplacement to be adjusted under the action of traction force;

step 3.6: repeating the steps 3.3 to 3.5 until the iteration times of the nodes are larger than the set maximum iteration times MI, and ending the cloth simulation process;

step 3.7: calculating the height difference delta H between the distribution node and the corresponding submarine laser point, and judging the size relation between delta H and the self-adaptive distance threshold H; if delta H is less than H, the node is reserved as a submarine point, otherwise, the node is judged to be a non-submarine point and deleted;

step 3.8: establishing a forward filtering surface F _F (P _T IMRI, H) and inverse filtering plane F _I (P _T ,IMRI,H)；

Defining a forward filtering surface as the simulated cloth formed after forward filtering, and defining a reverse filtering surface as the simulated cloth formed after reverse filtering; let the forward filtering surface be F _F (P _T IMRI, H), on which the distribution node is denoted P _F (x _F ,y _F ,z _F ) The method comprises the steps of carrying out a first treatment on the surface of the The reverse filtering surface is F _I (P _T IMRI, H), on which the distribution node is denoted P _I (x _I ,y _I ,z _I ) The method comprises the steps of carrying out a first treatment on the surface of the The two filter surfaces are represented as shown in equation (5):

wherein F is _F (P _T IMRI, H) and F _I (P _T IMRI, H) is the forward and reverse filtering face; p (P) _T (x _T ,y _T ,z _T )(z _T ＝f(x _T ,y _T ) A submarine laser spot); the IMRI is a rigidity parameter for controlling the softness of the cloth; h is a distance threshold value for judging that the distribution node is a submarine point and a non-submarine point; (x) _F ,y _F ,z _F ) And (x) _I ,y _I ,z _I ) For distributing the node P _F And P _I Is satisfied with z _F ＝f _F (x _F ,y _F ) And z _I ＝f _I (x _I ,y _I ) The method comprises the steps of carrying out a first treatment on the surface of the T represents the sea floor topography;

step 3.9:establishing a submarine filtering surface F (P) _T ,IMRI,H)；

Since some nodes of the forward and reverse filtering surfaces will be moved to the same subsea laser spot position, the subsea filtering surface F (P _T IMRI, H) can be expressed as forward filtering plane F _F (P _T IMRI, H) and inverse filtering plane F _I (P _T IMRI, H) as shown in equation (6):

F(P _T ,IMRI,H)＝F _F (P _T ,IMRI,H)∪F _I (P _T ,IMRI,H) (6)。

the invention has the beneficial technical effects that:

compared with the prior art, the invention utilizes a quadric surface LM algorithm to construct a transfer type iteration trend surface, eliminates negative abnormal data of the seabed, and constructs a bidirectional cloth simulation filtering model by combining a cloth simulation technology and bidirectional filtering so as to obtain a final seabed filtering surface, thereby solving the problem that excessive filtering is easy to generate in the prior art, effectively improving the filtering precision of the airborne laser sounding point cloud, and enhancing the robustness of sounding point cloud filtering.

Drawings

Fig. 1 is a flow chart of an airborne laser sounding point cloud filtering method based on bidirectional cloth simulation.

FIG. 2 is a schematic diagram of a transitive iterative trend surface construction in the present invention.

FIG. 3 is a flow chart of the construction of the simulation correction model of the bidirectional cloth in the invention.

Detailed Description

The invention is described in further detail below with reference to the attached drawings and detailed description:

the invention provides an airborne laser sounding point cloud filtering method based on bidirectional cloth simulation, and the flow of the filtering method is shown in figure 1.

The method comprises the following steps:

step 1: and acquiring airborne laser sounding data, and calculating the coordinates of the submarine laser points.

Acquiring sounding laser data by using an airborne LiDAR sounding system; the effective sounding laser point is the submarine laser point after attitude correction, multi-source data fusion and refraction correction.

Step 2: and (5) quadric surface fitting, and constructing a transitive iteration trend surface.

The natural topography of the sea floor is continuously variable, not abrupt, despite the gradient, as the sea dynamics erode the sea floor irregularities. Due to the fact that secondary reflection and the like of the submarine echo occur, negative abnormal points lower than the real submarine topography inevitably occur in actual submarine topography, however, in a cloth simulation filtering algorithm, filtering processing cannot be performed when negative abnormality exists. Therefore, it is necessary to establish a submarine trend surface, so that negative abnormal data can be removed, and continuity of submarine topography can be expressed. The construction of the transitive iterative trend surface is shown in fig. 2.

In a further embodiment, the step 2 specifically includes the following steps:

step 2.1: quadric LM algorithm fitting: assuming that a point P is taken as a center, r is taken as a radius, searching all submarine laser points in a given radius area, and fitting all laser points in the area to form a quadric surface model as follows:

wherein, (x) _j ,y _j ,z _j ) For the point set N _j Geographic coordinates of the laser point within (j=1, 2, …, k) in the local coordinate system; a. b, c, d, e and f are quadric fitting parameters, and fitting parameter values can be obtained based on an LM algorithm;

step 2.2: at point P ₁ Taking r as a radius for centering, searching all laser points in the area and performing quadric LM fitting. Calculating the water depth average value of laser points in curved surfaces

Standard deviation sigma of water depth _H Find and +.>

The difference is greater than 2 sigma _H The laser points are marked as negative abnormal points and stored into a suspicious point set, and the water depth values of all suspicious points are calculated into a fitted quadric surface to ensure a discrete point set N ₁ The method does not contain coarse differences and ensures the integrity of sounding data;

step 2.3: find and point P ₁ The nearest point, designated as point P ₂ Due to (N) ₁ ∩N ₂ ) The dot concentration of (2) no longer contains coarse differences, and is defined by N ₂ The quadric surface fitting accuracy is higher, and the (N ₂ -N ₁ ∩N ₂ ) Point concentration

The difference is greater than 2 sigma _H The laser points are marked as negative abnormal points and stored into a suspicious point set, and the water depth values of all suspicious points are calculated into a fitted quadric surface to ensure a discrete point set N ₂ The method does not contain coarse differences;

step 2.4: find and point P ₂ The nearest point (excluding the completed point) is denoted as point P ₃ And (2.3) repeating the step until all submarine laser points are traversed, storing no laser points into a suspicious point set, and obtaining the final submarine topography trend surface when the fitted topography trend surface reaches a stable state.

In the concrete implementation, the construction of the transmission type iteration trend surface not only effectively eliminates the negative abnormal points on the sea bottom, but also reduces the water depth values of the negative abnormal points to the fitted trend surface, thereby ensuring the integrity of sounding data and providing guarantee for smooth implementation of bidirectional cloth simulation filtering.

Step 3: and constructing a bidirectional cloth simulation correction model to obtain the submarine filtering surface.

The submarine topography is complex and changeable, when the submarine topography has more obvious positive topography (such as cliffs, submarine volcanoes, mud volcanic and reef) or negative topography (such as notching, notching and pit), the submarine laser point cloud filtering can generate more serious excessive filtering phenomenon by adopting a forward distribution simulation method, and the integrity of the submarine topography is destroyed. Therefore, the bidirectional cloth simulation correction model is adopted to carry out the filtering of the seabed laser point cloud, and the problem of excessive filtering of the point cloud is effectively avoided. The specific flow of the construction of the bidirectional cloth simulation correction model is shown in fig. 3.

In a further embodiment, the step 3 specifically includes the following steps:

step 3.1: setting filtering parameters: the grid resolution GR, namely the distance between two adjacent nodes; the stiffness parameter IMRI controls the softness of the simulated cloth; the maximum iteration number MI is used for judging whether the simulation process is terminated; judging classification thresholds of the sea bottom points and the non-sea bottom points according to the self-adaptive distance threshold H;

step 3.2: respectively performing forward filtering and backward filtering on the submarine laser points, namely respectively inverting the submarine laser points and not inverting the submarine laser points, and placing the simulated cloth right above the laser points and keeping a certain distance from the laser points; projecting the distribution node and the laser point to the same horizontal plane, performing grid-networking, determining the corresponding relation between the distribution node and the laser point, and recording the elevation values of the distribution node and the laser point;

step 3.3: the simulated cloth is enabled to do free falling movement, and the position and the stress effect of the cloth node meet the formula (2);

wherein m represents node Z _l Is the mass of (3); x represents Z _l Position at time t; g represents gravitational acceleration; c _e Representing the elastic deformation coefficient of the node spring; e represents Z _l A set of neighbors;

representing node Z _l And node Z _n Distance at time t>

Representing node Z _l And node Z _n Original distance between>

Representation ofZ _l Pointing to Z _n Is a unit vector of (a).

Step 3.4: in the process of simulating the descending of cloth, performing collision detection on cloth nodes;

assuming that the collision state is defined as C, the collision constraint is determined by the position and velocity constraints of the cloth node, and there are two time steps:

where dt is the time step; XC (t), XC (t+dt) and XC (t+2dt) are the positions of the collision at times t, t+dt and t+2dt, respectively; XC ' (t), XC ' (t+dt) and XC ' (t+2dt) are the movement velocities of the collisions at times t, t+dt and t+2dt, respectively; XC "(t) and XC" (t+dt) are the motion accelerations of the collision at times t and t+dt, respectively; Δxc "(t) and Δxc" (t+dt) are corrected accelerations of the collision at times t and t+dt, respectively, and thus equation (4) can be obtained:

wherein Δxc "(t) can be converted into additional impact forces acting on the distribution node;

step 3.5: based on collision detection, judging the position relationship between the distribution node and the submarine laser points:

if the height value H of the cloth node is satisfied _cp Less than the elevation value H of the submarine laser spot _sp Moving the node to the corresponding laser point location and marking the node as an immovable point; otherwise, calculating the displacement of the node which needs to be adjusted under the action of internal force (pulling force inside the node).

Step 3.6: repeating the steps 3.3-3.5 until the iteration times of the nodes are larger than the set maximum iteration times MI, and ending the cloth simulation process.

Step 3.7: calculating the height difference delta H between the distribution node and the corresponding submarine laser point, and judging the size relation between delta H and the self-adaptive distance threshold H; if ΔH < H, the node is reserved as a submarine point, otherwise, the node is judged as a non-submarine point and deleted.

Step 3.8: establishing a forward filtering surface F _F (P _T IMRI, H) and inverse filtering plane F _I (P _T ,IMRI,H)：

Defining a forward filtering surface as the simulated cloth formed after forward filtering, and defining a reverse filtering surface as the simulated cloth formed after reverse filtering; let the forward filtering surface be F _F (P _T IMRI, H), on which the distribution node is denoted P _F (x _F ,y _F ,z _F ) The method comprises the steps of carrying out a first treatment on the surface of the The reverse filtering surface is F _I (P _T IMRI, H), on which the distribution node is denoted P _I (x _I ,y _I ,z _I ) The method comprises the steps of carrying out a first treatment on the surface of the The two filtering surfaces can be expressed as

Wherein F is _F (P _T IMRI, H) and F _I (P _T IMRI, H) is the forward and reverse filtering face; p (P) _T (x _T ,y _T ,z _T )(z _T ＝f(x _T ,y _T ) A submarine laser spot); the IMRI is a rigidity parameter for controlling the softness of the cloth; h is a distance threshold value for judging that the distribution node is a submarine point and a non-submarine point; (x) _F ,y _F ,z _F ) And (x) _I ,y _I ,z _I ) For distributing the node P _F And P _I Is satisfied with z _F ＝f _F (x _F ,y _F ) And z _I ＝f _I (x _I ,y _I ) The method comprises the steps of carrying out a first treatment on the surface of the T represents the sea floor topography;

step 3.9: establishing a submarine filtering surface F (P) _T ,IMRI,H)：

Since some nodes of the forward and reverse filtering surfaces will be moved to the same subsea laser spot position, the subsea filtering surface F (P _T IMRI, H) can be expressed as forward filtering plane F _F (P _T IMRI, H) and inverse filtering plane F _I (P _T Union of IMRI, H), i.e

F(P _T ,IMRI,H)＝F _F (P _T ,IMRI,H)∪F _I (P _T ,IMRI,H) (6)。

In the concrete implementation, the submarine filtering surface is constructed according to the bidirectional cloth simulation correction model, so that non-submarine points are effectively filtered, the submarine micro-topography characteristics are fully reserved, and the integrity of submarine topography data is ensured. Therefore, the submarine filtering surface can accurately express the actual submarine topography, and the bidirectional cloth simulation filtering method has better applicability and robustness for airborne laser sounding data.

In summary, the invention provides an airborne laser sounding point cloud filtering method based on bidirectional cloth simulation, which comprises the following steps: on the premise of obtaining the submarine laser point coordinates, adopting a quadric LM algorithm to fit and construct a transfer type iteration trend surface, removing submarine negative abnormal data, and expressing the submarine continuity; and constructing a bidirectional distribution simulation correction model, and performing bidirectional filtering on the submarine laser point cloud according to the correction model to obtain a final submarine filtering surface. The method effectively solves the problems that the existing method can not identify negative abnormal data of airborne laser sounding and is easy to generate excessive filtering, and realizes high-precision and strong-robustness filtering of the airborne laser sounding point cloud.

It should be understood that the above description is not intended to limit the invention to the particular embodiments disclosed, but to limit the invention to the particular embodiments disclosed, and that the invention is not limited to the particular embodiments disclosed, but is intended to cover modifications, adaptations, additions and alternatives falling within the spirit and scope of the invention.

Claims (2)

1. An airborne laser sounding point cloud filtering method based on bidirectional cloth simulation is characterized in that: the method comprises the following steps:

step 1: acquiring airborne laser sounding point cloud data, and calculating the coordinates of the submarine laser points;

step 2: adopting quadric LM algorithm fitting, constructing a transfer type iteration trend surface, removing negative abnormal data of the seabed, and expressing the continuity of the seabed;

step 3: constructing a bidirectional distribution simulation correction model, and performing bidirectional filtering of the submarine laser point cloud according to the correction model to obtain a final filtering surface;

the specific construction method of the bidirectional cloth simulation correction model is as follows:

step 3.1: setting filtering parameters: the grid resolution GR, namely the distance between two adjacent nodes; the stiffness parameter IMRI controls the softness of the simulated cloth; the maximum iteration number MI is used for judging whether the simulation process is terminated; judging classification thresholds of the sea bottom points and the non-sea bottom points according to the self-adaptive distance threshold H;

step 3.2: respectively performing forward filtering and backward filtering on the submarine laser points, namely respectively inverting the submarine laser points and not inverting the submarine laser points, and placing the simulated cloth right above the laser points and keeping a certain distance from the laser points; projecting the distribution node and the laser point to the same horizontal plane, performing grid-networking, determining the corresponding relation between the distribution node and the laser point, and recording the elevation values of the distribution node and the laser point;

step 3.3: the simulated cloth is enabled to do free falling motion, and the position and the stress effect of the cloth node meet the formula (2):

wherein m represents node Z _l Is the mass of (3); x represents Z _l Position at time t; g represents gravitational acceleration; c _e Representing the elastic deformation coefficient of the node spring; e represents Z _l A set of neighbors;

representing node Z _l And node Z _n Distance at time t>

Representing node Z _l And node Z _n Original distance between>

Representing Z _l Pointing to Z _n Is a unit vector of (a);

step 3.4: in the process of simulating the descending of cloth, performing collision detection on cloth nodes;

assuming that the collision state is defined as C, the collision constraint is determined by the position and velocity constraints of the cloth node, and there are two time steps:

where dt is the time step; XC (t), XC (t+dt) and XC (t+2dt) are the positions of the collision at times t, t+dt and t+2dt, respectively; XC ' (t), XC ' (t+dt) and XC ' (t+2dt) are the movement velocities of the collisions at times t, t+dt and t+2dt, respectively; XC "(t) and XC" (t+dt) are the motion accelerations of the collision at times t and t+dt, respectively; Δxc "(t) and Δxc" (t+dt) are corrected accelerations of the collision at times t and t+dt, respectively, and thus equation (4) can be obtained:

wherein Δxc "(t) can be converted into additional impact forces acting on the distribution node;

step 3.5: judging the position relation between the distribution node and the submarine laser points on the basis of collision detection;

if the elevation value H of the material distribution node _cp Elevation value H smaller than the undersea laser spot _sp Moving the node to the corresponding laser point location and marking the node as an immovable point; otherwise, calculating the displacement of the node to be adjusted under the action of the traction force in the node;

step 3.6: repeating the steps 3.3 to 3.5 until the iteration times of the nodes are larger than the set maximum iteration times MI, and ending the cloth simulation process;

step 3.7: calculating the height difference delta H between the distribution node and the corresponding submarine laser point, and judging the size relation between delta H and the self-adaptive distance threshold H; if delta H is less than H, the node is reserved as a submarine point, otherwise, the node is judged to be a non-submarine point and deleted;

step 3.8: establishing a forward filtering surface F _F (P _T IMRI, H) and inverse filtering plane F _I (P _T ,IMRI,H)；

Defining a forward filtering surface as the simulated cloth formed after forward filtering, and defining a reverse filtering surface as the simulated cloth formed after reverse filtering; let the forward filtering surface be F _F (P _T IMRI, H), on which the distribution node is denoted P _F (x _F ,y _F ,z _F ) The method comprises the steps of carrying out a first treatment on the surface of the The reverse filtering surface is F _I (P _T IMRI, H), on which the distribution node is denoted P _I (x _I ,y _I ,z _I ) The method comprises the steps of carrying out a first treatment on the surface of the The two filter surfaces are represented as shown in equation (5):

wherein F is _F (P _T IMRI, H) and F _I (P _T IMRI, H) is the forward and reverse filtering face; p (P) _T (x _T ,y _T ,z _T ) Z is a submarine laser spot _T ＝f(x _T ,y _T ) The method comprises the steps of carrying out a first treatment on the surface of the The IMRI is a rigidity parameter for controlling the softness of the cloth; h is a distance threshold value for judging that the distribution node is a submarine point and a non-submarine point; (x) _F ,y _F ,z _F ) And (x) _I ,y _I ,z _I ) For distributing the node P _F And P _I Is satisfied with z _F ＝f _F (x _F ,y _F ) And z _I ＝f _I (x _I ,y _I ) The method comprises the steps of carrying out a first treatment on the surface of the T represents the sea floor topography;

step 3.9: establishing a submarine filtering surface F (P) _T ,IMRI,H)；

Since some nodes of the forward and reverse filtering surfaces will be moved to the same subsea laser spot position, the subsea filtering surface F (P _T IMRI, H) can be expressed as forward filtering plane F _F (P _T IMRI, H) and inverse filtering plane F _I (P _T IMRI, H) as shown in equation (6):

F(P _T ,IMRI,H)＝F _F (P _T ,IMRI,H)∪F _I (P _T ,IMRI,H) (6)。

2. the airborne laser sounding point cloud filtering method based on bidirectional cloth simulation according to claim 1, wherein the method is characterized in that: in step 2, the construction method of the transitive iterative trend surface is as follows:

step 2.1: quadric LM algorithm fitting: assuming that a point P is taken as a center, r is taken as a radius, searching all submarine laser points in a given radius area, and a quadric surface model formed by fitting all laser points in the area is shown as a formula (1):

wherein, (x) _j ,y _j ,z _j ) For the point set N _j Geographic coordinates of the laser points in the local coordinate system, j=1, 2, …, k; a. b, c, d, e and f are quadric fitting parameters;

step 2.2: at point P ₁ Taking r as a radius as a center, searching all laser points in the area and performing quadric surface LM fitting; calculating the water depth average value of laser points in curved surfaces

Standard deviation sigma of water depth _H Find and +.>

The difference is greater than 2 sigma _H The laser points are marked as negative abnormal points and stored into a suspicious point set, and the water depth values of all suspicious points are calculated into a fitted quadric surface to ensure a discrete point set N ₁ The method does not contain coarse differences and ensures the integrity of sounding data;

step 2.3: find and point P ₁ The nearest point, designated as point P ₂ Find (N) ₂ -N ₁ ∩N ₂ ) Point concentration

The difference is greater than 2 sigma _H The laser points are marked as negative abnormal points and stored into a suspicious point set, and the water depth values of all suspicious points are calculated into a fitted quadric surface to ensure a discrete point set N ₂ The method does not contain coarse differences;

step 2.4: find and point P ₂ The nearest point, excluding the completed point, is denoted as point P ₃ And (2.3) repeating the step until all submarine laser points are traversed, storing no laser points into a suspicious point set, and obtaining the final submarine topography trend surface when the fitted topography trend surface reaches a stable state.

Images