1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
37
38
39
40
41
42
43
44
45
46
47
48
49
50
51
52
53
54
55
56
57
58
59
60
61
62
63
64
65
66
67
68
69
70
71
72
73
74
75
76
77
78
79
80
81
82
83
84
85

/*! Hilbert transform.

[Wikipedia][wiki] has a bunch of math, but one use case for it is to
convert floating point values (think audio waveform) into upper
sideband.

Then again I guess you can do the same with a FloatToComplex plus
FftFilter.

This implementation is a pretty inefficient.

[wiki]: https://en.wikipedia.org/wiki/Hilbert_transform
*/

use crate::block::{Block, BlockRet};
use crate::fir::FIR;
use crate::stream::{new_streamp, Streamp};
use crate::{Complex, Error, Float};

/// Hilbert transformer block.
pub struct Hilbert {
    src: Streamp<Float>,
    dst: Streamp<Complex>,
    history: Vec<Float>,
    filter: FIR<Float>,
    ntaps: usize,
}

impl Hilbert {
    /// Create new hilber transformer with this many taps.
    pub fn new(src: Streamp<Float>, ntaps: usize) -> Self {
        assert!(ntaps & 1 == 1, "hilbert filter len must be odd");
        let taps = crate::fir::hilbert(ntaps); // TODO: provide window function.
        Self {
            src,
            ntaps,
            dst: new_streamp(),
            history: vec![0.0; ntaps],
            filter: FIR::new(&taps),
        }
    }
    /// Get the output stream.
    pub fn out(&self) -> Streamp<Complex> {
        self.dst.clone()
    }
}

impl Block for Hilbert {
    fn block_name(&self) -> &'static str {
        "Hilbert"
    }
    fn work(&mut self) -> Result<BlockRet, Error> {
        assert_eq!(self.ntaps, self.history.len());
        let (i, tags) = self.src.read_buf()?;
        if i.is_empty() {
            return Ok(BlockRet::Noop);
        }
        let mut o = self.dst.write_buf()?;
        if o.is_empty() {
            return Ok(BlockRet::Ok);
        }

        let inout = std::cmp::min(i.len(), o.len());
        let len = self.history.len() + inout;
        let n = len - self.ntaps;

        // TODO: Probably needless copy.
        let mut iv = Vec::with_capacity(len);
        iv.extend(&self.history);
        iv.extend(i.iter().take(inout).copied());

        // I tried a couple of variations of this loop, and this was
        // the fastest on my laptop.
        for i in 0..n {
            let t = &iv[i..(i + self.ntaps)];
            o.slice()[i] = Complex::new(iv[i + self.ntaps / 2], self.filter.filter(t));
        }

        o.produce(n, &tags);

        self.history[..self.ntaps].clone_from_slice(&iv[n..len]);
        i.consume(n);
        Ok(BlockRet::Ok)
    }
}