CN114726433A - 一种对流层散射传输损耗及传播延迟计算方法 - Google Patents
一种对流层散射传输损耗及传播延迟计算方法 Download PDFInfo
- Publication number
- CN114726433A CN114726433A CN202210220518.0A CN202210220518A CN114726433A CN 114726433 A CN114726433 A CN 114726433A CN 202210220518 A CN202210220518 A CN 202210220518A CN 114726433 A CN114726433 A CN 114726433A
- Authority
- CN
- China
- Prior art keywords
- sub
- scattering
- point
- receiving
- ray
- Prior art date
- Legal status (The legal status is an assumption and is not a legal conclusion. Google has not performed a legal analysis and makes no representation as to the accuracy of the status listed.)
- Pending
Links
- 230000005540 biological transmission Effects 0.000 title claims abstract description 90
- 239000005436 troposphere Substances 0.000 title claims abstract description 59
- 238000004364 calculation method Methods 0.000 title claims abstract description 49
- 238000000034 method Methods 0.000 claims description 108
- 238000001228 spectrum Methods 0.000 claims description 16
- 239000013598 vector Substances 0.000 claims description 12
- 238000013213 extrapolation Methods 0.000 claims description 10
- 238000010521 absorption reaction Methods 0.000 claims description 7
- 230000001133 acceleration Effects 0.000 claims description 4
- 230000010354 integration Effects 0.000 claims description 4
- 230000001174 ascending effect Effects 0.000 claims description 3
- 230000005484 gravity Effects 0.000 claims description 3
- HOWHQWFXSLOJEF-MGZLOUMQSA-N systemin Chemical compound NCCCC[C@H](N)C(=O)N[C@@H](CCSC)C(=O)N[C@@H](CCC(N)=O)C(=O)N[C@@H]([C@@H](C)O)C(=O)N[C@@H](CC(O)=O)C(=O)OC(=O)[C@@H]1CCCN1C(=O)[C@H]1N(C(=O)[C@H](CC(O)=O)NC(=O)[C@H](CCCN=C(N)N)NC(=O)[C@H](CCCCN)NC(=O)[C@H](CO)NC(=O)[C@H]2N(CCC2)C(=O)[C@H]2N(CCC2)C(=O)[C@H](CCCCN)NC(=O)[C@H](CO)NC(=O)[C@H](CCC(N)=O)NC(=O)[C@@H](NC(=O)[C@H](C)N)C(C)C)CCC1 HOWHQWFXSLOJEF-MGZLOUMQSA-N 0.000 claims description 2
- 108010050014 systemin Proteins 0.000 claims description 2
- 230000008859 change Effects 0.000 description 16
- 230000008569 process Effects 0.000 description 10
- 238000010586 diagram Methods 0.000 description 9
- 238000004458 analytical method Methods 0.000 description 7
- 238000004891 communication Methods 0.000 description 7
- 238000005452 bending Methods 0.000 description 4
- 230000000694 effects Effects 0.000 description 4
- 230000007246 mechanism Effects 0.000 description 3
- 238000006243 chemical reaction Methods 0.000 description 2
- 238000001514 detection method Methods 0.000 description 2
- 230000035772 mutation Effects 0.000 description 2
- 230000005855 radiation Effects 0.000 description 2
- 230000009466 transformation Effects 0.000 description 2
- XLYOFNOQVPJJNP-UHFFFAOYSA-N water Chemical compound O XLYOFNOQVPJJNP-UHFFFAOYSA-N 0.000 description 2
- 230000004075 alteration Effects 0.000 description 1
- 238000004422 calculation algorithm Methods 0.000 description 1
- 238000012512 characterization method Methods 0.000 description 1
- 238000010276 construction Methods 0.000 description 1
- 230000007547 defect Effects 0.000 description 1
- 239000006185 dispersion Substances 0.000 description 1
- 230000005284 excitation Effects 0.000 description 1
- 239000000284 extract Substances 0.000 description 1
- 230000006872 improvement Effects 0.000 description 1
- 230000007774 longterm Effects 0.000 description 1
- 230000003287 optical effect Effects 0.000 description 1
- 230000000737 periodic effect Effects 0.000 description 1
- 230000000630 rising effect Effects 0.000 description 1
- 238000004088 simulation Methods 0.000 description 1
- 238000006467 substitution reaction Methods 0.000 description 1
Images
Classifications
-
- H—ELECTRICITY
- H04—ELECTRIC COMMUNICATION TECHNIQUE
- H04B—TRANSMISSION
- H04B7/00—Radio transmission systems, i.e. using radiation field
- H04B7/22—Scatter propagation systems, e.g. ionospheric, tropospheric or meteor scatter
-
- G—PHYSICS
- G01—MEASURING; TESTING
- G01W—METEOROLOGY
- G01W1/00—Meteorology
-
- 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
- Y02E—REDUCTION OF GREENHOUSE GAS [GHG] EMISSIONS, RELATED TO ENERGY GENERATION, TRANSMISSION OR DISTRIBUTION
- Y02E60/00—Enabling technologies; Technologies with a potential or indirect contribution to GHG emissions mitigation
Landscapes
- Engineering & Computer Science (AREA)
- Environmental & Geological Engineering (AREA)
- Life Sciences & Earth Sciences (AREA)
- Atmospheric Sciences (AREA)
- Biodiversity & Conservation Biology (AREA)
- Ecology (AREA)
- Environmental Sciences (AREA)
- Computer Networks & Wireless Communication (AREA)
- Signal Processing (AREA)
- Radar Systems Or Details Thereof (AREA)
Abstract
一种对流层散射传输损耗及传播延迟计算方法，步骤如下：通过数值气象模型获取散射传播路径上的气象数据；基于气象数据计算发射天线主轴方向的电磁波传播路径和接收天线主轴方向的电磁波传播路径；对发射波束和接收波束进行剖分，从得到的发射子波束的射线路径和接收子波束的射线路径中截取从发射天线到散射点再到接收天线的路径，即为散射子波束对流层散射传播路径；计算各子波束公共散射体的体积；计算散射子波束的接收功率及散射子波束的传播时间；计算对流层散射传输损耗与传播延迟。本发明可以解决现有对流层散射模型无法准确反映不同区域不同时间大气环境变化对传输损耗和传播时延影响的问题，可以更精确地计算对流层散射损耗及传播延迟。
Description
技术领域
本发明属于对流层散射通信技术领域，尤其涉及一种对流层散射传输损耗及传输延迟的计算方法。
背景技术
对流层散射通信链路存在传输损耗大、多径效应严重、易受大气环境影响明显的不足，因此，为了提高应用对流层散射链路的系统性能，需要准确掌握对流层散射链路的传输损耗特性和传播延迟特性。除了半经验估计模型外，基于散射传播机制还产生了多种损耗预测模型，如张明高提出的基于广义散射截面的对流层散射传输损耗计算方法，李磊提出的基于抛物方程法的对流层散射传输损耗计算方法，Ergin Dinc提出的基于射线的对流层散射信道模型等。
对流层散射还存在明显的多径效应，电磁波通过不同路径到达接收机，接收信号存在传播延迟和时延扩展，现有模型气象特性表征能力不足，通常是假设电磁波沿直线光速传播，如Sunde从链路几何构形出发，基于对称链路和光滑球面推导的平均时延差的粗略计算公式；Bello在湍流非相干散射传播机制上，根据Boor-Gordon的散射横截面理论，以及二维平面假设从散射链路几何构形出发推导给出的散射时延功率谱计算方法；张明高在广义散射截面理论基础上结合散射链路几何构形推导的归一化延迟功率谱；Ergin Dinc基于射线的对流层散射信道群延迟模型等。
上述模型在辅助对流层散射通信链路的设计及散射系统设备参数的确定等方面发挥了重要作用，但随着对流层散射应用的不断拓展，基于对流层散射链路的时间同步系统、超视距探测系统对准确预测和分析不同气象条件下的传输损耗及传播延迟的预测和分析能力提出了更高要求。而现有传输损耗及传播延迟估计模型中气象参数一般使用气象参数经验拟合值或统计结果平均值来描述，不能给出不同区域不同时间的具体气象环境信息；而且，计算过程中均假设电磁波在大气中是沿直线以恒定光速传播，无法有效反映电磁波在复杂多变的大气环境下传播时发生的路径弯曲、传播速度变化以及大气吸收衰减等，在对流层散射传输损耗及传播延迟的计算精度上仍有改善的空间。
发明内容
本发明的目的是提供一种基于数值气象模型的对流层散射传输损耗及传播时延计算方法，可以更精确地计算对流层散射传输损耗及传播延迟。
为了实现上述目的，本发明采取如下的技术解决方案：
一种对流层散射传输损耗及传播延迟计算方法，包括以下步骤：
S1、通过数值气象模型获取散射传播路径上的气象数据；所述数值气象模型以网格形式存储气象数据，水平方向网格间隔为经纬度分辨率，竖直方向网格为等压面，基于数值气象模型采用插值法或外推法在获取指定位置的气象数据，所述气象数据包括温度、比湿度、气压、风速；
S2、基于通过数值气象模型获取的气象数据，计算发射天线主轴方向的电磁波传播路径和接收天线主轴方向的电磁波传播路径；
S3、对发射波束进行横纵方向的剖分及路径旋转，得到各发射子波束的射线路径，对接收波束进行横向剖分、路径旋转及射线对准，得到有效范围内的接收子波束的射线路径，从得到的发射子波束的射线路径和有效范围内的接收子波束的射线路径中截取从发射天线到散射点再到接收天线的路径，即为散射子波束对流层散射传播路径；
S4、计算各子波束公共散射体的体积，所述子波束公共散射体为发射子波束和接收子波束相交的部分；
S5、计算散射子波束的接收功率及散射子波束的传播时间；
S6、计算对流层散射传输损耗与传播延迟；
对流层散射传输损耗为：
其中，Pr_ij为序号下标为ij的散射子波束的接收功率，Pt为发射天线的发射功率；
将所有散射子波束的接收功率归一化，并按升序排列所有散射子波束的传播时间，得到对流层散射链路的归一化延迟功率谱，即为对流层散射传播延迟。
更具体的，所述步骤S1中基于数值气象模型获取指定位置的气象数据的步骤具体如下：
S101、将指定位置的大地高Hg转换为位势高
其中，
为指定位置的纬度，
为地球椭球面法向重力，γ45°为纬度45°的重力加速度的标准值，
为局部地球半径；
S102、通过比较插值点与数值气象模型的等压面的位势高相对位置关系确定相应的插值方法和外推方法，包括以下三种情况：
(1)插值点位势高在等压面内部时，通过插值点W的经纬度和位势高确定数值气象模型中与插值点W邻近的8个网格点m1～m8，分别对具有相同经纬度的4组上下相对的两个网格点向插值点位势高处插值：竖直方向上，温度、比湿度和风速采用线性插值，气压采用指数模型插值；水平方向上，温度、比湿度、风速、气压均采用双线性差值；
(2)插值点位势高在等压面最底层以下时，通过插值点W的经纬度和位势高确定与插值点W邻近的8个网格点m1～m8，对于8个网格点中的上层等压面的4个网格点，水平方向上使用双线性插值法获取插值点在等压面上对应位置处的气象参数，竖直方向上分别采用以下方法外推计算插值点在等压面上对应位置处的气象参数：气压使用指数模型；温度Tj＝Tv_j-0.0065hg_j，Tv_j为在等压面上通过水平双线性插值获取的温度；比湿度qj＝qv_j，qv_j为在等压面上通过水平双线性插值获取的比湿度；风速使用幂定律插值；
(3)插值点位势高位于等压面最高层以上时，通过插值点W的经纬度和位势高确定与插值点W邻近的8个网格点m1～m8，对于8个网格点中的下层等压面的4个网格点，水平方向上使用双线性插值法获取插值点在等压面上对应位置处的气象参数，竖直方向上分别采用以下方法外推计算插值点在等压面上对应位置处的气象参数：比湿度为零；气压在水平双线性插值的基础上使用指数模型外推；温度和风速分量使用CIRA86国际参考大气模型获取。
更具体的，所述步骤S2中，采用射线追踪法计算发射天线主轴方向的电磁波传播路径，步骤如下：
S201、以发射天线所在位置为坐标原点，基于右手螺旋准则建立局部直角坐标系xyz，坐标系x轴为射线的初始方向在地球表面切平面投影，z轴为地球球心到发射机所在位置的方向向量，α0和β0分别为发射天线主轴指向的方位角和俯仰角；
在该局部直角坐标系xyz中电磁波传播路径满足线性微分方程组
式中的n为大气折射指数，N为大气折射率；
上述线性微分方程组在初始边界条件
约束下以及对流层顶高限制下，基于数值气象模型得到的气象数据，通过数值法迭代求解即可获得发射天线主轴方向的电磁波传播路径；
对于接收天线，以接收天线作为起始点采用相同步骤计算，得到天线主轴方向的电磁波传播路径。
更具体的，所述步骤S3中，采用矩形截面剖分法对发射波束进行剖分，分别沿发射波束截面的横轴和纵轴以等角度间隔dw将发射波束剖分为K等份，将发射波束分解为多个紧密排列且互不重叠的正四棱锥形的发射子波束；
将发射天线主轴方向的电磁波传播路径ltransmit_ray_main从局部直角坐标系xyz 转换为球坐标系BLH，然后将球坐标系BLH下发射天线主轴方向的电磁波传播路径的初始指向旋转到各发射子波束射线的初始指向，序号下标为ij的发射子波束的射线路径ltransmit_ray_ij＝ltransmit_ray_main(B,L,H)+(dw·i,dw·j,0)，ltransmit_ray_main(B,L,H)为球坐标系BLH下的发射天线主轴方向的电磁波传播路径，i,j＝-K/2,-K/2+1,…,K/2-1,K/2。
更具体的，所述步骤S3中，采用横向剖分法对接收波束进行剖分，在接收波束截面以等角度间隔dw将接收波束横向剖分为K等份，得到方位角相同而俯仰角不同的若干接收子波束；将发射天线主轴方向的电磁波传播路径 ltransmit_ray_main从局部直角坐标系xyz转换为球坐标系BLH，然后将球坐标系 BLH下发射天线主轴方向的电磁波传播路径的初始指向旋转到各发射子波束射线的初始指向，序号下标为ij的发射子波束的射线路径ltransmit_ray_ij＝ltransmit_ray_main(B,L,H)+(dw·i,dw·j,0)，ltransmit_ray_main(B,L,H)为球坐标系BLH下的发射天线主轴方向的电磁波传播路径， i＝-K/2,-K/2+1,…,K/2-1,K/2；对于每一条发射子波束，计算接收子波束射线与该发射子波束射线之间的最小距离点以及这两个最小距离点之间的最小距离，获取这两个最小距离点关于接收天线发射点的夹角，按照夹角大小旋转接收子波束射线的方位角，直到两个最小距离点之间的最小距离小于设定阈值，计算此时接收子波束射线的初始方位角与接收天线主轴的方位角夹角，将方位角夹角超过波束宽度的接收子波束射线舍弃，得到有效范围内的接收子波束的射线路径。
更具体的，所述步骤S4中，通过以下步骤计算子波束公共散射体的体积：
在散射点S、发射点T和接收点R构成的平面内，以散射点S为起点，分别向发射点T、接收点R方向作子波束射线路径切向量，切向量长度分别设定为从发射点T到散射点S的实际射线路径长度以及从散射点S到接收点R的实际射线路径长度，切向量的终点T'和R'即为虚拟发射点和虚拟接收点；
在虚拟发射点T'、散射点S和虚拟接收点R'构成的平面内，以虚拟发射点 T'为原点o'，虚拟发射点T'到散射点S方向为x'轴，地球球心到虚拟发射点T' 方向为y'轴，基于右手螺旋准则构建三维直角坐标系x'y'z'，在该三维直角坐标系x'y'z'中，T'S和SR'即为电磁波的虚拟直线传播路径，将T'S和SR'分别绕原点o'和虚拟接收点R'在x'o'y'平面中旋转0.5dw和-0.5dw，获得发射子波束和接收子波束的上下边界，进一步获得发射子波束和接收子波束的上下边界交点P1～P4在x'y'z'中的坐标；通过P1～P4在x'y'z'中的坐标可获得子波束公共散射体各顶点坐标；将子波束公共散射体沿x'轴方向细分若干个梯形体，计算各梯形体的体积，所有梯形体的体积之和即为子波束公共散射体的体积。
更具体的，所述步骤S5中，采用双基地雷达方程计算散射子波束的接收功率的步骤如下：
计算各梯形体的接收功率，散射子波束中第l个梯形体的接收功率为：
其中，Pt_ij为序号下标为ij的散射子波束发射功率，Gt、Gr分别为发射天线增益和接收天线增益，gt、gr分别为发射天线的方向性函数和接收天线的方向性函数，ξ为电磁波波长，dVij_l为第l个梯形体的体积，Rij_l为第l个梯形体的中心点到接收天线的距离，Sij_l为第l个梯形体的中心点到发射天线的距离，σij_l为第l个梯形体的散射截面；
序号下标为ij的散射子波束的接收功率Pr_ij为该散射子波束中全部梯形体的接收功率之和，Pr_ij＝∑Pr_ij_l。
更具体的，所述步骤S5中，采用折射率积分法计算散射子波束的传播时间τij，
式中的St为从发射点到散射点的射线路径，Sr为从散射点到接收点的射线路径，n为射线路径上的大气折射指数，c为光速，ds为曲线微元。
优选的，所述步骤S6中，先基于子波束的传输损耗Lij对散射子波束的接收功率进行修正，修正后的散射子波束的接收功率为：
式中的Pr_ij为序号下标为ij的散射子波束的接收功率，Pt_ij为序号下标为ij的散射子波束的发射功率，
为序号下标为ij的散射子波束的大气吸收损耗，再根据修正后的散射子波束的接收功率Pr'_ij计算对流层散射传输损耗，对流层散射传输损耗
优选的，根据修正后的散射子波束的接收功率Pr'_ij计算对流层散射传输延时。
由以上技术方案可知，本发明将数值气象模型应用于对流层散射传输损耗与传播时延计算，使其具备了比现有模型更详细准确的气象分析能力；将射线追踪法应用于对流层散射路径计算，首次在对流层散射模型中考虑了路径弯曲和延迟对传输损耗和传播延迟的影响，使传播时延和传输损耗计算结果更加精确可靠；提出了基于湍流非相干散射机制的对流层散射波束剖分和射线对准的方法，实现了对对流层散射过程的准确模拟；给出了精确计算公共散射体体积的新方法，提高了传输损耗计算准确度。将本发明应用于对流层散射链路可分析不同地区不同气象环境下传输损耗与传播时延的变化特征，为对流层散射设备参数设计和散射链路性能分析提供参考。
附图说明
图1为本发明方法的流程图；
图2为数值气象模型气象数据插值计算示意图；
图3为局部直角坐标系下的三维射线追踪法示意图；
图4为发射波束剖分示意图；
图5为接收波束剖分示意图；
图6为射线对准法示意图；
图7为子波束散射公共体示意图；
图8为虚拟发射点和虚拟接受点示意图；
图9为散射公共体细分梯形体示意图；
图10为链路2063采用本发明方法和ITU-R P.617-2、ITU-R P.617-5方法在不同仰角条件下的传输损耗；
图11为链路2307采用本发明方法和ITU-R P.617-2、ITU-R P.617-5方法在不同仰角条件下的传输损耗；
图12为链路1441的地理位置示意图；
图13为链路2305的地理位置示意图；
图14为链路1441和2305传输损耗的逐小时变化图；
图15a至图15f分别为本发明方法和Bello方法计算的6条散射链路的归一化延迟功率谱；
图16a和图16b分别为1441链路和2305链路在2020年8月25日、26日和27日的逐小时归一化延迟功率谱热力图。
以下结合附图对本发明的具体实施方式作进一步详细地说明。
具体实施方式
下面将结合本发明实施例中的附图，对本发明实施例的技术方案进行清楚、完整地描述，显然，所描述的实施例仅仅是本发明一部分实施例，而不是全部的实施例。基于本发明中的实施例，本领域普通技术人员在没有做出创造性劳动前提下获得的所有其他实施例，都属于本发明保护的范围。
现有的对流层散射传输损耗及传播延迟的计算方法不能有效反映不同区域不同时间的气象环境对对流层散射电磁波传播的影响，无法准确计算电磁波在不同大气环境下传播时的路径变化、速度变化、大气吸收衰减和接收功率等，从而影响了对流层散射传输损耗及传播延迟计算的精度。为了解决以上问题，本发明的基本思路是：首先以数值气象模型为基础提取散射链路计算区域的气象参数，在气象参数支持下计算发射天线及接收天线主轴方向的电磁波传播路径；采用波束剖分和射线对准法模拟散射过程，通过旋转天线主轴的电磁波射线路径获得发射子波束及接收子波束的射线路径，再截取从发射端到散射点再到接收端路径作为各子波束对流层散射传播路径；然后计算子波束公共散射体的体积，基于子波束公共散射体的体积计算子波束接收功率及传播时间，最后根据各子波束接收功率和传播时间计算对流层散射传输损耗和传播延迟。
图1为本发明方法的流程图，下面结合图1，对本发明的方法作进一步的说明，计算过程中需要获取对流层散射链路的相关参数，这些对流层散射链路的相关参数包括发射天线及接收天线的位置及海拔高度、发射天线及接收天线的视平线夹角、发射天线及接收天线的仰角、电磁波频率、波束宽度及计算时间；在获取对流层散射链路的参数后，如图1所示，本发明方法包括以下步骤：
S1、通过数值气象模型获取散射传播路径(散射链路计算区域)上的气象数据；本发明使用的数值气象模型以网格形式存储气象数据，水平方向网格间隔为数值气象模型经纬度分辨率，竖直方向网格为等压面，基于数值气象模型采用插值法或外推法在获取指定位置(散射链路计算区域)的气象数据，这些气象数据包括温度、比湿度、气压、风速(风速包括水平风东向分量及水平风北向分量)。
基于数值气象模型获取指定位置的气象数据的方法具体如下：
S101、将指定位置的大地高Hg转换为位势高
数值气象模型一般使用位势高程系统，因此需要进行大地高和位势高的转换，转换公式为：
其中，
为指定位置的纬度，
为地球椭球面法向重力，γ45°为纬度45°的重力加速度的标准值，
为局部地球半径，
γ45o＝9.80665ms-2，
S102、在获取指定位置的位势高后，通过比较插值点与数值气象模型的等压面的位势高的相对位置关系来确定相应的插值方法和外推方法，包括以下三种情况：
(1)插值点位势高在等压面内部时，如图2所示，首先通过插值点W的经纬度和位势高确定数值气象模型中与插值点W邻近的8个网格点m1～m8，分别对具有相同经纬度的4组上下相对的两个网格点(m1和m5，m2和m6，m3和 m7，m4和m8)向插值点位势高处插值：
竖直方向上，温度、比湿度和风速采用线性插值，气压采用指数模型插值；
线性插值：
式中的m_valuej表示上、下网格点竖直方向的气象数据(温度、比湿度和风速)插值结果，m_valueg_j和m_valueg_j+4分别为上、下网格点的气象数据数值，hg_j为插值点位势高到上网格点的距离，hg_j+4为插值点位势高到下网格点的距离，此处j＝1，2，3，4，如图2中网格点下标序号；
指数模型插值：
式中的Pj为竖直方向的气压插值结果，Pc_j和Pc_j+4为分别为下网格点和上网格点通过指数外推公式外推计算得到的气压值，以Pc_j为例，其外推公式为：
Tv_g_j＝Tg_j(1+0.6077qg_j)，Pg_j为下网格点处的气压值，Tg_j为下网格点处的温度值，qg_j为下网格点处的比湿度值，g0为重力加速度常数，Rd为干空气气体常数，Rd＝287.054 J/K·kg；上网格点气压外推值Pc_j+4的计算方法与之相同；
水平方向上，所有的气象数据(温度、比湿度、风速、气压)均采用双线性差值；在通过线性插值和指数模型插值获得与插值点W同一位势高处的4个插值点w1,5、w2,6、w3,7、w4,8的气象数据差值结果m_value(w1,5)、 m_value(w2,6)、m_value(w3,7)、m_value(w4,8)后，使用双线性插值获取插值点W的最终差值结果m_value(W)：
(2)插值点位势高在等压面最底层以下时，同样通过插值点W的经纬度和位势高确定与插值点W邻近的8个网格点m1～m8，对于8个网格点中的上层等压面的4个网格点，水平方向上使用双线性插值法获取插值点在等压面上对应位置处的气象参数，竖直方向上分别采用以下方法外推计算插值点在等压面上对应位置处的气象参数：
气压使用指数模型；
温度以恒定的温度梯度-0.0065℃/m线性变化：Tj＝Tv_j-0.0065hg_j， Tv_j为在等压面上通过水平双线性插值获取的温度；
比湿度在最底层等压面以下保持恒定不变：qj＝qv_j，qv_j为在等压面上通过水平双线性插值获取的比湿度；
风速使用幂定律插值：
uv_j和vv_j分别为在等压面上通过水平双线性插值获取的水平风东向分量和水平风北向风量，hj为插值点的海拔高度，hv_j为网格点海拔高度，m为风切变系数，可通过查表获取；
(3)插值点位势高位于等压面最高层以上时，通过插值点W的经纬度和位势高确定与插值点W邻近的8个网格点m1～m8，对于8个网格点中的下层等压面的4个网格点，水平方向上使用双线性插值法获取插值点在等压面上对应位置处的气象参数，竖直方向上分别采用以下方法外推计算插值点在等压面上对应位置处的气象参数：
比湿度设定为零；
气压在水平双线性插值的基础上使用指数模型外推；
温度和风速分量使用CIRA86国际参考大气模型获取。
S2、基于通过数值气象模型获取的气象数据，计算发射天线主轴方向的电磁波传播路径和接收天线主轴方向的电磁波传播路径；
本实施例采用基于几何光学射线理论的射线追踪法，基于步骤S1获取的气象数据计算发射天线主轴方向的电磁波传播路径和接收天线主轴方向的电磁波传播路径，发射天线主轴方向的电磁波传播路径和接收天线主轴方向的电磁波传播路径的计算过程相同，下面以发射天线主轴方向的电磁波传播路径为例进行说明：
S201、建立射线追踪法局部直角坐标系，如图3所示，以发射天线所在位置为坐标原点，基于右手螺旋准则建立局部直角坐标系xyz，坐标系x轴为射线的初始方向在地球表面切平面投影，z轴为地球球心到发射机所在位置的方向向量，α0和β0分别为发射天线主轴指向的方位角和俯仰角；
在该局部直角坐标系xyz中电磁波传播路径满足线性微分方程组
式中的n为大气折射指数，n＝10-6N+1，N 为大气折射率，
P、T、q为通过步骤S1从数值气象模型中获取指定位置的气压、温度和比湿度，k1、k2、k3为大气折射系数，k1＝77.6890×10-2k/hPa，k2＝71.2952×10-2k/hPa，k3＝375.463×103 k2/hPa，Rd为干空气气体常数，Rv为水蒸气气体常数，Rv＝461.525J/(kg·K)，e 为水汽压；
上述线性微分方程组在初始边界条件
约束下以及对流层顶高(86km)限制下，基于数值气象模型的气象数据，通过数值法迭代求解即可获得电磁波从发射天线到对流层顶的传播路径，也就得到发射天线主轴方向的电磁波传播路径。对于接收天线，以接收天线作为起始点采用相同的方法计算即可。
S3、对发射波束进行横纵方向的剖分及路径旋转，得到各发射子波束的射线路径，对接收波束进行横向剖分、路径旋转及射线对准，得到有效范围内的接收子波束的射线路径，从得到的发射子波束的射线路径和有效范围内的接收子波束的射线路径中截取从发射天线到散射点再到接收天线的路径，即为散射子波束对流层散射传播路径；
散射波束包括发射波束和接收波束，对流层散射是典型的多径传输信道，单一的天线主轴射线不能反应信道的真实情况，本发明采用波束剖分法，将散射过程分解为多子波束传播过程，通过计算各发射子波束的接收功率和传播时间来确定对流层散射传输损耗和传播延迟；
在不考虑大气折射情况下，假定收发天线发出的是锥状波束(波束垂直宽度和水平宽度相等)，天线的3dB波束宽度是锥形波束的边界；如图4所示，对于发射波束采用矩形截面剖分法，分别沿发射波束截面(波束传播方向横截面)的横轴和纵轴以等角度间隔dw将发射波束剖分为K等份，将发射波束分解为多个紧密排列且互不重叠的正四棱锥形的发射子波束，每一个正四棱锥形子波束的中心对应一条射线，发射子波束射线对应的仰角β0_ij和方位角α0_ij为
α0和β0分别为发射天线主轴指向的方位角和俯仰角，i,j 为波束横向和纵向剖分的序号，i,j＝-K/2,-K/2+1,…,K/2-1,K/2；
由于散射波束宽度通常较窄，在波束宽度范围内各发射/接收子波束沿传播方向的大气折射率变化基本相同，可认为各发射/接收子波束弯曲程度也相同，因此采用旋转发射/接收天线主轴方向的电磁波传播路径的方式来获取各发射/ 接收子波束的射线路径；对于发射子波束，即将步骤S2得到的发射天线主轴方向的电磁波传播路径ltransmit_ray_main从局部直角坐标系xyz转换为球坐标系 BLH，坐标系转换的方法为本领域的公知技术，不是本发明的创新之处，在此不做赘述，然后将球坐标系BLH下发射天线主轴方向的电磁波传播路径的初始指向旋转到各发射子波束射线的初始指向即得到各发射子波束的射线路径：序号下标为ij的发射子波束的射线路径 ltransmit_ray_ij＝ltransmit_ray_main(B,L,H)+(dw·i,dw·j,0)，ltransmit_ray_main(B,L,H)为球坐标系BLH下的发射天线主轴方向的电磁波传播路径；
根据湍流非相干散射理论，接收波束能量来源于公共散射体中被发射波束照射而激发的偶极子的二次辐射，因此每个发射子波束在与接收天线波束的公共重叠交汇区都将激发二次辐射，即每条发射子波束对应多个接收子波束，对于接收波束，其剖分方式与发射波束的剖分方式不同，接收子波束俯仰角在接收波束宽度内等角度间隔分布，方位角使用射线对准算法调整，接收子波束具体通过横向剖分、路径旋转和射线对准三个过程构建：
横向剖分：如图5所示，在接收波束截面以等角度间隔dw将接收波束横向剖分为K等份，得到方位角相同而俯仰角不同的若干接收子波束，每一接收子波束的中心同样对应一条射线，接收子波束射线对应的方位角均为α0'，俯仰角β0_ij'为β0_ij'＝β0'+dw·i，β0'为接收天线主轴指向的俯仰角，i为波束横向剖分的序号，i＝-K/2,-K/2+1,…,K/2-1,K/2；
路径旋转：通过旋转接收天线主轴方向的电磁波传播路径获取接收子波束的射线路径，即将步骤S2得到的接收天线主轴方向的电磁波传播路径lr_ray_main从局部直角坐标系xyz转换为球坐标系BLH，然后将球坐标系BLH下接收天线主轴方向的电磁波传播路径的初始指向旋转到各接收子波束的初始指向即得到各接收子波束的射线路径：序号下标为ij的接收子波束的射线路径 lr_ray_i＝lr_ray_main(B,L,H)+(dw·i,0,0)，lr_ray_main(B,L,H)为球坐标系BLH下的接收天线主轴方向的电磁波传播路径；
射线对准：通常情况下通过路径旋转获得发射子波束射线与对应的接收子波束射线是三维空间中互不相交的异面曲线，需要调整接收子波束的方位角以实现收发射线对准。接收子波束的方位角与激发该接收子波束的发射子波束相关，如图6所示，针对每一条发射子波束，使用射线对准法来确定对应各接收子波束的方位角；具体的，对于每一条发射子波束，计算接收子波束射线与该发射子波束射线之间的最小距离点pt_min、pr_min以及最小距离点之间的最小距离dmin，获取pt_min和pr_min关于接收天线发射点的夹角Δα，对于每一条接收子波束射线，按照夹角Δα大小旋转该接收子波束射线的方位角，直到最小距离 dmin小于设定阈值(阈值设定为米级，本实施例的设定阈值为1m)，计算此时 (dmin小于设定阈值时)该接收子波束射线的初始方位角与接收天线主轴的方位角夹角，将方位角夹角超过波束宽度的接收子波束射线舍弃，得到有效范围内的接收子波束的射线路径，从得到的发射子波束的射线路径和有效范围内的接收子波束的射线路径中截取从发射点(发射天线)到散射点(在dmin小于设定阈值时，认定pt_min和pr_min即为散射点)再到接收点(接收天线)的路径，即为散射子波束对流层散射传播路径。发射射线和接收射线为空间异面曲线，两者不相交于一点，因此两者间的最小距离点有两个，一个是在发射射线上的最小距离点pt_min，另一个是在接收射线上的最小距离点pr_min，dmin为这两个点之间的距离，Δα是点pt_min及点pr_min分别和接收天线发射点的连线间的夹角。在接收射线旋转Δα后，dmin值小于设定阈值，其长度可忽略不计，pt_min和 pr_min在空间位置几乎重合，可认定为散射点。
S4、计算各子波束公共散射体的体积；子波束公共散射体就是一条发射子波束和一条接收子波束相交的部分，即发射子波束和接收子波束共有的部分，通过步骤S3得到发射子波束和接收子波束的射线路径，发射子波束为正四棱锥形，正四棱锥两正对侧面夹角已知为dw，接收子波束的俯仰方向边界夹角已知为dw，散射子波束公共散射体近似为图7中所示的六面体，其体积可通过数学几何方法进行计算；
电磁波在大气中受折射率变化影响，传播路径向地球方向发生弯曲，其作用效果等同于将收发天线抬高，为准确计算公共散射体体积，本实施例通过设定虚拟发射点和虚拟接收点计算子波束公共散射体的体积；如图8所示，在通过步骤S3得到的散射子波束对流层散射传播路径的基础上，在散射点S(即发射子波束射线和接收子波束射线的交点)、发射点T和接收点R构成的平面内，以散射点S为起点，分别向发射点T、接收点R方向作子波束射线路径切向量，切向量长度分别设定为从发射点T到散射点S的实际射线路径长度以及从散射点S到接收点R的实际射线路径长度，切向量的终点T'和R'即为虚拟发射点和虚拟接收点；
在虚拟发射点T'、散射点S和虚拟接收点R'构成的平面内，以虚拟发射点 T'为原点o'，虚拟发射点T'到散射点S方向为x'轴，地球球心到虚拟发射点T' 方向为y'轴，基于右手螺旋准则构建三维直角坐标系x'y'z'，在该三维直角坐标系x'y'z'中，T'S和SR'即为电磁波的虚拟直线传播路径，将T'S和SR'分别绕原点o'和虚拟接收点R'在x'o'y'平面中旋转0.5dw和-0.5dw，即获得发射子波束和接收子波束的上下边界，进一步获得发射子波束和接收子波束的上下边界交点P1～P4在x'y'z'中的坐标；发射子波束的波束夹角为dw，通过P1～P4在 x'y'z'中的坐标即可获得子波束公共散射体(六面体)各顶点坐标；如图9所示，将子波束公共散射体沿x'轴方向(发射波束传播方向)细分为长度为1km 的梯形体，计算各梯形体的体积，所有梯形体的体积之和即为子波束公共散射体的体积。
S5、计算散射子波束的接收功率和传播时间；
本实施例采用双基地雷达方程计算散射子波束的接收功率；首先分别计算各梯形体二次激发产生的接收功率，散射子波束中第l个梯形体的接收功率为：
其中，Pt_ij为序号下标为ij的散射子波束发射功率，Gt、Gr分别为发射天线增益和接收天线增益，gt、gr分别为发射天线的方向性函数和接收天线的方向性函数，σij_l为第l个梯形体的散射截面，ξ为电磁波波长，dVij_l为第l个梯形体的体积，Rij_l为第l个梯形体的中心点到接收天线的距离，Sij_l为第l个梯形体的中心点到发射天线的距离；
收发天线增益、方向性函数、电磁波波长等与气象环境无关，通过散射链路基本参数直接获取；梯形体的体积、梯形体的中心点到收发天线的距离与气象环境间接相关，通过基于数值气象模型获取的子波束射线路径获取；散射截面与气象环境密切相关，其中与气象环境相关的变量M、u、v和dT/dh均通过数值气象模型获取；更具体的，发射天线增益和接收天线增益可根据天线直径 D和电磁波波长ξ计算，即天线增益G＝10lg(4.5(D/ξ)2)；
天线的方向性函数为高斯型：
其中的ψ1、ψ2分别为发射子波束及接收子波束的3dB宽度，
β1、β2分别为发射天线的仰角和接收天线的仰角，α1、α2分别为发射天线的方位角和接收天线的方位角，β10、β20分别为发射天线的主轴仰角和接收天线的主轴仰角，α10、α20分别为发射天线的主轴方位角和接收天线的主轴方位角；
梯形体的散射截面σij_l根据柯尔莫哥洛夫理论和科尔莫哥洛夫-奥布霍夫定律计算，σij_l＝2πk4cos(Θij)2Φ(k)，其中的k为波数，k＝2π/ξ，Θij为通过散射子波束对流层散射传播路径确定的散射角，Φ(k)为空间谱函数，
为折射率结构常数，
M为折射率梯度，L0为湍流外尺度，可利用HMNSP99外尺度模式计算
dT/dh为温度梯度，u、v是水平风东向分量及水平风北向分量，h为海拔高度，M、u、v和dT/dh均通过步骤S1从数值气象模型中获取；
序号下标为ij的散射子波束的接收功率Pr_ij为散射子波束中全部梯形体的接收功率之和，Pr_ij＝∑Pr_ij_l；
序号下标为ij的散射子波束的传播时间τij采用折射率积分法计算：
式中的St为从发射点到散射点的射线路径，Sr为从散射点到接收点的射线路径，n为射线路径上的大气折射指数，可通过数值气象模型获取，c为光速，ds为沿射线路径的曲线微元。
S6、计算对流层散射传输损耗与传播延迟；
对流层散射传输损耗为：
其中，Pr_ij为序号下标为ij的散射子波束的接收功率，Pt为发射天线的发射功率；
传播延迟通过时延功率谱表示，根据步骤S5得到各散射子波束的传播时间以及各散射子波束的接收功率，将所有散射子波束的接收功率归一化，并按升序排列所有散射子波束的传播时间，得到对流层散射链路的归一化延迟功率谱，即为流层散射传播延迟。
大气衰减会对接收功率产生影响，优选的，本实施例基于子波束的传输损耗Lij对散射子波束的接收功率进行修正，修正后的散射子波束的接收功率为：
则对流层散射传输损耗为：
计算流层散射传播延迟时，各子波束的接收功率也优选采用修正后的散射子波束的接收功率。
序号下标为ij的散射子波束的传输损耗
Pr_ij为序号下标为ij的散射子波束的接收功率，Pt_ij为序号下标为ij的散射子波束的发射功率，
为序号下标为ij的散射子波束的大气吸收损耗；考虑到电磁波在大气中传播时，受大气中干空气和水汽影响产生大气吸收损耗，散射子波束的大气吸收损耗
采用ITU-R P.676-12建议书中的沿路径积分逐线求和法计算：
其中，ij_length为散射子波束对流层散射传播路径的分段长度(即路径数组长度)，aij_b为第b段路径的长度，γij_b为第b段路径的大气衰减比，γij_b与传播路径上的温度、大气压、比湿度相关，相应气象数据通过步骤S1获取，具体计算方法参考ITU-R P.676-12建议书，此处不再赘述。
为了验证本发明方法的传输损耗及传输延迟的计算结果，分别采用本发明方法和现有方法来计算传输损耗和传输延迟。
一、传输损耗的计算及分析
采用本发明方法和ITU-R P.617-2、ITU-R P.617-5两种方法进行传输损耗对比计算。本发明对比计算中使用的数值气象模型为欧洲中期天气预报中心 (European CentreFor Medium Range Weather Forecasts，ECMWF)提供的ERA5 气压分层产品，ERA5时间分辨率1h，空间水平分辨率0.25°，竖直方向按照等压面从1hPa到1000hPa分为37层，相应数据存储在网格模型中。
对比使用的10条对流层散射链路数据从CCIR1986年发布OT/TRER 16报告“对流层传输损失的测量和预测长期分布”中选择，10条对流层散射通信链路的参数分别如表1和表2所示。
表1
表2
ITU方法使用地平线角度，而不是实际的天线仰角，并且未体现波束宽度对传输损耗的影响。本发明方法使用实际的天线仰角和波束宽度计算损耗。在链路建设中，对流散射一般采用低仰角来降低损耗，在通信时会将天线调整到损耗最低的最佳仰角。然而，最佳仰角随气象环境的变化而变化，因此本发明方法根据2020年1月1日的ERA5气象数据，采用迭代计算出传输损耗最小的天线仰角，天线的3dB波束宽度是通过天线增益获得的。十条链路的天线增益是由中国无线电波传播研究所提供，仰角和总增益大小如表1中所示。
表3给出了十条链路相对于自由空间的传输损耗的测量中值和三种方法的传输损耗计算结果。本发明的气象数据是ERA5在链路观测时间范围内的月平均值。
表3三种方法计算结果与观测值对比
从表3的结果可以看出来，本发明方法的计算结果在平均误差和均方根误差方面优于两个ITU模型。
传输损耗与仰角密切相关。图10和图11分别展示了在不同仰角下LT2063 和LT2307两个链路的三种方法的计算结果。从图10和图11可以看出，传输损耗随着仰角的增加而增加，在低仰角时，不同方法之间的差异很小；随着仰角的增加，差异趋向于更加明显。ITU-R P.617-2模型在高仰角时给出的结果比本发明方法和ITU-R P.617-5大得多。本发明方法和ITU-R P.617-5之间的差异也随着仰角的增加而增加，但差异的大小比ITU-RP.617-2和本发明方法之间的差异小。
为体现本发明方法基于气象数据的小时级损耗分析能力，选择LT1441和 LT2305链路基于ERA5小时级气象数据分析对流散射传输损耗的变化特征。两条链路的位置分别如图12和图13所示，图12是位于加拿大纽芬兰岛的 LT1441链路，长277公里，散射链路部分穿过海湾；图13是位于日本东京和福岛之间的LT2305链路，长226公里，散射链路穿过一片陆地。图14给出了 2020年8月25日、26日和27日使用本发明计算的逐小时传输损耗。由图可知，本发明方法通过ERA5的气象数据成功获得了逐小时的传输损耗变化。两条链路都呈现出昼夜交替的损耗变化趋势，白天损耗增加，夜间损耗减少，但昼夜上升和下降的幅度并不一致。LT1441链路紧邻大西洋，穿过海湾，气象环境复杂多变。传输损耗在不同时期变化很大。最大小时变化可达6.3dB，最大昼夜损耗波动为20.4dB；LT2305链路经过陆域，气象条件相对稳定，损耗变化轻微，最大小时变化为2.4dB，最大昼夜损耗波动为5.3dB，损耗变化呈周期性规律。
二、传播延迟的计算及分析
将本发明方法的传播延迟计算能力与基于二维平面和电磁波直线光速传播假设的Bello模型进行对比分析。为了分析对比两者性能，设置了六个不同距离和波束宽度的链路，链路的参数见表4。图15a至图15f给出了两种方法计算出的归一化延迟功率谱。
表4
图5归一化延迟功率谱
从图15a至图15f可以看出，本发明方法和Bello方法的延迟功率分布基本相同，波束宽度相同的情况下通信距离越远，时延越大，时延扩展越宽；对比图15c和图15d，通信距离相同，波束越宽，时延扩展越大。但同时注意到，本发明方法和Bello方法的延迟谱分布存在一定偏移，即两者计算的散射传播延迟大小不同，本发明方法的传播时间均大于Bello方法，这是因为Bello方法是基于散射链路几何构形推导的传播延迟，在计算过程中设定电波沿直线传播，未考虑大气环境对电波传播的影响；而本发明方法通过引入气象参数，在散射波束路径生成和延迟统计中计算了大气环境引起的电波传播路径弯曲和延迟，因此传播延迟均大于Bello方法，显然本发明的计算结果比Bello方法更加准确，更适应于需要精确估计传输时延的情况。
气象环境变化会影响散射链路时延分布，选择1441链路和2305链路的在 2020年8月25日、26日和27日的ERA5气象数据计算传播时延的逐小时变化，为便于对比分析，两条链路除收发站点位置和海拔不同外，其余链路参数均相同：频率4.09GHz、收发天线高度5.2m、收发天线水平角20mrad、收发总增益97dB，其余链路参数参考表1和表2。
图16a和图16b出了采用本发明方法计算的两条链路逐小时归一化延迟功率谱热力图，由图可知，本发明方法提取到了对流层散射链路的逐小时时延变化特征，不同时刻变化的气象环境使散射链路的传播时间和群时延存在差异。 1441链路途径诺特达姆海湾，受海洋气候影响，气象环境复杂多变，散射时延规律性较弱，随机性较强，随气象参数变化而变化，如图16a中25日19时至 23时，时延功率谱在短时间内出现较大幅度突变，主要缘于该散射链路部分低仰角子波束的散射公共体所在位置垂直方向的温度梯度突变，由大约-5e-3℃/m 变化到-3e-3℃/m左右，进而引起子波束接收功率减小，而低仰角子波束传播时间小于高仰角子波束，在归一化延迟功率谱中则表现为最大接收功率波峰向传播时延较大区域移动。2305链路途径区间均为陆地，气象环境相对稳定，散射时延呈现较为明显的昼夜变化规律，白天传播时间明显高于夜间，昼夜传播时间变化幅度约为数十纳秒，因此对于需要精确估计传播时延和群延迟的时间同步系统和无源探测系统，气象环境带来的影响不容忽视。
以上所述，仅是本发明的较佳实施例而已，并非对本发明做任何形式上的限制，虽然本发明已以较佳实施例揭露如上，然而并非用以限定本发明，任何熟悉本专业的技术人员，在不脱离本发明技术方案范围内，当可利用上述揭示的技术内容做出些许更动或修饰为等同变化的等效实施例，但凡是未脱离本发明技术方案的内容，依据本发明的技术实质对以上实施例所作的任何简单修改、等同变化与修饰，均仍属于本发明技术方案的范围内。
Claims (10)
1.一种对流层散射传输损耗及传播延迟计算方法，其特征在于，包括以下步骤：
S1、通过数值气象模型获取散射传播路径上的气象数据；所述数值气象模型以网格形式存储气象数据，水平方向网格间隔为经纬度分辨率，竖直方向网格为等压面，基于数值气象模型采用插值法或外推法在获取指定位置的气象数据，所述气象数据包括温度、比湿度、气压、风速；
S2、基于通过数值气象模型获取的气象数据，计算发射天线主轴方向的电磁波传播路径和接收天线主轴方向的电磁波传播路径；
S3、对发射波束进行横纵方向的剖分及路径旋转，得到各发射子波束的射线路径，对接收波束进行横向剖分、路径旋转及射线对准，得到有效范围内的接收子波束的射线路径，从得到的发射子波束的射线路径和有效范围内的接收子波束的射线路径中截取从发射天线到散射点再到接收天线的路径，即为散射子波束对流层散射传播路径；
S4、计算各子波束公共散射体的体积，所述子波束公共散射体为一条发射子波束和一条接收子波束相交的部分；
S5、计算散射子波束的接收功率及散射子波束的传播时间；
S6、计算对流层散射传输损耗与传播延迟；
对流层散射传输损耗为：
其中，Pr_ij为序号下标为ij的散射子波束的接收功率，Pt为发射天线的发射功率；
将所有散射子波束的接收功率归一化，并按升序排列所有散射子波束的传播时间，得到对流层散射链路的归一化延迟功率谱，即为对流层散射传播延迟。
2.根据权利要求1所述的对流层散射传输损耗及传播延迟计算方法，其特征在于：所述步骤S1中基于数值气象模型获取指定位置的气象数据的步骤具体如下：
S101、将指定位置的大地高Hg转换为位势高
其中，
为指定位置的纬度，
为地球椭球面法向重力，γ45°为纬度45°的重力加速度的标准值，，
为局部地球半径；
S102、通过比较插值点与数值气象模型的等压面的位势高相对位置关系确定相应的插值方法和外推方法，包括以下三种情况：
(1)插值点位势高在等压面内部时，通过插值点的经纬度和位势高确定数值气象模型中与插值点邻近的8个网格点，分别对具有相同经纬度的4组上下相对的两个网格点向插值点位势高处插值：
竖直方向上，温度、比湿度和风速采用线性插值，气压采用指数模型插值；
水平方向上，温度、比湿度、风速、气压均采用双线性差值；
(2)插值点位势高在等压面最底层以下时，通过插值点的经纬度和位势高确定与插值点邻近的8个网格点，对于8个网格点中的上层等压面的4个网格点：
水平方向上，使用双线性插值法获取插值点在等压面上对应位置处的气象参数；
竖直方向上，分别采用以下方法外推计算插值点在等压面上对应位置处的气象参数：气压使用指数模型，温度Tj＝Tv_j-0.0065hg_j，Tv_j为在等压面上通过水平双线性插值获取的温度，比湿度qj＝qv_j，qv_j为在等压面上通过水平双线性插值获取的比湿度，风速使用幂定律插值；
(3)插值点位势高位于等压面最高层以上时，通过插值点的经纬度和位势高确定与插值点邻近的8个网格点，对于8个网格点中的下层等压面的4个网格点：
水平方向上，使用双线性插值法获取插值点在等压面上对应位置处的气象参数；
竖直方向上，分别采用以下方法外推计算插值点在等压面上对应位置处的气象参数：比湿度为零；气压在水平双线性插值的基础上使用指数模型外推；温度和风速分量使用CIRA86国际参考大气模型获取。
3.根据权利要求1所述的对流层散射传输损耗及传播延迟计算方法，其特征在于：所述步骤S2中，采用射线追踪法计算发射天线主轴方向的电磁波传播路径，步骤如下：
S201、以发射天线所在位置为坐标原点，基于右手螺旋准则建立局部直角坐标系xyz，坐标系x轴为射线的初始方向在地球表面切平面投影，z轴为地球球心到发射机所在位置的方向向量，α0和β0分别为发射天线主轴指向的方位角和俯仰角；
在该局部直角坐标系xyz中电磁波传播路径满足线性微分方程组
式中的n为大气折射指数，N为大气折射率；
上述线性微分方程组在初始边界条件
约束下以及对流层顶高限制下，基于数值气象模型得到的气象数据，通过数值法迭代求解即可获得发射天线主轴方向的电磁波传播路径；
对于接收天线，以接收天线作为起始点采用相同步骤计算，得到天线主轴方向的电磁波传播路径。
4.根据权利要求1所述的对流层散射传输损耗及传播延迟计算方法，其特征在于：所述步骤S3中，采用矩形截面剖分法对发射波束进行剖分，分别沿发射波束截面的横轴和纵轴以等角度间隔dw将发射波束剖分为K等份，将发射波束分解为多个紧密排列且互不重叠的正四棱锥形的发射子波束；
将发射天线主轴方向的电磁波传播路径ltransmit_ray_main从局部直角坐标系xyz转换为球坐标系BLH，然后将球坐标系BLH下发射天线主轴方向的电磁波传播路径的初始指向旋转到各发射子波束射线的初始指向，序号下标为ij的发射子波束的射线路径ltransmit_ray_ij＝ltransmit_ray_main(B,L,H)+(dw·i,dw·j,0)，ltransmit_ray_main(B,L,H)为球坐标系BLH下的发射天线主轴方向的电磁波传播路径，i,j＝-K/2,-K/2+1,…,K/2-1,K/2。
5.根据权利要求1所述的对流层散射传输损耗及传播延迟计算方法，其特征在于：所述步骤S3中，采用横向剖分法对接收波束进行剖分，在接收波束截面以等角度间隔dw将接收波束横向剖分为K等份，得到方位角相同而俯仰角不同的若干接收子波束；将发射天线主轴方向的电磁波传播路径ltransmit_ray_main从局部直角坐标系xyz转换为球坐标系BLH，然后将球坐标系BLH下发射天线主轴方向的电磁波传播路径的初始指向旋转到各发射子波束射线的初始指向，序号下标为ij的发射子波束的射线路径ltransmit_ray_ij＝ltransmit_ray_main(B,L,H)+(dw·i,dw·j,0)，ltransmit_ray_main(B,L,H)为球坐标系BLH下的发射天线主轴方向的电磁波传播路径，i＝-K/2,-K/2+1,…,K/2-1,K/2；对于每一条发射子波束，计算接收子波束射线与该发射子波束射线之间的最小距离点以及这两个最小距离点之间的最小距离，获取这两个最小距离点关于接收天线发射点的夹角，按照夹角大小旋转接收子波束射线的方位角，直到两个最小距离点之间的最小距离小于设定阈值，计算此时接收子波束射线的初始方位角与接收天线主轴的方位角夹角，将方位角夹角超过波束宽度的接收子波束射线舍弃，得到有效范围内的接收子波束的射线路径。
6.根据权利要求1所述的对流层散射传输损耗及传播延迟计算方法，其特征在于：所述步骤S4中，通过以下步骤计算子波束公共散射体的体积：
在散射点S、发射点T和接收点R构成的平面内，以散射点S为起点，分别向发射点T、接收点R方向作子波束射线路径切向量，切向量长度分别设定为从发射点T到散射点S的实际射线路径长度以及从散射点S到接收点R的实际射线路径长度，切向量的终点T'和R'即为虚拟发射点和虚拟接收点；
在虚拟发射点T'、散射点S和虚拟接收点R'构成的平面内，以虚拟发射点T'为原点o'，虚拟发射点T'到散射点S方向为x'轴，地球球心到虚拟发射点T'方向为y'轴，基于右手螺旋准则构建三维直角坐标系x'y'z'，在该三维直角坐标系x'y'z'中，T'S和SR'即为电磁波的虚拟直线传播路径，将T'S和SR'分别绕原点o'和虚拟接收点R'在x'o'y'平面中旋转0.5dw和-0.5dw，获得发射子波束和接收子波束的上下边界，进一步获得发射子波束和接收子波束的上下边界交点P1～P4在x'y'z'中的坐标；通过P1～P4在x'y'z'中的坐标可获得子波束公共散射体各顶点坐标；将子波束公共散射体沿x'轴方向细分若干个梯形体，计算各梯形体的体积，所有梯形体的体积之和即为子波束公共散射体的体积。
7.根据权利要求1所述的对流层散射传输损耗及传播延迟计算方法，其特征在于：所述步骤S5中，采用双基地雷达方程计算散射子波束的接收功率，步骤如下：
计算各梯形体的接收功率，散射子波束中第l个梯形体的接收功率为：
其中，Pt_ij为序号下标为ij的散射子波束发射功率，Gt、Gr分别为发射天线增益和接收天线增益，gt、gr分别为发射天线的方向性函数和接收天线的方向性函数，σij_l为第l个梯形体的散射截面，ξ为电磁波波长，dVij_l为第l个梯形体的体积，Rij_l为第l个梯形体的中心点到接收天线的距离，Sij_l为第l个梯形体的中心点到发射天线的距离；
序号下标为ij的散射子波束的接收功率Pr_ij为该散射子波束中全部梯形体的接收功率之和，Pr_ij＝∑Pr_ij_l。
8.根据权利要求1所述的对流层散射传输损耗及传播延迟计算方法，其特征在于：所述步骤S5中，采用折射率积分法计算散射子波束的传播时间τij，
式中的St为从发射点到散射点的射线路径，Sr为从散射点到接收点的射线路径，n为射线路径上的大气折射指数，c为光速，ds为曲线微元。
9.根据权利要求1所述的对流层散射传输损耗及传播延迟计算方法，其特征在于：所述步骤S6中，先基于子波束的传输损耗Lij对散射子波束的接收功率进行修正，修正后的散射子波束的接收功率为：
式中的Pr_ij为序号下标为ij的散射子波束的接收功率，Pt_ij为序号下标为ij的散射子波束的发射功率，
为序号下标为ij的散射子波束的大气吸收损耗，再根据修正后的散射子波束的接收功率P′r_ij计算对流层散射传输损耗，流层散射传输损耗
10.根据权利要求9所述的对流层散射传输损耗及传播延迟计算方法，其特征在于：根据修正后的散射子波束的接收功率P′r_ij计算对流层散射传输延时。
Priority Applications (1)
| Application Number | Priority Date | Filing Date | Title |
|---|---|---|---|
| CN202210220518.0A CN114726433A (zh) | 2022-03-08 | 2022-03-08 | 一种对流层散射传输损耗及传播延迟计算方法 |
Applications Claiming Priority (1)
| Application Number | Priority Date | Filing Date | Title |
|---|---|---|---|
| CN202210220518.0A CN114726433A (zh) | 2022-03-08 | 2022-03-08 | 一种对流层散射传输损耗及传播延迟计算方法 |
Publications (1)
| Publication Number | Publication Date |
|---|---|
| CN114726433A true CN114726433A (zh) | 2022-07-08 |
Family
ID=82238691
Family Applications (1)
| Application Number | Title | Priority Date | Filing Date |
|---|---|---|---|
| CN202210220518.0A Pending CN114726433A (zh) | 2022-03-08 | 2022-03-08 | 一种对流层散射传输损耗及传播延迟计算方法 |
Country Status (1)
| Country | Link |
|---|---|
| CN (1) | CN114726433A (zh) |
Cited By (2)
| Publication number | Priority date | Publication date | Assignee | Title |
|---|---|---|---|---|
| CN117194856A (zh) * | 2023-11-07 | 2023-12-08 | 天津云遥宇航科技有限公司 | 基于era5和cira的掩星弯曲角优化方法 |
| CN118013172A (zh) * | 2024-04-07 | 2024-05-10 | 山东科技大学 | 基于梯形公式求解积分方程的低频地波传播特性预测方法 |
-
2022
- 2022-03-08 CN CN202210220518.0A patent/CN114726433A/zh active Pending
Cited By (3)
| Publication number | Priority date | Publication date | Assignee | Title |
|---|---|---|---|---|
| CN117194856A (zh) * | 2023-11-07 | 2023-12-08 | 天津云遥宇航科技有限公司 | 基于era5和cira的掩星弯曲角优化方法 |
| CN117194856B (zh) * | 2023-11-07 | 2024-01-26 | 天津云遥宇航科技有限公司 | 基于era5和cira的掩星弯曲角优化方法 |
| CN118013172A (zh) * | 2024-04-07 | 2024-05-10 | 山东科技大学 | 基于梯形公式求解积分方程的低频地波传播特性预测方法 |
Similar Documents
| Publication | Publication Date | Title |
|---|---|---|
| Hitney et al. | Tropospheric radio propagation assessment | |
| Guan et al. | Deterministic propagation modeling for the realistic high-speed railway environment | |
| Zhang et al. | Applying the parabolic equation to tropospheric groundwave propagation: A review of recent achievements and significant milestones | |
| CN114726433A (zh) | 一种对流层散射传输损耗及传播延迟计算方法 | |
| Matolak et al. | Initial results for air-ground channel measurements & modeling for unmanned aircraft systems: Over-sea | |
| KR102175179B1 (ko) | 인공지능 기법을 활용한 공백 지역 레이더 합성장 정보 생성 시스템 및 방법 | |
| CN115236615B (zh) | 基于t矩阵法的机载极化气象雷达降水粒子回波仿真方法 | |
| Pederick et al. | A directional HF noise model: Calibration and validation in the Australian region | |
| CN113391315A (zh) | 基于电磁波抛物方程伴随模式的雷达回波资料反演大气波导的方法 | |
| CN112671453A (zh) | 低轨星座的电磁态势建模方法和装置 | |
| CN112383370A (zh) | 一种卫星时变信道的建模仿真方法及系统 | |
| Doebler et al. | Five-year simulation study of NASA's X-59 low-boom carpets across the contiguous United States of America | |
| Nilsson | Radio-wave propagation modelling over rough sea surfaces and inhomogeneous atmosphere | |
| Lindquist | Wave propagation models in the troposphere for long-range UHF/SHF radio connections | |
| Juaeni et al. | Retrieval of temperature profiles using radio acoustic sounding system (RASS) with the equatorial atmosphere radar (EAR) in West Sumatra, Indonesia | |
| Zhang et al. | A new model for estimating troposcatter loss and delays based on ray-tracing and beam splitting with ERA5 | |
| Grabner et al. | Clear-air propagation modeling using parabolic equation method | |
| CN115238737A (zh) | 一种基于合作式航海雷达信号的蒸发波导剖面反演方法 | |
| CN112149302B (zh) | 一种基于非均匀下垫面的台风建模方法及相关装置 | |
| Etokebe et al. | Determination of atmospheric effective earth radius factor (k-factor) under clear air in Lagos, Nigeria | |
| Durgin | Advanced site-specific propagation prediction techniques | |
| Sirkova | Anomalous tropospheric propagation: Usage possibilities and limitations in radar and wireless communications systems | |
| Zhang et al. | A New Model for Predicting Troposcatter Transmission Loss and Propagation Delay Based on GPT2w | |
| Grishin et al. | Use of ITU-R recommendations in calculating tropospheric signal attenuation in the simulation modeling problems of satellite systems | |
| Asiyo | Characterization and modelling of effects of clear air on multipath fading in terrestrial links |
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 |