WO2015101058A1 - 平面地表环境中地下建筑的红外成像探测定位方法 - Google Patents

平面地表环境中地下建筑的红外成像探测定位方法 Download PDF

Info

Publication number
WO2015101058A1
WO2015101058A1 PCT/CN2014/085713 CN2014085713W WO2015101058A1 WO 2015101058 A1 WO2015101058 A1 WO 2015101058A1 CN 2014085713 W CN2014085713 W CN 2014085713W WO 2015101058 A1 WO2015101058 A1 WO 2015101058A1
Authority
WO
WIPO (PCT)
Prior art keywords
infrared image
image
underground building
infrared
underground
Prior art date
Application number
PCT/CN2014/085713
Other languages
English (en)
French (fr)
Inventor
张天序
郝龙伟
鲁岑
马文绚
王岳环
桑农
杨卫东
Original Assignee
华中科技大学
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 华中科技大学 filed Critical 华中科技大学
Priority to US15/105,456 priority Critical patent/US9741120B2/en
Publication of WO2015101058A1 publication Critical patent/WO2015101058A1/zh

Links

Images

Classifications

    • G—PHYSICS
    • G06—COMPUTING; CALCULATING OR COUNTING
    • G06T—IMAGE DATA PROCESSING OR GENERATION, IN GENERAL
    • G06T7/00—Image analysis
    • G06T7/70—Determining position or orientation of objects or cameras
    • G—PHYSICS
    • G01—MEASURING; TESTING
    • G01J—MEASUREMENT OF INTENSITY, VELOCITY, SPECTRAL CONTENT, POLARISATION, PHASE OR PULSE CHARACTERISTICS OF INFRARED, VISIBLE OR ULTRAVIOLET LIGHT; COLORIMETRY; RADIATION PYROMETRY
    • G01J5/00—Radiation pyrometry, e.g. infrared or optical thermometry
    • G01J5/10—Radiation pyrometry, e.g. infrared or optical thermometry using electric radiation detectors
    • G—PHYSICS
    • G06—COMPUTING; CALCULATING OR COUNTING
    • G06T—IMAGE DATA PROCESSING OR GENERATION, IN GENERAL
    • G06T5/00—Image enhancement or restoration
    • G—PHYSICS
    • G01—MEASURING; TESTING
    • G01J—MEASUREMENT OF INTENSITY, VELOCITY, SPECTRAL CONTENT, POLARISATION, PHASE OR PULSE CHARACTERISTICS OF INFRARED, VISIBLE OR ULTRAVIOLET LIGHT; COLORIMETRY; RADIATION PYROMETRY
    • G01J5/00—Radiation pyrometry, e.g. infrared or optical thermometry
    • G01J2005/0077—Imaging
    • G—PHYSICS
    • G06—COMPUTING; CALCULATING OR COUNTING
    • G06T—IMAGE DATA PROCESSING OR GENERATION, IN GENERAL
    • G06T2207/00—Indexing scheme for image analysis or image enhancement
    • G06T2207/10—Image acquisition modality
    • G06T2207/10048—Infrared image
    • G—PHYSICS
    • G06—COMPUTING; CALCULATING OR COUNTING
    • G06T—IMAGE DATA PROCESSING OR GENERATION, IN GENERAL
    • G06T2207/00—Indexing scheme for image analysis or image enhancement
    • G06T2207/20—Special algorithmic details
    • G06T2207/20076—Probabilistic image processing

Definitions

  • the invention belongs to the field of interdisciplinary theory of geophysics and remote sensing technology, and more specifically relates to an infrared imaging detection and positioning method for underground buildings in a planar surface environment.
  • the physical basis is that a large amount of sunlight energy is irradiated to the soil to absorb heat.
  • the heated soil emits infrared radiation and is detected by a thermal infrared sensor.
  • the daily heating and cooling of natural solar energy has different effects on the buried object and the land surrounding it, resulting in a detectable temperature difference.
  • the presence of underground buildings can produce surface water/heat distribution anomalies that can be used to detect, detect, verify, and confirm underground structures.
  • Infrared imaging technology is used to detect underground buildings at home and abroad, but they all directly use infrared imaging sensors for imaging, and artificially interpret the infrared images obtained by them, which has great limitations.
  • the limitations are as follows: First, the thermal anomaly of underground buildings is modulated by the conduction of the buried layers, and the heat distribution to the surface has changed greatly from the underground buildings, showing heat. The diffusion is severe, the temperature difference is greatly reduced, and the heat signal is weak. Second, through the conduction modulation of the formation, the location of the underground building reflected by the surface thermal signal may change, difficult to find and locate. Third, the difficulty of manual interpretation is not conducive to the accurate detection and location of underground buildings.
  • the present invention proposes a new method for detecting underground buildings.
  • the ground layer conduction modulation in the infrared image is solved.
  • the problem of artificial detection and positioning is impossible to detect and locate the underground building in the plane environment.
  • the invention provides an infrared imaging detection and positioning method for an underground building in a plane ground environment, which is inversely modulated according to a Gaussian model of energy diffusion of an underground building, and an original infrared image formed by the underground building after being modulated by the ground layer.
  • Obtaining a target image of the underground building the method comprising the steps of:
  • step (3) It is judged whether the iteration termination condition is satisfied. If it is satisfied, the target image f n obtained by the current iteration is the final target image f; if not, the process returns to step (3), and the iterative calculation is continued.
  • the step (1) is specifically:
  • the iterative solution of the thermal expansion function h n and the target image f n in the step (3) is specifically calculated according to the following two equations:
  • n is the current number of iterations
  • (x, y) belongs to the thermal diffusion function support domain
  • (i, j) belongs to the image support domain
  • f(i, j) represents the target image
  • g(i, j) represents the local infrared image
  • h(x, y) is a thermal diffusion function
  • h(-x, -y) represents the conjugate of h(x, y)
  • f(-i, -j) represents the conjugate of f(i,j)
  • the iterative termination condition in the step (2) is the number of iteration terminations n>N 0 or the error ⁇ , and the step (4) determines whether the iteration condition is satisfied or not:
  • the method of the invention By inversely modulating the infrared image formed by the underground building after being modulated by the ground layer, the method of the invention not only makes the infrared image display of the original underground building clearer, but also can invert the real structure of the underground building. According to the known infrared information in the obtained infrared image, the key physical characteristics of the underground building can be inferred, and the purpose of "seeing” and further quantitatively measuring the key physical features of the underground building can be achieved.
  • Figure 1 is a schematic flow chart of the present invention
  • Figure 2 shows the surface heat distribution and classification results of underground buildings (depth 100m);
  • Figure 3 is a schematic diagram of a two-dimensional Gaussian function form
  • Figure 4(a) shows the infrared thermal imaging of the scaled model at 19:35;
  • Figure 4 (b) is an infrared thermal imaging of the scaled model at 20:00;
  • Figure 4(c) shows the infrared thermal imaging of the scaled model at 20:20;
  • Figure 4(d) shows the infrared thermal imaging of the scaled model at 20:45;
  • Figure 4(e) shows the infrared thermal imaging of the scaled model at 21:00
  • Figure 4(f) shows the infrared thermal imaging of the scaled model at 21:30;
  • Figure 5 is a flow chart of maximum likelihood estimation image recovery
  • Figure 6 shows the surface heat distribution (upper right) and the inverse transformation result (bottom right) of a typical underground structure in a flat, depth 100m environment;
  • Figure 7 (a) is the original infrared image obtained by the sand burial experiment 1;
  • Figure 7 (b) is the result of the anti-modulation treatment of the sand burial experiment
  • Figure 9 (a) is the original infrared image of the sand buried experiment 2;
  • Fig. 9(b) shows the results of the anti-modulation treatment of the sand burial experiment.
  • the invention preliminarily classifies infrared images of underground buildings and underground buildings without existing in a planar environment, and then performs inverse modulation processing on infrared images existing in underground buildings, and the detection problem is attributed to the inverse problem in mathematical physics. Detecting the structural information and location of the underground building, The detection and positioning of the underground building target is realized, and the same report as the present invention has not been seen in the existing domestic and foreign literatures.
  • the method for detecting and locating an underground building in a plane ground environment comprises the following steps:
  • the present invention has carried out theoretical analysis on the infrared characteristics of underground buildings and their relationship with the environment. Modeling and simulation of typical underground buildings using existing software. In order to display the infrared information of the ground surface after the underground building was modulated by the ground layer, the thermal distribution of the underground building with a depth of 100 meters in the plane ground environment was simulated by ANSYS. The underground building thermal radiation simulation is realized according to the thermal radiation analysis in the "ANSYS 12.0 Thermal Analysis Engineering Application Practical Manual" edited by Zhang Zhaohui.
  • the infrared image obtained includes the infrared information of the target area and the infrared information of the non-target area, and the input infrared image is divided by the nearest neighbor geometric constraint mean clustering algorithm.
  • the nearest neighbor geometric constraint mean clustering algorithm There are two types of infrared images of underground buildings and infrared images without underground buildings.
  • Target area formation modulation process Generally, underground buildings must maintain a constant temperature and humidity state due to various requirements and their state generally does not change. The thermal radiation existing in the ground itself is modulated by the formation, which leads to the detection area. Changes in the state of matter and energy transfer result in anomalies in the surface temperature distribution and produce unique infrared information fields that are different from other detection areas. Its modulation form appears as a form of modulation of Gaussian thermal diffusion.
  • the standard deviation parameter ⁇ of the Gaussian thermal diffusion function may have the following relationship with the buried depth h of the underground building and the thermal conductivity ⁇ of the formation material:
  • the method uses the surface thermal image detection to identify the underground building structure. It shows that the inverse inverse transformation is performed on the basis of the forward modulation and Gaussian thermal diffusion of the energy through the formation, because the infrared image obtained by the shooting is the underground building. After the formation of the stable infrared image, we can select a stable infrared image with better effect for processing, which can get better after processing.
  • the target image reveals the structural information of the underground building obscured by the stratum.
  • thermophysical equation Two 10cm ⁇ 10cm heating elements were buried in 10cm and 20cm thick sand and gravel and 1.5cm thick marble plates. The shooting distance was 4.5m, and the temperature of the heating device was 27°C. The temperature was 22 ° C and the surface was covered with a layer of grass. A medium-wave infrared camera with a medium-wave infrared window of 3.4 to 5 ⁇ m was used for shooting.
  • the target obtained by the test corresponds to the surface temperature difference with the environment and the underground building with 10m*10m, the target depth is 16m (concrete 2m, soil 14m) and buried depth 32m (concrete 2m, soil 30m)
  • the upper surface temperature is equal to the ambient temperature difference.
  • the average gradation difference ⁇ G of the target area is smaller than the average gradation difference ⁇ G of the target area and the non-target area, and the average gradation difference ⁇ G of the adjacent target area is also larger than the average of the target area and the non-target area.
  • the gradation difference ⁇ G is small, and the average gradation value of the target region is also larger than the average gradation value of the non-target region, so that a region having a large average gradation value and a small adjacent gradation difference value ⁇ G can be used as a region An infrared image g(i, j) that requires inverse modulation processing.
  • the infrared image model formed by the underground building after being modulated by the ground layer can generally be expressed as:
  • g(i,j) ⁇ h(i,j; ⁇ , ⁇ )f(i- ⁇ ,j- ⁇ )d ⁇ d ⁇ +n(i,j) (1)
  • g(i,j) is the infrared image at a certain time
  • f(i,j) is the target original image
  • n(i,j) is the sensor noise
  • h(i,j; ⁇ , ⁇ ) is Gaussian thermal diffusion function
  • the influence of formation modulation on the imaging of underground buildings can usually be assumed to be linear shift invariant, that is, the fuzzy operator (thermal diffusion function) has spatial shift invariance, which can be described as:
  • the target strength as a non-negative function: ⁇ f(x), x ⁇ X ⁇ , where X is the support domain for the target strength.
  • the thermal diffusion function that defines the influence of the formation on the target image is: ⁇ h(y
  • i(y) the intensity of the infrared image pair at coordinate y
  • the image data actually observed at a certain pixel position y of the infrared image is g(y).
  • the infrared image after the formation modulation has a Poisson distribution property. It can be seen that given the target intensity f and the thermal diffusion function h, it can be assumed that g(y) is an independent random variable obeying the Poisson distribution with i(y;f,h) as the mean.
  • the probability of taking the integer gray value g(y) at the element position y can be expressed as
  • equation (9) can be derived for each component f(x) and h(x) and its derivative is equal to zero, which can be derived.
  • equation (13) For the two-dimensional image in the present invention, the equation (13) should be:
  • n is the current number of iterations
  • (x, y) belongs to the thermal diffusion function support domain
  • (i, j) belongs to the image support domain
  • h(-x, -y) represents the conjugate of h(x, y)
  • f( -i, -j) denotes the conjugate of f(i,j).
  • the target image is within the support domain of the observed image. Since the original image of the target is normalized, the sum of the energy values of the target original image is 1, so there is
  • equation (17) For the two-dimensional image in the present invention, equation (17) should be:
  • n is the current number of iterations
  • (x, y) belongs to the thermal diffusion function support domain
  • (i, j) belongs to the image support domain
  • h(-x, -y) represents the conjugate of h(x, y)
  • f( -i, -j) denotes the conjugate of f(i,j).
  • the Gaussian thermal diffusion function h(x, y) and the processed infrared image f(x, y) can be obtained.
  • the infrared image of the underground building can be collected by the anti-modulation process to collect the heat energy after the Gaussian heat diffusion.
  • the infrared image after the inverse modulation process is the infrared image after collecting the energy, and the infrared image can reflect the heat distribution of the underground building.
  • the underground building structure information is restored as realistically as possible, so that people can accurately locate the underground building.
  • the invention is characterized in that: firstly, the method uses Gaussian thermal diffusion as a form of underground building modulation by the ground layer; secondly, the method finds that the infrared image of the underground building and its surrounding environment is carried out by the scaled model test. After the anti-modulation process, the location of the underground building can be clearly seen, and the intrinsic information of the underground building can be exposed, so that the underground building can be accurately detected. Third, the method is inversely modulated by the simulation result map. Not only does the infrared image of the original underground building show clearer, but also the invisible underground tubular facilities are also displayed, and the real structure of the underground building can be inverted. According to the infrared information known in the obtained infrared image To infer the key physical characteristics of underground buildings. Achieve "seeing" and further quantitatively measure the key physical characteristics of underground buildings.
  • the method of the present invention will be specifically described below by way of examples.
  • the flow of the present invention is shown in FIG. 1 , and the specific implementation method includes the following steps. Including: infrared characteristic analysis and modeling of underground buildings in the planar environment, stratigraphic modulation processing of the target area, anti-modulation processing method of the infrared image of the target area, restoration and positioning of the underground building.
  • infrared characteristic analysis and modeling of underground buildings in the planar environment stratigraphic modulation processing of the target area
  • anti-modulation processing method of the infrared image of the target area restoration and positioning of the underground building.
  • the invention firstly analyzes the infrared characteristics of the underground building and its relationship with the environment, in order to display the infrared information of the ground surface modulated by the underground building, and uses the ANSYS to the depth of the underground building in the plane ground environment at a depth of 100 meters.
  • the heat distribution was simulated.
  • the underground building thermal radiation simulation is realized according to the thermal radiation analysis in the "ANSYS 12.0 Thermal Analysis Engineering Application Practical Manual” edited by Zhang Zhaohui.
  • the specific modeling steps are as follows:
  • the obtained infrared image contains the infrared information of the target area and the infrared information of the non-target area.
  • the nearest infrared image is roughly classified into two types: the infrared image of the underground building and the infrared image without the underground building.
  • Figure 7(a) shows the thermal halo image formed by the heat treatment of the heat source through the formation. From the thermal infrared image, it is not easy to accurately locate the exact location of the heat source, from the shape of the heat source and related knowledge. It can be considered that its modulation form is a form of Gaussian thermal diffusion.
  • the heat conduction of underground buildings is spatial heat conduction. Since only the thermal diffusion information of underground buildings to the earth's surface can be detected, it can be equivalent to a two-dimensional Gaussian function model. The model is shown in Figure 3.
  • the infrared thermal image of the surface of the underground building is formed by Gaussian thermal diffusion modulation.
  • the general location of the underground building can be known, but the specific structure of the underground building cannot be determined.
  • the modulation process reveals the structural information of the underground building that is covered by the ground. This is verified by the sand squeezing ratio test.
  • the specific implementation process and results are as follows:
  • the average gray-scale difference ⁇ G G i+1 -G i of the adjacent infrared image regions, respectively, because the average gray-scale difference ⁇ G of the adjacent non-target regions is higher than the average gray of the target region and the non-target region.
  • the difference value ⁇ G is small, and the average gradation difference ⁇ G of the adjacent target area is also smaller than the average gradation difference ⁇ G of the target area and the non-target area, and the average gray value of the target area is also smaller than the non-target area. Since the average gradation value is large, a region where the average gradation value is large and the adjacent average gradation value difference ⁇ G is small can be used as the infrared image g(i, j) requiring the inverse modulation processing.
  • the infrared image model formed by the underground building after being modulated by the ground layer can generally be expressed as:
  • g(i,j) ⁇ h(i,j; ⁇ , ⁇ )f(i- ⁇ ,j- ⁇ )d ⁇ d ⁇ +n(i,j) (1)
  • g(i,j) is the infrared image at a certain time
  • f(i,j) is the target original image
  • n(i,j) is the sensor noise
  • h(i,j; ⁇ , ⁇ ) is Gaussian thermal diffusion function
  • the influence of formation modulation on the imaging of underground buildings can usually be assumed to be linear shift invariant, that is, the fuzzy operator (thermal diffusion function) has spatial shift invariance, which can be described as:
  • g(y) is an independent random variable obeying the Poisson distribution with i(y;f,h) as the mean, and therefore, at the pixel position y
  • the probability of taking the integer gray value g(y) can be expressed as
  • equation (9) can be derived for each component f(x) and h(x) and its derivative is equal to zero, which can be derived.
  • the target image is within the support domain of the observed image. Since the original image of the target is normalized, the sum of the energy values of the target original image is 1, so there is
  • the maximum likelihood estimation image recovery flowchart is shown in FIG. 5.
  • the Gaussian thermal diffusion function and the processed infrared target image can be obtained, which can be clearly seen from the processed infrared target image.
  • Figure 6 is a comparison of the effects of the computer simulation image before and after the anti-modulation process, and the actual underground structure and the inverse modulated infrared image are marked. It can be seen from Fig. 6 that the result map obtained by the anti-modulation process can more clearly detect and locate the target, improve the accuracy of the manual interpretation, and can invert the structure of the underground building.
  • Fig. 7(b) is a graph showing the result obtained after the inverse modulation processing of Fig. 7(a).
  • the infrared image of the underground building can be collected by the anti-modulation process to collect the heat energy after the Gaussian heat diffusion.
  • the infrared image after the inverse modulation process the infrared image after collecting the energy can be seen, and the infrared image can reflect the heat of the underground building.
  • the distribution situation so as to restore the underground building structure information as realistic as possible, so that people can accurately locate the underground building.
  • Fig. 9(b) is a view showing the result obtained after the inverse modulation processing of Fig. 9(a). Comparing Figure 9(a) with Figure 9(b), The position of the heat source is not visible at all from Fig. 9(a), and the position of the heat source can be clearly seen from Fig. 9(b). Realized the positioning of underground buildings. At the same time, it can be seen from Fig. 6 that the result map obtained by the inverse modulation process can more clearly detect and locate the target, improve the accuracy of the manual interpretation, and can invert the structure of the underground building.

Landscapes

  • Engineering & Computer Science (AREA)
  • Physics & Mathematics (AREA)
  • General Physics & Mathematics (AREA)
  • Theoretical Computer Science (AREA)
  • Computer Vision & Pattern Recognition (AREA)
  • Spectroscopy & Molecular Physics (AREA)
  • Geophysics And Detection Of Objects (AREA)
  • Image Processing (AREA)
  • Investigating Or Analyzing Materials Using Thermal Means (AREA)
  • Radiation Pyrometers (AREA)

Abstract

一种平面地表环境中地下建筑的红外成像探测定位方法，包括：获取地下建筑经地层调制后所形成的原始红外图像g0，并确定地下建筑在原始红外图像g0中的大体位置的局部红外图像g；设置迭代终止条件，并设定高斯热扩散函数的初始值h0；以所述局部红外图像g作为初始目标图像f0，根据所述高斯热扩散函数的初始值h0，利用最大似然估计算法迭代求解热扩展函数hn和目标图像fn；判断是否满足迭代终止条件，如果满足，则本次迭代求解得到的目标图像fn即为最终的目标图像f；若不满足，则继续迭代计算。该方法通过对地下建筑经地层调制后所形成的红外图像进行反调制处理，不仅使原来地下建筑的红外图像显示更清晰，还可以反演地下建筑的真实结构。

Description

平面地表环境中地下建筑的红外成像探测定位方法 [技术领域]
本发明属于地球物理学学科理论与遥感技术交叉的领域，更具体地，涉及一种平面地表环境中地下建筑的红外成像探测定位方法。
[背景技术]
目前，随着城市化的进程不断加快，越来越多的人口涌入城市，城市用地显得越来越紧张，这使得更多的建筑不得不选择修建在地下如：地下停车场，地下仓库，地下餐厅等；当然，也有很多设施是出于安全因素的考虑修建在地下，如：大型地下油库，军事设施等。这些都是地下建筑的典型例子，那么一旦这些地下建筑出现故障，就面临着勘察故障出现位置难度大的问题。当然不只是城市中的地下建筑存在勘察难的问题，考古、探矿、工程热物理、水坝探测等方面也面临着这样的问题。由此看来，对地下建筑的探测在民用上有着重要的意义。因此，有必要开展地下建筑探测识别的研究。
目前国内外主要是基于一种被动红外成像探测技术，其物理基础是大量的太阳光能量照射土壤被吸收产生热量，这些被加热的土壤发出红外辐射被热红外传感器探测。自然太阳能经过每日循环的加热和冷却对埋藏的物体和包围其周围的土地的影响是不同的，从而导致可探测的温差。地下建筑的存在可以产生地表水/热分布异常，可用于探测、发现、验证、确认地下建筑。
国内外使用红外成像技术探测地下建筑，但其均直接使用红外成像传感器进行成像，并对其得到的红外图像进行人工判读，这样就存在很大的局限性。其局限性在于：第一，地下建筑的热异常通过其埋入地层的传导调制，到达地表的热分布发生了与地下建筑有很大差异的变化，表现为热 扩散严重，温差大大降低，热信号微弱。第二，通过地层的传导调制，地表热信号所反映的地下建筑的位置可能发生变化，难于发现与定位。第三，人工判读困难，不利于地下建筑的准确探测与定位。
发明内容
针对现有技术的以上缺陷或改进需求，本发明提出了一种新的探测地下建筑的方法，通过对存在地下建筑的红外图像进行反调制处理，解决了地下建筑在红外图像中经地层传导调制后信号变弱，人工无法探测和定位的问题，从而对平面环境中地下建筑进行探测和定位。
本发明提供了一种平面地表环境中地下建筑的红外成像探测定位方法，所述方法根据地下建筑的能量扩散的高斯模型，对地下建筑经地层调制后所形成的原始红外图像进行反调制处理，得到地下建筑的目标图像，所述方法包括以下步骤：包括：
(1)获取地下建筑经地层调制后所形成的原始红外图像g0，并确定地下建筑在原始红外图像g0中的大体位置的局部红外图像g；
(2)设置迭代终止条件，并设定高斯热扩散函数的初始值h0；
(3)以所述局部红外图像g作为初始目标图像f0，根据所述高斯热扩散函数的初始值h0，利用最大似然估计算法迭代求解热扩展函数hn和目标图像fn；
(4)判断是否满足迭代终止条件，如果满足，则本次迭代求解得到的目标图像fn即为最终的目标图像f；若不满足，则返回步骤(3)，继续迭代计算。
优选地，所述步骤(1)具体为：
(1.1)将原始红外图像g0分成N个m×m像素大小的红外图像区域；
(1.2)计算得到N个红外图像区域的平均灰度值，其中第i幅红外图像区域的平均灰度值为
Figure PCTCN2014085713-appb-000001
i＝1□N，gij为第i幅红外图像区域中 的第j个像素的灰度值；
(1.3)分别计算出相邻红外图像区域的平均灰度差值□G＝Gi+1-Gi；
(1.4)将平均灰度值大且相邻平均灰度值差值□G小的区域作为需要反调制处理的局部红外图像g。
优选地，所述步骤(3)中迭代求解热扩展函数hn和目标图像fn具体根据下面两式迭代计算：
Figure PCTCN2014085713-appb-000002
Figure PCTCN2014085713-appb-000003
其中n为当前迭代次数,(x,y)属于热扩散函数支持域，(i,j)属于图像支持域，f(i,j)表示目标图像,g(i,j)表示局部红外图像,h(x,y)是热扩散函数，h(-x,-y)表示h(x,y)的共轭，f(-i,-j)表示f(i,j)的共轭，*为卷积运算符。
优选地，所述步骤(2)中的迭代终止条件为迭代终止次数n>N0或误差ε，所述步骤(4)中判断是否满足迭代条件具体为：
判断是否满足|g-hn+1*fn+1|<ε或者n>N0，如果二者中任一个满足，则满足迭代终止条件，否则不满足。
本发明方法通过对地下建筑经地层调制后所形成的红外图像进行反调制处理后，不仅使原来地下建筑的红外图像显示更清晰，还可以反演地下建筑的真实结构。根据所获得的红外图像中已知的红外信息可以推断出地下建筑关键物理特征，达到“看见”并进一步定量测量出地下建筑关键物理特征的目的。
附图说明
图1为本发明流程示意图；
图2为地下建筑(深度100m)的地表热分布及分类结果划分；
图3为二维高斯函数形式示意图；
图4(a)为缩比模型在19:35时刻的红外热成像；
图4(b)为缩比模型在20:00时刻的红外热成像；
图4(c)为缩比模型在20:20时刻的红外热成像；
图4(d)为缩比模型在20:45时刻的红外热成像；
图4(e)为缩比模型在21:00时刻的红外热成像；
图4(f)为缩比模型在21:30时刻的红外热成像；
图5为最大似然估计图象恢复流程图；
图6为典型地下结构在平地，深度100m环境下，地下建筑地表热分布(右上)及其逆变换结果(右下)；
图7(a)为沙埋实验一得到的原始红外图像；
图7(b)为沙埋实验一反调制处理后的结果；
图8(a)为等效16米的地下建筑估计的热扩散函数，热扩散函数标准差为σ＝6.1；
图8(b)为等效32米的地下建筑估计的热扩散函数，热扩散函数标准差为σ＝9.3；
图9(a)为沙埋实验二原始红外图像；
图9(b)为沙埋实验二反调制处理后的结果。
[具体实施方式]
为了使本发明的目的、技术方案及优点更加清楚明白，以下结合附图及实施例，对本发明进行进一步详细说明。应当理解，此处所描述的具体实施例仅仅用以解释本发明，并不用于限定本发明。此外，下面所描述的本发明各个实施方式中所涉及到的技术特征只要彼此之间未构成冲突就可以相互组合。
本发明针对平面环境中存在地下建筑和不存在地下建筑的红外图像进行初步分类，然后对存在地下建筑的红外图像进行反调制处理，把探测问题归结为数学物理中的求逆问题。将地下建筑的结构信息和位置探测出来， 实现对地下建筑目标的探测和定位，在现有的国内外文献中还没有看到与本发明相同的报道。
本发明所提供的一种平面地表环境中地下建筑的探测定位方法，包括下述几个步骤：
(1)平面环境中地下建筑红外特性分析与建模：以地球物理学、热学、遥感技术等学科的基本理论为指导，分析地下建筑与埋藏环境的物质/能量特征关系及相互作用规律，研究并建立目标区和背景区的红外信息映射模型，并对红外仿真的热图像进行分类学习。包括以下几个子步骤：
(1.1)本发明针对地下建筑红外特征及其与环境的关系先期开展了理论分析。利用现有的软件对典型的地下建筑进行建模仿真。为了将地下建筑经地层调制后的地表红外信息显示出来，利用ANSYS对地下建筑在平面地表环境下深度为100米时的地下建筑的热分布进行了仿真计算。地下建筑热辐射仿真是按照张朝晖编著的《ANSYS 12.0热分析工程应用实战手册》中的热辐射分析实现。
(1.2)根据上述地下建筑的地表红外图像仿真结果进行聚类分析，所得到的红外图像中包含目标区红外信息和非目标区红外信息，利用近邻几何约束均值聚类算法把输入红外图像分为下有地下建筑红外图像和没有地下建筑的红外图像两类。
(1.2.1)选择红外图像中下有地下建筑的样本区和没有地下建筑的样本区，确定样本区的中心位置，即样本区的类心；
(1.2.2)计算每个样本区的平均灰度值，其中第i个红外图像区域的平均灰度值为：
Figure PCTCN2014085713-appb-000004
(1.2.3)计算出不同红外图像目标区的平均灰度差值□G＝Gi+1-Gi，由于非目标区之间的类心距离和目标区之间的类心距离较小，而非目标区和 目标区的类间距离较大，且目标区的平均灰度值高于非目标去的平均灰度值。
(2)目标区域地层调制处理过程：一般的地下建筑由于各种要求必须保持恒温恒湿的状态且其状态一般不会发生改变，其本身存在的热辐射经过地层的调制，就会导致探测区域的物质、能量迁移状态发生改变，导致地表温度分布产生异常，并产生不同于其它探测区域的独特的红外信息场。其调制形式表现为一种高斯热扩散的调制形式。
我们假设拍摄环境中的地层是相对均匀的，拍摄的是在地下建筑热辐射稳定时的地表红外图像，我们做了试验说明在面积为150m2的单元，恒温25℃的条件下，在不同深度、不同地层介质时，通过能量交换在上覆地表体现出0.5K温差所需时间仿真计算结果如下表所示：
Figure PCTCN2014085713-appb-000005
由上面的表格可以知道，不同介质的地层挡不住目标能量交换。同时，高斯热扩散函数标准差参数σ与地下建筑的埋深h和地层材质导热率ε可能存在有如下关系：
(h/ε)＝kσ+l
(3)目标区域红外图像的反调制处理方法：
本方法利用地表热晕图像探测识别地下建筑结构的试验，表明在对能量经地层正向调制变换高斯热扩散的基础上，进行相应的反变换处理，因为拍摄所得到的红外图像是地下建筑经地层调制后稳定的红外图像，我们可以选取一幅效果较好的稳定红外图像进行处理，处理后可以得到较好的 目标图像，可揭示被地层遮蔽的地下建筑的结构信息。
(3.1)将两个10cm×10cm发热体分别埋于距正面10cm和20cm厚的沙土瓦砾子和1.5cm厚的大理石板中，拍摄距离为4.5m，烘鞋器发热体温度27℃，埋藏环境温度为22℃，且在表面覆盖一层草。利用中波红外窗口为3.4～5μm的中波红外相机进行拍摄。根据热物理学方程计算得知，该试验得到的目标对应地表与环境温差和10m*10m的地下建筑，埋深16m(混凝土2m，泥土14m)和埋深32m(混凝土2m，泥土30m)时目标上地表温度与环境温度差相当。
(3.2)目标区域的最大似然估计法的反卷积处理步骤：
(3.2.1)读入地下建筑经地层调制后所形成的原始红外图像g0(i,j)，然后确定地下建筑在红外图像g0(i,j)中的大体位置的局部红外图像。我们利用分类的方法对该局部目标红外图像进行提取和处理，首先将原始红外图像g0(i,j)中分成N个m×m像素大小的红外图像区域，计算得到N个红外图像区域的平均灰度值，其中第i幅红外图像区域的平均灰度值为：
Figure PCTCN2014085713-appb-000006
i＝1□N，
gij为第i幅红外图像区域中的第j个像素的灰度值；然后分别计算出相邻红外图像区域的平均灰度差值□G＝Gi+1-Gi，由于相邻非目标区域的平均灰度差值□G比目标区域与非目标区域的平均灰度差值□G小，且相邻目标区域的平均灰度差值□G也比目标区域与非目标区域的平均灰度差值□G小，目标区域的平均灰度值也比非目标区域的平均灰度值大，故可以将平均灰度值大且相邻平均灰度值差值□G小的区域作为需要反调制处理的红外图像g(i,j)。
(3.2.2)地下建筑经地层调制后所形成的红外图像模型一般可表示为：
g(i,j)＝∫∫h(i,j；α,β)f(i-α,j-β)dαdβ+n(i,j)        (1)
式中g(i,j)为某一时刻的红外图象，f(i,j)为目标原图象，n(i,j)为传感器噪声，h(i,j；α,β)为高斯热扩散函数。
地层调制对地下建筑成像的影响通常可假定为线性移位不变，即模糊算子(热扩散函数)具有空间移不变性，可描述为：
ht(i,j；α,β)＝ht(α,β) (i,j)∈Ω       (2)
将式(2)代入式(1)，则可得出如下卷积形式
gt(i,j)＝∫∫Dht(α,β)f(i-α,j-β)dαdβ+n(i,j)＝ht(i,j)*f(i,j)+n(i,j)   (3)
由式(3)可知，在忽略噪声的情况下，只有得到高斯热扩散的热扩散函数，才能从经地层调制后的红外图象中恢复出地下建筑的结构。对近邻几何约束均值聚类后下有地下建筑的红外图象开展反调制处理，通过此图像来估计热扩散函数和目标图像。
(3.2.3)由于高斯热扩散函数是未知的，首先设置迭代终止次数N和误差ε，通过误差ε来检查是否满足反调制处理结果，并初始化最初的高斯热扩散模型参数σ；
(3.2.4)为了简化记号，我们对图像使用一维描述。首先定义目标强度为非负函数：{f(x),x∈X}，X为目标强度的支持域。定义地层对目标图像影响的热扩散函数为：{h(y|x),y∈Y}，Y为热扩散函数的支持域。为了减小误差，提高探测准确性，我们可以选取拍摄得到的一张稳定后的地表红外图像，可以更好的恢复地下建筑的结构。热扩散函数可看成是空间不变的，我们定义i(y)为红外图像对在坐标y处的强度，则
Figure PCTCN2014085713-appb-000007
显然，由于h为空间不变的卷积算子，同时式(4)也可表示为
Figure PCTCN2014085713-appb-000008
在大多数实际情况下，i不可能被完美地检测到，它总是被一些噪声所污染，红外图象的某像元位置y处实际观测到的图像数据为g(y)。在许多 情况下，地层调制后的红外图像具有Poisson分布性质。由此可见，在给定目标强度f和热扩散函数h条件下，可以假定g(y)是一个以i(y；f,h)为均值的服从Poisson分布的独立随机变量，因此，在象元位置y处取整数灰度值g(y)的概率可以表达为
Figure PCTCN2014085713-appb-000009
假定观测图象各像元是相互独立的，则其联合概率分布
Figure PCTCN2014085713-appb-000010
对(6)式取对数，得到其对数似然函数为
Figure PCTCN2014085713-appb-000011
假定红外观测图象g(y1y2…yn)在统计上是互相独立的，则泊松联合概率分布的对数似然函数为：
Figure PCTCN2014085713-appb-000012
可以看到式(8)中的最后一项为常数，因此它不影响似然函数的变化，为了简化问题可以将其舍去。由式(4)和(8)，我们可以得到对数似然函数为
Figure PCTCN2014085713-appb-000013
我们感兴趣的是从红外图像数据g1(y1y2…yn)中估计出目标强度f。当热扩散函数h(y|x)已知时，我们可以通过一些普通的图象去模糊处理技术来获得目标强度的估计。然而地下建筑的热扩散函数一般都是未知的，这就 增加了去模糊的难度。
为了极大化对数似然函数，可将式(9)分别对各分量f(x)和h(x)求导并令其导数等于零，可以推导出
Figure PCTCN2014085713-appb-000014
由于高斯热扩散函数离散值之和为1，因此有
Figure PCTCN2014085713-appb-000015
由此，我们可以建立如下迭代关系：
Figure PCTCN2014085713-appb-000016
对于本发明中的二维图像，式(13)应该为：
Figure PCTCN2014085713-appb-000017
其中n为当前迭代次数,(x,y)属于热扩散函数支持域，(i,j)属于图像支持域，h(-x,-y)表示h(x,y)的共轭，f(-i,-j)表示f(i,j)的共轭。
为了便于对h(x)求导，将i用式(5)表示，代入式(9)，此时，对数似然函数可等效地表示为
Figure PCTCN2014085713-appb-000018
对h(x)求导并令其为0，得
Figure PCTCN2014085713-appb-000019
目标图象在观察图象的支持域之内，由于对目标原图像进行归一化处理，因此，目标原图像能量值之和为1，于是有
Figure PCTCN2014085713-appb-000020
同理，当目标图象fn+1(x)用式(13)估计出来后，我们可以建立如下求解新高斯热扩散函数hn+1(x)的迭代关系：
Figure PCTCN2014085713-appb-000021
对于本发明中的二维图像，式(17)应该为：
Figure PCTCN2014085713-appb-000022
其中n为当前迭代次数,(x,y)属于热扩散函数支持域，(i,j)属于图像支持域，h(-x,-y)表示h(x,y)的共轭，f(-i,-j)表示f(i,j)的共轭。
因此，我们可以通过最大似然估计法来恢复地下建筑的真实红外图像。
(3.2.5)利用最大似然估计的方法反复迭代红外图像g(x,y)，可以得到高斯热扩散函数h(x,y)和处理后的红外图像f(x,y)。
(4)地下建筑的恢复和定位：
存在地下建筑的红外图像经过反调制处理可以将高斯热扩散后的热能量收集起来，反调制处理后的红外图像就是收集能量之后的红外图像，该红外图像可以反映地下建筑的热分布情况，这样就尽可能真实的恢复了地下建筑结构信息，使人准确的对地下建筑进行定位。
本发明的特点在于：第一，本方法将高斯热扩散作为地下建筑经地层调制的形式；第二，本方法通过所做的缩比模型试验发现，对地下建筑和其周围环境的红外图像进行反调制处理后，可以很明显的看到地下建筑所处的位置，揭露地下建筑的本征信息，这样就可以准确的探测地下建筑；第三，本方法通过对仿真结果图反调制处理后，不仅使原来地下建筑的红外图像显示更清晰，而且将不可见的地下管状设施也表现出来了，同时可以反演地下建筑的真实结构。根据所获得的红外图像中已知的红外信息可 以推断出地下建筑关键物理特征。达到“看见”并进一步定量测量出地下建筑关键物理特征的目的。
以下以实例具体地对本发明方法进行说明，本发明流程如图1所示，具体实施方法包括以下步骤。包括：平面环境中地下建筑红外特性分析与建模、目标区域地层调制处理过程、目标区域红外图像的反调制处理方法、地下建筑的恢复和定位。以下结合附图和具体实例对本发明进一步说明。
(1)平面环境中地下建筑红外特性分析与建模步骤：
本发明针对地下建筑红外特征及其与环境的关系先期开展了理论分析，为了将地下建筑经地层调制后的地表红外信息显示出来，利用ANSYS对地下建筑在平面地表环境下深度为100米时的热分布进行了仿真计算。地下建筑热辐射仿真是按照张朝晖编著的《ANSYS 12.0热分析工程应用实战手册》中的热辐射分析实现的，具体的建模的步骤如下：
(1.1)选择表面效应单元SURF151，其利用实体表面的节点行程单元，并直接覆盖在实体单元的表面；其中，SURF151单元有主要实常数(角系数、Stefan-Boltzmann常数)、材料属性(密度、热辐射)、表面载荷(对流、热流密度)和体载荷(生热率)；
(1.2)设置表面效应单元对应的材料属性；设置密度为7800，比热为465，辐射率为1；
(1.3)设置表面效应单元的实常数；Stefan-Boltzmann常数为5.67×10-8，角系数的形状参数设置为1；
(1.4)创建几何模型、划分网格；Keypoint number设置为默认值，Global Element Sizes为3，即可生成有限元模型；
(1.5)利用Anasys加载求解；
(1.6)查看求解结果，并利用Anasys画出仿真图像，模型如图2所示。
(1.2)根据上述地下建筑的地表红外图像仿真结果进行聚类分析，所 得到的红外图像中包含目标区红外信息和非目标区红外信息，利用近邻几何约束均值聚类算法把输入红外图像粗分为下有地下建筑红外图像和没有地下建筑的红外图像两类。
(1.2.1)选择红外图像中下有地下建筑的样本区和没有地下建筑的样本区，确定样本区的中心位置，即样本区的类心；本实例中，其中数字1代表不存在地下建筑的非目标区域，数字2代表存在地下建筑的非目标区。区域示意图如图3。
(1.2.2)计算每个样本区的平均灰度值，其中第i个红外图像区域的平均灰度值为：
Figure PCTCN2014085713-appb-000023
由于是计算机仿真结果，可以直接读出各个点的温度值，然后计算平均温度：
Figure PCTCN2014085713-appb-000024
(1.2.3)计算出不同红外图像目标区的平均灰度差值□G＝Gi+1-Gi，由于非目标区之间的类心距离和目标区之间的类心距离较小，而非目标区和目标区的类间距离较大，且目标区的平均灰度值高于非目标区的平均灰度值。本实例中，区域2类内之间的平均温差为0.1K，区域1类内之间的平均温差为0K，区域1和区域2类间平均温差为0.3K。
(2)目标区域地层调制处理过程：
一般的地下建筑由于各种要求需要保持恒温恒湿的状态且其状态一般不会发生改变，其本身存在的热辐射经过地层的调制，就会导致探测区域的物质、能量迁移状态发生改变，导致地表温度分布产生异常，并产生不同于其它探测区域的独特的红外信息场。其调制形式表现为一种高斯热扩散的调制模型。我们通过试验来验证这一过程，试验具体步骤如下：
(2.1)将10cm×10cm发热体埋于距正面10cm厚沙土和1.5cm的大理石板中，拍摄距离4.5m，发热体温度27℃，埋藏环境的温度为22℃。通过对红外成像实验数据的时间变化情况分析发现，背景沙土区域的灰度变化缓慢，而目标所对应的沙土区域的灰度随着目标温度向表面的扩散而呈现较大的改变，图4(a)(b)(c)(d)(e)(f)是在不同时刻拍摄得到的红外图像。热辐射平衡之后拍摄得到的原始红外图像如图7(a)所示。
(2.2)图7(a)显示的是发热源经过地层的调制处理后形成的热晕图像，从该热晕红外图像中不容易准确定位发热源的准确位置，从热源的形状和相关的知识可以认为其调制形式为高斯热扩散的形式。地下建筑的热传导是空间上的热传导，由于只能探测到地下建筑向地表的热扩散信息，所以可以等效为二维高斯函数模型，模型如图3所示。
(3)目标区域红外图像的反调制处理方法：
地下建筑经地层高斯热扩散调制形成的地表热晕红外图像，经过步骤(1)可以知道地下建筑的大体位置，但是地下建筑的具体结构并无法确定，为了得到地下建筑的具体结构信息，开展反调制处理，可揭示被地层遮蔽的地下建筑的结构信息。通过沙埋缩比试验对此进行验证，具体实施过程和结果如下：
(3.1)我们所做的沙埋缩比试验是模拟实场拍摄。拍摄的器材和背景环境如下表：
表1.1拍摄试验器材及相关数据
Figure PCTCN2014085713-appb-000025
Figure PCTCN2014085713-appb-000026
对试验精心设计后，我们通过多次试验，建立比照与参考以使得数据更加的准确和具有可比性。同时我们对图像进行了定标测温，得到整幅图像场景的温度。
(3.2)将两个10cm×10cm发热体分别埋于距正面10cm和20cm厚的沙土和1.5cm厚的大理石板中，拍摄距离为4.5m，发热体温度27℃，埋藏环境温度为22℃，且在表面覆盖一层草。拍摄得到的原始红外图像如图9(a)所示。根据热物理学方程计算得知，该试验得到的目标对应地表与环境温差和10m*10m的地下建筑，埋深16m(混凝土2m，泥土14m)和埋深32m(混凝土2m，泥土30m)时目标上地表温度与环境温度差相当。
(3.3)目标区域的最大似然估计法的反卷积处理步骤：
(3.3.1)读入地下建筑经地层调制后所形成的原始红外图像g0(i,j)，然后确定地下建筑在红外图像g0(i,j)中的大体位置的局部红外图像。我们利用分类的方法对该局部目标红外图像进行提取和处理，首先将原始红外图像g0(i,j)中分成N个m*m像素大小的红外图像区域，计算得到N个红外图像区域的平均灰度值，其中第i幅红外图像区域的平均灰度值为：
Figure PCTCN2014085713-appb-000027
然后分别计算出相邻红外图像区域的平均灰度差值□G＝Gi+1-Gi，由于相邻非目标区域的平均灰度差值□G比目标区域与非目标区域的平均灰度差值□G小，且相邻目标区域的平均灰度差值□G也比目标区域与非目标区域的平均灰度差值□G小，目标区域的平均灰度值也比非目标区域的平均灰度值 大，故可以将平均灰度值大且相邻平均灰度值差值□G小的区域作为需要反调制处理的红外图像g(i,j)。
(3.3.2)地下建筑经地层调制后所形成的红外图像模型一般可表示为：
g(i,j)＝∫∫h(i,j；α,β)f(i-α,j-β)dαdβ+n(i,j)     (1)
式中g(i,j)为某一时刻的红外图象，f(i,j)为目标原图象，n(i,j)为传感器噪声，h(i,j；α,β)为高斯热扩散函数。
地层调制对地下建筑成像的影响通常可假定为线性移位不变，即模糊算子(热扩散函数)具有空间移不变性，可描述为：
ht(i,j；α,β)＝ht(α,β) (i,j)∈Ω     (2)
将式(2)代入式(1)，则可得出如下卷积形式
gt(i,j)＝∫∫Dht(α,β)f(i-α,j-β)dαdβ+n(i,j)＝ht(i,j)*f(i,j)+n(i,j)     (3)
由式(3)可知，在忽略噪声的情况下，只有得到高斯热扩散的热扩散函数，才能从经地层调制后的红外图象中恢复出地下建筑的结构。对近邻几何约束均值聚类后下有地下建筑的红外图象开展反调制处理，通过此图像来估计热扩散函数和目标图像。原始红外图像如图9(a)。
(3.3.3)设置迭代终止次数N和误差ε，通过误差ε来检查是否满足反调制处理结果，并设置最初的高斯热扩散模型参数σ；本实例中，N＝40，ε＝0.003，σ＝4.0。
(3.3.4)定义目标强度为非负函数：{f(x),x∈X}，其中X为目标强度的支持域，定义地层对目标图像影响的热扩散函数为：{h(y|x),y∈Y}，Y为观察图象的支持域。为了减小误差，提高探测准确性，我们可以选取拍摄得到的一张稳定后的地表红外图像，可以更好的恢复地下建筑的结构。热扩散函数可看成是空间不变的，我们定义i(y)为红外图像对在坐标y处的强度，则
Figure PCTCN2014085713-appb-000028
显然，由于h为空间不变的卷积算子，同时式(4)也可表示为
Figure PCTCN2014085713-appb-000029
在给定目标强度f和热扩散函数h条件下，可以假定g(y)是一个以i(y；f,h)为均值的服从Poisson分布的独立随机变量，因此，在象元位置y处取整数灰度值g(y)的概率可以表达为
Figure PCTCN2014085713-appb-000030
假定观测图象各像元是相互独立的，则其联合概率分布
Figure PCTCN2014085713-appb-000031
对(6)式取对数，得到其对数似然函数为
Figure PCTCN2014085713-appb-000032
假定红外观测图象g(y1y2…yn)在统计上是互相独立的，则泊松联合概率分布的对数似然函数为：
Figure PCTCN2014085713-appb-000033
可以看到式(8)中的最后一项为常数，因此它不影响似然函数的变化，为了简化问题可以将其舍去。由式(4)和(8)，我们可以得到对数似然函数为
Figure PCTCN2014085713-appb-000034
我们感兴趣的是从红外图像数据g1(y1y2…yn)中估计出目标强度f。当热扩散函数h(y|x)已知时，我们可以通过一些普通的图象去模糊处理技术来获得目标强度的估计。然而地下建筑的热扩散函数一般都是未知的，这就增加了去模糊的难度。
为了极大化对数似然函数，可将式(9)分别对各分量f(x)和h(x)求导并令其导数等于零，可以推导出
Figure PCTCN2014085713-appb-000035
由于高斯热扩散函数离散值之和为1，因此有
Figure PCTCN2014085713-appb-000036
由此，我们可以建立如下迭代关系：
Figure PCTCN2014085713-appb-000037
为了便于对h(x)求导，将i用式(5)表示，代入式(9)，此时，对数似然函数可等效地表示为
Figure PCTCN2014085713-appb-000038
对h(x)求导并令其为0，得
Figure PCTCN2014085713-appb-000039
目标图象在观察图象的支持域之内，由于对目标原图像进行归一化处理，因此，目标原图像能量值之和为1，于是有
Figure PCTCN2014085713-appb-000040
同理，当目标图象fn+1(x)用式(13)估计出来后，我们可以建立如下求解新热扩散函数hn+1(x)的迭代关系：
Figure PCTCN2014085713-appb-000041
上述实施例中，我们利用的是二维高斯函数，由于高斯函数具有可分离性，二维函数的卷积可以用两次一维函数卷积来代替，所以反调制的过程是可以使用一维来描述的，针对具体的二维图像，将公式对应的变换为二维形式即可。
因此，我们可以通过使用最大似然估计法来设法找到红外目标图像和高斯热扩散函数。最大似然估计图像恢复流程图如图5所示。
(3.3.5)利用最大似然估计的方法反复迭代红外图像g(x,y)，可以得到高斯热扩散函数和处理后的红外目标图像，从处理后的红外目标图像中可以清楚地看到地下建筑的具体结构。图6是计算机仿真图像反调制处理前后效果的对比图，同时对实际地下建筑和反调制后的红外图像进行了标记。从图6可以看出反调制处理得到的结果图可以更清晰的探测和定位目标，提高了人工判读的准确性，同时可以将看不到地下建筑的结构反演出来。利用所做的沙埋试验对本发明进行验证，得到图9(a)的等效16米的地下建筑的高斯热扩散函数如图8(a)所示，等效32米的地下建筑的高斯热扩散函数如图8(b)所示，其中σ为获得的高斯热扩散函数标准差。图7(b)是图7(a)经过反调制处理之后所得到的结果图。
(4)地下建筑的恢复和定位：
存在地下建筑的红外图像经过反调制处理可以将高斯热扩散后的热能量收集起来，在反调制处理后的红外图像中可以看到收集能量之后的红外图像，该红外图像可以反映地下建筑的热分布情况，这样就尽可能真实的恢复了地下建筑结构信息，使人准确的对地下建筑进行定位。图9(b)是图9(a)经过反调制处理之后所得到的结果图。对比图9(a)和图9(b)， 从图9(a)中根本无法看出发热源的位置，而从图9(b)中可以清楚地看到发热源所处的位置。实现了地下建筑的定位。同时从图6可以看出反调制处理得到的结果图可以更清晰的探测和定位目标，提高了人工判读的准确性，同时可以将看不到地下建筑的结构反演出来。

Claims (4)

  1. 一种平面地表环境中地下建筑的红外成像探测定位方法，其特征在于，所述方法根据地下建筑的能量扩散的高斯模型，对地下建筑经地层调制后所形成的原始红外图像进行反调制处理，得到地下建筑的目标图像，所述方法包括以下步骤：
    (1)获取地下建筑经地层调制后所形成的原始红外图像g0，并确定地下建筑在原始红外图像g0中的大体位置的局部红外图像g；
    (2)设置迭代终止条件，并设定高斯热扩散函数的初始值h0；
    (3)以所述局部红外图像g作为初始目标图像f0，根据所述高斯热扩散函数的初始值h0，利用最大似然估计算法迭代求解热扩展函数hn和目标图像fn；
    (4)判断是否满足迭代终止条件，如果满足，则本次迭代求解得到的目标图像fn即为最终的目标图像f；若不满足，则返回步骤(3)，继续迭代计算。
  2. 如权利要求1所述的方法，其特征在于，所述步骤(1)具体为：
    (1.1)将原始红外图像g0分成N个m×m像素大小的红外图像区域；
    (1.2)计算得到N个红外图像区域的平均灰度值，其中第i幅红外图像区域的平均灰度值为
    Figure PCTCN2014085713-appb-100001
    i＝1□N，gij为第i幅红外图像区域中的第j个像素的灰度值；
    (1.3)分别计算出相邻红外图像区域的平均灰度差值□G＝Gi+1-Gi；
    (1.4)将平均灰度值大且相邻平均灰度值差值□G小的区域作为需要反调制处理的局部红外图像g。
  3. 如权利要求1或2所述的方法，其特征在于，所述步骤(3)中迭代求解热扩展函数hn和目标图像fn具体根据下面两式迭代计算：
    Figure PCTCN2014085713-appb-100002
    Figure PCTCN2014085713-appb-100003
    其中n为当前迭代次数,(x,y)属于热扩散函数支持域，(i,j)属于图像支持域，f(i,j)表示目标图像,g(i,j)表示局部红外图像,h(x,y)是热扩散函数，h(-x,-y)表示h(x,y)的共轭，f(-i,-j)表示f(i,j)的共轭，*为卷积运算符。
  4. 如权利要求3所述的方法，其特征在于，所述步骤(2)中的迭代终止条件为迭代终止次数n>N0或误差ε，所述步骤(4)中判断是否满足迭代条件具体为：
    判断是否满足|g-hn+1*fn+1|<ε或者n>N0，如果二者中任一个满足，则满足迭代终止条件，否则不满足。
PCT/CN2014/085713 2013-12-31 2014-09-02 平面地表环境中地下建筑的红外成像探测定位方法 WO2015101058A1 (zh)

Priority Applications (1)

Application Number Priority Date Filing Date Title
US15/105,456 US9741120B2 (en) 2013-12-31 2014-09-02 Infrared imaging detection and positioning method for underground building in planar land surface environment

Applications Claiming Priority (2)

Application Number Priority Date Filing Date Title
CN201310754286.8A CN103884431B (zh) 2013-12-31 2013-12-31 平面地表环境中地下建筑的红外成像探测定位方法
CN201310754286.8 2013-12-31

Publications (1)

Publication Number Publication Date
WO2015101058A1 true WO2015101058A1 (zh) 2015-07-09

Family

ID=50953451

Family Applications (1)

Application Number Title Priority Date Filing Date
PCT/CN2014/085713 WO2015101058A1 (zh) 2013-12-31 2014-09-02 平面地表环境中地下建筑的红外成像探测定位方法

Country Status (3)

Country Link
US (1) US9741120B2 (zh)
CN (1) CN103884431B (zh)
WO (1) WO2015101058A1 (zh)

Cited By (1)

* Cited by examiner, † Cited by third party
Publication number Priority date Publication date Assignee Title
CN110852976A (zh) * 2019-11-22 2020-02-28 昆明物理研究所 一种红外图像明暗不均匀校正方法及计算机程序产品

Families Citing this family (11)

* Cited by examiner, † Cited by third party
Publication number Priority date Publication date Assignee Title
CN103884431B (zh) 2013-12-31 2015-09-09 华中科技大学 平面地表环境中地下建筑的红外成像探测定位方法
CN103744124B (zh) 2013-12-31 2015-07-22 华中科技大学 一种平面地形中地下管状设施红外成像探测定位方法
CN104104922A (zh) * 2014-07-24 2014-10-15 成都市晶林科技有限公司 一种考古探测系统及方法
CN104484577B (zh) * 2014-12-30 2017-06-16 华中科技大学 一种基于山脊能量校正的山地中带状地下目标的探测方法
CN105528586B (zh) * 2015-12-25 2017-04-19 华中科技大学 利用双时相探测定位山体中地下建筑的背景滤波方法
CN108053411B (zh) * 2017-12-21 2020-05-19 华中科技大学 一种边界热交换约束下的江底隧道遥感探测定位方法
CN108510453B (zh) * 2018-03-07 2021-06-29 西安工程大学 基于视觉注意机制的智能交通监控图像去模糊方法
US11243329B2 (en) 2019-03-01 2022-02-08 Nancy Kerr Del Grande Detecting subsurface objects and voids using thermal inertia
CN112183516A (zh) * 2020-09-22 2021-01-05 广东中鹏热能科技有限公司 基于图形图像识别技术对单照片多区域最高温的读取方法
CN113034469B (zh) * 2021-03-24 2023-10-24 东北大学 基于transformer的热成像检测管道内部缺陷的方法
CN114112069B (zh) * 2022-01-27 2022-04-26 华中科技大学 地质约束的城市深埋条带通道红外成像探测方法及系统

Citations (6)

* Cited by examiner, † Cited by third party
Publication number Priority date Publication date Assignee Title
JP2005283383A (ja) * 2004-03-30 2005-10-13 Mitsubishi Electric Corp 埋設物探査方法および埋設物探査装置
US20100025582A1 (en) * 2006-12-19 2010-02-04 Weil Gary J Remote sensing of subsurface artifacts by use of visual and thermal imagery
CN101699509A (zh) * 2009-11-11 2010-04-28 耿则勋 一种利用气象数据进行大气模糊遥感影像恢复的方法
CN103040889A (zh) * 2013-01-25 2013-04-17 四川科伦新光医药有限公司 一种人参刺五加口服液及其生产工艺
CN103279935A (zh) * 2013-06-09 2013-09-04 河海大学 基于map算法的热红外遥感图像超分辨率重建方法及系统
CN103884431A (zh) * 2013-12-31 2014-06-25 华中科技大学 平面地表环境中地下建筑的红外成像探测定位方法

Family Cites Families (4)

* Cited by examiner, † Cited by third party
Publication number Priority date Publication date Assignee Title
JP4997000B2 (ja) * 2007-07-05 2012-08-08 株式会社日立製作所 埋設配管経路探索システム、及び埋設配管経路探索方法
US8494220B2 (en) * 2010-10-15 2013-07-23 Nancy Kerr Del Grande Temporal thermal imaging method for detecting subsurface objects and voids
CN102776904A (zh) * 2012-07-19 2012-11-14 同济大学 一种利用红外热场快速预测地下连续墙侧壁渗漏的方法
CN103049889A (zh) 2012-12-06 2013-04-17 华中科技大学 一种强噪声气动光学效应退化图像的预处理方法

Patent Citations (6)

* Cited by examiner, † Cited by third party
Publication number Priority date Publication date Assignee Title
JP2005283383A (ja) * 2004-03-30 2005-10-13 Mitsubishi Electric Corp 埋設物探査方法および埋設物探査装置
US20100025582A1 (en) * 2006-12-19 2010-02-04 Weil Gary J Remote sensing of subsurface artifacts by use of visual and thermal imagery
CN101699509A (zh) * 2009-11-11 2010-04-28 耿则勋 一种利用气象数据进行大气模糊遥感影像恢复的方法
CN103040889A (zh) * 2013-01-25 2013-04-17 四川科伦新光医药有限公司 一种人参刺五加口服液及其生产工艺
CN103279935A (zh) * 2013-06-09 2013-09-04 河海大学 基于map算法的热红外遥感图像超分辨率重建方法及系统
CN103884431A (zh) * 2013-12-31 2014-06-25 华中科技大学 平面地表环境中地下建筑的红外成像探测定位方法

Non-Patent Citations (1)

* Cited by examiner, † Cited by third party
Title
LING, JIANGUO ET AL.: "INFRARED TARGET EXTRACTION METHOD BASED ON KERNEL DENSITY ESTIMATION", J. INFRARED MILLIM WAVE, vol. 25, no. 6, 31 December 2006 (2006-12-31), pages 434 - 438 *

Cited By (1)

* Cited by examiner, † Cited by third party
Publication number Priority date Publication date Assignee Title
CN110852976A (zh) * 2019-11-22 2020-02-28 昆明物理研究所 一种红外图像明暗不均匀校正方法及计算机程序产品

Also Published As

Publication number Publication date
US20160321795A1 (en) 2016-11-03
CN103884431B (zh) 2015-09-09
US9741120B2 (en) 2017-08-22
CN103884431A (zh) 2014-06-25

Similar Documents

Publication Publication Date Title
WO2015101058A1 (zh) 平面地表环境中地下建筑的红外成像探测定位方法
US10365399B2 (en) Method for infrared imaging detection and positioning of underground tubular facility in plane terrain
Quan et al. Multi-temporal trajectory of the urban heat island centroid in Beijing, China based on a Gaussian volume model
Wang et al. Effect of the sampling design of ground control points on the geometric correction of remotely sensed imagery
Lindberg et al. Continuous sky view factor maps from high resolution urban digital elevation models
Wang et al. Recovering land surface temperature under cloudy skies considering the solar‐cloud‐satellite geometry: Application to MODIS and Landsat‐8 data
Wang et al. Novel system for rapid investigation and damage detection in cultural heritage conservation based on deep learning
Jiao et al. Evaluation of four sky view factor algorithms using digital surface and elevation model data
Anders et al. Multitemporal terrestrial laser scanning point clouds for thaw subsidence observation at Arctic permafrost monitoring sites
Fadakar Alghalandis et al. The RANSAC method for generating fracture networks from micro-seismic event data
Zhang et al. Robust algorithm for registration of building point clouds using planar patches
Yu et al. Spatial interpolation-based analysis method targeting visualization of the indoor thermal environment
Statella et al. Image processing algorithm for the identification of Martian dust devil tracks in MOC and HiRISE images
Dogan et al. Surfer: A fast simulation algorithm to predict surface temperatures and mean radiant temperatures in large urban models
CN105426881A (zh) 山体背景热场模型约束的地下热源昼间遥感探测定位方法
Smith et al. Microtopography of bare peat: a conceptual model and objective classification from high‐resolution topographic survey data
Jia et al. Multi-scale remote sensing estimates of urban fractions and road widths for regional models
Papadopoulos et al. Electrical resistivity tomography for the modelling of cultural deposits and geomophological landscapes at Neolithic sites: a case study from Southeastern Hungary
Di Giuseppe et al. Electrical resistivity tomography imaging of the near-surface structure of the Solfatara crater, Campi Flegrei (Naples, Italy)
Su et al. Detect and identify earth rock embankment leakage based on UAV visible and infrared images
Hu et al. What you get is not always what you see—pitfalls in solar array assessment using overhead imagery
Scudero et al. Integrated Geophysical Investigations at the Greek Kamarina Site (Southern Sicily, Italy)
Feng et al. Robotic inspection of underground utilities for construction survey using a ground penetrating radar
Yu et al. Estimating exposure roughness based on Google earth
Danino et al. Quantifying shallow overland flow patterns under laboratory simulations using thermal and LiDAR imagery

Legal Events

Date Code Title Description
121 Ep: the epo has been informed by wipo that ep was designated in this application

Ref document number: 14876209

Country of ref document: EP

Kind code of ref document: A1

WWE Wipo information: entry into national phase

Ref document number: 15105456

Country of ref document: US

NENP Non-entry into the national phase

Ref country code: DE

122 Ep: pct application non-entry in european phase

Ref document number: 14876209

Country of ref document: EP

Kind code of ref document: A1