CMSIS DSP Library from CMSIS 2.0. See http://www.onarm.com/cmsis/ for full details

Dependents:   K22F_DSP_Matrix_least_square BNO055-ELEC3810 1BNO055 ECE4180Project--Slave2 ... more

Committer:
simon
Date:
Thu Mar 10 15:07:50 2011 +0000
Revision:
0:1014af42efd9

        

Who changed what in which revision?

UserRevisionLine numberNew contents of line
simon 0:1014af42efd9 1 /* ----------------------------------------------------------------------
simon 0:1014af42efd9 2 * Copyright (C) 2010 ARM Limited. All rights reserved.
simon 0:1014af42efd9 3 *
simon 0:1014af42efd9 4 * $Date: 29. November 2010
simon 0:1014af42efd9 5 * $Revision: V1.0.3
simon 0:1014af42efd9 6 *
simon 0:1014af42efd9 7 * Project: CMSIS DSP Library
simon 0:1014af42efd9 8 * Title: arm_biquad_cascade_df1_q15.c
simon 0:1014af42efd9 9 *
simon 0:1014af42efd9 10 * Description: Processing function for the
simon 0:1014af42efd9 11 * Q15 Biquad cascade DirectFormI(DF1) filter.
simon 0:1014af42efd9 12 *
simon 0:1014af42efd9 13 * Target Processor: Cortex-M4/Cortex-M3
simon 0:1014af42efd9 14 *
simon 0:1014af42efd9 15 * Version 1.0.3 2010/11/29
simon 0:1014af42efd9 16 * Re-organized the CMSIS folders and updated documentation.
simon 0:1014af42efd9 17 *
simon 0:1014af42efd9 18 * Version 1.0.2 2010/11/11
simon 0:1014af42efd9 19 * Documentation updated.
simon 0:1014af42efd9 20 *
simon 0:1014af42efd9 21 * Version 1.0.1 2010/10/05
simon 0:1014af42efd9 22 * Production release and review comments incorporated.
simon 0:1014af42efd9 23 *
simon 0:1014af42efd9 24 * Version 1.0.0 2010/09/20
simon 0:1014af42efd9 25 * Production release and review comments incorporated.
simon 0:1014af42efd9 26 *
simon 0:1014af42efd9 27 * Version 0.0.5 2010/04/26
simon 0:1014af42efd9 28 * incorporated review comments and updated with latest CMSIS layer
simon 0:1014af42efd9 29 *
simon 0:1014af42efd9 30 * Version 0.0.3 2010/03/10
simon 0:1014af42efd9 31 * Initial version
simon 0:1014af42efd9 32 * -------------------------------------------------------------------- */
simon 0:1014af42efd9 33
simon 0:1014af42efd9 34 #include "arm_math.h"
simon 0:1014af42efd9 35
simon 0:1014af42efd9 36 /**
simon 0:1014af42efd9 37 * @ingroup groupFilters
simon 0:1014af42efd9 38 */
simon 0:1014af42efd9 39
simon 0:1014af42efd9 40 /**
simon 0:1014af42efd9 41 * @addtogroup BiquadCascadeDF1
simon 0:1014af42efd9 42 * @{
simon 0:1014af42efd9 43 */
simon 0:1014af42efd9 44
simon 0:1014af42efd9 45 /**
simon 0:1014af42efd9 46 * @brief Processing function for the Q15 Biquad cascade filter.
simon 0:1014af42efd9 47 * @param[in] *S points to an instance of the Q15 Biquad cascade structure.
simon 0:1014af42efd9 48 * @param[in] *pSrc points to the block of input data.
simon 0:1014af42efd9 49 * @param[out] *pDst points to the location where the output result is written.
simon 0:1014af42efd9 50 * @param[in] blockSize number of samples to process per call.
simon 0:1014af42efd9 51 * @return none.
simon 0:1014af42efd9 52 *
simon 0:1014af42efd9 53 *
simon 0:1014af42efd9 54 * <b>Scaling and Overflow Behavior:</b>
simon 0:1014af42efd9 55 * \par
simon 0:1014af42efd9 56 * The function is implemented using a 64-bit internal accumulator.
simon 0:1014af42efd9 57 * Both coefficients and state variables are represented in 1.15 format and multiplications yield a 2.30 result.
simon 0:1014af42efd9 58 * The 2.30 intermediate results are accumulated in a 64-bit accumulator in 34.30 format.
simon 0:1014af42efd9 59 * There is no risk of internal overflow with this approach and the full precision of intermediate multiplications is preserved.
simon 0:1014af42efd9 60 * The accumulator is then shifted by <code>postShift</code> bits to truncate the result to 1.15 format by discarding the low 16 bits.
simon 0:1014af42efd9 61 * Finally, the result is saturated to 1.15 format.
simon 0:1014af42efd9 62 *
simon 0:1014af42efd9 63 * \par
simon 0:1014af42efd9 64 * Refer to the function <code>arm_biquad_cascade_df1_fast_q15()</code> for a faster but less precise implementation of this filter.
simon 0:1014af42efd9 65 */
simon 0:1014af42efd9 66
simon 0:1014af42efd9 67 void arm_biquad_cascade_df1_q15(
simon 0:1014af42efd9 68 const arm_biquad_casd_df1_inst_q15 * S,
simon 0:1014af42efd9 69 q15_t * pSrc,
simon 0:1014af42efd9 70 q15_t * pDst,
simon 0:1014af42efd9 71 uint32_t blockSize)
simon 0:1014af42efd9 72 {
simon 0:1014af42efd9 73 q15_t *pIn = pSrc; /* Source pointer */
simon 0:1014af42efd9 74 q15_t *pOut = pDst; /* Destination pointer */
simon 0:1014af42efd9 75 q31_t in; /* Temporary variable to hold input value */
simon 0:1014af42efd9 76 q31_t out; /* Temporary variable to hold output value */
simon 0:1014af42efd9 77 q15_t b0; /* Temporary variable to hold bo value */
simon 0:1014af42efd9 78 q31_t b1, a1; /* Filter coefficients */
simon 0:1014af42efd9 79 q31_t state_in, state_out; /* Filter state variables */
simon 0:1014af42efd9 80 q63_t acc; /* Accumulator */
simon 0:1014af42efd9 81 int32_t shift = (15 - (int32_t) S->postShift); /* Post shift */
simon 0:1014af42efd9 82 q15_t *pState = S->pState; /* State pointer */
simon 0:1014af42efd9 83 q15_t *pCoeffs = S->pCoeffs; /* Coefficient pointer */
simon 0:1014af42efd9 84 q31_t *pState_q31; /* 32-bit state pointer for SIMD implementation */
simon 0:1014af42efd9 85 uint32_t sample, stage = (uint32_t) S->numStages; /* Stage loop counter */
simon 0:1014af42efd9 86
simon 0:1014af42efd9 87
simon 0:1014af42efd9 88
simon 0:1014af42efd9 89 do
simon 0:1014af42efd9 90 {
simon 0:1014af42efd9 91 /* Initialize state pointer of type q31 */
simon 0:1014af42efd9 92 pState_q31 = (q31_t *) (pState);
simon 0:1014af42efd9 93
simon 0:1014af42efd9 94 /* Read the b0 and 0 coefficients using SIMD */
simon 0:1014af42efd9 95 b0 = *__SIMD32(pCoeffs)++;
simon 0:1014af42efd9 96
simon 0:1014af42efd9 97 /* Read the b1 and b2 coefficients using SIMD */
simon 0:1014af42efd9 98 b1 = *__SIMD32(pCoeffs)++;
simon 0:1014af42efd9 99
simon 0:1014af42efd9 100 /* Read the a1 and a2 coefficients using SIMD */
simon 0:1014af42efd9 101 a1 = *__SIMD32(pCoeffs)++;
simon 0:1014af42efd9 102
simon 0:1014af42efd9 103 /* Read the input state values from the state buffer: x[n-1], x[n-2] */
simon 0:1014af42efd9 104 state_in = (q31_t) (*pState_q31++);
simon 0:1014af42efd9 105
simon 0:1014af42efd9 106 /* Read the output state values from the state buffer: y[n-1], y[n-2] */
simon 0:1014af42efd9 107 state_out = (q31_t) (*pState_q31);
simon 0:1014af42efd9 108
simon 0:1014af42efd9 109 /* Apply loop unrolling and compute 2 output values simultaneously. */
simon 0:1014af42efd9 110 /* The variable acc hold output values that are being computed:
simon 0:1014af42efd9 111 *
simon 0:1014af42efd9 112 * acc = b0 * x[n] + b1 * x[n-1] + b2 * x[n-2] + a1 * y[n-1] + a2 * y[n-2]
simon 0:1014af42efd9 113 * acc = b0 * x[n] + b1 * x[n-1] + b2 * x[n-2] + a1 * y[n-1] + a2 * y[n-2]
simon 0:1014af42efd9 114 */
simon 0:1014af42efd9 115 sample = blockSize >> 1u;
simon 0:1014af42efd9 116
simon 0:1014af42efd9 117 /* First part of the processing with loop unrolling. Compute 2 outputs at a time.
simon 0:1014af42efd9 118 ** a second loop below computes the remaining 1 sample. */
simon 0:1014af42efd9 119 while(sample > 0u)
simon 0:1014af42efd9 120 {
simon 0:1014af42efd9 121
simon 0:1014af42efd9 122 /* Read the input */
simon 0:1014af42efd9 123 in = *__SIMD32(pIn)++;
simon 0:1014af42efd9 124
simon 0:1014af42efd9 125 /* out = b0 * x[n] + 0 * 0 */
simon 0:1014af42efd9 126 out = (q31_t) b0 * ((q15_t) in);
simon 0:1014af42efd9 127 /* acc += b1 * x[n-1] + b2 * x[n-2] + out */
simon 0:1014af42efd9 128 acc = __SMLALD(b1, state_in, out);
simon 0:1014af42efd9 129 /* acc += a1 * y[n-1] + a2 * y[n-2] */
simon 0:1014af42efd9 130 acc = __SMLALD(a1, state_out, acc);
simon 0:1014af42efd9 131
simon 0:1014af42efd9 132 /* The result is converted from 3.29 to 1.31 if postShift = 1, and then saturation is applied */
simon 0:1014af42efd9 133 out = __SSAT((acc >> shift), 16);
simon 0:1014af42efd9 134
simon 0:1014af42efd9 135 /* Every time after the output is computed state should be updated. */
simon 0:1014af42efd9 136 /* The states should be updated as: */
simon 0:1014af42efd9 137 /* Xn2 = Xn1 */
simon 0:1014af42efd9 138 /* Xn1 = Xn */
simon 0:1014af42efd9 139 /* Yn2 = Yn1 */
simon 0:1014af42efd9 140 /* Yn1 = acc */
simon 0:1014af42efd9 141 /* x[n-N], x[n-N-1] are packed together to make state_in of type q31 */
simon 0:1014af42efd9 142 /* y[n-N], y[n-N-1] are packed together to make state_out of type q31 */
simon 0:1014af42efd9 143 state_in = __PKHBT(in, state_in, 16);
simon 0:1014af42efd9 144 state_out = __PKHBT(out, state_out, 16);
simon 0:1014af42efd9 145
simon 0:1014af42efd9 146 /* out = b0 * x[n] + 0 * 0 */
simon 0:1014af42efd9 147 out = (q31_t) b0 * ((q15_t) (in >> 16));
simon 0:1014af42efd9 148 /* acc += b1 * x[n-1] + b2 * x[n-2] + out */
simon 0:1014af42efd9 149 acc = __SMLALD(b1, state_in, out);
simon 0:1014af42efd9 150 /* acc += a1 * y[n-1] + a2 * y[n-2] */
simon 0:1014af42efd9 151 acc = __SMLALD(a1, state_out, acc);
simon 0:1014af42efd9 152
simon 0:1014af42efd9 153 /* The result is converted from 3.29 to 1.31 if postShift = 1, and then saturation is applied */
simon 0:1014af42efd9 154 out = __SSAT((acc >> shift), 16);
simon 0:1014af42efd9 155
simon 0:1014af42efd9 156 /* Store the output in the destination buffer. */
simon 0:1014af42efd9 157 *__SIMD32(pOut)++ = __PKHBT(state_out, out, 16);
simon 0:1014af42efd9 158
simon 0:1014af42efd9 159 /* Every time after the output is computed state should be updated. */
simon 0:1014af42efd9 160 /* The states should be updated as: */
simon 0:1014af42efd9 161 /* Xn2 = Xn1 */
simon 0:1014af42efd9 162 /* Xn1 = Xn */
simon 0:1014af42efd9 163 /* Yn2 = Yn1 */
simon 0:1014af42efd9 164 /* Yn1 = acc */
simon 0:1014af42efd9 165 /* x[n-N], x[n-N-1] are packed together to make state_in of type q31 */
simon 0:1014af42efd9 166 /* y[n-N], y[n-N-1] are packed together to make state_out of type q31 */
simon 0:1014af42efd9 167 state_in = __PKHBT(in >> 16, state_in, 16);
simon 0:1014af42efd9 168 state_out = __PKHBT(out, state_out, 16);
simon 0:1014af42efd9 169
simon 0:1014af42efd9 170 /* Decrement the loop counter */
simon 0:1014af42efd9 171 sample--;
simon 0:1014af42efd9 172
simon 0:1014af42efd9 173 }
simon 0:1014af42efd9 174
simon 0:1014af42efd9 175 /* If the blockSize is not a multiple of 2, compute any remaining output samples here.
simon 0:1014af42efd9 176 ** No loop unrolling is used. */
simon 0:1014af42efd9 177
simon 0:1014af42efd9 178 if((blockSize & 0x1u) != 0u)
simon 0:1014af42efd9 179 {
simon 0:1014af42efd9 180 /* Read the input */
simon 0:1014af42efd9 181 in = *pIn++;
simon 0:1014af42efd9 182
simon 0:1014af42efd9 183 /* out = b0 * x[n] + 0 * 0 */
simon 0:1014af42efd9 184 out = (q31_t) in *b0;
simon 0:1014af42efd9 185 /* acc = b1 * x[n-1] + b2 * x[n-2] + out */
simon 0:1014af42efd9 186 acc = __SMLALD(b1, state_in, out);
simon 0:1014af42efd9 187 /* acc += a1 * y[n-1] + a2 * y[n-2] */
simon 0:1014af42efd9 188 acc = __SMLALD(a1, state_out, acc);
simon 0:1014af42efd9 189
simon 0:1014af42efd9 190 /* The result is converted from 3.29 to 1.31 if postShift = 1, and then saturation is applied */
simon 0:1014af42efd9 191 out = __SSAT((acc >> shift), 16);
simon 0:1014af42efd9 192
simon 0:1014af42efd9 193 /* Store the output in the destination buffer. */
simon 0:1014af42efd9 194 *pOut++ = (q15_t) out;
simon 0:1014af42efd9 195
simon 0:1014af42efd9 196 /* Every time after the output is computed state should be updated. */
simon 0:1014af42efd9 197 /* The states should be updated as: */
simon 0:1014af42efd9 198 /* Xn2 = Xn1 */
simon 0:1014af42efd9 199 /* Xn1 = Xn */
simon 0:1014af42efd9 200 /* Yn2 = Yn1 */
simon 0:1014af42efd9 201 /* Yn1 = acc */
simon 0:1014af42efd9 202 /* x[n-N], x[n-N-1] are packed together to make state_in of type q31 */
simon 0:1014af42efd9 203 /* y[n-N], y[n-N-1] are packed together to make state_out of type q31 */
simon 0:1014af42efd9 204 state_in = __PKHBT(in, state_in, 16);
simon 0:1014af42efd9 205 state_out = __PKHBT(out, state_out, 16);
simon 0:1014af42efd9 206
simon 0:1014af42efd9 207 }
simon 0:1014af42efd9 208
simon 0:1014af42efd9 209 /* The first stage goes from the input wire to the output wire. */
simon 0:1014af42efd9 210 /* Subsequent numStages occur in-place in the output wire */
simon 0:1014af42efd9 211 pIn = pDst;
simon 0:1014af42efd9 212
simon 0:1014af42efd9 213 /* Reset the output pointer */
simon 0:1014af42efd9 214 pOut = pDst;
simon 0:1014af42efd9 215
simon 0:1014af42efd9 216 /* Store the updated state variables back into the state array */
simon 0:1014af42efd9 217 *__SIMD32(pState)++ = __PKHBT(state_in, (state_in >> 16), 16);
simon 0:1014af42efd9 218 *__SIMD32(pState)++ = __PKHBT(state_out, (state_out >> 16), 16);
simon 0:1014af42efd9 219
simon 0:1014af42efd9 220 /* Decrement the loop counter */
simon 0:1014af42efd9 221 stage--;
simon 0:1014af42efd9 222
simon 0:1014af42efd9 223 } while(stage > 0u);
simon 0:1014af42efd9 224 }
simon 0:1014af42efd9 225
simon 0:1014af42efd9 226
simon 0:1014af42efd9 227 /**
simon 0:1014af42efd9 228 * @} end of BiquadCascadeDF1 group
simon 0:1014af42efd9 229 */