CN110246092A - A kind of three-dimensional laser point cloud denoising method for taking neighborhood point mean distance and slope into account - Google Patents
A kind of three-dimensional laser point cloud denoising method for taking neighborhood point mean distance and slope into account Download PDFInfo
- Publication number
- CN110246092A CN110246092A CN201910365632.0A CN201910365632A CN110246092A CN 110246092 A CN110246092 A CN 110246092A CN 201910365632 A CN201910365632 A CN 201910365632A CN 110246092 A CN110246092 A CN 110246092A
- Authority
- CN
- China
- Prior art keywords
- point
- slope
- gamma
- mean value
- formula
- Prior art date
- Legal status (The legal status is an assumption and is not a legal conclusion. Google has not performed a legal analysis and makes no representation as to the accuracy of the status listed.)
- Granted
Links
- 238000000034 method Methods 0.000 title claims description 23
- 230000001186 cumulative effect Effects 0.000 claims description 10
- 230000008569 process Effects 0.000 claims description 8
- 230000008859 change Effects 0.000 claims description 7
- 238000005315 distribution function Methods 0.000 claims description 7
- 238000013480 data collection Methods 0.000 claims description 4
- 238000013507 mapping Methods 0.000 claims description 3
- AYFVYJQAPQTCCC-GBXIJSLDSA-N L-threonine Chemical compound C[C@@H](O)[C@H](N)C(O)=O AYFVYJQAPQTCCC-GBXIJSLDSA-N 0.000 claims description 2
- 238000013461 design Methods 0.000 claims description 2
- 239000004744 fabric Substances 0.000 claims description 2
- 238000009825 accumulation Methods 0.000 claims 1
- 230000000694 effects Effects 0.000 abstract description 21
- 238000001914 filtration Methods 0.000 abstract description 17
- 230000006870 function Effects 0.000 description 12
- 238000002474 experimental method Methods 0.000 description 4
- 238000005516 engineering process Methods 0.000 description 3
- 238000007619 statistical method Methods 0.000 description 2
- 238000012360 testing method Methods 0.000 description 2
- 238000010146 3D printing Methods 0.000 description 1
- 238000004458 analytical method Methods 0.000 description 1
- 230000008901 benefit Effects 0.000 description 1
- 230000002146 bilateral effect Effects 0.000 description 1
- 238000005094 computer simulation Methods 0.000 description 1
- 238000005260 corrosion Methods 0.000 description 1
- 230000007797 corrosion Effects 0.000 description 1
- 238000007405 data analysis Methods 0.000 description 1
- 238000012217 deletion Methods 0.000 description 1
- 230000037430 deletion Effects 0.000 description 1
- 238000001514 detection method Methods 0.000 description 1
- 238000010586 diagram Methods 0.000 description 1
- 230000008030 elimination Effects 0.000 description 1
- 238000003379 elimination reaction Methods 0.000 description 1
- 239000000284 extract Substances 0.000 description 1
- 239000003999 initiator Substances 0.000 description 1
- 239000000203 mixture Substances 0.000 description 1
- VMXUWOKSQNHOCA-UKTHLTGXSA-N ranitidine Chemical compound [O-][N+](=O)\C=C(/NC)NCCSCC1=CC=C(CN(C)C)O1 VMXUWOKSQNHOCA-UKTHLTGXSA-N 0.000 description 1
- 230000009467 reduction Effects 0.000 description 1
- 238000011160 research Methods 0.000 description 1
- 238000006467 substitution reaction Methods 0.000 description 1
- 238000013519 translation Methods 0.000 description 1
- 239000002023 wood Substances 0.000 description 1
Classifications
-
- G—PHYSICS
- G06—COMPUTING; CALCULATING OR COUNTING
- G06T—IMAGE DATA PROCESSING OR GENERATION, IN GENERAL
- G06T5/00—Image enhancement or restoration
- G06T5/70—Denoising; Smoothing
-
- G—PHYSICS
- G06—COMPUTING; CALCULATING OR COUNTING
- G06T—IMAGE DATA PROCESSING OR GENERATION, IN GENERAL
- G06T2207/00—Indexing scheme for image analysis or image enhancement
- G06T2207/10—Image acquisition modality
- G06T2207/10028—Range image; Depth image; 3D point clouds
-
- G—PHYSICS
- G06—COMPUTING; CALCULATING OR COUNTING
- G06T—IMAGE DATA PROCESSING OR GENERATION, IN GENERAL
- G06T2207/00—Indexing scheme for image analysis or image enhancement
- G06T2207/20—Special algorithmic details
- G06T2207/20024—Filtering details
Landscapes
- Physics & Mathematics (AREA)
- General Physics & Mathematics (AREA)
- Engineering & Computer Science (AREA)
- Theoretical Computer Science (AREA)
- Image Processing (AREA)
- Complex Calculations (AREA)
Abstract
The present invention is arranged by rule of thumb for the threshold value of current neighbor point points does not ensure that applicability, and the capability of fitting of common filtering algorithm is limited, have some limitations equal series of problems, to ground three-dimensional laser point cloud (Terrestrial Laser Scanner, TLS) the denoising strategy of neighbor point distance and slope is taken in proposition into account, under the support of k-d tree index construct, the building for individually considering the one-dimensional Gamma denoising function of neighbor point distance and slope is realized, and takes the building of the two-dimension gamma denoising function of neighbor point distance and slope into account.Two constraint conditions of filter proposed by the present invention consider query point surrounding neighbors quantity, neighborhood the point mean value and situations such as distribution density at a distance from central point from different directions, for TLS data outliers, the noise types such as either sparse exceptional value or isolated noise point are effective, can reach the technical effect for filtering out noise spot automatically.
Description
Technical field
The present invention is one and belongs to mapping, Geographical Information Sciences field.
Background technique
In recent years, computer technology and three-dimensional digital information initiator are advanced by leaps and bounds, and pass through the modes such as three-dimensional silk profoundness
Object three dimensional point model data is obtained, and reverse modeling is carried out to it and has become current research hotspot.Three-dimensional point modeling is outstanding
It is used in the fields such as 3D printing, virtual reality and computer modeling.However point data is being obtained by spatial digitizer
In the process, the three-dimensional data measured always will appear random error, and the irregularity on testee surface itself is but also measure
Data information there are certain noise informations, therefore in modeling process, accuracy detection and denoising are carried out to the data obtained
Work plays the role of modeling result accuracy vital.Have much to the denoising method of point cloud model in the prior art,
Such as 1) mobile camber fitting point cloud denoising method, this method take full advantage of atural object or topographical surface in the company of micro- a small range
Continuous property, that is, there is second order can lead characteristic, can be approached with quadratic surface.2) the point cloud denoising based on Mathematical Morphology theory
The size of feature extracts characteristics of image in method " reduction " (corrosion) or " increase " (expansion) image.3) bilateral filtering point cloud
Picture noise not only can be effectively reduced in denoising method, but also can keep the algorithm of edge detail information to the maximum extent.4) a kind of base
In the isolated noise point elimination method of statistical analysis, the point cloud data of above ground structure is organized first with k-d tree first,
Secondly, threshold value is arranged using the method for statistical analysis in the principle based on normal distribution automatically, to filter out isolated noise point.
The present invention is arranged by rule of thumb for the threshold value of current neighbor point points does not ensure that applicability, and common filter
The capability of fitting of wave algorithm is limited, equal series of problems is had some limitations, to ground three-dimensional laser point cloud
(Terrestrial Laser Scanner, TLS) proposes to take the denoising strategy of neighbor point distance and slope into account, in k-d tree
Under index construct is supported, the building for individually considering the one-dimensional Gamma denoising function of neighbor point distance and slope is realized, and care for
And the building of the two-dimension gamma denoising function of neighbor point distance and slope.Two constraint conditions of filter proposed by the present invention
Query point surrounding neighbors quantity, neighborhood point feelings such as mean value and distribution density at a distance from central point are considered from different directions
Condition, for TLS data outliers, the noise types such as either sparse exceptional value or isolated noise point effectively be can reach certainly
The dynamic technical effect for filtering out noise spot.
Summary of the invention
The present invention relates to a kind of three-dimensional laser point cloud denoising methods for taking neighborhood point mean distance and slope into account, by following step
Rapid composition:
Step 1, original point cloud k-d tree index construct:
(1) original point cloud data collection P={ P is read1,P2,…,Pn, be deposited into three-dimensional array, to each dimension into
Row variance calculates, and selects the maximum dimension t of variance (t=1,2,3), then selects the intermediate value μ of all numerical value in t dimension as the
One comparison other, i.e. division axis, obtain two subsets, and store to a tree node N (t, μ)；
(2) (1) process in step 1 is repeated to two obtained subsets, until raw data set P={ P1,P2,…,Pn}
It can not divide again；If some subset can not be subdivided in the process, data in the subset are saved to leaf node.
Step 2, the nearest-neighbor lookup algorithm design of k-d tree index:
(1) by data point P to be checkedi(x, y, z) and the μ value pair in each tree node N (t, μ) since root node
Than entering right subtree and inquiring if the value in t dimension is greater than μ；Otherwise, it is inquired into left subtree, after reaching leaf node,
Calculate PiThe distance between the data point stored in (x, y, z) and this leaf node d, and record corresponding coordinate of ground point Po
(x,y,z)；
(2) back tracking operation judges in not visited data point with the presence or absence of distance Pi(x, y, z) closer point, sentences
Disconnected PiWhether the distance between not visited data point is greater than d under (x, y, z) and father node belonging to it, if it is, point
It is not present in branch and Pi(x, y, z) closer point；Otherwise, exist from Pi(x, y, z) closer point then accesses the point, repeats step
(1) in 2, the point P that (1) retains in replacement step 2o(x, y, z) is not present and P until tracing back in root nodei(x, y, z) more
Until close point；
Step 3, calculate query point around it between m neighborhood point at a distance from mean value and the overall situation vertex neighborhood apart from mean value:
(1) to original point cloud data collection P={ P1,P2,…,PnIn each data point Pi(x, y, z), with above-mentioned k-
D tree index carries out m neighborhood search, obtains nearest point set Q={ P1,P2,…,Pj…,Pm}；
(2) P is calculatediThe Euclidean distance d of (x, y, z) point each point into Qij, as formula 1. shown in, then calculate arrive all neighbours
Domain point apart from mean valueAs formula 2. shown in:
(3) point set P={ P is 1. 2. calculated with formula according to formula1,P2,…,PnIn each point m neighboring mean value
Step 4, to all neighborhoods of a point apart from mean valueCarry out gamma fitting of distribution:
(1) all neighborhoods of a point are constructed apart from mean valueGamma distribution probability density function；
To point set P={ P1,P2,…,PnMidpoint PiThe m neighboring mean value of (x, y, z) carries out gamma fitting of distribution, asks
Solve gamma profile parameter1And β1, Gamma distribution stochastic variable probability distributing density function such as formula 3. shown in:
α in formula1,β1, δ is respectively form parameter, scale parameter and the location parameter of gamma distribution, and wherein δ is model point
The lower limit of cloth is all larger than zero apart from mean value, stochastic variableIt is actually distributed in x-axis positive axis, then its probability density function
Real Equivalent Form is 4.:
(2) all the points neighboring mean value distance is establishedCumulative probability F1Function expression；
All the points neighboring mean value distanceCumulative probability F1Function expression such as formula 5. shown in:
(3) cumulative probability is selected to reach certain threshold valueAs the condition of preliminary judgement noise spot, determine in noise spot
It in the process, is more than cumulative probability F1, indicate query point PiThe mean value of m neighborhood point of (x, y, z) surrounding is larger, that is, is greater than
Such center query point can tentatively regard as noise spot；
Step 5, each query point P is solvediThe distance change slope k of m neighborhood point of (x, y, z) surroundingi, intercept biAnd
All the points slope ki':
(1) to the Euclidean distance d of m nearest-neighbor point around query pointijLinear fit, such as 6. formula are carried out, is solved every
A query point PiThe distance change slope k of m neighborhood point of (x, y, z) surroundingi:
yi=kidij+bi ⑥
(2) point set P={ P is 6. calculated according to formula1,P2,…,PnIn all the points slope ki′
Step 6, to all the points slope ki' carry out gamma fitting of distribution:
(1) to point set P={ P1,P2,…,PnIn point PiThe slope k of (x, y, z)i' progress gamma distribution is quasi-
Close, as formula 7. shown in, solve gamma be distributed α2And β2:
Similarly, α in formula2,β2The form parameter and scale parameter of respectively gamma distribution；
(2) all the points neighboring mean value slope k is establishedi' cumulative probability F2Function expression, cumulative probability F2Function
Expression formula is as shown in 8.:
(3) query point P is calculatediThe Euclidean distance of (x, y, z) surrounding neighbors point changes slope kiCumulative probability, if being more than
Cumulative probability F2, indicate query point PiThe distance change of m neighborhood point of (x, y, z) surrounding is larger, i.e., slope is greater than given threshold value
k′thre, show that the neighborhood point distributed pole around query point is uneven, such query point is identified as noise spot；
Step 7, the neighboring mean value distance of all the points is establishedWith neighboring mean value apart from slope ki' two-dimension gamma joint
Probability density function:
It (1), can be more by the two joint mapping probability density function since two above constraint condition is mutually indepedent
Noise spot is removed well, is based onAnd ki' two-dimension gamma joint probability density function such as the formula established is 9. shown:
Wherein, α1,β1,α2,β2It is all larger than zero, is claimedObedience parameter is α1,β1,α2,β2The two-dimension gamma of distribution joins
Distribution is closed, then 10. two-dimension gamma joint distribution function meets formula:
Separately there is the binary function in its distribution function to meet formula
(2) according to the property of two-dimension gamma Joint Distribution, setting meets the sky between a certain plane G and curved surface f (x, y)
Between in region point be in point, be otherwise exterior point, i.e. noise spot；
Step 8, all noise spots for not meeting step (7) are removed to get filtered point set is arrived.
Detailed description of the invention
Fig. 1 original point cloud data；
Fig. 2 inquires the joint apart from the histogram of mean value and slope k, matched curve and the two that vertex neighborhood number is 10
Probability density；
Fig. 3 inquires the joint apart from the histogram of mean value and slope k, matched curve and the two that vertex neighborhood number is 30
Probability density；
Fig. 4 inquires the joint apart from the histogram of mean value and slope k, matched curve and the two that vertex neighborhood number is 40
Probability density；
Fig. 5 inquires the joint apart from the histogram of mean value and slope k, matched curve and the two that vertex neighborhood number is 50
Probability density；
Fig. 6 inquires the joint apart from the histogram of mean value and slope k, matched curve and the two that vertex neighborhood number is 60
Probability density；
Each Parameters variation situation of the neighborhood apart from mean value when Fig. 7 difference neighborhood number；
Each Parameters variation situation of slope when Fig. 8 difference neighborhood number；
Number and variance variation are put when Fig. 9 difference neighborhood number after Gamma filter filtering；
When Figure 10 m=20, filtered result figure；
When Figure 11 m=30, filtered result figure；
When Figure 12 m=40, filtered result figure；
Figure 13 method flow diagram.
Wherein, (a) and (b) is schemed in Fig. 1 to be respectively as follows:
(a) original point cloud data front view；(b) original point cloud data side view.
Scheme (a) and (b) in Figure 10 to be respectively as follows:
(a) when m=20, effect picture is faced after filtering；
(b) when m=20, rear side viewing effect picture is filtered；
(a) is schemed in Figure 11 and figure (b) is respectively as follows:
(a) when m=30, effect picture is faced after filtering；
(b) when m=30, rear side viewing effect picture is filtered；
(a) is schemed in Figure 12 and figure (b) is respectively as follows:
(a) when m=40, effect picture is faced after filtering；
(b) when m=40, rear side viewing effect picture is filtered；
Specific embodiment
In order to verify the validity of this method, Institutes Of Technology Of Jiangxi's administrative building front point cloud data is acquired, such as
Include data point 820025 shown in Fig. 1 (a) and Fig. 1 (b), includes building, car, the short tree in periphery in original point cloud
The elements such as wood, scene is complex, and ground species are more, and the complicated performance for studying scene effectively verifies constructed filtering of the invention
The robustness of device.By being found to original point cloud data analysis, due to the factors such as dirt in air produce a large amount of non-ground and
Non- culture point cloud, these noise spots bring to geographic scenes three-dimensional modeling and greatly interfere with.The present invention tests hardware environment:
Intel (R) Xeon (R) CPU E3-1231v3 3.40GHz, interior 8GB is saved as；It is in MATLAB that Gamma filter, which tests algorithm,
Under 2017 translation and compiling environment of R2017 and MicrosoftVisual Studio, programs and realize in conjunction with PCL method base.
Determine ideal neighborhood point (m) number according to many experiments result, the present invention take respectively m=10,20,30,40,50,
When 60, experimental analysis is carried out, while comparing with the normal distribution experiment effect of the prior art four.Inventive algorithm and existing
The normal distribution algorithm of technology four under different neighborhood value conditions apart from mean value and neighborhood slope histogram, Gamma function
Gamma joint probability density, the normal distribution joint probability of matched curve, normal distyribution function matched curve and the two are close
Degree experiment effect is respectively as shown in Fig. 2~Fig. 6.
Normal distribution in order to preferably embody Gamma filter denoising effect proposed by the present invention, with the prior art four
Filter effect is compared from the joint probability density of matched curve and the two apart from mean value and the histogram of slope k respectively
Compared with when gamma distribution curve and normal distribution curve are that cumulative probability reaches 95%, corresponding neighborhood is apart from mean value and tiltedly
Rate value as shown in Fig. 2~Fig. 6,.
From Fig. 2~Fig. 6, when taking different neighborhood numbers, the neighboring mean value of query point and the threshold value value of slope all can
Change, when m takes 10,30,40,50,60, the value of neighboring mean value is respectively 0.13573,0.1979,0.2962,
0.32793,0.36866, the value of slope is respectively 0.0038785,0.0031125,0.0024233,0.0022197,
0.0020282.When being fitted respectively to neighboring mean value and slope histogram, the changing rule and histogram of gamma curve
Unanimously, and when accumulated probability reaches 95%, preferential normal distribution convergence, and the horizontal axis value of Gamma filter is all larger than zero,
Fitting effect more meets the actual conditions of neighboring mean value and slope greater than zero.As can be known from Fig. 7, when selection same vicinity number
When, constraint condition single for neighboring mean value filters out noise spot number than using when being filtered using gamma curve
Normal distribution curve filters out noise spot number when filtering is more, illustrates that Gamma filter can filter out more under the same conditions
Noise spot.As can be known from Fig. 8, when choosing same vicinity number, constraint condition single for slope, using gamma curve
The number for filtering out noise spot is filtered than using the number for filtering out noise spot when normal distribution curve filtering more.Cause
This is better than regardless of being the neighboring mean value or slope around query point when single condition using the effect that gamma curve filters
It is filtered using normal distribution curve.
The matched curve of Joint Distribution is it is found that when taking different neighborhood numbers from Fig. 2~Fig. 6, Gamma joint probability
Density map be it is convergent, in central area, this distribution situation is more conducive to selected threshold and is filtered place in high information point region
Reason, and two dimension joint normal distyribution function high information area do not converge on a certain central area, when setting neighboring mean value or
It can not include most of information content when the threshold value of person's slope.Meanwhile it is found that choosing different neighbours from Fig. 2~Fig. 6
When the number of domain, when to neighboring mean value and slope using two dimension joint Gamma fitting of distribution, center of curve is to belong to certain neighboring mean value
Query point number it is more, and the region that the query point number that belongs to certain slope is more, threshold range (i.e. red block inner region)
Outer point is noise spot.The neighboring mean value of query point is excessive, illustrates query point and other point distances farther out, it is believed that be noise
Point, as m=10, is observed along the horizontal axis of Gamma Joint Distribution in Fig. 2, works as query pointWhen, it is believed that this query point is noise spot；The slope of query point is bigger, illustrates query point week
It encloses neighborhood point to be unevenly distributed, in Fig. 2, as m=10, along the axial spotting of Gamma Joint Distribution, when query point
Slope k ' > k 'threWhen=0.0020282, it is believed that this query point is noise spot.
When Gamma filter being selected to carry out noise spot removal, different neighborhood point numbers is chosen, filter effect is different, is
Reach preferable filter effect, the preferable filter effect of effect is chosen in the neighborhood point number different from 6.From Fig. 7
The relationship between number and neighborhood point number is put after Gamma filtering it is found that when to choose neighboring mean value be constraint condition, when neighborhood point
When number is 30, it is more to filter out noise spot quantity.The relationship between number and neighborhood point number is put after Gamma filtering in Fig. 8
It is found that when selection slope is constraint condition, although the quantity for filtering out noise spot is more as m=20, as m=30, filter
Except the quantity of noise spot is not much different therewith, in addition, from β2Distribution from the point of view of, the β after m=302Variation is than more gentle.From Fig. 9
In it is found that when choosing two dimension joint Gamma probability density distribution, i.e. when two kinds of constraint condition collective effects, as m=30,
The quantity for filtering out noise spot is most.
Therefore, in conclusion when choosing m=30, in neighborhood under two kinds of constraint conditions of mean value and slope, joint
Gamma filter can effectively filter out noise spot.Experiment effect as shown in figs. 10,11 and 12, from the figure, it can be seen that when m takes
When value is 20, there is partial noise point to be not filtered out (red block), when m value is 40, part details is filtered out (blue on metope
Frame), when m value is 30, noise spot is filtered out completely, and can preferably retain metope detailed information, and filter effect is best.
The present invention can be used any Object-Oriented Programming Language and realize that now the implementation by taking C Plus Plus as an example is as follows:
(1) open source point cloud library PCL (Point Cloud Library, PCL) creation point cloud data is called to read and write class:
PCReadClass and PCWriteClass；
(2) creation point cloud k-d tree index creation class and neighborhood point search class: kdTreeClass and
NPFindClass；
(3) creation query point and neighborhood point mean distance calculate class: AvgDistCalcClass；
(4) creation query point and neighborhood point slope mean value computation class: AvgKCalcClass；
(5) creation all the points neighboring mean value calculates class apart from probability density Gamma distribution function:
AvgDistGammaClass；
(6) creation all the points neighboring mean value slope probability density Gamma distribution function calculates class: AvgKGammaClass；
(7) it creates all the points neighboring mean value distance and neighboring mean value slope two-dimensional Gamma distribution function calculates class:
DistAndKGammaClass；
(8) AvgDistGammaClass is called to judge whether current queries point mean distance is greater than threshold value, if it is greater than threshold
Value, then be stored into preliminary candidate point set Q；If current queries point mean distance is less than threshold value, call
AvgKGammaClass judges whether current queries vertex neighborhood slope mean value is greater than given threshold value, if it is greater than threshold value, then stores
Into preliminary candidate point set Q；
(9) DistAndKGammaClass is called further to carry out noise point deletion to candidate point set Q.
The foregoing is only a preferred embodiment of the present invention, but scope of protection of the present invention is not limited thereto,
In the technical scope disclosed by the present invention, any changes or substitutions that can be easily thought of by anyone skilled in the art,
All it is included within the scope of protection of the present invention.Therefore, protection scope of the present invention should be with the protection scope of claims
It is quasi-.
Claims (1)
1. a kind of three-dimensional laser point cloud denoising method for taking field point mean distance and slope into account, it is characterised in that have following step
It is rapid:
Step 1, original point cloud k-d tree index construct:
(1) original point cloud data collection P={ P is read1,P2,…,Pn, it is deposited into three-dimensional array, to each dimension progress side
Difference calculates, and selects the maximum dimension t of variance (t=1,2,3), then selects the intermediate value μ of all numerical value in t dimension as first
Comparison other, i.e. division axis, obtain two subsets, and store to a tree node N (t, μ)；
(2) (1) process in step 1 is repeated to two obtained subsets, until raw data set P={ P1,P2,…,PnCan not
Divide again；If some subset can not be subdivided in the process, data in the subset are saved to leaf node；
Step 2, the nearest-neighbor lookup algorithm design of k-d tree index:
(1) by data point P to be checkedi(x, y, z) is compared with the μ value in each tree node N (t, μ) since root node, if
Value in t dimension is greater than μ, then enters right subtree and inquire；Otherwise, it is inquired into left subtree, after reaching leaf node, calculates Pi
The distance between the data point stored in (x, y, z) and this leaf node d, and record corresponding coordinate of ground point Po(x,y,
z)；
(2) back tracking operation judges in not visited data point with the presence or absence of distance Pi(x, y, z) closer point, judges Pi
Whether the distance between not visited data point is greater than d under (x, y, z) and father node belonging to it, if it is, in branch
It is not present and Pi(x, y, z) closer point；Otherwise, exist from Pi(x, y, z) closer point then accesses the point, repeats in step 2
(1), the point P that (1) retains in replacement step 2o(x, y, z) is not present and P until tracing back in root nodei(x, y, z) is closer
Point until；
Step 3, calculate query point around it between m neighborhood point at a distance from mean value and the overall situation vertex neighborhood apart from mean value:
(1) to original point cloud data collection P={ P1,P2,…,PnIn each data point Pi(x, y, z), with above-mentioned k-d tree
Index carries out m neighborhood search, obtains nearest point set Q={ P1,P2,…,Pj…,Pm}；
(2) P is calculatediThe Euclidean distance d of (x, y, z) point each point into Qij, as formula 1. shown in, then calculate arrive all neighborhood points
Apart from mean valueAs formula 2. shown in:
(3) point set P={ P is 1. 2. calculated with formula according to formula1,P2,…,PnIn each point m neighboring mean value
Step 4, to all neighborhoods of a point apart from mean valueCarry out Gamma fitting of distribution:
(1) all neighborhoods of a point are constructed apart from mean valueGamma distribution probability density function；
To point set P={ P1,P2,…,PnMidpoint PiThe m neighboring mean value of (x, y, z) carries out Gamma fitting of distribution, solves
Gamma profile parameter1And β1, Gamma distribution stochastic variable probability distributing density function such as formula 3. shown in:
α in formula1,β1, δ is respectively form parameter, scale parameter and the location parameter of Gamma distribution, and wherein δ is model profile
Lower limit is all larger than zero apart from mean value, stochastic variableIt is actually distributed in x-axis positive axis, then its probability density function reality etc.
Valence formula is 4.:
(2) all the points neighboring mean value distance is establishedCumulative probability F1Function expression；
All the points neighboring mean value distanceCumulative probability F1Function expression such as formula 5. shown in:
(3) cumulative probability is selected to reach certain threshold valueAs the condition of preliminary judgement noise spot, in noise spot decision process
In, it is more than cumulative probability F1, indicate query point PiThe mean value of m neighborhood point of (x, y, z) surrounding is larger, that is, is greater thanIt is such
Center query point can tentatively regard as noise spot；
Step 5, each query point P is solvediThe distance change slope k of m neighborhood point of (x, y, z) surroundingi, intercept biAnd it is all
Point slope ki':
(1) to the Euclidean distance d of m nearest-neighbor point around query pointijLinear fit, such as 6. formula are carried out, each look into is solved
Ask point PiThe distance change slope k of m neighborhood point of (x, y, z) surroundingi:
yi=kidij+bi ⑥
(2) point set P={ P is 6. calculated according to formula1,P2,…,PnIn all the points slope ki'
Step 6, to all the points slope ki' carry out gamma curve fitting of distribution:
(1) to point set P={ P1,P2,…,PnIn point PiThe slope k of (x, y, z)i' gamma curve fitting of distribution is carried out, such as formula
Shown in 7., the α of Gamma distribution is solved2And β2:
Similarly, α in formula2,β2The form parameter and scale parameter of respectively Gamma distribution；
(2) all the points neighboring mean value slope k is establishedi' cumulative probability F2Function expression, cumulative probability F2Function representation
Formula is as shown in 8.:
(3) query point P is calculatediThe Euclidean distance of (x, y, z) surrounding neighbors point changes slope kiCumulative probability, if be more than accumulation
Probability F2, indicate query point PiThe distance change of m neighborhood point of (x, y, z) surrounding is larger, i.e., slope is greater than given threshold value k'thre,
Show that the neighborhood point distributed pole around query point is uneven, such query point is identified as noise spot；
Step 7, the neighboring mean value distance of all the points is establishedWith neighboring mean value apart from slope ki' two-dimension gamma joint probability
Density function:
It (1), can be preferably by the two joint mapping probability density function since two above constraint condition is mutually indepedent
Noise spot is removed, is based onAnd ki' establish two-dimension gamma joint probability density function such as formula 9. shown in:
Wherein, α1,β1,α2,β2It is all larger than zero, is claimedObedience parameter is α1,β1,α2,β2The two-dimension gamma joint point of distribution
Cloth, then 10. two-dimension gamma joint distribution function meets formula:
Separately there is the binary function in its distribution function to meet formula
F (X, Y)=P { X≤x, Y≤y }
(2) according to the property of two-dimension gamma Joint Distribution, setting meets the space region between a certain plane G and curved surface f (x, y)
Point is interior point in domain, is otherwise exterior point, i.e. noise spot；
Step 8, all noise spots for not meeting step (7) are removed to get filtered point set is arrived.
Priority Applications (1)
| Application Number | Priority Date | Filing Date | Title |
|---|---|---|---|
| CN201910365632.0A CN110246092B (en) | 2019-05-02 | 2019-05-02 | Three-dimensional laser point cloud denoising method considering neighborhood point mean distance and slope |
Applications Claiming Priority (1)
| Application Number | Priority Date | Filing Date | Title |
|---|---|---|---|
| CN201910365632.0A CN110246092B (en) | 2019-05-02 | 2019-05-02 | Three-dimensional laser point cloud denoising method considering neighborhood point mean distance and slope |
Publications (2)
| Publication Number | Publication Date |
|---|---|
| CN110246092A true CN110246092A (en) | 2019-09-17 |
| CN110246092B CN110246092B (en) | 2022-09-02 |
Family
ID=67883609
Family Applications (1)
| Application Number | Title | Priority Date | Filing Date |
|---|---|---|---|
| CN201910365632.0A Active CN110246092B (en) | 2019-05-02 | 2019-05-02 | Three-dimensional laser point cloud denoising method considering neighborhood point mean distance and slope |
Country Status (1)
| Country | Link |
|---|---|
| CN (1) | CN110246092B (en) |
Cited By (10)
| Publication number | Priority date | Publication date | Assignee | Title |
|---|---|---|---|---|
| CN111369610A (en) * | 2020-03-05 | 2020-07-03 | 山东交通学院 | Point cloud data gross error positioning and eliminating method based on credibility information |
| CN111524084A (en) * | 2020-05-13 | 2020-08-11 | 曹彬才 | Complex scene photon counting laser radar point cloud denoising algorithm based on multimodal Gaussian fitting |
| CN111612806A (en) * | 2020-01-10 | 2020-09-01 | 江西理工大学 | Building facade window extraction method and device |
| CN112986964A (en) * | 2021-02-26 | 2021-06-18 | 北京空间机电研究所 | Photon counting laser point cloud self-adaptive denoising method based on noise neighborhood density |
| CN113031524A (en) * | 2021-02-07 | 2021-06-25 | 南京航空航天大学 | Cubic spline-based press fitting force envelope generation method |
| CN113466826A (en) * | 2021-05-12 | 2021-10-01 | 武汉中仪物联技术股份有限公司 | Data denoising method, device, equipment and medium for range radar |
| CN114932562A (en) * | 2021-08-10 | 2022-08-23 | 南京航空航天大学 | Underground cable tunnel inspection robot based on laser radar and implementation method |
| CN116628429A (en) * | 2023-07-26 | 2023-08-22 | 青岛远度智能科技有限公司 | Intelligent control method for stable lifting of unmanned aerial vehicle |
| CN117828737A (en) * | 2024-01-04 | 2024-04-05 | 济南市园林绿化工程质量与安全中心 | Digital twin landscape design method |
| CN117974933A (en) * | 2024-03-28 | 2024-05-03 | 岐山县华强工贸有限责任公司 | 3D printing mold rapid scanning method for disc brake calipers |
Citations (6)
| Publication number | Priority date | Publication date | Assignee | Title |
|---|---|---|---|---|
| CN102629367A (en) * | 2012-01-17 | 2012-08-08 | 安徽建筑工业学院 | Bilateral filtering de-noising method of point cloud data based on KDTree |
| CN104573705A (en) * | 2014-10-13 | 2015-04-29 | 北京建筑大学 | Clustering method for building laser scan point cloud data |
| CN106918813A (en) * | 2017-03-08 | 2017-07-04 | 浙江大学 | A kind of three-dimensional sonar point cloud chart image intensifying method based on distance statistics |
| US20170193692A1 (en) * | 2015-12-30 | 2017-07-06 | Shenzhen Institutes Of Advanced Technology Chinese Academy Of Sciences | Three-dimensional point cloud model reconstruction method, computer readable storage medium and device |
| CN107169464A (en) * | 2017-05-25 | 2017-09-15 | 中国农业科学院农业资源与农业区划研究所 | A kind of Method for Road Boundary Detection based on laser point cloud |
| CN107392875A (en) * | 2017-08-01 | 2017-11-24 | 长安大学 | A kind of cloud data denoising method based on the division of k neighbours domain |
-
2019
- 2019-05-02 CN CN201910365632.0A patent/CN110246092B/en active Active
Patent Citations (6)
| Publication number | Priority date | Publication date | Assignee | Title |
|---|---|---|---|---|
| CN102629367A (en) * | 2012-01-17 | 2012-08-08 | 安徽建筑工业学院 | Bilateral filtering de-noising method of point cloud data based on KDTree |
| CN104573705A (en) * | 2014-10-13 | 2015-04-29 | 北京建筑大学 | Clustering method for building laser scan point cloud data |
| US20170193692A1 (en) * | 2015-12-30 | 2017-07-06 | Shenzhen Institutes Of Advanced Technology Chinese Academy Of Sciences | Three-dimensional point cloud model reconstruction method, computer readable storage medium and device |
| CN106918813A (en) * | 2017-03-08 | 2017-07-04 | 浙江大学 | A kind of three-dimensional sonar point cloud chart image intensifying method based on distance statistics |
| CN107169464A (en) * | 2017-05-25 | 2017-09-15 | 中国农业科学院农业资源与农业区划研究所 | A kind of Method for Road Boundary Detection based on laser point cloud |
| CN107392875A (en) * | 2017-08-01 | 2017-11-24 | 长安大学 | A kind of cloud data denoising method based on the division of k neighbours domain |
Non-Patent Citations (3)
| Title |
|---|
| CHINTHAKA DINESH 等: "local 3D point cloud denoising via bipartite graph approximation & totalvariation", 《2018 IEEE 20TH INTERNATIONAL WORKSHOP ON MULTIMEDIA SIGNAL PROCESSING》 * |
| 郭子选等: "基于频率直方图的K邻域稀疏离群点移除算法", 《计算机应用与软件》 * |
| 郭进等: "采用密度k-means和改进双边滤波的点云自适应去噪算法", 《传感器与微系统》 * |
Cited By (16)
| Publication number | Priority date | Publication date | Assignee | Title |
|---|---|---|---|---|
| CN111612806A (en) * | 2020-01-10 | 2020-09-01 | 江西理工大学 | Building facade window extraction method and device |
| CN111369610B (en) * | 2020-03-05 | 2022-09-06 | 山东交通学院 | Point cloud data gross error positioning and eliminating method based on credibility information |
| CN111369610A (en) * | 2020-03-05 | 2020-07-03 | 山东交通学院 | Point cloud data gross error positioning and eliminating method based on credibility information |
| CN111524084A (en) * | 2020-05-13 | 2020-08-11 | 曹彬才 | Complex scene photon counting laser radar point cloud denoising algorithm based on multimodal Gaussian fitting |
| CN111524084B (en) * | 2020-05-13 | 2023-08-22 | 曹彬才 | Photon counting laser radar point cloud denoising method based on multi-peak Gaussian fitting |
| CN113031524A (en) * | 2021-02-07 | 2021-06-25 | 南京航空航天大学 | Cubic spline-based press fitting force envelope generation method |
| CN112986964B (en) * | 2021-02-26 | 2023-03-31 | 北京空间机电研究所 | Photon counting laser point cloud self-adaptive denoising method based on noise neighborhood density |
| CN112986964A (en) * | 2021-02-26 | 2021-06-18 | 北京空间机电研究所 | Photon counting laser point cloud self-adaptive denoising method based on noise neighborhood density |
| CN113466826A (en) * | 2021-05-12 | 2021-10-01 | 武汉中仪物联技术股份有限公司 | Data denoising method, device, equipment and medium for range radar |
| CN114932562A (en) * | 2021-08-10 | 2022-08-23 | 南京航空航天大学 | Underground cable tunnel inspection robot based on laser radar and implementation method |
| CN114932562B (en) * | 2021-08-10 | 2024-04-19 | 南京航空航天大学 | Underground cable tunnel inspection robot based on laser radar and implementation method |
| CN116628429A (en) * | 2023-07-26 | 2023-08-22 | 青岛远度智能科技有限公司 | Intelligent control method for stable lifting of unmanned aerial vehicle |
| CN116628429B (en) * | 2023-07-26 | 2023-10-10 | 青岛远度智能科技有限公司 | Intelligent control method for stable lifting of unmanned aerial vehicle |
| CN117828737A (en) * | 2024-01-04 | 2024-04-05 | 济南市园林绿化工程质量与安全中心 | Digital twin landscape design method |
| CN117974933A (en) * | 2024-03-28 | 2024-05-03 | 岐山县华强工贸有限责任公司 | 3D printing mold rapid scanning method for disc brake calipers |
| CN117974933B (en) * | 2024-03-28 | 2024-06-11 | 岐山县华强工贸有限责任公司 | 3D printing mold rapid scanning method for disc brake calipers |
Also Published As
| Publication number | Publication date |
|---|---|
| CN110246092B (en) | 2022-09-02 |
Similar Documents
| Publication | Publication Date | Title |
|---|---|---|
| CN110246092A (en) | A kind of three-dimensional laser point cloud denoising method for taking neighborhood point mean distance and slope into account | |
| CN110363861B (en) | Laser radar point cloud-based field crop three-dimensional reconstruction method | |
| CN111932688B (en) | Indoor plane element extraction method, system and equipment based on three-dimensional point cloud | |
| Garland et al. | Fast polygonal approximation of terrains and height fields | |
| CN111986115A (en) | Accurate elimination method for laser point cloud noise and redundant data | |
| CN105844602B (en) | A kind of airborne LIDAR point cloud three-dimensional filtering method based on volume elements | |
| CN111598780B (en) | Terrain adaptive interpolation filtering method suitable for airborne LiDAR point cloud | |
| CN115222625A (en) | Laser radar point cloud denoising method based on multi-scale noise | |
| EP2887262B1 (en) | Point Cloud Simplification | |
| Pound et al. | A patch-based approach to 3D plant shoot phenotyping | |
| CN112462347B (en) | Laser radar point cloud rapid classification filtering algorithm based on density clustering | |
| CN111340723A (en) | Terrain-adaptive airborne LiDAR point cloud regularization thin plate spline interpolation filtering method | |
| CN115249330B (en) | Method and system for evaluating vegetation connectivity | |
| Zhu et al. | 3D reconstruction of plant leaves for high-throughput phenotyping | |
| CN116580048B (en) | Method and system for extracting contour line of right-angle house on house inclination model | |
| CN108986024A (en) | A kind of regularly arranged processing method of laser point cloud based on grid | |
| CN114898043A (en) | Laser point cloud data tile construction method | |
| CN108154516A (en) | Towards the point cloud topology dividing method and device of enclosure space | |
| CN114359437A (en) | Building structure two-dimensional plane map reconstruction method based on point cloud | |
| CN116523898A (en) | Tobacco phenotype character extraction method based on three-dimensional point cloud | |
| CN113744408A (en) | Grid generation method, device and storage medium | |
| Zhou et al. | Superpixels by bilateral geodesic distance | |
| Lei et al. | Extraction of maize leaf base and inclination angles using terrestrial laser scanning (TLS) data | |
| CN113129328B (en) | Target hotspot area fine analysis method | |
| CN113139982B (en) | Automatic segmentation method for indoor room point cloud |
Legal Events
| Date | Code | Title | Description |
|---|---|---|---|
| PB01 | Publication | ||
| PB01 | Publication | ||
| SE01 | Entry into force of request for substantive examination | ||
| SE01 | Entry into force of request for substantive examination | ||
| GR01 | Patent grant | ||
| GR01 | Patent grant |