Designing an FIR notch filter with python
Question: I am writing some code with Python using the scipy.signal library to filter
electromagnetic data that is mixed with various undesirable signatures that I
want to filter out. For example, I have power line harmonics at various
frequencies (i.e. 60, 120 Hz, etc....) with a width of only a few Hz that I
would like to remove from the data using a notch filter. Is there already an
existing function in python where I can merely inform the code how many data
points i wish to use for the filter, the center-line frequency that I wish to
remove and the width of the transition band or do I need to design a filter
from scratch? If it is the latter I would greatly appreciate an example of
notch filter design in Python to include window implementation to minimize
aliasing.


Answer: There are a few options for the solution on the scipy.signal website, but they
introduce allot of ringing, which will translate to artifacts in the convolved
signal. After trying many things I found the following function worked the
best for implementing an FIR notch filter.

    
    
    # Required input defintions are as follows;
    # time:   Time between samples
    # band:   The bandwidth around the centerline freqency that you wish to filter
    # freq:   The centerline frequency to be filtered
    # ripple: The maximum passband ripple that is allowed in db
    # order:  The filter order.  For FIR notch filters this is best set to 2 or 3,
    #         IIR filters are best suited for high values of order.  This algorithm
    #         is hard coded to FIR filters
    # filter_type: 'butter', 'bessel', 'cheby1', 'cheby2', 'ellip'
    # data:         the data to be filtered
    def Implement_Notch_Filter(time, band, freq, ripple, order, filter_type, data):
        from scipy.signal import iirfilter
        fs   = 1/time
        nyq  = fs/2.0
        low  = freq - band/2.0
        high = freq + band/2.0
        low  = low/nyq
        high = high/nyq
        b, a = iirfilter(order, [low, high], rp=ripple, btype='bandstop',
                         analog=False, ftype=filter_type)
        filtered_data = lfilter(b, a, data)
        return filtered_data
    

