CN107092020B - Merge the surface evenness monitoring method of unmanned plane LiDAR and high score image - Google Patents
Merge the surface evenness monitoring method of unmanned plane LiDAR and high score image Download PDFInfo
- Publication number
- CN107092020B CN107092020B CN201710257172.0A CN201710257172A CN107092020B CN 107092020 B CN107092020 B CN 107092020B CN 201710257172 A CN201710257172 A CN 201710257172A CN 107092020 B CN107092020 B CN 107092020B
- Authority
- CN
- China
- Prior art keywords
- point cloud
- flatness
- data
- road
- road surface
- 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.)
- Active
Links
- 238000000034 method Methods 0.000 title claims abstract description 127
- 238000012544 monitoring process Methods 0.000 title claims abstract description 60
- 230000003595 spectral effect Effects 0.000 claims abstract description 35
- 230000008569 process Effects 0.000 claims abstract description 31
- 238000001914 filtration Methods 0.000 claims abstract description 16
- 238000011156 evaluation Methods 0.000 claims description 46
- 238000004364 calculation method Methods 0.000 claims description 22
- 238000004422 calculation algorithm Methods 0.000 claims description 13
- PXFBZOLANLWPMH-UHFFFAOYSA-N 16-Epiaffinine Natural products C1C(C2=CC=CC=C2N2)=C2C(=O)CC2C(=CC)CN(C)C1C2CO PXFBZOLANLWPMH-UHFFFAOYSA-N 0.000 claims description 11
- 238000007637 random forest analysis Methods 0.000 claims description 11
- 238000013441 quality evaluation Methods 0.000 claims description 10
- 230000009466 transformation Effects 0.000 claims description 9
- 230000002159 abnormal effect Effects 0.000 claims description 8
- 201000010099 disease Diseases 0.000 claims description 8
- 208000037265 diseases, disorders, signs and symptoms Diseases 0.000 claims description 8
- 238000007781 pre-processing Methods 0.000 claims description 7
- 238000012732 spatial analysis Methods 0.000 claims description 7
- 230000036541 health Effects 0.000 claims description 6
- 230000009467 reduction Effects 0.000 claims description 6
- 238000012937 correction Methods 0.000 claims description 5
- 238000010801 machine learning Methods 0.000 claims description 5
- 238000005070 sampling Methods 0.000 claims description 5
- 230000003746 surface roughness Effects 0.000 claims description 5
- 238000000605 extraction Methods 0.000 claims description 4
- 238000012854 evaluation process Methods 0.000 claims description 2
- 238000007689 inspection Methods 0.000 claims description 2
- 238000002156 mixing Methods 0.000 claims description 2
- 238000012795 verification Methods 0.000 claims description 2
- NAWXUBYGYWOOIX-SFHVURJKSA-N (2s)-2-[[4-[2-(2,4-diaminoquinazolin-6-yl)ethyl]benzoyl]amino]-4-methylidenepentanedioic acid Chemical compound C1=CC2=NC(N)=NC(N)=C2C=C1CCC1=CC=C(C(=O)N[C@@H](CC(=C)C(O)=O)C(O)=O)C=C1 NAWXUBYGYWOOIX-SFHVURJKSA-N 0.000 claims 2
- 238000004458 analytical method Methods 0.000 abstract description 6
- 230000007547 defect Effects 0.000 abstract description 4
- 238000012549 training Methods 0.000 description 18
- 238000012545 processing Methods 0.000 description 17
- 238000005516 engineering process Methods 0.000 description 12
- 239000011159 matrix material Substances 0.000 description 11
- 238000010586 diagram Methods 0.000 description 10
- 238000012423 maintenance Methods 0.000 description 10
- 238000010276 construction Methods 0.000 description 9
- 230000000694 effects Effects 0.000 description 9
- 239000004575 stone Substances 0.000 description 9
- 238000005259 measurement Methods 0.000 description 8
- 238000006073 displacement reaction Methods 0.000 description 7
- 230000004927 fusion Effects 0.000 description 7
- 238000007726 management method Methods 0.000 description 6
- 239000004576 sand Substances 0.000 description 6
- 238000011161 development Methods 0.000 description 5
- 230000008901 benefit Effects 0.000 description 4
- 238000002474 experimental method Methods 0.000 description 4
- 230000000007 visual effect Effects 0.000 description 4
- ATJFFYVFTNAWJD-UHFFFAOYSA-N Tin Chemical compound [Sn] ATJFFYVFTNAWJD-UHFFFAOYSA-N 0.000 description 3
- 238000013016 damping Methods 0.000 description 3
- 238000001514 detection method Methods 0.000 description 3
- 230000001788 irregular Effects 0.000 description 3
- 238000013178 mathematical model Methods 0.000 description 3
- 239000013598 vector Substances 0.000 description 3
- 238000012800 visualization Methods 0.000 description 3
- 239000010426 asphalt Substances 0.000 description 2
- 238000006243 chemical reaction Methods 0.000 description 2
- 230000008878 coupling Effects 0.000 description 2
- 238000010168 coupling process Methods 0.000 description 2
- 238000005859 coupling reaction Methods 0.000 description 2
- 230000001186 cumulative effect Effects 0.000 description 2
- 238000007499 fusion processing Methods 0.000 description 2
- 230000002452 interceptive effect Effects 0.000 description 2
- 238000013507 mapping Methods 0.000 description 2
- 230000008439 repair process Effects 0.000 description 2
- 238000011160 research Methods 0.000 description 2
- 239000002689 soil Substances 0.000 description 2
- 238000007619 statistical method Methods 0.000 description 2
- 239000000725 suspension Substances 0.000 description 2
- 235000019733 Fish meal Nutrition 0.000 description 1
- 238000000342 Monte Carlo simulation Methods 0.000 description 1
- 238000012952 Resampling Methods 0.000 description 1
- 238000009825 accumulation Methods 0.000 description 1
- 230000009471 action Effects 0.000 description 1
- 238000013459 approach Methods 0.000 description 1
- 230000009286 beneficial effect Effects 0.000 description 1
- 230000008859 change Effects 0.000 description 1
- 238000007635 classification algorithm Methods 0.000 description 1
- 238000013145 classification model Methods 0.000 description 1
- 230000001419 dependent effect Effects 0.000 description 1
- 238000013461 design Methods 0.000 description 1
- 239000006185 dispersion Substances 0.000 description 1
- 230000005284 excitation Effects 0.000 description 1
- 238000013213 extrapolation Methods 0.000 description 1
- 239000004467 fishmeal Substances 0.000 description 1
- 230000010354 integration Effects 0.000 description 1
- 238000012986 modification Methods 0.000 description 1
- 230000004048 modification Effects 0.000 description 1
- 238000012805 post-processing Methods 0.000 description 1
- 238000002360 preparation method Methods 0.000 description 1
- 230000001737 promoting effect Effects 0.000 description 1
- 102000004169 proteins and genes Human genes 0.000 description 1
- 108090000623 proteins and genes Proteins 0.000 description 1
- 238000004445 quantitative analysis Methods 0.000 description 1
- 230000007115 recruitment Effects 0.000 description 1
- 230000000717 retained effect Effects 0.000 description 1
- 238000013077 scoring method Methods 0.000 description 1
- 238000012216 screening Methods 0.000 description 1
- 238000001228 spectrum Methods 0.000 description 1
- 238000006467 substitution reaction Methods 0.000 description 1
- 238000012360 testing method Methods 0.000 description 1
- 238000007794 visualization technique Methods 0.000 description 1
Classifications
-
- G—PHYSICS
- G01—MEASURING; TESTING
- G01S—RADIO DIRECTION-FINDING; RADIO NAVIGATION; DETERMINING DISTANCE OR VELOCITY BY USE OF RADIO WAVES; LOCATING OR PRESENCE-DETECTING BY USE OF THE REFLECTION OR RERADIATION OF RADIO WAVES; ANALOGOUS ARRANGEMENTS USING OTHER WAVES
- G01S17/00—Systems using the reflection or reradiation of electromagnetic waves other than radio waves, e.g. lidar systems
- G01S17/88—Lidar systems specially adapted for specific applications
-
- G—PHYSICS
- G01—MEASURING; TESTING
- G01C—MEASURING DISTANCES, LEVELS OR BEARINGS; SURVEYING; NAVIGATION; GYROSCOPIC INSTRUMENTS; PHOTOGRAMMETRY OR VIDEOGRAMMETRY
- G01C11/00—Photogrammetry or videogrammetry, e.g. stereogrammetry; Photographic surveying
-
- G—PHYSICS
- G01—MEASURING; TESTING
- G01C—MEASURING DISTANCES, LEVELS OR BEARINGS; SURVEYING; NAVIGATION; GYROSCOPIC INSTRUMENTS; PHOTOGRAMMETRY OR VIDEOGRAMMETRY
- G01C11/00—Photogrammetry or videogrammetry, e.g. stereogrammetry; Photographic surveying
- G01C11/04—Interpretation of pictures
-
- G—PHYSICS
- G01—MEASURING; TESTING
- G01S—RADIO DIRECTION-FINDING; RADIO NAVIGATION; DETERMINING DISTANCE OR VELOCITY BY USE OF RADIO WAVES; LOCATING OR PRESENCE-DETECTING BY USE OF THE REFLECTION OR RERADIATION OF RADIO WAVES; ANALOGOUS ARRANGEMENTS USING OTHER WAVES
- G01S17/00—Systems using the reflection or reradiation of electromagnetic waves other than radio waves, e.g. lidar systems
- G01S17/86—Combinations of lidar systems with systems other than lidar, radar or sonar, e.g. with direction finders
Landscapes
- Engineering & Computer Science (AREA)
- Physics & Mathematics (AREA)
- Radar, Positioning & Navigation (AREA)
- Remote Sensing (AREA)
- General Physics & Mathematics (AREA)
- Computer Networks & Wireless Communication (AREA)
- Electromagnetism (AREA)
- Multimedia (AREA)
- Length Measuring Devices By Optical Means (AREA)
- Image Processing (AREA)
Abstract
The invention discloses a kind of highroad pavement planeness monitoring methods for merging unmanned plane LiDAR and high-resolution remote sensing image, comprising: the LiDAR point cloud and high score image data obtain to unmanned plane pre-processes；It carries out the multiple dimensioned geometrical characteristic of the spectral signature of high score image and point cloud data to construct validity feature variables set, to realize the Optimum Classification of the road surface road surface point Yun Yufei point cloud；Further more accurate road surface point cloud data is obtained using Filtering Analysis；Construct fine road surface model；Surface evenness is calculated using IRI index.It is low that the present invention effectively overcomes existing surface evenness monitoring method inefficiency, the degree of automation, is unsuitable for a wide range of system-wide section monitoring, the defects of can only obtaining the road vertical section information of only a few, can not reflect road surfaces three-dimensional structure information comprehensively；It is practical, it can be widely applied to urban road and the flatness monitoring of other inferior grade roads.
Description
Technical Field
The invention relates to a remote sensing monitoring technology for the road surface evenness, which utilizes an unmanned aerial vehicle laser radar (LiDAR) And a high-resolution remote sensing technology to monitor And evaluate the road surface evenness in a large range, high efficiency And low cost aiming at the road surface.
Background
An Airborne laser radar (Airborne LiDAR) is an active earth observation system, integrates a global navigation satellite system, inertial navigation, laser ranging and computer processing technologies, can directly acquire a digital ground model with high precision and high resolution and three-dimensional space information of ground objects, has the superiority that a traditional photogrammetric method cannot replace, and is widely applied to the fields of topographic mapping, road survey design, vegetation parameter inversion, disaster monitoring, urban planning, three-dimensional modeling and the like. The highway is used as a modern transportation channel to play an increasingly important role in social economy, and plays a positive role in promoting logistics, resource development, recruitment and quotation, industrial structure adjustment, transverse economy union and the like along the line, so that the quality of the highway directly influences the economic development of regions and countries. In document [1] (Sun L, Zhang Z M, RuthJ. modeling index standards of surface roads [ J ]. Journal of transportation Engineering-ASCE,2001,127(2):105-111), road tests conducted by the United states Joint road and transportation administration in 1960 showed that: approximately 95% of the road service performance is dependent on the flatness of the road surface. The pavement evenness is used as one of main indexes for evaluating the pavement quality and is also an important decision basis for executing the pavement maintenance management of the highway. According to statistics of a statistical bulletin for development of transportation industry in 2015, the total highway mileage in 2015 reaches 457.73 kilometers, the total highway maintenance mileage reaches 446.56 kilometers, and the total highway mileage accounts for 97.6 percent, and the trend of rapid increase year by year is shown. With the rapid development of highway construction in China, how to realize the monitoring of the road pavement health condition with large range, high efficiency and low cost becomes a problem to be solved urgently.
Flatness monitoring is one of the important indexes for monitoring and evaluating the health condition of the road surface. The flatness evaluation Index is measured by a 3-meter ruler measurement maximum gap (h), an International flatness Index (IRI), a Running Quality Index (RQI), a bump cumulative Value (VBI), a flatness standard deviation (σ), a Power Spectral Density (PSD), and the like. The traditional road surface flatness monitoring method comprises the following steps: the characteristics of the 3m ruler, the continuous flatness meter, the hand-push type section meter, the vehicle-mounted jolt accumulation meter and the recent laser flatness detector (including a road detection vehicle integrated with the laser flatness meter) are summarized in table 1. One of the commonalities of these methods is that only a few longitudinal section information of roads can be obtained by ground measurement, the description of the road flatness is rough, and the three-dimensional geometric information of the road surface cannot be obtained quickly, so that the management department has very limited knowledge of the overall road condition.
TABLE 1. common road flatness monitoring method
The appearance of the low-altitude unmanned LiDAR technology not only provides possibility for acquiring three-dimensional information of a road surface, but also has the advantages of less terrain limitation, high precision, high efficiency, low cost (especially large-range monitoring), and capability of acquiring ground remote sensing images simultaneously, and provides a new effective solution for monitoring the flatness of the road surface. At present, the LiDAR technology is utilized to monitor the road surface flatness, and related cases and researches are rare. Document [2] (Chinese A, Olsen M J. Evaluation of technologies for Road Profile Capture, Analysis, and Evaluation [ J ]. Journal of surveying Engineering,2014,141(1): 4014011) describes that highway flatness monitoring is performed by using ground lidar technology, but this method does not attempt a spatial visualization method of Road flatness, nor does it explore a classification method of point clouds and a construction method of a Road surface model, and another fatal point of this technology is that point cloud data of a small range (radius <50M) can only be acquired near a certain fixed point position. If a large range of pavement information is to be acquired, a series of complex post-processing processes with continuously accumulated errors, such as point cloud splicing, coordinate correction, point cloud registration, coordinate conversion and the like, are required, so that the process is complex and very inefficient. Document [3] (Guojiao road flatness detection system research [ D ]. capital university, 2013) based on vehicle-mounted laser records that the vehicle-mounted laser radar technology is utilized to analyze the road section flatness characteristics and locate the position where the road damage exceeds the standard, but the method does not explore the point cloud classification method, and completes classification work by means of the traditional point cloud processing platform. In addition, the method only depends on a single numerical value for expression, cannot search the spatial distribution characteristics of the road surface quality, and requires stable running speed of a detection vehicle and stable GPS signals, which are difficult to apply to the complicated road sections with certain building shelters and bumpiness.
In summary, most of the conventional highway pavement evenness monitoring methods have low efficiency and low automation degree, cannot perform wide-range full-section monitoring, and can only obtain a small number of pieces of road longitudinal section information and cannot obtain three-dimensional structure information of the road surface. The flatness monitoring method of the existing laser radar technology has the defects of small measuring range, low working efficiency, visual expression of lack space and single data source, and is limited by terrain conditions. However, the monitoring of the road surface evenness needs to explore a new technical approach urgently, and the monitoring efficiency and accuracy of the road surface evenness in a large range are improved.
Disclosure of Invention
Aiming at the defects of the existing road flatness monitoring method, the invention provides a road flatness monitoring method integrating unmanned aerial vehicle LiDAR and high-resolution images and an implementation scheme thereof. Firstly, point cloud data and high-resolution image data acquired by an unmanned aerial vehicle LiDAR system are preprocessed, then spectral features are extracted from the high-resolution image data, and multi-scale geometric features are extracted from the point cloud data. Taking the obtained geometric features and spectral features as feature variables of point cloud classification, performing redundancy removal and dimensionality reduction on the feature variables, fusing the two types of data at a feature level, and performing point cloud classification by using a random forest classification algorithm in machine learning; and further obtaining a more accurate road point cloud data set through filtering operation and gross error correction. And then constructing a high-precision fine road surface Digital Surface Model (DSM) by using a three-dimensional interpolation method. And calculating the road surface evenness based on an IRI index principle, and performing space visual expression and evenness quality evaluation by using a GIS (geographic information system). Compared with the prior art that the whole section information can be expressed by only using a single IRI value, the method can provide the spatial distribution of the IRI values of different lanes of different sections of a certain road, improves the scientificity of the road maintenance management work decision, effectively improves the practicability of the unmanned aerial vehicle remote sensing technology in the road maintenance management application, and can be widely applied to the road surface flatness monitoring requirements of urban roads and other low-grade roads.
The principle of the invention is as follows: the point cloud data acquired by the unmanned aerial vehicle LiDAR remote sensing system comprises three-dimensional coordinate information of roads and surrounding ground objects, noise interference can be effectively removed by selecting a proper point cloud denoising means, and high-resolution image data acquired synchronously is rich in spectral information of the ground objects and is very easy to visually interpret. The point cloud data is not easy to visually interpret because of point set dispersion, so that spectral information and space geometric information can be integrated simultaneously by fusing the high-resolution image and the point cloud data, and more available information is provided for point cloud classification. At present, a challenge exists in classifying the ground objects of high-resolution images, namely that the ground object targets of foreign objects and spectra cannot be distinguished, and the space geometric information of point cloud data can be well supplemented. Similarly, the point cloud data itself is very discrete, which causes confusion at some edge boundaries and between ground objects with very similar geometric structural features, and the high-resolution image can assist this. Based on the method, on one hand, the spectral features of typical ground object types are extracted from high-resolution images, on the other hand, the geometric features of the typical ground object types are extracted from point clouds, the spatial form difference of the ground object types is strengthened by constructing multi-scale geometric features, then the point clouds are classified by adopting a high-efficiency and anti-noise strong machine learning classification algorithm, then rough difference correction and class filtering are carried out on the classified point clouds to obtain more accurate road point cloud data, and then the fine road DSM is constructed by combining a point cloud interpolation algorithm. On the basis, a series of longitudinal section curves are selected from the road surface, the international flatness value is calculated in a segmented mode and is expressed in a space visualization mode, and finally the flatness value is graded in an equal mode, so that the evaluation on the road flatness quality is achieved. The method can be used for carrying out rapid automatic flatness quality condition monitoring on the road surface of a large-range road through precision evaluation.
The technical scheme provided by the invention is as follows:
a road flatness monitoring method fusing unmanned aerial vehicle LiDAR and high-resolution images comprises the following steps: the method comprises the technical links of data acquisition and pretreatment, feature extraction and fusion, point cloud classification, pavement DSM model construction, pavement evenness evaluation and the like.
The data acquisition and preprocessing aim to de-noise point cloud data obtained by a low-altitude unmanned LiDAR system, and the high-resolution images obtained synchronously are registered and spliced to obtain ortho-image data of an experimental area, so that preparation is made for data fusion. The purpose of data fusion is to integrate geometric spatial information of point cloud data and spectral information of high-resolution image data, so that the point cloud data has both the spatial geometric information and the spectral information. The point cloud data are classified by using a random forest algorithm to obtain the road point cloud data, the available information is provided for point cloud classification in the data fusion process, typical ground object classification is carried out on an experimental area, a series of geometric characteristic parameters are constructed to express the difference of ground object types, the characteristic variables of a classifier are formed by combining spectral characteristic parameters, and then the point cloud data are classified by using the random forest algorithm to obtain the road surface point cloud data. The construction of the road surface DSM model can provide a basis for the expression of a three-dimensional space structure of the road surface, and each grid unit of the DSM model stores elevation information, so that the flatness calculation is facilitated, and conditions are provided for the space visualization expression of the flatness. The purpose of the road surface evenness evaluation is to quantitatively describe the evenness condition of the road surface by selecting a certain evenness index, and establish a grading system to objectively evaluate the road condition of an experimental area, so that prospective analysis and decision are given to the current overall condition of the road surface.
The method for monitoring the road flatness provided by the invention specifically comprises the following steps:
1) with a low altitude unmanned LiDAR system (including: integration of unmanned aerial vehicle platform, LiDAR scanner, multispectral camera, attitude positioning unit, autopilot, etc.), obtaining LiDAR point cloud data and high resolution image data of the experimental area;
2) converting the acquired LiDAR point cloud data into a current universal point cloud data file format (. las), and then preprocessing the point cloud data and high-resolution image data shot by an unmanned aerial vehicle, including point cloud denoising, image splicing and registering processing and orthoimage data generation;
3) analyzing the wave bands of the ortho-images obtained in the step 2), extracting spectral characteristic parameters reflecting the difference of the ground features, and forming a plurality of spectral characteristic sets containing spectral information together with the original wave bands of the high-resolution images;
4) performing feature fusion on the obtained spectral feature data and point cloud data, and performing geographic space registration by taking the point cloud as a space reference;
5) dividing main ground feature types of the experimental area, performing thinning on the point cloud data obtained in the step 4) to obtain thinned point cloud data, selecting a series of geometric characteristic parameters for the thinned point cloud data to calculate, obtaining point cloud data containing multi-scale geometric characteristic parameters and spectral characteristic parameters by changing the scale of spatial analysis, selecting part of point cloud data of corresponding ground feature types from the point cloud data as training samples for subsequent point cloud classification, and taking the rest thinned point cloud data as data to be classified. The point cloud thinning process reduces the operation amount in the geometric characteristic parameter calculation process, simultaneously reduces the data amount required to be processed by classification, and greatly improves the data processing efficiency;
6) the characteristic variables obtained in the step 5) are usually many, and strong correlation among some variables inevitably occurs, so that the characteristic variables are redundant. In order to find key characteristic variables, so that a dimensionality reduction effect is achieved on data, the complexity of a classifier construction process is further simplified, and spectral information of a multispectral image and spatial geometric information of laser point cloud data are fused on a characteristic level to form a new point cloud classification characteristic set;
7) and (3) performing model training on the training sample in the step 5) by utilizing a random forest model based on the feature set in the step 6), obtaining a satisfactory classification effect by continuously adjusting model parameters, applying the model training to sample point cloud data to be classified for classification, and assigning the class attribute of the sample point cloud to the denoised point cloud data according to the spatial nearest neighbor interpolation principle so as to complete the classification of all the point clouds. The classification result may also include point clouds not belonging to the road surface, and at this time, the classified point cloud data is subjected to manual inspection and gross error correction, and further subjected to category attribute filtering to obtain road point cloud data;
8) simulating the fluctuation change of the road surface by adopting an irregular triangular network for the elevation information of the road point cloud data obtained in the step 7), and then constructing a DSM (digital surface model) of the high-precision fine road surface of the regular grid by a natural neighborhood interpolation method;
9) selecting a series of longitudinal section lines from the DSM model data obtained in the step 8), operating an IRI index model to calculate IRI values corresponding to the section lines, and performing rasterization processing to obtain a spatial distribution map of the IRI values of the road surface;
10) classifying the IRI value obtained in the step 9) according to a grading standard of road flatness, thereby obtaining a flatness quality evaluation result corresponding to the road section to be tested.
Aiming at the method for monitoring the road surface evenness, further, LiDAR point cloud data related in the step 2) is specifically point cloud data meeting the following 4 conditions:
a) and (3) coordinate system: the coordinate system of the point cloud data is WGS84 coordinate system (longitude, latitude) or local plane projection coordinate system (X, Y) matched with the local area;
b) attribute information: the point cloud data attribute information comprises an elevation, a scanning angle and a reflection intensity besides a longitude and a latitude;
c) point cloud density: the density of the point cloud is more than 400 points/square meter (namely the radius of the laser foot point is within 50 mm);
aiming at the flatness monitoring method, further, the data preprocessing process related in the step 2) mainly refers to the following 4 types: the method comprises the steps of filtering noise points of LiDAR point clouds, filtering elevation abnormal points of the LiDAR point clouds, filtering scanning angles of the LiDAR point clouds and splicing and registering high-resolution images. Aiming at the flatness monitoring method, further, the high-resolution image data related in the step 2) is a high-definition digital photo shot by a multispectral digital camera with good performance on the market, and no person can provide a matched posture and Positioning (POS) information file. Aiming at the flatness monitoring method, further, the spectral characteristic parameters related in the step 3) are in the following 3 types:
a) the gray values of the red, green and blue bands of the high-resolution image.
b) The standard deviation of the gray values of the red, green and blue wave bands of the high-resolution image;
c) the reflected intensity value of the laser point cloud (the laser pulses emitted by the LiDAR scanner are also one of the electromagnetic waves, as well as having wavelength information).
Aiming at the flatness monitoring method, further, the characteristics involved in the step 4) are fused, specifically, the point cloud data and the spectral characteristic parameters of the high-resolution image data are fused, and the fusion process is based on the spatial reference of the point cloud data, and an affine inverse transformation model is adopted to match the spectral value of each pixel on the image to a point cloud data coordinate system;
aiming at the flatness monitoring method, further, the geometric characteristic parameters related in the step 5) are in the following 3 types:
a) local Roughness (LDR). The specific meaning of the characteristic variable is the distance from a certain reference point in the point cloud to a best fit plane formed by neighborhood points under a certain spatial scale of the reference point;
b) local Dimension Feature (LDF). Document [4] (Brodu N, lagued.3d specific data classification Of complex natural science using an exemplary-scale dimensional criterion: Applications in geographic [ J ]. isprssjoural Of geographic and Remote Sensing, 2012,68: 121-;
c) local Height Difference (LHD). The specific meaning of the characteristic variable is that the difference between the maximum distance and the minimum distance of the best fitting surface formed by reaching a certain reference point in the point cloud and the neighborhood point of the reference point under a certain spatial scale is the local height difference, and the value of the plane with larger fluctuation is larger.
For the flatness monitoring method, further, the specific meaning of the multi-scale geometric features related in step 5) is that when a certain geometric feature value is calculated, the spatial analysis scale is multiple (a corresponding spatial geometric feature parameter can be obtained by fixing a certain spatial analysis scale), and the geometric feature values of different surface feature types can be different under different spatial analysis scales.
Aiming at the flatness monitoring method, further, the point cloud data types related in the step 5) mainly comprise 6 categories of road point cloud (including traffic markings) and low vegetation (crops and grasslands), soil, trees, moving targets and power lines (including power poles).
In view of the above flatness monitoring method, further, the point cloud thinning process involved in step 5) is to sample in space at certain spatial intervals, and this process is to ensure that the distribution of the point cloud in space is uniform (rather than random) with respect to the original point cloud.
For the flatness monitoring method, further, the feature selection involved in step 6) mainly includes two processes of searching for the feature subset and evaluating the feature subset. The search of the feature subset adopts a greedy strategy, namely: the initial optimal subset is empty, all feature variables to be evaluated form an alternative feature set, the feature variables with the strongest interpretation capability are screened from the alternative feature set and added into the optimal feature subset, and only one variable is added in one cycle; after the addition of a certain feature variable is finished, the feature variable needs to be removed from the alternative feature set, and the next cycle is repeated. Evaluation criteria of the optimal feature subset: the correlation between the newly added feature variables and the selected feature variables is low.
For the flatness monitoring method, the resolution of the DSM model involved in step 8) is further as high as possible, and the elevation accuracy should be controlled to within 15 mm. The specific horizontal and elevation accuracy evaluation criterion is as follows: selecting some targets (such as pit slots) with obvious characteristics from a road surface on site, recording corresponding GPS coordinates, measuring the length, width and height attributes of the targets, finding the corresponding targets in a high-resolution image, easily finding the corresponding positions of the targets in a DSM model by utilizing the fusion of the image and LiDAR point clouds, further interactively measuring the length, width and height of the targets by utilizing a length measuring tool, and then performing statistical analysis to finish the precision evaluation on the level and elevation of the DSM model.
For the flatness monitoring method, further, for the road profile line involved in the step 9), the smaller the distance is, the finer the corresponding spatial distribution of the flatness of the road surface is. The theoretical basis for the IRI index model is a quarter car model, requiring strict control of the sampling interval (within 50 mm) and the length of the measurement unit (10m or 20m) during the calculation. The optimal sampling interval and length of the measuring unit can be determined by quantitatively analyzing the variation rule of the IRI error with the sampling interval and length of the measuring unit when the elevation error is constant by a Monte Carlo simulation method.
Aiming at the flatness monitoring method, further, the rasterization involved in the step 9) refers to assigning the IRI value of the pavement longitudinal section line to the surrounding grid units according to the nearest principle, so as to obtain the flatness spatial distribution map of the road section to be measured.
For the flatness monitoring method, further, a series of road surface samples containing different flatness quality grades are selected as reference data on the road section to be detected, so that the evaluation result of the method can be evaluated accurately, and the method comprises the following specific steps:
1) according to a set grading standard of the flatness quality, a series of road surface samples (including good and poor flatness quality, covering health and disease types, and selecting a certain number of ground actual measurement samples for each road surface) are selected at a road section to be detected at the same time to serve as reference data;
2) selecting 5 indexes of appearance, surface roughness, track depth, disease area and damage degree of the reference data, and using a professional scoring method to evaluate the flatness quality of the reference data as reference values;
3) calculating the flatness value of the road surface actual measurement sample by using the method, and obtaining the evaluation result of the road surface flatness quality of the method according to the grading standard of the flatness quality;
4) the evaluation of the flatness quality result can be similar to the evaluation of the remote sensing image classification result, and the confusion matrix is the most common form for evaluating the remote sensing image classification precision, so the method also adopts the confusion matrix to evaluate the road flatness quality result. And introducing an Overall Accuracy (OA) of an Accuracy evaluation index and a Kappa coefficient to analyze the result to obtain the Accuracy of the method.
Compared with the prior art, the invention has the beneficial effects that:
the method for constructing the point cloud classification feature set by fusing the spectral features of the image and the multi-scale geometric features of the LiDAR point cloud is firstly provided based on LiDAR data and a high-resolution remote sensing image acquired by a low-altitude unmanned aerial vehicle platform, and the problem of current road surface point cloud and non-road surface point cloud classification is solved; and then, constructing a fine three-dimensional road surface Digital Surface Model (DSM) by using the road surface point cloud, and rapidly calculating to obtain the flatness parameters of the road surface by introducing the calculation principle of the international flatness index so as to realize rapid monitoring and evaluation of the road surface quality. The invention provides a rapid and efficient road surface flatness monitoring method integrating low-altitude unmanned LiDAR and high-resolution images for the first time, and effectively overcomes the defects that the conventional road surface flatness monitoring method is low in efficiency, low in automation degree, not suitable for large-scale all-section monitoring, only few road longitudinal section flatness information can be obtained, and road three-dimensional structure information cannot be obtained. The method can quickly master the pavement evenness condition in a large range, thereby providing scientific basis for pavement maintenance and repair decision. Meanwhile, the point cloud data processing flow designed by the invention is simple and easy to realize, and the effect is obvious. In addition, the flatness quality spatial distribution map provided by the invention can comprehensively reflect the flatness condition of the road surface, compared with the prior flatness monitoring technology which can express the flatness information of the whole road section only by depending on a single numerical value, the invention effectively promotes the informatization and refinement development of road maintenance management and improves the scientificity of the work decision of road maintenance management.
Compared with the prior art, the method has the advantages that firstly, a low-altitude unmanned aerial vehicle LiDAR data source is adopted, fused spectral characteristic information and multi-scale geometric characteristic information are provided for the first time to optimize road point cloud classification, fine DSM is constructed on the basis, and road surface flatness information with spatial distribution is extracted quickly; secondly, the data processing process is simple and easy to realize, the used mathematical model is easy to understand, the effect is obvious, all operations can be completed on a personal common computer, and the requirement on hardware is not high; thirdly, a method for optimizing data processing is further provided on the aspect of providing feasible operation flows, verification is completed through experiments, and good repeatability is achieved.
Drawings
FIG. 1 is a flow chart of a method for monitoring road flatness provided by the present invention
Fig. 2 is a schematic diagram of a local roughness calculation method.
FIG. 3 is a simplified representation of a local dimensional feature;
wherein, (a) is a one-dimensional space distribution characteristic, (b) is a two-dimensional space distribution characteristic, and (c) is a three-dimensional space distribution characteristic.
Fig. 4 is a schematic diagram of local height difference calculation.
FIG. 5 is a schematic diagram of a natural neighborhood interpolation algorithm for performing point cloud interpolation in the method.
FIG. 6 is a schematic diagram of an international flatness index (IRI) algorithm used in the present method;
wherein Z issRepresenting the vertical displacement of the sprung mass, ZuRepresenting the vertical displacement of the unsprung mass, msRepresenting the magnitude of the sprung mass, muIndicating the magnitude of the unsprung mass, KsRepresenting the coefficient of stiffness of the spring connecting the sprung and unsprung masses, CsLinear damping coefficient, K, representing the coupling sprung and unsprung massessThe precision coefficient of a spring connected with the tire is shown, the damping effect produced by the tire can be simulated, y (x) represents the road surface elevation, and b represents the accommodating length (the part of the tire, which is in contact with the ground) of the tire to the ground.
FIG. 7 is a photograph of a pavement sample taken from a ground surface with different flatness characteristics according to an embodiment of the present invention
Wherein, (a1) the pit slot is characterized in that the stone is exposed and distributed in a massive form, the shape is approximately round or oval, and the flatness quality is evaluated as extremely poor (Failed); (a2) a pit characterized in that the stone is exposed but has a relatively small area, the shape is similar to a circle or an ellipse, and the flatness quality is evaluated as extremely poor (Failed); (b1) collapse, characterized by a strip-like distribution with a large depth (stone is not exposed), a large coverage area, and a very poor flatness quality (Failed) rating; (b2) collapse, characterized by a striped distribution with small depth (stone is not exposed), small coverage area, flatness quality rated as Poor (Poor) or extremely Poor (Failed), depending on the depth of the collapse and the coverage area; (c) cracks (including fine cracks and coarse cracks) are characterized by long and narrow distribution, the crack gap is within 10mm, and the flatness quality is evaluated as Fair (Fair) or Poor (Poor), which depends on the width and the area of the crack; (d) pitted, characterized by a very rough road surface or many fine pits, the flatness quality is rated as medium (Fair) or Poor (Poor), depending on the degree and area of surface damage; (e) covering with sand and stone, characterized in that the road surface is broken in block shape and covered with sand and stone, and the flatness quality is evaluated as Poor (Poor); (f) the old pavement with general flatness is characterized in that the texture of the pavement is rough but the whole surface is smooth, and the flatness quality is evaluated as medium (Fair); (g) the old pavement with Good flatness is characterized in that the pavement is bright in color, dense in texture, very smooth in surface and Good in flatness quality evaluation (Good).
FIG. 8 is LiDAR point cloud data acquired and format converted within an experimental area in an embodiment of the present invention.
Fig. 9 is a multi-spectral ortho image of an experimental region processed by image stitching and registration according to an embodiment of the present invention.
FIG. 10 is LiDAR point cloud data fused with image band information (blue band shown in the figure) in an embodiment of the present invention.
Fig. 11 is an attribute value distribution characteristic of 23 key feature variables and training samples under these feature variables, which are retained after feature selection is performed using 5 pieces of spectral information and 36 pieces of multi-scale geometric feature information in the embodiment of the present invention;
wherein, the reflection intensity, the R-G-B wave band gray standard deviation, the R-G-B wave band gray value and other 5 spectral characteristics of different ground objects are (a); (b) local roughness (LDR) and Local Dimensional Features (LDF) of different features; (c) local Dimensional Features (LDF) and local elevation differences (LHD) of different features; (d) local Height Difference (LHD) of different features. The length number in brackets represents the size of a space analysis scale, LDF _1 represents a one-dimensional local dimension characteristic, and LDF _1 represents a two-dimensional local dimension characteristic
FIG. 12 illustrates an embodiment of the present invention for generating a DSM expressed in TIN by interpolating a point cloud of a road surface.
FIG. 13 illustrates a road DSM obtained by a natural neighborhood interpolation algorithm in accordance with an embodiment of the present invention.
Fig. 14 is a spatial distribution diagram of the road flatness calculated by using the IRI index model according to the embodiment of the present invention.
Fig. 15 is a spatial distribution diagram of road flatness quality obtained by re-classifying the flatness values according to the embodiment of the present invention.
FIG. 16 is a statistical map of the spatial distribution of different pavement quality types in an embodiment of the present invention.
Fig. 17 is a schematic diagram of a confusion matrix for evaluating the flatness quality evaluation results established by the embodiment of the present invention.
Detailed Description
The following further describes the implementation of the present invention by way of specific examples, but not in any way limits the scope of the present invention.
Fig. 1 is a block flow diagram of a road flatness monitoring method fusing unmanned aerial vehicle LiDAR and high-resolution images, which includes the following steps:
step 1: the method comprises the steps that an unmanned aerial vehicle LiDAR system is used for simultaneously obtaining LiDAR point cloud data and high-resolution image data of a road section to be detected;
step 2: because point cloud data formats provided by different LiDAR system manufacturers are different (most manufacturers have a set of self-defined storage formats), in order to process the point cloud data on some general software, format conversion needs to be carried out on the LiDAR point cloud data, and the selection of a projection coordinate system and reserved field attributes are paid attention to during output;
and step 3: and (3) denoising the general point cloud file (. las) obtained in the step (2), and preliminarily removing interference noise information, thereby obtaining the dehumidified point cloud data. The specific denoising processing mode comprises the following steps:
a) and filtering noise points. The method mainly comprises the steps of eliminating abnormal points generated due to refraction of laser pulses or multipath effect in LiDAR point cloud;
b) and filtering elevation abnormal points. Mainly eliminates local elevation abnormal points (gross errors), and the specific judgment method comprises the following steps: setting a certain search radius threshold r, and if the difference between the elevation average value of a certain point and the elevation average value of a neighborhood point within the radius r of the point exceeds 3 times of standard deviation, determining the point to be an elevation abnormal point and removing the elevation abnormal point;
c) and (5) filtering the scanning angle. Because when the unmanned aerial vehicle track is planned, in order to comprehensively acquire the surface information of the road section, the flight route of the unmanned aerial vehicle is right above the road section, and therefore the scanning angle attribute of the point cloud is utilized to filter data, so that the data processing amount is reduced, and unnecessary data processing time is saved. Assuming that the horizontal distance from a certain point to the position right below the airplane is d, and the flying height of the airplane is h, the scanning angle θ of the point is defined as formula 1:
θ ═ arctan (d/h) (formula 1)
The maximum scanning angle theta of the road can be determined according to the boundary range of the roadmaxThus, setting the appropriate scan angle threshold separates the road regions. Considering that the flight path has certain deviation, the scanning angle threshold value needs to be set slightly larger as much as possible in the experiment;
and 4, step 4: carrying out image splicing and registration processing on high-resolution image data shot by an unmanned aerial vehicle in a road section to be detected by using unmanned aerial vehicle image splicing software (such as commercial software Pix4D Mapper developed by Pix4D company of Switzerland in the example), so as to obtain an ortho-remote sensing image of the road section to be detected;
and 5: the spectral feature extraction is carried out on the high-resolution image, and the spectral feature parameters used by the method comprise the following 4 parameters: the gray values of the red, green and blue bands and the standard deviations of the three bands, and the mathematical model of the standard deviations is formula 2:
wherein std. represents the standard deviation, gi(i ═ 1,2,3) are gray values in the red, green, and blue three bands.
Step 6: and performing feature fusion on the multi-band high-resolution image and the point cloud data, mapping pixel coordinates of the image to geographic coordinates of a pixel center point by affine transformation by taking the geographic space reference of the point cloud data as a reference, and matching corresponding spectral information on the image to the point cloud data. The mathematical model of affine transformation is formula 3:
however, in actual operation, the number of pixels of the image is much larger than the data amount of the point cloud, and the two are not completely overlapped, so in order to reduce the calculation amount, the process of converting geographic coordinates into image coordinates, namely, affine inverse transformation, is adopted by the invention and is expressed as formula 4:
in the formulas 3 and 4, G is an affine transformation matrix, ai,bi,ci(i is1, 2) is an affine transformation coefficient, mapX and mapY are X and Y coordinates of the point cloud, respectively, and row and col are line and column values of the image, respectively. The metadata of the GeoTiff or other image files with geographic coordinate information generally has affine coefficient information, and the affine coefficient can also be obtained by constructing an equation through 4 boundary point coordinates on the image and corresponding row and column number (pixel coordinates) information thereof to solve the affine coefficient.
And 7: and constructing an octree for storing the point cloud data, and rarefying the point cloud by taking the length of an octree grid unit as a minimum space interval to obtain rarefied point cloud data. And then, calculating geometric characteristic parameters of the rarefied point cloud data, wherein the geometric characteristic parameters used in the method comprise the following three parameters:
a) local Roughness (LDR), as shown in fig. 2, assuming that the whole three-dimensional point cloud point set is C ═ pi|pi＝(xi,yi,zi) I is1, 2, 3.., n }, the current calculation point is P ═ (x, y, z) ∈ C, and the neighborhood point set with radius R in the three-dimensional space is P ═ Pj|||pj-p||≤R,pjGiven that plane T, Ax + By + Cz + D ═ 0 is the best fit plane for the set of points Q ═ P ∪ P ≠ P, j ═ 1,2, 3.
A, B, C and D in the formula 5 respectively represent basic parameters of a plane equation, and x, y and z respectively represent x, y and z coordinates of the current search point. The optimal fitting plane T in the invention is a plane which enables the sum of squared errors S of z values to be minimum when all points in a point set are projected on an x-y plane, and the specific mathematical form is as shown in formula 6:
x in formula 6i,yi,ziX, Y and Z coordinates respectively representing the ith point, and the meanings of other parameters are the same as formula 5;
b) local Dimension Feature (LDF), as shown in fig. 3, the structures of the three-dimensional point cloud point set and the neighborhood point set are similar to Local roughness, and let the current search point be P ═ X, Y, Z ∈ C, and let the point set Q ═ P ∪ P ═ X Y Z ═ C]First, principal component transformation is performed on Q to obtain three principal component coefficients mu of a matrix Q1,μ2,μ3(μ1≥μ2≥μ3) The three principal component coefficients are further normalized by equation 7:
wherein λ is1、λ2、λ3Respectively corresponding to the degree of one-dimensional, two-dimensional and three-dimensional spatial distribution of the current search point in the neighborhood, namely one-dimensional, two-dimensional and three-dimensional characteristic values, due to lambda1+λ2+λ31, the local dimensional features of the current search point can be expressed only by the first two feature parameters, so that the geometric feature calculation efficiency of the point cloud data can be improved;
c) local Height Difference (LHD), as shown in fig. 4, wherein the definition of the best fit plane is consistent with that of the best fit plane in the Local roughness calculation, assuming that the current search point is P ═ x, y, z ∈ C, and the set of neighborhood points with radius R in the three-dimensional space is P ═ P { (P), the set of neighborhood points with radius R is assumed to be P ═ yj|||pj-p||≤R,pjNot equal to P, j 1,2,3,.., n }, and point set Q { P ∪ P } is a neighborhood set of points that contains the current search point, then the mathematical form of the LHD can be expressed as equation 8:
in formula 9Respectively representing the distances from the ith and the j neighborhood points to the optimal fitting plane, and the meanings of other parameters are the same as formula 5;
in the experimental process, in order to accelerate the speed of calculating the local height difference of the search point, 4 key parameters A, B, C and D of the plane T can be quickly obtained by carrying out SVD on Q, and then the point with the maximum distance from the plane T and the point with the minimum distance from the plane T are searched for distance difference calculation, so that the local height difference of the search point is obtained.
By continuously changing the size of the radius R of the search neighborhood and sequentially calculating the corresponding values of the geometric features (LDR, LDF and LHD) under different search space radii R, a group of column vectors of the three geometric features can be obtained, and the group of column vectors are multi-scale geometric feature parameters of the ground object target (the object where the search point is located).
And 8: after the high-resolution image and the point cloud data are fused, the point cloud data at the moment contain a large number of characteristic variables (geometric characteristic parameters and spectral characteristic parameters) to form a high-dimensional matrix (each point is a 3-dimensional vector), correlation possibly exists among the characteristic variables, information can generate certain redundancy, and the process of removing redundant information is the process of characteristic selection. The feature selection strategy used by the present invention is divided into two processes: feature subset search and feature subset evaluation. The strategy of feature subset search is forward search, and the evaluation of feature subsets mainly evaluates the degree of correlation inside the subsets by means of correlation coefficients. The specific implementation process is as follows:
1) firstly, establishing two data structures in a set form, and assuming that P and Q are respectively provided, P represents a candidate subset, Q represents an optimal subset, a variable of a category attribute is Y, and other attributes except the category attribute are X ═ X1,X2,X3,...,Xn}. Adding other characteristic variables except Y into the candidate subset, and setting a correlation threshold value inside the subset to be delta;
2) if it is notEntering step 4), otherwise, circularly reading each characteristic variable X in XiAnd calculating the correlation coefficient with Y according to equation 10, and finding the characteristic variable X in which the correlation coefficient with Y is the largestcAnd let δ be max { corr (X)i,Y)}(XiE.x), then go to 3);
in formula 10, XiY represents the ith characteristic and the classification attribute respectively,respectively represent the average of the ith feature, the classification attribute, andrespectively representing the ith characteristic and the category attribute value of the jth sample data.
3) The characteristic variable X to be selectedcPerforming relevance evaluation on each characteristic variable in Q if XcIf the correlation with any one characteristic variable in Q is less than delta, then X is addedcAdding the characteristic variable into Q, and deleting the characteristic variable from P, otherwise, directly deleting the characteristic variable from P, and continuing to the step 2);
4) and outputting Q, and finishing the feature subset search.
It should be noted that this step can be omitted if the number of feature variables is not large (e.g., within 30).
And step 9: after the point cloud data after rarefaction is obtained, a part of the point cloud data is selected as a training sample in an interactive mode, the training sample is required to cover all ground object types, the training sample cannot be too concentrated and is properly dispersed in a space where the point cloud is located, the number of the training samples of each ground object type is not too small, and otherwise, a classifier cannot construct a classification model. The random forest model has the advantages of strong adaptability to a data set (discrete type and numerical type data can be processed simultaneously), high training speed, strong noise resistance, difficulty in falling into overfitting, few model parameters, good classification effect and the like, and becomes one of the most popular algorithms in the field of machine learning.
The basic idea of the algorithm is as follows: assuming that the number of input sample data is N (each sample corresponds to one row), the number of feature variables is M (each feature variable corresponds to one column, that is, the size of original sample data is a matrix of N × M), initially, N samples are selected from the N samples to serve as a training set of a tree (each tree in a random forest is a binary tree), and a root node of the tree is the training set; randomly selecting M feature variables from the M features, and then selecting the feature variables from the M feature variables to split according to the principle of minimum purity, so that the tree can be divided into two nodes, wherein the two nodes respectively comprise a part of the training subset, and the splitting is performed recursively in such a way until the tree can not be split continuously (the maximum recursive depth of the binary tree is reached) or the M features are used completely. Repeating the process of just building trees results in a forest containing any trees (random forest), and the training samples of each tree are not identical, nor are the feature variables used. In the stage of prediction classification, when a sample to be classified is input, each tree in the forest can obtain a decision classification for the sample, the decision classification is counted to obtain the distribution frequency of the sample to be classified under each prediction classification, and the classification with the maximum frequency is output as the prediction classification, so that classification is finished.
The purity measurement index used in the invention is a Keyny coefficient, and a calculation method recorded in a random forest-based fish meal protein near infrared quantitative analysis [ J ] in agricultural machinery, 2015(05): 233-:
1) assuming that the training sample set contains P classes, the ratio of the number of samples in the ith class is PiThe Gini (P) value is:
dividing the training sample set into m subsets according to a certain characteristic variable, wherein the ith subsetThe number of the samples of the set is niIf the total number of samples in the training sample set is n, the Gini coefficient of the divided subset is Ginisplit(P) is:
and step 9: classifying the point cloud after rarefaction, filtering the point cloud according to the category attribute to obtain road surface point cloud data, and further obtaining DSM data of the road section to be measured through a related point cloud interpolation algorithm. In order to reduce the influence of elevation noise on the surface model effect and subsequent flatness calculation, the method firstly utilizes an irregular triangulation network (TIN) to generate three-dimensional terrain data, and then adopts a natural neighborhood interpolation method to construct a DSM model of the road. The natural neighborhood interpolation method described in document [6] (Sibson R.A briefdescription of natural neighbor interpolation [ J ]. interpolation MultivariateData,1981:21-36) obtains the elevation value of a point to be interpolated by calculation using formula 13:
w in formula 13iWeight coefficient, Z, representing the ith elevation point in the neighborhoodiThen, the elevation of the ith elevation point in the neighborhood is obtained, i represents the ith elevation point in the neighborhood, and Z is the elevation value of the point to be interpolated. The method determines the proximity relation by determining whether there is an intersection with the Thiessen (Voronoi) polygon of the point to be interpolated, and the weight of the proximity point is the ratio of the common area of the Thiessen polygon of the point to be interpolated and the Thiessen polygon of the proximity point (see FIG. 5). It should be noted that: all points are considered in the construction process of the Thiessen polygons of the points to be interpolated, and the points to be interpolated are excluded in the construction process of the Thiessen polygons of the adjacent points. The interpolation method is basically characterized in that the interpolation method has locality, only a subset of sample points around a point to be interpolated is used, and a weight coefficient is ensuredIt is verified that the result of the interpolation is within the range of the sample values used. The advantage is that no extrapolation tendency occurs and no peaks, valleys, ridges or valleys are generated which are not represented by the input sample point, and the resulting surface is smooth at all positions except the input sample position;
step 10: and selecting a series of longitudinal section lines with a certain distance (such as 0.2m) by using the obtained road DSM, and expressing the flatness value of the experimental road section by using an international flatness index IRI model. In the specific implementation of the invention, the road section flatness value is calculated by a method described in document [7] (Sayers M W, Gillespie T D, Queiroz C. the international road surface roughness experiment [ J ]. 1986). The core principle of the IRI index is a quarter-vehicle model (as shown in fig. 6), when a quarter-vehicle runs along the longitudinal section of a road at a certain speed, the system generates vibration under the excitation action of the road gradient, and the cumulative displacement value of the suspension system after the quarter-vehicle runs for a certain distance (for example, 1km) at a specified speed (for example, 80km/h) is calculated to be IRI, and the unit is m/km. In order to solve the relative displacement of the suspension system, a second-order vibration differential equation is established as shown in the formula 14-15:
wherein Y is the road elevation, ZsAnd ZuRespectively, representing sprung and unsprung mass displacements, mu, C, K1、K2Is a coefficient and is according to the document [7]]The following are described:
wherein m issRepresenting the magnitude of the sprung mass, muIndicating the magnitude of the unsprung mass, KsRepresenting the coefficient of stiffness of the spring connecting the sprung and unsprung masses, CsLinear damping coefficient, K, representing the coupling sprung and unsprung massessIndicating the coefficient of precision of the spring connecting the tire.
Construction state variablesThe original equation can be:
wherein:
equation 16 is a non-homogeneous linear differential equation, and the corresponding equation of state is:
wherein,
PR＝A-1(ST-I) B, (formula 19)
I in the formulae 18 and 19 is an identity matrix in combinationOfValue state and gradeSequences, which can be used to determine in turn any time by recursionAnd the following steps:to pairThe integral is carried out to obtain a key variable ZsAnd ZuIRI can be calculated from equation 20:
document [8](influencing factors of IRI evaluation index of Rong Jianfeng, Song hong Xun, Ma Rong Gui. road surface flatness [ J]The Chongqing university of traffic journal (Nature science edition), 2012,31(06): 1145-. In this example, to implement the calculation process of the IRI index, the DSM data needs to be sampled at equal intervals to obtain a road surface elevation point sequence Y ═ Y0,y1,...,ynAnd it is used as input data of the model. Then, state matrixes ST and PR are calculated according to the matrixes A and B, and then first-order difference operation is carried out on Y to obtainAnd let Z (0) be [ b 0]TAnd b is a speed value 11 meters from the starting point (if the link length is less than 11m, 0m is taken), so that the t is 0,1,2Value, and then pairAnd accumulating time intervals to obtain integral displacement of each time point, calculating the displacement of each time point according to a formula 20, taking absolute values, summing, and finally taking an average value to solve the IRI.
Step 11: the calculation result of step 6 is only to obtain the IRI value sequence of a single section line, and in order to simulate the flatness value of the whole road surface, the corresponding IRI value can be allocated to each grid unit on the road DSM according to the nearest neighbor principle, the error of the processing mode mainly depends on the space size of the longitudinal section line, and the smaller the space size, the smaller the error. The flatness spatial distribution of the whole road section to be measured can be obtained through the process.
Step 12: the obtained road surface IRI values were classified with reference to the following criteria (see table 2), i.e., the evaluation of the road surface flatness quality was completed. The grading setting standard needs to be set by combining with the relevant requirements of road surface maintenance quality of 'technical Specification for maintaining asphalt road surface of highway' in China, and the method is determined by combining with the specification of 'technical Specification for maintaining asphalt road surface of highway' (JTJ 073.2-2001). The corresponding distribution proportion can be calculated by counting the number of pixels of each quality grade of the road section to be detected, and a scientific decision can be made on what maintenance measures are needed for the road section by combining the relevant specifications of road maintenance and repair.
TABLE 2 grading Standard of the flatness quality of the road surface
In order to evaluate the reliability of the method, the following steps are adopted to evaluate the accuracy of the road surface flatness quality monitoring result of the method.
1) Selecting a series of road surface samples (including good and bad flatness quality, including health and disease types, and selecting a certain amount of ground actual measurement sample data for each road surface) as reference data on a road section to be detected according to a set flatness quality grading standard;
2) selecting 5 indexes of appearance, surface roughness, track depth, disease area and damage degree of the reference data, and using an expert evaluation method to evaluate the flatness quality of the reference data as a reference value;
3) calculating the flatness value of the pavement sample by using the method, and obtaining a corresponding evaluation result of the method according to the grading standard of the flatness quality;
4) the precision evaluation of the flatness quality result can be similar to that of the remote sensing image classification result, and the confusion matrix is the most common form for evaluating the remote sensing image classification precision. Therefore, the method also adopts the confusion matrix to evaluate the pavement evenness quality result. Analyzing the result by means of the overall precision (OA) of the classification precision evaluation index and the Kappa coefficient to obtain the precision of the method, wherein the specific precision evaluation index calculation method comprises the following steps:
in the formulas 21 and 22, N is the total number of samples, NiiNumber of samples representing the method in which the evaluation result and the reference evaluation result are both in the i-order, ni.Then the total number of samples with the reference evaluation result of i grade is represented, n.iThe total number of samples representing the evaluation result of the method is the i-grade.
The concrete implementation and the precision evaluation process of the method are described below by taking a county level highway in a certain city as an example.
(1) Data acquisition
Unmanned aerial vehicle LiDAR point cloud data
The county level highway is located near 44 degrees 24 '47' north latitude and 85 degrees 53 '47' east longitude, and the road width is about 8 m. The investigator performed experimental data acquisition directly above the experimental area on day 23/6/2016 with the aid of a laser scanner system (Rigel VUX-1LR) carried by Scout B1-100 unmanned helicopter manufactured by Switzerland, the basic parameter information of which is shown in Table 3 (source: http:// www.riegl.com/products/newriegl-VUX-1-series/newriegl-VUX-1LR /). In order to ensure that the density of the acquired point cloud data is high enough, the flying distance from the ground to the unmanned aerial vehicle in the experiment is set to be 30m, the flying speed is 5m/s, the scanning angle is 110 degrees, the scanning frequency is 550KHz, the density of the laser foot points obtained finally is 300 plus 600pts, the scanning width is 85.7m, and the data is acquired without (little) cloud on the sunny day.
TABLE 3 Key parameter information for laser scanner systems
Ground data acquisition
On the same day as LiDAR point cloud data is acquired, researchers collect feature point and road surface sample data of the ground features in the flight area. The feature point data of the ground features are mainly collected along the periphery of a road, geometric dimension measurement is carried out by selecting some key ground object targets (such as a vehicle body, a sign line, a pit slot and a large block crack), and real-time differential GPS equipment is utilized for positioning, writing and archiving; and then selecting different types of pavement samples (including good and poor flatness quality, coverage health and disease types) on the pavement surface of the flight area to perform sampling photographing and feature description and classified statistics, and performing position recording and photographing storage on all data acquisition works. In total, 10 key feature points of the ground feature were collected in this study, 9 different types of pavement samples (total number of sample points 52, approximately 5-7 for each pavement sample). The pavement sample types include: pit (divided into large and small blocks), collapse (divided stone is fully exposed and stone is not exposed), cracks (large and small cracks), pitted surface, sand and stone coverage, general road surface flatness and good flatness, as shown in fig. 7.
(2) Data processing
The specific processing steps of the data are shown in fig. 1, and the specific steps are executed as follows.
The first step is as follows: LiDAR point cloud data acquired by the Rigel VUX-1LR laser scanner contains attributes such as: three-dimensional coordinates, scanning angle, echo number, echo times and laser intensity. The original point cloud file is not in a universal point cloud file format, and the universal las format is obtained by converting the format of the point cloud data by means of matching software provided by Rigel company, as shown in FIG. 8;
the second step is that: denoising the LiDAR point cloud data, preliminarily removing noise points through visual identification in the experimental process, then setting a search radius to be 0.2m to remove local elevation abnormal points, calculating according to the road width and the unmanned aerial height to obtain a scanning angle range of the road surface point cloud within +/-8 degrees, and setting a scanning angle threshold to be +/-24 degrees in the experimental process in consideration of certain deviation of a route;
the third step: the orthophoto map of the road section to be measured can be obtained by performing image splicing and registration processing on high-resolution image data shot by an unmanned aerial vehicle in the road section to be measured by using commercial unmanned aerial vehicle data processing software Pix4D Mapper developed by Pix4D company of switzerland (fig. 9). Then, feature extraction is carried out on the ortho-image, the gray value standard deviation of the red, green and blue three wave bands is calculated to serve as a new feature wave band, and the new feature wave band is combined with the red, green and blue three wave bands to obtain 4 feature wave bands containing spectral information;
the fourth step: matching the point in the point cloud with the pixel in the image by using an affine transformation model, and then fusing the band value on the image into the point cloud data, wherein the result is shown in fig. 10 (taking blue band information as an example);
the fifth step: after feature fusion is performed on the point cloud data and the high-resolution image data, multi-scale geometric feature parameter calculation needs to be performed on the fused data. Because the amount of point cloud data is very large, in order to save time and not affect the reliability of results, the embodiment optimizes the storage structure of the point cloud data, namely constructs an octree storage structure, and sets a relatively proper numerical value (0.3 m is used in the embodiment) as the length of a space storage unit for resampling to obtain the point cloud data after thinning;
and a sixth step: and taking each point of the point cloud data as a center, setting a search radius within [0.2m,1.0m ] and taking 0.1m as an interval unit for searching, wherein the searched neighborhood points are point sets with the distance from the current search point to the search radius being less than or equal to the search radius in the denoised point cloud data, and then respectively calculating the local roughness (LDR), the local dimension characteristic (LDF) and the Local Height Difference (LHD) under the current scale. With the adjustment of the space search scale, point cloud data of geometric feature information of different scales can be obtained;
the seventh step: and (3) selecting partial point clouds covering all ground object types from the point cloud data obtained in the step six through visual observation to serve as training samples, then performing dimension reduction on the characteristic variables through characteristic selection (as shown in figure 11), and easily seeing that the multi-scale geometric characteristic parameters LHD after the characteristic selection have the strongest distinguishing capability on various ground objects relative to other characteristic variables. And further constructing a classifier by means of a random forest model. In this example, the number of trees is set to 200 and the depth of the trees is set to 40 during the construction process. Secondly, classifying the point cloud data after rarefaction by using the constructed random forest model so as to finish the classification process of the rarefied point cloud;
eighth step: since only the sparse point clouds are classified so far, all the point clouds after denoising need to be classified by a certain method. Because the attributes of two objects with closer spatial distance are more similar, the embodiment assigns the class attribute of the denoised point cloud data with the class attribute value of the classified point closest to the point according to the spatial nearest neighbor principle;
the ninth step: the classified point cloud data are subjected to category screening to obtain road point cloud data, a TIN model (figure 12) of the road surface is constructed through an irregular triangular network, a DSM model of the road surface is further obtained through natural neighborhood interpolation, in order to meet the requirement of IRI index calculation and simultaneously take data processing into consideration, the grid resolution is set to be 50mm, and the result is shown in figure 13;
the tenth step: in the Arcgis10, a series of longitudinal section lines with the interval of 0.125m are digitally generated in the road interior along the road direction by means of interactive vectorization, and then the IRI values on the section lines are calculated by using an IRI index model. To reflect the difference in the longitudinal direction of each link, an appropriate calculation unit length, 10m in this example, is set when calculating IRI. Rasterizing the computed IRI values to obtain a flatness spatial distribution result of the whole road section for a visualization effect, as shown in fig. 14;
the eleventh step: according to the flatness quality grading setting standard set by the method, reclassifying the obtained pavement flatness spatial distribution map to obtain a pavement flatness quality distribution map, as shown in fig. 15, determining each quality grade distribution proportion statistical map of the road section by performing statistical analysis on grid units of different quality grades, as shown in fig. 16;
(3) analysis of results
Generally, the longitudinal characteristic is that the quality evaluation result of the south road section is mainly Fair, and a small part of area is Failed; the quality evaluation of the middle part of the road is the worst, and the evaluation result is mainly 'Poor' or 'Failed'; the quality evaluation result of the north road section is more than "Fair" or "Failed", which is consistent with the spatial distribution trend of flatness. The transverse characteristic is that the quality grade near the center line is lower than that of the two-side roadway, namely 'Failed' as the main characteristic, the evaluation result near the roadway mainly adopts 'Fair' or 'Poor' as the main characteristic, the road surface quality grade of the road boundary zone is relatively lower and is 'Poor' or 'Failed', and the surface of the partial area is rough and uneven due to road surface subsidence, sand and soil mixing and the like. From the distribution diagram of each quality grade (fig. 15), it is known that the flatness quality of the whole road section is mainly "Fair", the distribution ratio thereof is 41.1%, the distribution ratio of the road section with the quality "Poor" or "Failed" is secondly, respectively 27.4% and 24.9%, and the ratio of the road section with the quality "Good" is minimum, only 6.6%. And the sum of the distribution proportion of the mass of the 'Poor' or 'Failed' accounts for more than half of the whole road section, which indicates that the whole quality of the road section deviates and needs to be repaired greatly.
(4) Evaluation of accuracy
52 road surface samples collected on the spot are taken as analysis data, and flatness quality evaluation is carried out on the road surface samples by selecting 5 indexes of appearance, surface roughness, track depth, disease area and damage degree by using an expert evaluation method, and the indexes are taken as reference values. Then map positioning is carried out on the road surface samples on the DSM model according to GPS point coordinates measured on the spot and assisted by high-resolution image data, then the flatness quality grades of the road surface samples are identified according to a road surface flatness spatial quality spatial distribution diagram, then the identification results of all the road surface samples are counted, and a confusion matrix for precision evaluation is obtained (as shown in figure 17). Finally, the overall classification accuracy of the method is 75% and the Kappa coefficient is 0.65, which are calculated by using the formula 21 and the formula 22 respectively, and the result shows that the overall evaluation result of the method is substantially the same as the actual result and the consistency degree is higher.
It is noted that the disclosed embodiments are intended to aid in further understanding of the invention, but those skilled in the art will appreciate that: various substitutions and modifications are possible without departing from the spirit and scope of the invention and appended claims. Therefore, the invention should not be limited to the embodiments disclosed, but the scope of the invention is defined by the appended claims.
Claims (9)
1. A road flatness monitoring method fusing unmanned aerial vehicle LiDAR and high-resolution remote sensing images is characterized in that a classification feature set is constructed by fusing spectral features of images and multi-scale geometric features of LiDAR point clouds based on LiDAR data and high-resolution remote sensing images acquired by a low-altitude unmanned aerial vehicle platform; then, a fine three-dimensional road surface digital surface model DSM is constructed by utilizing the road surface point cloud, the flatness parameters of the road surface are quickly calculated and obtained by an international flatness index calculation method, the road surface flatness grade is evaluated by combining with ground verification data, and the quick monitoring and evaluation of the road surface quality are realized; the method specifically comprises the following steps:
1.1) firstly, acquiring LiDAR point cloud data and high-resolution remote sensing image data by using a low-altitude unmanned remote sensing system integrating LiDAR and multispectral cameras; converting the LiDAR point cloud data into a.las file format; preprocessing LiDAR point cloud data and remote sensing image data, and extracting spectral features based on the image data;
the parameters of the spectral features are parameters reflecting the difference of the types of the ground features, and comprise: gray values of red, green and blue wave bands of the high-resolution image; the standard deviation of the gray values of the red, green and blue wave bands of the high-resolution image; the reflection intensity value of the laser point cloud;
1.2) fusing the spectral features into the attributes of the point cloud data, and performing multi-scale geometric feature extraction on the point cloud data; combining the obtained multi-scale geometric characteristic information and spectral characteristic information into a point cloud classification characteristic variable set;
the multi-scale geometric features comprise local roughness LDR, local dimension features LDF and local height difference LHD; the local roughness LDR refers to the distance from a certain reference point in the point cloud to a best fit plane formed by neighborhood points under a certain spatial scale of the reference point; the local dimension characteristic LDF refers to a characteristic metric value which is distributed in one dimension, two dimensions or three dimensions between a certain reference point in the point cloud and a neighborhood point of the reference point in space; the local elevation difference LHD is the difference between the maximum distance and the minimum distance of a best fit plane formed by reaching a point set in the whole point set formed by a certain reference point and a neighborhood point of the reference point under a certain spatial scale;
1.3) carrying out redundancy removal and dimension reduction on the point cloud classification characteristic variables, including a searching process of a characteristic subset and an evaluation process of the characteristic subset, so as to obtain an optimal characteristic subset;
the strategy of the feature subset search is forward search, and the feature subset search is carried out through dimension reduction; for the feature subset, when the newly added feature variable has high correlation with the category variable and has low correlation with the feature variable in the current subset, obtaining an optimal feature subset; the realization process is as follows:
1) establishing two data structures P and Q in a set form, wherein P represents a candidate subset, Q represents an optimal subset, a variable of a category attribute is Y, and other attributes except the category attribute are X { X ═ X1,X2,X3,...,Xn}; adding other characteristic variables except Y into the candidate subset, and setting a correlation threshold value inside the subset to be delta;
2) if it is notEntering step 4), otherwise, circularly reading each characteristic variable X in XiAnd calculating the correlation coefficient of the characteristic variable X with Y according to the formula (10), and finding out the characteristic variable X with the maximum correlation coefficient with YcAnd let δ be max { corr (X)i,Y)}，XiE.g. X, and then enter step 3);
in the formula (10), XiY respectively represent the ith feature variable value and the corresponding category attribute value of the feature set X,respectively representing the mean value of the ith characteristic variable and the mean value of the corresponding class attribute of the characteristic set X, XijAnd YjRespectively representing the ith characteristic variable value and the corresponding category attribute value of the jth sample data;
3) mixing XcPerforming relevance evaluation on each characteristic variable in Q if XcIf the correlation with any one characteristic variable in Q is less than delta, then X is addedcAdding the characteristic variable into Q, and deleting the characteristic variable from P, otherwise, directly deleting the characteristic variable from P, and continuing to the step 2);
4) outputting Q to obtain an optimal feature subset, and ending the search;
then selecting a proper machine learning classification method for point cloud classification, and further filtering and correcting gross errors to obtain more accurate road point cloud data;
1.4) constructing high-precision DSM through three-dimensional interpolation;
1.5) selecting a series of longitudinal section lines from the obtained DSM model data, calculating by using an international flatness index IRI model to obtain IRI values corresponding to the section lines, and rasterizing to obtain a spatial distribution map of the IRI values of the road surface;
1.6) the calibration is carried out by combining with ground measured data, different grades of the pavement evenness are distinguished, the type of the pavement quality from poor to good is represented, and the remote sensing monitoring and evaluation of the pavement quality are realized.
2. The method for monitoring the flatness of the pavement according to claim 1, wherein, in step 1.1), data are acquired by flying over the near ground by carrying LiDAR scanner equipment by an unmanned aerial vehicle, and the acquired data comprise LiDAR point cloud data and contemporaneous high-resolution image data; the high-resolution image data is a high-definition digital photo obtained by shooting by a multispectral camera and is provided with a corresponding posture and positioning information file provided by an unmanned aerial vehicle; the LiDAR point cloud data at least satisfies the following conditions:
2.1) the coordinate system of the LiDAR point cloud data is a WGS84 coordinate system or a local plane projection coordinate system matched with the local;
2.2) the LiDAR point cloud data contains attribute information including longitude, latitude, elevation, scan angle, reflection intensity;
2.3) the point cloud density of the LiDAR point cloud data is above 400 points per square meter.
3. The method for monitoring the pavement evenness according to claim 1, wherein the preprocessing in the step 1.1) comprises the preprocessing of point cloud data and the preprocessing of high-resolution image data; the main pretreatment process comprises the following steps: filtering noise points of the LiDAR point cloud; filtering elevation abnormal points of LiDAR point cloud; filtering the scanning angle of the LiDAR point cloud; and (4) splicing and registering the high-resolution images.
4. The method for monitoring the pavement evenness according to claim 1, wherein the step 1.2) of fusing the spectral features into the attributes of the point cloud data is to fuse the spectral features of the image into the attributes of the point cloud based on a spatial reference of the point cloud data; specifically, an affine inverse transformation model is adopted to match the wave band value information of each pixel on the image to the point cloud data coordinates of the corresponding position.
5. The method for monitoring the flatness of the road surface according to claim 1, wherein the multi-scale geometric features are specifically: setting a plurality of spatial analysis scales when calculating a certain geometric characteristic value; fixing a certain spatial analysis scale for calculation to obtain a corresponding spatial geometric characteristic parameter; the geometric characteristics of different ground object types under different spatial analysis scales are different; thereby obtaining a multi-scale geometric feature set.
6. The method for monitoring the pavement evenness as claimed in claim 1, wherein in step 1.3), specifically, a greedy strategy is adopted for the dimension reduction and feature subset search; the machine learning classification method comprises the steps of classifying sample point clouds by adopting a random forest classification algorithm, and then interpolating the category attributes of the sample point clouds by adopting a point cloud interpolation algorithm to obtain the classification of overall point cloud data; and when the classification result of the overall point cloud data contains point clouds which do not belong to the road surface, carrying out manual inspection and gross error correction on the classified point cloud data, and further filtering through category attributes to obtain relatively accurate road point cloud data.
7. The method for monitoring the flatness of the road surface according to claim 6, wherein the point cloud interpolation algorithm adopts a natural neighborhood method.
8. The method of claim 1, wherein the DSM of step 1.4) has a grid resolution controlled to within 50mm and an elevation accuracy within 15 mm; and step 1.5) when the international flatness index IRI model is adopted to calculate the road surface flatness value, the sampling interval is controlled within 50mm, and the length of a measuring unit is 10-20 m.
9. The method for monitoring the road flatness according to claim 1, wherein after the step 1.6) obtains the road flatness, namely the IRI value, a flatness quality grading standard is set, and the road is classified according to the road flatness quality grading standard; further adopting a remote sensing image classification precision evaluation method to carry out precision evaluation, thereby obtaining a flatness quality evaluation result corresponding to the road section to be measured; the precision evaluation method comprises the following steps:
10.1) selecting a contemporaneous road surface sample as reference data on a road section to be detected according to the flatness quality grading standard; the selected contemporaneous pavement sample comprises pavement samples with good and poor flatness quality, and health coverage and disease types; randomly selecting a plurality of samples from each road surface sample as reference data;
10.2) selecting indexes of appearance, surface roughness, rutting depth, disease area and damage degree for the reference data, carrying out flatness quality evaluation by using an Analytic Hierarchy Process (AHP) expert evaluation method, and taking an evaluation score as a reference value;
10.3) grading the pavement flatness value of the pavement sample according to the set flatness quality grading standard to obtain a flatness grading result corresponding to the pavement;
10.4) calculating the grading result of the flatness by adopting the precision evaluation index to obtain the monitoring precision of the pavement; the accuracy evaluation index includes an Overall Accuracy (OA) or/and a Kappa coefficient, and a larger value of both indicates a higher accuracy.
Priority Applications (1)
| Application Number | Priority Date | Filing Date | Title |
|---|---|---|---|
| CN201710257172.0A CN107092020B (en) | 2017-04-19 | 2017-04-19 | Merge the surface evenness monitoring method of unmanned plane LiDAR and high score image |
Applications Claiming Priority (1)
| Application Number | Priority Date | Filing Date | Title |
|---|---|---|---|
| CN201710257172.0A CN107092020B (en) | 2017-04-19 | 2017-04-19 | Merge the surface evenness monitoring method of unmanned plane LiDAR and high score image |
Publications (2)
| Publication Number | Publication Date |
|---|---|
| CN107092020A CN107092020A (en) | 2017-08-25 |
| CN107092020B true CN107092020B (en) | 2019-09-13 |
Family
ID=59636957
Family Applications (1)
| Application Number | Title | Priority Date | Filing Date |
|---|---|---|---|
| CN201710257172.0A Active CN107092020B (en) | 2017-04-19 | 2017-04-19 | Merge the surface evenness monitoring method of unmanned plane LiDAR and high score image |
Country Status (1)
| Country | Link |
|---|---|
| CN (1) | CN107092020B (en) |
Families Citing this family (54)
| Publication number | Priority date | Publication date | Assignee | Title |
|---|---|---|---|---|
| CN107576960B (en) * | 2017-09-04 | 2021-03-16 | 赵建辉 | Target detection method and system for visual radar space-time information fusion |
| CN107657636B (en) * | 2017-10-16 | 2021-07-30 | 南京市测绘勘察研究院股份有限公司 | Method for automatically extracting road topographic map elevation points based on vehicle-mounted laser radar data |
| CN109073744A (en) * | 2017-12-18 | 2018-12-21 | 深圳市大疆创新科技有限公司 | Landform prediction technique, equipment, system and unmanned plane |
| CN108198230A (en) * | 2018-02-05 | 2018-06-22 | 西北农林科技大学 | A kind of crop and fruit three-dimensional point cloud extraction system based on image at random |
| CN108492329B (en) * | 2018-03-19 | 2022-04-12 | 北京航空航天大学 | Three-dimensional reconstruction point cloud precision and integrity evaluation method |
| CN108595553B (en) * | 2018-04-10 | 2022-02-08 | 红云红河烟草(集团)有限责任公司 | Industrial data acquisition time sequence data compression storage and decompression query method based on relational database |
| CN108562913B (en) * | 2018-04-19 | 2021-12-17 | 武汉大学 | Unmanned ship false target detection method based on three-dimensional laser radar |
| CN108564096B (en) * | 2018-04-26 | 2019-11-19 | 电子科技大学 | A kind of neighborhood fitting RCS sequence characteristic extracting method |
| CN109190673B (en) * | 2018-08-03 | 2021-09-03 | 西安电子科技大学 | Ground target classification method based on random forest and data rejection |
| CN109446983A (en) * | 2018-10-26 | 2019-03-08 | 福州大学 | A kind of coniferous forest felling accumulation evaluation method based on two phase unmanned plane images |
| CN109459759B (en) * | 2018-11-13 | 2020-06-30 | 中国科学院合肥物质科学研究院 | Urban terrain three-dimensional reconstruction method based on quad-rotor unmanned aerial vehicle laser radar system |
| CN109684929A (en) * | 2018-11-23 | 2019-04-26 | 中国电建集团成都勘测设计研究院有限公司 | Terrestrial plant ECOLOGICAL ENVIRONMENTAL MONITORING method based on multi-sources RS data fusion |
| CN109508508B (en) * | 2018-12-08 | 2024-03-22 | 河北省地质矿产勘查开发局国土资源勘查中心(河北省矿山和地质灾害应急救援中心) | Surface mine governance investigation design method |
| US11675083B2 (en) * | 2019-01-03 | 2023-06-13 | Nvidia Corporation | Removal of ephemeral points from point cloud of a high-definition map for navigating autonomous vehicles |
| CN109829425B (en) * | 2019-01-31 | 2020-12-22 | 沈阳农业大学 | Farmland landscape small-scale ground feature classification method and system |
| CN110084116B (en) * | 2019-03-22 | 2022-02-01 | 深圳市速腾聚创科技有限公司 | Road surface detection method, road surface detection device, computer equipment and storage medium |
| CN109945853B (en) * | 2019-03-26 | 2023-08-15 | 西安因诺航空科技有限公司 | Geographic coordinate positioning system and method based on 3D point cloud aerial image |
| US11618502B2 (en) * | 2019-03-28 | 2023-04-04 | Volkswagen Aktiengesellschaft | On-road localization methodologies and equipment utilizing road surface characteristics |
| CN110009054B (en) * | 2019-04-12 | 2021-01-29 | 南京大学 | Hierarchical airborne LiDAR point cloud classification method utilizing geometric and intensity features |
| CN110197223B (en) * | 2019-05-29 | 2021-02-09 | 北方民族大学 | Point cloud data classification method based on deep learning |
| CN110210500A (en) * | 2019-06-06 | 2019-09-06 | 上海黑塞智能科技有限公司 | A kind of point cloud classifications method based on the insertion of multiple dimensioned local feature |
| CN110070544A (en) * | 2019-06-06 | 2019-07-30 | 江苏省农业科学院 | One planting fruit-trees target three-dimensional data compensation method and compensation system |
| CN110264572B (en) * | 2019-06-21 | 2021-07-30 | 哈尔滨工业大学 | Terrain modeling method and system integrating geometric characteristics and mechanical characteristics |
| CN110490888B (en) * | 2019-07-29 | 2022-08-30 | 武汉大学 | Highway geometric feature vectorization extraction method based on airborne laser point cloud |
| CN110516653A (en) * | 2019-09-03 | 2019-11-29 | 武汉天擎空间信息技术有限公司 | A kind of method for extracting roads based on multispectral airborne laser radar point cloud data |
| CN110726998B (en) * | 2019-10-24 | 2020-08-07 | 西安科技大学 | Method for measuring mining subsidence basin in mining area through laser radar scanning |
| CN110763223B (en) * | 2019-10-31 | 2022-03-18 | 苏州大学 | Sliding window based indoor three-dimensional grid map feature point extraction method |
| CN112470032A (en) * | 2019-11-04 | 2021-03-09 | 深圳市大疆创新科技有限公司 | Topographic prediction method and device for undulating ground, radar, unmanned aerial vehicle and operation control method |
| CN110658850B (en) * | 2019-11-12 | 2022-07-12 | 重庆大学 | Greedy strategy-based flight path planning method for unmanned aerial vehicle |
| CN110670461A (en) * | 2019-11-14 | 2020-01-10 | 上海宝冶建筑工程有限公司 | Method for detecting flatness of airport pavement |
| CN110880202B (en) * | 2019-12-02 | 2023-03-21 | 中电科特种飞机***工程有限公司 | Three-dimensional terrain model creating method, device, equipment and storage medium |
| CN111381248B (en) * | 2020-03-23 | 2023-06-09 | 湖南大学 | Obstacle detection method and system considering vehicle bump |
| CN111638185B (en) * | 2020-05-09 | 2022-05-17 | 哈尔滨工业大学 | Remote sensing detection method based on unmanned aerial vehicle platform |
| US11769238B2 (en) | 2020-05-18 | 2023-09-26 | Roadbotics, Inc. | Systems and methods for creating and/or analyzing three-dimensional models of infrastructure assets |
| CN111784966A (en) * | 2020-06-15 | 2020-10-16 | 武汉烽火众智数字技术有限责任公司 | Personnel management and control method and system based on machine learning |
| CN112446343B (en) * | 2020-12-07 | 2024-03-15 | 园测信息科技股份有限公司 | Vehicle-mounted point cloud road shaft-shaped object machine learning automatic extraction method integrating multi-scale features |
| CN112633092B (en) * | 2020-12-09 | 2022-06-14 | 西南交通大学 | Road information extraction method based on vehicle-mounted laser scanning point cloud |
| CN112686396B (en) * | 2020-12-15 | 2024-02-23 | 中公高科养护科技股份有限公司 | Pavement maintenance property selection method, medium and system based on disease quantity |
| CN112819273B (en) * | 2020-12-31 | 2023-11-28 | 东北大学 | Regional highway network quality evaluation method and system |
| CN113326865A (en) * | 2021-04-15 | 2021-08-31 | 东南大学 | Highway road surface disease three-dimensional information detecting system based on deep learning |
| CN113325440B (en) * | 2021-05-06 | 2022-08-16 | 武汉大学 | Polarization laser radar data inversion method and system based on image recognition and signal characteristic decomposition |
| CN113551636B (en) * | 2021-07-02 | 2023-06-06 | 武汉光谷卓越科技股份有限公司 | Flatness detection method based on abnormal data correction |
| CN113962301B (en) * | 2021-10-20 | 2022-06-17 | 北京理工大学 | Multi-source input signal fused pavement quality detection method and system |
| CN113689565B (en) * | 2021-10-21 | 2022-03-18 | 北京中科慧眼科技有限公司 | Road flatness grade detection method and system based on binocular stereo vision and intelligent terminal |
| CN114155447B (en) * | 2021-12-02 | 2022-06-24 | 北京中科智易科技有限公司 | Artificial intelligence big data acquisition system |
| CN114646305B (en) * | 2022-03-03 | 2024-04-02 | 湖南省测绘科技研究所 | Intelligent recognition method for unmanned aerial vehicle mapping behavior |
| CN114283070B (en) * | 2022-03-07 | 2022-05-03 | 中国铁路设计集团有限公司 | Method for manufacturing terrain section by fusing unmanned aerial vehicle image and laser point cloud |
| CN114858214B (en) * | 2022-04-27 | 2023-08-25 | 中徽建技术有限公司 | Urban road performance monitoring system |
| CN114863695B (en) * | 2022-05-30 | 2023-04-18 | 中邮建技术有限公司 | Overproof vehicle detection system and method based on vehicle-mounted laser and camera |
| CN115323876A (en) * | 2022-07-15 | 2022-11-11 | 东南大学 | Airport cement concrete pavement flatness detection system |
| CN115290012B (en) * | 2022-08-09 | 2023-03-31 | 南京市计量监督检测院 | Road surface flatness standard field magnitude tracing method based on laser point cloud data |
| CN115168975B (en) * | 2022-08-17 | 2024-02-09 | 中国建筑第二工程局有限公司 | BIM technology-based pavement evenness quality control method and device |
| CN116612400B (en) * | 2023-05-30 | 2024-03-19 | 衡水金湖交通发展集团有限公司 | Road management method and system based on road flatness |
| CN117437368B (en) * | 2023-12-22 | 2024-04-26 | 深圳大学 | Unmanned plane-based pavement evenness measuring method, system, terminal and medium |
Citations (13)
| Publication number | Priority date | Publication date | Assignee | Title |
|---|---|---|---|---|
| JPH09264743A (en) * | 1996-03-29 | 1997-10-07 | Pub Works Res Inst Ministry Of Constr | Control method for finished form of pavement |
| CN102103202A (en) * | 2010-12-01 | 2011-06-22 | 武汉大学 | Semi-supervised classification method for airborne laser radar data fusing images |
| CN102163342A (en) * | 2011-04-07 | 2011-08-24 | 北京农业信息技术研究中心 | Fruit three morphological structure three-dimensional (3D) reconstruction method based on multi-scale measurement data |
| CN102902883A (en) * | 2012-09-24 | 2013-01-30 | 北京师范大学 | Method for establishing bidirectional reflectance distribution function (BRDF) prototype library based on multi-angular measurement |
| CN103017739A (en) * | 2012-11-20 | 2013-04-03 | 武汉大学 | Manufacturing method of true digital ortho map (TDOM) based on light detection and ranging (LiDAR) point cloud and aerial image |
| CN104233935A (en) * | 2014-08-28 | 2014-12-24 | 吉林大学 | Identification method for pavement quality grade on basis of information of longitudinal section of road |
| CN104867180A (en) * | 2015-05-28 | 2015-08-26 | 南京林业大学 | UAV and LiDAR integrated forest stand characteristic inversion method |
| CN104933708A (en) * | 2015-06-07 | 2015-09-23 | 浙江大学 | Barrier detection method in vegetation environment based on multispectral and 3D feature fusion |
| CN105488770A (en) * | 2015-12-11 | 2016-04-13 | 中国测绘科学研究院 | Object-oriented airborne laser radar point cloud filtering method |
| CN106056591A (en) * | 2016-05-25 | 2016-10-26 | 哈尔滨工业大学 | Method for estimating urban density through fusion of optical spectrum image and laser radar data |
| CN106127771A (en) * | 2016-06-28 | 2016-11-16 | 上海数联空间科技有限公司 | Tunnel orthography system and method is obtained based on laser radar LIDAR cloud data |
| CN106296814A (en) * | 2015-05-26 | 2017-01-04 | 中国公路工程咨询集团有限公司 | Highway maintenance detection and virtual interactive interface method and system |
| CN106355643A (en) * | 2016-08-31 | 2017-01-25 | 武汉理工大学 | Method for generating three-dimensional real scene road model of highway |
-
2017
- 2017-04-19 CN CN201710257172.0A patent/CN107092020B/en active Active
Patent Citations (13)
| Publication number | Priority date | Publication date | Assignee | Title |
|---|---|---|---|---|
| JPH09264743A (en) * | 1996-03-29 | 1997-10-07 | Pub Works Res Inst Ministry Of Constr | Control method for finished form of pavement |
| CN102103202A (en) * | 2010-12-01 | 2011-06-22 | 武汉大学 | Semi-supervised classification method for airborne laser radar data fusing images |
| CN102163342A (en) * | 2011-04-07 | 2011-08-24 | 北京农业信息技术研究中心 | Fruit three morphological structure three-dimensional (3D) reconstruction method based on multi-scale measurement data |
| CN102902883A (en) * | 2012-09-24 | 2013-01-30 | 北京师范大学 | Method for establishing bidirectional reflectance distribution function (BRDF) prototype library based on multi-angular measurement |
| CN103017739A (en) * | 2012-11-20 | 2013-04-03 | 武汉大学 | Manufacturing method of true digital ortho map (TDOM) based on light detection and ranging (LiDAR) point cloud and aerial image |
| CN104233935A (en) * | 2014-08-28 | 2014-12-24 | 吉林大学 | Identification method for pavement quality grade on basis of information of longitudinal section of road |
| CN106296814A (en) * | 2015-05-26 | 2017-01-04 | 中国公路工程咨询集团有限公司 | Highway maintenance detection and virtual interactive interface method and system |
| CN104867180A (en) * | 2015-05-28 | 2015-08-26 | 南京林业大学 | UAV and LiDAR integrated forest stand characteristic inversion method |
| CN104933708A (en) * | 2015-06-07 | 2015-09-23 | 浙江大学 | Barrier detection method in vegetation environment based on multispectral and 3D feature fusion |
| CN105488770A (en) * | 2015-12-11 | 2016-04-13 | 中国测绘科学研究院 | Object-oriented airborne laser radar point cloud filtering method |
| CN106056591A (en) * | 2016-05-25 | 2016-10-26 | 哈尔滨工业大学 | Method for estimating urban density through fusion of optical spectrum image and laser radar data |
| CN106127771A (en) * | 2016-06-28 | 2016-11-16 | 上海数联空间科技有限公司 | Tunnel orthography system and method is obtained based on laser radar LIDAR cloud data |
| CN106355643A (en) * | 2016-08-31 | 2017-01-25 | 武汉理工大学 | Method for generating three-dimensional real scene road model of highway |
Non-Patent Citations (1)
| Title |
|---|
| 机载LiDAR点云数据与影像数据融合;尚大帅等;《测绘工程》;20120229;第21卷(第1期);18-20,24 * |
Also Published As
| Publication number | Publication date |
|---|---|
| CN107092020A (en) | 2017-08-25 |
Similar Documents
| Publication | Publication Date | Title |
|---|---|---|
| CN107092020B (en) | Merge the surface evenness monitoring method of unmanned plane LiDAR and high score image | |
| KR101933216B1 (en) | River topography information generation method using drone and geospatial information | |
| CN104376595B (en) | A kind of three-dimensional road generation method cooperateed with based on airborne LiDAR and GIS | |
| CN102506824B (en) | Method for generating digital orthophoto map (DOM) by urban low altitude unmanned aerial vehicle | |
| CN113034689A (en) | Laser point cloud-based terrain three-dimensional model, terrain map construction method and system, and storage medium | |
| Aplin | 2 Comparison of simulated IKONOS and SPOT HRV imagery for | |
| Karantanellis et al. | 3D hazard analysis and object-based characterization of landslide motion mechanism using UAV imagery | |
| Atencio et al. | Calibration of UAV flight parameters for pavement pothole detection using orthogonal arrays | |
| Demir | Automated detection of 3D roof planes from Lidar data | |
| CN110580468B (en) | Single wood structure parameter extraction method based on image matching point cloud | |
| Zeybek et al. | Road surface and inventory extraction from mobile LiDAR point cloud using iterative piecewise linear model | |
| Suleymanoglu et al. | Road infrastructure mapping by using iPhone 14 Pro: An accuracy assessment | |
| CN114004740B (en) | Building wall line extraction method based on unmanned aerial vehicle laser radar point cloud | |
| Pahlavani et al. | 3D reconstruction of buildings from LiDAR data considering various types of roof structures | |
| Saponaro et al. | Influence of co-alignment procedures on the co-registration accuracy of multi-epoch SFM points clouds | |
| Gonçalves et al. | Planar projection of mobile laser scanning data in tunnels | |
| Sefercik et al. | Point-based and model-based geolocation analysis of airborne laser scanning data | |
| Polat | LIDAR Derived 3d City Modelling | |
| Susetyo et al. | Automatic building model extraction using LiDAR data | |
| Huang et al. | Ground filtering algorithm for mobile LIDAR using order and neighborhood point information | |
| Di Stefano et al. | Point cloud classification of an urban environment using a semi-automatic approach | |
| DOWAJY et al. | Comparative Analysis of Road Scanning Techniques | |
| Yu et al. | A cue line based method for building modeling from LiDAR and satellite imagery | |
| Wang | Semi-automated generation of high-accuracy digital terrain models along roads using mobile laser scanning data | |
| Kersting et al. | Object oriented transmission line corridors classification using lidar technology and a non-metric digital camera |
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 |