Wideband autonomous SDR analysis engine forked from sdr-visual-suite
You can not select more than 25 topics Topics must start with a letter or number, can include dashes ('-') and can be up to 35 characters long.

303 line
7.3KB

  1. package detector
  2. import (
  3. "math"
  4. "sort"
  5. "time"
  6. "sdr-visual-suite/internal/classifier"
  7. )
  8. type Event struct {
  9. ID int64 `json:"id"`
  10. Start time.Time `json:"start"`
  11. End time.Time `json:"end"`
  12. CenterHz float64 `json:"center_hz"`
  13. Bandwidth float64 `json:"bandwidth_hz"`
  14. PeakDb float64 `json:"peak_db"`
  15. SNRDb float64 `json:"snr_db"`
  16. FirstBin int `json:"first_bin"`
  17. LastBin int `json:"last_bin"`
  18. Class *classifier.Classification `json:"class,omitempty"`
  19. }
  20. type Detector struct {
  21. ThresholdDb float64
  22. MinDuration time.Duration
  23. Hold time.Duration
  24. EmaAlpha float64
  25. HysteresisDb float64
  26. MinStableFrames int
  27. GapTolerance time.Duration
  28. binWidth float64
  29. nbins int
  30. sampleRate int
  31. ema []float64
  32. active map[int64]*activeEvent
  33. nextID int64
  34. }
  35. type activeEvent struct {
  36. id int64
  37. start time.Time
  38. lastSeen time.Time
  39. centerHz float64
  40. bwHz float64
  41. peakDb float64
  42. snrDb float64
  43. firstBin int
  44. lastBin int
  45. class *classifier.Classification
  46. stableHits int
  47. }
  48. type Signal struct {
  49. FirstBin int `json:"first_bin"`
  50. LastBin int `json:"last_bin"`
  51. CenterHz float64 `json:"center_hz"`
  52. BWHz float64 `json:"bw_hz"`
  53. PeakDb float64 `json:"peak_db"`
  54. SNRDb float64 `json:"snr_db"`
  55. Class *classifier.Classification `json:"class,omitempty"`
  56. }
  57. func New(thresholdDb float64, sampleRate int, fftSize int, minDur, hold time.Duration, emaAlpha, hysteresis float64, minStable int, gapTolerance time.Duration) *Detector {
  58. if minDur <= 0 {
  59. minDur = 250 * time.Millisecond
  60. }
  61. if hold <= 0 {
  62. hold = 500 * time.Millisecond
  63. }
  64. if emaAlpha <= 0 || emaAlpha > 1 {
  65. emaAlpha = 0.2
  66. }
  67. if hysteresis <= 0 {
  68. hysteresis = 3
  69. }
  70. if minStable <= 0 {
  71. minStable = 3
  72. }
  73. if gapTolerance <= 0 {
  74. gapTolerance = hold
  75. }
  76. return &Detector{
  77. ThresholdDb: thresholdDb,
  78. MinDuration: minDur,
  79. Hold: hold,
  80. EmaAlpha: emaAlpha,
  81. HysteresisDb: hysteresis,
  82. MinStableFrames: minStable,
  83. GapTolerance: gapTolerance,
  84. binWidth: float64(sampleRate) / float64(fftSize),
  85. nbins: fftSize,
  86. sampleRate: sampleRate,
  87. ema: make([]float64, fftSize),
  88. active: map[int64]*activeEvent{},
  89. nextID: 1,
  90. }
  91. }
  92. func (d *Detector) Process(now time.Time, spectrum []float64, centerHz float64) ([]Event, []Signal) {
  93. signals := d.detectSignals(spectrum, centerHz)
  94. finished := d.matchSignals(now, signals)
  95. return finished, signals
  96. }
  97. // UpdateClasses refreshes active event classes from current signals.
  98. func (d *Detector) UpdateClasses(signals []Signal) {
  99. for _, s := range signals {
  100. for _, ev := range d.active {
  101. if overlapHz(s.CenterHz, s.BWHz, ev.centerHz, ev.bwHz) && math.Abs(s.CenterHz-ev.centerHz) < (s.BWHz+ev.bwHz)/2.0 {
  102. if s.Class != nil {
  103. if ev.class == nil || s.Class.Confidence >= ev.class.Confidence {
  104. ev.class = s.Class
  105. }
  106. }
  107. }
  108. }
  109. }
  110. }
  111. func (d *Detector) detectSignals(spectrum []float64, centerHz float64) []Signal {
  112. n := len(spectrum)
  113. if n == 0 {
  114. return nil
  115. }
  116. smooth := d.smoothSpectrum(spectrum)
  117. thresholdOn := d.ThresholdDb
  118. thresholdOff := d.ThresholdDb - d.HysteresisDb
  119. noise := median(smooth)
  120. var signals []Signal
  121. in := false
  122. start := 0
  123. peak := -1e9
  124. peakBin := 0
  125. for i := 0; i < n; i++ {
  126. v := smooth[i]
  127. if v >= thresholdOn {
  128. if !in {
  129. in = true
  130. start = i
  131. peak = v
  132. peakBin = i
  133. } else if v > peak {
  134. peak = v
  135. peakBin = i
  136. }
  137. } else if in && v < thresholdOff {
  138. signals = append(signals, d.makeSignal(start, i-1, peak, peakBin, noise, centerHz, smooth))
  139. in = false
  140. }
  141. }
  142. if in {
  143. signals = append(signals, d.makeSignal(start, n-1, peak, peakBin, noise, centerHz, smooth))
  144. }
  145. return signals
  146. }
  147. func (d *Detector) makeSignal(first, last int, peak float64, peakBin int, noise float64, centerHz float64, spectrum []float64) Signal {
  148. centerBin := float64(first+last) / 2.0
  149. centerFreq := centerHz + (centerBin-float64(d.nbins)/2.0)*d.binWidth
  150. bw := float64(last-first+1) * d.binWidth
  151. snr := peak - noise
  152. return Signal{
  153. FirstBin: first,
  154. LastBin: last,
  155. CenterHz: centerFreq,
  156. BWHz: bw,
  157. PeakDb: peak,
  158. SNRDb: snr,
  159. }
  160. }
  161. func (d *Detector) smoothSpectrum(spectrum []float64) []float64 {
  162. if d.ema == nil || len(d.ema) != len(spectrum) {
  163. d.ema = make([]float64, len(spectrum))
  164. copy(d.ema, spectrum)
  165. return d.ema
  166. }
  167. alpha := d.EmaAlpha
  168. for i := range spectrum {
  169. v := spectrum[i]
  170. d.ema[i] = alpha*v + (1-alpha)*d.ema[i]
  171. }
  172. return d.ema
  173. }
  174. func (d *Detector) matchSignals(now time.Time, signals []Signal) []Event {
  175. used := make(map[int64]bool, len(d.active))
  176. for _, s := range signals {
  177. var best *activeEvent
  178. var candidates []struct {
  179. ev *activeEvent
  180. dist float64
  181. }
  182. for _, ev := range d.active {
  183. if overlapHz(s.CenterHz, s.BWHz, ev.centerHz, ev.bwHz) && math.Abs(s.CenterHz-ev.centerHz) < (s.BWHz+ev.bwHz)/2.0 {
  184. candidates = append(candidates, struct {
  185. ev *activeEvent
  186. dist float64
  187. }{ev: ev, dist: math.Abs(s.CenterHz - ev.centerHz)})
  188. }
  189. }
  190. if len(candidates) > 0 {
  191. sort.Slice(candidates, func(i, j int) bool { return candidates[i].dist < candidates[j].dist })
  192. best = candidates[0].ev
  193. }
  194. if best == nil {
  195. id := d.nextID
  196. d.nextID++
  197. d.active[id] = &activeEvent{
  198. id: id,
  199. start: now,
  200. lastSeen: now,
  201. centerHz: s.CenterHz,
  202. bwHz: s.BWHz,
  203. peakDb: s.PeakDb,
  204. snrDb: s.SNRDb,
  205. firstBin: s.FirstBin,
  206. lastBin: s.LastBin,
  207. class: s.Class,
  208. stableHits: 1,
  209. }
  210. continue
  211. }
  212. used[best.id] = true
  213. best.lastSeen = now
  214. best.stableHits++
  215. best.centerHz = (best.centerHz + s.CenterHz) / 2.0
  216. if s.BWHz > best.bwHz {
  217. best.bwHz = s.BWHz
  218. }
  219. if s.PeakDb > best.peakDb {
  220. best.peakDb = s.PeakDb
  221. }
  222. if s.SNRDb > best.snrDb {
  223. best.snrDb = s.SNRDb
  224. }
  225. if s.FirstBin < best.firstBin {
  226. best.firstBin = s.FirstBin
  227. }
  228. if s.LastBin > best.lastBin {
  229. best.lastBin = s.LastBin
  230. }
  231. if s.Class != nil {
  232. if best.class == nil || s.Class.Confidence >= best.class.Confidence {
  233. best.class = s.Class
  234. }
  235. }
  236. }
  237. var finished []Event
  238. for id, ev := range d.active {
  239. if used[id] {
  240. continue
  241. }
  242. if now.Sub(ev.lastSeen) < d.GapTolerance {
  243. continue
  244. }
  245. duration := ev.lastSeen.Sub(ev.start)
  246. if duration < d.MinDuration || ev.stableHits < d.MinStableFrames {
  247. delete(d.active, id)
  248. continue
  249. }
  250. finished = append(finished, Event{
  251. ID: ev.id,
  252. Start: ev.start,
  253. End: ev.lastSeen,
  254. CenterHz: ev.centerHz,
  255. Bandwidth: ev.bwHz,
  256. PeakDb: ev.peakDb,
  257. SNRDb: ev.snrDb,
  258. FirstBin: ev.firstBin,
  259. LastBin: ev.lastBin,
  260. Class: ev.class,
  261. })
  262. delete(d.active, id)
  263. }
  264. return finished
  265. }
  266. func overlapHz(c1, b1, c2, b2 float64) bool {
  267. l1 := c1 - b1/2.0
  268. r1 := c1 + b1/2.0
  269. l2 := c2 - b2/2.0
  270. r2 := c2 + b2/2.0
  271. return l1 <= r2 && l2 <= r1
  272. }
  273. func median(vals []float64) float64 {
  274. if len(vals) == 0 {
  275. return 0
  276. }
  277. cpy := append([]float64(nil), vals...)
  278. sort.Float64s(cpy)
  279. mid := len(cpy) / 2
  280. if len(cpy)%2 == 0 {
  281. return (cpy[mid-1] + cpy[mid]) / 2.0
  282. }
  283. return cpy[mid]
  284. }