CN114895346B - 适用于全地形的航空γ能谱测量系统及其测量方法 - Google Patents
适用于全地形的航空γ能谱测量系统及其测量方法 Download PDFInfo
- Publication number
- CN114895346B CN114895346B CN202210481171.5A CN202210481171A CN114895346B CN 114895346 B CN114895346 B CN 114895346B CN 202210481171 A CN202210481171 A CN 202210481171A CN 114895346 B CN114895346 B CN 114895346B
- Authority
- CN
- China
- Prior art keywords
- coordinate
- detector
- aviation
- point
- laser
- 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
Links
- 238000001228 spectrum Methods 0.000 title claims abstract description 101
- 238000000034 method Methods 0.000 title claims abstract description 29
- 238000005259 measurement Methods 0.000 claims abstract description 128
- 238000012876 topography Methods 0.000 claims abstract description 41
- 238000004364 calculation method Methods 0.000 claims description 40
- 230000005251 gamma ray Effects 0.000 claims description 27
- ZSLUVFAKFWKJRC-IGMARMGPSA-N 232Th Chemical compound [232Th] ZSLUVFAKFWKJRC-IGMARMGPSA-N 0.000 claims description 15
- 229910052776 Thorium Inorganic materials 0.000 claims description 15
- 229910052770 Uranium Inorganic materials 0.000 claims description 15
- JFALSRSLKYAFGM-UHFFFAOYSA-N uranium(0) Chemical compound [U] JFALSRSLKYAFGM-UHFFFAOYSA-N 0.000 claims description 15
- 108091006146 Channels Proteins 0.000 claims description 13
- 238000005070 sampling Methods 0.000 claims description 11
- ZLMJMSJWJFRBEC-UHFFFAOYSA-N Potassium Chemical compound [K] ZLMJMSJWJFRBEC-UHFFFAOYSA-N 0.000 claims description 9
- 229910052700 potassium Inorganic materials 0.000 claims description 9
- 239000011591 potassium Substances 0.000 claims description 9
- 102000004257 Potassium Channel Human genes 0.000 claims description 6
- 108020001213 potassium channel Proteins 0.000 claims description 6
- 238000003705 background correction Methods 0.000 claims description 4
- 230000002285 radioactive effect Effects 0.000 claims description 4
- 230000000630 rising effect Effects 0.000 claims description 4
- JJWKPURADFRFRB-UHFFFAOYSA-N carbonyl sulfide Chemical compound O=C=S JJWKPURADFRFRB-UHFFFAOYSA-N 0.000 claims description 3
- 229910052704 radon Inorganic materials 0.000 claims description 3
- SYUHGPGVQRZVTB-UHFFFAOYSA-N radon atom Chemical compound [Rn] SYUHGPGVQRZVTB-UHFFFAOYSA-N 0.000 claims description 3
- 238000007670 refining Methods 0.000 claims description 2
- 230000005856 abnormality Effects 0.000 abstract description 5
- 238000005516 engineering process Methods 0.000 abstract description 3
- 230000000873 masking effect Effects 0.000 abstract description 2
- 230000002159 abnormal effect Effects 0.000 abstract 1
- 239000002245 particle Substances 0.000 description 33
- 230000000875 corresponding effect Effects 0.000 description 26
- 230000005855 radiation Effects 0.000 description 21
- 238000001514 detection method Methods 0.000 description 10
- 238000006243 chemical reaction Methods 0.000 description 6
- 238000002347 injection Methods 0.000 description 5
- 239000007924 injection Substances 0.000 description 5
- 239000000523 sample Substances 0.000 description 4
- 230000003068 static effect Effects 0.000 description 3
- 238000010586 diagram Methods 0.000 description 2
- 238000001730 gamma-ray spectroscopy Methods 0.000 description 2
- 230000001788 irregular Effects 0.000 description 2
- 239000012857 radioactive material Substances 0.000 description 2
- 239000000941 radioactive substance Substances 0.000 description 2
- 230000000284 resting effect Effects 0.000 description 2
- 239000000243 solution Substances 0.000 description 2
- 230000001133 acceleration Effects 0.000 description 1
- 238000009825 accumulation Methods 0.000 description 1
- 238000004458 analytical method Methods 0.000 description 1
- 238000010276 construction Methods 0.000 description 1
- 239000013078 crystal Substances 0.000 description 1
- 238000001914 filtration Methods 0.000 description 1
- 238000000691 measurement method Methods 0.000 description 1
Classifications
-
- G—PHYSICS
- G01—MEASURING; TESTING
- G01T—MEASUREMENT OF NUCLEAR OR X-RADIATION
- G01T1/00—Measuring X-radiation, gamma radiation, corpuscular radiation, or cosmic radiation
- G01T1/36—Measuring spectral distribution of X-rays or of nuclear radiation spectrometry
-
- G—PHYSICS
- G01—MEASURING; TESTING
- G01C—MEASURING DISTANCES, LEVELS OR BEARINGS; SURVEYING; NAVIGATION; GYROSCOPIC INSTRUMENTS; PHOTOGRAMMETRY OR VIDEOGRAMMETRY
- G01C21/00—Navigation; Navigational instruments not provided for in groups G01C1/00 - G01C19/00
- G01C21/10—Navigation; Navigational instruments not provided for in groups G01C1/00 - G01C19/00 by using measurements of speed or acceleration
- G01C21/12—Navigation; Navigational instruments not provided for in groups G01C1/00 - G01C19/00 by using measurements of speed or acceleration executed aboard the object being navigated; Dead reckoning
- G01C21/16—Navigation; Navigational instruments not provided for in groups G01C1/00 - G01C19/00 by using measurements of speed or acceleration executed aboard the object being navigated; Dead reckoning by integrating acceleration or speed, i.e. inertial navigation
- G01C21/165—Navigation; Navigational instruments not provided for in groups G01C1/00 - G01C19/00 by using measurements of speed or acceleration executed aboard the object being navigated; Dead reckoning by integrating acceleration or speed, i.e. inertial navigation combined with non-inertial navigation instruments
- G01C21/1652—Navigation; Navigational instruments not provided for in groups G01C1/00 - G01C19/00 by using measurements of speed or acceleration executed aboard the object being navigated; Dead reckoning by integrating acceleration or speed, i.e. inertial navigation combined with non-inertial navigation instruments with ranging devices, e.g. LIDAR or RADAR
-
- 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
- G01S17/00—Systems using the reflection or reradiation of electromagnetic waves other than radio waves, e.g. lidar systems
- G01S17/88—Lidar systems specially adapted for specific applications
-
- 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
- G01S19/00—Satellite radio beacon positioning systems; Determining position, velocity or attitude using signals transmitted by such systems
- G01S19/38—Determining a navigation solution using signals transmitted by a satellite radio beacon positioning system
- G01S19/39—Determining a navigation solution using signals transmitted by a satellite radio beacon positioning system the satellite radio beacon positioning system transmitting time-stamped messages, e.g. GPS [Global Positioning System], GLONASS [Global Orbiting Navigation Satellite System] or GALILEO
- G01S19/42—Determining position
- G01S19/48—Determining position by combining or switching between position solutions derived from the satellite radio beacon positioning system and position solutions derived from a further system
-
- 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
- G01S19/00—Satellite radio beacon positioning systems; Determining position, velocity or attitude using signals transmitted by such systems
- G01S19/38—Determining a navigation solution using signals transmitted by a satellite radio beacon positioning system
- G01S19/39—Determining a navigation solution using signals transmitted by a satellite radio beacon positioning system the satellite radio beacon positioning system transmitting time-stamped messages, e.g. GPS [Global Positioning System], GLONASS [Global Orbiting Navigation Satellite System] or GALILEO
- G01S19/42—Determining position
- G01S19/48—Determining position by combining or switching between position solutions derived from the satellite radio beacon positioning system and position solutions derived from a further system
- G01S19/49—Determining position by combining or switching between position solutions derived from the satellite radio beacon positioning system and position solutions derived from a further system whereby the further system is an inertial position system, e.g. loosely-coupled
Landscapes
- Engineering & Computer Science (AREA)
- Remote Sensing (AREA)
- Radar, Positioning & Navigation (AREA)
- Physics & Mathematics (AREA)
- General Physics & Mathematics (AREA)
- Computer Networks & Wireless Communication (AREA)
- Spectroscopy & Molecular Physics (AREA)
- Life Sciences & Earth Sciences (AREA)
- Molecular Biology (AREA)
- High Energy & Nuclear Physics (AREA)
- Automation & Control Theory (AREA)
- Health & Medical Sciences (AREA)
- Electromagnetism (AREA)
- Measurement Of Radiation (AREA)
Abstract
本发明提供一种适用于全地形的航空γ能谱测量系统及其测量方法，包括航空γ能谱仪、机载GPS系统、IMU、激光雷达扫描仪、数据采集存储计算机、飞机平台。本发明可消除航空γ能谱测量地形因素对测量结果的影响，避免地形因素给航空γ能谱测量造成异常掩盖或出现假异常，提高航空γ能谱测量的质量和精度，提高后期地质的推断解释的准确性，保证航空γ能谱测量技术在在地形平缓测区和地形复杂测区均可有效应用。
Description
技术领域
本发明属于航空γ能谱测量技术领域，具体为一种适用于全地形的航空γ能谱测量系统及其测量方法。
背景技术
现行应用的航空γ能谱测量系统主要由配置大体积NaI(Tl)晶体探测器的航空γ能谱仪、导航定位系统、雷达高度计和数据收录仪等构成。航空γ能谱测量过程在记录每一测点γ能谱各道计数率数据的同时，还记录该测点的位置坐标和飞机离地高度信息。
在现行的常规测量中为便于能谱数据高度修正的程序自动化，采用负指数函数来模拟γ照射量率(计数率、计数)随测量高度的变化规律，采用如下公式进行高度衰减修正：
NS＝Nm·exp[μ·(h0-h)]
式中，NS是换算到规定测量高度h0上的计数率；
Nm是STP等等效高度上经剥离和本底修正的计数率；
h为测量系统雷达高度计测量的离地高度值；
μ为衰减系数。
航空γ能谱各测点的能谱数据是来自该测点所对应整个测量周期作用带的γ辐射，而雷达高度计仅记录飞机所在测点位置的离地高度。因此进行高度衰减改正的前提是地面相对于待改正点为近似平坦且无限大的辐射体，即只适用于平原、戈壁和丘陵等地形比较平缓的测区。当在地形起伏变化的山区进行测量时，只通过高度改正必然产生误差，会导致凸地形上空改正不足，测量结果较真实值偏低；凹地形上空改正过大，测量结果偏高，往往造成异常的掩盖或者出现假异常，给后期的推断解释带来困难。因此在地形切割严重的测区只通过高度改正不能消除地形因素对航空γ能谱测量的影响。
另一种航空γ能谱测量针对地形对航空γ能谱测量的影响，数据处理过程中地形被分割成若干个辐射体，根据辐射体γ场的公式计算探测点处所有辐射体的总γ场与半无限宽辐射体在该探测点同样测量条件下的比值，通过比值对该测点进行地形校正：其中方法一是高度改正过程中除了取测点本身的高度外,还沿测线方向前后各取若干点的高度值,并赋予它们由作用带所决定的不同系数，使高度改正考虑地形因素的影响获得的测量结果更准确；方法二是根据矩形辐射体航空γ辐射场理论，可以把空中任一测点的照射量率认为是地面上n个带状辐射体的累积贡献。如果这n个带状辐射体按照改正点对地面辐射体张角来划分，且改正点相对每个辐射体的张角皆为π/n，分别计算每个辐射体对改正点的影响，即按角度地形改正方法；方法三是作用带中垂直测线的方向将其视为处于不同地形高度的多个线状辐射源，线状辐射源的照射量率进行积分可得作用带在探测器处的照射量。
上述针对地形较复杂区域的航空γ能谱测量，仍采用的是传统测量系统，对地形因素影响修正时，采用的是二维或准二维地形改正方法。当测区为复杂的地形时，二维地形改正有许多局限性，这种改正方法只考虑沿测线剖面方向的地形变化影响，未考虑垂直测线方面以及与测线斜交方向的地形影响，而航空γ能谱仪探测器接收的信号来自整个360°方位作用带，故仅通过二维校正不能完全消除地形因素的影响。
同时上述方法中沿测线剖面方向的地形数据一般采用系统测量过程中采集的GPS高度值与雷达高度作差获得，而GPS高度值与雷达高度的采集频率有限，因此获得的沿测线剖面方向的地形数据的精度也受到很大限制。
随着GIS和数字高程模型(DEM)的普及，DEM应用于航空γ能谱测量地形改正方法中，通过测区内的DEM数据快速划分网格或采用有限元等方法进行三维地形影响校正。地形影响改正所采用的DEM数据不是航空γ能谱测量过程中的实时地形数据，由于测量时间、测量环境等的差别，实时性较差，可能具有误差。同时通常获取的大图幅的DEM数据是由不同扫描日期和不同数据源的数据通过拼接及滤波处理形成，该类型DEM数据并不一定适合作为航空γ能谱单点地形改正所需较小范围作用带内的地形数据。
发明内容
本发明涉及适用于全地形测区的航空γ能谱测量系统，可实现对航空γ能谱测量数据三维地形影响改正，消除地形因素对测量结果的干扰，避免地形因素给航空γ能谱测量结果造成异常掩盖或出现假异常的影响，提高航空γ能谱测量的质量和精度，保证航空γ能谱测量技术在在地形平缓测区和地形复杂测区均可有效应用。
具体的技术方案为：
适用于全地形的航空γ能谱测量系统包括航空γ能谱仪、机载GPS系统、IMU、激光雷达扫描仪、数据采集存储计算机和飞机平台。
机载GPS系统分别向航空γ能谱仪、IMU、激光雷达扫描仪提供PPS信号和GPS信息，以GPS接收机输出的PPS秒脉冲信号作为时间同步基准，对航空γ能谱仪、IMU、激光雷达扫描仪三个子系统采样触发信号进行时间标定，实现各子系统输出结果的时间同步，通过GPS输出信息中的UTC时间，递推出不同测量频率下各子系统测量结果所对应的精确UTC时间，为各子系统测量数据添加相应的时间戳，将不同测量子系统的结果统一到同一个时间轴上，实现在同一时间维度下数据精确采集。航空γ能谱仪、IMU、激光雷达扫描仪将标有时间戳的测量数据包传输到数据采集存储计算机上。
航空γ能谱测量单个测点的测量周期为1秒，机载GPS系统输出的PPS秒脉冲信号中的上升沿作为航空γ能谱仪上一测量周期结束控制信号，并作为下一测量周期的开始控制信号。航空测量为动态测量过程，在单个测量周期中航空γ能谱仪探测器位置是动态变化，将单个测量周期中航空γ能谱仪探测器移动范围的中间点位置作为该测量周期能谱数据对应的测点位置，航空γ能谱仪接收的机载GPS系统提供更新率为2Hz的GPS信息，选用每秒的第2个GPS信息中的坐标位置信息为该周期测量结果对应的测点坐标。
将航空γ能谱测量作用带中所有激光雷达扫描仪获取的激光脚点视为点状放射源，将作用带表面细化为离散点源，航空γ能谱各测量周期的结果为作用带内多个激光脚点对应点源共同作用的结果。
适用于全地形的航空γ能谱测量方法，包括以下步骤：
S1系统获取数据处理
S1.1POS数据处理
(1)GPS数据和IMU数据的融合，经解算获得与IMU相同采样频率的200Hz GPS天线坐标和IMU姿态数据。
(2)以激光雷达扫描仪数据采样频率为基准，对(1)中解算获得200Hz的GPS天线WGS84坐标数据进行内插，得到与激光点采集频率一致的20kHz GPS天线WGS-84坐标数据。
(3)以激光雷达扫描仪数据采样频率为基准，对IMU的200Hz的姿态数据进行内插，得到与激光点采集频率一致的20kHz IMU平台姿态数据。
S1.2激光扫描数据处理
(1)通过激光雷达扫描仪获得的原始数据建立每个激光脚点在激光发射瞬时坐标系下坐标(0,0,ρ)。
(2)通过激光雷达扫描仪瞬时扫描角θ，将每个激光脚点在激光发射瞬时坐标系下坐标(0,0,ρ)转化为激光脚点在扫描仪平台坐标系下坐标(XUA,YUA,ZUA)。
(3)通过偏心改正将激光脚点在扫描仪平台坐标系下坐标(XUA,YUA,ZUA)转化为激光脚点在以GPS天线为原点的IMU平台坐标系下坐标(XGA,YGA,ZGA)。
(4)通过IMU姿态角，将激光脚点在以GPS天线为原点的IMU平台坐标系下坐标(XGA,YGA,ZGA)转化为激光脚点在以GPS天线为原点在当地水平参考坐标系下的坐标(XGIA,YGIA,ZGIA)。
(5)通过GPS天线的WG-S84坐标，将激光脚点在以GPS天线为原点的当地水平参考坐标系下的坐标(XGIA,YGIA,ZGIA)转化为激光脚点在WGS-84坐标系的坐标(X84,Y84,Z84)。
(6)激光脚点在WGS-84坐标系的坐标(X84,Y84,Z84)转化为激光脚点的平面坐标(XN,YE)。
(7)激光脚点以平面坐标(XN,YE)和海拔高度ZH的形式来表示其三维坐标，形成测区范围内的激光点云三维坐标，即第j个激光脚点的三维坐标位置(xfpj,yfpj,zfpj)采用平面坐标系中坐标(xfpj,yfpj)和海拔高度zfpj表示。
S1.3能谱窗口数据地形及高度综合改正公式的建立
(1)建立航空γ能谱仪探测器在以GPS天线为原点的平台坐标系下的坐标
通过偏心测量建立与IMU数据更新频率一致的200Hz航空γ能谱仪探测器以GPS天线为原点在平台坐标系下的坐标(XGD,YGD,ZGD)。
(2)建立航空γ能谱仪探测器以GPS天线为原点的当地水平参考坐标系的坐标
通过S1.1中POS数据解算获得的200Hz IMU姿态角将200Hz的航空γ能谱仪探测器以GPS天线为原点在平台坐标系下坐标(XGD,YGD,ZGD)转化为以GPS天线为原点的当地水平参考坐标系的坐标(XGID,YGID,ZGID)。
(3)建立航空γ能谱仪探测器在WGS-84坐标系的坐标
通过S1.1中POS数据解算获得的200Hz GPS天线WGS-84坐标将200Hz的航空γ能谱仪探测器以GPS天线为原点的当地水平参考坐标系的坐标(XGID,YGID,ZGID)转化为航空γ能谱仪探测器在WGS-84坐标系的坐标(X84D,Y84D,Z84D)。
(4)建立航空γ能谱仪探测器的平面坐标
航空γ能谱仪探测器在WGS-84坐标系下的坐标(X84D,Y84D,Z84D)转化为航空γ能谱仪探测器的平面坐标(XND,YED)。
(5)建立航空γ能谱仪探测器以平面坐标和海拔高度表示的三维坐标
采用航空γ能谱仪探测器平面坐标系下坐标(XNd,YEd)和海拔高度ZHd来表示其三维坐标(XNd,YEd,ZHd)。
(6)构建航空γ能谱仪探测器γ射线信号接收边界曲面
在平坦地区假设航空γ能谱仪探测器静止于某一飞行高度H处，经验总结航空γ能谱探测作用带(地表探测范围)为位于探测器正下方半径为H的圆面积,即航空γ能谱测量作用带的直径约为飞行高度的2倍，因此航空γ能谱仪探测器收γ射线信号的视场角为90°。
根据航空γ能谱仪接收γ射线信号的视场角为90°，构建基于三维坐标的航空γ能谱仪探测器接收射线信号的边界曲面，是为以探测器为顶点、顶角为90°的圆锥曲面，表达式为：
其中(x,y,z)为曲面上任意一点的三维坐标，(XNd,YEd,ZHd)为航空γ能谱仪探测器的三维坐标，采用其平面坐标系下坐标(XNd,YEd)和海拔高度ZHd表示。
(7)实际地形下航空γ能谱仪探测器作用带范围内激光脚点的确定
根据辐射场的可叠加原理，任何形状的辐射源都可视为若干点源的叠加。航空γ能谱测量过程中空中探测器处的粒子注射率可视为作用带内无数个点源叠加作用的结果。将航空γ能谱测量作用带中所有激光雷达扫描仪获取得的激光脚点视为点状放射源，即可将作用带表面细化为离散点源，航空γ能谱各测量周期的结果为作用带内多个激光脚点对应点源共同作用的结果，即探测器处的粒子注射率为作用带内多个点源叠加作用的结果。
单能不带电粒子的辐射场中，γ照射量率与粒子注量率成正比。单能γ辐射点源处于均匀介质中，距点源R处的γ射线粒子注量率Φ点的计算公式为：
其中m为放射性物质的质量；μ为介质对γ射线的衰减系数；K是γ常数[C·m2/(kg·s·g)]；f为γ粒子注量率与γ照射量率的换算系数。
在实际测区地形条件下探测器接收射线信号的边界曲面与地表曲面的交线为航空γ能谱测量作用带的边界，为不规则图形。激光雷达扫描仪获得的地表激光脚点视为地面点状辐射源，第j个激光脚点(点状辐射源)的三维坐标(xfpj,yfpj,zfpj)采用平面坐标系中坐标(xfpj,yfpj)和海拔高度zfpj表示。当激光脚点位置满足式1-3，该激光脚点在作用带范围内，将其视为作用带内的点源，参与航空γ能谱测量地形改正的计算。
其中(XNd,YEd,ZHd)为采用航空γ能谱仪探测器的平面坐标系下坐标(XNd,YEd)和探测器的海拔高度ZHd表示的三维坐标。
(8)理想平坦地形条件下航空γ能谱仪探测器作用带范围内激光脚点的确定
在理想平坦测区航空γ能谱仪探测器静止于某一飞行高度H处，根据探测器接收γ射线信号的90°视场角计算，航空γ能谱探测作用带为位于探测器正下方半径为H的圆面积，即航空γ能谱测量作用带的直径约为飞行高度的2倍；
系统中的激光雷达扫描仪获得的地表激光脚点视为地面点状辐射源，第j个激光脚点(点状辐射源)的三维坐标(xfpj,yfpj,zfpj)采用平面坐标系中坐标(xfpj,yfpj)和海拔高度zfpj表示。
航空γ能谱仪探测器的三维坐标(XNd,YEd,ZHd)采用探测器的平面坐标系下坐标(XNd,YEd)和探测器的海拔高度ZHd来表示。
航空γ能谱仪探测器与获得的第j个激光脚点的平面坐标分别为(XNd,YEd)和(xfpj,yfpj)，取(XNd-xfpj)2+(YEd-yfpj)2为最小值的激光脚点作为探测器在正下方地表上的投影，设该激光脚点的平面坐标为(x′fp y′fp)，海拔高度为z′fp。则假设理想平坦地形条件下各激光脚点海拔相同均为z′fp，各激光脚点的三维坐标为(xfpj,yfpj,z′fp)。探测器的离地高度hi的计算公式为式1-4：
h＝ZHd-z′fp 1-4
当理想平坦地形条件下激光脚点位置满足式1-5，视该激光脚点为作用带内的点源，参与航空γ能谱测量地形改正的计算。
(9)实际地形条件下航空γ能谱仪探测器静止于空中时的粒子注量率计算
实际地形条件下航空γ能谱仪探测器静止于空中时的粒子注量率计算公式如下：
ΦC为实际地形条件下探测器静止于空中时的粒子注量率；
ΦC点j为实际测量地形条件下γ能谱测量作用带内的第j个点源在探测器处产生的γ射线粒子注量率；
k为实际测量地形条件下航空γ能谱测量作用带内点源数量；
m为放射性物质的质量；
μ为介质对γ射线的衰减系数；
K是γ常数[C·m2/(kg·s·g)]；f为γ粒子注量率与γ照射量率的换算系数。
RCj为实际测量地形条件下航空γ能谱测量作用带内第j个点源距探测器的距离，计算公式如下:
式中(XNd,YEd)为探测器的平面坐标系下坐标，ZHd为探测器的海拔高度；(xfpj,yfpj)为作用带内第j个激光脚点(点源)的平面坐标，zfpj为作用带内第j个激光脚点(点源)的海拔高度。
(10)理想平坦地形条件下航空γ能谱仪探测器静止于空中时的粒子注量率计算
理想平坦地形条件下航空γ能谱仪探测器静止于空中时的粒子注量率ΦPi计算公式如下：
ΦP为理想平坦地形条件下探测器静止于空中时的粒子注量率；
ΦP点j为理想平坦地形条件下γ能谱测量作用带内的第j个点源在探测器处产生的γ射线粒子注量率；
n为理想平坦地形条件下航空γ能谱测量作用带内点源数量；
m为放射性物质的质量；
μ为介质对γ射线的衰减系数；
K是γ常数[C·m2/(kg·s·g)]；f为γ粒子注量率与γ照射量率的换算系数。
RPj为理想平坦地形条件下航空γ能谱测量作用带内第j个点源距探测器的距离，计算公式如下:
1-9式中(XNd,YEd)为探测器在平面坐标系中坐标，探测器海拔高度为ZHd，理想平坦地形条件下作用带内第j个激光脚点(视为点源)的平面坐标为(xfpj,yfpj)，海拔高度统一为z′fp。
(11)地形改正公式建立
航空γ能谱测量为动态测量过程，在1秒的测量周期内探测器处于运动状态，沿测线方向进行连续数据采集，每个测点处的计数率为相应测量周期中计数累加的结果。根据POS数据解算的输出频率，将每个测量周期内航空γ能谱仪探测器移动路径视为由探测器在空中200个不同位置组成，位置对应POS数据解算的探测器坐标；每个测量周期测量结果视为探测器处在空中200个不同位置测量共同作用的结果，即每个测量周期的测量对象为探测器静止在空中200个不同位置时所对应的作用带；
每个测量周期为1秒的航空γ能谱测点的各能量窗口计数率地形及高度综合改正公式：
其中N为地形及高度改正后的计数率；N0为地形及高度改正前的计数率。
POS数据解算的输出频率为200Hz，即每个测量周期内探测器移动路径视为由探测器在空中200个不同坐标位置组成。
ΦPiΦCi分别为理想平坦地形条件下和实际地形条件下航空γ能谱仪探测器位于每个测量周期内200个不同坐标位置中的第i个位置时的粒子注量率。
hi为航空γ能谱仪探测器位于每个测量周期内200个不同坐标位置中的第i个位置时的离地高度，h0为能谱地形和高度改正的统一目标高度。
μ为γ射线高度衰减系数。
将式1-6、1-8代入式1-10，每个测点航空γ能谱测量计数率地形及高度综合改正公式：
式中ki和ni分别为实际测量地形条件下和理想平坦地形条件下航空γ能谱仪探测器位于每个测量周期内200个不同坐标位置中的第i个位置时对应作用带内点源数量；
RCij和RPij分别为实际测量地形条件下和理想平坦地形条件下航空γ能谱仪探测器位于每个测量周期内200个不同坐标位置中的第i个位置时对应作用带内第j个点源距探测器的距离，计算公式分别为公式1-7和1-9。
S2航空γ能谱数据处理
(1)分别对航空γ能谱测量获取的总窗计数率(Tc)、钾窗计数率(K)、铀窗计数率(U)、钍窗计数率(Th)进行飞机本底及宇宙射线本底改正；
(2)总窗计数率(Tc)、钾窗计数率(K)、铀窗计数率(U)、钍窗计数率(Th)大气氡气本底改正；
(3)钾道计数率(K)、铀道计数率(U)、钍道计数率(Th)进行康普顿散射改正。
(4)采用S1.3中建立的地形及高度综合改正公式，分别对总窗计数率(Tc)、钾窗计数率(K)、铀窗计数率(U)、钍窗计数率(Th)进行地形及高度综合改正。
(5)将地形及高度综合改正后的钾道计数率(K)、铀道计数率(U)、钍道计数率(Th)换算为对应的核素含量。
本发明可消除地形因素对航空γ能谱测量结果的影响，避免地形因素给航空γ能谱测量造成异常掩盖或出现假异常的影响，提高航空γ能谱测量的质量和精度，提高后期地质推断解释的准确性，保证航空γ能谱测量技术在地形平缓测区和地形复杂测区均可有效应用。
附图说明
图1为本发明的系统结构示意图；
图2为本发明的PPS信号；
图3(1)为本发明空中测量数据处理流程之一；
图3(2)为本发明空中测量数据处理流程之二；
图4为本发明实际地形下γ射线探测器作用带范围示意图；
图5为本发明理想平坦地形条件下γ射线探测器作用带范围示意图；
图6为本发明航空γ能谱测量为动态测量过程探测器作用带范围示意图。
具体实施方式
结合附图说明本发明的具体技术方案。
适用于全地形的航空γ能谱测量系统由航空γ能谱仪、机载GPS系统、IMU、激光雷达扫描仪、数据采集存储计算机和飞机平台构成(如图1)。
机载GPS系统由机载GPS天线及机载GPS接收机构成，分别向航空γ能谱仪、IMU、激光雷达扫描仪提供秒脉冲(PPS)信号和GPS信息(UTC时间、GPS天线位置、速度数据)；其中机载GPS系统向航空γ能谱仪提供2Hz的GPS信息，向IMU提供20Hz的GPS信息，向激光雷达扫描仪提供20Hz的GPS信息；以GPS接收机输出的PPS秒脉冲信号作为时间同步基准，对航空γ能谱仪、IMU、激光雷达扫描仪三个子系统采样触发信号进行时间标定，通过GPS输出信息中的UTC时间，递推出不同测量频率下各子系统测量结果所对应的精确UTC时间，为各子系统测量数据添加相应的时间戳，实现不同的测量子系统在同一时间维度下进行数据精确采集，将结果数据统一到同一个时间轴上。
航空γ能谱仪包括γ射线探测器和能谱仪主机，将航空γ能谱数据和相应的时间戳打包形成带时间戳的航空γ能谱数据包传输到数据采集存储计算机，更新频率为1Hz；IMU将IMU数据和相应的时间戳打包形成带时间戳的IMU数据包传输到数据采集存储计算机，更新频率为200Hz；激光雷达扫描仪将所测激光测距数据和相应的时间戳打包形成带时间戳的激光测距数据包传输到数据采集存储计算机，更新频率为20kHz。
在航空γ能谱仪中GPS输出的秒脉冲信号(PPS)的上升沿作为每个测点测量周期的开始/停止测量控制信号。如图2所示，PPS秒脉冲信号中的上升沿即作为上一测量周期结束控制信号，也作为下一测量周期的开始控制信号。
航空γ能谱测量结果以每个测点的位置坐标和所对应能谱数据的形式表示，航空γ能谱测量是γ射线信号分析累加过程，单个测点能谱的测量周期为1秒。而航空测量为动态测量过程，飞机及γ射线探测器处于移动状态，位置是动态变化的，因此将单个测量周期中航空γ能谱仪探测器的移动范围的中间点位置作为该测量周期能谱数据对应的测点位置。航空γ能谱仪接收的GPS系统提供的GPS信息(GPGGA格式数据)为2Hz更新率，因此选用每秒的第2个GPS信息中的坐标位置信息为该周期测量结果对应的测点坐标。
适用于全地形的航空γ能谱测量方法，包括以下步骤：
S1系统获取数据处理
适用于全地形的航空γ能谱测量系统获取的空中测量数据经过如图3(1)和图3(2)所示数据处理流程，获得最终测量数据。具体的：
S1.1 POS数据处理
(1)GPS数据(GPS天线在WGS84坐标系下坐标)和IMU数据(三轴加速度、三周角速度)的融合，经解算获得与IMU相同采样频率的200Hz GPS天线坐标和IMU姿态数据。
(2)以激光雷达扫描仪数据采样频率为基准，对(1)中解算获得200Hz的GPS天线WGS84坐标数据进行内插，得到与激光点采集频率一致的20kHz GPS天线WGS84坐标数据。
(3)以激光雷达扫描仪数据采样频率为基准，对IMU的200Hz的姿态数据进行内插，得到与激光点采集频率一致的20kHz IMU平台姿态数据。
S1.2激光扫描数据处理
(1)通过激光雷达扫描仪获得的原始数据建立每个激光脚点在激光发射瞬时坐标系下坐标(0,0,ρ)。
(2)通过激光雷达扫描仪瞬时扫描角θ，将每个激光脚点在激光发射瞬时坐标系下坐标(0,0,ρ)转化为激光脚点在扫描仪平台坐标系下坐标(XUA,YUA,ZUA)。
(3)通过偏心改正将激光脚点在扫描仪平台坐标系下坐标(XUA,YUA,ZUA)转化为激光脚点在以GPS天线为原点的IMU平台坐标系下坐标(XGA,YGA,ZGA)。
(4)通过IMU姿态角，将激光脚点在以GPS天线为原点的IMU平台坐标系下坐标(XGA,YGA,ZGA)转化为激光脚点在以GPS天线为原点在当地水平参考坐标系下的坐标(XGIA,YGIA,ZGIA)。
(5)通过GPS天线的WGS84坐标，将激光脚点在以GPS天线为原点的当地水平参考坐标系下的坐标(XGIA,YGIA,ZGIA)转化为激光脚点在WGS-84坐标系的坐标(X84,Y84,Z84)。
(6)激光脚点在WGS-84坐标系的坐标(X84,Y84,Z84)转化为激光脚点的平面坐标(XN,YE)。
(7)激光脚点以平面坐标(XN,YE)和海拔高度ZH的形式来表示其三维坐标，形成测区范围内的激光点云三维坐标，即第j个激光脚点的三维坐标位置(xfpj,yfpj,zfpj)采用平面坐标系中坐标(xfpj,yfpj)和海拔高度zfpj表示。
S1.3能谱窗口数据地形及高度综合改正公式的建立
(1)建立航空γ能谱仪探测器在以GPS天线为原点的平台坐标系下的坐标
通过偏心测量建立与IMU数据更新频率一致的200Hz航空γ能谱仪探测器以GPS天线为原点在平台坐标系下的坐标(XGD,YGD,ZGD)。
(2)建立航空γ能谱仪探测器以GPS天线为原点的当地水平参考坐标系的坐标
通过S1.1中POS数据解算获得的200Hz IMU姿态角将200Hz的航空γ能谱仪探测器以GPS天线为原点在平台坐标系下坐(XGD,YGD,ZGD)转化为以GPS天线为原点的当地水平参考坐标系的坐标(XGID,YGID,ZGID)。
(3)建立航空γ能谱仪探测器在WGS-84坐标系的坐标
通过S1.1中POS数据解算获得的200Hz GPS天线WGS-84坐标将200Hz的航空γ能谱仪探测器为以GPS天线为原点的当地水平参考坐标系的坐标(XGID,YGID,ZGID)转化为探测器在WGS-84坐标系的坐标(X84D,Y84D,Z84D)。
(4)建立航空γ能谱仪探测器的平面坐标
航空γ能谱仪探测器在WGS-84坐标系的坐标(X84D,Y84D,Z84D)转化为探测器的平面坐标(XND,YED)。
(5)建立航空γ能谱仪探测器以平面坐标和海拔高度的形式表示的三维坐标
γ射线探测器的三维坐标(XNd,YEd,ZHd)，采用航空γ能谱仪探测器的平面坐标系下坐标(XNd,YEd)和探测器的海拔高度ZHd来表示。
(6)构建航空γ能谱仪探测器γ射线信号接收边界曲面
在平坦地区假设航空γ能谱仪探测器静止于某一飞行高度H处，经验总结航空γ能谱探测作用带(地表探测范围)为位于探测器正下方半径为H的圆面积,即航空γ能谱测量作用带的直径约为飞行高度的2倍，因此航空γ能谱探测器收γ射线信号的视场角为90°。
根据航空γ能谱仪接收γ射线信号的视场角为90°，构建基于三维坐标的航空γ能谱探测器接收射线信号的边界曲面，是为以探测器为顶点、顶角为90°的圆锥曲面，表达式为：
其中(x,y,z)为曲面上任意一点的三维坐标，(XNd,YEd,ZHd)为航空γ能谱仪探测器的三维坐标，采用其平面坐标系坐标(XNd,YEd)和海拔高度ZHd表示。
(7)实际地形下航空γ能谱仪探测器作用带范围内激光脚点的确定(保护点)
根据辐射场的可叠加原理，任何形状的辐射源都可视为若干点源的叠加。航空γ能谱测量过程中空中探测器处的粒子注射率可视为作用带内无数个点源叠加作用的结果。将航空γ能谱测量作用带中所有激光雷达扫描仪获取得的激光脚点视为点状放射源，即可将作用带表面细化为离散点源，航空γ能谱各测量周期的结果为作用带内多个激光脚点对应点源共同作用的结果，即探测器处的粒子注射率为作用带内多个点源叠加作用的结果。
单能不带电粒子的辐射场中，γ照射量率与粒子注量率成正比。单能γ辐射点源处于均匀介质中，距点源R处的γ射线粒子注量率Φ点的计算公式为：
其中m为放射性物质的质量；μ为介质对γ射线的衰减系数；K是γ常数[C·m2/(kg·s·g)]；f为γ粒子注量率与γ照射量率的换算系数。
√在实际测区地形条件下探测器接收射线信号的边界曲面与地表曲面的交线为航空γ能谱测量作用带的边界，为不规则图形。作用带内的激光脚点视为点源，探测器处的粒子注射率为作用带内多个点源叠加作用的结果，如图4。
激光雷达扫描仪获得的地表激光脚点视为地面点状辐射源，第j个激光脚点(点状辐射源)的三维坐标位置(xfpj,yfpj,zfpj),采用平面坐标系中坐标(xfpj,yfpj)和海拔高度zfpj表示。当激光脚点位置满足式1-3，该激光脚点在作用带范围内，将其视为作用带内的点源，参与航空γ能谱测量地形改正的计算。
其中(XNd,YEd,ZHd)为采用探测器的平面坐标系下坐标(XNd,YEd)和探测器的海拔高度ZHd表示的三维坐标。
(8)理想平坦地形条件下构建航空γ能谱仪探测器作用带范围内的激光脚点的确定(保护点)
在理想平坦测区航空γ能谱仪探测器静止于某一飞行高度H处，根据探测器接收γ射线信号的90°视场角计算，航空γ能谱探测作用带(地表探测范围)为位于探测器正下方半径为H的圆面积，即航空γ能谱测量作用带的直径约为飞行高度的2倍，如图5。
系统中的激光雷达扫描仪获得的地表激光脚点视为地面点状辐射源，第j个激光脚点(点状辐射源)的三维坐标(xfpj,yfpj,zfpj)采用平面坐标系中坐标(xfpj,yfpj)和海拔高度zfpj表示。
航空γ能谱仪探测器的三维坐标(XNd,YEd,ZHd)采用探测器的平面坐标系下坐标(XNd,YEd)和探测器的海拔高度ZHd来表示。
航空γ能谱仪探测器与获得的第j个激光脚点的平面坐标分别为(XNd,YEd)和(xfpj,yfpj)，取(XNd-xfpj)2+(YEd-yfpj)2为最小值的激光脚点作为探测器在正下方地表上的投影，设该激光脚点的平面坐标为(x′fp y′fp)，海拔高度为z′fp。则假设理想平坦地形条件下各激光脚点海拔相同均为z′fp，各激光脚点的三维坐标为(xfpj,yfpj,z′fp)。探测器的离地高度hi的计算公式为式1-4：
h＝ZHd-z′fp 1-4
当理想平坦地形条件下激光脚点位置满足式1-5，视该激光脚点为作用带内的点源，参与航空γ能谱测量地形改正的计算。
(9)实际地形条件下航空γ能谱仪探测器静止于空中时的粒子注量率计算
实际地形条件下航空γ能谱仪探测器静止于空中时的粒子注量率ΦCi计算公式如下：
ΦC为实际地形条件下探测器静止于空中时的粒子注量率；
ΦC点j为实际测量地形条件下γ能谱测量作用带内的第j个点源在探测器处产生的γ射线粒子注量率；
k为实际测量地形条件下航空γ能谱测量作用带内点源数量；
m为放射性物质的质量；
μ为介质对γ射线的衰减系数；
K是γ常数[C·m2/(kg·s·g)]；f为γ粒子注量率与γ照射量率的换算系数。
RCj为实际测量地形条件下航空γ能谱测量作用带内第j个点源距探测器的距离，计算公式如下:
式中(XNd,YEd)为探测器的平面坐标系下坐标，ZHd为探测器的海拔高度；(xfpj,yfpj)为作用带内第j个激光脚点(点源)的平面坐标，zfpj为第j个激光脚点(点源)的海拔高度。
(10)理想平坦地形条件下航空γ能谱仪探测器静止于空中时的粒子注量率计算
理想平坦地形条件下航空γ能谱仪探测器静止于空中时的粒子注量率计算公式如下：
ΦP为理想平坦地形条件下探测器静止于空中时的粒子注量率；
ΦP点j为理想平坦地形条件下γ能谱测量作用带内的第j个点源在探测器处产生的γ射线粒子注量率；
n为理想平坦地形条件下航空γ能谱测量作用带内点源数量；
m为放射性物质的质量；
μ为介质对γ射线的衰减系数；
K是γ常数[C·m2/(kg·s·g)]；f为γ粒子注量率与γ照射量率的换算系数。
RPj为理想平坦地形条件下航空γ能谱测量作用带内第j个点源距探测器的距离，计算公式如下:
1-9式中(XNd,YEd)为探测器在平面坐标系中坐标，探测器海拔高度为ZHd，理想平坦地形条件下作用带内第j个激光脚点(视为点源)的平面坐标为(xfpj,yfpj)，海拔高度统一为z′fp。
(11)地形改正公式建立
航空γ能谱测量为动态测量过程，在1秒的测量周期内探测器处于运动状态，沿测线方向进行连续数据采集，每个测点处的计数率为相应测量周期中计数累加的结果。根据POS数据解算的输出频率(200Hz)，将每个测量周期内航空γ能谱仪探测器移动路径视为由探测器在空中200个不同位置组成，位置对应POS数据解算的探测器坐标；每个测量周期测量结果视为探测器处在空中200个不同位置测量共同作用的结果，即每个测量周期的测量对象为探测器静止在空中200个不同位置时所对应的作用带(如图6)。
每个测量周期为1秒的航空γ能谱测点的各能量窗口计数率地形及高度综合改正公式：
其中N为地形及高度改正后的计数率；N0为地形及高度改正前的计数率。
POS数据解算的输出频率为200Hz，即每个测量周期内探测器移动路径视为由探测器在空中200个不同坐标位置组成。
ΦPiΦCi分别为理想平坦地形条件下和实际地形条件下航空γ能谱仪探测器位于每个测量周期内200个不同坐标位置中的第i个位置时的粒子注量率。
hi为航空γ能谱仪探测器位于每个测量周期内200个不同坐标位置中的第i个位置时的离地高度，h0为能谱地形和高度改正的统一高度。
μ为γ射线高度衰减系数。
将式1-6、1-8代入式1-10，每个测点航空γ能谱测量计数率地形及高度综合改正公式：
式中ki和ni分别为实际测量地形条件下和理想平坦地形条件下航空γ能谱仪探测器位于每个测量周期内200个不同坐标位置中的第i个位置时对应作用带内点源数量；
RCij和RPij分别为实际测量地形条件下和理想平坦地形条件下航空γ能谱仪探测器位于每个测量周期内200个不同坐标位置中的第i个位置时对应作用带内第j个点源距探测器的距离，计算公式分别为公式1-7和1-9。
S2航空γ能谱数据处理
(1)分别对航空γ能谱测量获取的总窗计数率(Tc)、钾窗计数率(K)、铀窗计数率(U)、钍窗计数率(Th)进行飞机本底及宇宙射线本底改正；
(2)总窗计数率(Tc)、钾窗计数率(K)、铀窗计数率(U)、钍窗计数率(Th)进行大气氡气本底改正；
(3)钾道计数率(K)、铀道计数率(U)、钍道计数率(Th)进行康普顿散射改正。
(4)采用S1.3中建立的地形及高度综合改正公式，分别对总窗计数率(Tc)、钾窗计数率(K)、铀窗计数率(U)、钍窗计数率(Th)进行地形及高度综合改正。
(5)将地形及高度综合改正后的钾道计数率(K)、铀道计数率(U)、钍道计数率(Th)换算为对应的核素含量。
Claims (2)
1.适用于全地形的航空γ能谱测量系统，其特征在于，包括航空γ能谱仪、机载GPS系统、IMU、激光雷达扫描仪、数据采集存储计算机和飞机平台；
机载GPS系统分别向航空γ能谱仪、IMU、激光雷达扫描仪提供PPS信号和GPS信息，以GPS接收机输出的PPS秒脉冲信号作为时间同步基准，对航空γ能谱仪、IMU、激光雷达扫描仪三个子系统的采样触发信号进行时间标定，实现各子系统输出结果的时间同步，通过GPS输出信息中的UTC时间，递推出不同测量频率下各子系统测量结果所对应的精确UTC时间，为各子系统的测量数据添加相应的时间戳，将不同测量子系统的结果统一到同一个时间轴上，实现在同一时间维度下数据精确采集；航空γ能谱仪、IMU、激光雷达扫描仪分别将标有时间戳的测量数据包传输到数据采集存储计算机上；
航空γ能谱测量单个测点的测量周期为1秒，机载GPS系统输出的PPS秒脉冲信号中的上升沿作为航空γ能谱仪上一测量周期结束控制信号，并作为下一测量周期的开始控制信号；将单个测量周期中航空γ能谱仪探测器移动范围的中间点位置作为该测量周期能谱数据对应的测点位置，航空γ能谱仪接收机载GPS系统提供更新率为2Hz的GPS信息，选用每秒的第2个GPS信息中的坐标位置信息为该周期测量结果对应的测点坐标；
将航空γ能谱测量作用带中所有激光雷达扫描仪获取的激光脚点视为点状放射源，将作用带表面细化为离散点源，航空γ能谱各测量周期的结果为作用带内多个激光脚点对应点源共同作用的结果。
2.适用于全地形的航空γ能谱测量方法，其特征在于，包括以下步骤：
S1系统获取数据处理
S1.1 POS数据处理
(1)GPS数据和IMU数据的融合，经解算获得与IMU相同采样频率的200Hz GPS天线坐标和IMU姿态数据；
(2)以激光雷达扫描仪数据采样频率为基准，对(1)中解算获得200Hz的GPS天线WGS-84坐标数据进行内插，得到与激光点采集频率一致的20kHz GPS天线WGS-84坐标数据；
(3)以激光雷达扫描仪数据采样频率为基准，对IMU的200Hz的姿态数据进行内插，得到与激光点采集频率一致的20kHz IMU平台姿态数据；
S1.2激光扫描数据处理
(1)通过激光雷达扫描仪获得的原始数据建立每个激光脚点在激光发射瞬时坐标系下坐标(0,0,ρ)；
(2)通过激光雷达扫描仪瞬时扫描角θ，将每个激光脚点在激光发射瞬时坐标系下坐标(0,0,ρ)转化为激光脚点在扫描仪平台坐标系下坐标(XUA,YUA,ZUA)；
(3)通过偏心改正将激光脚点在扫描仪平台坐标系下坐标(XUA,YUA,ZUA)转化为激光脚点在以GPS天线为原点的IMU平台坐标系下坐标(XGA,YGA,ZGA)；
(4)通过IMU姿态角，将激光脚点在以GPS天线为原点的IMU平台坐标系下坐标(XGA,YGA,ZGA)转化为激光脚点在以GPS天线为原点在当地水平参考坐标系下的坐标(XGIA,YGIA,ZGIA)；
(5)通过GPS天线的WG-S84坐标，将激光脚点在以GPS天线为原点的当地水平参考坐标系下的坐标(XGIA,YGIA,ZGIA)转化为激光脚点在WGS-84坐标系的坐标(X84,Y84,Z84)；
(6)激光脚点在WGS-84坐标系的坐标(X84,Y84,Z84)转化为激光脚点的平面坐标(XN,YE)；
(7)激光脚点以平面坐标(XN,YE)和海拔高度ZH的形式来表示其三维坐标，形成测区范围内的激光点云三维坐标，即第j个激光脚点的三维坐标位置(xfpj,yfpj,zfpj)采用平面坐标系中坐标(xfpj,yfpj)和海拔高度zfpj表示；
S1.3能谱窗口数据地形及高度综合改正公式的建立
(1)建立航空γ能谱仪探测器在以GPS天线为原点的平台坐标系下的坐标
通过偏心测量建立与IMU数据更新频率一致的200Hz航空γ能谱仪探测器以GPS天线为原点在平台坐标系下的坐标(XGD,YGD,ZGD)；
(2)建立航空γ能谱仪探测器以GPS天线为原点的当地水平参考坐标系的坐标
通过S1.1中POS数据解算获得的200Hz IMU姿态角将200Hz的航空γ能谱仪探测器以GPS天线为原点在平台坐标系下坐(XGD,YGD,ZGD)转化为以GPS天线为原点的当地水平参考坐标系的坐标(XGID,YGID,ZGID)；
(3)建立航空γ能谱仪探测器在WGS-84坐标系的坐标
通过S1.1中POS数据解算获得的200Hz GPS天线WGS-84坐标将200Hz的航空γ能谱仪探测器以GPS天线为原点的当地水平参考坐标系的坐标(XGID,YGID,ZGID)转化为探测器在WGS-84坐标系的坐标(X84D,Y84D,Z84D)；
(4)建立航空γ能谱仪探测器的平面坐标
航空γ能谱仪探测器在WGS-84坐标系下的坐标(X84D,Y84D,Z84D)转化为探测器的平面坐标(XND,YED)；
(5)建立航空γ能谱仪探测器以平面坐标和海拔高度表示的三维坐标
采用航空γ能谱仪探测器平面坐标系下坐标(XNd,YEd)和海拔高度ZHd来表示其三维位置坐标(XNd,YEd,ZHd)；
(6)构建航空γ能谱仪探测器γ射线信号接收边界曲面
基于三维坐标的航空γ能谱仪探测器接收射线信号的边界曲面表达式为：
其中(x,y,z)为曲面上任意一点的三维坐标；(XNd,YEd,ZHd)为航空γ能谱仪探测器的三维坐标，采用其平面坐标系下坐标(XNd,YEd)和海拔高度ZHd表示；
(7)实际地形下航空γ能谱仪探测器作用带范围内激光脚点的确定
实际地形测区内第j个激光脚点的三维坐标(xfpj,yfpj,zfpj)采用平面坐标系坐标(xfpj,yfpj)和海拔高度zfpj表示；当激光脚点位置满足式1-2，视该激光脚点为作用带内的点源，参与航空γ能谱测量地形改正的计算；
其中(XNd,YEd,ZHd)为采用航空γ能谱仪探测器的平面坐标系下坐标(XNd,YEd)和海拔高度ZHd表示的三维坐标；
(8)理想平坦地形条件下航空γ能谱仪探测器作用带范围内激光脚点的确定
理想平坦地形测区内第j个激光脚点的三维坐标(xfpj,yfpj,zfpj)采用平面坐标系中坐标(xfpj,yfpj)和海拔高度zfpj表示；
航空γ能谱仪探测器的三维坐标(XNd,YEd,ZHd)采用其平面坐标系下坐标(XNd,YEd)和海拔高度ZHd来表示；
取(XNd-xfpj)2+(YEd-yfpj)2为最小值的激光脚点作为航空γ能谱仪探测器在正下方地表上的投影，设该激光脚点的平面坐标为(x′fp y′fp)，海拔高度为z′fp；则理想平坦地形条件下各激光脚点海拔相同均为z′fp，各激光脚点的三维坐标为(xfpj,yfpj,z′fp)；航空γ能谱仪探测器的离地高度hi的计算公式为式1-3：
h＝ZHd-z′fp 1-3
当理想平坦地形条件下激光脚点位置满足式1-4，视该激光脚点为作用带内的点源，参与航空γ能谱测量地形改正的计算；
(9)建立地形改正公式
每个测量周期为1秒的航空γ能谱测点的各能量窗口计数率地形及高度综合改正公式1-5：
其中N为地形及高度改正后的计数率；N0为地形及高度改正前的计数率；
200为由POS数据解算确定的每个测量周期内航空γ能谱仪探测器不同坐标位置个数；
hi为航空γ能谱仪探测器位于每个测量周期内200个不同坐标位置中的第i个位置时的离地高度，h0为高度改正的统一目标高度；
μ为γ射线高度衰减系数；
ki为实际测量地形条件下航空γ能谱仪探测器位于每个测量周期内200个不同坐标位置中的第i个位置时测量作用带内点源数量；
ni为理想平坦地形条件下航空γ能谱仪探测器位于每个测量周期内200个不同坐标位置中的第i个位置时测量作用带内点源数量；
RCij为实际测量地形条件下航空γ能谱仪探测器位于每个测量周期内200个不同坐标位置中的第i个位置时作用带内第j个点源距探测器的距离，计算公式如下1-6:
式1-6中(XNdi,YEdi)为探测器位于第i个位置时平面坐标系下坐标，ZHdi为探测器位于第i个位置时的海拔高度；(xfpij,yfpij)为实际测量地形条件下探测器位于第i个位置时作用带内第j个激光脚点的平面坐标，zfpij为探测器位于第i个位置时作用带内第j个激光脚点的海拔高度；
RPij为理想平坦地形条件下航空γ能谱仪探测器位于每个测量周期内200个不同坐标位置中的第i个位置时作用带内第j个点源距探测器的距离，计算公式如下1-7:
式中(XNdi,YEdi)为探测器位于第i个位置时在平面坐标系中坐标，ZHdi为探测器位于第i个位置时的海拔高度；(xfpij,yfpij)为理想平坦地形条件下探测器位于第i个位置时作用带内第j个激光脚点的平面坐标，z′fpi为探测器位于第i个位置时作用带内各激光脚点的统一海拔高度；
S2航空γ能谱数据处理
(1)分别对航空γ能谱测量获取的总窗计数率Tc、钾窗计数率K、铀窗计数率U、钍窗计数率Th进行飞机本底及宇宙射线本底改正；
(2)分别对总窗计数率Tc、钾窗计数率K、铀窗计数率U、钍窗计数率Th进行大气氡气本底改正；
(3)分别对钾道计数率K、铀道计数率U、钍道计数率Th进行康普顿散射改正；
(4)采用S1.3中建立的地形及高度综合改正公式，分别对总窗计数率Tc、钾窗计数率K、铀窗计数率U、钍窗计数率Th进行地形及高度综合改正；
(5)将地形及高度综合改正后的钾道计数率K、铀道计数率U、钍道计数率Th换算为对应的核素含量。
Priority Applications (1)
| Application Number | Priority Date | Filing Date | Title |
|---|---|---|---|
| CN202210481171.5A CN114895346B (zh) | 2022-05-05 | 2022-05-05 | 适用于全地形的航空γ能谱测量系统及其测量方法 |
Applications Claiming Priority (1)
| Application Number | Priority Date | Filing Date | Title |
|---|---|---|---|
| CN202210481171.5A CN114895346B (zh) | 2022-05-05 | 2022-05-05 | 适用于全地形的航空γ能谱测量系统及其测量方法 |
Publications (2)
| Publication Number | Publication Date |
|---|---|
| CN114895346A CN114895346A (zh) | 2022-08-12 |
| CN114895346B true CN114895346B (zh) | 2024-05-17 |
Family
ID=82719107
Family Applications (1)
| Application Number | Title | Priority Date | Filing Date |
|---|---|---|---|
| CN202210481171.5A Active CN114895346B (zh) | 2022-05-05 | 2022-05-05 | 适用于全地形的航空γ能谱测量系统及其测量方法 |
Country Status (1)
| Country | Link |
|---|---|
| CN (1) | CN114895346B (zh) |
Families Citing this family (1)
| Publication number | Priority date | Publication date | Assignee | Title |
|---|---|---|---|---|
| CN116611026B (zh) * | 2023-05-25 | 2024-01-09 | 中国自然资源航空物探遥感中心 | 一种航空γ能谱数据融合处理方法和系统 |
Citations (6)
| Publication number | Priority date | Publication date | Assignee | Title |
|---|---|---|---|---|
| CN201173973Y (zh) * | 2008-03-10 | 2008-12-31 | 成都理工大学 | 手持式一体化多功能γ能谱仪 |
| RU2388018C1 (ru) * | 2009-06-26 | 2010-04-27 | Александр Прокопьевич Елохин | Способ дистанционного измерения загрязнения радионуклидами подстилающей поверхности в следе радиоактивного выброса радиационно-опасных предприятий и система для его осуществления |
| CN105783875A (zh) * | 2016-05-10 | 2016-07-20 | 绵阳天眼激光科技有限公司 | 一种集成非扫描激光雷达的航空摄影测量系统及方法 |
| CN108508473A (zh) * | 2018-05-30 | 2018-09-07 | 国家海洋技术中心 | 基于水下自主航行器的水体γ辐射立体监测系统及方法 |
| CN112904402A (zh) * | 2021-02-03 | 2021-06-04 | 中国自然资源航空物探遥感中心 | 航空γ能谱测量质量控制方法 |
| CN113189636A (zh) * | 2021-05-10 | 2021-07-30 | 核工业航测遥感中心 | 基于多能量响应的航空γ能谱分析方法 |
-
2022
- 2022-05-05 CN CN202210481171.5A patent/CN114895346B/zh active Active
Patent Citations (6)
| Publication number | Priority date | Publication date | Assignee | Title |
|---|---|---|---|---|
| CN201173973Y (zh) * | 2008-03-10 | 2008-12-31 | 成都理工大学 | 手持式一体化多功能γ能谱仪 |
| RU2388018C1 (ru) * | 2009-06-26 | 2010-04-27 | Александр Прокопьевич Елохин | Способ дистанционного измерения загрязнения радионуклидами подстилающей поверхности в следе радиоактивного выброса радиационно-опасных предприятий и система для его осуществления |
| CN105783875A (zh) * | 2016-05-10 | 2016-07-20 | 绵阳天眼激光科技有限公司 | 一种集成非扫描激光雷达的航空摄影测量系统及方法 |
| CN108508473A (zh) * | 2018-05-30 | 2018-09-07 | 国家海洋技术中心 | 基于水下自主航行器的水体γ辐射立体监测系统及方法 |
| CN112904402A (zh) * | 2021-02-03 | 2021-06-04 | 中国自然资源航空物探遥感中心 | 航空γ能谱测量质量控制方法 |
| CN113189636A (zh) * | 2021-05-10 | 2021-07-30 | 核工业航测遥感中心 | 基于多能量响应的航空γ能谱分析方法 |
Non-Patent Citations (2)
| Title |
|---|
| Gamma Ray and Neutron Sensors for Remote Monitoring Using Aerial Robotic Platforms;Monia KAZEMEINI et al.;《Sensors & Transducers》;20190131;第229卷;第47-54页 * |
| 航空γ能谱测量地形校正方法研究;刘秋实 等;《航空γ 能谱测量地形校正方法研究》;20210930;第44卷(第9期);第1-8页 * |
Also Published As
| Publication number | Publication date |
|---|---|
| CN114895346A (zh) | 2022-08-12 |
Similar Documents
| Publication | Publication Date | Title |
|---|---|---|
| Connor et al. | Airborne radiation mapping: overview and application of current and future aerial systems | |
| US10473794B2 (en) | Radiation detection device and method | |
| Sanada et al. | Temporal variation of dose rate distribution around the Fukushima Daiichi nuclear power station using unmanned helicopter | |
| JP2016515204A5 (zh) | ||
| WO2006078310A2 (en) | Method and apparatus for determining a location associated with an image | |
| CN114895346B (zh) | 适用于全地形的航空γ能谱测量系统及其测量方法 | |
| US20160259043A1 (en) | Method for determining trajectories of moving physical objects in a space on the basis of sensor data of a plurality of sensors | |
| CN111505608B (zh) | 一种基于星载激光单片足印影像的激光指向在轨检校方法 | |
| Chen et al. | In-flight performance of the Advanced Radiation Detector for UAV Operations (ARDUO) | |
| Girdwood et al. | Design and field campaign validation of a multi-rotor unmanned aerial vehicle and optical particle counter | |
| Moore et al. | Autonomous Inspection of Electrical Transmission Structures with Airborne UV Sensors-NASA Report on Dominion Virginia Power Flights of November 2016 | |
| CN103644907A (zh) | 一种基于双卫星平台的脉冲星角位置测量系统及方法 | |
| RU2654883C2 (ru) | Способ определения временной привязки производимых с космического аппарата снимков земной поверхности | |
| JP6143219B2 (ja) | 原子力施設近隣の放射性セシウム沈着量の測定方法 | |
| WO2022126085A1 (en) | Camera triggering and multi-camera photogrammetry | |
| CN113640848A (zh) | 无人机的地面激光足印数据采集方法、系统、介质及设备 | |
| Boyns | Aerial Radiological Measuring System (ARMS): systems, procedures and sensitivity (1976) | |
| CN110873890A (zh) | 一种航空放射性谱数据大气氡修正方法 | |
| CN115268491A (zh) | 一种基于北斗+无人机网格码高速场的数据采集、处理系统 | |
| CN111623802B (zh) | 一种基于足印影像的星载激光测高仪在轨指向检校方法 | |
| Cardarelli II et al. | Environmental Protection Agency (EPA) airborne gamma spectrometry system for environmental and emergency response surveys | |
| Gjøystdal et al. | One-array and two-array location capabilities | |
| Consiglio et al. | Nuclear Emulsion technique for volcanoes radiography with cosmic ray muons: status of art and future perspectives | |
| Wasiolek et al. | AMS/NRCan Joint Survey Report: Aerial Campaign | |
| Haber et al. | AMS/DSA Joint Survey Report |
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 |