1/*
2 * Copyright (C) 2001-2002 Jonas Beskow <beskow@speech.kth.se>
3 *
4 * This file is part of the Snack Sound Toolkit.
5 * The latest version can be found at http://www.speech.kth.se/snack/
6 *
7 * This program is free software; you can redistribute it and/or modify
8 * it under the terms of the GNU General Public License as published by
9 * the Free Software Foundation; either version 2 of the License, or
10 * (at your option) any later version.
11 *
12 * This program is distributed in the hope that it will be useful,
13 * but WITHOUT ANY WARRANTY; without even the implied warranty of
14 * MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE.  See the
15 * GNU General Public License for more details.
16 *
17 * You should have received a copy of the GNU General Public License
18 * along with this program; if not, write to the Free Software
19 * Foundation, Inc., 675 Mass Ave, Cambridge, MA 02139, USA.
20 */
21
22#include "snack.h"
23#include <math.h>
24#include <stdlib.h>
25#include <string.h>
26
27#if defined(WIN) || defined(MAC)
28#  define M_PI 3.14159265358979323846
29#endif
30
31/* linear interpolation macro */
32#define LIN(X,FRAC) (X[0]+(X[1]-X[0])*FRAC)
33
34struct formantFilter {
35  configProc *configProc;
36  startProc  *startProc;
37  flowProc   *flowProc;
38  freeProc   *freeProc;
39  Tcl_Interp *interp;
40  Snack_Filter prev, next;
41  Snack_StreamInfo si;
42  double dataRatio;
43  int    reserved[4];
44  /* formant filter params */
45  double bw;
46  double freq;
47  double a0,b0,c0; /* saved coeffs */
48  float mem[2];    /* last two samples */
49} formantFilter;
50
51typedef struct formantFilter *formantFilter_t;
52
53static int formantConfigProc(Snack_Filter f, Tcl_Interp *interp, int objc,
54			     Tcl_Obj *CONST objv[]) {
55  formantFilter_t mf = (formantFilter_t) f;
56  /*  fprintf(stderr,"formantConfigProc\n");*/
57
58  switch (objc) {
59  case 1:
60    if (Tcl_GetDoubleFromObj(interp, objv[0], &(mf->freq)) != TCL_OK) {
61      return TCL_ERROR;
62    }
63    break;
64  case 2:
65    if (Tcl_GetDoubleFromObj(interp, objv[0], &(mf->freq)) != TCL_OK) {
66      return TCL_ERROR;
67    }
68    if (Tcl_GetDoubleFromObj(interp, objv[1], &(mf->bw)) != TCL_OK) {
69      return TCL_ERROR;
70    }
71    break;
72  default:
73    Tcl_SetResult(interp,"wrong # args. should be \"filter configure freq ?bandwidth?\"",TCL_STATIC);
74    return TCL_ERROR;
75  }
76
77  return TCL_OK;
78}
79
80static void calcCoeffs(double f,double bw, double *a, double *b, double *c) {
81  double tmp = exp(-M_PI*bw);
82  *c = -tmp*tmp;
83  *b = 2.0*tmp*cos(2*M_PI*f);
84  *a = 1.0 - *b - *c;
85}
86
87static Snack_Filter
88formantCreateProc(Tcl_Interp *interp, int objc, Tcl_Obj *CONST objv[])
89{
90  formantFilter_t mf;
91  /*  fprintf(stderr,"formantCreateProc\n");*/
92
93  mf = (formantFilter_t) ckalloc(sizeof(formantFilter));
94
95  mf->freq = 0;
96  mf->bw = 1;
97
98  if (formantConfigProc((Snack_Filter) mf, interp, objc, objv) != TCL_OK) {
99    return (Snack_Filter) NULL;
100  }
101  return (Snack_Filter) mf;
102}
103
104static int
105formantStartProc(Snack_Filter f, Snack_StreamInfo si)
106{
107  formantFilter_t mf = (formantFilter_t) f;
108
109  /* only mono sounds... */
110  if (si->outWidth != 1) return TCL_ERROR;
111
112  calcCoeffs(1.0*mf->freq/si->rate,
113	     1.0*mf->bw/si->rate,
114	     &mf->a0,&mf->b0,&mf->c0);
115  mf->mem[0] = 0.0;
116  mf->mem[1] = 0.0;
117  /*  fprintf(stderr,"formantStartProc, width = %d\n",si->outWidth);*/
118
119  return TCL_OK;
120}
121
122static int
123formantFlowProc(Snack_Filter f, Snack_StreamInfo si, float *in, float *out,
124		int *inFrames, int *outFrames) {
125  int i,nf = 0;
126  double a[2],b[2],c[2],frac,z;
127  formantFilter_t mf = (formantFilter_t) f;
128
129  /*fprintf(stderr,"entering   formantFlowProc, *inFrames = %d\t*outFrames = %d, in = %p, out = %p\n",*inFrames,*outFrames,in,out); */
130
131  calcCoeffs(1.0*mf->freq/si->rate,
132	     1.0*mf->bw/si->rate,
133	     &a[1],&b[1],&c[1]);
134  a[0]=mf->a0;b[0]=mf->b0;c[0]=mf->c0;
135  /* only mono sounds... */
136  if (si->outWidth != 1) {
137    *outFrames = 0;
138    *inFrames = 0;
139    return TCL_ERROR;
140  }
141
142  nf = (*inFrames)<(*outFrames)?(*inFrames):(*outFrames);
143  if (nf != 0) {
144    z = 1.0/nf;
145    if (nf>=1) {
146      out[0] = (float)(LIN(a,0*z)*in[0]+LIN(b,0*z)*mf->mem[0]+LIN(c,0*z)*mf->mem[1]);
147    }
148    if (nf>=2) {
149      out[1] = (float)(LIN(a,1*z)*in[1]+LIN(b,1*z)*out[0]+LIN(c,1*z)*mf->mem[0]);
150    }
151    for (i=2;i<nf;i++) {
152      frac = i*z;
153      out[i] = (float)(LIN(a,frac)*in[i] + LIN(b,frac)*out[i-1] + LIN(c,frac)*out[i-2]);
154    }
155  }
156  if (nf>=1) mf->mem[0]=out[nf-1];
157  if (nf>=2) mf->mem[1]=out[nf-2];
158  mf->a0=a[1];
159  mf->b0=b[1];
160  mf->c0=c[1];
161
162  *outFrames = nf;
163  *inFrames = nf;
164  /*fprintf(stderr,"leaving    formantFlowProc, *inFrames = %d\t*outFrames = %d\n",*inFrames,*outFrames); */
165  return TCL_OK;
166}
167
168static void
169formantFreeProc(Snack_Filter f) {
170  ckfree((char*)f);
171}
172
173static Snack_FilterType snackFormantType = {
174  "formant",
175  formantCreateProc,
176  formantConfigProc,
177  formantStartProc,
178  formantFlowProc,
179  formantFreeProc,
180  (Snack_FilterType *) NULL
181};
182
183
184#define RECTANGLE 1
185#define TRIANGLE 2
186#define SINE 3
187#define NOISE 4
188#define SAMPLED 5
189
190
191#define MAX_SAMPLES 1600
192
193struct generatorFilter {
194  /* general filter stuff */
195  configProc *configProc;
196  startProc  *startProc;
197  flowProc   *flowProc;
198  freeProc   *freeProc;
199  Tcl_Interp *interp;
200  Snack_Filter prev, next;
201  Snack_StreamInfo si;
202  double dataRatio;
203  int    reserved[4];
204  /* generator filter stuff */
205  double freq[2];
206  double ampl[2];
207  double shape[2];     /* pulse width or peak */
208  int type;            /* TRIANGLE, RECTANGLE, SINE, NOISE or SAMPLED*/
209  double _phase;
210  float samples[MAX_SAMPLES];
211  float maxval;
212  int nSamples;
213  int ntot;
214  int ngen;
215} generatorFilter;
216
217typedef struct generatorFilter *generatorFilter_t;
218
219
220
221static int generatorConfigProc(Snack_Filter f, Tcl_Interp *interp, int objc,
222			     Tcl_Obj *CONST objv[]) {
223  generatorFilter_t mf = (generatorFilter_t) f;
224  char *typestr;
225
226  /*  fprintf(stderr,"generatorConfigProc, objc = %d\n",objc);*/
227
228  /*
229    syntax: generator configure freq ?ampl? ?shape? ?type? ?ntot?
230  */
231  switch (objc) {
232  case 5:
233    if (Tcl_GetIntFromObj(interp,objv[4],&mf->ntot)==TCL_ERROR) {
234      return TCL_ERROR;
235    }
236
237  case 4:
238    typestr = Tcl_GetStringFromObj(objv[3], NULL);
239    if (strncmp(typestr,"rec",3)==0) {
240      mf->type = RECTANGLE;
241    } else if (strncmp(typestr,"tri",3)==0) {
242      mf->type = TRIANGLE;
243    } else if (strncmp(typestr,"sin",3)==0) {
244      mf->type = SINE;
245    } else if (strncmp(typestr,"noi",3)==0) {
246      mf->type = NOISE;
247    } else if (strncmp(typestr,"sam",3)==0) {
248      mf->type = SAMPLED;
249    } else {
250      Tcl_SetResult(interp,
251		    "bad waveform type, must be rectangle, triangle, sine, noise or sampled",
252		    TCL_STATIC);
253      return TCL_ERROR;
254    }
255
256  case 3:
257    if (Tcl_GetDoubleFromObj(interp,objv[2],&mf->shape[1])==TCL_ERROR) {
258      return TCL_ERROR;
259    }
260
261  case 2:
262    if (Tcl_GetDoubleFromObj(interp,objv[1],&mf->ampl[1])==TCL_ERROR) {
263      return TCL_ERROR;
264    }
265
266  case 1:
267    if (Tcl_GetDoubleFromObj(interp,objv[0],&mf->freq[1])==TCL_ERROR) {
268      return TCL_ERROR;
269    }
270    break;
271
272  default:
273    Tcl_SetResult(interp, "wrong # args, should be \"generator configure freq ?ampl? ?shape? ?type?\"", TCL_STATIC);
274    return TCL_ERROR;
275  }
276
277  return TCL_OK;
278}
279
280static Snack_Filter
281generatorCreateProc(Tcl_Interp *interp, int objc, Tcl_Obj *CONST objv[]) {
282  generatorFilter_t mf;
283  /*  fprintf(stderr,"generatorCreateProc, objc = %d\n",objc);*/
284  mf = (generatorFilter_t) ckalloc(sizeof(generatorFilter));
285  mf->freq[1] = 0;
286  mf->ampl[1] = 0;
287  mf->type = SINE;
288  mf->shape[1] = 0.5;
289  mf->_phase = 0.0;
290  mf->ntot = -1;
291  if (generatorConfigProc((Snack_Filter) mf, interp, objc, objv) != TCL_OK) {
292    return (Snack_Filter) NULL;
293  }
294return (Snack_Filter) mf;
295}
296
297
298static int
299generatorStartProc(Snack_Filter f, Snack_StreamInfo si) {
300
301  generatorFilter_t mf = (generatorFilter_t) f;
302  /*  fprintf(stderr,"entering generatorStartProc\n"); */
303
304  mf->freq[0]  = mf->freq[1];
305  mf->ampl[0]  = mf->ampl[1];
306  mf->shape[0] = mf->shape[1];
307  mf->_phase   = 0.0;
308  mf->nSamples = 0;
309  mf->ngen     = 0;
310  mf->maxval   = 1.0;
311
312  return TCL_OK;
313}
314
315static int
316generatorFlowProc(Snack_Filter f, Snack_StreamInfo si, float *in, float *out,
317		  int *inFrames, int *outFrames) {
318
319  generatorFilter_t mf = (generatorFilter_t) f;
320  int i, i0, i1, ii, fr, wi;
321  double y = 1.0/RAND_MAX, z = 1.0/(*outFrames), ph = mf->_phase;
322  double frac,a,b,v[2],tmp;
323
324  /*  fprintf(stderr,"entering generatorFlowProc, *inFrames = %d\t*outFrames = %d, in = %p, out = %p\n",*inFrames,*outFrames,in,out); */
325
326  if (mf->ntot>0 && mf->ngen + (*outFrames) > mf->ntot) {
327    *outFrames = mf->ntot - mf->ngen;
328  }
329
330  /* read input samples (if any) - this is used by the SAMPLED wave type*/
331  for (i=0;i<*inFrames;i++) {
332    ii = i+mf->nSamples;
333    if (ii>=MAX_SAMPLES) break;
334    tmp = mf->samples[ii] = in[i];
335    if (fabs(tmp) > mf->maxval) { mf->maxval = (float) fabs(tmp);}
336  }
337  mf->nSamples += i;
338  *inFrames = i;
339
340  i = 0;
341  switch (mf->type) {
342  case NOISE:
343
344    for (fr = 0; fr < *outFrames; fr++) {
345      frac = fr*z;
346      for (wi = 0; wi < si->outWidth; wi++) {
347	out[i++]=(float)(LIN(mf->ampl,frac)*2*(y*rand()-.5));
348      }
349      i += (si->streamWidth - si->outWidth);
350    }
351    *inFrames = 0;
352    break;
353
354  case RECTANGLE:
355    for (fr = 0; fr < *outFrames; fr++) {
356      frac = fr*z;
357      ph = fmod (ph + LIN(mf->freq,frac)/si->rate,1.0);
358      for (wi = 0; wi < si->outWidth; wi++) {
359	out[i++] = (float)(LIN(mf->ampl,frac)*(ph<LIN(mf->shape,frac)?-1:1));
360      }
361      i += (si->streamWidth - si->outWidth);
362    }
363    *inFrames = 0;
364    break;
365
366  case TRIANGLE:
367    for (fr = 0; fr < *outFrames; fr++) {
368      frac = fr*z;
369      ph = fmod (ph + LIN(mf->freq,frac)/si->rate,1.0);
370      for (wi = 0; wi < si->outWidth; wi++) {
371	if (ph < LIN(mf->shape,frac)) {
372	  out[i++] = (float)(LIN(mf->ampl,frac)*(-1+2*(ph)/LIN(mf->shape,frac)));
373	} else if (ph > LIN(mf->shape,frac)) {
374	  out[i++] = (float)(LIN(mf->ampl,frac)*(1-2*(ph-LIN(mf->shape,frac))/(1-LIN(mf->shape,frac))));
375	} else {
376	  out[i++] = (float)LIN(mf->ampl,frac);
377	}
378      }
379      i += (si->streamWidth - si->outWidth);
380    }
381    *inFrames = 0;
382    break;
383
384  case SINE:
385    for (fr = 0; fr < *outFrames; fr++) {
386      frac = fr*z;
387      ph = fmod (ph + LIN(mf->freq,frac)/si->rate,1.0);
388      a = sin(ph * 2 * M_PI);
389      b = 2*LIN(mf->shape,frac)-1;
390      a = a>b?a:b;
391      for (wi = 0; wi < si->outWidth; wi++) {
392	if (1-b==0.0) {
393	  out[i++] = 0.0;
394	} else {
395	  out[i++] = (float)(LIN(mf->ampl,frac)*(a-.5-.5*b)/(1-b));
396	}
397      }
398      i += (si->streamWidth - si->outWidth);
399      *inFrames = 0;
400    }
401    break;
402
403  case SAMPLED:
404    if (mf->nSamples > 0) {
405      for (fr = 0; fr < *outFrames; fr++) {
406	frac = fr*z;
407	ph = fmod (ph + LIN(mf->freq,frac)/si->rate,1.0);
408	a = ph*(mf->nSamples);
409	i0 = (int)floor(a);
410	i1 = (int)ceil(a)%mf->nSamples;
411	v[0] = mf->samples[i0];
412	v[1] = mf->samples[i1];
413	frac = (a-i0);
414	for (wi = 0; wi < si->outWidth; wi++) {
415	  out[i++] = (float)(LIN(v,frac)*LIN(mf->ampl,frac)/mf->maxval);
416	}
417	i += (si->streamWidth - si->outWidth);
418      }
419    } else {
420      for (fr = 0; fr < *outFrames; fr++) {
421	for (wi = 0; wi < si->outWidth; wi++) {
422	  out[i++] = 0;
423	}
424	i += (si->streamWidth - si->outWidth);
425      }
426    }
427  }
428  mf->_phase = ph; /* save current phase value */
429  mf->freq[0] = mf->freq[1];
430  mf->ampl[0] = mf->ampl[1];
431  mf->shape[0] = mf->shape[1];
432
433  mf->ngen += *outFrames;
434
435  /*  fprintf(stderr,"leaving  generatorFlowProc, *inFrames = %d\t*outFrames = %d (i=%d)\n-----------------------------------------\n",*inFrames,*outFrames,i); */
436  return TCL_OK;
437}
438
439
440static void
441generatorFreeProc(Snack_Filter f) {
442  ckfree((char*)f);
443}
444
445static Snack_FilterType snackGeneratorType = {
446  "generator",
447  generatorCreateProc,
448  generatorConfigProc,
449  generatorStartProc,
450  generatorFlowProc,
451  generatorFreeProc,
452  (Snack_FilterType *) NULL
453};
454
455
456void
457createSynthesisFilters() {
458  Snack_CreateFilterType(&snackGeneratorType);
459  Snack_CreateFilterType(&snackFormantType);
460}
461