Prev Next pwlfir

Piecewise Linear FIR Filter Design
Syntax y = pwlfir(fd,Hd,fs,Ntaps)
Syntax y = pwlfir(fd,Hd,fs,Ntaps,wintype)
Syntax y = pwlfir(fd,Hd,fs,Ntaps,wintype,winparam)
Include: include spt\pwlfir.oms
See Also linfir , hilbert , pbfir , rcfir , srrcfir

ARGUMENTS:
   INPUTS:
      fd       = INPUT, COLUMN VECTOR, any numerical type. Specifies the
                 N+1 frequency segment endpoints for the filter. Values
                 in [Hz], 0<=f<=fs/2. Coerced to DOUBLE before local
                 processing.
      Hd       = INPUT, COLUMN VECTOR, any numerical type. Specifies
                 the 2*N frequency magnitude values for the filter.
      fs       = INPUT, scalar, any numerical type. Represents sampling
                 rate in [Samples/second]. Coerced to DOUBLE.
      Ntaps    = INPUT, scalar, any numerical type. Filter length in
                 [taps]. Coerced to INTEGER.
      wintype  = INPUT, string. One of ["RECTANGULAR"|
                 "TRIANGULAR"|"HAMMING"|"HANNING"|"NUTTALL"|
                 "BLACKMAN"|"GAUSSIAN"|"KAISER"].
      winparam = INPUT, scalar, any numerical type. Parameter
                 used by GAUSSIAN or KAISER window functions.
                 Coerced to DOUBLE.
   RETURN: COLUMN VECTOR, DOUBLE, piece-wise linear FIR filter taps.

Description

Create a linear phase FIR filter with a frequency response that is constructed from a piece-wise linear specification of N contiguous segments. With this function any arbitrary frequency response may be obtained that may be approximated in a piece-wise linear fashion.

This function returns a linear phase FIR based on the specification of N piece-wise linear segments given in the frequency domain. The linear segments need not be continuous; discontinuities may be included. The piece-wise linear frequency domain segment definitions are given by the user in column vector 'fd' which contains N+1 frequency points and in a matching frequency response magnitude definition column vector 'Hd' that contains 2*N magnitude values.

'fd' lists, in increasing order, the N+1 frequency axis endpoints of the N frequency segments comprising the filter's desired frequency response. The frequency domain over which the filter response is specified is 0 to fs/2 (the Nyquist frequency). The first 'fd' entry must be 0.0(zero) and the last 'fd' entry must be equal to fs/2 while all intervening entries must be such that 0<fn<(fs/2) (strict inequality). For example: if you divide the domain 0 to fs/2 into 4 segments, there would be 5 endpoints, 0 and fs/s being the first and last endpoints. They would be listed in 'fd' as: fd = {0d0,f1,f2,f3,fs/2}. 'fd' may be any numerical type but is coerced to DOUBLE before local processing.

'Hd' lists, in increasing order, the 2*N frequency response magnitude values corresponding to the endpoints given in 'fd'. It is constructed as N successive pairs of magnitude specifications, one pair for each segment, where the first value is the frequency response magnitude as the lower endpoint is approached from above, and the second value is the magnitude as the upper endpoint is approached from below. Thus:

   Hd(1)     = mag at fd(1)=0.0  , approached from above, (Seg 1)
   Hd(2)     = mag at fd(2)=f1   , approached from below, (Seg 1)
   Hd(3)     = mag at fd(2)=f1   , approached from above, (Seg 2)
   Hd(4)     = mag at fd(3)=f2   , approached from below, (Seg 2)
   Hd(5)     = mag at fd(3)=f2   , approached from above, (Seg 3)
   ... etc.
   Hd(2*N-1) = mag at fd(3)=f2   , approached from above, (Seg N)
   Hd(2*N)   = mag at fd(N)=fs/2 , approached from below. (Seg N)

'Ntaps' is the length of the requested filter. Must be >=2.

'wintype' and 'winparam' specify a window function that may be applied to the filter. See function 'winddata()' for details.

Example

#######################################################################
# Arbitrary frequency response: 7 'random' segment definitions
#   Sampling Rate  = 10KHz
#   Segment #1 = fd=[0.0KHz,0.5KHz] Hd=[0.5, 0.5]
#   Segment #2 = fd=[0.5KHz,1.0KHz] Hd=[0.5, 0.7]
#   Segment #3 = fd=[1.0KHz,2.0KHz] Hd=[0.0, 0.0]
#   Segment #4 = fd=[2.0KHz,2.5KHz] Hd=[0.0, 0.3]
#   Segment #5 = fd=[2.5KHz,3.5KHz] Hd=[0.3, 1.0]
#   Segment #6 = fd=[3.5KHz,4.5KHz] Hd=[1.0, 0.6]
#   Segment #7 = fd=[4.5KHz,5.0KHz] Hd=[0.6, 0.6]
#   Note that fd 'covers' the whole 0 to fs/2 domain.
#######################################################################

#--- A 'clear' way to construct fd and Hd ---
#--- fd just lists the successive endpoints ---

fd = {0.0d3,0.5d3,1.0d3,2.0d3,2.5d3,3.5d3,4.5d3,5.0d3}

#--- Hd can be constructed from the concatenation of the ---
#--- successive pairs of segment endpoint magnitudes     ---

Hd1 = {0.5, 0.5}
Hd2 = {0.5, 0.7}
Hd3 = {0.0, 0.0}
Hd4 = {0.0, 0.3}
Hd5 = {0.3, 1.0}
Hd6 = {1.0, 0.6}
Hd7 = {0.6, 0.6}
Hd  = {Hd1,Hd2,Hd3,Hd4,Hd5,Hd6,Hd7}

#--- Make a template to show specifications in plot ---
ft1  = {0.0d3,0.5d3}
ft2  = {0.5d3,1.0d3}
ft2a = {1.0d3,1.0d3} # a discontinuity here
ft3  = {1.0d3,2.0d3}
ft4  = {2.0d3,2.5d3}
ft5  = {2.5d3,3.5d3}
ft6  = {3.5d3,4.5d3}
ft7  = {4.5d3,5.0d3}

Ht1  = {0.5d0,0.5d0}
Ht2  = {0.5d0,0.7d0}
Ht2a = {0.7d0,0.0d0} # a discontinuity here
Ht3  = {0.0d0,0.0d0}
Ht4  = {0.0d0,0.3d0}
Ht5  = {0.3d0,1.0d0}
Ht6  = {1.0d0,0.6d0}
Ht7  = {0.6d0,0.6d0}

ft = {ft1,ft2,ft2a,ft3,ft4,ft5,ft6,ft7}; # final spec templates
Ht = {Ht1,Ht2,Ht2a,Ht3,Ht4,Ht5,Ht6,Ht7};

#--- Set up and make the filter ---

fs       = 10d3; # Sampling Rate [Samples/Sec]
Ntaps    = 65;   # Filter Tap length
wintype  = "GAUSSIAN"; # Window Type
winparam = 4d0;        # Window Parameter
h        = pwlfir( fd, Hd, fs, Ntaps, wintype, winparam);

 Plot the filter along with it's specification template.



Reference

Chan et. al.