1 /*
2 * Talkbox
3 *
4 * This module is ported from the mdaTalkbox plugin by Paul Kellet
5 * (maxim digital audio).
6 *
7 * More information on his plugins and the original code can be found here:
8 *
9 * http://mda.smartelectronix.com/
10 *
11 */
12
13 #include <math.h>
14 #include <stdlib.h>
15 #include <string.h>
16 #include "soundpipe.h"
17
18 #ifndef TWO_PI
19 #define TWO_PI 6.28318530717958647692528676655901
20 #endif
21
22 #define ORD_MAX 50
23
lpc_durbin(SPFLOAT * r,int p,float * k,float * g)24 static void lpc_durbin(SPFLOAT *r, int p, float *k, float *g)
25 {
26 int i, j;
27 SPFLOAT a[ORD_MAX], at[ORD_MAX], e=r[0];
28
29 for(i=0; i<=p; i++) a[i] = at[i] = 0.0f;
30
31 for(i=1; i<=p; i++)
32 {
33 k[i] = -r[i];
34
35 for(j=1; j<i; j++)
36 {
37 at[j] = a[j];
38 k[i] -= a[j] * r[i-j];
39 }
40 if(fabs(e) < 1.0e-20f) { e = 0.0f; break; }
41 k[i] /= e;
42
43 a[i] = k[i];
44 for(j=1; j<i; j++) a[j] = at[j] + k[i] * at[i-j];
45
46 e *= 1.0f - k[i] * k[i];
47 }
48
49 if(e < 1.0e-20f) e = 0.0f;
50 *g = (float)sqrt(e);
51 }
52
lpc(float * buf,float * car,uint32_t n,uint32_t o)53 static void lpc(float *buf, float *car, uint32_t n, uint32_t o)
54 {
55 SPFLOAT z[ORD_MAX], r[ORD_MAX], k[ORD_MAX], G, x;
56 uint32_t i, j, nn=n;
57 SPFLOAT min;
58
59 /* buf[] is already emphasized and windowed */
60 for(j=0; j<=o; j++, nn--) {
61 z[j] = r[j] = 0.0f;
62 for(i=0; i<nn; i++) r[j] += buf[i] * buf[i+j]; /* autocorrelation */
63 }
64
65 r[0] *= 1.001f; /* stability fix */
66
67 min = 0.00001f;
68 if(r[0] < min) { for(i=0; i<n; i++) buf[i] = 0.0f; return; }
69
70 lpc_durbin(r, o, k, &G); /* calc reflection coeffs */
71
72 for(i=1; i<=o; i++) {
73 if(k[i] > 0.995f) k[i] = 0.995f; else if(k[i] < -0.995f) k[i] = -.995f;
74 }
75
76 for(i=0; i<n; i++) {
77 x = G * car[i];
78 /* lattice filter */
79 for(j=o; j>0; j--) {
80 x -= k[j] * z[j-1];
81 z[j] = z[j-1] + k[j] * x;
82 }
83 buf[i] = z[0] = x; /* output buf[] will be windowed elsewhere */
84 }
85 }
86
sp_talkbox_create(sp_talkbox ** p)87 int sp_talkbox_create(sp_talkbox **p)
88 {
89 *p = malloc(sizeof(sp_talkbox));
90 return SP_OK;
91 }
92
sp_talkbox_destroy(sp_talkbox ** p)93 int sp_talkbox_destroy(sp_talkbox **p)
94 {
95 free(*p);
96 return SP_OK;
97 }
98
sp_talkbox_init(sp_data * sp,sp_talkbox * p)99 int sp_talkbox_init(sp_data *sp, sp_talkbox *p)
100 {
101 uint32_t n;
102 p->quality = 1.f;
103 p->N = 1;
104 p->K = 0;
105
106 n = (uint32_t)(0.01633f * sp->sr);
107 if(n > SP_TALKBOX_BUFMAX) n = SP_TALKBOX_BUFMAX;
108
109 /* calculate hanning window */
110 if(n != p->N) {
111 p->N = n;
112 SPFLOAT dp = TWO_PI / (SPFLOAT)p->N;
113 SPFLOAT pos = 0.0f;
114 for(n=0; n < p->N; n++) {
115 p->window[n] = 0.5f - 0.5f * (SPFLOAT)cos(pos);
116 pos += dp;
117 }
118 }
119
120 /* zero out variables and buffers */
121 p->pos = p->K = 0;
122 p->emphasis = 0.0f;
123 p->FX = 0;
124
125 p->u0 = p->u1 = p->u2 = p->u3 = p->u4 = 0.0f;
126 p->d0 = p->d1 = p->d2 = p->d3 = p->d4 = 0.0f;
127
128 memset(p->buf0, 0, SP_TALKBOX_BUFMAX * sizeof(SPFLOAT));
129 memset(p->buf1, 0, SP_TALKBOX_BUFMAX * sizeof(SPFLOAT));
130 memset(p->car0, 0, SP_TALKBOX_BUFMAX * sizeof(SPFLOAT));
131 memset(p->car1, 0, SP_TALKBOX_BUFMAX * sizeof(SPFLOAT));
132 return SP_OK;
133 }
134
sp_talkbox_compute(sp_data * sp,sp_talkbox * t,SPFLOAT * src,SPFLOAT * exc,SPFLOAT * out)135 int sp_talkbox_compute(sp_data *sp, sp_talkbox *t, SPFLOAT *src, SPFLOAT *exc, SPFLOAT *out)
136 {
137 int32_t p0=t->pos, p1 = (t->pos + t->N/2) % t->N;
138 SPFLOAT e=t->emphasis, w, o, x, fx=t->FX;
139 SPFLOAT p, q, h0=0.3f, h1=0.77f;
140 SPFLOAT den;
141
142 t->O = (uint32_t)((0.0001f + 0.0004f * t->quality) * sp->sr);
143
144 o = *src;
145 x = *exc;
146
147 p = t->d0 + h0 * x;
148 t->d0 = t->d1;
149 t->d1 = x - h0 * p;
150
151 q = t->d2 + h1 * t->d4;
152 t->d2 = t->d3;
153 t->d3 = t->d4 - h1 * q;
154
155 t->d4 = x;
156
157 x = p + q;
158
159 if(t->K++) {
160 t->K = 0;
161
162 /* carrier input */
163 t->car0[p0] = t->car1[p1] = x;
164
165 /* 6dB/oct pre-emphasis */
166 x = o - e; e = o;
167
168 /* 50% overlapping hanning windows */
169 w = t->window[p0]; fx = t->buf0[p0] * w; t->buf0[p0] = x * w;
170 if(++p0 >= t->N) { lpc(t->buf0, t->car0, t->N, t->O); p0 = 0; }
171
172 w = 1.0f - w; fx += t->buf1[p1] * w; t->buf1[p1] = x * w;
173 if(++p1 >= t->N) { lpc(t->buf1, t->car1, t->N, t->O); p1 = 0; }
174 }
175
176 p = t->u0 + h0 * fx;
177 t->u0 = t->u1;
178 t->u1 = fx - h0 * p;
179
180 q = t->u2 + h1 * t->u4;
181 t->u2 = t->u3;
182 t->u3 = t->u4 - h1 * q;
183
184 t->u4 = fx;
185 x = p + q;
186
187 o = x * 0.5;
188 *out = o;
189 t->emphasis = e;
190 t->pos = p0;
191 t->FX = fx;
192
193 den = 1.0e-10f;
194 /* anti-denormal */
195 if(fabs(t->d0) < den) t->d0 = 0.0f;
196 if(fabs(t->d1) < den) t->d1 = 0.0f;
197 if(fabs(t->d2) < den) t->d2 = 0.0f;
198 if(fabs(t->d3) < den) t->d3 = 0.0f;
199 if(fabs(t->u0) < den) t->u0 = 0.0f;
200 if(fabs(t->u1) < den) t->u1 = 0.0f;
201 if(fabs(t->u2) < den) t->u2 = 0.0f;
202 if(fabs(t->u3) < den) t->u3 = 0.0f;
203 return SP_OK;
204 }
205