Movatterモバイル変換


[0]ホーム

URL:


CN116576827A - Remote sensing terrain three-dimensional mapping method based on polarization three-dimensional imaging - Google Patents

Remote sensing terrain three-dimensional mapping method based on polarization three-dimensional imaging
Download PDF

Info

Publication number
CN116576827A
CN116576827ACN202310484598.5ACN202310484598ACN116576827ACN 116576827 ACN116576827 ACN 116576827ACN 202310484598 ACN202310484598 ACN 202310484598ACN 116576827 ACN116576827 ACN 116576827A
Authority
CN
China
Prior art keywords
target
area
polarization
angle
normal vector
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
CN202310484598.5A
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.)
Hangzhou Research Institute Of Xi'an University Of Electronic Science And Technology
Xidian University
Original Assignee
Hangzhou Research Institute Of Xi'an University Of Electronic Science And Technology
Xidian University
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 Hangzhou Research Institute Of Xi'an University Of Electronic Science And Technology, Xidian UniversityfiledCriticalHangzhou Research Institute Of Xi'an University Of Electronic Science And Technology
Priority to CN202310484598.5ApriorityCriticalpatent/CN116576827A/en
Publication of CN116576827ApublicationCriticalpatent/CN116576827A/en
Pendinglegal-statusCriticalCurrent

Links

Classifications

Landscapes

Abstract

Translated fromChinese

本发明公开了一种基于偏振三维成像的遥感地形三维测绘方法,包括:采集不同目标地形区域的偏振图像;计算目标地形区域表面的偏振度以及入射光的入射角、方位角;构建法向量与入射光的入射角、方位角之间的第一映射关系,获得第一法向量信息;计算太阳的天顶角与方位角,构建光源矢量与天顶角、方位角之间的第二映射关系,获得目标地形区域表面的光源矢量;构建目标地形区域表面的辐射模型,并计算先验法向量信息;基于先验法向量信息对第一法向量信息进行修正,得到修正法向量信息;根据修正法向量信息和预设角度下不同目标地形区域的偏振图像,进行目标区域的三维场景重建。本发明可以实现大场景下的卫星遥感地形三维重建,提高了运算速度。

The invention discloses a remote sensing terrain three-dimensional mapping method based on polarization three-dimensional imaging, which includes: collecting polarization images of different target terrain areas; calculating the degree of polarization of the surface of the target terrain area and the incident angle and azimuth angle of incident light; constructing a normal vector and The first mapping relationship between the incident angle and azimuth angle of the incident light to obtain the first normal vector information; calculate the zenith angle and azimuth angle of the sun, and construct the second mapping relationship between the light source vector and the zenith angle and azimuth angle , to obtain the light source vector on the surface of the target terrain area; construct the radiation model of the surface of the target terrain area, and calculate the prior normal vector information; based on the prior normal vector information, correct the first normal vector information to obtain the corrected normal vector information; Normal vector information and polarization images of different target terrain areas at preset angles are used to reconstruct the 3D scene of the target area. The invention can realize the three-dimensional reconstruction of the satellite remote sensing topography in a large scene, and improves the operation speed.

Description

Translated fromChinese
一种基于偏振三维成像的遥感地形三维测绘方法A remote sensing terrain 3D mapping method based on polarization 3D imaging

技术领域Technical Field

本发明属于遥感成像技术领域,具体涉及一种基于偏振三维成像的遥感地形三维测绘方法。The invention belongs to the technical field of remote sensing imaging, and in particular relates to a remote sensing terrain three-dimensional mapping method based on polarization three-dimensional imaging.

背景技术Background Art

航天光学遥感技术通常指利用卫星、航天飞机等在高空对地面上的目标进行远距离的探测,通常卫星平台搭载着空间相机或成像光谱仪等设备来收集目标所反射或发射出的辐射信息,利用所捕获到的信息进行光电转换、信息处理等来实现目标的监测。同时,利用卫星遥感影像的优势,在国土资源统计与监测、灾情监测与城市变迁等领域有着重要应用。随着科技水平的飞速发展,原始的二维遥感影像已无法满足应用需求,因此,需要寻求一种三维成像技术来对地形建立三维重建模型,而三维重建的主要目标是通过计算获取目标区域地形的数字高程信息。Aerospace optical remote sensing technology usually refers to the use of satellites, space shuttles, etc. to detect targets on the ground from a distance at high altitudes. Usually, the satellite platform is equipped with space cameras or imaging spectrometers and other equipment to collect radiation information reflected or emitted by the target, and use the captured information for photoelectric conversion and information processing to achieve target monitoring. At the same time, the advantages of satellite remote sensing images have important applications in the fields of land and resources statistics and monitoring, disaster monitoring, and urban changes. With the rapid development of science and technology, the original two-dimensional remote sensing images can no longer meet application requirements. Therefore, it is necessary to seek a three-dimensional imaging technology to establish a three-dimensional reconstruction model for the terrain. The main goal of three-dimensional reconstruction is to obtain digital elevation information of the terrain of the target area through calculation.

目前,常见的三维重建方法是利用航空航天影像结合摄影测量理论,通过获取立体像对进行三维重建,其中,遥感卫星的立体像对是指对相同目标区域从不同的角度进行拍摄,获取到两幅或者多幅图像,再利用模拟人眼的立体视觉,最终实现目标地形区域的三维重建。通过确定立体像对的内、外方位元素,结合共线条件以及两幅或多幅匹配图像(通常是同一卫星对同一目标区域不同角度所拍摄得到的多幅,具有一定重叠区域的图像)上的像点坐标进行影像匹配,得到立体像对上的同名像点,再进一步构建出两条同名射线,通过两条同名射线的交点即可计算出对应地面点的三维坐标,最终实现基于立体像对的三维重建。At present, the common 3D reconstruction method is to use aerospace images combined with photogrammetry theory to obtain stereo image pairs for 3D reconstruction. Among them, the stereo image pairs of remote sensing satellites refer to the same target area taken from different angles to obtain two or more images, and then use the stereoscopic vision of the simulated human eye to finally achieve 3D reconstruction of the target terrain area. By determining the internal and external orientation elements of the stereo image pair, combining the collinearity conditions and the image point coordinates on two or more matching images (usually multiple images with a certain overlapping area taken by the same satellite at different angles of the same target area), image matching is performed to obtain the same-name image points on the stereo image pair, and then further construct two same-name rays. The three-dimensional coordinates of the corresponding ground point can be calculated through the intersection of the two same-name rays, and finally 3D reconstruction based on stereo image pairs is achieved.

由于卫星遥感图像地形多为山地且起伏较大,影像匹配得到的同名点数量往往会存在极大的限制,导致最终重建效果受到影响;另外,虽然影像匹配方法较多,但其往往存在对噪声与光照影响精度与算法运行速度上无法同时兼顾的缺陷;再者,由于基于立体像对遥感影像三维重建技术匹配时往往需要多幅图片中包含相同的重叠区域,因此对于大范围大场景的目标区域拼接往往较为繁琐复杂,无法高效精确实现。Since the terrain of satellite remote sensing images is mostly mountainous and undulating, the number of same-name points obtained by image matching is often greatly limited, which affects the final reconstruction effect. In addition, although there are many image matching methods, they often have the defects of being unable to take into account both the accuracy affected by noise and illumination and the algorithm running speed. Furthermore, since the three-dimensional reconstruction technology based on stereo pair remote sensing images often requires multiple pictures to contain the same overlapping areas, the stitching of target areas in large areas and large scenes is often cumbersome and complicated, and cannot be achieved efficiently and accurately.

发明内容Summary of the invention

为了解决现有技术中存在的上述问题,本发明提供了一种基于偏振三维成像的遥感地形三维测绘方法。本发明要解决的技术问题通过以下技术方案实现:In order to solve the above problems existing in the prior art, the present invention provides a remote sensing terrain three-dimensional mapping method based on polarization three-dimensional imaging. The technical problem to be solved by the present invention is achieved through the following technical solutions:

本发明提供一种基于偏振三维成像的遥感地形三维测绘方法,包括:The present invention provides a remote sensing terrain three-dimensional mapping method based on polarization three-dimensional imaging, comprising:

采集预设角度下不同目标地形区域的偏振图像;Collect polarization images of different target terrain areas at preset angles;

针对每个目标地形区域,计算目标地形区域表面的偏振度,并根据所述偏振度计算所述目标地形区域表面的入射光的入射角θ和方位角For each target terrain area, the polarization degree of the surface of the target terrain area is calculated, and the incident angle θ and azimuth angle of the incident light on the surface of the target terrain area are calculated according to the polarization degree.

根据所述入射角θ和所述方位角构建所述目标地形区域表面的法向量与所述入射角θ及所述方位角之间的第一映射关系,获得第一法向量信息;According to the incident angle θ and the azimuth angle Construct the normal vector of the target terrain area surface and the incident angle θ and the azimuth angle A first mapping relationship between the two, obtaining first normal vector information;

计算太阳的天顶角θz与方位角θs,构建光源矢量与所述天顶角θz及所述方位角θs之间的第二映射关系,获得所述目标地形区域表面的光源矢量;Calculating the sun's zenith angleθz and azimuth angleθs , constructing a second mapping relationship between the light source vector and the zenith angleθz and the azimuth angleθs , and obtaining the light source vector on the surface of the target terrain area;

构建所述目标地形区域表面的辐射模型,并结合所述目标地形区域表面的光源矢量计算先验法向量信息;Constructing a radiation model of the surface of the target terrain area, and calculating prior normal vector information in combination with the light source vector of the surface of the target terrain area;

基于所述先验法向量信息对所述第一法向量信息进行修正,得到所述目标地形区域表面的修正法向量信息;Correcting the first normal vector information based on the prior normal vector information to obtain corrected normal vector information of the surface of the target terrain area;

根据所述修正法向量信息和所述预设角度下不同目标地形区域的偏振图像,进行目标区域的三维场景重建。The three-dimensional scene of the target area is reconstructed according to the corrected normal vector information and the polarization images of different target terrain areas at the preset angles.

在本发明的一个实施例中,所述预设角度包括0°、45°、90°和135°;In one embodiment of the present invention, the preset angles include 0°, 45°, 90° and 135°;

针对每个目标地形区域,计算目标地形区域表面的偏振度,并根据所述偏振度计算所述目标地形区域表面的入射光的入射角θ和方位角的步骤,包括:For each target terrain area, the polarization degree of the surface of the target terrain area is calculated, and the incident angle θ and azimuth angle of the incident light on the surface of the target terrain area are calculated according to the polarization degree. The steps include:

针对每个目标地形区域,获取采集预设角度下的偏振图像时通过0°偏光片的光强I00、通过45°偏光片的光强I+450、通过90°偏光片的光强I900和通过135°方向线偏振器的光强I-45°,并计算斯托克斯参数S0、S1和S2For each target terrain area, the light intensity I00 passing through the 0° polarizer, the light intensity I+450 passing through the 45° polarizer, the light intensity I900 passing through the 90° polarizer, and the light intensity I-45° passing through the 135° linear polarizer are obtained when collecting polarization images at a preset angle, and the Stokes parameters S0 , S1 and S2 are calculated:

S0=I+I90°S0 =I +I90° ;

S1=I-I90°S1 =I -I90° ;

S2=I+45°-I-45°S2 =I+45° -I-45° ;

根据所述斯托克斯参数S0、S1和S2,计算目标地形区域表面的偏振度:According to the Stokes parameters S0 , S1 and S2 , the polarization degree of the surface of the target terrain area is calculated:

根据所述斯托克斯参数S0、S1及所述偏振度,计算目标地形区域表面的入射光的入射角θ和方位角According to the Stokes parameters S0 , S1 and the polarization degree, the incident angle θ and the azimuth angle of the incident light on the surface of the target terrain area are calculated.

在本发明的一个实施例中,根据所述斯托克斯参数S0、S1及所述偏振度,计算目标地形区域表面的入射光的入射角θ和方位角的步骤,包括:In one embodiment of the present invention, the incident angle θ and the azimuth angle of the incident light on the surface of the target terrain area are calculated according to the Stokes parameters S0 , S1 and the polarization degree. The steps include:

根据所述斯托克斯参数S0、S1计算目标地形区域表面的入射光的方位角Calculate the azimuth of the incident light on the surface of the target terrain area according to the Stokes parameters S0 and S1

根据所述偏振度DOLP计算入射光的入射角θ:The incident angle θ of the incident light is calculated according to the degree of polarization DOLP:

式中,n表示目标地形区域表面的折射率,n=1.5。Wherein, n represents the refractive index of the surface of the target terrain area, n=1.5.

在本发明的一个实施例中,根据所述入射角θ和所述方位角构建所述目标地形区域表面的法向量与所述入射角θ及所述方位角之间的第一映射关系,获得第一法向量信息的步骤,包括:In one embodiment of the present invention, according to the incident angle θ and the azimuth angle Construct the normal vector of the target terrain area surface and the incident angle θ and the azimuth angle The step of obtaining first normal vector information according to the first mapping relationship between the two comprises:

根据所述入射角θ和所述方位角构建所述目标地形区域表面的法向量与所述入射角θ及所述方位角之间的第一映射关系:According to the incident angle θ and the azimuth angle Construct the normal vector of the target terrain area surface and the incident angle θ and the azimuth angle The first mapping relationship between:

式中,nx、ny、nz分别表示所述目标地形区域表面的法向量在x方向、y方向和z方向上的分量,其中,z方向与所述目标地形区域表面垂直,x方向与y方向垂直且二者所在平面与所述z方向垂直;Wherenx ,ny andnz represent the surface normal vectors of the target terrain area respectively. Components in the x-direction, y-direction and z-direction, wherein the z-direction is perpendicular to the surface of the target terrain area, the x-direction is perpendicular to the y-direction and the plane where the two directions are located is perpendicular to the z-direction;

根据所述第一映射关系计算第一法向量信息D。The first normal vector information D is calculated according to the first mapping relationship.

在本发明的一个实施例中,计算太阳的天顶角θz与方位角θs,构建光源矢量与所述天顶角θz及所述方位角θs之间的第二映射关系,获得所述目标地形区域表面的光源矢量的步骤,包括:In one embodiment of the present invention, the steps of calculating the sun's zenith angleθz and azimuth angleθs , constructing a second mapping relationship between a light source vector and the zenith angleθz and the azimuth angleθs , and obtaining the light source vector on the surface of the target terrain area include:

根据所述目标地形区域的被拍摄时间,计算太阳赤纬角:According to the time when the target terrain area was photographed, the solar declination angle is calculated:

式中,N表示积日数;In the formula, N represents the cumulative number of days;

根据所述太阳赤纬角δ、目标地形区域的纬度φ和太阳时角ω计算太阳的高度角:The solar altitude angle is calculated according to the solar declination angle δ, the latitude φ of the target terrain area and the solar hour angle ω:

sinhs=sinδ·sinφ+cosδ·cosφ·cosω;sinhs = sinδ·sinφ+cosδ·cosφ·cosω;

式中,ω=(12-t)×15°,t表示目标地形区域被拍摄时的当地时间;Where, ω = (12-t) × 15°, t represents the local time when the target terrain area is photographed;

根据太阳赤纬角δ、目标地形区域的经纬度信息和太阳高度角hs计算太阳的方位角:The azimuth of the sun is calculated based on the solar declination angle δ, the longitude and latitude information of the target terrain area, and the solar altitude anglehs :

式中,φ表示目标地形区域的纬度;Where φ represents the latitude of the target terrain area;

根据太阳的高度角hs计算太阳的天顶角θzCalculate the solar zenith angleθz based on the solar altitude anglehs :

θz=90°-hsθz = 90°-hs ;

构建光源矢量与所述天顶角θz及所述方位角θs之间的第二映射关系:Construct a second mapping relationship between the light source vector and the zenith angleθz and the azimuth angleθs :

根据所述第二映射关系,获得所述目标地形区域表面的光源矢量(ps,qs,-1)。According to the second mapping relationship, the light source vector (ps ,qs , -1) of the surface of the target terrain area is obtained.

在本发明的一个实施例中,构建所述目标地形区域表面的辐射模型,并结合所述目标地形区域表面的光源矢量计算先验法向量信息的步骤,包括:In one embodiment of the present invention, the step of constructing a radiation model of the surface of the target terrain area and calculating the prior normal vector information in combination with the light source vector of the surface of the target terrain area includes:

根据朗伯反射模型、所述目标地形区域表面的光源矢量(ps,qs,-1)和光滑性约束,构建所述目标地形区域表面的辐射模型:According to the Lambertian reflection model, the light source vector (ps ,qs , -1) on the surface of the target terrain area and the smoothness constraint, a radiation model of the surface of the target terrain area is constructed:

所述光滑性约束为:The smoothness constraint is:

式中,I(x',y')表示目标地形区域(x',y')处在0°和90°下的偏振图像的光源强度之和,ρ表示预设的偏振图像反照率,px,py,qx,qy分别为p(x',y')、q(x',y')对应于x方向和y方向上的导数;Where I(x',y') represents the sum of the light source intensities of the polarization images of the target terrain area (x',y') at 0° and 90°, ρ represents the preset polarization image albedo,px ,py ,qx ,qy are the derivatives of p(x',y') and q(x',y') in the x and y directions respectively;

根据所述辐射模型及光滑性约束,计算得到所述目标地形区域表面(x',y')处的先验法向量信息Y=(p(x',y'),q(x',y'),-1)。According to the radiation model and the smoothness constraint, the prior normal vector information Y=(p(x',y'),q(x',y'),-1) at the surface (x',y') of the target terrain area is calculated.

在本发明的一个实施例中,基于所述先验法向量信息对所述第一法向量信息进行修正,得到所述目标地形区域表面修正后的法向量的步骤,包括:In one embodiment of the present invention, the step of correcting the first normal vector information based on the prior normal vector information to obtain the corrected normal vector of the surface of the target terrain area includes:

按照如下公式计算二元算子The binary operator is calculated according to the following formula

式中,T表示非线性算子;In the formula, T represents the nonlinear operator;

基于所述二元算子对第一法向量信息D进行修正:Based on the binary operator Correct the first normal vector information D:

式中,Ncor表示所述目标地形区域表面的修正法向量信息。Wherein, Ncor represents the corrected normal vector information of the surface of the target terrain area.

在本发明的一个实施例中,所述不同目标区域包括第一目标区域和第二目标区域,所述第一目标区域与所述第二目标区域至少存在部分重叠区域;In one embodiment of the present invention, the different target areas include a first target area and a second target area, and the first target area and the second target area have at least a partial overlapping area;

根据所述修正法向量信息和所述预设角度下不同目标地形区域的偏振图像,进行目标区域的三维场景重建的步骤,包括:The step of reconstructing a three-dimensional scene of the target area according to the corrected normal vector information and the polarization images of different target terrain areas at the preset angles includes:

从所述第一目标区域的任意一幅偏振图像中选取一特征点,所述特征点位于重叠区域内;Selecting a feature point from any one of the polarization images of the first target area, wherein the feature point is located in the overlapping area;

获取所述特征点在所述第二目标区域中的位置信息,并建立所述第一目标区域的偏振图像与所述第二目标区域的偏振图像之间的变换关系;Acquire position information of the feature point in the second target area, and establish a transformation relationship between the polarization image of the first target area and the polarization image of the second target area;

根据所述变换关系拼接所述第一目标区域的偏振图像与所述第二目标区域的偏振图像;splicing the polarization image of the first target area and the polarization image of the second target area according to the transformation relationship;

根据拼接后的图像、第一目标区域表面的修正法向量信息和第二目标区域表面的修正法向量信息进行三维重建。Three-dimensional reconstruction is performed based on the stitched image, the corrected normal vector information of the surface of the first target area, and the corrected normal vector information of the surface of the second target area.

在本发明的一个实施例中,根据拼接后的图像、第一目标区域表面的修正法向量信息和第二目标区域表面的修正法向量信息,按照如下公式进行三维重建:In one embodiment of the present invention, three-dimensional reconstruction is performed according to the following formula based on the stitched image, the corrected normal vector information of the surface of the first target area, and the corrected normal vector information of the surface of the second target area:

式中,F{}和F-1{}分别表示离散傅里叶变换和傅里叶逆变换,M、N分别表示拼接后的图像中像素点的行数及列数,u、ν分别表示离散微算子在x和y方向上的傅里叶系数。Where F{} and F-1 {} represent discrete Fourier transform and inverse Fourier transform, respectively; M and N represent the number of rows and columns of pixels in the spliced image, respectively; u and ν represent the Fourier coefficients of the discrete differential operator in the x and y directions, respectively.

与现有技术相比,本发明的有益效果在于:Compared with the prior art, the present invention has the following beneficial effects:

本发明提供一种基于偏振三维成像的遥感地形三维测绘方法,通过利用偏振图像中的信息进行三维重建,避免了现有技术中基于立体像对遥感影像三维重建时无法避开的图像匹配的过程,提高了运算速度;与此同时,本发明利用偏振三维成像进行逐像素拼接,可以快速高效地实现大场景下的卫星遥感地形三维重建,解决了现有基于立体像对遥感影像三维重建过程中需要多幅图像且只能对重叠区域进行三维重建的问题。另外,基于偏振三维重建的遥感影像地形恢复技术具有被动式、高精度、远距离、干扰较小的优势,使得利用卫星遥感图像进行三维重建的精度进一步提升成为可能。The present invention provides a remote sensing terrain three-dimensional mapping method based on polarization three-dimensional imaging. By using the information in the polarization image for three-dimensional reconstruction, the image matching process that cannot be avoided in the three-dimensional reconstruction of remote sensing images based on stereo pairs in the prior art is avoided, and the operation speed is improved; at the same time, the present invention uses polarization three-dimensional imaging for pixel-by-pixel splicing, which can quickly and efficiently realize the three-dimensional reconstruction of satellite remote sensing terrain in large scenes, and solves the problem that multiple images are required in the three-dimensional reconstruction of remote sensing images based on stereo pairs and only overlapping areas can be reconstructed in three dimensions. In addition, the remote sensing image terrain restoration technology based on polarization three-dimensional reconstruction has the advantages of being passive, high-precision, long-distance, and less interference, making it possible to further improve the accuracy of three-dimensional reconstruction using satellite remote sensing images.

以下将结合附图及实施例对本发明做进一步详细说明。The present invention will be further described in detail below with reference to the accompanying drawings and embodiments.

附图说明BRIEF DESCRIPTION OF THE DRAWINGS

图1是本发明实施例提供的基于偏振三维成像的遥感地形三维测绘方法的一种流程图;FIG1 is a flow chart of a method for three-dimensional surveying and mapping of remote sensing terrain based on polarization three-dimensional imaging provided by an embodiment of the present invention;

图2是本发明实施例提供的基于偏振三维成像的遥感地形三维测绘方法的另一种流程图;2 is another flow chart of a method for three-dimensional remote sensing terrain mapping based on polarization three-dimensional imaging provided by an embodiment of the present invention;

图3是本发明实施例提供的目标地形区域表面法向量与入射光的方位角和入射角的示意图;3 is a schematic diagram of a surface normal vector of a target terrain area and an azimuth and incident angle of incident light provided by an embodiment of the present invention;

图4是本发明实施例提供的目标地形区域表面法向量与太阳天顶角和方位角的示意图;4 is a schematic diagram of the surface normal vector of the target terrain area and the solar zenith angle and azimuth angle provided by an embodiment of the present invention;

图5是本发明实施例提供的逐像素拼接的示意图。FIG. 5 is a schematic diagram of pixel-by-pixel stitching provided by an embodiment of the present invention.

具体实施方式DETAILED DESCRIPTION

下面结合具体实施例对本发明做进一步详细的描述,但本发明的实施方式不限于此。The present invention is further described in detail below with reference to specific embodiments, but the embodiments of the present invention are not limited thereto.

图1是本发明实施例提供的基于偏振三维成像的遥感地形三维测绘方法的一种流程图,图2是本发明实施例提供的基于偏振三维成像的遥感地形三维测绘方法的另一种流程图。如图1-2所示,本发明提供一种基于偏振三维成像的遥感地形三维测绘方法,包括:FIG1 is a flow chart of a method for three-dimensional surveying and mapping of remote sensing terrain based on polarization three-dimensional imaging provided by an embodiment of the present invention, and FIG2 is another flow chart of a method for three-dimensional surveying and mapping of remote sensing terrain based on polarization three-dimensional imaging provided by an embodiment of the present invention. As shown in FIG1-2, the present invention provides a method for three-dimensional surveying and mapping of remote sensing terrain based on polarization three-dimensional imaging, comprising:

S1、采集预设角度下不同目标地形区域的偏振图像;S1, collecting polarization images of different target terrain areas at preset angles;

S2、针对每个目标地形区域,计算目标地形区域表面的偏振度,并根据偏振度计算目标地形区域表面的入射光的入射角θ和方位角S2. For each target terrain area, calculate the polarization degree of the surface of the target terrain area, and calculate the incident angle θ and azimuth angle of the incident light on the surface of the target terrain area according to the polarization degree.

S3、根据入射角θ和方位角构建目标地形区域表面的法向量与入射角θ及方位角之间的第一映射关系,获得第一法向量信息;S3, according to the incident angle θ and azimuth Construct the normal vector, incident angle θ and azimuth of the target terrain area surface A first mapping relationship between the two, obtaining first normal vector information;

S4、计算太阳的天顶角θz与方位角θs,构建光源矢量与天顶角θz及方位角θs之间的第二映射关系,获得目标地形区域表面的光源矢量;S4, calculating the sun's zenith angleθz and azimuth angleθs , constructing a second mapping relationship between the light source vector and the zenith angleθz and azimuth angleθs , and obtaining the light source vector on the surface of the target terrain area;

S5、构建目标地形区域表面的辐射模型,并结合目标地形区域表面的光源矢量计算先验法向量信息;S5, constructing a radiation model of the surface of the target terrain area, and calculating prior normal vector information in combination with the light source vector of the surface of the target terrain area;

S6、基于先验法向量信息对第一法向量信息进行修正,得到目标地形区域表面的修正法向量信息;S6. Correcting the first normal vector information based on the prior normal vector information to obtain corrected normal vector information of the surface of the target terrain area;

S7、根据修正法向量信息和预设角度下不同目标地形区域的偏振图像,进行目标区域的三维场景重建。S7. Reconstruct a three-dimensional scene of the target area according to the corrected normal vector information and polarization images of different target terrain areas at preset angles.

本实施例中,预设角度包括0°、45°、90°和135°。具体地,可以利用卫星上所搭载的成像探测器CMOS相机采集目标地形区域的反射光,得到每个目标地形区域在0°、45°、90°和135°的偏振图像。In this embodiment, the preset angles include 0°, 45°, 90° and 135°. Specifically, the imaging detector CMOS camera carried by the satellite can be used to collect the reflected light of the target terrain area to obtain the polarization images of each target terrain area at 0°, 45°, 90° and 135°.

上述步骤S2中,针对每个目标地形区域,计算目标地形区域表面的偏振度,并根据偏振度计算目标地形区域表面的入射光的入射角θ和方位角的步骤,包括:In the above step S2, for each target terrain area, the polarization degree of the surface of the target terrain area is calculated, and the incident angle θ and azimuth angle of the incident light on the surface of the target terrain area are calculated according to the polarization degree. The steps include:

S201、针对每个目标地形区域,获取采集预设角度下的偏振图像时通过0°偏光片的光强I、通过45°偏光片的光强I+45°、通过90°偏光片的光强I90°和通过135°方向线偏振器的光强I-45°,并计算斯托克斯参数S0、S1和S2S201, for each target terrain area, obtain the light intensity I passing through the 0° polarizer, the light intensity I+45° passing through the 45° polarizer, the light intensity I90° passing through the 90° polarizer, and the light intensity I-45° passing through the 135° linear polarizer when collecting the polarization image at a preset angle, and calculate the Stokes parameters S0 , S1 and S2 :

S0=I+I90°S0 =I +I90° ;

S1=I-I90°S1 =I -I90° ;

S2=I+45°-I-45°S2 =I+45° -I-45° ;

S202、根据斯托克斯参数S0、S1和S2,计算目标地形区域表面的偏振度:S202. Calculate the polarization degree of the target terrain area surface according to Stokes parameters S0 , S1 and S2 :

S203、根据斯托克斯参数S0、S1及所述偏振度,计算目标地形区域表面的入射光的入射角θ和方位角S203, calculating the incident angle θ and azimuth angle of the incident light on the surface of the target terrain area according to the Stokes parameters S0 , S1 and the polarization degree

应当理解,利用两个或两个以上具有不同偏振态的光子可以描述光的偏振状态,偏振状态及其和由斯托克参数描述,其分量为四个斯托克斯参数:It should be understood that the polarization state of light can be described using two or more photons with different polarization states, and the polarization state and its sum are described by Stokes parameters, whose components are four Stokes parameters:

其中:in:

在实际应用中,斯托克斯参数可以通过偏振器测量获得:In practical applications, the Stokes parameters can be obtained by polarizer measurement:

上式中,I表示通过0°偏光片的光强,I90°是通过90°偏光片的光强,I+45°是与水平正方向成45°偏光片的光强,I-45°是135°方向的线偏振器得到的光强,IRH和ILH分别表示右旋圆偏振光分量和左旋圆偏振光分量。In the above formula, I represents the light intensity passing through the 0° polarizer, I90° is the light intensity passing through the 90° polarizer, I+45° is the light intensity of the polarizer at an angle of 45° to the horizontal positive direction, I-45° is the light intensity obtained by the linear polarizer in the 135° direction, and IRH and ILH represent the right-handed circularly polarized light component and the left-handed circularly polarized light component, respectively.

进一步地,利用斯托克斯参数可以直接确定任何光信号的偏振度(DoP),偏振度表征偏振分量的比例:Furthermore, the Stokes parameters allow the degree of polarization (DoP) of any optical signal to be determined directly. The degree of polarization characterizes the ratio of the polarization components:

由于本实施例中采集偏振图像时的光束是线性偏振,因而右旋圆偏振光分量与左旋圆偏振光分量均为零,则目标地形区域表面的偏振度可表示为:Since the light beam when collecting the polarization image in this embodiment is linearly polarized, the right-handed circularly polarized light component and the left-handed circularly polarized light component are both zero, and the polarization degree of the surface of the target terrain area can be expressed as:

式中,DOLP一般称为线性偏振度,表征了光束的线性偏振度。In the formula, DOLP is generally referred to as the degree of linear polarization, which characterizes the linear polarization degree of the light beam.

上述步骤S203中,根据斯托克斯参数S0、S1及偏振度,计算目标地形区域表面的入射光的入射角θ和方位角的步骤,包括:In the above step S203, the incident angle θ and the azimuth angle of the incident light on the surface of the target terrain area are calculated according to the Stokes parameters S0 , S1 and the polarization degree. The steps include:

根据斯托克斯参数S0、S1计算目标地形区域表面的入射光的方位角Calculate the azimuth of the incident light on the surface of the target terrain area based on the Stokes parameters S0 and S1

根据偏振度DOLP计算入射光的入射角θ:Calculate the incident angle θ of the incident light according to the degree of polarization DOLP:

式中,n表示目标地形区域表面的折射率,n=1.5。Wherein, n represents the refractive index of the surface of the target terrain area, n=1.5.

图3是本发明实施例提供的目标地形区域表面法向量与入射光的方位角和入射角的示意图。请参见图3,上述步骤S3中,根据入射角θ和所述方位角构建目标地形区域表面的法向量与入射角θ及方位角之间的第一映射关系,获得第一法向量信息的步骤,包括:FIG3 is a schematic diagram of the surface normal vector of the target terrain area and the azimuth and incident angle of the incident light provided by an embodiment of the present invention. Referring to FIG3, in the above step S3, according to the incident angle θ and the azimuth Construct the normal vector, incident angle θ and azimuth of the target terrain area surface The step of obtaining first normal vector information according to the first mapping relationship between the two comprises:

根据入射角θ和方位角构建目标地形区域表面的法向量与入射角θ及方位角之间的第一映射关系:According to the incident angle θ and azimuth Construct the normal vector, incident angle θ and azimuth of the target terrain area surface The first mapping relationship between:

式中,nx、ny、nz分别表示目标地形区域表面的法向量在x方向、y方向和z方向上的分量,其中,z方向与目标地形区域表面垂直,x方向与y方向垂直且二者所在平面与z方向垂直;Wherenx ,ny , andnz represent the surface normal vectors of the target terrain area. Components in the x-, y- and z-directions, where the z-direction is perpendicular to the surface of the target terrain area, the x-direction is perpendicular to the y-direction and the plane where the two directions are located is perpendicular to the z-direction;

根据第一映射关系计算第一法向量信息D。The first normal vector information D is calculated according to the first mapping relationship.

需要说明的是,在光强一定时存在对应多个偏振器旋转角度的情况,因此通过光强来计算入射光的方位角时会存在一个180°的多值性问题,也就无法准确获取目标地形区域表面的法向量信息。It should be noted that when the light intensity is constant, there are multiple corresponding polarizer rotation angles, so the azimuth of the incident light is calculated by the light intensity. There will be a 180° multi-value problem, and it is impossible to accurately obtain the normal vector information of the surface of the target terrain area.

基于上述分析,本实施例还需要进一步对根据第一映射关系计算获得的第一法向量信息D进行修正。Based on the above analysis, this embodiment also needs to further correct the first normal vector information D calculated according to the first mapping relationship.

图4是本发明实施例提供的目标地形区域表面法向量与太阳天顶角和方位角的示意图。如图4所示,上述步骤S4中,计算太阳的天顶角θz与方位角θs,构建光源矢量与天顶角θz及方位角θs之间的第二映射关系,获得目标地形区域表面的光源矢量的步骤,包括:FIG4 is a schematic diagram of the surface normal vector of the target terrain area and the solar zenith angle and azimuth angle provided by an embodiment of the present invention. As shown in FIG4, in the above step S4, the step of calculating the solar zenith angleθz and azimuth angleθs , constructing a second mapping relationship between the light source vector and the zenith angleθz and azimuth angleθs , and obtaining the light source vector on the surface of the target terrain area includes:

S401、根据目标地形区域的被拍摄时间,计算太阳赤纬角:S401, calculating the solar declination angle according to the time when the target terrain area was photographed:

式中,N表示积日数;In the formula, N represents the cumulative number of days;

S402、根据太阳赤纬角δ、目标地形区域的纬度φ和太阳时角ω计算太阳的高度角:S402, calculating the solar altitude angle according to the solar declination angle δ, the latitude φ of the target terrain area and the solar hour angle ω:

sinhs=sinδ·sinφ+cosδ·cosφ·cosω;sinhs = sinδ·sinφ+cosδ·cosφ·cosω;

式中,ω=(12-t)×15°,t表示目标地形区域被拍摄时的当地时间;Where, ω = (12-t) × 15°, t represents the local time when the target terrain area is photographed;

S403、根据太阳赤纬角δ、目标地形区域的经纬度信息和太阳高度角hs计算太阳的方位角:S403, calculating the azimuth of the sun according to the solar declination angle δ, the longitude and latitude information of the target terrain area and the solar altitude anglehs :

式中,φ表示目标地形区域的纬度;Where φ represents the latitude of the target terrain area;

S404、根据太阳的高度角hs计算太阳的天顶角θzS404. Calculate the solar zenith angleθz according to the solar altitude anglehs :

θz=90°-hsθz = 90°-hs ;

S405、构建光源矢量与天顶角θz及方位角θs之间的第二映射关系:S405, construct a second mapping relationship between the light source vector and the zenith angleθz and the azimuth angleθs :

S406、根据第二映射关系,获得目标地形区域表面的光源矢量(ps,qs,-1)。S406. Obtain the light source vector (ps ,qs , -1) on the surface of the target terrain area according to the second mapping relationship.

具体而言,首先根据目标地形区域的经纬度信息与被拍摄时间,计算太阳的方位角θs与高度角hs为:Specifically, firstly, according to the latitude and longitude information of the target terrain area and the shooting time, the azimuth angleθs and the altitude anglehs of the sun are calculated as:

sinhs=sinδ·sinφ+cosδ·cosφ·cosω;sinhs = sinδ·sinφ+cosδ·cosφ·cosω;

式中,δ为太阳中心点和地球中心点的连线与赤道所在平面之间所夹的锐角,又称为太阳赤纬角,φ表示目标地形区域的纬度,ω表示太阳时角。Where δ is the acute angle between the line connecting the center of the sun and the center of the earth and the plane where the equator is located, also known as the solar declination angle, φ represents the latitude of the target terrain area, and ω represents the solar hour angle.

其中,赤纬角δ可由下式解得:The declination angle δ can be solved by the following formula:

N表示从拍摄时间所在年份的1月1日到拍摄当天的日期之间的天数,又称积日数。N represents the number of days from January 1st of the year in which the shooting was conducted to the date of the shooting, also known as the accumulated days.

太阳时角ω的计算公式为:The calculation formula of solar hour angle ω is:

ω=(12-t)×15°ω=(12-t)×15°

t表示目标地形区域被拍摄时的当地时间。t represents the local time when the target terrain area is photographed.

请参见图4,太阳的高度角hs与天顶角θz互为余角,因此太阳的天顶角可由高度角hs确定:Please refer to Figure 4. The sun's altitude anglehs and zenith angleθz are complementary angles, so the sun's zenith angle can be determined by the altitude anglehs :

θz=90°-hsθz =90°-hs .

上述步骤S5中,构建目标地形区域表面的辐射模型,并结合目标地形区域表面的光源矢量计算先验法向量信息的步骤,包括:In the above step S5, the step of constructing a radiation model of the surface of the target terrain area and calculating the prior normal vector information in combination with the light source vector of the surface of the target terrain area includes:

S501、根据朗伯反射模型、目标地形区域表面的光源矢量(ps,qs,-1)和光滑性约束,构建目标地形区域表面的辐射模型:S501, constructing a radiation model of the surface of the target terrain area according to the Lambertian reflection model, the light source vector (ps ,qs , -1) of the surface of the target terrain area and the smoothness constraint:

光滑性约束为:The smoothness constraint is:

式中,光滑性约束表示相邻像素点法向量方向接近,I(x',y')表示目标地形区域(x',y')处在0°和90°下的偏振图像的光源强度之和,ρ表示预设的偏振图像反照率,px,py,qx,qy分别为p(x',y')、q(x',y')对应于x方向和y方向上的导数;In the formula, the smoothness constraint Indicates that the directions of the normal vectors of adjacent pixels are close, I(x',y') represents the sum of the light source intensities of the polarization images of the target terrain area (x',y') at 0° and 90°, ρ represents the preset polarization image albedo,px ,py ,qx ,qy are the derivatives of p(x',y') and q(x',y') in the x and y directions respectively;

S502、根据辐射模型及光滑性约束,计算得到目标地形区域表面(x',y')处的先验法向量信息Y=(p(x',y'),q(x',y'),-1)。S502, according to the radiation model and the smoothness constraint, calculate and obtain the prior normal vector information Y=(p(x',y'),q(x',y'),-1) at the surface (x',y') of the target terrain area.

上述步骤S6中,基于先验法向量信息对第一法向量信息进行修正,得到目标地形区域表面的修正法向量信息的步骤,包括:In the above step S6, the step of correcting the first normal vector information based on the prior normal vector information to obtain the corrected normal vector information of the surface of the target terrain area includes:

按照如下公式计算二元算子The binary operator is calculated according to the following formula

式中,T表示非线性算子;In the formula, T represents the nonlinear operator;

基于所述二元算子对第一法向量信息D进行修正:Based on the binary operator Correct the first normal vector information D:

式中,Ncor表示目标地形区域表面的修正法向量信息。Where Ncor represents the corrected normal vector information of the surface of the target terrain area.

本实施例利用法向梯度场的先验信息所提供的目标地形区域表面在校正过程中的变化趋势,解决了方位角模糊所导致表面法向量不准确的问题,得到修正后的目标表面法线,进一步获得目标表面梯度场,然后可利用表面梯度场重建方法可以恢复目标表面形状。This embodiment uses the change trend of the target terrain area surface during the correction process provided by the prior information of the normal gradient field to solve the problem of inaccurate surface normal vector caused by azimuth ambiguity, obtains the corrected target surface normal, and further obtains the target surface gradient field. Then, the surface gradient field reconstruction method can be used to restore the target surface shape.

可选地,不同目标区域包括第一目标区域和第二目标区域,第一目标区域与第二目标区域至少存在部分重叠区域;Optionally, the different target areas include a first target area and a second target area, and the first target area and the second target area have at least a partial overlapping area;

上述步骤S7中,根据修正法向量信息和预设角度下不同目标地形区域的偏振图像,进行目标区域的三维场景重建的步骤,包括:In the above step S7, the step of reconstructing the three-dimensional scene of the target area according to the corrected normal vector information and the polarization images of different target terrain areas at preset angles includes:

S701、从第一目标区域的任意一幅偏振图像中选取一特征点,特征点位于重叠区域内;S701, selecting a feature point from any polarization image of the first target area, where the feature point is located in the overlapping area;

S702、获取特征点在第二目标区域中的位置信息,并建立第一目标区域的偏振图像与第二目标区域的偏振图像之间的变换关系;S702, obtaining position information of the feature point in the second target area, and establishing a transformation relationship between the polarization image of the first target area and the polarization image of the second target area;

S703、根据变换关系拼接第一目标区域的偏振图像与第二目标区域的偏振图像;S703, splicing the polarization image of the first target area and the polarization image of the second target area according to the transformation relationship;

S704、根据拼接后的图像、第一目标区域表面的修正法向量信息和第二目标区域表面的修正法向量信息进行三维重建。S704 , performing three-dimensional reconstruction according to the stitched image, the corrected normal vector information of the surface of the first target area, and the corrected normal vector information of the surface of the second target area.

对于实际中大场景大区域的目标区域地形三维重建时,由于相机视场有限,无法使一幅图像中包含所有目标区域,因此需要对多区域的图像进行拼接。步骤S701~S702中,首先从第一目标区域的任意一幅偏振图像中选取中一个特征点,然后结合这一特征点在第二目标区域的任意一幅偏振图像中的位置信息,建立两幅偏振图像之间的变换关系,从而实现图像的无缝拼接。同时,由于偏振三维成像是基于每一个独立的像素点进行的,因此可以同时实现大场景的逐像素拼接与重建。When reconstructing the terrain of a target area in a large scene or large area in practice, due to the limited field of view of the camera, it is impossible to include all target areas in one image, so it is necessary to stitch images of multiple areas. In steps S701 to S702, first select a feature point from any polarization image of the first target area, and then combine the position information of this feature point in any polarization image of the second target area to establish a transformation relationship between the two polarization images, thereby achieving seamless stitching of images. At the same time, since polarization three-dimensional imaging is based on each independent pixel, pixel-by-pixel stitching and reconstruction of large scenes can be achieved simultaneously.

图5是本发明实施例提供的逐像素拼接的示意图。如图5所示,以方格代表偏振图像中的像素、实心方格表示选取的特征点,其中,I(m,n)表示在第一目标区域的偏振图像中位于(m,n)处的像素点,I(i,j)表示第二目标区域的偏振图像中位于(i,j)处的像素点,结合两幅图像之间的变换关系以及修正法向量信息即可实现大场景下遥感地形的逐像素重建。Figure 5 is a schematic diagram of pixel-by-pixel stitching provided by an embodiment of the present invention. As shown in Figure 5, squares represent pixels in the polarization image, and solid squares represent selected feature points, where I(m,n) represents the pixel at (m,n) in the polarization image of the first target area, and I(i,j) represents the pixel at (i,j) in the polarization image of the second target area. By combining the transformation relationship between the two images and the corrected normal vector information, pixel-by-pixel reconstruction of remote sensing terrain in a large scene can be achieved.

具体而言,将目标曲面记作Z(x,y),利用可积函数的梯度表达形式来表示目标区域的表面法线:Specifically, the target surface is recorded as Z(x,y), and the surface normal of the target area is expressed using the gradient expression of the integrable function:

结合法向量与入射角θ及方位角之间的映射关系,联立可得:Combined normal vector With the incident angle θ and azimuth The mapping relationship between them can be obtained by combining:

定义函数u(Z)表示利用梯度场求解目标曲面形状与实际形状之间的误差:The function u(Z) is defined to represent the error between the target surface shape and the actual shape solved by the gradient field:

u(Z)=(Zx-p)2+(Zy-q)2u(Z)=(Zx -p)2+ (Zy -q)2

其中,Zx与Zy分别表示目标曲面Z(x,y)在x方向和y方向上的偏导数。WhereZx andZy represent the partial derivatives of the target surface Z(x,y) in the x-direction and y-direction respectively.

为了得到最真实的目标曲面Z(x,y),需要使得目标曲面上所有点的误差最小化,因此可定义目标代价函数f(Z)来进行优化处理:In order to obtain the most realistic target surface Z(x,y), it is necessary to minimize the errors of all points on the target surface. Therefore, the target cost function f(Z) can be defined for optimization:

f(Z)=∫∫u(Z)dxdy。f(Z)=∫∫u(Z)dxdy.

对于全局连续的法线梯度场,通常采用全局或局部积分的方法。局部积分方法容易被待处理数据中的噪声所干扰,而实际采集初始数据时往往会因为环境杂散光以及本身探测器的影响,存在不可避免的噪声影响,而全局积分方法对噪声具有更好的鲁棒性,重建的表面更光滑。因此,本实施例采用全局积分方法来求解目标区域曲面使得误差最小化,并利用Frankot Chellappa算法获得最终重建结果。For a globally continuous normal gradient field, a global or local integration method is usually used. The local integration method is easily disturbed by the noise in the data to be processed, and when the initial data is actually collected, there is often an inevitable noise influence due to the influence of ambient stray light and the detector itself, while the global integration method has better robustness to noise and the reconstructed surface is smoother. Therefore, this embodiment uses a global integration method to solve the target area surface to minimize the error, and uses the Frankot Chellappa algorithm to obtain the final reconstruction result.

Frankot Chellappa算法主要是通过在频率域空间中来重建出目标区域的表面函数,再通过傅里叶逆变换来实现目标的三维重建模型。可选地,按照如下公式进行三维重建:The Frankot Chellappa algorithm mainly reconstructs the surface function of the target area in the frequency domain space, and then realizes the three-dimensional reconstruction model of the target through the inverse Fourier transform. Optionally, the three-dimensional reconstruction is performed according to the following formula:

式中,F{}和F-1{}分别表示离散傅里叶变换和傅里叶逆变换,M、N分别表示拼接后的图像中像素点的行数及列数,u、ν分别表示离散微算子在x和y方向上的傅里叶系数。Where F{} and F-1 {} represent discrete Fourier transform and inverse Fourier transform, respectively; M and N represent the number of rows and columns of pixels in the spliced image, respectively; u and ν represent the Fourier coefficients of the discrete differential operator in the x and y directions, respectively.

通过上述各实施例可知,本发明的有益效果在于:It can be seen from the above embodiments that the beneficial effects of the present invention are:

本发明提供一种基于偏振三维成像的遥感地形三维测绘方法,通过利用偏振图像中的信息进行三维重建,避免了现有技术中基于立体像对遥感影像三维重建时无法避开的图像匹配的过程,提高了运算速度;与此同时,本发明利用偏振三维成像进行逐像素拼接,可以快速高效地实现大场景下的卫星遥感地形三维重建,解决了现有基于立体像对遥感影像三维重建过程中需要多幅图像且只能对重叠区域进行三维重建的问题。另外,基于偏振三维重建的遥感影像地形恢复技术具有被动式、高精度、远距离、干扰较小的优势,使得利用卫星遥感图像进行三维重建的精度进一步提升成为可能。The present invention provides a three-dimensional surveying and mapping method of remote sensing terrain based on polarization three-dimensional imaging. By using the information in the polarization image for three-dimensional reconstruction, the image matching process that cannot be avoided in the three-dimensional reconstruction of remote sensing images based on stereo pairs in the prior art is avoided, and the operation speed is improved; at the same time, the present invention uses polarization three-dimensional imaging for pixel-by-pixel splicing, which can quickly and efficiently realize the three-dimensional reconstruction of satellite remote sensing terrain in large scenes, and solves the problem that multiple images are required in the three-dimensional reconstruction process of remote sensing images based on stereo pairs, and only overlapping areas can be reconstructed in three dimensions. In addition, the remote sensing image terrain restoration technology based on polarization three-dimensional reconstruction has the advantages of being passive, high-precision, long-distance, and less interference, making it possible to further improve the accuracy of three-dimensional reconstruction using satellite remote sensing images.

在本发明的描述中,术语“第一”、“第二”仅用于描述目的,而不能理解为指示或暗示相对重要性或者隐含指明所指示的技术特征的数量。由此,限定有“第一”、“第二”的特征可以明示或者隐含地包括一个或者更多个该特征。在本发明的描述中,“多个”的含义是两个或两个以上,除非另有明确具体的限定。In the description of the present invention, the terms "first" and "second" are used for descriptive purposes only and should not be understood as indicating or implying relative importance or implicitly indicating the number of the indicated technical features. Thus, a feature defined as "first" or "second" may explicitly or implicitly include one or more of the features. In the description of the present invention, the meaning of "plurality" is two or more, unless otherwise clearly and specifically defined.

参考术语“一个实施例”、“一些实施例”、“示例”、“具体示例”、或“一些示例”等的描述意指结合该实施例或示例描述的具体特征、结构、材料或者特点包含于本发明的至少一个实施例或示例中。在本说明书中,对上述术语的示意性表述不必须针对的是相同的实施例或示例。而且,描述的具体特征、结构、材料或者特点可以在任何的一个或多个实施例或示例中以合适的方式结合。此外,本领域的技术人员可以将本说明书中描述的不同实施例或示例进行接合和组合。The description with reference to the terms "one embodiment", "some embodiments", "example", "specific example", or "some examples" etc. means that the specific features, structures, materials or characteristics described in conjunction with the embodiment or example are included in at least one embodiment or example of the present invention. In this specification, the schematic representations of the above terms do not necessarily refer to the same embodiment or example. Moreover, the specific features, structures, materials or characteristics described may be combined in any one or more embodiments or examples in a suitable manner. In addition, those skilled in the art may combine and combine different embodiments or examples described in this specification.

尽管在此结合各实施例对本申请进行了描述,然而,在实施所要求保护的本申请过程中,本领域技术人员通过查看所述附图、公开内容、以及所附权利要求书,可理解并实现所述公开实施例的其他变化。Although the present application is described herein in conjunction with various embodiments, in the process of implementing the present application for which protection is sought, a person skilled in the art may understand and implement other variations of the disclosed embodiments by reviewing the drawings, the disclosure, and the appended claims.

以上内容是结合具体的优选实施方式对本发明所作的进一步详细说明,不能认定本发明的具体实施只局限于这些说明。对于本发明所属技术领域的普通技术人员来说,在不脱离本发明构思的前提下,还可以做出若干简单推演或替换,都应当视为属于本发明的保护范围。The above contents are further detailed descriptions of the present invention in combination with specific preferred embodiments, and it cannot be determined that the specific implementation of the present invention is limited to these descriptions. For ordinary technicians in the technical field to which the present invention belongs, several simple deductions or substitutions can be made without departing from the concept of the present invention, which should be regarded as falling within the scope of protection of the present invention.

Claims (9)

Translated fromChinese
1.一种基于偏振三维成像的遥感地形三维测绘方法,其特征在于,包括:1. A method for three-dimensional surveying and mapping of remote sensing terrain based on polarization three-dimensional imaging, characterized by comprising:采集预设角度下不同目标地形区域的偏振图像;Collect polarization images of different target terrain areas at preset angles;针对每个目标地形区域,计算目标地形区域表面的偏振度,并根据所述偏振度计算所述目标地形区域表面的入射光的入射角θ和方位角For each target terrain area, the polarization degree of the surface of the target terrain area is calculated, and the incident angle θ and azimuth angle of the incident light on the surface of the target terrain area are calculated according to the polarization degree.根据所述入射角θ和所述方位角构建所述目标地形区域表面的法向量与所述入射角θ及所述方位角之间的第一映射关系,获得第一法向量信息;According to the incident angle θ and the azimuth angle Construct the normal vector of the target terrain area surface and the incident angle θ and the azimuth angle A first mapping relationship between the two, obtaining first normal vector information;计算太阳的天顶角θz与方位角θs,构建光源矢量与所述天顶角θz及所述方位角θs之间的第二映射关系,获得所述目标地形区域表面的光源矢量;Calculating the sun's zenith angleθz and azimuth angleθs , constructing a second mapping relationship between the light source vector and the zenith angleθz and the azimuth angleθs , and obtaining the light source vector on the surface of the target terrain area;构建所述目标地形区域表面的辐射模型,并结合所述目标地形区域表面的光源矢量计算先验法向量信息;Constructing a radiation model of the surface of the target terrain area, and calculating prior normal vector information in combination with the light source vector of the surface of the target terrain area;基于所述先验法向量信息对所述第一法向量信息进行修正,得到所述目标地形区域表面的修正法向量信息;Correcting the first normal vector information based on the prior normal vector information to obtain corrected normal vector information of the surface of the target terrain area;根据所述修正法向量信息和所述预设角度下不同目标地形区域的偏振图像,进行目标区域的三维场景重建。The three-dimensional scene of the target area is reconstructed according to the corrected normal vector information and the polarization images of different target terrain areas at the preset angles.2.根据权利要求1所述的基于偏振三维成像的遥感地形三维测绘方法,其特征在于,所述预设角度包括0°、45°、90°和135°;2. The remote sensing terrain three-dimensional mapping method based on polarization three-dimensional imaging according to claim 1, characterized in that the preset angles include 0°, 45°, 90° and 135°;针对每个目标地形区域,计算目标地形区域表面的偏振度,并根据所述偏振度计算所述目标地形区域表面的入射光的入射角θ和方位角的步骤,包括:For each target terrain area, the polarization degree of the surface of the target terrain area is calculated, and the incident angle θ and azimuth angle of the incident light on the surface of the target terrain area are calculated according to the polarization degree. The steps include:针对每个目标地形区域,获取采集预设角度下的偏振图像时通过0°偏光片的光强I、通过45°偏光片的光强I+45°、通过90°偏光片的光强I90°和通过135°方向线偏振器的光强I-45°,并计算斯托克斯参数S0、S1和S2For each target terrain area, the light intensity I passing through the 0° polarizer, the light intensity I+45° passing through the 45° polarizer, the light intensity I90° passing through the 90° polarizer, and the light intensity I-45° passing through the 135° linear polarizer are obtained when collecting polarization images at a preset angle, and the Stokes parameters S0 , S1 and S2 are calculated:S0=I+I90°S0 =I +I90° ;S1=I-I90°S1 =I -I90° ;S2=I+45°-I-45°S2 =I+45° -I-45° ;根据所述斯托克斯参数S0、S1和S2,计算目标地形区域表面的偏振度:According to the Stokes parameters S0 , S1 and S2 , the polarization degree of the surface of the target terrain area is calculated:根据所述斯托克斯参数S0、S1及所述偏振度,计算目标地形区域表面的入射光的入射角θ和方位角According to the Stokes parameters S0 , S1 and the polarization degree, the incident angle θ and the azimuth angle of the incident light on the surface of the target terrain area are calculated.3.根据权利要求2所述的基于偏振三维成像的遥感地形三维测绘方法,其特征在于,根据所述斯托克斯参数S0、S1及所述偏振度,计算目标地形区域表面的入射光的入射角θ和方位角的步骤,包括:3. The remote sensing terrain three-dimensional mapping method based on polarization three-dimensional imaging according to claim 2 is characterized in that the incident angle θ and the azimuth angle of the incident light on the surface of the target terrain area are calculated according to the Stokes parameters S0 , S1 and the polarization degree. The steps include:根据所述斯托克斯参数S0、S1计算目标地形区域表面的入射光的方位角Calculate the azimuth of the incident light on the surface of the target terrain area according to the Stokes parameters S0 and S1根据所述偏振度DOLP计算入射光的入射角θ:The incident angle θ of the incident light is calculated according to the degree of polarization DOLP:式中,n表示目标地形区域表面的折射率,n=1.5。Wherein, n represents the refractive index of the surface of the target terrain area, n=1.5.4.根据权利要求2所述的基于偏振三维成像的遥感地形三维测绘方法,其特征在于,根据所述入射角θ和所述方位角构建所述目标地形区域表面的法向量与所述入射角θ及所述方位角之间的第一映射关系,获得第一法向量信息的步骤,包括:4. The remote sensing terrain 3D mapping method based on polarization 3D imaging according to claim 2, characterized in that according to the incident angle θ and the azimuth angle Construct the normal vector of the target terrain area surface and the incident angle θ and the azimuth angle The step of obtaining first normal vector information according to the first mapping relationship between the two comprises:根据所述入射角θ和所述方位角构建所述目标地形区域表面的法向量与所述入射角θ及所述方位角之间的第一映射关系:According to the incident angle θ and the azimuth angle Construct the normal vector of the target terrain area surface and the incident angle θ and the azimuth angle The first mapping relationship between:式中,nx、ny、nz分别表示所述目标地形区域表面的法向量在x方向、y方向和z方向上的分量,其中,z方向与所述目标地形区域表面垂直,x方向与y方向垂直且二者所在平面与所述z方向垂直;Wherenx ,ny andnz represent the surface normal vectors of the target terrain area respectively. Components in the x-direction, y-direction and z-direction, wherein the z-direction is perpendicular to the surface of the target terrain area, the x-direction is perpendicular to the y-direction and the plane where the two directions are located is perpendicular to the z-direction;根据所述第一映射关系计算第一法向量信息D。The first normal vector information D is calculated according to the first mapping relationship.5.根据权利要求4所述的基于偏振三维成像的遥感地形三维测绘方法,其特征在于,计算太阳的天顶角θz与方位角θs,构建光源矢量与所述天顶角θz及所述方位角θs之间的第二映射关系,获得所述目标地形区域表面的光源矢量的步骤,包括:5. The method for three-dimensional surveying and mapping of remote sensing terrain based on polarization three-dimensional imaging according to claim 4, characterized in that the step of calculating the sun's zenith angleθz and azimuth angleθs , constructing a second mapping relationship between the light source vector and the zenith angleθz and the azimuth angleθs , and obtaining the light source vector on the surface of the target terrain area comprises:根据所述目标地形区域的被拍摄时间,计算太阳赤纬角:According to the time when the target terrain area was photographed, the solar declination angle is calculated:式中,N表示积日数;In the formula, N represents the cumulative number of days;根据所述太阳赤纬角δ、目标地形区域的纬度φ和太阳时角ω计算太阳的高度角:The solar altitude angle is calculated according to the solar declination angle δ, the latitude φ of the target terrain area and the solar hour angle ω:sinhs=sinδ·sinφ+cosδ·cosφ·cosω;sinhs = sinδ·sinφ+cosδ·cosφ·cosω;式中,ω=(12-t)×15°,t表示目标地形区域被拍摄时的当地时间;Where, ω = (12-t) × 15°, t represents the local time when the target terrain area is photographed;根据太阳赤纬角δ、目标地形区域的经纬度信息和太阳高度角hs计算太阳的方位角:The azimuth of the sun is calculated based on the solar declination angle δ, the longitude and latitude information of the target terrain area, and the solar altitude anglehs :式中,φ表示目标地形区域的纬度;Where φ represents the latitude of the target terrain area;根据太阳的高度角hs计算太阳的天顶角θzCalculate the solar zenith angleθz based on the solar altitude anglehs :θz=90°-hsθz = 90°-hs ;构建光源矢量与所述天顶角θz及所述方位角θs之间的第二映射关系:Construct a second mapping relationship between the light source vector and the zenith angleθz and the azimuth angleθs :根据所述第二映射关系,获得所述目标地形区域表面的光源矢量(ps,qs,-1)。According to the second mapping relationship, the light source vector (ps ,qs , -1) of the surface of the target terrain area is obtained.6.根据权利要求5所述的基于偏振三维成像的遥感地形三维测绘方法,其特征在于,构建所述目标地形区域表面的辐射模型,并结合所述目标地形区域表面的光源矢量计算先验法向量信息的步骤,包括:6. The method for remote sensing terrain three-dimensional mapping based on polarization three-dimensional imaging according to claim 5 is characterized in that the step of constructing a radiation model of the surface of the target terrain area and calculating the prior normal vector information in combination with the light source vector of the surface of the target terrain area comprises:根据朗伯反射模型、所述目标地形区域表面的光源矢量(ps,qs,-1)和光滑性约束,构建所述目标地形区域表面的辐射模型:According to the Lambertian reflection model, the light source vector (ps ,qs , -1) on the surface of the target terrain area and the smoothness constraint, a radiation model of the surface of the target terrain area is constructed:所述光滑性约束为:The smoothness constraint is:式中,I(x',y')表示目标地形区域(x',y')处在0°和90°下的偏振图像的光源强度之和,ρ表示预设的偏振图像反照率,px,py,qx,qy分别为p(x',y')、q(x',y')对应于x方向和y方向上的导数;Where I(x',y') represents the sum of the light source intensities of the polarization images of the target terrain area (x',y') at 0° and 90°, ρ represents the preset polarization image albedo,px ,py ,qx ,qy are the derivatives of p(x',y') and q(x',y') in the x and y directions respectively;根据所述辐射模型及光滑性约束,计算得到所述目标地形区域表面(x',y')处的先验法向量信息Y=(p(x',y'),q(x',y'),-1)。According to the radiation model and the smoothness constraint, the prior normal vector information Y=(p(x',y'),q(x',y'),-1) at the surface (x',y') of the target terrain area is calculated.7.根据权利要求6所述的基于偏振三维成像的遥感地形三维测绘方法,其特征在于,基于所述先验法向量信息对所述第一法向量信息进行修正,得到所述目标地形区域表面修正后的法向量的步骤,包括:7. The method for remote sensing terrain three-dimensional mapping based on polarization three-dimensional imaging according to claim 6 is characterized in that the step of correcting the first normal vector information based on the prior normal vector information to obtain the corrected normal vector of the surface of the target terrain area comprises:按照如下公式计算二元算子The binary operator is calculated according to the following formula式中,T表示非线性算子;In the formula, T represents the nonlinear operator;基于所述二元算子对第一法向量信息D进行修正:Based on the binary operator Correct the first normal vector information D:式中,Ncor表示所述目标地形区域表面的修正法向量信息。Wherein, Ncor represents the corrected normal vector information of the surface of the target terrain area.8.根据权利要求1所述的基于偏振三维成像的遥感地形三维测绘方法,其特征在于,所述不同目标区域包括第一目标区域和第二目标区域,所述第一目标区域与所述第二目标区域至少存在部分重叠区域;8. The remote sensing terrain 3D mapping method based on polarization 3D imaging according to claim 1, characterized in that the different target areas include a first target area and a second target area, and the first target area and the second target area have at least a partial overlapping area;根据所述修正法向量信息和所述预设角度下不同目标地形区域的偏振图像,进行目标区域的三维场景重建的步骤,包括:The step of reconstructing a three-dimensional scene of the target area according to the corrected normal vector information and the polarization images of different target terrain areas at the preset angles includes:从所述第一目标区域的任意一幅偏振图像中选取一特征点,所述特征点位于重叠区域内;Selecting a feature point from any one of the polarization images of the first target area, wherein the feature point is located in the overlapping area;获取所述特征点在所述第二目标区域中的位置信息,并建立所述第一目标区域的偏振图像与所述第二目标区域的偏振图像之间的变换关系;Acquire position information of the feature point in the second target area, and establish a transformation relationship between the polarization image of the first target area and the polarization image of the second target area;根据所述变换关系拼接所述第一目标区域的偏振图像与所述第二目标区域的偏振图像;splicing the polarization image of the first target area and the polarization image of the second target area according to the transformation relationship;根据拼接后的图像、第一目标区域表面的修正法向量信息和第二目标区域表面的修正法向量信息进行三维重建。Three-dimensional reconstruction is performed based on the stitched image, the corrected normal vector information of the surface of the first target area, and the corrected normal vector information of the surface of the second target area.9.根据权利要求8所述的基于偏振三维成像的遥感地形三维测绘方法,其特征在于,根据拼接后的图像、第一目标区域表面的修正法向量信息和第二目标区域表面的修正法向量信息,按照如下公式进行三维重建:9. The method for three-dimensional surveying and mapping of remote sensing terrain based on polarization three-dimensional imaging according to claim 8 is characterized in that three-dimensional reconstruction is performed according to the following formula based on the spliced image, the corrected normal vector information of the surface of the first target area, and the corrected normal vector information of the surface of the second target area:式中,F{}和F-1{}分别表示离散傅里叶变换和傅里叶逆变换,M、N分别表示拼接后的图像中像素点的行数及列数,u、ν分别表示离散微算子在x和y方向上的傅里叶系数。Where F{} and F-1 {} represent discrete Fourier transform and inverse Fourier transform, respectively; M and N represent the number of rows and columns of pixels in the spliced image, respectively; u and ν represent the Fourier coefficients of the discrete differential operator in the x and y directions, respectively.
CN202310484598.5A2023-04-282023-04-28Remote sensing terrain three-dimensional mapping method based on polarization three-dimensional imagingPendingCN116576827A (en)

Priority Applications (1)

Application NumberPriority DateFiling DateTitle
CN202310484598.5ACN116576827A (en)2023-04-282023-04-28Remote sensing terrain three-dimensional mapping method based on polarization three-dimensional imaging

Applications Claiming Priority (1)

Application NumberPriority DateFiling DateTitle
CN202310484598.5ACN116576827A (en)2023-04-282023-04-28Remote sensing terrain three-dimensional mapping method based on polarization three-dimensional imaging

Publications (1)

Publication NumberPublication Date
CN116576827Atrue CN116576827A (en)2023-08-11

Family

ID=87544591

Family Applications (1)

Application NumberTitlePriority DateFiling Date
CN202310484598.5APendingCN116576827A (en)2023-04-282023-04-28Remote sensing terrain three-dimensional mapping method based on polarization three-dimensional imaging

Country Status (1)

CountryLink
CN (1)CN116576827A (en)

Cited By (2)

* Cited by examiner, † Cited by third party
Publication numberPriority datePublication dateAssigneeTitle
CN119085616A (en)*2024-11-072024-12-06山东省国土测绘院 Three-dimensional terrain mapping system based on remote sensing technology
CN119672208A (en)*2024-11-062025-03-21西安电子科技大学杭州研究院 A high-precision stereo reconstruction method for complex scenes

Citations (2)

* Cited by examiner, † Cited by third party
Publication numberPriority datePublication dateAssigneeTitle
CN109147029A (en)*2018-06-252019-01-04西安电子科技大学A kind of monocular polarization three-dimensional rebuilding method
CN115235369A (en)*2022-03-142022-10-25西安电子科技大学 A polarization three-dimensional imaging method, system and application of a color target

Patent Citations (2)

* Cited by examiner, † Cited by third party
Publication numberPriority datePublication dateAssigneeTitle
CN109147029A (en)*2018-06-252019-01-04西安电子科技大学A kind of monocular polarization three-dimensional rebuilding method
CN115235369A (en)*2022-03-142022-10-25西安电子科技大学 A polarization three-dimensional imaging method, system and application of a color target

Non-Patent Citations (1)

* Cited by examiner, † Cited by third party
Title
李轩: "复杂场景偏振三维成像关键技术研究", 《中国博士学位论文全文数据库(电子期刊)基础科学辑》, no. 2, 15 February 2023 (2023-02-15), pages 005 - 63*

Cited By (3)

* Cited by examiner, † Cited by third party
Publication numberPriority datePublication dateAssigneeTitle
CN119672208A (en)*2024-11-062025-03-21西安电子科技大学杭州研究院 A high-precision stereo reconstruction method for complex scenes
CN119085616A (en)*2024-11-072024-12-06山东省国土测绘院 Three-dimensional terrain mapping system based on remote sensing technology
CN119085616B (en)*2024-11-072025-03-04山东省国土测绘院Three-dimensional topography mapping system based on remote sensing technology

Similar Documents

PublicationPublication DateTitle
Teller et al.Calibrated, registered images of an extended urban area
CN107767440B (en)Cultural relic sequence image fine three-dimensional reconstruction method based on triangulation network interpolation and constraint
CN116051659B (en) A joint calibration method of line scan camera and 2D laser scanner
Haala et al.Quality of 3D point clouds from highly overlapping UAV imagery
AU2011312140C1 (en)Rapid 3D modeling
CN106643669B (en)A kind of more camera lens multi-detector aerial camera single centre projection transform methods
CN116576827A (en)Remote sensing terrain three-dimensional mapping method based on polarization three-dimensional imaging
CN108168521A (en)One kind realizes landscape three-dimensional visualization method based on unmanned plane
Verykokou et al.Oblique aerial images: a review focusing on georeferencing procedures
CN105300362B (en)A kind of photogrammetric survey method applied to RTK receiver
CN104913780B (en)The high-precision deviation of plumb line method for fast measuring of integrated GNSS and CCD zenith telescopes
Parente et al.Optimising the quality of an SfM‐MVS slope monitoring system using fixed cameras
CN102735216B (en)CCD stereoscopic camera three-line imagery data adjustment processing method
CN110986888A (en)Aerial photography integrated method
CN108269234A (en)A kind of lens of panoramic camera Attitude estimation method and panorama camera
CN108955642B (en) A seamless stitching method for large-format equivalent central projection images
CN119963737A (en) A passive three-dimensional reconstruction method and device based on single polarization information
Zhang et al.High-precision photogrammetric 3d modeling technology based on multi-source data fusion and deep learning-enhanced feature learning using internet of things big data
Elhalawani et al.Implementation of close range photogrammetry using modern non-metric digital cameras for architectural documentation
CN109029379A (en)A kind of high-precision stereo mapping with low base-height ratio method
WuPhotogrammetry: 3-D from imagery
CN107941201B (en)The zero intersection optical satellite image simultaneous adjustment method and system that light is constrained with appearance
Swaminathan et al.Polycameras: Camera clusters for wide angle imaging
Zhou et al.Calibration of fisheye camera with colinear constraint of the main optical axis and epipolar line orthogonality-transverse axis
Zhang et al.Tests and performance evaluation of DMC images and new methods for their processing

Legal Events

DateCodeTitleDescription
PB01Publication
PB01Publication
SE01Entry into force of request for substantive examination
SE01Entry into force of request for substantive examination

[8]ページ先頭

©2009-2025 Movatter.jp