[Commit-gnuradio] [gnuradio] 03/07: gr-analog: Add safety and default for FM preemphasis filter

[git]

This is an automated email from the git hooks/post-receive script.

jcorgan pushed a commit to branch master
in repository gnuradio.

commit 515e40039e614ac86fb4c7be8dc4d067aa2bd9a3
Author: Andy Walls <address@hidden>
Date:   Wed Mar 30 09:15:45 2016 -0400

    gr-analog: Add safety and default for FM preemphasis filter
    
    The FM preemphasis filter design now precludes the user from
    inducing a pole on the unit circle at z = -1.0 and z = 1.0.
    A pole at either of these locations makes the filter unstable and
    useless: feeding back "+/-inf" into an IIR filter has no good
    recovery.
    
    Also provide a reasonable, maximally safe default of 0.925*fs/2.0
    for the high frequency corner, fh.  This keeps the slope of the
    preemphasis filter looking reasonable sane in the whole band; at
    least for tau=75e-6 and fs=48000.
---
 docs/exploring-gnuradio/fm_tx.grc   |  2 +-
 gr-analog/examples/fmtest.py        |  4 +-
 gr-analog/grc/analog_fm_preemph.xml |  2 +-
 gr-analog/grc/analog_nbfm_tx.xml    |  2 +-
 gr-analog/grc/analog_wfm_tx.xml     |  2 +-
 gr-analog/python/analog/fm_emph.py  | 79 +++++++++++++++----------------------
 gr-analog/python/analog/nbfm_tx.py  |  4 +-
 gr-analog/python/analog/wfm_tx.py   |  4 +-
 gr-filter/examples/synth_to_chan.py |  2 +-
 gr-uhd/examples/python/fm_tx4.py    |  4 +-
 10 files changed, 44 insertions(+), 61 deletions(-)

diff --git a/docs/exploring-gnuradio/fm_tx.grc 
b/docs/exploring-gnuradio/fm_tx.grc
index be8fe77..bb13417 100644
--- a/docs/exploring-gnuradio/fm_tx.grc
+++ b/docs/exploring-gnuradio/fm_tx.grc
@@ -793,7 +793,7 @@
     </param>
     <param>
       <key>fh</key>
-      <value>0.0</value>
+      <value>0.925 * tx_rate/2.0</value>
     </param>
     <param>
       <key>affinity</key>
diff --git a/gr-analog/examples/fmtest.py b/gr-analog/examples/fmtest.py
index 22448e6..7ed08ca 100755
--- a/gr-analog/examples/fmtest.py
+++ b/gr-analog/examples/fmtest.py
@@ -48,8 +48,8 @@ class fmtx(gr.hier_block2):
                                 gr.io_signature(1, 1, gr.sizeof_float),
                                 gr.io_signature(1, 1, gr.sizeof_gr_complex))
 
-        fmtx = analog.nbfm_tx(audio_rate, if_rate, max_dev=5e3, tau=75e-6,
-                                                               fh=0.0)
+        fmtx = analog.nbfm_tx(audio_rate, if_rate, max_dev=5e3,
+                             tau=75e-6, fh=0.925*if_rate/2.0)
 
         # Local oscillator
         lo = analog.sig_source_c(if_rate,            # sample rate
diff --git a/gr-analog/grc/analog_fm_preemph.xml 
b/gr-analog/grc/analog_fm_preemph.xml
index e5a8c35..0e8928c 100644
--- a/gr-analog/grc/analog_fm_preemph.xml
+++ b/gr-analog/grc/analog_fm_preemph.xml
@@ -23,7 +23,7 @@
        <param>
                <name>High Corner Freq</name>
                <key>fh</key>
-               <value>0.0</value>
+               <value>0.925*samp_rate/2.0</value>
                <type>real</type>
        </param>
        <sink>
diff --git a/gr-analog/grc/analog_nbfm_tx.xml b/gr-analog/grc/analog_nbfm_tx.xml
index 3e91e9c..d496a18 100644
--- a/gr-analog/grc/analog_nbfm_tx.xml
+++ b/gr-analog/grc/analog_nbfm_tx.xml
@@ -46,7 +46,7 @@
   <param>
     <name>Preemphasis High Corner Freq</name>
     <key>fh</key>
-    <value>0.0</value>
+    <value>0.925*float(quad_rate)/2.0</value>
     <type>real</type>
   </param>
 
diff --git a/gr-analog/grc/analog_wfm_tx.xml b/gr-analog/grc/analog_wfm_tx.xml
index cd67d6a..dccfeb0 100644
--- a/gr-analog/grc/analog_wfm_tx.xml
+++ b/gr-analog/grc/analog_wfm_tx.xml
@@ -40,7 +40,7 @@
        <param>
                <name>Preemphasis High Corner Freq</name>
                <key>fh</key>
-               <value>0.0</value>
+               <value>0.925*float(quad_rate)/2.0</value>
                <type>real</type>
        </param>
        <check>($quad_rate)%($audio_rate) == 0</check>
diff --git a/gr-analog/python/analog/fm_emph.py 
b/gr-analog/python/analog/fm_emph.py
index 7637743..bfa4742 100644
--- a/gr-analog/python/analog/fm_emph.py
+++ b/gr-analog/python/analog/fm_emph.py
@@ -251,69 +251,52 @@ class fm_preemph(gr.hier_block2):
     """
     FM Preemphasis IIR filter.
     """
-    def __init__(self, fs, tau=75e-6, fh=0.0):
+    def __init__(self, fs, tau=75e-6, fh=-1.0):
         """
 
         Args:
             fs: sampling frequency in Hz (float)
             tau: Time constant in seconds (75us in US, 50us in EUR) (float)
-            fh: High frequency at which to flatten out; 0.0 means none (float)
+            fh: High frequency at which to flatten out (< 0 means default of 
0.925*fs/2.0) (float)
         """
         gr.hier_block2.__init__(self, "fm_preemph",
                                 gr.io_signature(1, 1, gr.sizeof_float),  # 
Input signature
                                 gr.io_signature(1, 1, gr.sizeof_float))  # 
Output signature
 
-       if fh > 0.0 and fh < (fs / 2.0):
-               # Digital corner frequencies
-               w_cl = 1.0 / tau
-               w_ch = 2.0 * math.pi * fh
+       # Set fh to something sensible, if needed.
+       # N.B. fh == fs/2.0 or fh == 0.0 results in a pole on the unit circle
+       # at z = -1.0 or z = 1.0 respectively.  That makes the filter unstable
+       # and useless.
+       if fh <= 0.0 or fh >= fs/2.0:
+               fh = 0.925 * fs/2.0
 
-               # Prewarped analog corner frequencies
-               w_cla = 2.0 * fs * math.tan(w_cl / (2.0 * fs))
-               w_cha = 2.0 * fs * math.tan(w_ch / (2.0 * fs))
+       # Digital corner frequencies
+       w_cl = 1.0 / tau
+       w_ch = 2.0 * math.pi * fh
 
-               # Resulting digital pole, zero, and gain term from the bilinear
-               # transformation of H(s) = (s + w_cla) / (s + w_cha) to
-               # H(z) = b0 (1 - z1 z^-1)/(1 - p1 z^-1)
-               kl = -w_cla / (2.0 * fs)
-               kh = -w_cha / (2.0 * fs)
-               z1 = (1.0 + kl) / (1.0 - kl)
-               p1 = (1.0 + kh) / (1.0 - kh)
-               b0 = (1.0 - kl) / (1.0 - kh)
+       # Prewarped analog corner frequencies
+       w_cla = 2.0 * fs * math.tan(w_cl / (2.0 * fs))
+       w_cha = 2.0 * fs * math.tan(w_ch / (2.0 * fs))
 
-               # Since H(s = infinity) = 1.0, then H(z = -1) = 1.0 and
-               # this filter  has 0 dB gain at fs/2.0.
-               # That isn't what users are going to expect, so adjust with a
-               # gain, g, so that H(z = 1) = 1.0 for 0 dB gain at DC.
-               w_0dB = 2.0 * math.pi * 0.0
-               g =        abs(1.0 - p1 * cmath.rect(1.0, -w_0dB))  \
-                  / (b0 * abs(1.0 - z1 * cmath.rect(1.0, -w_0dB)))
+       # Resulting digital pole, zero, and gain term from the bilinear
+       # transformation of H(s) = (s + w_cla) / (s + w_cha) to
+       # H(z) = b0 (1 - z1 z^-1)/(1 - p1 z^-1)
+       kl = -w_cla / (2.0 * fs)
+       kh = -w_cha / (2.0 * fs)
+       z1 = (1.0 + kl) / (1.0 - kl)
+       p1 = (1.0 + kh) / (1.0 - kh)
+       b0 = (1.0 - kl) / (1.0 - kh)
 
-               btaps = [ g * b0 * 1.0, g * b0 * -z1 ]
-               ataps = [          1.0,          -p1 ]
+       # Since H(s = infinity) = 1.0, then H(z = -1) = 1.0 and
+       # this filter  has 0 dB gain at fs/2.0.
+       # That isn't what users are going to expect, so adjust with a
+       # gain, g, so that H(z = 1) = 1.0 for 0 dB gain at DC.
+       w_0dB = 2.0 * math.pi * 0.0
+       g =        abs(1.0 - p1 * cmath.rect(1.0, -w_0dB))  \
+          / (b0 * abs(1.0 - z1 * cmath.rect(1.0, -w_0dB)))
 
-       else:
-               # Just use H(s) = (s + 1/RC)/(1/RC) as the transfer function
-
-               # Digital corner frequencies
-               w_cl = 1.0 / tau
-
-               # Prewarped analog corner frequencies
-               w_cla = 2.0 * fs * math.tan(w_cl / (2.0 * fs))
-
-               # Resulting digital pole, zero, and gain term from the bilinear
-               # transformation of H(s) = (s + w_cl)/w_cl to
-               # H(z) = b0 (1 - z1 z^-1)/(1 - p1 z^-1)
-               kl = -w_cla / (2.0 * fs)
-               z1 = (1.0 + kl) / (1.0 - kl)
-               p1 = -1.0
-               b0 = (1.0 - kl) / -kl
-
-               # Since H(s = 0) = 1.0, then H(z = 1) = 1.0 and
-               # has 0 dB gain at DC
-
-               btaps = [ b0 * 1.0, b0 * -z1 ]
-               ataps = [      1.0,      -p1 ]
+       btaps = [ g * b0 * 1.0, g * b0 * -z1 ]
+       ataps = [          1.0,          -p1 ]
 
         if 0:
             print "btaps =", btaps
diff --git a/gr-analog/python/analog/nbfm_tx.py 
b/gr-analog/python/analog/nbfm_tx.py
index cd11c2f..aa6c1ec 100644
--- a/gr-analog/python/analog/nbfm_tx.py
+++ b/gr-analog/python/analog/nbfm_tx.py
@@ -29,7 +29,7 @@ except ImportError:
     import analog_swig as analog
 
 class nbfm_tx(gr.hier_block2):
-    def __init__(self, audio_rate, quad_rate, tau=75e-6, max_dev=5e3, fh=0.0):
+    def __init__(self, audio_rate, quad_rate, tau=75e-6, max_dev=5e3, fh=-1.0):
         """
         Narrow Band FM Transmitter.
 
@@ -41,7 +41,7 @@ class nbfm_tx(gr.hier_block2):
             quad_rate: sample rate of output stream (integer)
             tau: preemphasis time constant (default 75e-6) (float)
             max_dev: maximum deviation in Hz (default 5e3) (float)
-            fh: high frequency at which to flatten preemphasis; 0.0 means none 
(float)
+            fh: high frequency at which to flatten preemphasis; < 0 means 
default of 0.925*quad_rate/2.0 (float)
 
         quad_rate must be an integer multiple of audio_rate.
         """
diff --git a/gr-analog/python/analog/wfm_tx.py 
b/gr-analog/python/analog/wfm_tx.py
index 7363ccd..a1b5893 100644
--- a/gr-analog/python/analog/wfm_tx.py
+++ b/gr-analog/python/analog/wfm_tx.py
@@ -30,7 +30,7 @@ except ImportError:
     import analog_swig as analog
 
 class wfm_tx(gr.hier_block2):
-    def __init__(self, audio_rate, quad_rate, tau=75e-6, max_dev=75e3, fh=0.0):
+    def __init__(self, audio_rate, quad_rate, tau=75e-6, max_dev=75e3, 
fh=-1.0):
         """
         Wide Band FM Transmitter.
 
@@ -42,7 +42,7 @@ class wfm_tx(gr.hier_block2):
             quad_rate: sample rate of output stream (integer)
             tau: preemphasis time constant (default 75e-6) (float)
             max_dev: maximum deviation in Hz (default 75e3) (float)
-            fh: high frequency at which to flatten preemphasis; 0.0 means none 
(float)
+            fh: high frequency at which to flatten preemphasis; < 0 means 
default of 0.925*quad_rate/2.0 (float)
 
         quad_rate must be an integer multiple of audio_rate.
         """
diff --git a/gr-filter/examples/synth_to_chan.py 
b/gr-filter/examples/synth_to_chan.py
index 7a295ea..88fb080 100755
--- a/gr-filter/examples/synth_to_chan.py
+++ b/gr-filter/examples/synth_to_chan.py
@@ -54,7 +54,7 @@ def main():
     fmtx = list()
     for fi in freqs:
         s = analog.sig_source_f(fs, analog.GR_SIN_WAVE, fi, 1)
-        fm = analog.nbfm_tx(fs, 4*fs, max_dev=10000, tau=75e-6, fh=0.0)
+        fm = analog.nbfm_tx(fs, 4*fs, max_dev=10000, tau=75e-6, 
fh=0.925*(4*fs)/2.0)
         sigs.append(s)
         fmtx.append(fm)
 
diff --git a/gr-uhd/examples/python/fm_tx4.py b/gr-uhd/examples/python/fm_tx4.py
index ffc74ab..516033d 100755
--- a/gr-uhd/examples/python/fm_tx4.py
+++ b/gr-uhd/examples/python/fm_tx4.py
@@ -63,8 +63,8 @@ class pipeline(gr.hier_block2):
             sys.exit(1)
 
         print audio_rate, if_rate
-        fmtx = analog.nbfm_tx(audio_rate, if_rate, max_dev=5e3, tau=75e-6,
-                                                               fh=0.0)
+        fmtx = analog.nbfm_tx(audio_rate, if_rate, max_dev=5e3,
+                             tau=75e-6, fh=0.925*if_rate/2.0)
 
         # Local oscillator
         lo = analog.sig_source_c(if_rate,            # sample rate