[go: up one dir, main page]

CN115222914A - Aggregate three-dimensional morphology spherical harmonic reconstruction method based on three-dimensional point cloud data - Google Patents

Aggregate three-dimensional morphology spherical harmonic reconstruction method based on three-dimensional point cloud data Download PDF

Info

Publication number
CN115222914A
CN115222914A CN202210722175.8A CN202210722175A CN115222914A CN 115222914 A CN115222914 A CN 115222914A CN 202210722175 A CN202210722175 A CN 202210722175A CN 115222914 A CN115222914 A CN 115222914A
Authority
CN
China
Prior art keywords
dimensional
aggregate
point cloud
reconstruction
spherical harmonic
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.)
Pending
Application number
CN202210722175.8A
Other languages
Chinese (zh)
Inventor
仇培云
陈昊
韦江雄
高鹏
吴娟
王大海
崔力波
王江卡
余大鹏
王鹏
Current Assignee (The listed assignees may be inaccurate. Google has not performed a legal analysis and makes no representation or warranty as to the accuracy of the list.)
South China University of Technology SCUT
Guangzhou Metro Group Co Ltd
China Railway Construction South China Construction Co Ltd
China Railway Guangzhou Investment and Development Co Ltd
Original Assignee
South China University of Technology SCUT
Guangzhou Metro Group Co Ltd
China Railway Construction South China Construction Co Ltd
China Railway Guangzhou Investment and Development Co Ltd
Priority date (The priority date 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 date listed.)
Filing date
Publication date
Application filed by South China University of Technology SCUT, Guangzhou Metro Group Co Ltd, China Railway Construction South China Construction Co Ltd, China Railway Guangzhou Investment and Development Co Ltd filed Critical South China University of Technology SCUT
Priority to CN202210722175.8A priority Critical patent/CN115222914A/en
Publication of CN115222914A publication Critical patent/CN115222914A/en
Pending legal-status Critical Current

Links

Images

Classifications

    • GPHYSICS
    • G06COMPUTING; CALCULATING OR COUNTING
    • G06TIMAGE DATA PROCESSING OR GENERATION, IN GENERAL
    • G06T17/00Three dimensional [3D] modelling, e.g. data description of 3D objects
    • G06T17/20Finite element generation, e.g. wire-frame surface description, tesselation
    • GPHYSICS
    • G06COMPUTING; CALCULATING OR COUNTING
    • G06FELECTRIC DIGITAL DATA PROCESSING
    • G06F3/00Input arrangements for transferring data to be processed into a form capable of being handled by the computer; Output arrangements for transferring data from processing unit to output unit, e.g. interface arrangements
    • G06F3/06Digital input from, or digital output to, record carriers, e.g. RAID, emulated record carriers or networked record carriers
    • G06F3/0601Interfaces specially adapted for storage systems
    • G06F3/0602Interfaces specially adapted for storage systems specifically adapted to achieve a particular effect
    • G06F3/0608Saving storage space on storage systems
    • GPHYSICS
    • G06COMPUTING; CALCULATING OR COUNTING
    • G06TIMAGE DATA PROCESSING OR GENERATION, IN GENERAL
    • G06T7/00Image analysis
    • G06T7/70Determining position or orientation of objects or cameras
    • G06T7/73Determining position or orientation of objects or cameras using feature-based methods
    • GPHYSICS
    • G06COMPUTING; CALCULATING OR COUNTING
    • G06TIMAGE DATA PROCESSING OR GENERATION, IN GENERAL
    • G06T2207/00Indexing scheme for image analysis or image enhancement
    • G06T2207/10Image acquisition modality
    • G06T2207/10028Range image; Depth image; 3D point clouds
    • GPHYSICS
    • G06COMPUTING; CALCULATING OR COUNTING
    • G06TIMAGE DATA PROCESSING OR GENERATION, IN GENERAL
    • G06T2207/00Indexing scheme for image analysis or image enhancement
    • G06T2207/10Image acquisition modality
    • G06T2207/10072Tomographic images
    • G06T2207/10081Computed x-ray tomography [CT]

Landscapes

  • Engineering & Computer Science (AREA)
  • Physics & Mathematics (AREA)
  • Theoretical Computer Science (AREA)
  • General Physics & Mathematics (AREA)
  • Computer Graphics (AREA)
  • Geometry (AREA)
  • Software Systems (AREA)
  • Computer Vision & Pattern Recognition (AREA)
  • Human Computer Interaction (AREA)
  • General Engineering & Computer Science (AREA)
  • Length Measuring Devices With Unspecified Measuring Means (AREA)

Abstract

The invention discloses an aggregate three-dimensional shape spherical harmonic reconstruction method based on three-dimensional point cloud data, which comprises the following steps: extracting three-dimensional point cloud coordinates of aggregate particles; surface parameterization, namely converting the aggregate three-dimensional point cloud data into a spherical coordinate system from a Cartesian coordinate system, and obtaining the radius of a reconstruction point by applying an angle approximation algorithm according to the regular increment of a polar angle and an azimuth angle of a Gaussian orthogonal point; calculating a reconstructed point-sphere function based on a Legendre polynomial in combination with a recursive algorithm, and solving a spherical harmonic coefficient by using the sphere function; and summing the spherical harmonic coefficients and the generalized Fourier series on the spherical surface to obtain a surface function reflecting the real morphology of the aggregate, and carrying out visual processing on the reconstructed morphology. According to the invention, the three-dimensional point cloud data of the aggregate particles are converted into spherical harmonic coefficient description, so that the calculation of surface morphology parameters such as sphericity and edge angle of the aggregate particles is greatly facilitated; meanwhile, the storage mode of the three-dimensional point cloud data is converted into the storage mode of the spherical harmonic coefficient, the storage compression ratio of the three-dimensional point cloud data can reach 1%, and a large amount of storage space is saved.

Description

一种基于三维点云数据的骨料三维形貌球谐重构方法A spherical harmonic reconstruction method of aggregate 3D morphology based on 3D point cloud data

技术领域technical field

本发明属于骨料三维形貌重构技术领域,具体涉及一种基于三维点云数据的骨料三维形貌球谐重构方法。The invention belongs to the technical field of three-dimensional topography reconstruction of aggregates, and in particular relates to a spherical harmonic reconstruction method of three-dimensional topography of aggregates based on three-dimensional point cloud data.

背景技术Background technique

骨料在混凝土中主要起骨架和填充作用,同时可以减少混凝土硬化过程中因湿胀干缩引起的体积变化,骨料对砂浆的约束作用也直接影响着混凝土的强度以及微裂缝的形成,具有不同形貌的骨料,其约束作用的方向不同,引起裂缝形成的部位也不同。因此,开展骨料形貌的研究是十分有必要的。Aggregate mainly plays the role of skeleton and filling in concrete. At the same time, it can reduce the volume change caused by wet expansion and dry shrinkage in the process of concrete hardening. The restraint effect of aggregate on mortar also directly affects the strength of concrete and the formation of micro-cracks. Aggregates with different morphologies have different confinement directions and different positions that cause cracks to form. Therefore, it is very necessary to carry out the research of aggregate morphology.

骨料形状特征主要包括颗粒尺寸、球形度、棱角性、表面纹理等,是影响混凝土工作性、强度和收缩变形的重要因素之一。近年来不少学者开展了骨料形状特征对混凝土性能影响的研究,主要研究方向可分为骨料粒径及级配对混凝土性能的影响与骨料粒形对混凝土性能的影响两类。过去很长一段时间,相关科研工作者更多关注粒径与级配对混凝土性能的影响,在粒形对混凝土性能的影响方面研究不多;随着高强高性能混凝土研究的深入,粒形对混凝土形貌的影响才逐渐引起科研工作者的重视。对于粗骨料颗粒的粒形,目前还没有一个较全面统一的评定方法,实际工程中大多将针片状颗粒含量作为控制粗骨料颗粒形状的指标,而这无法从本质上精准表征粗骨料颗粒形状的不规则程度,不能满足当前混凝土材料研究细微观的要求。因此,需开展骨料的三维形貌研究。The shape characteristics of aggregate mainly include particle size, sphericity, angularity, surface texture, etc., and are one of the important factors affecting the workability, strength and shrinkage deformation of concrete. In recent years, many scholars have carried out research on the influence of aggregate shape characteristics on concrete performance. For a long time in the past, relevant scientific researchers paid more attention to the influence of particle size and gradation on concrete performance, and there was not much research on the influence of particle shape on concrete performance. The influence of morphology has gradually attracted the attention of researchers. For the particle shape of coarse aggregate particles, there is currently no comprehensive and unified evaluation method. In practical projects, the content of needle-like particles is mostly used as an index to control the shape of coarse aggregate particles, which cannot essentially accurately characterize coarse aggregate particles. The irregularity of the shape of the material particles cannot meet the requirements of the current research on concrete materials. Therefore, it is necessary to carry out research on the three-dimensional morphology of aggregates.

目前,骨料的三维重建方法有基于CT、激光扫描和SfM这三种方式获得骨料的三维数据。但是这三种重建方法存在生成数据量大、形貌参数计算困难等问题。球谐重构通过对三维数据进行处理并重构三维形貌,可大大降低数据的存储空间,方便后续对骨料球形度、棱角度、表面粗糙度等相关参数的计算。因此,球谐重构具有较大的应用潜力。At present, three-dimensional reconstruction methods of aggregates are based on CT, laser scanning and SfM to obtain three-dimensional data of aggregates. However, these three reconstruction methods have problems such as large amount of generated data and difficulty in calculating topographic parameters. By processing the 3D data and reconstructing the 3D topography, the spherical harmonic reconstruction can greatly reduce the storage space of the data and facilitate the subsequent calculation of related parameters such as aggregate sphericity, angularity, and surface roughness. Therefore, spherical harmonic reconstruction has great application potential.

此外,美国国家标准与技术研究院,荷兰代尔夫特理工大学等已基于CT扫描的骨料三维体素数据,建立球谐重构算法。同时,现有技术实施方案为直接通过三维扫描数据进行计算或基于体素数据的球谐重构,直接通过三维扫描数据计算,数据分布不规律,计算量大,且计算效果不好;基于体素数据的球谐重构的数据获取通常为CT,该方法对设备要求高。In addition, the National Institute of Standards and Technology of the United States, Delft University of Technology in the Netherlands, etc. have established spherical harmonic reconstruction algorithms based on the three-dimensional voxel data of aggregates from CT scans. At the same time, the prior art implementation is to directly calculate through 3D scanning data or perform spherical harmonic reconstruction based on voxel data, and directly calculate through 3D scanning data, the data distribution is irregular, the amount of calculation is large, and the calculation effect is not good; The data acquisition of spherical harmonic reconstruction of pixel data is usually CT, and this method requires high equipment.

发明内容SUMMARY OF THE INVENTION

为了克服现有技术存在的缺陷与不足,本发明提供了一种基于三维点云数据的骨料三维形貌球谐重构方法,该方法将三维点云数据的存储方式转换为球谐系数的存储方式,实现了低存储量压缩率,从而节约大量的存储空间,以达到更快速、批量可视化处理与形貌参数计算的目的。In order to overcome the defects and deficiencies of the prior art, the present invention provides a spherical harmonic reconstruction method of aggregate three-dimensional shape based on three-dimensional point cloud data, which converts the storage mode of three-dimensional point cloud data into spherical harmonic coefficients. The storage method realizes a low storage volume compression rate, thereby saving a large amount of storage space, so as to achieve the purpose of faster, batch visualization processing and shape parameter calculation.

为了达到上述目的,本发明采用以下技术方案:In order to achieve the above object, the present invention adopts the following technical solutions:

一种基于三维点云数据的骨料三维形貌球谐重构计算方法,包括步骤:A method for calculating spherical harmonic reconstruction of three-dimensional topography of aggregates based on three-dimensional point cloud data, comprising the steps of:

步骤(1):从三维数据中提取骨料的三维点云坐标;Step (1): extract the three-dimensional point cloud coordinates of the aggregate from the three-dimensional data;

步骤(2):表面参数化,将骨料三维点云数据由笛卡尔坐标系转换为球坐标系,根据高斯正交点取极角θ与方位角

Figure BDA0003711927450000021
的正则增量,基于角度逼近算法计算得到重构点半径值
Figure BDA0003711927450000022
Step (2): Surface parameterization, convert the three-dimensional point cloud data of the aggregate from the Cartesian coordinate system to the spherical coordinate system, and take the polar angle θ and azimuth angle according to the Gaussian orthogonal point
Figure BDA0003711927450000021
The regular increment of , and the radius value of the reconstructed point is calculated based on the angle approximation algorithm
Figure BDA0003711927450000022

步骤(3):基于勒让德多项式结合递归算法求取球谐系数anmStep (3): obtain spherical harmonic coefficient a nm based on Legendre polynomial combined with recursive algorithm;

步骤(4):运用步骤(3)求取的球谐系数anm与球函数

Figure BDA0003711927450000023
进行傅里叶级数求和,求取表征骨料形貌的曲面函数
Figure BDA0003711927450000024
运用曲面函数
Figure BDA0003711927450000025
可视化球谐重构结果。Step (4): Use the spherical harmonic coefficient a nm and spherical function obtained in step (3)
Figure BDA0003711927450000023
Perform a Fourier series summation to obtain the surface function that characterizes the aggregate morphology
Figure BDA0003711927450000024
Use surface functions
Figure BDA0003711927450000025
Visualize spherical harmonic reconstruction results.

进一步地,步骤(1)具体包括:Further, step (1) specifically includes:

步骤(1.1):对三维数据进行格式转化为预设格式;Step (1.1): convert the three-dimensional data into a preset format;

步骤(1.2):根据预设格式的结构判断存储骨料三维坐标数组所在位置;Step (1.2): determine the location of the stored aggregate three-dimensional coordinate array according to the structure of the preset format;

步骤(1.3):预设三维坐标数组提取的初始条件与结束条件,提取骨料的三维点云坐标。Step (1.3): Preset the initial conditions and end conditions of the three-dimensional coordinate array extraction, and extract the three-dimensional point cloud coordinates of the aggregate.

进一步地,预设格式具体采用Meshlab中的wrl格式。Further, the preset format specifically adopts the wrl format in Meshlab.

进一步地,三维坐标数组提取的初始条件为:分割字符串的最后一个元素等于point,且下一行存在提取起始字符;Further, the initial conditions for the extraction of the three-dimensional coordinate array are: the last element of the split string is equal to point, and the next line has an extraction start character;

三维坐标数组提取的结束条件为:从开始提取三维坐标往下遍历,首次遇见提取结束字符则完成三维坐标数组提取。The end condition of the extraction of the three-dimensional coordinate array is: traverse down from the beginning of extracting the three-dimensional coordinates, and complete the extraction of the three-dimensional coordinate array when the extraction end character is encountered for the first time.

进一步地,步骤(2)具体包括:Further, step (2) specifically includes:

步骤(2.1):将坐标原点移动到颗粒内部,选定z轴的正方向的单位向量为极角计算的初始向量,x轴正方向的单位向量为方位角计算的初始向量,根据点云坐标计算对应点云坐标点的极角、方位角与半径值;Step (2.1): Move the coordinate origin to the inside of the particle, select the unit vector in the positive direction of the z-axis as the initial vector for polar angle calculation, and the unit vector in the positive direction of the x-axis is the initial vector for azimuth calculation, according to the point cloud coordinates Calculate the polar angle, azimuth angle and radius value of the corresponding point cloud coordinate point;

步骤(2.2):设定极角θ取值范围为[0,π],方位角

Figure BDA0003711927450000031
的取值范围为[0,2π];结合高斯-勒让德求积公式,运用牛顿迭代算法计算勒让德多项式零点,取勒让德多项式的零点作为高斯正交点,并结合勒让德多项式计算高斯正交点权重;将极角与方位角按照求取的高斯正交点分别分为间距不等的g份;得到重构选取g2个点的极角θ(i)与方位角
Figure BDA0003711927450000032
Step (2.2): Set the value range of the polar angle θ to [0, π], the azimuth angle
Figure BDA0003711927450000031
The value range of is [0,2π]; Combined with the Gauss-Legendre quadrature formula, use the Newton iteration algorithm to calculate the Legendre polynomial zeros, take the Legendre polynomial zeros as the Gaussian orthogonal points, and combine Legendre polynomial zeros Calculate the weight of the Gaussian orthogonal point by polynomial; divide the polar angle and the azimuth angle into g parts with unequal spacing according to the obtained Gaussian orthogonal points; obtain the polar angle θ(i) and azimuth angle of the reconstructed and selected g 2 points
Figure BDA0003711927450000032

步骤(2.3):运用角度逼近算法,将步骤(2.1)计算所得极角、方位角与步骤(2.2)重构选取的点的极角θ(i)与方位角

Figure BDA0003711927450000033
分别做差,取差值绝对值之和最小点的半径为重构点半径,遍历所有重构点,得到描述骨料颗粒表面形貌的所有重构选取点的半径值
Figure BDA0003711927450000034
Step (2.3): Using the angle approximation algorithm, reconstruct the polar angle and azimuth angle calculated in step (2.1) with the polar angle θ(i) and azimuth angle of the selected point in step (2.2).
Figure BDA0003711927450000033
Make the difference respectively, take the radius of the minimum point of the sum of the absolute values of the difference as the radius of the reconstruction point, traverse all the reconstruction points, and get the radius value of all the reconstruction points that describe the surface topography of the aggregate particles
Figure BDA0003711927450000034

进一步地,根据点云坐标计算对应点云坐标的极角θ′、方位角

Figure BDA0003711927450000035
与半径值r′,具体采用第一坐标转换公式,所述第一坐标转换公式表示为:Further, calculate the polar angle θ′ and azimuth angle corresponding to the point cloud coordinates according to the point cloud coordinates
Figure BDA0003711927450000035
and the radius value r', the first coordinate conversion formula is specifically adopted, and the first coordinate conversion formula is expressed as:

Figure BDA0003711927450000041
Figure BDA0003711927450000041

Figure BDA0003711927450000042
Figure BDA0003711927450000042

Figure BDA0003711927450000043
Figure BDA0003711927450000043

其中,a为以坐标原点为起点指向x轴正向的单位向量,c为以坐标原点为起点指向z轴正向的单位向量,d为坐标原点到骨料表面点M的向量,b为向量d在x0y平面的投影,x、y、z为骨料点云在笛卡尔坐标系的坐标。Among them, a is the unit vector starting from the coordinate origin and pointing to the positive direction of the x-axis, c is the unit vector starting from the coordinate origin and pointing to the positive direction of the z-axis, d is the vector from the coordinate origin to the aggregate surface point M, and b is the vector d is the projection of the x0y plane, and x, y, and z are the coordinates of the aggregate point cloud in the Cartesian coordinate system.

进一步地,步骤(3)包括:Further, step (3) includes:

步骤(3.1):设定重构级数n,运用连带勒让德函数与勒让德多项式计算重构点球函数

Figure BDA0003711927450000044
重构点球函数
Figure BDA0003711927450000045
具体为:Step (3.1): Set the reconstruction series n, and use the associated Legendre function and Legendre polynomial to calculate the reconstructed penalty function
Figure BDA0003711927450000044
Refactoring the penalty function
Figure BDA0003711927450000045
Specifically:

Figure BDA0003711927450000046
Figure BDA0003711927450000046

Figure BDA0003711927450000047
Figure BDA0003711927450000047

Figure BDA0003711927450000048
Figure BDA0003711927450000048

其中,x=cos(θ),

Figure BDA0003711927450000049
为m阶n次连带勒让德函数,Pn(x)为n次勒让德多项式;where x=cos(θ),
Figure BDA0003711927450000049
is the m-order n-order associated Legendre function, and P n (x) is the n-order Legendre polynomial;

步骤(3.2):基于步骤(2.3)计算得到的重构点的半径

Figure BDA00037119274500000410
与步骤(3.1)计算得到的重构点的球函数
Figure BDA00037119274500000411
采用球谐系数计算公式求解球谐系数anm,具体表示为:Step (3.2): Based on the radius of the reconstructed point calculated in step (2.3)
Figure BDA00037119274500000410
with the spherical function of the reconstructed point calculated in step (3.1)
Figure BDA00037119274500000411
Use the spherical harmonic coefficient calculation formula to solve the spherical harmonic coefficient a nm , which is specifically expressed as:

Figure BDA0003711927450000051
Figure BDA0003711927450000051

进一步地,步骤(4)包括:Further, step (4) includes:

步骤(4.1):运用球面上的广义傅里叶级数,求取展开级数为n时的曲面函数

Figure BDA0003711927450000052
Step (4.1): Use the generalized Fourier series on the sphere to find the surface function when the expansion series is n
Figure BDA0003711927450000052

Figure BDA0003711927450000053
Figure BDA0003711927450000053

其中,

Figure BDA0003711927450000054
为描述颗粒表面形貌的曲面函数;anm为展开级数为n时球谐系数;
Figure BDA0003711927450000055
为展开级数为n时的球函数;in,
Figure BDA0003711927450000054
is the surface function describing the particle surface morphology; a nm is the spherical harmonic coefficient when the expansion series is n;
Figure BDA0003711927450000055
is the spherical function when the expansion series is n;

步骤(4.2):运用曲面函数

Figure BDA0003711927450000056
计算骨料的三维重构坐标;Step (4.2): Apply the surface function
Figure BDA0003711927450000056
Calculate the three-dimensional reconstruction coordinates of the aggregate;

步骤(4.3):根据wrl格式的结构,生成可视化所需的参数。Step (4.3): Generate parameters required for visualization according to the structure of the wrl format.

进一步地,在步骤(4.2)中,基于表面曲面函数

Figure BDA0003711927450000057
采用第二坐标转换公式,计算骨料的三维重构坐标(x’,y’,z’),其中第二坐标转换公式具体为:Further, in step (4.2), based on the surface surface function
Figure BDA0003711927450000057
The second coordinate conversion formula is used to calculate the three-dimensional reconstructed coordinates (x', y', z') of the aggregate, where the second coordinate conversion formula is specifically:

Figure BDA0003711927450000058
Figure BDA0003711927450000058

Figure BDA0003711927450000059
Figure BDA0003711927450000059

Figure BDA00037119274500000510
Figure BDA00037119274500000510

其中x’、y’、z’构成骨料的三维重构坐标(x’,y’,z’)。Where x', y', z' constitute the three-dimensional reconstruction coordinates (x', y', z') of the aggregate.

本发明与现有技术相比,具有如下优点和有益效果:Compared with the prior art, the present invention has the following advantages and beneficial effects:

相比现有技术实施方案中采用直接对三维扫描数据进行计算或运用CT技术基于体素数据的球谐重构,本发明提出的基于三维点云数据的骨料三维形貌球谐重构方法,采用基于点云数据的球谐重构,首先在数据的获取方面更简单,其次本发明通过将三维点云数据的存储方式转换为球谐系数的存储方式,实现了低存储量压缩率,从而在节约大量的存储空间的同时简化了形貌参数的计算。Compared with the prior art implementation that directly calculates the three-dimensional scan data or uses CT technology to perform spherical harmonic reconstruction based on voxel data, the method for spherical harmonic reconstruction of aggregate three-dimensional topography based on three-dimensional point cloud data proposed by the present invention , using spherical harmonic reconstruction based on point cloud data, firstly, the data acquisition is simpler, and secondly, the present invention achieves a low storage volume compression rate by converting the storage method of three-dimensional point cloud data into the storage method of spherical harmonic coefficients, Therefore, the calculation of the topography parameters is simplified while saving a large amount of storage space.

附图说明Description of drawings

图1为本发明基于三维点云数据的骨料三维形貌球谐重构方法的步骤流程图。FIG. 1 is a flow chart of the steps of the method for spherical harmonic reconstruction of three-dimensional topography of aggregates based on three-dimensional point cloud data according to the present invention.

图2为本发明wrl格式中提取三维坐标数组的初始条件与截止条件示意图。FIG. 2 is a schematic diagram of initial conditions and cut-off conditions for extracting a three-dimensional coordinate array in the wrl format of the present invention.

图3为本发明运用向量方法将笛卡尔坐标系转换为球坐标系示意图。FIG. 3 is a schematic diagram of converting a Cartesian coordinate system into a spherical coordinate system using a vector method in the present invention.

图4(a)为本发明三维扫描仪获取骨料在球谐重构前的结构示意图;Figure 4(a) is a schematic structural diagram of the aggregate obtained by the three-dimensional scanner of the present invention before spherical harmonic reconstruction;

图4(b)为本发明三维扫描仪获取骨料在球谐重构后的结构示意图。Figure 4(b) is a schematic structural diagram of the aggregate obtained by the three-dimensional scanner of the present invention after spherical harmonic reconstruction.

图5为本发明中骨料在球谐重构前后内存占用量示意图。FIG. 5 is a schematic diagram of the memory occupancy of aggregates before and after spherical harmonic reconstruction in the present invention.

具体实施方式Detailed ways

在本公开的描述中,需要说明的是,术语“中心”、“上”、“下”、“左”、“右”、“竖直”、“水平”、“内”、“外”等指示的方位或位置关系为基于附图所示的方位或位置关系,仅是为了便于描述本公开和简化描述,而不是指示或暗示所指的装置或元件必须具有特定的方位、以特定的方位构造和操作,因此不能理解为对本公开的限制。In the description of the present disclosure, it should be noted that the terms "center", "upper", "lower", "left", "right", "vertical", "horizontal", "inner", "outer", etc. The indicated orientation or positional relationship is based on the orientation or positional relationship shown in the drawings, and is only for the convenience of describing the present disclosure and simplifying the description, rather than indicating or implying that the indicated device or element must have a specific orientation or a specific orientation. construction and operation, and therefore should not be construed as limiting the present disclosure.

此外,术语“第一”、“第二”、“第三”仅用于描述目的,而不能理解为指示或暗示相对重要性。同样,“一个”、“一”或者“该”等类似词语也不表示数量限制,而是表示存在至少一个。“包括”或者“包含”等类似的词语意指出现在该词前面的元素或者物件涵盖出现在该词后面列举的元素或者物件及其等同,而不排除其他元素或者物件。“连接”或者“相连”等类似的词语并非限定于物理的或者机械的连接,而是可以包括电性的连接,不管是直接的还是间接的。Furthermore, the terms "first", "second", and "third" are used for descriptive purposes only and should not be construed to indicate or imply relative importance. Likewise, words such as "a," "an," or "the" do not denote a limitation of quantity, but rather denote the presence of at least one. "Comprises" or "comprising" and similar words mean that the elements or things appearing before the word encompass the elements or things recited after the word and their equivalents, but do not exclude other elements or things. Words like "connected" or "connected" are not limited to physical or mechanical connections, but may include electrical connections, whether direct or indirect.

在本公开的描述中,需要说明的是,除非另有明确的规定和限定,否则术语“安装”、“相连”、“连接”应做广义理解。例如,可以是固定连接,也可以是可拆卸连接,或一体地连接;可以是机械连接,也可以是电连接;可以是直接相连,也可以通过中间媒介间接相连,可以是两个元件内部的连通。对于本领域的普通技术人员而言,可以具体情况理解上述术语在本公开中的具体含义。此外,下面所描述的本公开不同实施方式中所涉及的技术特征只要彼此之间未构成冲突就可以相互结合。In the description of the present disclosure, it should be noted that the terms "installed", "connected" and "connected" should be construed in a broad sense unless otherwise expressly specified and limited. For example, it can be a fixed connection, a detachable connection, or an integral connection; it can be a mechanical connection or an electrical connection; it can be a direct connection, or an indirect connection through an intermediate medium, or an internal connection between two components. Connected. For those of ordinary skill in the art, the specific meanings of the above terms in the present disclosure can be understood in specific situations. In addition, the technical features involved in the different embodiments of the present disclosure described below can be combined with each other as long as they do not conflict with each other.

为了使本发明的目的、技术方案及优点更加清楚明白,以下结合附图及实施例,对本发明进行进一步详细说明。应当理解,此处所描述的具体实施例仅仅用以解释本发明,并不用于限定本发明。In order to make the objectives, technical solutions and advantages of the present invention clearer, the present invention will be further described in detail below with reference to the accompanying drawings and embodiments. It should be understood that the specific embodiments described herein are only used to explain the present invention, but not to limit the present invention.

实施例Example

如图1所示,本实施例提供了一种基于三维点云数据的骨料三维形貌球谐重构方法,该方法包括步骤:As shown in FIG. 1 , this embodiment provides a method for spherical harmonic reconstruction of three-dimensional topography of aggregates based on three-dimensional point cloud data. The method includes the following steps:

步骤(1):从三维扫描、SfM等方法生成的三维数据中提取骨料的三维点云坐标;Step (1): extract the three-dimensional point cloud coordinates of the aggregate from the three-dimensional data generated by three-dimensional scanning, SfM and other methods;

在本实施例中,步骤(1)具体包括:In this embodiment, step (1) specifically includes:

步骤(1.1):对三维数据进行格式转化为预设格式;Step (1.1): convert the three-dimensional data into a preset format;

步骤(1.2):根据预设格式的结构判断存储骨料三维坐标数组所在位置;Step (1.2): determine the location of the stored aggregate three-dimensional coordinate array according to the structure of the preset format;

步骤(1.3):根据预设的三维坐标数组提取的初始条件与结束条件进行提取骨料的三维点云坐标。Step (1.3): Extract the three-dimensional point cloud coordinates of the aggregate according to the initial conditions and the end conditions of the preset three-dimensional coordinate array extraction.

实际应用时,首先将三维扫描获得的三维数据运用格式转换软件,例如Meshlab等转换为wrl格式。根据wrl格式的结构,判断存储骨料三维坐标数组所在位置,进而利用预设的三维坐标数组提取的初始条件与结束条件进行提取骨料的三维点云坐标。In practical application, the 3D data obtained by 3D scanning is first converted into wrl format using format conversion software, such as Meshlab. According to the structure of the wrl format, determine the location of the three-dimensional coordinate array of the stored aggregate, and then extract the three-dimensional point cloud coordinates of the aggregate by using the initial conditions and end conditions of the preset three-dimensional coordinate array extraction.

步骤(2):表面参数化,将骨料三维点云数据由笛卡尔坐标系转换为球坐标系,根据高斯正交点取极角θ与方位角

Figure BDA0003711927450000071
的正则增量,基于角度逼近算法得到重构点半径值
Figure BDA0003711927450000081
Step (2): Surface parameterization, convert the three-dimensional point cloud data of aggregates from Cartesian coordinate system to spherical coordinate system, and take polar angle θ and azimuth angle according to Gaussian orthogonal points
Figure BDA0003711927450000071
The regular increment of , and the radius value of the reconstructed point is obtained based on the angle approximation algorithm
Figure BDA0003711927450000081

在本实施例中,步骤(2)具体包括:In this embodiment, step (2) specifically includes:

步骤(2.1):将坐标原点移动到颗粒内部,选定z轴的正方向的单位向量为极角计算的初始向量,x轴正方向的单位向量为方位角计算的初始向量,根据点云坐标计算对应点云坐标点的极角、方位角与半径值;Step (2.1): Move the coordinate origin to the inside of the particle, select the unit vector in the positive direction of the z-axis as the initial vector for polar angle calculation, and the unit vector in the positive direction of the x-axis is the initial vector for azimuth calculation, according to the point cloud coordinates Calculate the polar angle, azimuth angle and radius value of the corresponding point cloud coordinate point;

步骤(2.2):设定极角θ取值范围为[0,π],方位角

Figure BDA0003711927450000082
的取值范围为[0,2π];结合高斯-勒让德求积公式,运用牛顿迭代算法计算勒让德多项式零点,取勒让德多项式的零点作为高斯正交点,并结合勒让德多项式计算高斯正交点权重;将极角与方位角按照求取的高斯正交点分别分为间距不等的g份;得到重构选取g2个点的极角θ(i)与方位角
Figure BDA0003711927450000083
Step (2.2): Set the value range of the polar angle θ to [0, π], the azimuth angle
Figure BDA0003711927450000082
The value range of is [0,2π]; combined with the Gauss-Legendre quadrature formula, use the Newton iterative algorithm to calculate the zeros of the Legendre polynomial, take the zeros of the Legendre polynomial as the Gauss orthogonal point, and combine the Legendre polynomial Calculate the weight of the Gaussian orthogonal point by polynomial; divide the polar angle and the azimuth angle into g parts with unequal spacing according to the obtained Gaussian orthogonal points; obtain the polar angle θ(i) and azimuth angle of the reconstructed and selected g 2 points
Figure BDA0003711927450000083

Figure BDA0003711927450000084
Figure BDA0003711927450000084

nPn(x)=(2n-1)xPn-1(x)-(n-1)Pn-2(x) (式2)nP n (x)=(2n-1)xP n-1 (x)-(n-1)P n-2 (x) (Equation 2)

(1-x2)Pn'(x)=nPn-1(x)-nxPn(x) (式3)(1-x 2 )P n '(x)=nP n-1 (x)-nxP n (x) (Equation 3)

Figure BDA0003711927450000085
Figure BDA0003711927450000085

xi=xk+1,(lim(xk+1-xk)=0) (式5)x i =x k+1 , (lim(x k+1 -x k )=0) (Equation 5)

Figure BDA0003711927450000086
Figure BDA0003711927450000086

Figure BDA0003711927450000087
Figure BDA0003711927450000087

Figure BDA0003711927450000091
Figure BDA0003711927450000091

其中Pn(x)为n次勒让德多项式,P’n(x)为n次勒让德多项式的倒数,xk为进行牛顿迭代时第k-1次迭代时对应的x的值,xi为勒让德多项式的零点,wi为高斯正交点xi的权重,极角θ(i)与方位角

Figure BDA0003711927450000092
表示对应xi,xj处的极角与方位角值,其中i,j∈(1,2,3,4......g);式1为勒让德多项式的级数与微分表示,式2、式3为勒让德多项式及其导数的递推公式,式4、式5为运用牛顿迭代算法计算勒让德多项式的零点,式6为高斯勒让德正交点权重算法,式7、式8为极角与方位角的划分方法。where P n (x) is the n-th Legendre polynomial, P' n (x) is the reciprocal of the n-th Legendre polynomial, and x k is the value of x corresponding to the k-1th iteration during Newton iteration, xi is the zero point of Legendre polynomial, wi is the weight of the Gaussian orthogonal point xi , the polar angle θ(i) is related to the azimuth angle
Figure BDA0003711927450000092
Represents the polar angle and azimuth angle value corresponding to x i , x j , where i,j∈(1,2,3,4...g); Equation 1 is the series and differential of Legendre polynomials Representation, Equation 2 and Equation 3 are the recursive formulas of Legendre polynomial and its derivatives, Equation 4 and Equation 5 are the use of Newton iteration algorithm to calculate the zero point of Legendre polynomial, Equation 6 is the Gauss Legendre orthogonal point weight algorithm , Equation 7 and Equation 8 are the division methods of polar angle and azimuth angle.

步骤(2.3):运用角度逼近算法,将步骤(2.1)计算所得极角、方位角与步骤(2.2)重构选取的点的极角θ(i)与方位角

Figure BDA0003711927450000093
分别做差,取差值绝对值之和最小点的半径为重构点半径,遍历所有重构点,得到描述骨料颗粒表面形貌的所有重构选取点的半径值
Figure BDA0003711927450000094
Step (2.3): Using the angle approximation algorithm, reconstruct the polar angle and azimuth angle calculated in step (2.1) with the polar angle θ(i) and azimuth angle of the selected point in step (2.2).
Figure BDA0003711927450000093
Make the difference respectively, take the radius of the minimum point of the sum of the absolute values of the difference as the radius of the reconstructed point, traverse all the reconstructed points, and obtain the radius value of all the reconstructed selected points describing the surface morphology of the aggregate particles
Figure BDA0003711927450000094

步骤(3):基于勒让德多项式结合递归算法求取球谐系数;Step (3): obtain spherical harmonic coefficient based on Legendre polynomial combined with recursive algorithm;

在本实施例中,步骤(3)具体包括:In this embodiment, step (3) specifically includes:

步骤(3.1):设定重构级数n,运用连带勒让德函数与勒让德多项式计算步骤(2.2)中重构点的球函数

Figure BDA0003711927450000095
Step (3.1): Set the reconstruction series n, and use the associated Legendre function and Legendre polynomial to calculate the spherical function of the reconstructed point in step (2.2)
Figure BDA0003711927450000095

Figure BDA0003711927450000096
Figure BDA0003711927450000096

Figure BDA0003711927450000097
Figure BDA0003711927450000097

其中x表示重构选取点的极角θ的余弦值,即x=cos(θ),

Figure BDA0003711927450000098
为m阶n次连带勒让德函数m∈(-n,-n+1,....,0,.....n-1,n),Pn(x)为n次勒让德多项式。where x represents the cosine value of the polar angle θ of the reconstructed selected point, that is, x=cos(θ),
Figure BDA0003711927450000098
is the m-order n-time associated Legendre function m∈(-n,-n+1,....,0,.....n-1,n), and P n (x) is the n-time Legendre function German polynomial.

步骤(3.2):基于步骤(2.3)计算得到的重构点的半径

Figure BDA0003711927450000099
与步骤(3.1)计算得到的重构点的球函数
Figure BDA0003711927450000101
采用球谐系数计算公式求解球谐系数anm;Step (3.2): Based on the radius of the reconstructed point calculated in step (2.3)
Figure BDA0003711927450000099
with the spherical function of the reconstructed point calculated in step (3.1)
Figure BDA0003711927450000101
Use the spherical harmonic coefficient calculation formula to solve the spherical harmonic coefficient a nm ;

Figure BDA0003711927450000102
Figure BDA0003711927450000102

其中上式为球谐系数anm计算公式,

Figure BDA0003711927450000103
的值计算时可直接选用
Figure BDA0003711927450000104
中极角θ(i)与方位角
Figure BDA0003711927450000105
在该高斯正交点处的权重w(i)w(j)。The above formula is the calculation formula of spherical harmonic coefficient a nm ,
Figure BDA0003711927450000103
can be directly selected when calculating the value of
Figure BDA0003711927450000104
Mid-polar angle θ(i) and azimuth
Figure BDA0003711927450000105
Weight w(i)w(j) at this Gaussian quadrature point.

步骤(4):运用步骤(3)求取的球谐系数anm与球函数

Figure BDA0003711927450000106
进行傅里叶级数求和,求取表征骨料形貌的曲面函数
Figure BDA0003711927450000107
运用曲面函数
Figure BDA0003711927450000108
可视化球谐重构结果。Step (4): Use the spherical harmonic coefficient a nm and spherical function obtained in step (3)
Figure BDA0003711927450000106
Perform a Fourier series summation to obtain the surface function that characterizes the aggregate morphology
Figure BDA0003711927450000107
Use surface functions
Figure BDA0003711927450000108
Visualize spherical harmonic reconstruction results.

在本实施例中,步骤(4)具体包括:In this embodiment, step (4) specifically includes:

步骤(4.1):运用球面上的广义傅里叶级数,求取展开级数为n时的曲面函数

Figure BDA0003711927450000109
Step (4.1): Use the generalized Fourier series on the sphere to find the surface function when the expansion series is n
Figure BDA0003711927450000109

Figure BDA00037119274500001010
Figure BDA00037119274500001010

其中

Figure BDA00037119274500001011
为描述颗粒表面形貌的曲面函数;anm为展开级数为n时球谐系数;
Figure BDA00037119274500001012
为展开级数为n时的球函数。in
Figure BDA00037119274500001011
is the surface function describing the particle surface morphology; a nm is the spherical harmonic coefficient when the expansion series is n;
Figure BDA00037119274500001012
is the spherical function when the expansion series is n.

步骤(4.2):运用曲面函数

Figure BDA00037119274500001013
计算骨料的三维重构坐标;Step (4.2): Apply the surface function
Figure BDA00037119274500001013
Calculate the three-dimensional reconstruction coordinates of the aggregate;

步骤(4.3):根据wrl格式的结构,生成可视化所需的参数。Step (4.3): Generate parameters required for visualization according to the structure of the wrl format.

下面结合附图对本发明做进一步说明:The present invention will be further described below in conjunction with the accompanying drawings:

结合图1所示,一种基于三维点云数据的骨料三维形貌球谐重构方法,包括步骤:With reference to Figure 1, a spherical harmonic reconstruction method of aggregate 3D morphology based on 3D point cloud data, including steps:

(1)执行前处理工作:(1) Pre-processing work:

①将三维扫描获取的stl格式数据运用MeshLab软件将存储格式转换为wrl格式;①Use MeshLab software to convert the stl format data obtained by 3D scanning to wrl format;

②根据wrl格式的结构,判断骨料三维坐标数组所在位置;如图2所示,即为wrl格式的基本结构,point数组中所包含的数据即为骨料的三维坐标。② According to the structure of the wrl format, determine the location of the three-dimensional coordinate array of the aggregate; as shown in Figure 2, it is the basic structure of the wrl format, and the data contained in the point array is the three-dimensional coordinate of the aggregate.

③遍历阅读wrl文件的每一行,给出三维坐标数组提取的初始条件为:分割字符串的最后一个元素等于point,且下一行存在“[”字符;结束条件为:从开始提取三维坐标往下遍历,首次遇见“]”字符则完成三维坐标数组提取。③ Traverse and read each line of the wrl file, and give the initial conditions for the extraction of the three-dimensional coordinate array: the last element of the split string is equal to point, and there is a "[" character in the next line; the end condition is: from the beginning to extract the three-dimensional coordinates down Traverse, the first time the "]" character is encountered, the three-dimensional coordinate array extraction is completed.

(2)球谐重构:(2) Spherical harmonic reconstruction:

①如图3所示,将坐标原点移动到颗粒内部,选定z轴的正方向的单位向量为极角θ角度计算的初始向量,x轴正方向的单位向量为方位角

Figure BDA0003711927450000111
角度计算的初始向量,根据点云坐标计算对应点云坐标的极角θ′、方位角
Figure BDA0003711927450000112
与半径值r′,其中根据点云坐标计算对应点云坐标的极角θ′、方位角
Figure BDA0003711927450000113
与半径值r′具体采用第一坐标转换公式,该第一坐标转换公式表示为:①As shown in Figure 3, move the coordinate origin to the inside of the particle, select the unit vector in the positive direction of the z-axis as the initial vector calculated by the polar angle θ, and the unit vector in the positive direction of the x-axis is the azimuth angle
Figure BDA0003711927450000111
The initial vector for angle calculation, and the polar angle θ′ and azimuth angle of the corresponding point cloud coordinates are calculated according to the point cloud coordinates.
Figure BDA0003711927450000112
and the radius value r′, in which the polar angle θ′ and azimuth angle of the corresponding point cloud coordinates are calculated according to the point cloud coordinates
Figure BDA0003711927450000113
The first coordinate conversion formula is specifically used for the radius value r', and the first coordinate conversion formula is expressed as:

Figure BDA0003711927450000114
Figure BDA0003711927450000114

Figure BDA0003711927450000115
Figure BDA0003711927450000115

Figure BDA0003711927450000116
Figure BDA0003711927450000116

其中,a为以坐标原点为起点指向x轴正向的单位向量,c为以坐标原点为起点指向z轴正向的单位向量,d为坐标原点O到骨料表面点M的向量,b为向量d在x0y平面的投影,x、y、z为骨料点云在笛卡尔坐标系的坐标。Among them, a is the unit vector starting from the coordinate origin and pointing to the positive direction of the x-axis, c is the unit vector starting from the coordinate origin and pointing to the positive direction of the z-axis, d is the vector from the coordinate origin O to the aggregate surface point M, and b is The projection of the vector d on the x0y plane, x, y, and z are the coordinates of the aggregate point cloud in the Cartesian coordinate system.

②设定极角θ取值范围为[0,π],方位角

Figure BDA0003711927450000117
的取值范围为[0,2π],结合高斯-勒让德求积公式,运用牛顿迭代算法计算勒让德多项式的零点,取勒让德多项式的零点作为高斯正交点,并结合勒让德多项式计算高斯正交点的权重;将极角与方位角按照求取的高斯正交点分别分为间距不等的g(g>=2n)份;得到重构选取的极角θ(i)与方位角
Figure BDA0003711927450000121
i,j∈(1,2,3,4......g);②Set the value range of the polar angle θ to [0, π], and the azimuth angle
Figure BDA0003711927450000117
The value range of [0,2π], combined with the Gauss-Legendre quadrature formula, use the Newton iterative algorithm to calculate the zero point of the Legendre polynomial, take the zero point of the Legendre polynomial as the Gaussian orthogonal point, and combine the Legendre polynomial with the zero point. De polynomial calculates the weight of the Gaussian orthogonal point; divides the polar angle and the azimuth angle into g (g>=2n) parts with unequal spacing according to the obtained Gaussian orthogonal point; obtains the polar angle θ (i ) and azimuth
Figure BDA0003711927450000121
i,j∈(1,2,3,4...g);

③运用角度逼近算法,对点云坐标计算所得极角θ′、方位角

Figure BDA0003711927450000122
与重构选取的点的极角θ与方位角
Figure BDA0003711927450000123
分别做差,取差值的绝对值和最小点的半径为重构点半径
Figure BDA0003711927450000124
遍历所有重构点,得到所有重构点的半径
Figure BDA0003711927450000125
③Using the angle approximation algorithm to calculate the polar angle θ′ and azimuth angle from the point cloud coordinates
Figure BDA0003711927450000122
Polar angle θ and azimuth angle with the reconstructed selected point
Figure BDA0003711927450000123
Make the difference respectively, take the absolute value of the difference and the radius of the minimum point as the radius of the reconstructed point
Figure BDA0003711927450000124
Traverse all reconstructed points to get the radius of all reconstructed points
Figure BDA0003711927450000125

Figure BDA0003711927450000126
Figure BDA0003711927450000126

Figure BDA0003711927450000127
Figure BDA0003711927450000127

其中θ′(l)为骨料点云坐标计算出的极角,

Figure BDA0003711927450000128
为骨料点云坐标计算出的方位角(l∈(1,L),L为扫描得到的点云总数),θ′min为所有点云中与该重构点极角和方位角之差的绝对值之和最小的点的极角,
Figure BDA0003711927450000129
为所有点云中与该重构点极角和方位角之差的绝对值之和最小的点的方位角。where θ′(l) is the polar angle calculated from the aggregate point cloud coordinates,
Figure BDA0003711927450000128
The azimuth angle calculated for the aggregate point cloud coordinates (l∈(1,L), L is the total number of scanned point clouds), θ′ min is the difference between the polar and azimuth angles of all point clouds and the reconstructed point The polar angle of the point with the smallest sum of absolute values,
Figure BDA0003711927450000129
is the azimuth of the point with the smallest sum of absolute values of the difference between the polar and azimuth angles of the reconstructed point in all point clouds.

④设定重构级数为n(n>=10),运用式9-式10计算重构点的球函数

Figure BDA00037119274500001210
④ Set the reconstruction series to n (n>=10), and use equations 9 to 10 to calculate the spherical function of the reconstruction point
Figure BDA00037119274500001210

⑤运用式11求解球谐系数anm⑤ Use Equation 11 to solve the spherical harmonic coefficient a nm .

(3)后处理工作:(3) Post-processing work:

①运用式12,求取展开级数为n时的曲面函数

Figure BDA00037119274500001211
①Using Equation 12, find the surface function when the expansion series is n
Figure BDA00037119274500001211

②基于表面曲面函数

Figure BDA00037119274500001212
采用第二坐标转换公式,计算骨料的三维重构坐标(x',y',z'),其中第二坐标转换公式具体为:② Based on surface surface function
Figure BDA00037119274500001212
The second coordinate conversion formula is used to calculate the three-dimensional reconstructed coordinates (x', y', z') of the aggregate, where the second coordinate conversion formula is specifically:

Figure BDA00037119274500001213
Figure BDA00037119274500001213

Figure BDA00037119274500001214
Figure BDA00037119274500001214

Figure BDA00037119274500001215
Figure BDA00037119274500001215

其中x'、y'、z'构成骨料的三维重构坐标(x',y',z')。Among them, x', y', z' constitute the three-dimensional reconstruction coordinates (x', y', z') of the aggregate.

③结合图2所示,根据wrl格式结构,生成可视化所需参数,并将计算得到的三维坐标与点的连接方式输出到wrl文本中。③ Combined with Figure 2, according to the wrl format structure, generate the parameters required for visualization, and output the connection mode of the calculated three-dimensional coordinates and points to the wrl text.

如图4(a)和图4(b)所示,为三维扫描仪获取的骨料与球谐重构可视化结果。如图5所示,本实施例将三维点云数据的存储方式转换为球谐系数的存储方式,文件压缩率小于1%,极大地节约了存储空间,且重构获取表面函数,更便于求取形貌参数;As shown in Fig. 4(a) and Fig. 4(b), the visualization results of aggregate and spherical harmonic reconstruction obtained by the 3D scanner are shown. As shown in Figure 5, in this embodiment, the storage method of 3D point cloud data is converted to the storage method of spherical harmonic coefficients, and the file compression rate is less than 1%, which greatly saves storage space, and reconstructs and obtains surface functions, which is more convenient to find Take shape parameters;

上述实施例为本发明较佳的实施方式,但本发明的实施方式并不受上述实施例的限制,其他的任何未背离本发明的精神实质与原理下所作的改变、修饰、替代、组合、简化,均应为等效的置换方式,都包含在本发明的保护范围之内。The above-mentioned embodiments are preferred embodiments of the present invention, but the embodiments of the present invention are not limited by the above-mentioned embodiments, and any other changes, modifications, substitutions, combinations, The simplification should be equivalent replacement manners, which are all included in the protection scope of the present invention.

Claims (9)

1. An aggregate three-dimensional shape spherical resonance reconstruction method based on three-dimensional point cloud data is characterized by comprising the following steps:
step (1): extracting three-dimensional point cloud coordinates of aggregate from the three-dimensional data;
step (2): surface parameterization, namely converting the aggregate three-dimensional point cloud data into a spherical coordinate system from a Cartesian coordinate system, and taking a polar angle theta and an azimuth angle according to a Gaussian orthogonal point
Figure FDA0003711927440000011
Based on angle approximation algorithm to obtain radius value of reconstruction point
Figure FDA0003711927440000012
And (3): method for solving spherical harmonic coefficient a based on Legendre polynomial combined with recursive algorithm nm
And (4): using the spherical harmonic coefficient a nm Function of and sphere
Figure FDA0003711927440000013
Fourier series summation is carried out to obtain a curved surface function representing the shape of the aggregate
Figure FDA0003711927440000014
Using curved surface functions
Figure FDA0003711927440000015
And visualizing the spherical harmonic reconstruction result.
2. The method for reconstructing the spherical harmonic of the three-dimensional morphology of the aggregate based on the three-dimensional point cloud data according to claim 1, wherein the step (1) specifically comprises the following steps:
step (1.1): converting the format of the three-dimensional data into a preset format;
step (1.2): judging the position of a three-dimensional coordinate array of the stored aggregate according to a structure in a preset format;
step (1.3): presetting initial conditions and finishing conditions for extracting the three-dimensional coordinate array, and extracting three-dimensional point cloud coordinates of the aggregate.
3. The method for reconstructing the spherical harmonic of the three-dimensional morphology of the aggregate based on the three-dimensional point cloud data as claimed in claim 2, wherein in the step (1.1), the preset format is wrl format in Meshlab.
4. The method for reconstructing the spherical resonance of the three-dimensional morphology of the aggregate based on the three-dimensional point cloud data as claimed in claim 2, wherein the initial conditions for extracting the three-dimensional coordinate array are as follows: the last element of the segmentation character string is equal to point, and the next line has an extraction starting character;
the three-dimensional coordinate array extraction end condition is as follows: and traversing from the beginning of extracting the three-dimensional coordinates to the bottom, and finishing the extraction of the three-dimensional coordinate array when the character is extracted and ended in the first time.
5. The method for reconstructing the spherical resonance of the three-dimensional morphology of the aggregate based on the three-dimensional point cloud data according to claim 1, wherein the step (2) specifically comprises the following steps:
step (2.1): moving the origin of coordinates to the inside of the particles, selecting a unit vector in the positive direction of the z axis as an initial vector for polar angle calculation, selecting a unit vector in the positive direction of the x axis as an initial vector for azimuth angle calculation, and calculating polar angles, azimuth angles and radius values of corresponding point cloud coordinate points according to point cloud coordinates;
step (2.2): setting the polar angle theta to be in the value range of [0, pi]Azimuth angle
Figure FDA0003711927440000021
Has a value range of [0,2 pi](ii) a Calculating Legendre polynomial zero points by combining a Gaussian-Legendre product formula and applying a Newton iterative algorithm, taking the zero points of Legendre polynomials as Gaussian orthogonal points, and calculating the weight of the Gaussian orthogonal points by combining Legendre polynomials; dividing the polar angle and the azimuth angle into g parts with unequal intervals according to the solved Gaussian orthogonal point; obtaining reconstructed selection g 2 Polar angle θ (i) and azimuth angle of point
Figure FDA0003711927440000022
Step (2.3): applying an angle approximation algorithm to reconstruct the polar angle and the azimuth angle calculated in the step (2.1) and the polar angle theta (i) and the azimuth angle of the point selected by reconstruction
Figure FDA0003711927440000023
Respectively making difference, taking the radius of the minimum point of the sum of the absolute values of the difference values as the radius of the reconstruction point, traversing all the reconstruction points to obtain the radius values of all the reconstruction selection points describing the surface appearance of the aggregate particles
Figure FDA0003711927440000024
6. The method for reconstructing the spherical harmonic of the three-dimensional morphology of the aggregate based on the three-dimensional point cloud data as claimed in claim 5, wherein the polar angle θ' and the azimuth angle of the corresponding point cloud coordinate are calculated according to the point cloud coordinate
Figure FDA0003711927440000025
And the radius value r' specifically adopts a first coordinate conversion formula, wherein the first coordinate conversion formula is expressed as follows:
Figure FDA0003711927440000026
Figure FDA0003711927440000027
Figure FDA0003711927440000028
wherein a is a unit vector pointing to the positive direction of the x axis by taking the origin of coordinates as a starting point, c is a unit vector pointing to the positive direction of the z axis by taking the origin of coordinates as a starting point, d is a vector from the origin of coordinates to an aggregate surface point M, b is a projection of the vector d on an x0y plane, and x, y and z are coordinates of the aggregate point cloud in a Cartesian coordinate system.
7. The method for reconstructing the spherical resonance of the three-dimensional morphology of the aggregate based on the three-dimensional point cloud data according to claim 1, wherein the step (3) specifically comprises the following steps:
step (3.1): setting reconstruction series n, calculating reconstruction point-sphere function by using legendre function and legendre polynomial
Figure FDA0003711927440000031
Reconstructing a point-sphere function
Figure FDA0003711927440000032
The method specifically comprises the following steps:
Figure FDA0003711927440000033
Figure FDA0003711927440000034
Figure FDA0003711927440000035
wherein x = cos (θ),
Figure FDA0003711927440000036
is an m-th order n-th order associated Legendre function, P n (x) Is an n-th Legendre polynomial;
step (3.2): radius based on calculated reconstruction points
Figure FDA0003711927440000037
And (3.1) calculating the spherical function of the reconstructed point
Figure FDA0003711927440000038
Solving the spherical harmonic coefficient a nm Specifically, it is represented as:
Figure FDA0003711927440000039
8. the method for reconstructing the spherical harmonic of the three-dimensional morphology of the aggregate based on the three-dimensional point cloud data as claimed in claim 1, wherein the step (4) specifically comprises:
step (4.1): the generalized Fourier series on the spherical surface is used to solve the curved surface function with the expansion series of n
Figure FDA00037119274400000310
Figure FDA00037119274400000311
Wherein,
Figure FDA00037119274400000312
is a curved surface function for describing the surface topography of the particles; a is nm The spherical harmonic coefficient when the expansion series is n;
Figure FDA00037119274400000313
the spherical function is a spherical function when the expansion series is n;
step (4.2): using curved surface functions
Figure FDA00037119274400000314
Calculating three-dimensional reconstruction coordinates of the aggregate;
step (4.3): and generating parameters required by visualization according to the structure of the wrl format.
9. The method for reconstructing the spherical harmonic of the three-dimensional morphology of the aggregate based on the three-dimensional point cloud data as claimed in claim 8, wherein in the step (4.2), the method is based on a surface function
Figure FDA0003711927440000041
Polar angle theta and azimuth angle
Figure FDA0003711927440000042
Calculating the three-dimensional reconstruction coordinates (x ', y ', z ') of the aggregate by adopting a second coordinate conversion formula, wherein the second coordinate conversion formula specifically comprises the following steps:
Figure FDA0003711927440000043
Figure FDA0003711927440000044
Figure FDA0003711927440000045
wherein x ', y', z 'constitute the three-dimensional reconstruction coordinates (x', y ', z') of the aggregate.
CN202210722175.8A 2022-06-24 2022-06-24 Aggregate three-dimensional morphology spherical harmonic reconstruction method based on three-dimensional point cloud data Pending CN115222914A (en)

Priority Applications (1)

Application Number Priority Date Filing Date Title
CN202210722175.8A CN115222914A (en) 2022-06-24 2022-06-24 Aggregate three-dimensional morphology spherical harmonic reconstruction method based on three-dimensional point cloud data

Applications Claiming Priority (1)

Application Number Priority Date Filing Date Title
CN202210722175.8A CN115222914A (en) 2022-06-24 2022-06-24 Aggregate three-dimensional morphology spherical harmonic reconstruction method based on three-dimensional point cloud data

Publications (1)

Publication Number Publication Date
CN115222914A true CN115222914A (en) 2022-10-21

Family

ID=83610556

Family Applications (1)

Application Number Title Priority Date Filing Date
CN202210722175.8A Pending CN115222914A (en) 2022-06-24 2022-06-24 Aggregate three-dimensional morphology spherical harmonic reconstruction method based on three-dimensional point cloud data

Country Status (1)

Country Link
CN (1) CN115222914A (en)

Cited By (1)

* Cited by examiner, † Cited by third party
Publication number Priority date Publication date Assignee Title
CN116385642A (en) * 2023-03-31 2023-07-04 浙江大学 Three-dimensional scalar information compression reconstruction method based on spherical Shearlet

Cited By (2)

* Cited by examiner, † Cited by third party
Publication number Priority date Publication date Assignee Title
CN116385642A (en) * 2023-03-31 2023-07-04 浙江大学 Three-dimensional scalar information compression reconstruction method based on spherical Shearlet
CN116385642B (en) * 2023-03-31 2023-09-12 浙江大学 Three-dimensional scalar information compression and reconstruction method based on spherical Shearlet

Similar Documents

Publication Publication Date Title
CN100561523C (en) A 3D Model Mesh Reconstruction Method
CN110095060A (en) Steel construction rapid quality detection method based on 3-D scanning technology
Sun et al. Quantitative evaluation for shape characteristics of aggregate particles based on 3D point cloud data
CN105469404B (en) A kind of rotary body approximating method and device based on three dimensional point cloud
CN101533529A (en) Range image-based 3D spatial data processing method and device
CN109389596B (en) A method for evaluating the overall roughness of three-dimensional irregular particle surfaces
CN112819962B (en) Non-uniform grid division and local grid density method in digital image correlation
CN114117702A (en) Point cloud-based automatic reverse modeling method for power transmission line
CN117274527B (en) Method for constructing three-dimensional visualization model data set of generator equipment
CN112150430A (en) Numerical analysis method utilizing rock microscopic structure digital image
CN103077559A (en) Cluster three-dimensional rebuilding method based on sequence image
CN117237569A (en) A three-dimensional digital method for soil microstructure characteristics
CN115222914A (en) Aggregate three-dimensional morphology spherical harmonic reconstruction method based on three-dimensional point cloud data
CN118262065A (en) A 3D modeling method for subway tunnels based on laser point cloud
Pfeifer A subdivision algorithm for smooth 3D terrain models
CN115690316A (en) Aggregate three-dimensional reconstruction and random generation method based on spherical DOG wavelet
Kang et al. Point cloud smooth sampling and surface reconstruction based on moving least squares
CN115578524A (en) Infrared three-dimensional reconstruction method, system, storage medium and computer equipment
CN117974887A (en) Tunnel wall modeling method and system based on three-dimensional laser point cloud
Sun et al. Datum feature extraction and deformation analysis method based on normal vector of point cloud
CN115797573A (en) Method, device and medium for measuring point cloud twinning geometric accuracy
CN114636403A (en) Power transmission tower inclination measuring method based on three-dimensional laser point cloud
Pan et al. Improved QEM simplification algorithm based on local area feature information constraint
CN115082527A (en) Reference point cloud registration method and device, and high-precision point cloud registration evaluation method and device
Zou et al. Research on optimization of rendering efficiency of point cloud data of transmission lines in three-dimensional GIS

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