1// Package reception estimates over-the-air signal strength from transmitter2// parameters and a terrain profile.3//4// The model is free-space loss plus Deygout multiple knife-edge diffraction5// over the terrain profile on a 4/3 effective earth. It captures line of6// sight and ridge blockage, which dominate reception in mountainous areas,7// but not multipath, foliage, buildings or directional transmit patterns.8// Treat results as a ranking and a rough margin, not a guarantee.9package reception1011import (12 "math"1314 "airwaves/internal/fcc"15 "airwaves/internal/terrain"16)1718// Tier buckets a noise margin into a reception expectation.19type Tier string2021// Tiers from best to worst.22const (23 Strong Tier = "strong"24 Good Tier = "good"25 Fair Tier = "fair"26 Weak Tier = "weak"27 Unlikely Tier = "unlikely"28 Unknown Tier = "unknown"29)3031// Receivable reports whether a tier should decode reliably.32func (t Tier) Receivable() bool { return t == Strong || t == Good || t == Fair }3334const (35 // effectiveEarthRadiusM is 4/3 of the earth's radius, the standard36 // allowance for atmospheric refraction.37 effectiveEarthRadiusM = 8_495_000.038 // thermalNoiseDBm is kTB for a 6 MHz channel plus a 5 dB tuner noise39 // figure.40 thermalNoiseDBm = -106.2 + 541 // requiredCNRDB is the ATSC 1.0 (8-VSB) decoding threshold. Robust42 // ATSC 3.0 modes decode a few dB lower.43 requiredCNRDB = 15.244 // dipoleDBi converts ERP (relative to a dipole) to EIRP.45 dipoleDBi = 2.1546)4748// Preset describes a receive antenna installation.49type Preset struct {50 Name string `json:"name"`51 Label string `json:"label"`52 HeightM float64 `json:"heightM"`53 GainDBi map[string]float64 `json:"gainDbi"` // by band, see fcc.Band54 LossDB float64 `json:"lossDb"` // walls, roof, cable55}5657// Presets are the installations every transmitter is evaluated against.58var Presets = []Preset{59 {60 Name: "indoor", Label: "Indoor antenna", HeightM: 3, LossDB: 12,61 GainDBi: map[string]float64{"UHF": 2, "VHF-Hi": -2, "VHF-Lo": -6},62 },63 {64 Name: "attic", Label: "Attic antenna", HeightM: 6, LossDB: 7,65 GainDBi: map[string]float64{"UHF": 6, "VHF-Hi": 3, "VHF-Lo": -2},66 },67 {68 Name: "rooftop", Label: "Rooftop antenna, 30 ft", HeightM: 9.1, LossDB: 2,69 GainDBi: map[string]float64{"UHF": 9, "VHF-Hi": 6, "VHF-Lo": 2},70 },71}7273// manMadeNoiseDB is ambient electrical noise above thermal, by band. Low VHF74// suffers most from household interference.75var manMadeNoiseDB = map[string]float64{"UHF": 0, "VHF-Hi": 4, "VHF-Lo": 12}7677// Transmitter is the subset of facility data the model needs.78type Transmitter struct {79 RFChannel int80 ERPkW float6481 RCAMSLm float6482 RCAGLm float6483 HAATm float6484}8586// FromFacility extracts model inputs from an FCC record.87func FromFacility(f fcc.Facility) Transmitter {88 return Transmitter{RFChannel: f.RFChannel, ERPkW: f.ERPkW, RCAMSLm: f.RCAMSLm, RCAGLm: f.RCAGLm, HAATm: f.HAATm}89}9091// Estimate is the predicted signal for one transmitter and preset.92type Estimate struct {93 Preset string `json:"preset"`94 PowerDBm float64 `json:"powerDbm"` // at the tuner input95 NoiseMarginDB float64 `json:"noiseMarginDb"`96 FieldDBuVm float64 `json:"fieldDbuVm"` // field strength at the antenna97 PathLossDB float64 `json:"pathLossDb"`98 DiffractionDB float64 `json:"diffractionDb"`99 LineOfSight bool `json:"lineOfSight"`100 Tier Tier `json:"tier"`101}102103// Predict estimates reception of tx over prof with antenna preset p.104func Predict(tx Transmitter, prof terrain.Profile, p Preset) Estimate {105 n := len(prof.Elevations)106 if n < 3 || tx.ERPkW <= 0 || prof.DistanceKm <= 0 {107 return Estimate{Preset: p.Name, Tier: Unknown}108 }109 fMHz := fcc.CenterMHz(tx.RFChannel)110 band := fcc.Band(tx.RFChannel)111 lambda := 299.792458 / fMHz112 dM := prof.DistanceKm * 1000113114 pos := make([]float64, n)115 z := make([]float64, n)116 for i, e := range prof.Elevations {117 pos[i] = dM * float64(i) / float64(n-1)118 // Raise terrain by the earth bulge relative to the straight chord.119 z[i] = e + pos[i]*(dM-pos[i])/(2*effectiveEarthRadiusM)120 }121 rxH := prof.Elevations[0] + p.HeightM122 txH := txHeightAMSL(tx, prof.Elevations[n-1])123124 _, vMain := worstEdge(pos, z, 0, n-1, rxH, txH, lambda)125 diff := deygout(pos, z, 0, n-1, rxH, txH, lambda, 2)126 fspl := 32.45 + 20*math.Log10(prof.DistanceKm) + 20*math.Log10(fMHz)127 isoDBm := 60 + 10*math.Log10(tx.ERPkW) + dipoleDBi - fspl - diff128 power := isoDBm + p.GainDBi[band] - p.LossDB129 margin := power - (thermalNoiseDBm + manMadeNoiseDB[band] + requiredCNRDB)130131 return Estimate{132 Preset: p.Name,133 PowerDBm: round1(power),134 NoiseMarginDB: round1(margin),135 FieldDBuVm: round1(isoDBm + 20*math.Log10(fMHz) + 77.2),136 PathLossDB: round1(fspl + diff),137 DiffractionDB: round1(diff),138 LineOfSight: vMain < 0,139 Tier: tierFor(margin),140 }141}142143func tierFor(margin float64) Tier {144 switch {145 case margin >= 20:146 return Strong147 case margin >= 10:148 return Good149 case margin >= 3:150 return Fair151 case margin >= -5:152 return Weak153 default:154 return Unlikely155 }156}157158func txHeightAMSL(tx Transmitter, ground float64) float64 {159 h := tx.RCAMSLm160 if h <= 0 {161 agl := tx.RCAGLm162 if agl <= 0 {163 agl = max(tx.HAATm, 30)164 }165 h = ground + agl166 }167 // Guard against records whose height sits below the sampled ground.168 return max(h, ground+10)169}170171// worstEdge finds the terrain point with the largest Fresnel-Kirchhoff172// parameter between indexes i and j.173func worstEdge(pos, z []float64, i, j int, hi, hj, lambda float64) (int, float64) {174 k, v := -1, math.Inf(-1)175 span := pos[j] - pos[i]176 for m := i + 1; m < j; m++ {177 d1, d2 := pos[m]-pos[i], pos[j]-pos[m]178 los := hi + (hj-hi)*d1/span179 vm := (z[m] - los) * math.Sqrt(2*span/(lambda*d1*d2))180 if vm > v {181 k, v = m, vm182 }183 }184 return k, v185}186187// deygout sums knife-edge losses: the dominant edge, then recursively the188// dominant edges of the sub-paths on either side of it.189func deygout(pos, z []float64, i, j int, hi, hj, lambda float64, depth int) float64 {190 if j-i < 2 {191 return 0192 }193 k, v := worstEdge(pos, z, i, j, hi, hj, lambda)194 if k < 0 || v <= -0.78 {195 return 0196 }197 loss := knifeEdgeDB(v)198 if depth > 0 {199 loss += deygout(pos, z, i, k, hi, z[k], lambda, depth-1)200 loss += deygout(pos, z, k, j, z[k], hj, lambda, depth-1)201 }202 return loss203}204205// knifeEdgeDB is the ITU-R P.526 approximation of single knife-edge loss.206func knifeEdgeDB(v float64) float64 {207 if v <= -0.78 {208 return 0209 }210 return 6.9 + 20*math.Log10(math.Sqrt((v-0.1)*(v-0.1)+1)+v-0.1)211}212213func round1(v float64) float64 { return math.Round(v*10) / 10 }