CN115061138B - 探地雷达自聚焦并行掩膜后向投影成像方法、设备及介质 - Google Patents

探地雷达自聚焦并行掩膜后向投影成像方法、设备及介质 Download PDF

Info

Publication number
CN115061138B
CN115061138B CN202210992032.9A CN202210992032A CN115061138B CN 115061138 B CN115061138 B CN 115061138B CN 202210992032 A CN202210992032 A CN 202210992032A CN 115061138 B CN115061138 B CN 115061138B
Authority
CN
China
Prior art keywords
value
background medium
dielectric constant
mask
relative dielectric
Prior art date
Legal status (The legal status is an assumption and is not a legal conclusion. Google has not performed a legal analysis and makes no representation as to the accuracy of the status listed.)
Active
Application number
CN202210992032.9A
Other languages
English (en)
Other versions
CN115061138A (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.)
Central South University
Original Assignee
Central South 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 Central South University filed Critical Central South University
Priority to CN202210992032.9A priority Critical patent/CN115061138B/zh
Publication of CN115061138A publication Critical patent/CN115061138A/zh
Application granted granted Critical
Publication of CN115061138B publication Critical patent/CN115061138B/zh
Active legal-status Critical Current
Anticipated expiration legal-status Critical

Links

Images

Classifications

    • G—PHYSICS
    • G01—MEASURING; TESTING
    • G01S—RADIO DIRECTION-FINDING; RADIO NAVIGATION; DETERMINING DISTANCE OR VELOCITY BY USE OF RADIO WAVES; LOCATING OR PRESENCE-DETECTING BY USE OF THE REFLECTION OR RERADIATION OF RADIO WAVES; ANALOGOUS ARRANGEMENTS USING OTHER WAVES
    • G01S13/00—Systems using the reflection or reradiation of radio waves, e.g. radar systems; Analogous systems using reflection or reradiation of waves whose nature or wavelength is irrelevant or unspecified
    • G01S13/88—Radar or analogous systems specially adapted for specific applications
    • G01S13/89—Radar or analogous systems specially adapted for specific applications for mapping or imaging
    • G01S13/90—Radar or analogous systems specially adapted for specific applications for mapping or imaging using synthetic aperture techniques, e.g. synthetic aperture radar [SAR] techniques
    • G01S13/9004—SAR image acquisition techniques
    • G01S13/9019—Auto-focussing of the SAR signals
    • G—PHYSICS
    • G01—MEASURING; TESTING
    • G01S—RADIO DIRECTION-FINDING; RADIO NAVIGATION; DETERMINING DISTANCE OR VELOCITY BY USE OF RADIO WAVES; LOCATING OR PRESENCE-DETECTING BY USE OF THE REFLECTION OR RERADIATION OF RADIO WAVES; ANALOGOUS ARRANGEMENTS USING OTHER WAVES
    • G01S7/00—Details of systems according to groups G01S13/00, G01S15/00, G01S17/00
    • G01S7/02—Details of systems according to groups G01S13/00, G01S15/00, G01S17/00 of systems according to group G01S13/00
    • G01S7/41—Details of systems according to groups G01S13/00, G01S15/00, G01S17/00 of systems according to group G01S13/00 using analysis of echo signal for target characterisation; Target signature; Target cross-section
    • G—PHYSICS
    • G06—COMPUTING; CALCULATING OR COUNTING
    • G06T—IMAGE DATA PROCESSING OR GENERATION, IN GENERAL
    • G06T11/00—2D [Two Dimensional] image generation
    • G06T11/003—Reconstruction from projections, e.g. tomography
    • G06T11/005—Specific pre-processing for tomographic reconstruction, e.g. calibration, source positioning, rebinning, scatter correction, retrospective gating
    • G—PHYSICS
    • G06—COMPUTING; CALCULATING OR COUNTING
    • G06T—IMAGE DATA PROCESSING OR GENERATION, IN GENERAL
    • G06T11/00—2D [Two Dimensional] image generation
    • G06T11/003—Reconstruction from projections, e.g. tomography
    • G06T11/006—Inverse problem, transformation from projection-space into object-space, e.g. transform methods, back-projection, algebraic methods
    • G—PHYSICS
    • G06—COMPUTING; CALCULATING OR COUNTING
    • G06T—IMAGE DATA PROCESSING OR GENERATION, IN GENERAL
    • G06T11/00—2D [Two Dimensional] image generation
    • G06T11/003—Reconstruction from projections, e.g. tomography
    • G06T11/008—Specific post-processing after tomographic reconstruction, e.g. voxelisation, metal artifact correction
    • G—PHYSICS
    • G06—COMPUTING; CALCULATING OR COUNTING
    • G06T—IMAGE DATA PROCESSING OR GENERATION, IN GENERAL
    • G06T5/00—Image enhancement or restoration
    • G06T5/20—Image enhancement or restoration using local operators
    • G06T5/30—Erosion or dilatation, e.g. thinning
    • G—PHYSICS
    • G06—COMPUTING; CALCULATING OR COUNTING
    • G06T—IMAGE DATA PROCESSING OR GENERATION, IN GENERAL
    • G06T5/00—Image enhancement or restoration
    • G06T5/73—Deblurring; Sharpening
    • G—PHYSICS
    • G06—COMPUTING; CALCULATING OR COUNTING
    • G06T—IMAGE DATA PROCESSING OR GENERATION, IN GENERAL
    • G06T7/00—Image analysis
    • G06T7/10—Segmentation; Edge detection
    • G06T7/11—Region-based segmentation
    • Y—GENERAL TAGGING OF NEW TECHNOLOGICAL DEVELOPMENTS; GENERAL TAGGING OF CROSS-SECTIONAL TECHNOLOGIES SPANNING OVER SEVERAL SECTIONS OF THE IPC; TECHNICAL SUBJECTS COVERED BY FORMER USPC CROSS-REFERENCE ART COLLECTIONS [XRACs] AND DIGESTS
    • Y02—TECHNOLOGIES OR APPLICATIONS FOR MITIGATION OR ADAPTATION AGAINST CLIMATE CHANGE
    • Y02A—TECHNOLOGIES FOR ADAPTATION TO CLIMATE CHANGE
    • Y02A90/00—Technologies having an indirect contribution to adaptation to climate change
    • Y02A90/30—Assessment of water resources

Landscapes

  • Engineering & Computer Science (AREA)
  • Physics & Mathematics (AREA)
  • General Physics & Mathematics (AREA)
  • Theoretical Computer Science (AREA)
  • Remote Sensing (AREA)
  • Radar, Positioning & Navigation (AREA)
  • Computer Networks & Wireless Communication (AREA)
  • Mathematical Physics (AREA)
  • Pure & Applied Mathematics (AREA)
  • Computer Vision & Pattern Recognition (AREA)
  • Mathematical Optimization (AREA)
  • Mathematical Analysis (AREA)
  • Algebra (AREA)
  • Electromagnetism (AREA)
  • Radar Systems Or Details Thereof (AREA)

Abstract

本发明提供一种探地雷达自聚焦并行掩膜后向投影成像方法、设备及介质，其中方法包括：获得B‑scan回波数据；对B‑scan回波数据归一化成像，进行二值化处理与形态学滤波，得到预掩膜矩阵；设定背景介质相对介电常数估计区间，令估计区间中值为背景介质相对介电常数估计值；根据该估计值和预掩膜矩阵，并行计算获得并行掩膜后向投影成像结果；根据二分法和并行掩膜后向投影成像结果更新估计区间，迭代运算，直至估计区间满足设定条件，获得背景介质最终估计值与自聚焦并行掩膜后向投影成像结果。本发明能够提高成像质量，缩短计算时间，在未知背景介质相对介电常数条件下实现背景介质相对介电常数估计与自聚焦并行掩膜后向投影成像。

Description

探地雷达自聚焦并行掩膜后向投影成像方法、设备及介质
技术领域
本发明属于探地雷达成像的技术领域，具体涉及一种探地雷达自聚焦并行掩膜后向投影成像方法、一种计算机设备和一种计算机可读存储介质。
背景技术
探地雷达是一种利用天线发射和接收高频电磁波来探测介质内部物质特性和分布的设备，该设备根据电磁波在传播过程中遇到与背景介质电磁性质不同的物体或结构时会发生散射的特性，通过对回波信号分析以对地下目标的空间位置、几何结构以及物质特性等关键参数的估计，能够实现对地下目标体的无损探测。由于探地雷达能够进行非接触式的无损探测，因此被广泛应用于考古、矿产勘查与灾害地质调查等众多领域。
探地雷达有多种数据记录方式，例如，可以采用固定雷达的发射与接收天线完成一次探测以获得的单道回波数据，该单道回波数据称为A-scan回波数据。A-scan回波数据是探地雷达最基本的数据格式，其中记录了电磁波在探测区域传播过程中的瞬时振幅、瞬时相位、瞬时频率与双程传播时间等信息，但是仅通过单道的A-scan回波数据，难以实现对探测区域的准确有效分析。将发射天线与接收天线在地面上沿测线等间距移动探测，按空间测点顺序组合A-scan回波数据的数据格式为B-scan回波数据。在B-scan回波数据所成的B-scan图像中，由于不同测点的天线位置同地下异常目标的空间距离不同，令地下异常目标在不同A-scan回波数据中对应幅值具有不同时延，体现为B-scan图像中开口向下的双曲线特征。
聚焦成像是探地雷达信号处理与分析的重要方法，在众多的聚焦成像方法中，由于基于“时延-求和”的后向投影方法适用于对复杂有耗介质中点散射目标的成像处理，因此其被广泛应用于探地雷达信号处理中。现有技术中的后向投影聚焦成像方法的流程如图1所示，包括：
在步骤S01中，获取B-scan回波数据，对其进行去直耦波、零点校正等预处理，并对B-scan回波数据成像获得B-scan图像；
在步骤S02中，获取探测时的工作参数，该工作参数包括：天线步进距离、时窗长度和时窗采样间隔，并根据先验知识，获得背景介质相对介电常数；
在步骤S03中，设定成像区域，将天线步进间隔作为成像区域单元格划分单位，获得一个用于存储散射强度的累加值的二维矩阵；
在步骤 S04 中，计算电磁波在成像区域中各成像单元和每一测点之间传播的双程时
延，并根据该时延值提取对应测点处 A-scan 回波数据中的幅值，生成成像单元的一维信
号序列，该序列的和即为成像单元的成像结果。
现有技术中的后向投影聚焦成像方法进行成像处理前需要提前获得背景介质的相对介电常数估计值等参数后才能进行迭代求解。在实际工作中，实验记录缺失、无法直接获得背景介质性质等问题严重影响着该算法的应用。在B-scan回波数据中，双曲线特征区域面积较小，对成像区域所有单元格的遍历求解产生了大量冗余计算，影响了算法的时间经济性，且会造成较强的伪影，影响成像质量。
综上所述可知，现有技术中在对探地雷达的回波数据进行后向投影聚焦成像时，其存在成像质量差、计算时间长及无法在缺失背景介质相关信息的条件下使用的问题。
发明内容
本发明的目的在于提供一种探地雷达自聚焦并行掩膜后向投影成像方法、设备及介质，为至少解决上述现有技术中在对探地雷达的回波数据进行后向投影聚焦成像时，其存在成像质量差、计算时间长及无法在缺失背景介质相对介电常数的条件下使用的问题。
为至少解决上述技术问题，第一方面，本发明提供了一种探地雷达自聚焦并行掩膜后向投影成像方法，包括：
获取B-scan回波数据，以及探测时的工作参数，对B-scan回波数据进行去直耦波和校正零点预处理，并进行归一化处理并成像，获得B-scan图像；
根据优化函数获取最优阈值，根据该最优阈值对所述B-scan图像进行二值化处理，并对二值化处理后的B-scan图像进行形态学滤波，以得到预掩膜矩阵；
根据背景介质性质，设定背景介质相对介电常数估计区间上限值与下限值，并将该估计区间中值作为背景介质相对介电常数估计值；
根据所述背景介质相对介电常数估计值、预掩膜矩阵和探测时的工作参数，通过并行计算获得并行掩膜后向投影成像结果；
根据二分法和所述并行掩膜后向投影成像结果，更新背景介质相对介电常数估计区间，计算背景介质相对介电常数估计值，再次通过并行计算获得并行掩膜后向投影成像结果并更新背景介质相对介电常数估计区间，迭代运算，直至背景介质相对介电常数估计区间满足设定条件，获得背景介质相对介电常数最终估计值和自聚焦并行掩膜后向投影成像结果。
根据本发明的一个实施例，所述对B-scan回波数据进行去直耦波和校正零点预处理，包括：
对每行B-scan回波数据减去该行均值去除B-scan回波数据中的直耦波；
将目标按一定深度埋在已知相对介电常数的均匀介质中，在其正上方进行探测，提取回波数据中对应该目标幅值的响应时延值，计算回波数据中对应该目标幅值的理想时延值，校正使响应时延值与理想时延值相等。
根据本发明的一个实施例，所述根据优化函数获取最优阈值，包括：
获取所有待选阈值；
获取所述各待选阈值所对应的优化函数输出值；
根据各优化函数的输出值，从所述待选阈值中确定所述最优阈值。
根据本发明的一个实施例，所述获取所述各待选阈值所对应的优化函数输出值，包括：
根据待选阈值对所述B-scan回波数据进行分类，以得到第一类数据和第二类数据；
获取所述第一类数据的特征值和第二类数据的特征值，并根据所述第一类数据的特征值和第二类数据的特征值，得到该阈值对应的优化函数输出值；
所述根据各优化函数的输出值，从所述待选阈值中确定所述最优阈值包括：获取优化函数最大输出值所对应的阈值，并将该阈值作为所述最优阈值。
根据本发明的一个实施例，所述对二值化处理后的B-scan图像进行形态学滤波包括：
采用矩形滤波器对所述二值化处理后的B-scan图像进行第一设定次数的膨胀-腐蚀滤波，以得到初次滤波后的B-scan图像；
分别采用上斜滤波器和下斜滤波器对初次滤波后的B-scan图像进行第二设定次数的腐蚀处理和第三设定次数的膨胀处理。
根据本发明的一个实施例，所述根据所述背景介质相对介电常数估计值、预掩膜矩阵和探测时的工作参数，通过并行计算获得并行掩膜后向投影成像结果，包括：
根据时延-深度对应关系和所述背景介质相对介电常数估计值、预掩膜矩阵和探测时的工作参数，获得成像掩膜矩阵，根据成像掩膜矩阵获得掩膜扩展矩阵，设置与掩膜扩展矩阵同维度的成像区域，对成像区域中对应掩膜扩展矩阵前景元素相同位置的成像单元进行后向投影成像，以获得并行掩膜后向投影成像结果。
根据本发明的一个实施例，所述根据二分法和所述并行掩膜后向投影成像结果，更新背景介质相对介电常数估计区间，计算背景介质相对介电常数估计值，包括；
根据所述并行掩膜后向投影成像结果，计算出地下目标上方区域聚焦度和地下目标下方区域聚焦度；
根据地下目标上方区域聚焦度和地下目标下方区域聚焦度对所述背景介质相对介电常数估计区间上限值或下限值进行更新，将更新后的估计区间中值作为背景介质相对介电常数新的估计值，包括增加所述背景介质相对介电常数的估计值和减小所述背景介质相对介电常数的估计值；
所述更新背景介质相对介电常数估计区间包括：根据二分法，若所述地下目标上方区域的聚焦度大于地下目标下方区域的聚焦度，则更新所述背景介质相对介电常数估计区间下限值为当前背景介质相对介电常数估计值，否则更新所述背景介质相对介电常数估计范围上限值为当前背景介质相对介电常数估计值；
所述计算背景介质相对介电常数的估计值包括：计算背景介质相对介电常数估计区间中值，该值即为背景介质相对介电常数估计值。
根据本发明的一个实施例，所述根据二分法和所述并行掩膜后向投影成像结果，迭代运算直至背景介质相对介电常数估计区间满足设定条件，获得背景介质相对介电常数最终估计值和自聚焦并行掩膜后向投影成像结果，包括：
所述迭代运算包括：根据所述并行掩膜后向投影成像结果，更新背景介质相对介电常数估计区间，将更新后的估计区间中值作为背景介质相对介电常数估计值，根据所述背景介质相对介电常数估计值与预掩膜矩阵，再次并行计算获得并行掩膜后向投影结果，并再次更新背景介质相对介电常数估计区间，直至背景介质相对介电常数估计区间满足设定条件；
所述设定条件包括：更新后的背景介质相对介电常数估计区间的上限值与下限值之差小于设定差值。
所述背景介质相对介电常数最终估计值为背景介质相对介电常数估计区间满足设定条件时的区间中值；
所述获得自聚焦并行掩膜后向投影结果包括：根据时延-深度关系与所述背景介质相对介电常数最终估计值与预掩膜矩阵，获得成像掩膜矩阵，设置与成像掩膜矩阵同维度的成像区域，通过并行计算，对成像区域中对应成像掩膜矩阵前景元素相同位置的成像单元进行后向投影成像，获得自聚焦并行掩膜后向投影成像结果。
第二方面，本发明还提供了一种计算机设备，包括处理器和存储器，所述存储器上存储有用于在所述处理器上执行的计算机程序指令，所述处理器执行该计算机程序指令时，实现如上述任意一项实施例所述的探地雷达自聚焦并行掩膜后向投影成像方法。
第三方面，本发明还提供了一种计算机可读存储介质，其存储有计算机指令，所述计算机指令由处理器运行时实现如上述任意一项实施例所述的探地雷达自聚焦并行掩膜后向投影成像方法。
本发明的技术方案，采用优化函数获取最优阈值与对B-scan图像进行二值化处理，并对该二值化处理后的B-scan图像进行形态学滤波以得到预掩膜矩阵；设定背景介质相对介电常数估计区间上限值与下限值，并将估计区间中值作为背景介质相对介电常数估计值；通过并行计算，根据时延-深度对应关系、背景介质相对介电常数估计值与预掩膜矩阵，获得并行掩膜后向投影成像结果；根据二分法与并行掩膜后向投影成像结果更新背景介质相对介电常数估计区间，迭代运算直至估计区间满足设定条件，获得背景介质相对介电常数的最终估计值，并根据该最终估计值和成像掩膜矩阵计算获得自聚焦掩膜的并行成像结果。与现有技术相比，本发明的技术方案能够实现对背景介质相对介电常数的估计，在未知背景介质相对介电常数的条件下实现自聚焦并行掩膜后向投影成像，通过并行计算缩短成像所用的时间，同时提高了成像质量。
附图说明
通过参考附图阅读下文的详细描述，本发明示例性实施方式的上述以及其他目的、特征和优点将变得易于理解。在附图中，以示例性而非限制性的方式示出了本发明的若干实施方式，并且相同或对应的标号表示相同或对应的部分，其中：
图1为根据现有技术中的一种获取后向投影聚焦方法的流程图；
图2为根据本发明实施例的一种探地雷达自聚焦并行掩膜后向投影成像方法的流程图；
图3为根据本发明实施例的一种B-scan回波数据的成像示意图；
图4为根据本发明实施例的一种预处理后的B-scan回波数据的成像示意图；
图5为根据现有成像方法的后向投影成像的三维示意图；
图6为根据本发明实施例的掩膜矩阵成像的示意图；
图7为根据本发明实施例的背景介质相对介电常数估计过程中的并行掩膜后向投影成像示意图，其中a、b、c、d、e、f、g、h、i、j表示背景相对介电常数取值上限和下限值不同时并行掩膜后向投影成像示意图；
图8为根据本发明实施例的自聚焦并行掩膜后向投影成像的三维示意图；
图9为根据本发明实施例的一种计算机设备的结构示意图。
具体实施方式
下面结合本申请实施例中的附图，对本申请实施例中的技术方案进行清楚、完整地描述，本领域技术人员应知，下文所描述的实施例是本申请一部分实施例，而不是全部的实施例。基于本申请中的实施例，本领域技术人员在没有做出创造性劳动前提下所获得的所有其他实施例，都属于本申请保护的范围。
请参阅图2，图2所示出的是本申请所提供的一种探地雷达自聚焦并行掩膜后向投影成像方法的流程，该方法用于对探地雷达的B-scan回波数据进行处理，以得到对应的高质量后向投影成像。下面结合图2所示出的流程，对本申请的探地雷达自聚焦并行掩膜后向投影成像的方法做详细的介绍。
如图2所示，本申请的探地雷达自聚焦并行掩膜后向投影成像方法，其包括：
在步骤S11中，获取探地雷达探测得到的B-scan回波数据，以及探测雷达在探测时的工作参数，并对B-scan回波数据成像，获得B-scan图像。在本步骤S11中，记测线方向为
Figure 321278DEST_PATH_IMAGE001
方向，测线范围为
Figure 679579DEST_PATH_IMAGE002
，共有
Figure 824252DEST_PATH_IMAGE003
个测点，其中第
Figure 344226DEST_PATH_IMAGE004
个测点的坐标为
Figure 207140DEST_PATH_IMAGE005
，记在第
Figure 369131DEST_PATH_IMAGE004
个测点处获得的A-scan回波数据为
Figure 368311DEST_PATH_IMAGE006
，则可将得到的B-scan回波数据记为
Figure 855924DEST_PATH_IMAGE007
，即一个维度为N×W的二维矩阵。
探测雷达探测时的工作参数包括天线步进间隔
Figure 206134DEST_PATH_IMAGE008
、时窗长度
Figure 640658DEST_PATH_IMAGE009
、时间维采样点数量
Figure 759923DEST_PATH_IMAGE010
和天线零功率波瓣宽度
Figure 152859DEST_PATH_IMAGE011
，并通过如下公式计算出时窗采样间隔
Figure 255944DEST_PATH_IMAGE012
：
Figure 25317DEST_PATH_IMAGE013
（1）
在步骤S12中，采用优化函数获取最优阈值，并根据该最优阈值对上述B-scan图像进行二值化处理，并对该二值化处理后的B-scan图像进行形态学滤波，以得到预掩膜矩阵。在本步骤S12中，可以根据优化函数的输出值确定出最优阈值，上述形态学滤波包括膨胀-腐蚀滤波、腐蚀处理和膨胀处理。
在步骤S13中，根据背景介质性质设定背景介质相对介电常数估计区间上限值与下限值，将所述区间中值作为背景介质相对介电常数的估计值；根据时延-深度对应关系与所述背景介质相对介电常数估计值、预掩膜矩阵与探测时的工作参数，获得成像掩膜矩阵与掩膜扩展矩阵，设置与掩膜扩展矩阵同维度的成像区域，通过并行计算，对成像区域中对应掩膜扩展矩阵为前景元素相同位置的成像单元进行后向投影成像，获得并行掩膜后向投影成像结果。
在步骤S14中，根据二分法和所述并行掩膜后向投影成像结果，更新背景介质相对介电常数估计区间，计算背景介质相对介电常数估计值，再次通过并行计算获得并行掩膜后向投影成像结果并更新背景介质相对介电常数估计区间，迭代运算，直至背景介质相对介电常数估计区间满足设定条件，获得背景介质相对介电常数最终估计值。根据时延-深度对应关系与所述背景介质相对介电常数最终估计值与预掩膜矩阵，获得成像掩膜矩阵并设置与成像掩膜矩阵相同维度成像区域，对成像区域中对应成像掩膜矩阵前景元素相同位置的成像单元进行后向投影成像，获得自聚焦并行掩膜后向投影成像结果。
综上所述可知，本申请的技术方案中，采用优化函数获取最优阈值与对B-scan图像进行二值化处理，并对该二值化处理后的B-scan图像进行形态学滤波以得到预掩膜矩阵；根据背景介质性质设定背景介质相对介电常数估计区间上限值与下限值，将估计区间中值作为背景介质相对介电常数估计值，通过并行计算获得并行掩膜后向投影成像结果；根据二分法与并行掩膜后向投影成像结果更新背景介质相对介电常数估计区间，迭代运算直至背景介质相对介电常数估计区间满足设定条件，获得背景介质相对介电常数最终估计值，并根据该值与预掩膜矩阵获得自聚焦并行掩膜后向投影成像结果。因此与现有技术相比，本申请的技术方案能够在未知背景介质相对介电常数的条件下实现自聚焦并行掩膜后向投影成像，缩短了计算时间，提高了成像质量。
上文中对本申请的技术方案做了详细的介绍，下面结合具体应用场景，对本申请的探地雷达自聚焦并行掩膜后向投影成像方法的步骤S11进行详细的说明。
在一个实施例中，上述步骤S11中获取到探地雷达探测得到的B-scan回波数据后，首先对其进行预处理以消除其中的噪音数据。
本实施例中，对B-scan回波数据的预处理包括去直耦波处理和零点校正处理。对每行B-scan回波数据减去该行均值去除B-scan回波数据中的直耦波；将目标按一定深度埋在已知相对介电常数的均匀介质中，在其正上方进行探测，提取回波数据中对应该目标幅值的响应时延值，计算回波数据中对应该目标幅值的理想时延值，校正使响应时延值与理想时延值相等。
在一个实施例中，包括对得到的B-scan回波数据进行归一化处理并成像获得B-scan图像的步骤。在本实施例中，假设上述经过预处理后的B-scan回波数据为
Figure 993230DEST_PATH_IMAGE014
，则对该B-scan回波数据
Figure 557066DEST_PATH_IMAGE015
进行归一化处理时，所采用的计算公式为：
Figure 881868DEST_PATH_IMAGE016
（2）
其中
Figure 923774DEST_PATH_IMAGE017
为B-scan回波数据
Figure 283211DEST_PATH_IMAGE018
的最小值，
Figure 17949DEST_PATH_IMAGE019
为B-scan回波数据
Figure 564467DEST_PATH_IMAGE020
中的最大值，S为颜色深度取值255；
Figure 675643DEST_PATH_IMAGE021
为向上取整函数，
Figure 624007DEST_PATH_IMAGE022
为归一化处理并成像的B-scan图像。
上文中对本申请技术方案中的B-scan回波数据的预处理和归一化处理做了详细的介绍，下面结合具体应用场景，对上述步骤S12中获取最优阈值的方法做详细的说明。
在一个实施例中，上述步骤S12中获取最优阈值的方法包括：首先获取待选阈值的集合，在本实施例中，待选阈值的集合为0到S之间的整数值，即待选阈值的取值为0到255之间的整数；然后获取上述集合中各待选阈值所对应的优化函数输出值，即根据各待选阈值对优化函数进行求解，得到各待选阈值所对应的优化函数输出值；最后根据各优化函数输出值确定出最优阈值。下面以一种具体应用场景为例，对该获取最优阈值的方法做详细的介绍，可以理解的是，下文中对该方法的介绍是示例性的而非限制性的，因此上文中对该方法的描述，同样也适用于下文中对该方法的介绍。
在一个实施例中，设其中一个待选阈值为
Figure 529647DEST_PATH_IMAGE023
，则根据该阈值
Figure 563462DEST_PATH_IMAGE023
将上述B-scan图像
Figure 947170DEST_PATH_IMAGE024
中的值分为第一类数据和第二类数据，设其中第一类数据为P，该类数据中的值均小于等于上述待选阈值
Figure 750040DEST_PATH_IMAGE023
；上述第二类数据为Q，该类数据中值均大于上述待选阈值
Figure 92160DEST_PATH_IMAGE023
，则首先计算上述第一类数据的特征值
Figure 613271DEST_PATH_IMAGE025
和第二类数据的特征值
Figure 800670DEST_PATH_IMAGE026
，计算公式分别为：
Figure 989206DEST_PATH_IMAGE027
（3）
Figure 254226DEST_PATH_IMAGE028
（4）
其中
Figure 262633DEST_PATH_IMAGE029
为
Figure 519302DEST_PATH_IMAGE030
中值为j的元素个数。
然后根据上述第一类数据的特征值
Figure 296765DEST_PATH_IMAGE031
和第二类数据的特征值
Figure 980687DEST_PATH_IMAGE032
，计算出待选阈值
Figure 210812DEST_PATH_IMAGE023
所对应的优化函数输出值，计算时所采用的公式为：
Figure 5592DEST_PATH_IMAGE033
（5）
通过上述步骤可以计算出集合内各待选阈值所对应的优化函数输出值，然后获取其中最大值所对应的待选阈值，该待选阈值即为上述最优阈值。设所得到的最优阈值为
Figure 903141DEST_PATH_IMAGE034
，则根据该最优阈值
Figure 226806DEST_PATH_IMAGE034
对B-scan图像
Figure 475385DEST_PATH_IMAGE035
进行二值化处理时，所采用的计算公式为：
Figure 808277DEST_PATH_IMAGE036
（6）
其中
Figure 294753DEST_PATH_IMAGE037
为二值化处理后的B-scan图像，其中值为1的元素为前景元素，值为0的元素为背景元素。在得到二值化处理后的B-scan图像后，需要对其进行形态学滤波，以得到预掩膜矩阵。下面结合具体应用场景，对形态学滤波的方法做详细的介绍。
在一个实施例中，上述步骤S12中对二值化处理后的B-scan图像进行形态学滤波包括：膨胀-腐蚀滤波，腐蚀处理，膨胀处理。
腐蚀处理和膨胀处理是形态学滤波的基本方法之一。腐蚀处理包括：设定被腐蚀图像
Figure 586058DEST_PATH_IMAGE038
和滤波器
Figure 790774DEST_PATH_IMAGE039
，记
Figure 661778DEST_PATH_IMAGE038
中所有前景元素的集合为
Figure 533919DEST_PATH_IMAGE040
。令
Figure 996124DEST_PATH_IMAGE039
在
Figure 416698DEST_PATH_IMAGE040
中所有元素上依次平移，当
Figure 91393DEST_PATH_IMAGE039
移动到集合
Figure 818041DEST_PATH_IMAGE041
中元素
Figure 388831DEST_PATH_IMAGE042
上时，若此时
Figure 568139DEST_PATH_IMAGE039
中的前景元素对应
Figure 312104DEST_PATH_IMAGE038
上的元素同样为前景元素，则仍将
Figure 627679DEST_PATH_IMAGE043
置为前景元素，否则将
Figure 431687DEST_PATH_IMAGE043
置为背景元素。膨胀处理包括：设定被膨胀图像
Figure 832712DEST_PATH_IMAGE038
和滤波器
Figure 114789DEST_PATH_IMAGE039
，记
Figure 550450DEST_PATH_IMAGE038
中所有前景元素的集合为
Figure 994201DEST_PATH_IMAGE044
。令
Figure 413681DEST_PATH_IMAGE039
在
Figure 233869DEST_PATH_IMAGE044
中所有元素上依次平移，当
Figure 524036DEST_PATH_IMAGE039
移动到集合
Figure 138688DEST_PATH_IMAGE044
中元素
Figure 508447DEST_PATH_IMAGE043
上时，即将
Figure 132326DEST_PATH_IMAGE038
中被
Figure 542579DEST_PATH_IMAGE039
前景元素覆盖的位置为前景元素。膨胀-腐蚀滤波是使用同一滤波器，首先进行一次膨胀处理，然后进行一次腐蚀处理的操作。
所述形态学滤波包括：首先采用矩形滤波器对二值化处理后的B-scan图像进行第一设定次数的膨胀-腐蚀滤波，以得到初次滤波后的B-scan图像；然后分别采用上斜滤波器和下斜滤波器对初次滤波后的B-scan图像进行第二设定次数的腐蚀处理和第三设定次数的膨胀处理，以得到预掩膜矩阵。
在本实施例中，上述矩形滤波器、上斜滤波器和下斜滤波器的结构分别为：
Figure 859290DEST_PATH_IMAGE045
、
Figure 722204DEST_PATH_IMAGE046
、
Figure 884195DEST_PATH_IMAGE047
通过本实施例的设置方式对二值化处理后的B-scan图像进行形态学滤波，可以得到预掩膜矩阵
Figure 414534DEST_PATH_IMAGE048
。
上文中对形态学滤波的方法做了示例性的说明，下面结合具体应用场景，对上述步骤S13中设定背景介质相对介电常数估计区间上限值、下限值、估计值及其更新方法，通过并行计算获得并行掩膜后向投影成像结果做详细的介绍。
在一个实施例中，上述步骤S13中获取背景介质相对介电常数的估计值包括：根据背景介质性质设定相对介电常数估计值的估计区间上限值和下限值，将区间中值作为背景介质相对介电常数估计值。
以获取初始的取背景介质相对介电常数的估计值为例，假设背景介质常数的初始估计区间为
Figure 636568DEST_PATH_IMAGE049
其中
Figure 986778DEST_PATH_IMAGE050
为该初始估计区间的下限值，
Figure 952460DEST_PATH_IMAGE051
为该初始估计区间的上限值，则可以通过如下公式获取初始的背景介质相对介电常数的估计值：
Figure 71725DEST_PATH_IMAGE052
（7）
在一个实施例中，上述步骤S13中通过并行计算获得并行掩膜后向投影成像结果包括：
在步骤S21中，以初始背景介质相对介电常数的估计值
Figure 730240DEST_PATH_IMAGE053
为背景介质相对介电常数、以天线步进间隔
Figure 567746DEST_PATH_IMAGE054
为成像区单元格划分单位对成像区域进行划分，得到一个维度为N×M的用于存储成像结果的二维矩阵，并设该二维矩阵为
Figure 71539DEST_PATH_IMAGE055
，其中M的计算公式为：
Figure 45312DEST_PATH_IMAGE056
（8）
在步骤S22中，生成B-scan回波数据的深度索引序列和成像区域的深度索引序列，设B-scan回波数据的深度索引序列为
Figure 874727DEST_PATH_IMAGE057
，成像区域的深度索引序列为
Figure 933950DEST_PATH_IMAGE058
，则：
Figure 507014DEST_PATH_IMAGE059
（9）
Figure 595013DEST_PATH_IMAGE060
（10）
在步骤S23中，根据B-scan回波数据和成像区域的一维深度索引序列确定时延-深度对应关系。对B-scan回波数据的深度索引序列
Figure 595330DEST_PATH_IMAGE061
中每一元素
Figure 141849DEST_PATH_IMAGE062
，都对应于成像区域的深度索引序列
Figure 987445DEST_PATH_IMAGE063
中的区间
Figure 935809DEST_PATH_IMAGE064
，其中low的值满足如下条件：
Figure 841448DEST_PATH_IMAGE065
（11）
Figure 875264DEST_PATH_IMAGE066
（12）
high的值满足如下条件：
Figure 790130DEST_PATH_IMAGE067
（13）
Figure 124159DEST_PATH_IMAGE068
（14）
在步骤S24中，根据时延-深度对应关系、预掩膜矩阵
Figure 200700DEST_PATH_IMAGE069
和B-scan回波数据与成像区域的一维深度索引序列，对预掩膜矩阵
Figure 721811DEST_PATH_IMAGE070
进行映射获得成像掩膜矩阵，本实施例中将成像掩膜矩阵设为
Figure 174789DEST_PATH_IMAGE071
，则该映射方法为：
Figure 97746DEST_PATH_IMAGE072
（15）
在步骤S25中，生成同成像掩膜矩阵
Figure 345187DEST_PATH_IMAGE073
相同维度的成像掩膜扩展矩阵
Figure 353595DEST_PATH_IMAGE074
，并根据成像掩膜矩阵设定成像掩膜扩展矩阵各单元格的值，设定的方法为：
当满足如下条件时，将
Figure 344684DEST_PATH_IMAGE075
的值设置为1，否则设置为0；
Figure 122148DEST_PATH_IMAGE076
（16）
在步骤S26中，计算成像区域中每个成像单元的成像结果，从而得到整个成像区域的成像结果
Figure 558069DEST_PATH_IMAGE077
。以成像区域的一个成像单元
Figure 319351DEST_PATH_IMAGE078
为例，其成像结果计算方法为：如果
Figure 848553DEST_PATH_IMAGE079
，则并行计算（如可以通过numba开源工具库进行计算）得出该成像单元对每一测点
Figure 11681DEST_PATH_IMAGE080
对应的电磁波传播双程时延
Figure 600925DEST_PATH_IMAGE081
，根据
Figure 583925DEST_PATH_IMAGE081
提取测点
Figure 182396DEST_PATH_IMAGE080
处A-scan回波数据的对应幅值并生成一维信号序列
Figure 934452DEST_PATH_IMAGE082
，将
Figure 694597DEST_PATH_IMAGE083
的和作为成像单元
Figure 899314DEST_PATH_IMAGE084
的成像结果。
上文中对设定背景介质相对介电常数估计区间上限值、下限值、估计值及其更新方法，通过并行计算获得并行掩膜后向投影成像结果做了示例性的说明，下面结合具体应用场景，对上述步骤S14中根据二分法和所述并行掩膜后向投影成像结果，更新背景介质相对介电常数估计区间，计算背景介质相对介电常数估计值；迭代运算，直至背景介质相对介电常数估计区间满足设定条件，获得背景介质相对介电常数最终估计值；获得自聚焦并行掩膜后向投影成像结果做详细的介绍。
在一个实施例中，上述步骤S14中根据二分法和所述并行掩膜后向投影结果，更新背景介质相对介电常数估计区间，计算背景介质相对介电常数估计值的方法包括：
在步骤S31中，根据背景介质相对介电常数估计值
Figure 35897DEST_PATH_IMAGE085
与B-scan图像双曲线特征顶点坐标序列L，得出成像中地下目标坐标序列为：
Figure 908038DEST_PATH_IMAGE086
（17）
其中
Figure 839085DEST_PATH_IMAGE087
（18）
其中，
Figure 531097DEST_PATH_IMAGE088
是第i个地下目标在B-scan图像中双曲线特征顶点的时延值。
在步骤S32中，根据背景介质相对介电常数估计值
Figure 471371DEST_PATH_IMAGE085
、预掩膜矩阵
Figure 932440DEST_PATH_IMAGE089
及探地雷达探测时的工作参数获得并行掩膜后向投影成像结果
Figure 34388DEST_PATH_IMAGE090
，计算
Figure 479276DEST_PATH_IMAGE090
中对成像中地下目标上方区域聚焦度
Figure 951802DEST_PATH_IMAGE091
和下方区域聚焦度
Figure 532956DEST_PATH_IMAGE092
，计算方法为：
Figure 71385DEST_PATH_IMAGE093
（19）
Figure 472411DEST_PATH_IMAGE094
（20）
在步骤S33中，根据二分法，如果
Figure 20066DEST_PATH_IMAGE095
，则更新相对介电常数估计区间的下限值为
Figure 190148DEST_PATH_IMAGE085
；如果
Figure 899478DEST_PATH_IMAGE096
，则更新相对介电常数估计区间的上限值为
Figure 787799DEST_PATH_IMAGE085
，从而实现对相对介电常数估计区间的更新，然后根据更新后的相对介电常数估计区间的上限值和下限值，计算背景介质相对介电常数的估计值为更新后的背景介质相对介电常数估计区间中值。
在一个实施例中，上述步骤S14中迭代运算，直至背景介质相对介电常数估计区间满足设定条件，获得背景介质相对介电常数最终估计值的方法包括：
根据所述并行掩膜后向投影结果，更新背景介质相对介电常数估计区间，计算背景介质相对介电常数估计值；
根据时延-深度对应关系和所述预掩膜矩阵、探测时的工作参数与更新的背景介质相对介电常数估计值，通过并行计算，获得新的并行掩膜后向投影成像结果；
根据新的并行掩膜后向投影成像结果，更新背景介质相对介电常数估计区间，重复上述运算，直至背景介质相对介电常数估计区间满足设定条件；
所述设定条件包括：更新后的背景介质相对介电常数估计区间的上限值与下限值之差小于设定差值；
所述背景介质相对介电常数最终估计值为背景介质相对介电常数估计区间满足设定条件时的区间中值；
在一个实施例中，上述步骤S14中获得自聚焦并行掩膜后向投影成像结果的方法包括：
根据时延-深度对应关系与所述背景介质相对介电常数最终估计值、预掩膜矩阵与探地雷达探测时的工作参数，获得成像掩膜矩阵
Figure 873567DEST_PATH_IMAGE097
，设置与
Figure 163734DEST_PATH_IMAGE097
相同维度的二维矩阵
Figure 43965DEST_PATH_IMAGE098
为成像区域；
通过并行计算，对
Figure 419583DEST_PATH_IMAGE099
中对应
Figure 43463DEST_PATH_IMAGE097
前景元素相同位置的成像单元进行后向投影成像，获得自聚焦并行掩膜后向投影成像结果。
下面以一个实验为例，对本申请技术方案的技术效果做详细的介绍。
在本实验中使用设备为配备2GHz天线的SIR4000探地雷达对探测区域进行合成孔径扫描，数据处理平台配置Intel i7-7700中央处理器和NVIDIA GeForce GTX 1060图形处理器，主存容量为16GB。探测区域背景介质为碳酸钙细砂，相对介电常数
Figure 188136DEST_PATH_IMAGE100
，
Figure 504848DEST_PATH_IMAGE101
为空气-地面分界线，上方为空气，下方为探测区域。在探测区域（0.318m，0.162m）处放置一根直径为0.008m的钢筋，在（0.530m，0.105m）处放置一根直径为0.010m的钢筋，两根钢筋长度约为0.3m。接收天线与发射天线沿
Figure 367762DEST_PATH_IMAGE101
的测线对探测区域进行合成孔径扫描，扫描孔径从
Figure 795332DEST_PATH_IMAGE102
开始，到
Figure 54232DEST_PATH_IMAGE103
终止，各孔径间距为
Figure 276266DEST_PATH_IMAGE104
，共有887个孔径点。在每个孔径点，发射天线向探测区域发射电磁波，接收天线接收来自探测区域的散射回波，共887道A-scan回波数据，将所有A-scan回波数据沿扫描顺序组合为B-scan回波数据，成像如图3所示。对该B-scan回波数据进行预处理，预处理结果如图4所示。观测可得B-scan回波数据中存在2个双曲线特征，顶点位置为（385，318）与（256，530）。
由于已验知背景介质相对介电常数，对该B-scan回波数据使用现有的方法成像，成像的三维结果如图5所示，共耗时5663秒。
对预处理后的B-scan回波数据进行二值化处理和形态学滤波，得到预掩膜矩阵，对预掩膜矩阵成像如图6所示，该步骤耗时小于1秒。
现假设未验知背景介质相对介电常数，根据上述获取以背景介质相对介电常数的精准估计值，考虑无任何背景介质的信息的极端情况，设定背景相对介电常数取值上限
Figure DEST_PATH_IMAGE105
为100，下限
Figure 564159DEST_PATH_IMAGE106
为1，记录迭代求解过程如表1所示。
表1
Figure DEST_PATH_IMAGE107
在第10次迭代后
Figure 467524DEST_PATH_IMAGE108
时退出迭代，迭代过程中计算
Figure DEST_PATH_IMAGE109
和
Figure 790052DEST_PATH_IMAGE110
如图7中的框选区域所示，其中的a部分-j部分的横坐标均为测线维度，纵坐标均为深度维度。通过本发明技术方案得到的背景介质相对介电常数估计值为2.26，与真实值偏离度为1.7%，该步骤约耗时约15秒。
采用本申请的技术方案，根据对背景介质相对介电常数估计过程中最后一次迭代中的成像掩膜矩阵与成像结果获得探测区域并行掩膜后向投影成像，成像的三维结果如图8所示，从图8中可知，成像质量较传统方法明显提高，成像效果显著优于传统方法，并且算进过程的总耗时为15秒，较传统方法提高约400倍。
根据本发明的第二方面，本发明还提供了计算机设备，如图9所示，该设备包括处理器、存储器、通信接口和通信总线，处理器、存储器和通信接口通过通信总线完成相互间的通信。处理器用于提供计算和控制能力。存储器包括非易失性存储介质、内存储器。该非易失性存储介质存储有操作系统和计算机程序指令。该内存储器为非易失性存储介质中的操作系统和计算机程序指令的运行提供环境。上述装置的通信接口用于与外部的终端进行有线或无线方式的通信，无线方式可通过WIFI、运营商网络、NFC(近场通信)或其他技术实现。存储器用于存储计算机程序指令，该计算机程序指令使处理器执行上述本发明如图2所示的探地雷达自聚焦并行掩膜后向投影成像方法。
根据本发明的第三方面，本发明还提供了一种计算机可读存储介质。本领域普通技术人员可以理解实现上述实施例方法中的全部或部分流程，是可以通过计算机程序指令来指令相关的硬件来完成，所述的计算机程序指令可存储于一非易失性计算机可读取存储介质中，该计算机程序指令在执行时，可实现上述本发明如图2所示的探地雷达自聚焦并行掩膜后向投影成像方法。其中，本发明所提供的各实施例中所使用的对存储器、存储、数据库或其它介质的任何引用，均可包括非易失性和/或易失性存储器。非易失性存储器可包括只读存储器(ROM)、可编程ROM(PROM)、电可编程ROM(EPROM)、电可擦除可编程ROM(EEPROM)或闪存。易失性存储器可包括随机存取存储器(RAM)或者外部高速缓冲存储器。作为说明而非局限，RAM以多种形式可得，诸如静态RAM(SRAM)、动态RAM(DRAM)、同步DRAM(SDRAM)、双数据率SDRAM(DDRSDRAM)、增强型SDRAM(ESDRAM)、同步链路(Synchlink)DRAM(SLDRAM)、存储器总线(Rambus)直接RAM(RDRAM)、直接存储器总线动态RAM(DRDRAM)、以及存储器总线动态RAM(RDRAM)等。
本说明书中所使用的术语“第一”或“第二”等用于指代编号或序数的术语仅用于描述目的，而不能理解为明示或暗示相对重要性或者隐含指明所指示的技术特征的数量。由此，限定有“第一”或“第二”的特征可以明示或者隐含地包括至少一个该特征。在本说明书的描述中，“多个”的含义是至少两个，例如两个，三个或更多个等，除非另有明确具体地限定。
虽然本说明书已经示出和描述了本发明的多个实施例，但对于本领域技术人员显而易见的是，这样的实施例只是以示例的方式提供的。本领域技术人员会在不偏离本发明思想和精神的情况下想到许多更改、改变和替代的方式。应当理解的是在实践本发明的过程中，可以采用本文所描述的本发明实施例的各种替代方案。所附权利要求书旨在限定本发明的保护范围，并因此覆盖这些权利要求保护范围内的模块组成、等同或替代方案。

Claims (10)

1.一种探地雷达自聚焦并行掩膜后向投影成像方法，其特征在于，包括：
获取B-scan回波数据，以及探测时的工作参数，对B-scan回波数据进行去直耦波和校正零点预处理，并进行归一化处理并成像，获得B-scan图像；
根据优化函数获取最优阈值，根据该最优阈值对所述B-scan图像进行二值化处理，并对二值化处理后的B-scan图像进行形态学滤波，以得到预掩膜矩阵；
根据背景介质性质，设定背景介质相对介电常数估计区间上限值与下限值，并将该估计区间中值作为背景介质相对介电常数估计值；
根据所述背景介质相对介电常数估计值、预掩膜矩阵和探测时的工作参数，通过并行计算获得并行掩膜后向投影成像结果；
根据二分法和所述并行掩膜后向投影成像结果，更新背景介质相对介电常数估计区间，计算背景介质相对介电常数估计值，再次通过并行计算获得并行掩膜后向投影成像结果并更新背景介质相对介电常数估计区间，迭代运算直至背景介质相对介电常数估计区间满足设定条件，获得背景介质相对介电常数最终估计值和自聚焦并行掩膜后向投影成像结果。
2.根据权利要求1所述的探地雷达自聚焦并行掩膜后向投影成像方法，其特征在于，所述对B-scan回波数据进行去直耦波和校正零点预处理，包括：
对每行B-scan回波数据减去该行均值去除B-scan回波数据中的直耦波；
将目标按一定深度埋在已知相对介电常数的均匀介质中，在其正上方进行探测，提取回波数据中对应该目标幅值的响应时延值，计算回波数据中对应该目标幅值的理想时延值，校正使响应时延值与理想时延值相等。
3.根据权利要求1所述的探地雷达自聚焦并行掩膜后向投影成像方法，其特征在于，所述根据优化函数获取最优阈值，包括：
获取所有待选阈值；
获取所述各待选阈值所对应的优化函数输出值；
根据各优化函数的输出值，从所述待选阈值中确定所述最优阈值。
4.根据权利要求3所述的探地雷达自聚焦并行掩膜后向投影成像方法，其特征在于，所述获取所述各待选阈值所对应的优化函数输出值，包括：
根据待选阈值对所述B-scan回波数据进行分类，以得到第一类数据和第二类数据；
获取所述第一类数据的特征值和第二类数据的特征值，并根据所述第一类数据的特征值和第二类数据的特征值，得到该阈值对应的优化函数输出值；
所述根据各优化函数的输出值，从所述待选阈值中确定所述最优阈值，包括：获取优化函数最大输出值所对应的阈值，并将该阈值作为所述最优阈值。
5.根据权利要求1所述的探地雷达自聚焦并行掩膜后向投影成像方法，其特征在于，所述对二值化处理后的B-scan图像进行形态学滤波包括：
采用矩形滤波器对所述二值化处理后的B-scan图像进行第一设定次数的膨胀-腐蚀滤波，以得到初次滤波后的B-scan图像；
分别采用上斜滤波器和下斜滤波器对初次滤波后的B-scan图像进行第二设定次数的腐蚀处理和第三设定次数的膨胀处理，以得到预掩膜矩阵。
6.根据权利要求1所述的探地雷达自聚焦并行掩膜后向投影成像方法，其特征在于，根据所述背景介质相对介电常数估计值、预掩膜矩阵和探测时的工作参数，通过并行计算获得并行掩膜后向投影成像结果，包括：
根据时延-深度对应关系和所述背景介质相对介电常数估计值、预掩膜矩阵和探测时的工作参数，获得成像掩膜矩阵，根据成像掩膜矩阵获得掩膜扩展矩阵，设置与掩膜扩展矩阵同维度的成像区域，通过并行计算，对成像区域中对应掩膜扩展矩阵前景元素相同位置的成像单元进行后向投影成像，以获得并行掩膜后向投影成像结果。
7.根据权利要求1所述的探地雷达自聚焦并行掩膜后向投影成像方法，其特征在于，根据二分法和所述并行掩膜后向投影成像结果，更新背景介质相对介电常数估计区间，计算背景介质相对介电常数估计值；
根据所述并行掩膜后向投影成像结果，计算出地下目标上方区域聚焦度和地下目标下方区域聚焦度；
根据地下目标上方区域聚焦度和地下目标下方区域聚焦度对所述背景介质相对介电常数估计区间上限值或下限值进行更新，将更新后的估计区间中值作为背景介质相对介电常数新的估计值，包括增加所述背景介质相对介电常数的估计值和减小所述背景介质相对介电常数的估计值；
所述更新背景介质相对介电常数估计区间包括：根据二分法，若所述地下目标上方区域的聚焦度大于地下目标下方区域的聚焦度，则更新所述背景介质相对介电常数估计区间下限值为当前背景介质相对介电常数估计值，否则更新所述背景介质相对介电常数估计范围上限值为当前背景介质相对介电常数估计值；
所述计算背景介质相对介电常数的估计值包括：计算背景介质相对介电常数估计区间中值，该值即为背景介质相对介电常数估计值。
8.根据权利要求1所述的探地雷达自聚焦并行掩膜后向投影成像方法，其特征在于，根据二分法和所述并行掩膜后向投影成像结果，迭代运算直至背景介质相对介电常数估计区间满足设定条件，获得背景介质相对介电常数最终估计值和自聚焦并行掩膜后向投影成像结果，包括：
所述迭代运算包括：根据所述并行掩膜后向投影成像结果，更新背景介质相对介电常数估计区间，将更新后的估计区间中值作为背景介质相对介电常数估计值，根据所述背景介质相对介电常数估计值与预掩膜矩阵，再次并行计算获得并行掩膜后向投影结果，并再次更新背景介质相对介电常数估计区间，直至背景介质相对介电常数估计区间满足设定条件；
所述设定条件包括：更新后的背景介质相对介电常数估计区间的上限值与下限值之差小于设定差值；
所述背景介质相对介电常数最终估计值为背景介质相对介电常数估计区间满足设定条件时的区间中值；
所述获得自聚焦并行掩膜后向投影结果包括：根据时延-深度关系与所述背景介质相对介电常数最终估计值与预掩膜矩阵，获得成像掩膜矩阵，设置与成像掩膜矩阵同维度的成像区域，通过并行计算，对成像区域中对应成像掩膜矩阵前景元素相同位置的成像单元进行后向投影成像，获得自聚焦并行掩膜后向投影成像结果。
9.一种计算机设备，其特征在于，包括处理器和存储器，所述存储器上存储有用于在所述处理器上执行的计算机程序指令，所述处理器执行该计算机程序指令时，实现如权利要求1-8任意一项所述的探地雷达自聚焦并行掩膜后向投影成像方法。
10.一种计算机可读存储介质，其特征在于，其存储有计算机指令，所述计算机指令由处理器运行时，实现根据权利要求1-8的任意一项所述的探地雷达自聚焦并行掩膜后向投影成像方法。
CN202210992032.9A 2022-08-18 2022-08-18 探地雷达自聚焦并行掩膜后向投影成像方法、设备及介质 Active CN115061138B (zh)

Priority Applications (1)

Application Number Priority Date Filing Date Title
CN202210992032.9A CN115061138B (zh) 2022-08-18 2022-08-18 探地雷达自聚焦并行掩膜后向投影成像方法、设备及介质

Applications Claiming Priority (1)

Application Number Priority Date Filing Date Title
CN202210992032.9A CN115061138B (zh) 2022-08-18 2022-08-18 探地雷达自聚焦并行掩膜后向投影成像方法、设备及介质

Publications (2)

Publication Number Publication Date
CN115061138A CN115061138A (zh) 2022-09-16
CN115061138B true CN115061138B (zh) 2022-10-28

Family

ID=83208097

Family Applications (1)

Application Number Title Priority Date Filing Date
CN202210992032.9A Active CN115061138B (zh) 2022-08-18 2022-08-18 探地雷达自聚焦并行掩膜后向投影成像方法、设备及介质

Country Status (1)

Country Link
CN (1) CN115061138B (zh)

Families Citing this family (2)

* Cited by examiner, † Cited by third party
Publication number Priority date Publication date Assignee Title
CN116819482B (zh) * 2023-08-28 2023-11-10 四川省石棉县恒达粉体材料有限责任公司 一种基于雷达数据的方解石探测方法
CN117310696A (zh) * 2023-09-26 2023-12-29 中南大学 一种探地雷达自聚焦后向投影成像方法及装置

Citations (6)

* Cited by examiner, † Cited by third party
Publication number Priority date Publication date Assignee Title
US6198424B1 (en) * 1999-01-21 2001-03-06 Rosemount Inc. Multiple process product interface detection for a low power radar level transmitter
CN102830401A (zh) * 2012-08-27 2012-12-19 中南大学 一种探地雷达加窗加权后向投影成像方法
CN103913741A (zh) * 2014-03-18 2014-07-09 电子科技大学 一种合成孔径雷达高效自聚焦后向投影bp方法
CN104833973A (zh) * 2015-05-08 2015-08-12 电子科技大学 基于半正定规划的线阵sar后向投影自聚焦成像方法
CN113361367A (zh) * 2021-06-01 2021-09-07 中南大学 基于深度学习的地下目标电磁反演方法及系统
CN113671496A (zh) * 2021-08-27 2021-11-19 中国人民解放军国防科技大学 一种基于图像滤波器的mimo雷达自聚焦三维成像方法

Family Cites Families (1)

* Cited by examiner, † Cited by third party
Publication number Priority date Publication date Assignee Title
US6927698B2 (en) * 2001-08-27 2005-08-09 Larry G. Stolarczyk Shuttle-in receiver for radio-imaging underground geologic structures

Patent Citations (6)

* Cited by examiner, † Cited by third party
Publication number Priority date Publication date Assignee Title
US6198424B1 (en) * 1999-01-21 2001-03-06 Rosemount Inc. Multiple process product interface detection for a low power radar level transmitter
CN102830401A (zh) * 2012-08-27 2012-12-19 中南大学 一种探地雷达加窗加权后向投影成像方法
CN103913741A (zh) * 2014-03-18 2014-07-09 电子科技大学 一种合成孔径雷达高效自聚焦后向投影bp方法
CN104833973A (zh) * 2015-05-08 2015-08-12 电子科技大学 基于半正定规划的线阵sar后向投影自聚焦成像方法
CN113361367A (zh) * 2021-06-01 2021-09-07 中南大学 基于深度学习的地下目标电磁反演方法及系统
CN113671496A (zh) * 2021-08-27 2021-11-19 中国人民解放军国防科技大学 一种基于图像滤波器的mimo雷达自聚焦三维成像方法

Non-Patent Citations (2)

* Cited by examiner, † Cited by third party
Title
A Multi-Scale Weighted Back Projection Imaging Technique for Ground Penetrating Radar Applications;Wentai Lei 等;《remote sensing》;20140605;全文 *
面向探地雷达B-scan图像的目标检测算法综述;侯斐斐 等;《电子与信息学报》;20200131;第42卷(第1期);全文 *

Also Published As

Publication number Publication date
CN115061138A (zh) 2022-09-16

Similar Documents

Publication Publication Date Title
CN115061138B (zh) 探地雷达自聚焦并行掩膜后向投影成像方法、设备及介质
Merazi-Meksen et al. Mathematical morphology for TOFD image analysis and automatic crack detection
Oliveri et al. A nested multi-scaling inexact-Newton iterative approach for microwave imaging
Chen et al. Probabilistic robust hyperbola mixture model for interpreting ground penetrating radar data
EP2681589A1 (en) Methods for deblending of seismic shot gathers
CN111290021A (zh) 基于梯度结构张量谱分解的碳酸盐岩缝洞增强识别方法
CN110031898B (zh) 数据优化方法及积分法叠前深度偏移方法
CN113533867B (zh) 基于Fourier插值的远场方向图快速测量方法
CN110837079B (zh) 一种基于雷达的目标检测方法及装置
CN114003981A (zh) 一种基于时空一体数字地球的电磁频谱可视化分析方法
Maisto et al. Transverse resolution in microwave imaging for strip objects buried in a half-space medium
CN116931101A (zh) 地下浅层病害体探测方法、装置、电子设备及存储介质
CN116908853B (zh) 高相干点选取方法、装置和设备
CN110319797B (zh) 一种基于双通道信息的图像域富钴结壳厚度提取方法
CN116503281A (zh) 图像滤波方法、装置、电子设备及存储介质
CN115311532A (zh) 一种基于ResNet网络模型的探地雷达地下空洞目标自动识别方法
CN111399038B (zh) 斜率参数提取方法、装置及计算机可读存储介质
CN112630840B (zh) 基于统计特征参数的随机反演方法及处理器
Bitenc et al. Potential and limitations of terrestrial laser scanning for discontinuity roughness estimation
CN110764062A (zh) 基于分数阶傅里叶域滤波的多分量线性调频信号参数估计方法
Zhang et al. Non-local Phase Linking for Distributed Scatterer Interferometry
CN111177886B (zh) 基于物探数据分析的海上测距规划及土厚预测的方法
CN117079147B (zh) 一种道路内部病害识别方法、电子设备及存储介质
CN115236415B (zh) 天线增益获取方法、装置、电子设备及存储介质
Milisavljevic et al. Improving Mine Recognition through Processing and Dempster–Shafer Fusion of Multisensor Data

Legal Events

Date Code Title Description
PB01 Publication
PB01 Publication
SE01 Entry into force of request for substantive examination
SE01 Entry into force of request for substantive examination
GR01 Patent grant
GR01 Patent grant