feat: add interactive exploration of Shannon's capacity formula with Plotly graphs
All checks were successful
Build & Deploy Shannon / 🏗️ Build & Deploy Shannon (push) Successful in 3m1s
All checks were successful
Build & Deploy Shannon / 🏗️ Build & Deploy Shannon (push) Successful in 3m1s
- Implemented bandwidth sensitivity and power sensitivity plots. - Created a contour map for bit rate multiplying factors. - Added input parameters for C/N and bandwidth with validation. - Displayed computed results and sensitivity analysis metrics. - Integrated interactive graphs for user exploration. - Included background information section for user guidance.
This commit is contained in:
280
core/calculations.py
Normal file
280
core/calculations.py
Normal file
@@ -0,0 +1,280 @@
|
||||
"""
|
||||
Shannon Equation - Core Calculations Module
|
||||
|
||||
All scientific computation functions extracted from the original Shannon.py.
|
||||
These are pure functions with no UI dependency.
|
||||
"""
|
||||
|
||||
from math import log, pi, sqrt, cos, acos, atan
|
||||
import numpy as np
|
||||
|
||||
# ──────────────────────────────────────────────────────────────────────────────
|
||||
# Fundamental Shannon Functions
|
||||
# ──────────────────────────────────────────────────────────────────────────────
|
||||
|
||||
def combine_cnr(*cnr_values: float) -> float:
|
||||
"""Combine multiple Carrier-to-Noise Ratios (in dB) into one equivalent C/N.
|
||||
|
||||
Uses the summation of normalized noise variances:
|
||||
1/CNR_total = sum(1/CNR_i)
|
||||
"""
|
||||
ncr_linear = 0.0
|
||||
for cnr_db in cnr_values:
|
||||
ncr_linear += 10 ** (-cnr_db / 10)
|
||||
return -10 * log(ncr_linear, 10)
|
||||
|
||||
|
||||
def shannon_capacity(bw: float = 36.0, cnr: float = 10.0, penalty: float = 0.0) -> float:
|
||||
"""Shannon channel capacity (bit rate in Mbps).
|
||||
|
||||
Args:
|
||||
bw: Bandwidth in MHz.
|
||||
cnr: Carrier-to-Noise Ratio in dB.
|
||||
penalty: Implementation penalty in dB.
|
||||
"""
|
||||
cnr_linear = 10 ** ((cnr - penalty) / 10)
|
||||
return bw * log(1 + cnr_linear, 2)
|
||||
|
||||
|
||||
def br_multiplier(bw_mul: float = 1.0, p_mul: float = 2.0, cnr: float = 10.0) -> float:
|
||||
"""Bit Rate multiplying factor when BW and Power are scaled."""
|
||||
cnr_linear = 10 ** (cnr / 10)
|
||||
return bw_mul * log(1 + cnr_linear * p_mul / bw_mul, 2) / log(1 + cnr_linear, 2)
|
||||
|
||||
|
||||
def shannon_points(bw: float = 36.0, cnr: float = 10.0):
|
||||
"""Compute key Shannon operating points.
|
||||
|
||||
Returns:
|
||||
(cnr_linear, br_infinity, c_n0, br_constrained)
|
||||
"""
|
||||
cnr_linear = 10 ** (cnr / 10)
|
||||
c_n0 = cnr_linear * bw
|
||||
br_infinity = c_n0 / log(2)
|
||||
br_constrained = shannon_capacity(bw, cnr)
|
||||
return cnr_linear, br_infinity, c_n0, br_constrained
|
||||
|
||||
|
||||
# ──────────────────────────────────────────────────────────────────────────────
|
||||
# Satellite Link Budget Calculations
|
||||
# ──────────────────────────────────────────────────────────────────────────────
|
||||
|
||||
def compute_satellite_link(
|
||||
freq_ghz: float,
|
||||
hpa_power_w: float,
|
||||
sat_loss_db: float,
|
||||
sat_cir_list: list[float],
|
||||
sat_beam_deg: float,
|
||||
gain_offset_db: float,
|
||||
sat_alt_km: float,
|
||||
sat_lat: float,
|
||||
sat_lon: float,
|
||||
gs_lat: float,
|
||||
gs_lon: float,
|
||||
availability_pct: float,
|
||||
) -> dict:
|
||||
"""Compute all satellite link parameters.
|
||||
|
||||
Returns a dict with all computed values.
|
||||
"""
|
||||
import itur
|
||||
|
||||
R_EARTH = 6378 # km
|
||||
SAT_ANT_EFF = 0.65
|
||||
|
||||
# Signal power after losses
|
||||
sig_power = hpa_power_w * 10 ** (-sat_loss_db / 10)
|
||||
|
||||
# Satellite antenna
|
||||
wavelength = 300e6 / freq_ghz / 1e9 # meters
|
||||
sat_gain_linear = SAT_ANT_EFF * (pi * 70 / sat_beam_deg) ** 2
|
||||
sat_gain_linear *= 10 ** (-gain_offset_db / 10)
|
||||
|
||||
# EIRP
|
||||
eirp_linear = sig_power * sat_gain_linear
|
||||
|
||||
# Satellite C/I
|
||||
sat_cir = combine_cnr(*sat_cir_list)
|
||||
|
||||
# Path geometry
|
||||
path_length = sqrt(
|
||||
sat_alt_km ** 2
|
||||
+ 2 * R_EARTH * (R_EARTH + sat_alt_km)
|
||||
* (1 - cos(np.radians(sat_lat - gs_lat)) * cos(np.radians(sat_lon - gs_lon)))
|
||||
)
|
||||
|
||||
phi = acos(cos(np.radians(sat_lat - gs_lat)) * cos(np.radians(sat_lon - gs_lon)))
|
||||
if phi > 0:
|
||||
elevation = float(np.degrees(
|
||||
atan((cos(phi) - R_EARTH / (R_EARTH + sat_alt_km)) / sqrt(1 - cos(phi) ** 2))
|
||||
))
|
||||
else:
|
||||
elevation = 90.0
|
||||
|
||||
# Atmospheric attenuation
|
||||
if elevation <= 0:
|
||||
atm_loss_db = 999.0
|
||||
else:
|
||||
atm_loss_db = float(
|
||||
itur.atmospheric_attenuation_slant_path(
|
||||
gs_lat, gs_lon, freq_ghz, elevation, 100 - availability_pct, 1
|
||||
).value
|
||||
)
|
||||
|
||||
# Path dispersion
|
||||
path_loss_linear = 4 * pi * (path_length * 1000) ** 2
|
||||
free_space_loss_linear = (4 * pi * path_length * 1000 / wavelength) ** 2
|
||||
|
||||
# PFD
|
||||
pfd_linear = eirp_linear / path_loss_linear * 10 ** (-atm_loss_db / 10)
|
||||
|
||||
return {
|
||||
"sig_power": sig_power,
|
||||
"wavelength": wavelength,
|
||||
"sat_gain_linear": sat_gain_linear,
|
||||
"eirp_linear": eirp_linear,
|
||||
"sat_cir": sat_cir,
|
||||
"path_length": path_length,
|
||||
"elevation": elevation,
|
||||
"atm_loss_db": atm_loss_db,
|
||||
"path_loss_linear": path_loss_linear,
|
||||
"free_space_loss_linear": free_space_loss_linear,
|
||||
"pfd_linear": pfd_linear,
|
||||
}
|
||||
|
||||
|
||||
def compute_receiver(
|
||||
pfd_linear: float,
|
||||
atm_loss_db: float,
|
||||
wavelength: float,
|
||||
cpe_ant_d: float,
|
||||
cpe_t_clear: float,
|
||||
) -> dict:
|
||||
"""Compute receiver-side parameters."""
|
||||
CPE_ANT_EFF = 0.6
|
||||
K_BOLTZ = 1.38e-23 # J/K
|
||||
|
||||
cpe_t_att = (cpe_t_clear - 40) + 40 * 10 ** (-atm_loss_db / 10) + 290 * (1 - 10 ** (-atm_loss_db / 10))
|
||||
|
||||
cpe_ae = pi * cpe_ant_d ** 2 / 4 * CPE_ANT_EFF
|
||||
cpe_gain_linear = (pi * cpe_ant_d / wavelength) ** 2 * CPE_ANT_EFF
|
||||
cpe_g_t = 10 * log(cpe_gain_linear / cpe_t_att, 10)
|
||||
|
||||
rx_power = pfd_linear * cpe_ae
|
||||
n0 = K_BOLTZ * cpe_t_att
|
||||
c_n0_hz = rx_power / n0
|
||||
c_n0_mhz = c_n0_hz / 1e6
|
||||
|
||||
br_infinity = c_n0_mhz / log(2)
|
||||
|
||||
# Spectral efficiency points
|
||||
bw_spe_1 = c_n0_mhz
|
||||
bw_spe_double = c_n0_mhz / (2 ** 2 - 1)
|
||||
|
||||
br_spe_1 = bw_spe_1
|
||||
br_spe_double = bw_spe_double * 2
|
||||
|
||||
return {
|
||||
"cpe_ae": cpe_ae,
|
||||
"cpe_gain_linear": cpe_gain_linear,
|
||||
"cpe_g_t": cpe_g_t,
|
||||
"cpe_t_att": cpe_t_att,
|
||||
"rx_power": rx_power,
|
||||
"n0": n0,
|
||||
"c_n0_hz": c_n0_hz,
|
||||
"c_n0_mhz": c_n0_mhz,
|
||||
"br_infinity": br_infinity,
|
||||
"bw_spe_1": bw_spe_1,
|
||||
"br_spe_1": br_spe_1,
|
||||
"br_spe_double": br_spe_double,
|
||||
}
|
||||
|
||||
|
||||
def compute_baseband(
|
||||
c_n0_mhz: float,
|
||||
br_infinity: float,
|
||||
bw_spe_1: float,
|
||||
sat_cir: float,
|
||||
bandwidth: float,
|
||||
rolloff: float,
|
||||
overheads: float,
|
||||
cnr_imp_list: list[float],
|
||||
penalties: float,
|
||||
) -> dict:
|
||||
"""Compute baseband processing results."""
|
||||
cnr_imp = combine_cnr(*cnr_imp_list)
|
||||
|
||||
cnr_spe_1 = 0.0 # dB
|
||||
cnr_bw = cnr_spe_1 + 10 * log(bw_spe_1 / bandwidth, 10)
|
||||
bw_nyq = bandwidth / (1 + rolloff / 100)
|
||||
cnr_nyq = cnr_spe_1 + 10 * log(bw_spe_1 / bw_nyq, 10)
|
||||
cnr_rcv = combine_cnr(cnr_nyq, cnr_imp, sat_cir)
|
||||
|
||||
br_nyq = shannon_capacity(bw_nyq, cnr_nyq)
|
||||
br_rcv = shannon_capacity(bw_nyq, cnr_rcv, penalties)
|
||||
br_rcv_higher = br_rcv / (1 + overheads / 100)
|
||||
|
||||
spe_nyq = br_nyq / bandwidth
|
||||
bits_per_symbol = br_nyq / bw_nyq
|
||||
spe_rcv = br_rcv / bandwidth
|
||||
spe_higher = br_rcv_higher / bandwidth
|
||||
|
||||
return {
|
||||
"cnr_bw": cnr_bw,
|
||||
"cnr_nyq": cnr_nyq,
|
||||
"cnr_rcv": cnr_rcv,
|
||||
"cnr_imp": cnr_imp,
|
||||
"bw_nyq": bw_nyq,
|
||||
"br_nyq": br_nyq,
|
||||
"br_rcv": br_rcv,
|
||||
"br_rcv_higher": br_rcv_higher,
|
||||
"br_nyq_norm": br_nyq / br_infinity,
|
||||
"br_rcv_norm": br_rcv / br_infinity,
|
||||
"br_rcv_h_norm": br_rcv_higher / br_infinity,
|
||||
"spe_nyq": spe_nyq,
|
||||
"bits_per_symbol": bits_per_symbol,
|
||||
"spe_rcv": spe_rcv,
|
||||
"spe_higher": spe_higher,
|
||||
"bandwidth": bandwidth,
|
||||
"rolloff": rolloff,
|
||||
"overheads": overheads,
|
||||
"penalties": penalties,
|
||||
}
|
||||
|
||||
|
||||
# ──────────────────────────────────────────────────────────────────────────────
|
||||
# Formatting Helpers
|
||||
# ──────────────────────────────────────────────────────────────────────────────
|
||||
|
||||
def fmt_br(br: float) -> str:
|
||||
return f"{br:.1f} Mbps"
|
||||
|
||||
|
||||
def fmt_power(p: float) -> str:
|
||||
p_db = 10 * log(p, 10)
|
||||
if 1 < p < 1e4:
|
||||
return f"{p:.1f} W · {p_db:.1f} dBW"
|
||||
elif 1e-3 < p <= 1:
|
||||
return f"{p:.4f} W · {p_db:.1f} dBW"
|
||||
else:
|
||||
return f"{p:.1e} W · {p_db:.1f} dBW"
|
||||
|
||||
|
||||
def fmt_pfd(p: float) -> str:
|
||||
p_db = 10 * log(p, 10)
|
||||
return f"{p:.1e} W/m² · {p_db:.1f} dBW/m²"
|
||||
|
||||
|
||||
def fmt_psd(p: float) -> str:
|
||||
p_db = 10 * log(p, 10)
|
||||
return f"{p:.1e} W/MHz · {p_db:.1f} dBW/MHz"
|
||||
|
||||
|
||||
def fmt_gain(g: float) -> str:
|
||||
g_db = 10 * log(g, 10)
|
||||
return f"{g:.1f} · {g_db:.1f} dBi"
|
||||
|
||||
|
||||
def fmt_ploss(loss: float) -> str:
|
||||
loss_db = 10 * log(loss, 10)
|
||||
return f"{loss:.2e} m² · {loss_db:.1f} dBm²"
|
||||
Reference in New Issue
Block a user