Further WIP

python/.gitignore

@@ -0,0 +1,5 @@
+.ipynb_checkpoints/
+__pycache__/
+build/
+fir.py
+pyfir_wrap.cxx

python/fir++.h

@@ -0,0 +1 @@
+../obj/fir++.h

python/pyfir.i

@@ -0,0 +1,24 @@
+%module fir
+%include <complex.i>
+
+/*
+ * List conversion
+ *
+ *  SWIG knows how to turn std::vector<T> func() into a function that returns a
+ *  Python list. We do need to tell it about all <T> types with %template,
+ *  however.
+ */
+%include <std_vector.i>
+%template(vectord) std::vector<double>;
+
+
+
+%{
+#include "fir++.h"
+%}
+
+
+
+%include "fir++.h"
+
+// vim: ts=4:sw=4

python/setup.py

@@ -0,0 +1,16 @@
+#!/usr/bin/env python3
+
+from distutils.core import setup, Extension
+
+fir_module = Extension('_fir',
+	sources=['pyfir_wrap.cxx'],
+	libraries=['fir++'],
+	)
+
+setup(name = 'fir',
+	version = '0.0.0',
+	author = 'Laurence Withers',
+	description = 'FIR filter',
+	ext_modules = [fir_module],
+	py_modules = ['fir'],
+	)

src/libfir++/000_TopHeader.h

@@ -10,6 +10,7 @@
 
 // standard includes, or includes needed for type declarations
 #include <vector>
+#include <complex>
 
 // forward declares, allowing us to avoid including the C library header
 extern "C" {

src/libfir++/100_Filter.cc

@@ -32,6 +32,13 @@ Filter::Filter(const std::vector<double>& coeff, Symmetry symm)
 
 
 
+Filter::Filter(struct fir_filter_t* fi)
+{
+    fi_ = fi;
+}
+
+
+
 Filter::~Filter()
 {
     fir_filter_free(fi_);
@@ -71,6 +78,14 @@ Filter::coeff(int idx) const
 
 
 
+std::complex<double>
+Filter::response(double freq) const
+{
+    return std::complex<double>(fir_filter_response(fi_, freq));
+}
+
+
+
 /* options for text editors
 vim: expandtab:ts=4:sw=4
 */

src/libfir++/100_Filter.h

@@ -39,7 +39,7 @@ enum Symmetry {
 
 
 
-/*! \brief Basic FIR filter object.
+/*! \brief Basic %FIR filter object.
 
 \ingroup libfirpp_fir
 
@@ -88,12 +88,25 @@ public:
 
 
 
+    /*! \brief Take ownership of underlying C filter object.
+
+    \param fi Filter.
+
+    This constructor may be used to build a filter object using an underlying C
+    library filter object. The new instance will take ownership of \a fi and
+    free it when the instance is deleted.
+
+    */
+    explicit Filter(struct fir_filter_t* fi);
+
+
+
     /*! \brief Destructor (frees underlying memory). */
     virtual ~Filter();
 
 
 
-    /*! \brief Filter a sample.
+    /*! \brief %Filter a sample.
 
     \param x Input sample.
     \returns Filter output.
@@ -141,6 +154,20 @@ public:
 
 
 
+    /*! \brief Evaluate response of filter at given frequency.
+
+    \param freq Frequency as a fraction of sampling rate.
+    \returns Complex response of filter at given frequency.
+
+    Performs Z-domain analysis of the filter to determine its complex response
+    (phase and magnitude) to a signal of the given frequency \a freq, which is
+    expressed as a fraction of the sampling rate being used with the filter.
+
+    */
+    std::complex<double> response(double freq) const;
+
+
+
 private:
     struct fir_filter_t* fi_;
 

src/libfir++/200_CommonFilters.cc

@@ -0,0 +1,28 @@
+/* libfir/src/libfir++/200_CommonFilters.cc
+ * 
+ *  Copyright: ©2014, Laurence Withers.
+ *  Author: Laurence Withers <l@lwithers.me.uk>
+ *  License: GPLv3
+ */
+
+
+
+Filter
+SincLowpass(int npoints, double corner, double gain)
+{
+    return Filter(fir_sinc_lowpass(npoints, corner, gain));
+}
+
+
+
+Filter
+LanczosLowpass(int npoints, double corner, double gain, int a)
+{
+    return Filter(fir_lanczos_lowpass(npoints, corner, gain, a));
+}
+
+
+
+/* options for text editors
+vim: expandtab:ts=4:sw=4
+*/

src/libfir++/200_CommonFilters.h

@@ -0,0 +1,48 @@
+/* libfir/src/libfir++/200_CommonFilters.h
+ * 
+ *  Copyright: ©2014, Laurence Withers.
+ *  Author: Laurence Withers <l@lwithers.me.uk>
+ *  License: GPLv3
+ */
+
+
+
+/*! \defgroup libfirpp_common Common filter types
+
+\ingroup libfirpp
+
+*/
+/*!@{*/
+
+
+
+/*! \brief Brickwall low-pass linear phase: sinc filter.
+
+\param npoints Number of points required (≥1).
+\param corner Corner frequency as fraction of sample rate.
+\param gain Gain of filter.
+\returns Newly-allocated %FIR filter.
+
+Builds a truncated-sinc filter which has a low-pass, brickwall effect (cutoff
+at \a corner, which is a fraction of sampling rate, i.e. 0 &lt; \a corner &lt;
+0.5).
+
+The number of points \a npoints is a trade-off between accuracy of filter (an
+ideal sinc filter would have an infinite number of points) and computation time.
+\a npoints may be odd or even; symmetry will be exploited in either case.
+
+The DC gain of the filter is normalised to \a gain, which is often 1.0.
+
+*/
+Filter SincLowpass(int npoints, double corner, double gain);
+
+
+
+Filter LanczosLowpass(int npoints, double corner, double gain, int a);
+
+
+
+/*!@}*/
+/* options for text editors
+vim: expandtab:ts=4:sw=4:syntax=cpp.doxygen
+*/

src/libfir/000_TopHeader.h

@@ -8,14 +8,15 @@
 #ifndef HEADER_libfir
 #define HEADER_libfir
 
+/* standard includes, or includes needed for type declarations */
+#include <complex.h>
+
 #ifdef __cplusplus
 extern "C" {
 #endif
 
 /*! \defgroup libfir C library */
 
-/* standard includes, or includes needed for type declarations */
-
 /* options for text editors
 vim: expandtab:ts=4:sw=4:syntax=ch.doxygen
 */

src/libfir/200_common_filters.c

@@ -29,7 +29,7 @@ fir_sinc_lowpass(int npoints, double corner, double gain)
             sum += 1; /* not 2, as we don't repeat this coeff */
         } else {
             t = i + 0.5 - npoints * 0.5;
-            t *= corner * 2;
+            t *= corner * 2 * M_PI;
             coeff[i] = sin(t) / t;
             sum += coeff[i] * 2;
         }
@@ -47,6 +47,46 @@ fir_sinc_lowpass(int npoints, double corner, double gain)
 
 
 
+struct fir_filter_t*
+fir_lanczos_lowpass(int npoints, double corner, double gain, int a)
+{
+    double* coeff, t, sum = 0;
+    int i, spoints;
+
+    if(npoints < 1) {
+        errno = EINVAL;
+        return 0;
+    }
+
+    /* find number of points after taking symmetry into account */
+    spoints = (npoints + 1) / 2;
+    coeff = alloca(sizeof(double) * spoints);
+
+    /* only fill out for t ≤ 0 */
+    for(i = 0; i < spoints; ++i) {
+        if((i * 2 + 1) == npoints) {
+            coeff[i] = 1;
+            sum += 1; /* not 2, as we don't repeat this coeff */
+        } else {
+            t = i + 0.5 - npoints * 0.5;
+            t *= corner * 2 * M_PI;
+            coeff[i] = (sin(t) / t) * sin(t / a) / (t / a);
+            sum += coeff[i] * 2;
+        }
+    }
+
+    /* normalise gain */
+    gain /= sum;
+    for(i = 0; i < spoints; ++i) {
+        coeff[i] *= gain;
+    }
+
+    return fir_filter_new(spoints, coeff,
+        (npoints & 1) ? fir_symmetry_odd : fir_symmetry_even);
+}
+
+
+
 /* options for text editors
 vim: expandtab:ts=4:sw=4
 */

src/libfir/200_common_filters.h

@@ -41,6 +41,10 @@ struct fir_filter_t* fir_sinc_lowpass(int npoints, double corner, double gain);
 
 
 
+struct fir_filter_t* fir_lanczos_lowpass(int npoints, double corner, double gain, int a);
+
+
+
 /*!@}*/
 /* options for text editors
 vim: expandtab:ts=4:sw=4:syntax=ch.doxygen

src/libfir/300_response.c

@@ -0,0 +1,28 @@
+/* libfir/src/libfir/300_response.c
+ *
+ *  Copyright: ©2014, Laurence Withers.
+ *  Author: Laurence Withers <l@lwithers.me.uk>
+ *  License: GPLv3
+ */
+
+
+
+double complex
+fir_filter_response(const struct fir_filter_t* fi, double freq)
+{
+    int i, ncoeff = fir_filter_get_ncoeff(fi);
+    double omega = 2 * M_PI * freq;
+    double complex resp = 0;
+
+    for(i = 0; i < ncoeff; ++i) {
+        resp += fir_filter_get_coeff(fi, i) * cexp(I * omega * i);
+    }
+
+    return resp;
+}
+
+
+
+/* options for text editors
+vim: expandtab:ts=4:sw=4:syntax=ch.doxygen
+*/

src/libfir/300_response.h

@@ -0,0 +1,39 @@
+/* libfir/src/libfir/300_response.h
+ *
+ *  Copyright: ©2014, Laurence Withers.
+ *  Author: Laurence Withers <l@lwithers.me.uk>
+ *  License: GPLv3
+ */
+
+
+
+/*! \defgroup libfir_response Response analysis
+
+\ingroup libfir
+
+Routines in this section allow the evaluation of the filter's response.
+
+*/
+/*!@{*/
+
+
+
+/*! \brief Measure phase and frequency response of filter at frequency.
+
+\param fi Filter to evaluate.
+\param freq Frequency (as a ratio of sample rate).
+\returns Complex response of filter at given frequency.
+
+Performs Z-domain analysis of the filter to determine the magnitude and phase of
+its response to a signal of the given frequency. \a freq is expressed as a
+fraction of the sample rate to be used with the filter.
+
+*/
+double complex fir_filter_response(const struct fir_filter_t* fi, double freq);
+
+
+
+/*!@}*/
+/* options for text editors
+vim: expandtab:ts=4:sw=4:syntax=ch.doxygen
+*/