CN104050473A - Road data extraction method based on rectangular neighborhood analysis - Google Patents
Road data extraction method based on rectangular neighborhood analysis Download PDFInfo
- Publication number
- CN104050473A CN104050473A CN201410213465.5A CN201410213465A CN104050473A CN 104050473 A CN104050473 A CN 104050473A CN 201410213465 A CN201410213465 A CN 201410213465A CN 104050473 A CN104050473 A CN 104050473A
- Authority
- CN
- China
- Prior art keywords
- point
- grid
- road
- elevation
- points
- Prior art date
- Legal status (The legal status is an assumption and is not a legal conclusion. Google has not performed a legal analysis and makes no representation as to the accuracy of the status listed.)
- Granted
Links
- 238000000034 method Methods 0.000 title claims abstract description 29
- 238000004458 analytical method Methods 0.000 title claims abstract description 14
- 238000013075 data extraction Methods 0.000 title claims description 10
- 238000007781 pre-processing Methods 0.000 claims abstract description 4
- 238000005260 corrosion Methods 0.000 claims description 14
- 230000007797 corrosion Effects 0.000 claims description 14
- 238000001914 filtration Methods 0.000 claims description 12
- 230000003628 erosive effect Effects 0.000 claims description 5
- 230000008520 organization Effects 0.000 claims description 4
- 238000004364 calculation method Methods 0.000 claims description 3
- 238000012216 screening Methods 0.000 abstract 1
- 238000000605 extraction Methods 0.000 description 27
- 238000010586 diagram Methods 0.000 description 7
- 238000005516 engineering process Methods 0.000 description 3
- 238000013507 mapping Methods 0.000 description 3
- 238000002310 reflectometry Methods 0.000 description 3
- 230000000877 morphologic effect Effects 0.000 description 2
- 238000012545 processing Methods 0.000 description 2
- 238000013138 pruning Methods 0.000 description 2
- 238000011160 research Methods 0.000 description 2
- ATJFFYVFTNAWJD-UHFFFAOYSA-N Tin Chemical compound [Sn] ATJFFYVFTNAWJD-UHFFFAOYSA-N 0.000 description 1
- 230000001427 coherent effect Effects 0.000 description 1
- 238000012937 correction Methods 0.000 description 1
- 238000001514 detection method Methods 0.000 description 1
- 230000010339 dilation Effects 0.000 description 1
- 230000000694 effects Effects 0.000 description 1
- 230000008030 elimination Effects 0.000 description 1
- 238000003379 elimination reaction Methods 0.000 description 1
- 238000012986 modification Methods 0.000 description 1
- 230000004048 modification Effects 0.000 description 1
- 230000011218 segmentation Effects 0.000 description 1
Landscapes
- Image Processing (AREA)
- Traffic Control Systems (AREA)
Abstract
本发明公开了一种基于矩形邻域分析的道路数据提取方法,其包括数据预处理步骤、筛选道路候选点步骤、提取道路中心线和剔除平坦地域点集步骤、道路连接步骤。本发明主要解决在没有反射强度或者反射强度不能较好反应地物属性,或者没有航空影像等情况,根据道路拓扑结构特点与周围环境的不同进行道路数据提取。
The invention discloses a method for extracting road data based on rectangular neighborhood analysis, which includes a data preprocessing step, a step of screening road candidate points, a step of extracting a road center line and eliminating flat area point sets, and a road connection step. The invention mainly solves the problem of extracting road data according to the difference between road topological structure characteristics and surrounding environment when there is no reflection intensity or the reflection intensity cannot better reflect the ground object attributes, or there is no aerial image.
Description
技术领域technical field
本发明涉及一种基于矩形邻域分析的道路数据提取方法。The invention relates to a road data extraction method based on rectangular neighborhood analysis.
背景技术Background technique
建筑物和道路提取一直是机载LiDAR(Light Detection AndRanging)数据特征提取的重点和难点,在道路提取方面,主要围绕城区道路网提取研究较多,针对野外地域中道路提取方法研究较少。其中机载LiDAR是一种集激光测距、全球定位系统(GPS)和惯性导航系统(INS)三种技术于一体的系统,并用于获取数据并生成精确三维地形(DEM);机载LiDAR(或称点云数据)是由三维激光扫描仪采集的物体表面点数据,每个点有空间三维坐标,有的还带有颜色和反色强度等信息;点云滤波是指去除点云中的噪声、树木、房屋等非地面点;点云组织是指将点云按照一定的数据格式进行管理,便于数据的查询和显示;道路提取是指从点云数据中自动提取道路点。现有的道路提取方法主要有基于反射强信息的道路提取方法、基于机载LiDAR数据和遥感影像的道路提取方法等,常用的基于点云、反射强度和影像的道路数据提取方法如下:The extraction of buildings and roads has always been the focus and difficulty of airborne LiDAR (Light Detection And Ranging) data feature extraction. In terms of road extraction, there are more researches on the extraction of road networks in urban areas, and less research on road extraction methods in wild areas. Among them, airborne LiDAR is a system that integrates three technologies of laser ranging, global positioning system (GPS) and inertial navigation system (INS), and is used to acquire data and generate accurate three-dimensional terrain (DEM); airborne LiDAR ( Or point cloud data) is the point data on the surface of an object collected by a 3D laser scanner. Each point has three-dimensional coordinates in space, and some also have information such as color and anti-color intensity; Non-ground points such as noise, trees, and houses; point cloud organization refers to managing point clouds according to a certain data format, which is convenient for data query and display; road extraction refers to automatically extracting road points from point cloud data. Existing road extraction methods mainly include road extraction methods based on strong reflection information, road extraction methods based on airborne LiDAR data and remote sensing images, etc. Commonly used road data extraction methods based on point clouds, reflection intensity, and images are as follows:
(1)车载激光扫描数据的结构化道路自动提取方法(1) Automatic extraction method of structured roads from vehicle-mounted laser scanning data
车载激光扫描是指将激光扫描仪架设在汽车等移动平台上进行动态车辆,与机载激光扫描不同点在于机载激光扫描的范围更大、速度更快,而车载激光扫描获取的数据都是沿着道路的,因此提取较为容易。Vehicle-mounted laser scanning refers to setting up a laser scanner on a mobile platform such as a car for dynamic vehicles. The difference from airborne laser scanning is that the range of airborne laser scanning is larger and the speed is faster, while the data obtained by vehicle-mounted laser scanning are all Along the road, so extraction is easier.
结构化道路是指城区道路或者高速公路等,这些道路的两旁都有高出道路点的特征,通过扫描线提取、地面点云滤波、道路边界自动提取、道路边界跟踪等步骤完成道路的提取。具体方法详见《测绘学报》2013年42卷2期《车载激光扫描数据的结构化道路自动提取方法》。Structured roads refer to urban roads or expressways, etc. Both sides of these roads have features higher than road points. Road extraction is completed through steps such as scanning line extraction, ground point cloud filtering, automatic road boundary extraction, and road boundary tracking. For the specific method, please refer to "Automatic Extraction Method of Structured Roads from Vehicular Laser Scanning Data" in "Journal of Surveying and Mapping", Volume 42, Issue 2, 2013.
(2)基于LiDAR回波信息的道路提取(2) Road extraction based on LiDAR echo information
首先采用一种方法(如:逐层加密TIN方法)提取或者生成数字地形模型(DTM),根据点云的回波信息从DTM中提取道路信息,最后通过搜索孤立点的滤波算法删除其中的噪声点。具体方法详见《测绘科学》2011年36卷2期《基于LiDAR回波信息的道路提取》。First, a method (such as layer-by-layer encryption TIN method) is used to extract or generate a digital terrain model (DTM), and road information is extracted from the DTM according to the echo information of the point cloud, and finally the noise is deleted through the filtering algorithm of searching isolated points point. The specific method is detailed in "Surveying and Mapping Science", Volume 36, Issue 2, "Road Extraction Based on LiDAR Echo Information" in 2011.
(3)基于机载LiDAR粗糙度指数和回波强度的道路提取(3) Road extraction based on airborne LiDAR roughness index and echo intensity
首先,由点云数据衍生归一化数字表面模型和粗糙度指数,然后对配准后的多源数据进行多分辨率分割,进而使用粗糙度指数和回波强度、道路描述因子等特征进行分类,最后,去除道路噪声,并获取准确的道路架网。详见《测绘科学技术学报》2013年30卷1期《基于机载LiDAR粗糙度指数和回波强度的道路提取》。First, the normalized digital surface model and roughness index are derived from the point cloud data, and then multi-resolution segmentation is performed on the registered multi-source data, and then the roughness index, echo intensity, road description factor and other features are used for classification , and finally, remove road noise and obtain an accurate road network. For details, see "Road Extraction Based on Airborne LiDAR Roughness Index and Echo Intensity" in Journal of Surveying and Mapping Science and Technology, Volume 30, Issue 1, 2013.
(4)基于机载LiDAR和高分辨率遥感影像的城市道路网提取(4) Urban road network extraction based on airborne LiDAR and high-resolution remote sensing images
首先将两者进行精确配准然后利用伪道路信息去除的方法分别将植被信息和建筑物信息等去除得到基本的道路轮廓,再利用形态细化算法提取道路的中心线,最后在ArcGIS和Matlab编程环境下实现了改进的道路修剪算法(IRT)利用该算法进行道路修剪得到了平滑和连贯的城市道路网。详见《遥感技术与应用》2013年28卷4期《基于机载LiDAR和高分辨率遥感影像的城市道路网提取》。First, the two are accurately registered, and then the vegetation information and building information are removed by the method of pseudo road information removal to obtain the basic road outline, and then the centerline of the road is extracted by using the morphological refinement algorithm, and finally programmed in ArcGIS and Matlab Under the environment, the improved road pruning algorithm (IRT) is realized, and the smooth and coherent urban road network is obtained by using this algorithm for road pruning. For details, see "Urban Road Network Extraction Based on Airborne LiDAR and High Resolution Remote Sensing Images" in "Remote Sensing Technology and Application", Volume 28, Issue 4, 2013.
然而激光回波强度与地物目标的反射率成正比关系,即反射率越大激光强度值越高,但是激光雷达系统记录的回波强度是尚未经过校正的值,会受大气衰减、激光入射角度等多种因素的影响,而且激光回波强度的校正至今还是一个比较困难的课题,因此激光回波强度很难真实地反映地物的反射率信息。基于遥感影像的道路提取法的提取效果主要取决于影像的清晰程度,但是在多云或者能见度不够的情况下采集的遥感影像往往并不能很好地反应地物属性,因此分类也较为困难。However, the laser echo intensity is directly proportional to the reflectivity of ground objects, that is, the greater the reflectivity, the higher the laser intensity value, but the echo intensity recorded by the laser radar system is an uncorrected value, which will be affected by atmospheric attenuation, laser incident Influenced by many factors such as angle, and the correction of laser echo intensity is still a relatively difficult subject, it is difficult for laser echo intensity to truly reflect the reflectivity information of ground objects. The extraction effect of the road extraction method based on remote sensing images mainly depends on the clarity of the images, but the remote sensing images collected under cloudy or insufficient visibility often cannot reflect the attributes of ground features well, so the classification is also more difficult.
发明内容Contents of the invention
本发明目的是:提供一种基于矩形邻域分析的道路数据提取方法,在不采用反射强度和遥感影像的条件下,通过分析丘林地域道路的结构特征,采用矩形邻域判断分析法进行道路的提取。The purpose of the present invention is to provide a method for extracting road data based on rectangular neighborhood analysis. Under the condition of not using reflection intensity and remote sensing images, by analyzing the structural characteristics of roads in hilly forest areas, the rectangular neighborhood judgment analysis method is used for road extraction. extraction.
本发明的技术方案是:一种基于矩形邻域分析的道路数据提取方法,其包括:The technical scheme of the present invention is: a kind of road data extraction method based on rectangular neighborhood analysis, it comprises:
第一步:数据预处理,计算点云间距,建立规则网格组织,然后采用基于规则网格的快速数学形态学滤波方法进行滤波,剔除非地面点;The first step: data preprocessing, calculate point cloud spacing, establish regular grid organization, and then use fast mathematical morphology filtering method based on regular grid to filter, and eliminate non-ground points;
第二步:筛选道路候选点,初步确定道路点集,缩小判断范围,对每个点搜索其邻域范围内的点,并比较其高程值,若都在阈值范围内,则作为道路候选点;Step 2: Screen candidate road points, preliminarily determine the road point set, narrow the judgment range, search for points within its neighborhood for each point, and compare their elevation values. If they are all within the threshold range, they will be used as road candidate points ;
第三步:提取道路中心线,剔除平坦地域点集,对每个道路候选点,计算不同方位角的矩形范围内高程方差值,根据方差值间的差异判断当前点是否为道路中心点;Step 3: Extract the road centerline, remove the flat area point set, calculate the elevation variance value within the rectangular range of different azimuths for each road candidate point, and judge whether the current point is the road center point according to the difference between the variance values ;
第四步:道路连接,对间断的道路进行连接,平滑。The fourth step: road connection, connect and smooth the discontinuous road.
在上述技术方案的基础上,进一步包括如下附属技术方案:On the basis of the above-mentioned technical solutions, the following subsidiary technical solutions are further included:
所述第一步包括:The first step includes:
(1)首先读取点云数据范围(xmin,xmax,ymin,ymax),统计点的数量pt_Num;(1) First read the point cloud data range (xmin, xmax, ymin, ymax), and count the number of points pt_Num;
(2)计算点云间距:(2) Calculate point cloud spacing:
(3)根据点云间距pt_dis和计算需要,建立覆盖点云范围(xmin,xmax,ymin,ymax)的规则网格;(3) According to the point cloud spacing pt_dis and calculation needs, establish a regular grid covering the point cloud range (xmin, xmax, ymin, ymax);
(4)根据点的坐标将机载LiDAR数据分配到各个网格中,(4) distribute the airborne LiDAR data into individual grids according to the coordinates of the points,
GridNum=XAxisNum*YAxisNum (公式3)GridNum=XAxisNum*YAxisNum (Formula 3)
XAxisNum、YAxisNum为x、y方向的网格数,GridNum为网格总数;XAxisNum and YAxisNum are the number of grids in the x and y directions, and GridNum is the total number of grids;
机载LiDAR数据中的任意一点pi(xi,yi)所在的网格为:The grid of any point p i ( xi , y i ) in the airborne LiDAR data is:
其中pi为pi所在的X方向网格数、Yi为Pi所在的Y方向网格数、Gi为Pi所在的一维数组记录的网格数;Among them, p i is the number of grids in the X direction where p i is located, Yi is the number of grids in the Y direction where Pi is located, and Gi is the number of grids recorded in the one-dimensional array where Pi is located;
对机载LiDAR数据进行非地面滤波,保留包含道路点在内的地面点数据,以规则网格为单元的数学形态学滤波,网格腐蚀结果为w×w网格内的最小高程值,膨胀结果为w×w网格内最大高程值,其中w代表窗口的大小,并进行公式5的运算:Carry out non-ground filtering on the airborne LiDAR data, retain ground point data including road points, and use regular grid as the unit of mathematical morphology filtering. The result of grid erosion is the minimum elevation value in the w×w grid, and the expansion The result is the maximum elevation value in the w×w grid, where w represents the size of the window, and the operation of formula 5 is performed:
其中A代表目标集合,B代表结构元素,AΘB表示腐蚀,表示膨胀,式中zi′和zi′′分别为网格单元腐蚀和膨胀后的高程值,Gridi(zmin)为网格内高程最小值。where A represents the target set, B represents the structural element, AΘB represents corrosion, In the formula, zi ′ and zi ′ ′ are the elevation values after grid cell erosion and expansion respectively, and Grid i (z min ) is the minimum elevation value in the grid.
所述第一步进一步包括:The first step further includes:
(1)由网格间距计算网格数量,并为网格分配内存空间,将点逐个按平面坐标存储到相应的网格单元中,同时与该网格单元中已保存的最小高程值进行比较,记录比较后的最小高程值;(1) Calculate the number of grids from the grid spacing, and allocate memory space for the grid, store the points one by one in the corresponding grid unit according to the plane coordinates, and compare with the minimum elevation value saved in the grid unit , record the minimum elevation value after comparison;
(2)遍历每个网格,比较周围w×w个网格内最小高程值,采用公式5的第一个公式,取最小值为当前网格腐蚀后的高程值,并将该网格最大高程值设为腐蚀后的高程值;(2) Traverse each grid, compare the minimum elevation value in the surrounding w×w grids, use the first formula of formula 5, take the minimum value as the elevation value of the current grid after corrosion, and set the maximum value of the grid The elevation value is set to the elevation value after corrosion;
(3)遍历每个网格,比较周围w×w个网格内腐蚀高程值,采用公式5的第二个公式,取最大值为当前网格膨胀后的高程值,并将该网格最小高程值设为膨胀后的高程值;(3) Traversing through each grid, comparing the corrosion elevation values in the surrounding w×w grids, using the second formula of formula 5, taking the maximum value as the elevation value after the expansion of the current grid, and making the grid minimum The elevation value is set to the inflated elevation value;
(4)遍历每个网格,求取网格内每个点高程值与膨胀后的高程值之差,当其绝对值小于设定阈值时,该点为地面点,否则为非地面点。(4) Traversing each grid, calculating the difference between the elevation value of each point in the grid and the inflated elevation value, when its absolute value is less than the set threshold, the point is a ground point, otherwise it is a non-ground point.
所述第二步包括:The second step includes:
(1)针对每个点pi(xi,yi)由公式4计算其所在网格;(1) For each point p i (x i , y i ), calculate its grid by formula 4;
(2)计算搜索半径,
(3)设点pi(xi,yi)的矩形邻域内的点集为(3) Let the point set in the rectangular neighborhood of point p i (xi , y i ) be
Pi_near={pj(xj,yj)|xj∈(xi-r,xi+r),yj∈(yi-r,yi+r)}(其中j是点的点号或者说是下标,并与当前点pi相区别),由此以当前点pi为中心,搜索周围网格内的点集Pi_near;P i_near ={p j (x j ,y j )|x j ∈(xi -r , xi +r),y j ∈(y i -r,y i +r)} (where j is the point The point number or subscript is different from the current point pi), so the point set P i_near in the surrounding grid is searched with the current point pi as the center;
(4)计算Pi_near中每个点与pi(xi,yi)的高程差绝对值之和ΣHi,若ΣHi小于阈值,则当前点为道路候选点;(4) Calculate the sum ΣH i of the absolute value of the elevation difference between each point in P i_near and p i (xi , y i ), if ΣH i is less than the threshold, the current point is a road candidate point;
zj是Pi_near中第j个点的高程值,zi是pi(xi,yi)的高程值,n为Pi_near的总点数;zj is the elevation value of the jth point in P i_near , zi is the elevation value of p i ( xi , y i ), and n is the total number of points in P i_near ;
(5)重复以上步骤,直至完成对所有点的判断,且将保留的道路候选点集为Pc。(5) Repeat the above steps until the judgment of all points is completed, and the reserved road candidate points are set as P c .
所述第三步包括:The third step includes:
(1)对每个道路候选点pi,根据公式7计算矩形范围点pi,r1,pi,r2,pi,r3,pi,r4,即矩形的四个顶点,pi,r1与当前点重合,其余点按逆时针排列,其中dθ为矩形的旋转角度,将360度划分为n个间隔;(1) For each road candidate point p i , calculate the rectangular range points p i,r1 , p i,r2 , p i,r3 , p i,r4 according to formula 7, that is, the four vertices of the rectangle, p i,r1 It coincides with the current point, and the other points are arranged counterclockwise, where dθ is the rotation angle of the rectangle, and divides 360 degrees into n intervals;
k=0,1,2…18 (公式7)k=0,1,2...18 (Formula 7)
其中pi,r1(x)为pi,r1的x坐标值,pi,r1(y)为pi,r1的y坐标值,其余的同理,l为矩形长度,w为矩形宽度,dθ为矩形的旋转角度间隔,k为旋转的次数或称为迭代次数;Among them, p i, r1 (x) is the x coordinate value of p i, r1 , p i, r1 (y) is the y coordinate value of p i, r1 , and the rest are the same, l is the length of the rectangle, w is the width of the rectangle, dθ is the rotation angle interval of the rectangle, and k is the number of rotations or the number of iterations;
(2)根据第4步的邻域判断方法,搜索pi点的邻近点Pi_near,然后判断Pi_near是否在pi,r1,pi,r2,pi,r3,pi,r4构成的矩形范围内,将在范围内的保存为点集Pi_r,k,并计算其高程方差 (2) According to the neighborhood judgment method in step 4, search for the neighboring point P i_near of point p i , and then judge whether P i_near is formed by p i,r1 , p i,r2 , p i,r3 , p i,r4 Within the rectangular range, save the point set P i_r,k within the range, and calculate its elevation variance
高程均值:其中zj为点集Pi_r,k中第j个点的高程值,m为点集Pi_r,k点个数;Elevation mean: Where z j is the elevation value of the jth point in the point set P i_r,k , and m is the number of points in the point set P i_r,k ;
高程方差值:
(3)根据以上两步,计算点pi的所有矩形(k=0,1,2…18)邻域内点的高程方差值,得到高程方差值集合分别比较集合Ψ中的高程方差值,记录高程方差值最小、第二小和第三小的三个值和和对应的矩形编号k1、k2、k3,若同时满足以下两个方程则当前点pi记为道路中心点;(3) According to the above two steps, calculate the elevation variance value of points in the neighborhood of all rectangles (k=0,1,2...18) of point p i , and obtain the elevation variance value set Compare the elevation variance values in the set Ψ respectively, and record the three values with the smallest, second smallest, and third smallest elevation variance values and and the corresponding rectangle numbers k 1 , k 2 , k 3 , if the following two equations are satisfied at the same time, the current point p i is recorded as the road center point;
式中ψ为用户设定道路粗糙阈值,μ为用户设定道路转弯阈值;where ψ is the road roughness threshold set by the user, and μ is the road turning threshold set by the user;
(4)重复以上步骤,直至完成对所有点的判断,并记道路中心点集为Pk。(4) Repeat the above steps until the judgment of all points is completed, and record the road center point set as P k .
所述第四步包括:The fourth step includes:
根据计算得到的道路候选点集Pc和道路中心点集Pk,道路候选点集中可能含有大量停车场、运动场、屋顶等平坦地域的点,通过以下的方法将这些点剔除:According to the calculated road candidate point set P c and road center point set P k , the road candidate point set may contain a large number of points in flat areas such as parking lots, sports fields, and roofs. These points are eliminated by the following methods:
(1)将点集Pc和Pk进行规格网格划分,且用编码1和2表示点的种类为道路候选点还是道路中心点;(1) Divide the point sets P c and P k into standard grids, and use codes 1 and 2 to indicate whether the point type is a road candidate point or a road center point;
(2)对点集Pc中的每个点,搜索其周围r/2范围内属性为2的点,若有则保留该点,否则就判断为非道路点,参数r/2表示道路宽度的二分之一与网格间距的比值。(2) For each point in the point set Pc , search for a point with an attribute of 2 within the r/2 range around it, if there is, keep the point, otherwise it will be judged as a non-road point, and the parameter r/2 represents the width of the road The ratio of one-half of to the grid spacing.
本发明优点是:Advantage of the present invention is:
本发明主要解决在没有反射强度或者反射强度不能较好反应地物属性,或者没有航空影像等情况,根据道路拓扑结构特点与周围环境的不同进行道路提取。只采用点云数据进行丘林地域的道路提取,较少了对原始数据的依赖性。The present invention mainly solves the problem of extracting roads according to the differences between road topological structure characteristics and the surrounding environment when there is no reflection intensity or the reflection intensity cannot better reflect the ground object attributes, or there is no aerial image. Only point cloud data is used for road extraction in hilly areas, which reduces the dependence on original data.
附图说明Description of drawings
下面结合附图及实施例对本发明作进一步描述:The present invention will be further described below in conjunction with accompanying drawing and embodiment:
图1为本发明将机载LiDAR数据分配到各个网格的划分示意图;Fig. 1 is the division schematic diagram that airborne LiDAR data is distributed to each grid by the present invention;
图2为本发明的w×w个网格示意图;Fig. 2 is a schematic diagram of w×w grids of the present invention;
图3为本发明的规则网格划分示意图;Fig. 3 is a schematic diagram of regular grid division of the present invention;
图4为本发明的邻域点集示意图;Fig. 4 is a schematic diagram of a neighborhood point set of the present invention;
图5为本发明的矩形范围点示意图;Fig. 5 is a schematic diagram of rectangular range points of the present invention;
图6为本发明的旋转矩形示意图。Fig. 6 is a schematic diagram of a rotating rectangle in the present invention.
具体实施方式Detailed ways
实施例:参考图1-6所示,本发明提供了一种基于矩形邻域分析的道路数据提取方法,尤其是一种基于矩形邻域分析的机载LiDAR丘林地域道路提取方法,总体框架如下:Embodiment: Referring to Figures 1-6, the present invention provides a method for extracting road data based on rectangular neighborhood analysis, especially an airborne LiDAR road extraction method based on rectangular neighborhood analysis, the overall framework as follows:
第一步:数据预处理,先计算点云间距,建立规则网格组织,然后采用基于规则网格的快速数学形态学滤波方法进行滤波,剔除非地面点;The first step: data preprocessing, first calculate the point cloud spacing, establish a regular grid organization, and then use a fast mathematical morphology filtering method based on a regular grid to filter out non-ground points;
第二步:筛选道路候选点,初步确定道路点集,缩小判断范围,对每个点搜索其邻域范围内的点,并比较其高程值,若都在阈值范围内,则作为道路候选点。Step 2: Screen candidate road points, preliminarily determine the road point set, narrow the judgment range, search for points within its neighborhood for each point, and compare their elevation values. If they are all within the threshold range, they will be used as road candidate points .
第三步:提取道路中心线,对每个道路候选点,计算不同方位角的矩形范围内高程方差值,根据方差值间的差异判断当前点是否为道路中心点。Step 3: extract the road centerline, calculate the elevation variance value within the rectangular range of different azimuths for each road candidate point, and judge whether the current point is the road center point according to the difference between the variance values.
第四步:非道路点剔除。Step 4: Eliminate non-road points.
具体详细步骤如下:The detailed steps are as follows:
1、根据机载LiDAR数据范围和数量估算点云间距pt_dis,其中机载LiDAR数据范围是指点云所覆盖的范围,一般用四个角的坐标表示;1. Estimate the point cloud distance pt_dis according to the range and quantity of the airborne LiDAR data, where the airborne LiDAR data range refers to the range covered by the point cloud, which is generally represented by the coordinates of the four corners;
(1)首先通过机载LiDAR读取机载LiDAR数据范围(xmin,xmax,ymin,ymax),统计点的数量pt_num;(1) First read the airborne LiDAR data range (xmin, xmax, ymin, ymax) through the airborne LiDAR, and the number of statistical points pt_num;
(2)计算点云间距pt_dis:(2) Calculate the point cloud distance pt_dis:
2、对机载LiDAR数据进行规则网格划分,便于搜索邻域点,具体步骤及说明如下:2. Carry out regular grid division on the airborne LiDAR data to facilitate the search for neighboring points. The specific steps and instructions are as follows:
(1)根据点云间距pt_dis和计算需要,建立覆盖点云范围(xmin,xmax,ymin,ymax)的规则网格,本发明取网格间距d优选为点云间距pt_dis的2~3倍。(1) According to the point cloud spacing pt_dis and calculation requirements, a regular grid covering the point cloud range (xmin, xmax, ymin, ymax) is established. The grid spacing d is preferably 2 to 3 times the point cloud spacing pt_dis in the present invention.
(2)根据点的坐标将机载LiDAR数据分配到各个网格中。(2) Distribute the airborne LiDAR data into individual grids according to the coordinates of the points.
GridNum=XAxisNum*YAxisNum(公式3)GridNum=XAxisNum*YAxisNum (Formula 3)
XAxisNum、YAxisNum为x、y方向的网格数,GridNum为网格总数。XAxisNum and YAxisNum are the number of grids in the x and y directions, and GridNum is the total number of grids.
机载LiDAR数据中的任意一点pi(xi,yi)所在的网格为:The grid of any point p i ( xi , y i ) in the airborne LiDAR data is:
其中Xi为Pi所在的X方向网格数、Yi为Pi所在的Y方向网格数、Gi为Pi所在的一维数组记录的网格数,如图1所示,图中的数据计算为XAxisNum=11,YAxisNum=7,GridNum=77,Xi=7,Yi=4,Gi=4*11+7=51,表示Pi点在第51个网格内。Among them, Xi is the number of grids in the X direction where Pi is located, Yi is the number of grids in the Y direction where Pi is located, and Gi is the number of grids recorded in the one-dimensional array where Pi is located. As shown in Figure 1, the data in the figure is calculated as XAxisNum =11, YAxisNum=7, GridNum=77, Xi=7, Yi=4, Gi=4*11+7=51, indicating that the Pi point is in the 51st grid.
3、通过数学形态学的基本运算及其扩展运算对机载LiDAR数据进行非地面滤波,保留包含道路点在内的地面点数据。3. Perform non-ground filtering on the airborne LiDAR data through the basic operation of mathematical morphology and its extended operation, and retain the ground point data including road points.
以规则网格为单元的数学形态学滤波,网格腐蚀结果为w×w网格内的最小高程值,膨胀结果为w×w网格内最大高程值,其中w代表窗口的大小,也就是网格大小,如图2所示,其为w=2的示意图。Mathematical morphological filtering with a regular grid as the unit, the grid erosion result is the minimum elevation value in the w×w grid, and the dilation result is the maximum elevation value in the w×w grid, where w represents the size of the window, that is, Grid size, as shown in Figure 2, which is a schematic diagram of w=2.
其中公式5是典型的规则数学形态学公式,A代表目标集合,B代表结构元素,AΘB表示腐蚀,表示膨胀,式中zi′和zi′′分别为网格单元腐蚀和膨胀后的高程值,Gridi(zmin)为网格内高程最小值,w为窗口的大小。Among them, Formula 5 is a typical regular mathematical morphology formula, A represents the target set, B represents the structural element, AΘB represents corrosion, In the formula, zi ′ and zi ′ ′ are the elevation values of the grid unit after corrosion and expansion respectively, Grid i (z min ) is the minimum elevation value in the grid, and w is the size of the window.
具体步骤如下:Specific steps are as follows:
(1)由网格间距计算网格数量,并为网格分配内存空间,将点逐个按平面坐标存储到相应的网格单元中,同时与该网格单元中已保存的最小高程值进行比较(初始默认最小高程值为一个很大的值),记录比较后的最小高程值;(1) Calculate the number of grids from the grid spacing, and allocate memory space for the grid, store the points one by one in the corresponding grid unit according to the plane coordinates, and compare with the minimum elevation value saved in the grid unit (The initial default minimum elevation value is a very large value), record the minimum elevation value after comparison;
(2)遍历每个网格,比较周围w×w个网格内最小高程值,采用公式5的第一个公式,取最小值为当前网格腐蚀后的高程值,并将该网格最大高程值设为腐蚀后的高程值,具体是将周围一定范围内网格内的最小高程值进行比较,保留最小的高程值,并将其赋值给当前网格的最大高程值这一变量;(2) Traverse each grid, compare the minimum elevation value in the surrounding w×w grids, use the first formula of formula 5, take the minimum value as the elevation value of the current grid after corrosion, and set the maximum value of the grid The elevation value is set to the elevation value after corrosion, specifically comparing the minimum elevation value in the surrounding grid within a certain range, retaining the minimum elevation value, and assigning it to the variable of the maximum elevation value of the current grid;
(3)遍历每个网格,比较周围w×w个网格内腐蚀高程值,采用公式5的第二个公式,取最大值为当前网格膨胀后的高程值,并将该网格最小高程值设为膨胀后的高程值,具体是将周围一定范围内网格内的腐蚀高程值(上一步计算的结果)进行比较,保留最小的高程值,并将其赋值给当前网格的最小高程值这一变量;(3) Traversing through each grid, comparing the corrosion elevation values in the surrounding w×w grids, using the second formula of formula 5, taking the maximum value as the elevation value after the expansion of the current grid, and making the grid minimum The elevation value is set to the elevation value after expansion. Specifically, compare the corrosion elevation values (results calculated in the previous step) in the grid within a certain range around, keep the minimum elevation value, and assign it to the minimum value of the current grid. The variable of elevation value;
(4)遍历每个网格,求取网格内每个点高程值与膨胀后的高程值(即该网格内最小高程值)之差,当其绝对值小于设定阈值时,该点为地面点,否则为非地面点。(4) Traversing each grid, calculating the difference between the elevation value of each point in the grid and the height value after expansion (that is, the minimum elevation value in the grid), when its absolute value is less than the set threshold, the point is a ground point, otherwise it is a non-ground point.
经过以上四步处理,即可快速完成对机载LiDAR数据的滤波处理。After the above four steps of processing, the filtering processing of the airborne LiDAR data can be quickly completed.
4、估计道路的宽度road_width,选取邻域内道路候选点。4. Estimate the width road_width of the road, and select the road candidate points in the neighborhood.
(1)针对每个点pi(xi,yi)由公式4计算其所在网格;(1) For each point p i (x i , y i ), calculate its grid by formula 4;
(2)计算搜索半径,
(3)设点pi(xi,yi)的矩形邻域内的点集为(3) Let the point set in the rectangular neighborhood of point p i (xi , y i ) be
Pi_near={pj(xj,yj)|xj∈(xi-r,xi+r),yj∈(yi-r,yi+r)}(其中j是点的点号或者说是下标,为了区别于当前点pi的下标i),由此以当前点pi为中心,搜索周围网格内的点集Pi_near,如图4所示,带有箭头指向的点为pi(xi,yi),R=3时,Pi_near点集为标出的点。P i_near ={p j (x j ,y j )|x j ∈(xi -r , xi +r),y j ∈(y i -r,y i +r)} (where j is the point The point number or subscript, in order to distinguish it from the subscript i) of the current point pi, so that the current point pi is used as the center to search for the point set P i_near in the surrounding grid, as shown in Figure 4, with arrows pointing to The point of is p i ( xi , y i ), when R=3, the point set of P i_near is the marked point.
(4)计算Pi_near中每个点与pi(xi,yi)的高程差绝对值之和ΣHi,若ΣHi小于阈值,则当前点为道路候选点。(4) Calculate the sum ΣH i of the absolute value of the elevation difference between each point in P i_near and p i (xi , y i ), if ΣH i is smaller than the threshold, the current point is a road candidate point.
zj是Pi_near中第j个点的高程值,zi是pi(xi,yi)的高程值,n为Pi_near的总点数。zj is the elevation value of the jth point in P i_near , zi is the elevation value of p i ( xi , y i ), and n is the total number of points in P i_near .
(5)重复以上步骤,直至完成对所有点的判断,且将保留的道路候选点集记为Pc。(5) Repeat the above steps until the judgment of all points is completed, and record the reserved road candidate point set as P c .
5、搜索道路中心线5. Search for road centerline
对每个点云,建立矩形搜索邻域,矩形的长度l优选为宽度w的2~3倍。For each point cloud, a rectangular search neighborhood is established, and the length l of the rectangle is preferably 2 to 3 times the width w.
(1)对每个道路候选点pi,根据公式7计算矩形范围点pi,r1,pi,r2,pi,r3,pi,r4,即矩形的四个顶点,pi,r1与当前点重合,其余点按逆时针排列,其中dθ为矩形的旋转角度,将360度划分为n个间隔,本发明取dθ=20度,即n=18,如图4所示。(1) For each road candidate point p i , calculate the rectangular range points p i,r1 , p i,r2 , p i,r3 , p i,r4 according to formula 7, that is, the four vertices of the rectangle, p i,r1 Coincident with the current point, all the other points are arranged counterclockwise, wherein dθ is the rotation angle of the rectangle, and 360 degrees are divided into n intervals. The present invention takes dθ=20 degrees, i.e. n=18, as shown in Figure 4.
k=0,1,2…18 (公式7)k=0,1,2...18 (Formula 7)
其中pi,r1(x)为pi,r1的x坐标值,pi,r1(y)为pi,r1的y坐标值,其余的同理,l为矩形长度,w为矩形宽度,dθ为矩形的旋转角度间隔,k为旋转的次数或称为迭代次数。Among them, p i, r1 (x) is the x coordinate value of p i, r1 , p i, r1 (y) is the y coordinate value of p i, r1 , and the rest are the same, l is the length of the rectangle, w is the width of the rectangle, dθ is the rotation angle interval of the rectangle, and k is the number of rotations or the number of iterations.
(2)根据第4步的邻域判断方法,搜索pi点的邻近点Pi_near,然后判断Pi_near是否在pi,r1,pi,r2,pi,r3,pi,r4构成的矩形范围内,将在范围内的保存为点集Pi_r,k,并计算其高程方差 (2) According to the neighborhood judgment method in step 4, search for the neighboring point P i_near of point p i , and then judge whether P i_near is formed by p i,r1 , p i,r2 , p i,r3 , p i,r4 Within the rectangular range, save the point set P i_r,k within the range, and calculate its elevation variance
高程均值:其中zj为点集Pi_r,k中第j个点的高程值,m为点集Pi_r,k点个数;Elevation mean: Where z j is the elevation value of the jth point in the point set P i_r,k , and m is the number of points in the point set P i_r,k ;
高程方差值:
(3)根据以上两步,计算点pi的所有矩形(k=0,1,2…18)邻域内点的高程方差值,得到高程方差值集合分别比较集合Ψ中的高程方差值,记录高程方差值最小、第二小和第三小的三个值和和对应的矩形编号k1、k2、k3,若同时满足以下两个方程则当前点pi记为道路中心点。(3) According to the above two steps, calculate the elevation variance value of points in the neighborhood of all rectangles (k=0,1,2...18) of point p i , and obtain the elevation variance value set Compare the elevation variance values in the set Ψ respectively, and record the three values with the smallest, second smallest, and third smallest elevation variance values and and the corresponding rectangle numbers k 1 , k 2 , k 3 , if the following two equations are satisfied at the same time, the current point p i is recorded as the road center point.
式中ψ为用户设定道路粗糙阈值,本发明取为0.5,μ为用户设定道路转弯阈值,本发明取为1,即转弯角度最大为20度。In the formula, ψ is the road roughness threshold set by the user, which is taken as 0.5 in the present invention, and μ is the road turning threshold set by the user, which is taken as 1 in the present invention, that is, the maximum turning angle is 20 degrees.
(4)重复以上步骤,直至完成对所有点的判断,并记道路中心点集为Pk。(4) Repeat the above steps until the judgment of all points is completed, and record the road center point set as P k .
6、非道路点剔除6. Elimination of non-road points
根据第4、5步分别计算得到的道路候选点集Pc和道路中心点集Pk,道路候选点集中可能含有大量停车场、运动场、屋顶等平坦地域的点,通过以下的方法将这些点剔除。According to the road candidate point set P c and the road center point set P k calculated in steps 4 and 5 respectively, the road candidate point set may contain a large number of points in flat areas such as parking lots, sports fields, and roofs. These points are divided by the following method remove.
(1)采用与第2步类似的方法,将点集Pc和Pk进行规格网格划分,且用编码1和2表示点的种类为道路候选点还是道路中心点。(1) Using a method similar to step 2, divide the point sets P c and P k into standard grids, and use codes 1 and 2 to indicate whether the point type is a road candidate point or a road center point.
(2)对点集点集Pc中的每个点,搜索其周围r/2范围内属性为2的点,若有则保留该点,否则就判断为非道路点。参数r/2表示道路宽度的二分之一与网格间距的比值,本发明中r/2取为2。(2) For each point in the point set Pc , search for points with an attribute of 2 within the r/2 range around it, if there is, keep the point, otherwise it is judged as a non-road point. The parameter r/2 represents the ratio of 1/2 of the road width to the grid spacing, and r/2 is taken as 2 in the present invention.
当然上述实施例只为说明本发明的技术构思及特点,其目的在于让熟悉此项技术的人能够了解本发明的内容并据以实施,并不能以此限制本发明的保护范围。凡根据本发明主要技术方案的精神实质所做的等效变换或修饰,都应涵盖在本发明的保护范围之内。Of course, the above-mentioned embodiments are only to illustrate the technical conception and characteristics of the present invention, and its purpose is to enable those skilled in the art to understand the content of the present invention and implement it accordingly, and not to limit the protection scope of the present invention. All equivalent changes or modifications made according to the spirit of the main technical solutions of the present invention shall fall within the protection scope of the present invention.
Claims (6)
Priority Applications (1)
Application Number | Priority Date | Filing Date | Title |
---|---|---|---|
CN201410213465.5A CN104050473B (en) | 2014-05-20 | 2014-05-20 | A Road Data Extraction Method Based on Rectangular Neighborhood Analysis |
Applications Claiming Priority (1)
Application Number | Priority Date | Filing Date | Title |
---|---|---|---|
CN201410213465.5A CN104050473B (en) | 2014-05-20 | 2014-05-20 | A Road Data Extraction Method Based on Rectangular Neighborhood Analysis |
Publications (2)
Publication Number | Publication Date |
---|---|
CN104050473A true CN104050473A (en) | 2014-09-17 |
CN104050473B CN104050473B (en) | 2019-07-19 |
Family
ID=51503285
Family Applications (1)
Application Number | Title | Priority Date | Filing Date |
---|---|---|---|
CN201410213465.5A Expired - Fee Related CN104050473B (en) | 2014-05-20 | 2014-05-20 | A Road Data Extraction Method Based on Rectangular Neighborhood Analysis |
Country Status (1)
Country | Link |
---|---|
CN (1) | CN104050473B (en) |
Cited By (8)
Publication number | Priority date | Publication date | Assignee | Title |
---|---|---|---|---|
CN105893961A (en) * | 2016-03-30 | 2016-08-24 | 广东中冶地理信息股份有限公司 | Method for extracting road center line |
CN106709892A (en) * | 2016-11-15 | 2017-05-24 | 昂纳自动化技术(深圳)有限公司 | Rapid region expansion algorithm and device of any structure element based on stroke coding |
CN107657636A (en) * | 2017-10-16 | 2018-02-02 | 南京市测绘勘察研究院股份有限公司 | A kind of method that route topography figure elevational point is automatically extracted based on mobile lidar data |
CN108021844A (en) * | 2016-10-31 | 2018-05-11 | 高德软件有限公司 | A kind of road edge recognition methods and device |
CN111158015A (en) * | 2019-12-31 | 2020-05-15 | 飞燕航空遥感技术有限公司 | Detection method and system for point cloud data of airborne laser radar to be wrongly divided into ground points |
CN111783721A (en) * | 2020-07-13 | 2020-10-16 | 湖北亿咖通科技有限公司 | Lane line extraction method of laser point cloud and electronic equipment |
CN111783722A (en) * | 2020-07-13 | 2020-10-16 | 湖北亿咖通科技有限公司 | Lane line extraction method of laser point cloud and electronic equipment |
CN111932477A (en) * | 2020-08-07 | 2020-11-13 | 武汉中海庭数据技术有限公司 | Noise removal method and device based on single line laser radar point cloud |
Citations (5)
Publication number | Priority date | Publication date | Assignee | Title |
---|---|---|---|---|
US20030103649A1 (en) * | 2001-11-30 | 2003-06-05 | Nissan Motor Co., Ltd. | Road white line recognition apparatus and method |
CN101901343A (en) * | 2010-07-20 | 2010-12-01 | 同济大学 | Road Extraction Method Based on Stereo Constraint from Remote Sensing Image |
CN103500338A (en) * | 2013-10-16 | 2014-01-08 | 厦门大学 | Road zebra crossing automatic extraction method based on vehicle-mounted laser scanning point cloud |
CN103605135A (en) * | 2013-11-26 | 2014-02-26 | 中交第二公路勘察设计研究院有限公司 | Road feature extracting method based on fracture surface subdivision |
CN103714339A (en) * | 2013-12-30 | 2014-04-09 | 武汉大学 | SAR image road damaging information extracting method based on vector data |
-
2014
- 2014-05-20 CN CN201410213465.5A patent/CN104050473B/en not_active Expired - Fee Related
Patent Citations (5)
Publication number | Priority date | Publication date | Assignee | Title |
---|---|---|---|---|
US20030103649A1 (en) * | 2001-11-30 | 2003-06-05 | Nissan Motor Co., Ltd. | Road white line recognition apparatus and method |
CN101901343A (en) * | 2010-07-20 | 2010-12-01 | 同济大学 | Road Extraction Method Based on Stereo Constraint from Remote Sensing Image |
CN103500338A (en) * | 2013-10-16 | 2014-01-08 | 厦门大学 | Road zebra crossing automatic extraction method based on vehicle-mounted laser scanning point cloud |
CN103605135A (en) * | 2013-11-26 | 2014-02-26 | 中交第二公路勘察设计研究院有限公司 | Road feature extracting method based on fracture surface subdivision |
CN103714339A (en) * | 2013-12-30 | 2014-04-09 | 武汉大学 | SAR image road damaging information extracting method based on vector data |
Non-Patent Citations (2)
Title |
---|
王勇: ""地面激光扫描数据滤波研究"", 《中国优秀硕士学位论文全文数据库 基础科学辑》 * |
陈卓,马洪超,李云帆: ""结合角度纹理信息和Snake方法从LiDAR点云数据中提取道路交叉口"", 《国土资源遥感》 * |
Cited By (10)
Publication number | Priority date | Publication date | Assignee | Title |
---|---|---|---|---|
CN105893961A (en) * | 2016-03-30 | 2016-08-24 | 广东中冶地理信息股份有限公司 | Method for extracting road center line |
CN108021844A (en) * | 2016-10-31 | 2018-05-11 | 高德软件有限公司 | A kind of road edge recognition methods and device |
CN108021844B (en) * | 2016-10-31 | 2020-06-02 | 阿里巴巴(中国)有限公司 | Road edge identification method and device |
CN106709892A (en) * | 2016-11-15 | 2017-05-24 | 昂纳自动化技术(深圳)有限公司 | Rapid region expansion algorithm and device of any structure element based on stroke coding |
CN107657636A (en) * | 2017-10-16 | 2018-02-02 | 南京市测绘勘察研究院股份有限公司 | A kind of method that route topography figure elevational point is automatically extracted based on mobile lidar data |
CN111158015A (en) * | 2019-12-31 | 2020-05-15 | 飞燕航空遥感技术有限公司 | Detection method and system for point cloud data of airborne laser radar to be wrongly divided into ground points |
CN111158015B (en) * | 2019-12-31 | 2020-11-24 | 飞燕航空遥感技术有限公司 | Detection method and system for point cloud data of airborne laser radar to be wrongly divided into ground points |
CN111783721A (en) * | 2020-07-13 | 2020-10-16 | 湖北亿咖通科技有限公司 | Lane line extraction method of laser point cloud and electronic equipment |
CN111783722A (en) * | 2020-07-13 | 2020-10-16 | 湖北亿咖通科技有限公司 | Lane line extraction method of laser point cloud and electronic equipment |
CN111932477A (en) * | 2020-08-07 | 2020-11-13 | 武汉中海庭数据技术有限公司 | Noise removal method and device based on single line laser radar point cloud |
Also Published As
Publication number | Publication date |
---|---|
CN104050473B (en) | 2019-07-19 |
Similar Documents
Publication | Publication Date | Title |
---|---|---|
CN104050473A (en) | Road data extraction method based on rectangular neighborhood analysis | |
Yan et al. | Urban land cover classification using airborne LiDAR data: A review | |
Höfle et al. | Urban vegetation detection using radiometrically calibrated small-footprint full-waveform airborne LiDAR data | |
Goodwin et al. | Characterizing urban surface cover and structure with airborne lidar technology | |
CN114877838A (en) | A road geometric feature detection method based on vehicle laser scanning system | |
CN107798657A (en) | A kind of vehicle-mounted laser point cloud filtering method based on circular cylindrical coordinate | |
Yadav et al. | Identification of trees and their trunks from mobile laser scanning data of roadway scenes | |
Karsli et al. | Automatic building extraction from very high-resolution image and LiDAR data with SVM algorithm | |
CN111611900A (en) | Target point cloud identification method and device, electronic equipment and storage medium | |
CN110207676A (en) | The acquisition methods and device of a kind of field ditch pool parameter | |
Jiangui et al. | A method for main road extraction from airborne LiDAR data in urban area | |
RU2638638C1 (en) | Method and system of automatic constructing three-dimensional models of cities | |
Rau et al. | Semi-automatic shallow landslide detection by the integration of airborne imagery and laser scanning data | |
Meng et al. | Canopy structure attributes extraction from LiDAR data based on tree morphology and crown height proportion | |
CN105701856B (en) | A kind of vegetation extracting method and system | |
CN109522787B (en) | A method of small road recognition based on remote sensing data | |
CN117523522B (en) | Tunnel scene identification method, device, equipment and storage medium | |
CN110618145B (en) | A method for rapid determination of the location of springs in the Loess Plateau based on unmanned aerial vehicles | |
CN114004740B (en) | Building wall line extraction method based on unmanned aerial vehicle laser radar point cloud | |
Ma et al. | Road Curbs Extraction from Mobile Laser Scanning Point Clouds with Multidimensional Rotation‐Invariant Version of the Local Binary Pattern Features | |
CN114092419A (en) | Automatic inspection method for point cloud spatial position quality based on earth surface point location | |
Gao et al. | Experimental Study on Precise Recognition of Settlements in Mountainous Areas Based on UAV Image and LIDAR Point Cloud | |
Hujebri et al. | Automatic building extraction from lidar point cloud data in the fusion of orthoimage | |
Vizireanu et al. | The potential of airborne LiDAR for detection of new archaeological site in Romania | |
CN117437537B (en) | Building target level change detection method and system based on airborne LiDAR point cloud data |
Legal Events
Date | Code | Title | Description |
---|---|---|---|
C06 | Publication | ||
PB01 | Publication | ||
C10 | Entry into substantive examination | ||
SE01 | Entry into force of request for substantive examination | ||
GR01 | Patent grant | ||
GR01 | Patent grant | ||
CF01 | Termination of patent right due to non-payment of annual fee | ||
CF01 | Termination of patent right due to non-payment of annual fee |
Granted publication date: 20190719 Termination date: 20200520 |