Diff: BilinearDesignLpfHpf/BilinearDesignLH.cpp

Revision:

0:33908268d9ea

diff -r 000000000000 -r 33908268d9ea BilinearDesignLpfHpf/BilinearDesignLH.cpp
--- /dev/null	Thu Jan 01 00:00:00 1970 +0000
+++ b/BilinearDesignLpfHpf/BilinearDesignLH.cpp	Fri Sep 11 09:54:45 2015 +0000
@@ -0,0 +1,81 @@
+//------------------------------------------------------------------------------
+//  Design of Butterworth LPF and HPF using bilinear transform
+//
+//   2014/06/29, Copyright (c) 2014 MIKAMI, Naoki
+//------------------------------------------------------------------------------
+
+#include "BilinearDesignLH.hpp"
+
+namespace Mikami
+{
+    // Execute design
+    //      input
+    //          fc: Cutoff frequency
+    //      output
+    //          c : Coefficients for cascade structure
+    //          g : Gain factor for cascade structure
+    void BilinearDesign::Execute(float fc, Coefs c[], float& g)
+    {
+        Butterworth();
+        Bilinear(fc);
+        ToCascade();
+        GetGain();
+        GetCoefs(c, g);
+    }
+
+    // Get poles for Butterworth characteristics
+    void BilinearDesign::Butterworth()
+    {
+        float pi_2order = PI_/(2.0f*ORDER_);
+        for (int j=0; j<ORDER_/2; j++)  // Pole with imaginary part >= 0
+        {
+            float theta = (2.0f*j + 1.0f)*pi_2order;
+            sP_[j] = Complex(-cosf(theta), sinf(theta));
+        }
+    }
+    
+    // Bilinear transform
+    //      fc: Cutoff frequency
+    void BilinearDesign::Bilinear(float fc)
+    {
+        float wc = tanf(fc*PI_FS_);
+        for (int k=0; k<ORDER_/2; k++)
+            zP_[k] = (1.0f + wc*sP_[k])/(1.0f - wc*sP_[k]);
+    }
+    
+    // Convert to coefficients for cascade structure
+    void BilinearDesign::ToCascade()
+    {
+        for (int j=0; j<ORDER_/2; j++)
+        {
+            ck_[j].a1 = 2.0f*real(zP_[j]);              // a1m
+            ck_[j].a2 = -norm(zP_[j]);                  // a2m
+            ck_[j].b1 = (PB_ == LPF) ? 2.0f : -2.0f;    // b1m
+        }
+    }
+
+    // Calculate gain factor
+    void BilinearDesign::GetGain(){
+        float u = (PB_ == LPF) ? 1.0f : -1.0f;
+        float g0 = 1.0f;
+        for (int k=0; k<ORDER_/2; k++)
+            g0 = g0*(1.0f - (ck_[k].a1 + ck_[k].a2*u)*u)/
+                    (1.0f + (ck_[k].b1 + u)*u);
+        gain_ = g0;
+    }
+    
+    // Get coefficients
+    void BilinearDesign::GetCoefs(Coefs c[], float& gain)
+    {
+        for (int k=0; k<ORDER_/2; k++)
+        {
+            c[k].a1 = ck_[k].a1;
+            c[k].a2 = ck_[k].a2;
+            c[k].b1 = ck_[k].b1;
+        }
+        gain = gain_;
+    }
+}
+
+
+