comparison projects/heavy/envelopeTrigger/SignalBiquad.h @ 162:c3e8226a5651 heavy-updated

- added additional flags to C rules (-DNDEBUG, -mfpu=neon) - sample-accurate envelope triggering pd/heavy example
author chnrx <chris.heinrichs@gmail.com>
date Thu, 12 Nov 2015 14:59:46 +0000
parents
children
comparison
equal deleted inserted replaced
161:07735c9d95c8 162:c3e8226a5651
1 /**
2 * Copyright (c) 2014, 2015, Enzien Audio Ltd.
3 *
4 * Permission to use, copy, modify, and/or distribute this software for any
5 * purpose with or without fee is hereby granted, provided that the above
6 * copyright notice and this permission notice appear in all copies.
7 *
8 * THE SOFTWARE IS PROVIDED "AS IS" AND THE AUTHOR DISCLAIMS ALL WARRANTIES WITH
9 * REGARD TO THIS SOFTWARE INCLUDING ALL IMPLIED WARRANTIES OF MERCHANTABILITY
10 * AND FITNESS. IN NO EVENT SHALL THE AUTHOR BE LIABLE FOR ANY SPECIAL, DIRECT,
11 * INDIRECT, OR CONSEQUENTIAL DAMAGES OR ANY DAMAGES WHATSOEVER RESULTING FROM
12 * LOSS OF USE, DATA OR PROFITS, WHETHER IN AN ACTION OF CONTRACT, NEGLIGENCE OR
13 * OTHER TORTIOUS ACTION, ARISING OUT OF OR IN CONNECTION WITH THE USE OR
14 * PERFORMANCE OF THIS SOFTWARE.
15 */
16
17 #ifndef _HEAVY_SIGNAL_BIQUAD_H_
18 #define _HEAVY_SIGNAL_BIQUAD_H_
19
20 #include "HvBase.h"
21
22 // http://en.wikipedia.org/wiki/Digital_biquad_filter
23 typedef struct SignalBiquad {
24 #if HV_SIMD_AVX
25 __m256 xm1;
26 __m256 xm2;
27 #elif HV_SIMD_SSE
28 __m128 xm1;
29 __m128 xm2;
30 #elif HV_SIMD_NEON
31 float32x4_t xm1;
32 float32x4_t xm2;
33 #else // HV_SIMD_NONE
34 float x1;
35 float x2;
36 #endif
37 float y1;
38 float y2;
39 } SignalBiquad;
40
41 hv_size_t sBiquad_init(SignalBiquad *o);
42
43 void __hv_biquad_f(SignalBiquad *o,
44 hv_bInf_t bIn, hv_bInf_t bX0, hv_bInf_t bX1, hv_bInf_t bX2, hv_bInf_t bY1, hv_bInf_t bY2,
45 hv_bOutf_t bOut);
46
47 typedef struct SignalBiquad_k {
48 #if HV_SIMD_AVX || HV_SIMD_SSE
49 // preprocessed filter coefficients
50 __m128 coeff_xp3;
51 __m128 coeff_xp2;
52 __m128 coeff_xp1;
53 __m128 coeff_x0;
54 __m128 coeff_xm1;
55 __m128 coeff_xm2;
56 __m128 coeff_ym1;
57 __m128 coeff_ym2;
58
59 // filter state
60 __m128 xm1;
61 __m128 xm2;
62 __m128 ym1;
63 __m128 ym2;
64 #elif HV_SIMD_NEON
65 float32x4_t coeff_xp3;
66 float32x4_t coeff_xp2;
67 float32x4_t coeff_xp1;
68 float32x4_t coeff_x0;
69 float32x4_t coeff_xm1;
70 float32x4_t coeff_xm2;
71 float32x4_t coeff_ym1;
72 float32x4_t coeff_ym2;
73 float32x4_t xm1;
74 float32x4_t xm2;
75 float32x4_t ym1;
76 float32x4_t ym2;
77 #else // HV_SIMD_NONE
78 float xm1;
79 float xm2;
80 float ym1;
81 float ym2;
82 #endif
83 // original filter coefficients
84 float b0; // x[0]
85 float b1; // x[-1]
86 float b2; // x[-2]
87 float a1; // y[-1]
88 float a2; // y[-2]
89 } SignalBiquad_k;
90
91 hv_size_t sBiquad_k_init(SignalBiquad_k *o, float x0, float x1, float x2, float y1, float y2);
92
93 void sBiquad_k_onMessage(SignalBiquad_k *o, int letIn, const HvMessage *const m);
94
95 static inline void __hv_biquad_k_f(SignalBiquad_k *o, hv_bInf_t bIn, hv_bOutf_t bOut) {
96 #if HV_SIMD_AVX
97 const __m128 c_xp3 = o->coeff_xp3;
98 const __m128 c_xp2 = o->coeff_xp2;
99 const __m128 c_xp1 = o->coeff_xp1;
100 const __m128 c_x0 = o->coeff_x0;
101 const __m128 c_xm1 = o->coeff_xm1;
102 const __m128 c_xm2 = o->coeff_xm2;
103 const __m128 c_ym1 = o->coeff_ym1;
104 const __m128 c_ym2 = o->coeff_ym2;
105
106 // lower half
107 __m128 x3 = _mm_set1_ps(bIn[3]);
108 __m128 x2 = _mm_set1_ps(bIn[2]);
109 __m128 x1 = _mm_set1_ps(bIn[1]);
110 __m128 x0 = _mm_set1_ps(bIn[0]);
111 __m128 xm1 = o->xm1;
112 __m128 xm2 = o->xm2;
113 __m128 ym1 = o->ym1;
114 __m128 ym2 = o->ym2;
115
116 __m128 a = _mm_mul_ps(c_xp3, x3);
117 __m128 b = _mm_mul_ps(c_xp2, x2);
118 __m128 c = _mm_mul_ps(c_xp1, x1);
119 __m128 d = _mm_mul_ps(c_x0, x0);
120 __m128 e = _mm_mul_ps(c_xm1, xm1);
121 __m128 f = _mm_mul_ps(c_xm2, xm2);
122 __m128 g = _mm_mul_ps(c_ym1, ym1);
123 __m128 h = _mm_mul_ps(c_ym2, ym2);
124
125 __m128 i = _mm_add_ps(a, b);
126 __m128 j = _mm_add_ps(c, d);
127 __m128 k = _mm_add_ps(e, f);
128 __m128 l = _mm_add_ps(g, h);
129 __m128 m = _mm_add_ps(i, j);
130 __m128 n = _mm_add_ps(k, l);
131
132 __m128 lo_y = _mm_add_ps(m, n); // lower part of output buffer
133
134 // upper half
135 xm1 = x3;
136 xm2 = x2;
137 x3 = _mm_set1_ps(bIn[7]);
138 x2 = _mm_set1_ps(bIn[6]);
139 x1 = _mm_set1_ps(bIn[5]);
140 x0 = _mm_set1_ps(bIn[4]);
141 ym1 = _mm_set1_ps(lo_y[3]);
142 ym2 = _mm_set1_ps(lo_y[2]);
143
144 a = _mm_mul_ps(c_xp3, x3);
145 b = _mm_mul_ps(c_xp2, x2);
146 c = _mm_mul_ps(c_xp1, x1);
147 d = _mm_mul_ps(c_x0, x0);
148 e = _mm_mul_ps(c_xm1, xm1);
149 f = _mm_mul_ps(c_xm2, xm2);
150 g = _mm_mul_ps(c_ym1, ym1);
151 h = _mm_mul_ps(c_ym2, ym2);
152
153 i = _mm_add_ps(a, b);
154 j = _mm_add_ps(c, d);
155 k = _mm_add_ps(e, f);
156 l = _mm_add_ps(g, h);
157 m = _mm_add_ps(i, j);
158 n = _mm_add_ps(k, l);
159
160 __m128 up_y = _mm_add_ps(m, n); // upper part of output buffer
161
162 o->xm1 = x3;
163 o->xm2 = x2;
164 o->ym1 = _mm_set1_ps(up_y[3]);
165 o->ym2 = _mm_set1_ps(up_y[2]);
166
167 *bOut = _mm256_insertf128_ps(_mm256_castps128_ps256(lo_y), up_y, 1);
168 #elif HV_SIMD_SSE
169 __m128 x3 = _mm_set1_ps(bIn[3]);
170 __m128 x2 = _mm_set1_ps(bIn[2]);
171 __m128 x1 = _mm_set1_ps(bIn[1]);
172 __m128 x0 = _mm_set1_ps(bIn[0]);
173
174 __m128 a = _mm_mul_ps(o->coeff_xp3, x3);
175 __m128 b = _mm_mul_ps(o->coeff_xp2, x2);
176 __m128 c = _mm_mul_ps(o->coeff_xp1, x1);
177 __m128 d = _mm_mul_ps(o->coeff_x0, x0);
178 __m128 e = _mm_mul_ps(o->coeff_xm1, o->xm1);
179 __m128 f = _mm_mul_ps(o->coeff_xm2, o->xm2);
180 __m128 g = _mm_mul_ps(o->coeff_ym1, o->ym1);
181 __m128 h = _mm_mul_ps(o->coeff_ym2, o->ym2);
182 __m128 i = _mm_add_ps(a, b);
183 __m128 j = _mm_add_ps(c, d);
184 __m128 k = _mm_add_ps(e, f);
185 __m128 l = _mm_add_ps(g, h);
186 __m128 m = _mm_add_ps(i, j);
187 __m128 n = _mm_add_ps(k, l);
188
189 __m128 y = _mm_add_ps(m, n);
190
191 o->xm1 = x3;
192 o->xm2 = x2;
193 o->ym1 = _mm_set1_ps(y[3]);
194 o->ym2 = _mm_set1_ps(y[2]);
195
196 *bOut = y;
197 #elif HV_SIMD_NEON
198 float32x4_t x3 = vdupq_n_f32(bIn[3]);
199 float32x4_t x2 = vdupq_n_f32(bIn[2]);
200 float32x4_t x1 = vdupq_n_f32(bIn[1]);
201 float32x4_t x0 = vdupq_n_f32(bIn[0]);
202
203 float32x4_t a = vmulq_f32(o->coeff_xp3, x3);
204 float32x4_t b = vmulq_f32(o->coeff_xp2, x2);
205 float32x4_t c = vmulq_f32(o->coeff_xp1, x1);
206 float32x4_t d = vmulq_f32(o->coeff_x0, x0);
207 float32x4_t e = vmulq_f32(o->coeff_xm1, o->xm1);
208 float32x4_t f = vmulq_f32(o->coeff_xm2, o->xm2);
209 float32x4_t g = vmulq_f32(o->coeff_ym1, o->ym1);
210 float32x4_t h = vmulq_f32(o->coeff_ym2, o->ym2);
211 float32x4_t i = vaddq_f32(a, b);
212 float32x4_t j = vaddq_f32(c, d);
213 float32x4_t k = vaddq_f32(e, f);
214 float32x4_t l = vaddq_f32(g, h);
215 float32x4_t m = vaddq_f32(i, j);
216 float32x4_t n = vaddq_f32(k, l);
217 float32x4_t y = vaddq_f32(m, n);
218
219 o->xm1 = x3;
220 o->xm2 = x2;
221 o->ym1 = vdupq_n_f32(y[3]);
222 o->ym2 = vdupq_n_f32(y[2]);
223
224 *bOut = y;
225 #else // HV_SIMD_NONE
226 float y = o->b0*bIn + o->b1*o->xm1 + o->b2*o->xm2 - o->a1*o->ym1 - o->a2*o->ym2;
227 o->xm2 = o->xm1;
228 o->xm1 = bIn;
229 o->ym2 = o->ym1;
230 o->ym1 = y;
231 *bOut = y;
232 #endif
233 }
234
235 #endif // _HEAVY_SIGNAL_BIQUAD_H_