衛星回線設計 — リンクバジェットの計算と設計マージン

import numpy as np
import matplotlib.pyplot as plt

# 物理定数
Re = 6371         # 地球半径 [km]
c = 3e8           # 光速 [m/s]
kB = 1.38e-23     # ボルツマン定数 [W/(Hz*K)]
kB_dB = -228.6    # dB(W/Hz/K)

def antenna_gain(diameter_m, freq_ghz, efficiency=0.6):
    """パラボラアンテナの利得 [dBi]"""
    lam = c / (freq_ghz * 1e9)
    G = efficiency * (np.pi * diameter_m / lam)**2
    return 10 * np.log10(G)

def slant_range(altitude_km, elev_deg):
    """衛星までのスラント距離 [km]"""
    elev = np.radians(elev_deg)
    h = altitude_km
    d = -Re * np.sin(elev) + np.sqrt((Re * np.sin(elev))**2 +
         2 * Re * h + h**2)
    return d

def fspl(distance_km, freq_ghz):
    """自由空間伝搬損失 [dB]"""
    return 20 * np.log10(distance_km * 1e3) + 20 * np.log10(freq_ghz * 1e9) + \
           20 * np.log10(4 * np.pi / c)

def link_budget(eirp_dbw, freq_ghz, distance_km, gt_dbk, bandwidth_hz,
                atm_loss_db=0.3, rain_loss_db=0, point_loss_db=0.3,
                feeder_loss_db=0.5, misc_loss_db=0.5):
    """
    リンクバジェットを計算する
    返り値: dict with all parameters
    """
    # 自由空間損失
    L_fspl = fspl(distance_km, freq_ghz)

    # 全損失
    L_total = L_fspl + atm_loss_db + rain_loss_db + point_loss_db + \
              feeder_loss_db + misc_loss_db

    # C/N計算
    cn = eirp_dbw + kB_dB - L_total + gt_dbk - 10 * np.log10(bandwidth_hz)

    # C/N0 (帯域幅正規化)
    cn0 = cn + 10 * np.log10(bandwidth_hz)

    return {
        'EIRP': eirp_dbw,
        'FSPL': L_fspl,
        'Atm_loss': atm_loss_db,
        'Rain_loss': rain_loss_db,
        'Point_loss': point_loss_db,
        'Feeder_loss': feeder_loss_db,
        'Misc_loss': misc_loss_db,
        'L_total': L_total,
        'G/T': gt_dbk,
        'Bandwidth': bandwidth_hz,
        'C/N': cn,
        'C/N0': cn0,
    }

# === GEO Ku帯ダウンリンクの例 ===
print("=" * 60)
print("GEO Ku-Band Downlink Budget (BSat-like)")
print("=" * 60)

# パラメータ
sat_alt = 35786       # km
elev = 38             # 度（東京からBSat方向）
freq_dl = 12.0        # GHz
sat_tx_power = 20     # dBW (100 W)
sat_ant_gain = 36     # dBi
eirp = sat_tx_power + sat_ant_gain  # = 56 dBW

# 受信側
rx_ant_dia = 0.45     # m
rx_ant_gain = antenna_gain(rx_ant_dia, freq_dl)
sys_noise_temp = 120  # K
gt = rx_ant_gain - 10 * np.log10(sys_noise_temp)

# 距離
dist = slant_range(sat_alt, elev)

# 帯域幅（1トランスポンダ 36 MHz）
bw = 36e6

# 計算
result = link_budget(
    eirp_dbw=eirp,
    freq_ghz=freq_dl,
    distance_km=dist,
    gt_dbk=gt,
    bandwidth_hz=bw,
    atm_loss_db=0.3,
    rain_loss_db=0,     # 晴天時
    point_loss_db=0.2,
    feeder_loss_db=0.3,
    misc_loss_db=0.5
)

print(f"\n--- 送信側 ---")
print(f"  衛星送信電力: {sat_tx_power} dBW ({10**(sat_tx_power/10):.0f} W)")
print(f"  衛星アンテナ利得: {sat_ant_gain} dBi")
print(f"  EIRP: {eirp} dBW")
print(f"\n--- 伝搬 ---")
print(f"  スラント距離: {dist:.0f} km")
print(f"  仰角: {elev}°")
print(f"  FSPL: {result['FSPL']:.1f} dB")
print(f"  大気損失: {result['Atm_loss']:.1f} dB")
print(f"  降雨損失: {result['Rain_loss']:.1f} dB (晴天)")
print(f"  ポインティング損失: {result['Point_loss']:.1f} dB")
print(f"  フィーダー損失: {result['Feeder_loss']:.1f} dB")
print(f"  その他損失: {result['Misc_loss']:.1f} dB")
print(f"  全損失: {result['L_total']:.1f} dB")
print(f"\n--- 受信側 ---")
print(f"  受信アンテナ直径: {rx_ant_dia} m")
print(f"  受信アンテナ利得: {rx_ant_gain:.1f} dBi")
print(f"  システム雑音温度: {sys_noise_temp} K ({10*np.log10(sys_noise_temp):.1f} dBK)")
print(f"  G/T: {gt:.1f} dB/K")
print(f"\n--- 結果 ---")
print(f"  帯域幅: {bw/1e6:.0f} MHz")
print(f"  C/N: {result['C/N']:.1f} dB")
print(f"  C/N0: {result['C/N0']:.1f} dB·Hz")

# DVB-S2の必要C/N（QPSK 3/4で約4.7 dB）
cn_required = 4.7
margin = result['C/N'] - cn_required
print(f"\n  必要C/N (QPSK 3/4): {cn_required} dB")
print(f"  マージン: {margin:.1f} dB")

import numpy as np
import matplotlib.pyplot as plt

# GEO vs LEO の比較
print("\n" + "=" * 70)
print("GEO vs LEO Link Budget Comparison (Ka-Band Downlink)")
print("=" * 70)

configs = {
    'GEO Ka-band': {
        'altitude': 35786, 'elev': 38, 'freq': 20.0,
        'sat_power_dbw': 20, 'sat_gain_dbi': 42,
        'rx_diameter': 0.60, 'sys_temp': 200,
        'atm': 0.5, 'rain': 0, 'point': 0.3,
        'feeder': 0.3, 'misc': 0.5
    },
    'LEO Ka-band\n(Starlink-like)': {
        'altitude': 550, 'elev': 40, 'freq': 20.0,
        'sat_power_dbw': 7, 'sat_gain_dbi': 38,
        'rx_diameter': 0.50, 'sys_temp': 250,
        'atm': 0.5, 'rain': 0, 'point': 0.8,
        'feeder': 0.3, 'misc': 0.5
    }
}

fig, axes = plt.subplots(1, 2, figsize=(14, 7))

for idx, (name, cfg) in enumerate(configs.items()):
    dist = slant_range(cfg['altitude'], cfg['elev'])
    eirp = cfg['sat_power_dbw'] + cfg['sat_gain_dbi']
    rx_gain = antenna_gain(cfg['rx_diameter'], cfg['freq'])
    gt = rx_gain - 10 * np.log10(cfg['sys_temp'])

    result = link_budget(
        eirp_dbw=eirp, freq_ghz=cfg['freq'], distance_km=dist,
        gt_dbk=gt, bandwidth_hz=500e6,
        atm_loss_db=cfg['atm'], rain_loss_db=cfg['rain'],
        point_loss_db=cfg['point'], feeder_loss_db=cfg['feeder'],
        misc_loss_db=cfg['misc']
    )

    # ウォーターフォールチャート（段階的な積み上げ）
    items = [
        ('EIRP', eirp, '#4CAF50'),
        ('FSPL', -result['FSPL'], '#F44336'),
        ('Atm Loss', -cfg['atm'], '#FF9800'),
        ('Rain Loss', -cfg['rain'], '#FFEB3B'),
        ('Point Loss', -cfg['point'], '#FF5722'),
        ('Feeder Loss', -cfg['feeder'], '#795548'),
        ('Misc Loss', -cfg['misc'], '#9E9E9E'),
        ('G/T', gt, '#2196F3'),
        ('-kB', 228.6, '#00BCD4'),
        ('-10log(B)', -10*np.log10(500e6), '#673AB7'),
    ]

    cumulative = 0
    labels = []
    values = []
    colors = []
    bottoms = []

    for label, val, color in items:
        labels.append(label)
        values.append(val)
        colors.append(color)
        bottoms.append(cumulative if val > 0 else cumulative + val)
        cumulative += val

    bars = axes[idx].barh(labels, [abs(v) for v in values],
                            left=bottoms, color=colors, alpha=0.8, height=0.6)

    # C/N値を表示
    axes[idx].axvline(x=cumulative, color='red', linewidth=2, linestyle='--')
    axes[idx].text(cumulative + 1, -0.5, f'C/N = {cumulative:.1f} dB',
                    fontsize=12, color='red', fontweight='bold')

    axes[idx].set_xlabel('Cumulative [dB]', fontsize=12)
    axes[idx].set_title(name.replace('\n', ' '), fontsize=13)
    axes[idx].grid(True, alpha=0.3, axis='x')
    axes[idx].invert_yaxis()

plt.tight_layout()
plt.savefig('link_budget_waterfall.png', dpi=150, bbox_inches='tight')
plt.show()

import numpy as np
import matplotlib.pyplot as plt

# ベースラインパラメータ（GEO Ku帯ダウンリンク）
baseline = {
    'eirp': 56, 'freq': 12.0, 'altitude': 35786, 'elev': 38,
    'rx_dia': 0.45, 'sys_temp': 120, 'bandwidth': 36e6,
    'atm': 0.3, 'rain': 0, 'point': 0.2, 'feeder': 0.3, 'misc': 0.5
}

def compute_cn(params):
    """ベースラインからC/Nを計算"""
    dist = slant_range(params['altitude'], params['elev'])
    rx_gain = antenna_gain(params['rx_dia'], params['freq'])
    gt = rx_gain - 10 * np.log10(params['sys_temp'])
    result = link_budget(
        eirp_dbw=params['eirp'], freq_ghz=params['freq'],
        distance_km=dist, gt_dbk=gt, bandwidth_hz=params['bandwidth'],
        atm_loss_db=params['atm'], rain_loss_db=params['rain'],
        point_loss_db=params['point'], feeder_loss_db=params['feeder'],
        misc_loss_db=params['misc']
    )
    return result['C/N']

cn_baseline = compute_cn(baseline)

fig, axes = plt.subplots(2, 3, figsize=(16, 10))

# 1. 受信アンテナ直径の影響
dias = np.linspace(0.3, 1.5, 50)
cn_dia = []
for d in dias:
    p = baseline.copy()
    p['rx_dia'] = d
    cn_dia.append(compute_cn(p))

axes[0, 0].plot(dias * 100, cn_dia, 'b-', linewidth=2)
axes[0, 0].axhline(y=cn_baseline, color='gray', linestyle=':', alpha=0.5)
axes[0, 0].axhline(y=4.7, color='red', linestyle='--', label='Required C/N (QPSK 3/4)')
axes[0, 0].axvline(x=45, color='green', linestyle='--', alpha=0.5, label='BS antenna (45cm)')
axes[0, 0].set_xlabel('Antenna Diameter [cm]', fontsize=11)
axes[0, 0].set_ylabel('C/N [dB]', fontsize=11)
axes[0, 0].set_title('Effect of Rx Antenna Size', fontsize=12)
axes[0, 0].legend(fontsize=9)
axes[0, 0].grid(True, alpha=0.3)

# 2. 降雨減衰の影響
rain = np.linspace(0, 15, 50)
cn_rain = []
for r in rain:
    p = baseline.copy()
    p['rain'] = r
    cn_rain.append(compute_cn(p))

axes[0, 1].plot(rain, cn_rain, 'r-', linewidth=2)
axes[0, 1].axhline(y=4.7, color='red', linestyle='--', label='Required C/N')
axes[0, 1].set_xlabel('Rain Attenuation [dB]', fontsize=11)
axes[0, 1].set_ylabel('C/N [dB]', fontsize=11)
axes[0, 1].set_title('Effect of Rain Attenuation', fontsize=12)
axes[0, 1].legend(fontsize=9)
axes[0, 1].grid(True, alpha=0.3)

# 降雨で通信断に至る閾値
cn_threshold = 4.7
rain_max = rain[np.argmin(np.abs(np.array(cn_rain) - cn_threshold))]
axes[0, 1].axvline(x=rain_max, color='orange', linestyle=':', alpha=0.7)
axes[0, 1].text(rain_max + 0.5, cn_baseline - 2,
                 f'Max rain: {rain_max:.1f} dB', fontsize=10, color='orange')

# 3. システム雑音温度の影響
temps = np.linspace(50, 500, 50)
cn_temp = []
for t in temps:
    p = baseline.copy()
    p['sys_temp'] = t
    cn_temp.append(compute_cn(p))

axes[0, 2].plot(temps, cn_temp, 'g-', linewidth=2)
axes[0, 2].axhline(y=4.7, color='red', linestyle='--', label='Required C/N')
axes[0, 2].axvline(x=120, color='green', linestyle='--', alpha=0.5, label='Typical LNB')
axes[0, 2].set_xlabel('System Noise Temperature [K]', fontsize=11)
axes[0, 2].set_ylabel('C/N [dB]', fontsize=11)
axes[0, 2].set_title('Effect of System Noise Temperature', fontsize=12)
axes[0, 2].legend(fontsize=9)
axes[0, 2].grid(True, alpha=0.3)

# 4. 仰角の影響
elevs = np.linspace(5, 90, 50)
cn_elev = []
for e in elevs:
    p = baseline.copy()
    p['elev'] = e
    # 仰角が低いと大気損失も増加
    p['atm'] = 0.1 / np.sin(np.radians(e))
    cn_elev.append(compute_cn(p))

axes[1, 0].plot(elevs, cn_elev, 'm-', linewidth=2)
axes[1, 0].axhline(y=4.7, color='red', linestyle='--', label='Required C/N')
axes[1, 0].set_xlabel('Elevation Angle [deg]', fontsize=11)
axes[1, 0].set_ylabel('C/N [dB]', fontsize=11)
axes[1, 0].set_title('Effect of Elevation Angle', fontsize=12)
axes[1, 0].legend(fontsize=9)
axes[1, 0].grid(True, alpha=0.3)

# 5. EIRP の影響
eirps = np.linspace(40, 70, 50)
cn_eirp = []
for e in eirps:
    p = baseline.copy()
    p['eirp'] = e
    cn_eirp.append(compute_cn(p))

axes[1, 1].plot(eirps, cn_eirp, 'c-', linewidth=2)
axes[1, 1].axhline(y=4.7, color='red', linestyle='--', label='Required C/N')
axes[1, 1].set_xlabel('EIRP [dBW]', fontsize=11)
axes[1, 1].set_ylabel('C/N [dB]', fontsize=11)
axes[1, 1].set_title('Effect of EIRP', fontsize=12)
axes[1, 1].legend(fontsize=9)
axes[1, 1].grid(True, alpha=0.3)

# 6. 帯域幅 vs C/N（スループットとのトレードオフ）
bws = np.logspace(6, 9, 50)  # 1 MHz〜1 GHz
cn_bw = []
for b in bws:
    p = baseline.copy()
    p['bandwidth'] = b
    cn_bw.append(compute_cn(p))

axes[1, 2].semilogx(bws / 1e6, cn_bw, 'k-', linewidth=2)
axes[1, 2].axhline(y=4.7, color='red', linestyle='--', label='Required C/N')
axes[1, 2].axvline(x=36, color='green', linestyle='--', alpha=0.5,
                     label='Standard (36 MHz)')
axes[1, 2].set_xlabel('Bandwidth [MHz]', fontsize=11)
axes[1, 2].set_ylabel('C/N [dB]', fontsize=11)
axes[1, 2].set_title('Effect of Bandwidth', fontsize=12)
axes[1, 2].legend(fontsize=9)
axes[1, 2].grid(True, alpha=0.3)

plt.suptitle('Link Budget Sensitivity Analysis\n(GEO Ku-Band Downlink)',
              fontsize=14, fontweight='bold', y=1.02)
plt.tight_layout()
plt.savefig('link_budget_sensitivity.png', dpi=150, bbox_inches='tight')
plt.show()

import numpy as np
import matplotlib.pyplot as plt

def combined_cn(cn_up_db, cn_down_db, cn_im_db=25):
    """アップリンクとダウンリンクの合成C/N [dB]"""
    cn_up = 10**(cn_up_db / 10)
    cn_down = 10**(cn_down_db / 10)
    cn_im = 10**(cn_im_db / 10)
    cn_total = 1 / (1/cn_up + 1/cn_down + 1/cn_im)
    return 10 * np.log10(cn_total)

# 合成C/Nの特性
cn_up_range = np.linspace(5, 35, 200)
cn_down_values = [8, 12, 16, 20]

fig, axes = plt.subplots(1, 2, figsize=(14, 6))

for cn_down in cn_down_values:
    cn_total = [combined_cn(cu, cn_down) for cu in cn_up_range]
    axes[0].plot(cn_up_range, cn_total, linewidth=2,
                  label=f'$C/N_{{down}}$ = {cn_down} dB')

# 対角線（上限）
axes[0].plot(cn_up_range, cn_up_range, 'k:', alpha=0.3, label='C/N = C/N_up')
axes[0].set_xlabel('Uplink C/N [dB]', fontsize=12)
axes[0].set_ylabel('Total C/N [dB]', fontsize=12)
axes[0].set_title('Combined C/N vs Uplink C/N', fontsize=13)
axes[0].legend(fontsize=10)
axes[0].grid(True, alpha=0.3)
axes[0].set_xlim([5, 35])
axes[0].set_ylim([5, 22])

# 降雨時のリンクバジェット変化
rain_atten_values = np.linspace(0, 15, 100)

# GEO Ku帯の例
cn_up_clear = 25   # dB (晴天時アップリンクC/N)
cn_down_clear = 14  # dB (晴天時ダウンリンクC/N)

# 降雨がダウンリンクに影響
cn_down_rain = cn_down_clear - rain_atten_values
cn_total_rain = [combined_cn(cn_up_clear, cd) for cd in cn_down_rain]

axes[1].plot(rain_atten_values, cn_down_rain, 'b-', linewidth=2,
              label='Downlink C/N')
axes[1].plot(rain_atten_values, cn_total_rain, 'r-', linewidth=2.5,
              label='Total C/N')
axes[1].axhline(y=cn_up_clear, color='green', linestyle='--', alpha=0.5,
                 label=f'Uplink C/N = {cn_up_clear} dB')

# 必要C/Nと降雨マージン
cn_req = 4.7
axes[1].axhline(y=cn_req, color='red', linestyle=':', linewidth=2,
                 label=f'Required C/N = {cn_req} dB')

# 降雨マージンを計算
cn_total_arr = np.array(cn_total_rain)
rain_max = rain_atten_values[np.argmin(np.abs(cn_total_arr - cn_req))]
axes[1].axvline(x=rain_max, color='orange', linestyle='--', alpha=0.7)
axes[1].fill_between(rain_atten_values, cn_req,
                       cn_total_arr, where=cn_total_arr >= cn_req,
                       alpha=0.15, color='green', label='Available margin')

axes[1].set_xlabel('Rain Attenuation on Downlink [dB]', fontsize=12)
axes[1].set_ylabel('C/N [dB]', fontsize=12)
axes[1].set_title(f'C/N Degradation with Rain\n'
                    f'(Rain margin: {rain_max:.1f} dB)', fontsize=13)
axes[1].legend(fontsize=9)
axes[1].grid(True, alpha=0.3)
axes[1].set_xlim([0, 15])
axes[1].set_ylim([-5, 30])

plt.tight_layout()
plt.savefig('combined_cn_analysis.png', dpi=150, bbox_inches='tight')
plt.show()

print(f"\n--- 合成C/N解析 ---")
print(f"晴天時アップリンクC/N: {cn_up_clear} dB")
print(f"晴天時ダウンリンクC/N: {cn_down_clear} dB")
cn_total_clear = combined_cn(cn_up_clear, cn_down_clear)
print(f"晴天時合成C/N: {cn_total_clear:.1f} dB")
print(f"必要C/N: {cn_req} dB")
print(f"降雨マージン: {rain_max:.1f} dB")
print(f"（この降雨減衰までダウンリンクが耐えられる）")