tomwalters@0
|
1 /*
|
tomwalters@0
|
2 Copyright (c) Applied Psychology Unit, Medical Research Council. 1988, 1989
|
tomwalters@0
|
3 ===========================================================================
|
tomwalters@0
|
4
|
tomwalters@0
|
5 Permission to use, copy, modify, and distribute this software without fee
|
tomwalters@0
|
6 is hereby granted for research purposes, provided that this copyright
|
tomwalters@0
|
7 notice appears in all copies and in all supporting documentation, and that
|
tomwalters@0
|
8 the software is not redistributed for any fee (except for a nominal shipping
|
tomwalters@0
|
9 charge). Anyone wanting to incorporate all or part of this software in a
|
tomwalters@0
|
10 commercial product must obtain a license from the Medical Research Council.
|
tomwalters@0
|
11
|
tomwalters@0
|
12 The MRC makes no representations about the suitability of this
|
tomwalters@0
|
13 software for any purpose. It is provided "as is" without express or implied
|
tomwalters@0
|
14 warranty.
|
tomwalters@0
|
15
|
tomwalters@0
|
16 THE MRC DISCLAIMS ALL WARRANTIES WITH REGARD TO THIS SOFTWARE, INCLUDING
|
tomwalters@0
|
17 ALL IMPLIED WARRANTIES OF MERCHANTABILITY AND FITNESS, IN NO EVENT SHALL THE
|
tomwalters@0
|
18 A.P.U. BE LIABLE FOR ANY SPECIAL, INDIRECT OR CONSEQUENTIAL DAMAGES OR ANY
|
tomwalters@0
|
19 DAMAGES WHATSOEVER RESULTING FROM LOSS OF USE, DATA OR PROFITS, WHETHER IN
|
tomwalters@0
|
20 AN ACTION OF CONTRACT, NEGLIGENCE OR OTHER TORTIOUS ACTION, ARISING OUT OF
|
tomwalters@0
|
21 OR IN CONNECTION WITH THE USE OR PERFORMANCE OF THIS SOFTWARE.
|
tomwalters@0
|
22 */
|
tomwalters@0
|
23
|
tomwalters@0
|
24 /*
|
tomwalters@0
|
25 Acknowledgment:
|
tomwalters@0
|
26 ==============
|
tomwalters@0
|
27
|
tomwalters@0
|
28 The source code provided in this file was originally developed by
|
tomwalters@0
|
29 Christian Giguere as part of a Ph.D degree at the Department of
|
tomwalters@0
|
30 Engineering of the University of Cambridge from April 1990 to
|
tomwalters@0
|
31 November 1993. The code was subsequently adapted under a grant
|
tomwalters@0
|
32 from the Hearing Research Trust for full compatibility with
|
tomwalters@0
|
33 AIM Release 6.15.
|
tomwalters@0
|
34
|
tomwalters@0
|
35 Christian Giguere 25/03/94
|
tomwalters@0
|
36
|
tomwalters@0
|
37 */
|
tomwalters@0
|
38
|
tomwalters@0
|
39 /*
|
tomwalters@0
|
40 ===========================================================
|
tomwalters@0
|
41 fir.c
|
tomwalters@0
|
42 ===========================================================
|
tomwalters@0
|
43
|
tomwalters@0
|
44 Finite Impulse Response (FIR) filtering module.
|
tomwalters@0
|
45
|
tomwalters@0
|
46 Author : Christian Giguere
|
tomwalters@0
|
47 First written : 12th November, 1990
|
tomwalters@0
|
48 Last edited : 22nd September, 1991
|
tomwalters@0
|
49
|
tomwalters@0
|
50 Reference:
|
tomwalters@0
|
51 (1) A.V.Oppenheim and R.W.Schafer (1975). Digital Signal
|
tomwalters@0
|
52 Processing (Prentice-Hall), Sections 4.5.5 and 5.5.
|
tomwalters@0
|
53 ===========================================================
|
tomwalters@0
|
54 */
|
tomwalters@0
|
55
|
tomwalters@0
|
56
|
tomwalters@0
|
57 /***** includes *****/
|
tomwalters@0
|
58
|
tomwalters@0
|
59 #include <math.h>
|
tomwalters@0
|
60 #include <stdio.h>
|
tomwalters@0
|
61 #include "stitch.h"
|
tomwalters@0
|
62 #include "calc.h"
|
tomwalters@0
|
63 #include "calc_tl.h"
|
tomwalters@0
|
64 #include "fir.h"
|
tomwalters@0
|
65
|
tomwalters@0
|
66
|
tomwalters@0
|
67 /***** defines *****/
|
tomwalters@0
|
68
|
tomwalters@0
|
69 #define NUMBER_OF_COEFFS ( 51 ) /* odd order */
|
tomwalters@0
|
70 #define NUMBER_OF_STATES ( NUMBER_OF_COEFFS )
|
tomwalters@0
|
71
|
tomwalters@0
|
72
|
tomwalters@0
|
73 /***** functions *****/
|
tomwalters@0
|
74
|
tomwalters@0
|
75 /*******************************************************************************
|
tomwalters@0
|
76 * name: function:
|
tomwalters@0
|
77 *
|
tomwalters@0
|
78 * NewLpFIRfilter() Design of LowPass FIR filter using windowing technique
|
tomwalters@0
|
79 * (Hamming window). It generates the multiplier coefficients
|
tomwalters@0
|
80 * and initializes the states of the filter. It returns a
|
tomwalters@0
|
81 * pointer to a structure containing all relevant information
|
tomwalters@0
|
82 * for the realisation of the filter.
|
tomwalters@0
|
83 *
|
tomwalters@0
|
84 * DoFIRfilter() Direct form realisation for FIR filter of odd order linear
|
tomwalters@0
|
85 * phase. It is called by function ``upsample_callback()''
|
tomwalters@0
|
86 * in file ``upsample.c'' once for a specified number of input
|
tomwalters@0
|
87 * points. It computes the filter output for each input point
|
tomwalters@0
|
88 * and keeps track of the filter states.
|
tomwalters@0
|
89 *
|
tomwalters@0
|
90 * CloseFIRfilter() Deletes all private data structures and arrays of the FIR
|
tomwalters@0
|
91 * filter upon comletion of filtering. It is called by function
|
tomwalters@0
|
92 * ``upsample_callback()'' in file ``upsample.c''.
|
tomwalters@0
|
93 ********************************************************************************/
|
tomwalters@0
|
94
|
tomwalters@0
|
95
|
tomwalters@0
|
96 /************************* NewLpFIRfilter() *******************************/
|
tomwalters@0
|
97
|
tomwalters@0
|
98 FIRfilterState *NewLpFIRfilter( wc, upsample_factor, delay )
|
tomwalters@0
|
99 double wc ;
|
tomwalters@0
|
100 int upsample_factor, *delay ;
|
tomwalters@0
|
101 {
|
tomwalters@0
|
102 DeclareNew( FIRfilterState *, filter_state ) ;
|
tomwalters@0
|
103 CoeffType *hn ;
|
tomwalters@0
|
104 double filterGain, hamming() ;
|
tomwalters@0
|
105 int n ;
|
tomwalters@0
|
106
|
tomwalters@0
|
107 /***** generate and store FIR filter coefficients (linear phase, odd order, window desing technique) *****/
|
tomwalters@0
|
108 filter_state->number_of_coeffs = ( NUMBER_OF_COEFFS + 1 ) / 2 ;
|
tomwalters@0
|
109 filter_state->coeffs = ( char * ) NewArray( CoeffType, filter_state->number_of_coeffs, "fir.c for coefficients\n" ) ;
|
tomwalters@0
|
110
|
tomwalters@0
|
111 *delay = ( NUMBER_OF_COEFFS - 1 ) / 2 ;
|
tomwalters@0
|
112 hn = ( CoeffType * ) filter_state->coeffs ;
|
tomwalters@0
|
113 filterGain = ( double ) upsample_factor ;
|
tomwalters@0
|
114
|
tomwalters@0
|
115 for( n = 0 ; n < filter_state->number_of_coeffs - 1 ; n++ )
|
tomwalters@0
|
116 *hn++ = sin( wc * ( n - *delay ) ) / ( Pi * ( n - *delay ) ) * hamming( n, NUMBER_OF_COEFFS ) * filterGain ;
|
tomwalters@0
|
117 *hn = wc / Pi * filterGain ;
|
tomwalters@0
|
118
|
tomwalters@0
|
119
|
tomwalters@0
|
120 /***** initialise filter states *****/
|
tomwalters@0
|
121 filter_state->number_of_states = NUMBER_OF_STATES ;
|
tomwalters@0
|
122 filter_state->states = ( char * ) NewZeroedArray( StateType, NUMBER_OF_STATES, "fir.c for filter states\n" ) ;
|
tomwalters@0
|
123
|
tomwalters@0
|
124
|
tomwalters@0
|
125 /***** return *****/
|
tomwalters@0
|
126 return( filter_state ) ;
|
tomwalters@0
|
127 }
|
tomwalters@0
|
128
|
tomwalters@0
|
129
|
tomwalters@0
|
130 /****************************** DoFIRfilter() *******************************/
|
tomwalters@0
|
131
|
tomwalters@0
|
132 void DoFIRfilter( filter_state, inoutput, points )
|
tomwalters@0
|
133 FIRfilterState *filter_state ;
|
tomwalters@0
|
134 DataType *inoutput ;
|
tomwalters@0
|
135 int points ;
|
tomwalters@0
|
136 {
|
tomwalters@0
|
137 register StateType *st = ( StateType * ) filter_state->states ;
|
tomwalters@0
|
138 register CoeffType *cf = ( CoeffType * ) filter_state->coeffs ;
|
tomwalters@0
|
139 register DataType *inout_ptr = inoutput ;
|
tomwalters@0
|
140 register double acc ;
|
tomwalters@0
|
141 register int n, Ns = filter_state->number_of_states ;
|
tomwalters@0
|
142 register int Nc = filter_state->number_of_coeffs ;
|
tomwalters@0
|
143
|
tomwalters@0
|
144 /***** for all points *****/
|
tomwalters@0
|
145 for( ; points > 0 ; points-- ) {
|
tomwalters@0
|
146 acc = 0. ;
|
tomwalters@0
|
147
|
tomwalters@0
|
148 /*** shift the states ***/
|
tomwalters@0
|
149 for( n = Ns - 1 ; n > 0 ; n-- )
|
tomwalters@0
|
150 st[n] = st[n-1] ;
|
tomwalters@0
|
151 st[0] = *inout_ptr ;
|
tomwalters@0
|
152
|
tomwalters@0
|
153 /*** filter current point ***/
|
tomwalters@0
|
154 for( n = 0 ; n < Nc - 1 ; n++ )
|
tomwalters@0
|
155 acc += cf[n] * ( st[n] + st[Ns-1-n] ) ;
|
tomwalters@0
|
156 acc += cf[n] * st[n] ;
|
tomwalters@0
|
157
|
tomwalters@0
|
158 /*** store output point ***/
|
tomwalters@0
|
159 *inout_ptr++ = acc ;
|
tomwalters@0
|
160 }
|
tomwalters@0
|
161
|
tomwalters@0
|
162 /***** return *****/
|
tomwalters@0
|
163 return ;
|
tomwalters@0
|
164
|
tomwalters@0
|
165 }
|
tomwalters@0
|
166
|
tomwalters@0
|
167 /************************* CloseFIRfilter() *******************************/
|
tomwalters@0
|
168
|
tomwalters@0
|
169 void CloseFIRfilter( filter_state )
|
tomwalters@0
|
170 FIRfilterState *filter_state ;
|
tomwalters@0
|
171 {
|
tomwalters@0
|
172 Delete( filter_state->states ) ;
|
tomwalters@0
|
173 Delete( filter_state->coeffs ) ;
|
tomwalters@0
|
174 Delete( filter_state ) ;
|
tomwalters@0
|
175
|
tomwalters@0
|
176 return ;
|
tomwalters@0
|
177
|
tomwalters@0
|
178 }
|
tomwalters@0
|
179
|
tomwalters@0
|
180
|
tomwalters@0
|
181 /*********************** low level functions ***********************/
|
tomwalters@0
|
182
|
tomwalters@0
|
183 double hamming( n, N )
|
tomwalters@0
|
184 int n, N ;
|
tomwalters@0
|
185 {
|
tomwalters@0
|
186 return ( 0.54 - 0.46 * cos( TwoPi * n / ( N - 1. ) ) ) ;
|
tomwalters@0
|
187 }
|
tomwalters@0
|
188
|