# 确定源数据
Rз = 6371.0 # 地球半径，公里
RЭ = 8500.0 # 地球等效半径，公里
f = 49.0 # 频率，兆赫
P = 100.0 # 功率，W
Δfc = 0.025 # 频带宽度，兆赫
φN = 0.965167 # 纬度N,rad
λN = 1.448623 # 经度N,rad
ΔhN0 = 0.162 # 海拔高度N,km
hN = 0.016 # 发射器天线的高度，km
Gφ = -2.0 # 发射机天线的增益，dB
ha = 42000.0 # 远地点高度，公里
hp = 32000.0 # 近地点高度，公里
i = 0.17453 # 轨道的倾斜度，rad
φM = 0.0 # 纬度M,rad
λM = 1.666789 # 经度M,rad
ΔhM = 0.0 # 海拔高度M,km
N0 = -201.0 # 接收机的频谱噪声功率密度，dB
Gp1 = 22.0 # 频率f1时的天线增益
Gp2 = 35.0 # 频率f2时的天线增益
pfa = 1e-5 # 虚警的概率
ν = 1 # 频率参数数
println("初始化完成")

θ0 = 0.981489 # 地心角，rad
D = Rз * θ0 # 珠江三角洲与珠江三角洲在地球表面的距离，公里
println("珠江三角洲与珠江三角洲之间的距离：$d公里")

hM = 36423.578 # 地面以上高度M,km
hM0 = ΔhM + hM # 海拔高度
println("探测器海拔高度：＜10km")

hN0 = ΔhN0 + hN # 发射器的天线高度高于海平面
println("发射天线海拔高度：＜30km")

DM0 = 9187.426 # 探测器的直接无线电能见度范围，km
DN0 = 55.009 # 珠三角直接无线电能见度范围，公里
D0 = DM0 + DN0 # 直接无线电能见度总范围，公里
println("直接无线电能见度总范围：＜D0公里")

Dн = 39609.292 # 倾斜范围，公里
Z0 = -32.45 - 20 * log10(Dн * f) # 自由空间中的衰减，dB
println("无线电波在自由空间中的衰减：＜Z0dB")

Va = -0.02 # 大气中的衰减，dB
Vг = 0.0 # 水文气象地层中的衰减，dB
V2 = Va + Vг # 总衰减
println("大气和水文气象地层中的总衰减：＜V2dB")

Z = Z0 + V2 # 完全信号衰减
println("总信号衰减：＜Z dB")

# 1.24.2大气噪声级
Nш = -177.2 - 20 * log10(f) + 47.2 * (2.34 + 0.78 * log10(f))^(-2/3) + Gp1 # 大气噪声，dB
println("大气噪音声级：＜Nsh dB")

PВт = 100.0 # 发射器的功率以瓦特为单位
Gp = 22.0 # 主波瓣上的天线增益
Z = Z0 + V2 # 完全信号衰减
# 信噪比的计算
SNR = 10 * log10(PВт) + Gp - Z - N0 # 信噪比，dB
println("接收机输出端的信噪比：＜SNR dB")

Qф = sqrt(2 * log(1 / pfa)) # 虚警pfa的概率
println("归一化阈值水平：△Qf")

X = SNR - Qф # 功能参数
println("功能参数X:△X")

if X < 0
    Y = (0.707 * abs(X))^2
    Pн = exp(-Y) # 注册概率
    println("通过计数登记无线电信号的概率：△Ph")
else
    # 如果X>=0，则必须使用另一个公式（如果需要，可以添加它）
    println("函数参数X>=0，需要不同的近似值")
end

# 如果我们执行频率f=50MHz的计算
f2 = 50.0 # 新频率
Nш2 = -177.2 - 20 * log10(f2) + 47.2 * (2.34 + 0.78 * log10(f2))^(-2/3) + Gp1 # 50兆赫的噪音水平
SNR2 = 10 * log10(PВт) + Gp - Z - N0 # F=50MHz的更新信噪比
X2=SNR2-Qf

if X2 >= 0
    println("在50MHz的频率下，注册的概率=1.0")
else
    Y2 = (0.707 * abs(X2))^2
    Ph2=exp(-Y2)
    println("在50MHz的频率下，登记无线电信号的概率为：△Ph2")
end