AGC example

Author: PythonCHB

Examples/week-08-extensions/cython/AGC_example/README.rst

@@ -3,8 +3,24 @@ Automatic Gain Control Example
 
 Automatic Gain Control is a signal processing procedure:
 
-http://en.wikipedia.org/wiki/Automatic_gain_control
+http://seismicreflections.globeclaritas.com/2013/04/agc-equaliser.html
 
 This is an example of using fortran, C, and Cython to produce a fast AGC
 filter for Python/numpy.
 
+Thanks for the folks at IRIS (http://www.iris.edu) for this example and
+the Fortran code.
+
+Using these examples:
+=====================
+
+The c and cython versions should be usable with any pyton set up to
+compile extensions. They use setuptools to call the compiler.
+
+The Fortran / f2c version requires a compatible Fortran compiler. g77
+should work fine, but getting a compatible build on anything but Linux
+is a trick, and out of scope for this README (sorry). Google for it --
+you should find what you need.
+
+
+

Examples/week-08-extensions/cython/AGC_example/agc_c.c

@@ -1,6 +1,7 @@
 /* C version of AGC code */
 
 #include <stdio.h>
+#include <math.h>
 
 void AGC(int nAGC, int npts, float *amp, float *ampAGC);
 

Examples/week-08-extensions/cython/AGC_example/setup.py

@@ -1,24 +1,17 @@
 #!/usr/bin/env python
 
-from distutils.core import setup
+from setuptools import setup, Extension
 from Cython.Build import cythonize
-from distutils.extension import Extension
 
 import numpy
 
 cy = Extension("agc_cython",
-	           ["agc_cython.pyx",],
-		       include_dirs = [numpy.get_include(),],
-		       )
+               ["agc_cython.pyx"],
+               include_dirs=[numpy.get_include()],
+               )
 c_wrap = Extension("agc_c_cy",
-	               ["agc_c_cy.pyx", "agc_c.c"],
-	               include_dirs = [numpy.get_include(),],
-	               )
-
-setup(
-    # ext_modules = [cythonize("agc_cython.pyx"),
-    #                cythonize([c_wrap]),
-    #               ]
-    ext_modules = cythonize( [cy, c_wrap] ),
-)
+                   ["agc_c_cy.pyx", "agc_c.c"],
+                   include_dirs=[numpy.get_include()]
+                   )
 
+setup(ext_modules=cythonize([cy, c_wrap]))

Examples/week-08-extensions/cython/AGC_example/test_agc.py

@@ -28,45 +28,46 @@
 import agc_cython
 import agc_python
 import agc_c_cy
-import agc_subroutine
+
+# The Fortran / f2py version
+# import agc_subroutine
+
 
 def test_cython():
     # just make sure it runs.
     signal = np.arange(20, dtype=np.float32)
 
     result = agc_cython.agc(4, signal)
 
+
 def test_c_wrap():
     # just make sure it runs.
     signal = np.arange(20, dtype=np.float32)
 
     result = agc_c_cy.agc(4, signal)
 
-def test_subroutine():
-    # just make sure it runs.
-    signal = np.arange(20, dtype=np.float32)
-
-    result = agc_subroutine.agc(4, signal)
 
+# def test_subroutine():
+#     # the Fortran / f2py version
+#     # just make sure it runs.
+#     signal = np.arange(20, dtype=np.float32)
 
+#     result = agc_subroutine.agc(4, signal)
 
 
 def test_cy_py_same():
     signal = np.arange(20, dtype=np.float32)
-    
+
     cy_result = agc_cython.agc(4, signal)
     py_result = agc_python.agc(4, signal)
     c_cy_result = agc_c_cy.agc(4, signal)
-    sub_result = agc_subroutine.agc(4, signal)
-    
+#    sub_result = agc_subroutine.agc(4, signal)
+
     print "cy:", cy_result
     print "py:", py_result
     print "c_cy", c_cy_result
-    print "subroutine", sub_result
+#    print "subroutine", sub_result
 
     assert np.array_equal(cy_result, py_result)
     assert np.array_equal(cy_result, c_cy_result)
-    assert np.array_equal(cy_result, sub_result)
-
-
-
+#    assert np.array_equal(cy_result, sub_result)

Examples/week-08-extensions/cython/integrate/cy_integrate1.pyx

@@ -18,9 +18,7 @@ def integrate_f(a, b, N):
     """
 
     s = 0
-    dx = (b-a)/N
+    dx = (b - a) / N
     for i in range(N):
-        s += f(a+i*dx)
+        s += f(a + i * dx)
     return s * dx
-
-

Examples/week-08-extensions/cython/integrate/test_integrate.py

@@ -1,12 +1,18 @@
 #!/usr/bin/env python
 
+"""
+Simple test script to get gross timings
+
+Not unit tests
+"""
+
 from timer_context import Timer
 
 
 def run_test(msg):
     print msg
     with Timer():
-        result = integrate_f(0.0,10.0,1000000)
+        result = integrate_f(0.0, 10.0, 1000000)
     print "result:", result
     print
 
@@ -34,4 +40,3 @@ def run_test(msg):
 from cy_integrate7 import integrate_f
 run_test("Seventh Cython version:")
 
-

slides_sources/source/extensions.rst

@@ -780,6 +780,25 @@ At each step, we'll time and look at the output from::
 
   $cython -a cy_integrate1.pyx
 
+AGC Example
+-----------
+
+Another useful example is in:
+
+Examples/week-08-extensions/AGC_example
+
+This one impliments an Automatic Gain Control Signal processing filter.
+
+It turns out that you can use some advanced numpy tricks to get pretty
+good performancew with this filter, but you can't get full-on speed
+without some compiled code.
+
+This example uses:
+ * Pure Cython
+ * C called from Cython
+ * f2py and Fortran
+
+
 Auto-generated wrappers
 =======================
 
@@ -793,7 +812,7 @@ XDress
 
 [also Boost-Python -- not really a wrapper generator]
 
-f2Py -- for Fortran
+f2py -- for Fortran
 
 
 SWIG