Você não pode selecionar mais de 25 tópicos Os tópicos devem começar com uma letra ou um número, podem incluir traços ('-') e podem ter até 35 caracteres.

621 linhas
15KB

  1. package detector
  2. import (
  3. "math"
  4. "sort"
  5. "time"
  6. "sdr-visual-suite/internal/cfar"
  7. "sdr-visual-suite/internal/classifier"
  8. "sdr-visual-suite/internal/config"
  9. )
  10. type Event struct {
  11. ID int64 `json:"id"`
  12. Start time.Time `json:"start"`
  13. End time.Time `json:"end"`
  14. CenterHz float64 `json:"center_hz"`
  15. Bandwidth float64 `json:"bandwidth_hz"`
  16. PeakDb float64 `json:"peak_db"`
  17. SNRDb float64 `json:"snr_db"`
  18. FirstBin int `json:"first_bin"`
  19. LastBin int `json:"last_bin"`
  20. Class *classifier.Classification `json:"class,omitempty"`
  21. }
  22. type Detector struct {
  23. ThresholdDb float64
  24. MinDuration time.Duration
  25. Hold time.Duration
  26. EmaAlpha float64
  27. HysteresisDb float64
  28. MinStableFrames int
  29. GapTolerance time.Duration
  30. CFARScaleDb float64
  31. EdgeMarginDb float64
  32. MaxSignalBwHz float64
  33. MergeGapHz float64
  34. classHistorySize int
  35. classSwitchRatio float64
  36. binWidth float64
  37. nbins int
  38. sampleRate int
  39. ema []float64
  40. active map[int64]*activeEvent
  41. nextID int64
  42. cfarEngine cfar.CFAR
  43. lastThresholds []float64
  44. lastNoiseFloor float64
  45. }
  46. type activeEvent struct {
  47. id int64
  48. start time.Time
  49. lastSeen time.Time
  50. centerHz float64
  51. bwHz float64
  52. peakDb float64
  53. snrDb float64
  54. firstBin int
  55. lastBin int
  56. class *classifier.Classification
  57. stableHits int
  58. classHistory []classifier.SignalClass
  59. classIdx int
  60. }
  61. type Signal struct {
  62. FirstBin int `json:"first_bin"`
  63. LastBin int `json:"last_bin"`
  64. CenterHz float64 `json:"center_hz"`
  65. BWHz float64 `json:"bw_hz"`
  66. PeakDb float64 `json:"peak_db"`
  67. SNRDb float64 `json:"snr_db"`
  68. NoiseDb float64 `json:"noise_db,omitempty"`
  69. Class *classifier.Classification `json:"class,omitempty"`
  70. }
  71. func New(detCfg config.DetectorConfig, sampleRate int, fftSize int) *Detector {
  72. minDur := time.Duration(detCfg.MinDurationMs) * time.Millisecond
  73. hold := time.Duration(detCfg.HoldMs) * time.Millisecond
  74. gapTolerance := time.Duration(detCfg.GapToleranceMs) * time.Millisecond
  75. emaAlpha := detCfg.EmaAlpha
  76. hysteresis := detCfg.HysteresisDb
  77. minStable := detCfg.MinStableFrames
  78. cfarMode := detCfg.CFARMode
  79. binWidth := float64(sampleRate) / float64(fftSize)
  80. cfarGuard := int(math.Ceil(detCfg.CFARGuardHz / binWidth))
  81. if cfarGuard < 0 {
  82. cfarGuard = 0
  83. }
  84. cfarTrain := int(math.Ceil(detCfg.CFARTrainHz / binWidth))
  85. if cfarTrain < 1 {
  86. cfarTrain = 1
  87. }
  88. cfarRank := detCfg.CFARRank
  89. cfarScaleDb := detCfg.CFARScaleDb
  90. cfarWrap := detCfg.CFARWrapAround
  91. thresholdDb := detCfg.ThresholdDb
  92. edgeMarginDb := detCfg.EdgeMarginDb
  93. maxSignalBwHz := detCfg.MaxSignalBwHz
  94. mergeGapHz := detCfg.MergeGapHz
  95. classHistorySize := detCfg.ClassHistorySize
  96. classSwitchRatio := detCfg.ClassSwitchRatio
  97. if minDur <= 0 {
  98. minDur = 250 * time.Millisecond
  99. }
  100. if hold <= 0 {
  101. hold = 500 * time.Millisecond
  102. }
  103. if emaAlpha <= 0 || emaAlpha > 1 {
  104. emaAlpha = 0.2
  105. }
  106. if hysteresis <= 0 {
  107. hysteresis = 3
  108. }
  109. if minStable <= 0 {
  110. minStable = 3
  111. }
  112. if gapTolerance <= 0 {
  113. gapTolerance = hold
  114. }
  115. if cfarScaleDb <= 0 {
  116. cfarScaleDb = 6
  117. }
  118. if edgeMarginDb <= 0 {
  119. edgeMarginDb = 3.0
  120. }
  121. if maxSignalBwHz <= 0 {
  122. maxSignalBwHz = 150000
  123. }
  124. if mergeGapHz <= 0 {
  125. mergeGapHz = 5000
  126. }
  127. if classHistorySize <= 0 {
  128. classHistorySize = 10
  129. }
  130. if classSwitchRatio <= 0 || classSwitchRatio > 1 {
  131. classSwitchRatio = 0.6
  132. }
  133. if cfarRank <= 0 || cfarRank > 2*cfarTrain {
  134. cfarRank = int(math.Round(0.75 * float64(2*cfarTrain)))
  135. if cfarRank <= 0 {
  136. cfarRank = 1
  137. }
  138. }
  139. var cfarEngine cfar.CFAR
  140. if cfarMode != "" && cfarMode != "OFF" {
  141. cfarEngine = cfar.New(cfar.Config{
  142. Mode: cfar.Mode(cfarMode),
  143. GuardCells: cfarGuard,
  144. TrainCells: cfarTrain,
  145. Rank: cfarRank,
  146. ScaleDb: cfarScaleDb,
  147. WrapAround: cfarWrap,
  148. })
  149. }
  150. return &Detector{
  151. ThresholdDb: thresholdDb,
  152. MinDuration: minDur,
  153. Hold: hold,
  154. EmaAlpha: emaAlpha,
  155. HysteresisDb: hysteresis,
  156. MinStableFrames: minStable,
  157. GapTolerance: gapTolerance,
  158. CFARScaleDb: cfarScaleDb,
  159. EdgeMarginDb: edgeMarginDb,
  160. MaxSignalBwHz: maxSignalBwHz,
  161. MergeGapHz: mergeGapHz,
  162. classHistorySize: classHistorySize,
  163. classSwitchRatio: classSwitchRatio,
  164. binWidth: float64(sampleRate) / float64(fftSize),
  165. nbins: fftSize,
  166. sampleRate: sampleRate,
  167. ema: make([]float64, fftSize),
  168. active: map[int64]*activeEvent{},
  169. nextID: 1,
  170. cfarEngine: cfarEngine,
  171. }
  172. }
  173. func (d *Detector) Process(now time.Time, spectrum []float64, centerHz float64) ([]Event, []Signal) {
  174. signals := d.detectSignals(spectrum, centerHz)
  175. finished := d.matchSignals(now, signals)
  176. return finished, signals
  177. }
  178. func (d *Detector) LastThresholds() []float64 {
  179. if len(d.lastThresholds) == 0 {
  180. return nil
  181. }
  182. return append([]float64(nil), d.lastThresholds...)
  183. }
  184. func (d *Detector) LastNoiseFloor() float64 {
  185. return d.lastNoiseFloor
  186. }
  187. func (ev *activeEvent) updateClass(newCls *classifier.Classification, historySize int, switchRatio float64) {
  188. if newCls == nil {
  189. return
  190. }
  191. if historySize <= 0 {
  192. historySize = 10
  193. }
  194. if switchRatio <= 0 || switchRatio > 1 {
  195. switchRatio = 0.6
  196. }
  197. if len(ev.classHistory) != historySize {
  198. ev.classHistory = make([]classifier.SignalClass, historySize)
  199. ev.classIdx = 0
  200. }
  201. ev.classHistory[ev.classIdx%len(ev.classHistory)] = newCls.ModType
  202. ev.classIdx++
  203. if ev.class == nil {
  204. clone := *newCls
  205. ev.class = &clone
  206. return
  207. }
  208. counts := map[classifier.SignalClass]int{}
  209. filled := ev.classIdx
  210. if filled > len(ev.classHistory) {
  211. filled = len(ev.classHistory)
  212. }
  213. for i := 0; i < filled; i++ {
  214. c := ev.classHistory[i]
  215. if c != "" {
  216. counts[c]++
  217. }
  218. }
  219. var majority classifier.SignalClass
  220. majorityCount := 0
  221. for c, n := range counts {
  222. if n > majorityCount {
  223. majority = c
  224. majorityCount = n
  225. }
  226. }
  227. threshold := int(math.Ceil(float64(filled) * switchRatio))
  228. if threshold < 1 {
  229. threshold = 1
  230. }
  231. if majorityCount >= threshold && majority != ev.class.ModType {
  232. clone := *newCls
  233. clone.ModType = majority
  234. ev.class = &clone
  235. } else if majority == ev.class.ModType && newCls.Confidence > ev.class.Confidence {
  236. ev.class.Confidence = newCls.Confidence
  237. ev.class.Features = newCls.Features
  238. ev.class.SecondBest = newCls.SecondBest
  239. ev.class.Scores = newCls.Scores
  240. }
  241. }
  242. func (d *Detector) UpdateClasses(signals []Signal) {
  243. for _, s := range signals {
  244. for _, ev := range d.active {
  245. if overlapHz(s.CenterHz, s.BWHz, ev.centerHz, ev.bwHz) && math.Abs(s.CenterHz-ev.centerHz) < (s.BWHz+ev.bwHz)/2.0 {
  246. if s.Class != nil {
  247. ev.updateClass(s.Class, d.classHistorySize, d.classSwitchRatio)
  248. }
  249. }
  250. }
  251. }
  252. }
  253. func (d *Detector) detectSignals(spectrum []float64, centerHz float64) []Signal {
  254. n := len(spectrum)
  255. if n == 0 {
  256. return nil
  257. }
  258. smooth := d.smoothSpectrum(spectrum)
  259. var thresholds []float64
  260. if d.cfarEngine != nil {
  261. thresholds = d.cfarEngine.Thresholds(smooth)
  262. }
  263. d.lastThresholds = append(d.lastThresholds[:0], thresholds...)
  264. noiseGlobal := median(smooth)
  265. d.lastNoiseFloor = noiseGlobal
  266. var signals []Signal
  267. in := false
  268. start := 0
  269. peak := -1e9
  270. peakBin := 0
  271. for i := 0; i < n; i++ {
  272. v := smooth[i]
  273. thresholdOn := d.ThresholdDb
  274. if thresholds != nil && !math.IsNaN(thresholds[i]) {
  275. thresholdOn = thresholds[i]
  276. }
  277. thresholdOff := thresholdOn - d.HysteresisDb
  278. if v >= thresholdOn {
  279. if !in {
  280. in = true
  281. start = i
  282. peak = v
  283. peakBin = i
  284. } else if v > peak {
  285. peak = v
  286. peakBin = i
  287. }
  288. } else if in && v < thresholdOff {
  289. noise := noiseGlobal
  290. if thresholds != nil && peakBin >= 0 && peakBin < len(thresholds) && !math.IsNaN(thresholds[peakBin]) {
  291. noise = thresholds[peakBin] - d.CFARScaleDb
  292. }
  293. signals = append(signals, d.makeSignal(start, i-1, peak, peakBin, noise, centerHz, smooth))
  294. in = false
  295. }
  296. }
  297. if in {
  298. noise := noiseGlobal
  299. if thresholds != nil && peakBin >= 0 && peakBin < len(thresholds) && !math.IsNaN(thresholds[peakBin]) {
  300. noise = thresholds[peakBin] - d.CFARScaleDb
  301. }
  302. signals = append(signals, d.makeSignal(start, n-1, peak, peakBin, noise, centerHz, smooth))
  303. }
  304. signals = d.expandSignalEdges(signals, smooth, noiseGlobal, centerHz)
  305. for i := range signals {
  306. centerBin := float64(signals[i].FirstBin+signals[i].LastBin) / 2.0
  307. signals[i].CenterHz = d.centerFreqForBin(centerBin, centerHz)
  308. signals[i].BWHz = float64(signals[i].LastBin-signals[i].FirstBin+1) * d.binWidth
  309. }
  310. return signals
  311. }
  312. func (d *Detector) expandSignalEdges(signals []Signal, smooth []float64, noiseFloor float64, centerHz float64) []Signal {
  313. n := len(smooth)
  314. if n == 0 || len(signals) == 0 {
  315. return signals
  316. }
  317. margin := d.EdgeMarginDb
  318. if margin <= 0 {
  319. margin = 3.0
  320. }
  321. maxExpansionBins := int(d.MaxSignalBwHz / d.binWidth)
  322. if maxExpansionBins < 10 {
  323. maxExpansionBins = 10
  324. }
  325. for i := range signals {
  326. seed := signals[i]
  327. peakDb := seed.PeakDb
  328. localNoise := noiseFloor
  329. leftProbe := seed.FirstBin - 50
  330. rightProbe := seed.LastBin + 50
  331. if leftProbe >= 0 && rightProbe < n {
  332. leftNoise := minInRange(smooth, maxInt(0, leftProbe), maxInt(0, seed.FirstBin-5))
  333. rightNoise := minInRange(smooth, minInt(n-1, seed.LastBin+5), minInt(n-1, rightProbe))
  334. localNoise = math.Min(leftNoise, rightNoise)
  335. }
  336. edgeThreshold := localNoise + margin
  337. newFirst := seed.FirstBin
  338. prevVal := smooth[newFirst]
  339. for j := 0; j < maxExpansionBins; j++ {
  340. next := newFirst - 1
  341. if next < 0 {
  342. break
  343. }
  344. val := smooth[next]
  345. if val <= edgeThreshold {
  346. break
  347. }
  348. if val > prevVal+1.0 && val < peakDb-6.0 {
  349. break
  350. }
  351. prevVal = val
  352. newFirst = next
  353. }
  354. newLast := seed.LastBin
  355. prevVal = smooth[newLast]
  356. for j := 0; j < maxExpansionBins; j++ {
  357. next := newLast + 1
  358. if next >= n {
  359. break
  360. }
  361. val := smooth[next]
  362. if val <= edgeThreshold {
  363. break
  364. }
  365. if val > prevVal+1.0 && val < peakDb-6.0 {
  366. break
  367. }
  368. prevVal = val
  369. newLast = next
  370. }
  371. signals[i].FirstBin = newFirst
  372. signals[i].LastBin = newLast
  373. centerBin := float64(newFirst+newLast) / 2.0
  374. signals[i].CenterHz = d.centerFreqForBin(centerBin, centerHz)
  375. signals[i].BWHz = float64(newLast-newFirst+1) * d.binWidth
  376. }
  377. signals = d.mergeOverlapping(signals, centerHz)
  378. return signals
  379. }
  380. func (d *Detector) mergeOverlapping(signals []Signal, centerHz float64) []Signal {
  381. if len(signals) <= 1 {
  382. return signals
  383. }
  384. gapBins := 0
  385. if d.MergeGapHz > 0 && d.binWidth > 0 {
  386. gapBins = int(math.Ceil(d.MergeGapHz / d.binWidth))
  387. }
  388. sort.Slice(signals, func(i, j int) bool {
  389. return signals[i].FirstBin < signals[j].FirstBin
  390. })
  391. merged := []Signal{signals[0]}
  392. for i := 1; i < len(signals); i++ {
  393. last := &merged[len(merged)-1]
  394. cur := signals[i]
  395. gap := cur.FirstBin - last.LastBin - 1
  396. if gap <= gapBins {
  397. if cur.LastBin > last.LastBin {
  398. last.LastBin = cur.LastBin
  399. }
  400. if cur.PeakDb > last.PeakDb {
  401. last.PeakDb = cur.PeakDb
  402. }
  403. if cur.SNRDb > last.SNRDb {
  404. last.SNRDb = cur.SNRDb
  405. }
  406. centerBin := float64(last.FirstBin+last.LastBin) / 2.0
  407. last.BWHz = float64(last.LastBin-last.FirstBin+1) * d.binWidth
  408. last.CenterHz = d.centerFreqForBin(centerBin, centerHz)
  409. if cur.NoiseDb < last.NoiseDb || last.NoiseDb == 0 {
  410. last.NoiseDb = cur.NoiseDb
  411. }
  412. } else {
  413. merged = append(merged, cur)
  414. }
  415. }
  416. return merged
  417. }
  418. func (d *Detector) centerFreqForBin(bin float64, centerHz float64) float64 {
  419. return centerHz + (bin-float64(d.nbins)/2.0)*d.binWidth
  420. }
  421. func minInRange(s []float64, from, to int) float64 {
  422. if len(s) == 0 {
  423. return 0
  424. }
  425. if from < 0 {
  426. from = 0
  427. }
  428. if to >= len(s) {
  429. to = len(s) - 1
  430. }
  431. if from > to {
  432. return s[minInt(maxInt(from, 0), len(s)-1)]
  433. }
  434. m := s[from]
  435. for i := from + 1; i <= to; i++ {
  436. if s[i] < m {
  437. m = s[i]
  438. }
  439. }
  440. return m
  441. }
  442. func maxInt(a, b int) int {
  443. if a > b {
  444. return a
  445. }
  446. return b
  447. }
  448. func minInt(a, b int) int {
  449. if a < b {
  450. return a
  451. }
  452. return b
  453. }
  454. func overlapHz(center1, bw1, center2, bw2 float64) bool {
  455. left1 := center1 - bw1/2
  456. right1 := center1 + bw1/2
  457. left2 := center2 - bw2/2
  458. right2 := center2 + bw2/2
  459. return left1 <= right2 && left2 <= right1
  460. }
  461. func median(vals []float64) float64 {
  462. if len(vals) == 0 {
  463. return 0
  464. }
  465. cp := append([]float64(nil), vals...)
  466. sort.Float64s(cp)
  467. mid := len(cp) / 2
  468. if len(cp)%2 == 0 {
  469. return (cp[mid-1] + cp[mid]) / 2
  470. }
  471. return cp[mid]
  472. }
  473. func (d *Detector) makeSignal(first, last int, peak float64, peakBin int, noise float64, centerHz float64, spectrum []float64) Signal {
  474. centerBin := float64(first+last) / 2.0
  475. centerFreq := centerHz + (centerBin-float64(d.nbins)/2.0)*d.binWidth
  476. bw := float64(last-first+1) * d.binWidth
  477. snr := peak - noise
  478. return Signal{
  479. FirstBin: first,
  480. LastBin: last,
  481. CenterHz: centerFreq,
  482. BWHz: bw,
  483. PeakDb: peak,
  484. SNRDb: snr,
  485. NoiseDb: noise,
  486. }
  487. }
  488. func (d *Detector) smoothSpectrum(spectrum []float64) []float64 {
  489. if d.ema == nil || len(d.ema) != len(spectrum) {
  490. d.ema = make([]float64, len(spectrum))
  491. copy(d.ema, spectrum)
  492. return d.ema
  493. }
  494. alpha := d.EmaAlpha
  495. for i := range spectrum {
  496. v := spectrum[i]
  497. d.ema[i] = alpha*v + (1-alpha)*d.ema[i]
  498. }
  499. return d.ema
  500. }
  501. func (d *Detector) matchSignals(now time.Time, signals []Signal) []Event {
  502. used := make(map[int64]bool, len(d.active))
  503. for _, s := range signals {
  504. var best *activeEvent
  505. var candidates []struct {
  506. ev *activeEvent
  507. dist float64
  508. }
  509. for _, ev := range d.active {
  510. if overlapHz(s.CenterHz, s.BWHz, ev.centerHz, ev.bwHz) && math.Abs(s.CenterHz-ev.centerHz) < (s.BWHz+ev.bwHz)/2.0 {
  511. candidates = append(candidates, struct {
  512. ev *activeEvent
  513. dist float64
  514. }{ev: ev, dist: math.Abs(s.CenterHz - ev.centerHz)})
  515. }
  516. }
  517. if len(candidates) > 0 {
  518. sort.Slice(candidates, func(i, j int) bool { return candidates[i].dist < candidates[j].dist })
  519. best = candidates[0].ev
  520. }
  521. if best == nil {
  522. id := d.nextID
  523. d.nextID++
  524. d.active[id] = &activeEvent{
  525. id: id,
  526. start: now,
  527. lastSeen: now,
  528. centerHz: s.CenterHz,
  529. bwHz: s.BWHz,
  530. peakDb: s.PeakDb,
  531. snrDb: s.SNRDb,
  532. firstBin: s.FirstBin,
  533. lastBin: s.LastBin,
  534. class: s.Class,
  535. stableHits: 1,
  536. }
  537. continue
  538. }
  539. used[best.id] = true
  540. best.lastSeen = now
  541. best.stableHits++
  542. best.centerHz = (best.centerHz + s.CenterHz) / 2.0
  543. if best.bwHz <= 0 {
  544. best.bwHz = s.BWHz
  545. } else {
  546. const alpha = 0.15
  547. best.bwHz = alpha*s.BWHz + (1-alpha)*best.bwHz
  548. }
  549. if s.PeakDb > best.peakDb {
  550. best.peakDb = s.PeakDb
  551. }
  552. if s.SNRDb > best.snrDb {
  553. best.snrDb = s.SNRDb
  554. }
  555. if s.FirstBin < best.firstBin {
  556. best.firstBin = s.FirstBin
  557. }
  558. if s.LastBin > best.lastBin {
  559. best.lastBin = s.LastBin
  560. }
  561. if s.Class != nil {
  562. best.updateClass(s.Class, d.classHistorySize, d.classSwitchRatio)
  563. }
  564. }
  565. var finished []Event
  566. for id, ev := range d.active {
  567. if used[id] {
  568. continue
  569. }
  570. if now.Sub(ev.lastSeen) < d.GapTolerance {
  571. continue
  572. }
  573. duration := ev.lastSeen.Sub(ev.start)
  574. if duration < d.MinDuration || ev.stableHits < d.MinStableFrames {
  575. delete(d.active, id)
  576. continue
  577. }
  578. finished = append(finished, Event{
  579. ID: ev.id,
  580. Start: ev.start,
  581. End: ev.lastSeen,
  582. CenterHz: ev.centerHz,
  583. Bandwidth: ev.bwHz,
  584. PeakDb: ev.peakDb,
  585. SNRDb: ev.snrDb,
  586. FirstBin: ev.firstBin,
  587. LastBin: ev.lastBin,
  588. Class: ev.class,
  589. })
  590. delete(d.active, id)
  591. }
  592. return finished
  593. }