-/* One can do this the long way by generating the transfer function in
- the time domain and taking the forward FFT of the result. The
- results from direct calculation are cleaner and faster.
-
- This version does a linear curve generation and then later
- interpolates the log curve from the linear curve. This could stand
- optimization; it could both be more precise as well as not compute
- quite a few unused values */
-
-void _vlpc_de_helper(double *curve,double *lpc,double amp,
- lpc_lookup *l){
- int i;
- memset(curve,0,sizeof(double)*l->ln*2);
- if(amp==0)return;
-
- for(i=0;i<l->m;i++){
- curve[i*2+1]=lpc[i]/4/amp;
- curve[i*2+2]=-lpc[i]/4/amp;
- }
-
- drft_backward(&l->fft,curve); /* reappropriated ;-) */
-
- {
- int l2=l->ln*2;
- double unit=1./amp;
- curve[0]=(1./(curve[0]*2+unit));
- for(i=1;i<l->ln;i++){
- double real=(curve[i]+curve[l2-i]);
- double imag=(curve[i]-curve[l2-i]);
- curve[i]=(1./hypot(real+unit,imag));
- }
- }
-}
-
-
-/* generate the whole freq response curve of an LPC IIR filter */
-
-void vorbis_lpc_to_curve(double *curve,double *lpc,double amp,lpc_lookup *l){
- double *lcurve=alloca(sizeof(double)*(l->ln*2));
- int i;
-
- if(amp==0){
- memset(curve,0,sizeof(double)*l->n);
- return;
- }
- _vlpc_de_helper(lcurve,lpc,amp,l);
-
-#ifdef ANALYSIS
- {
- static int frameno=0;
- int j;
- FILE *out;
- char buffer[80];
-
- sprintf(buffer,"loglpc%d.m",frameno++);
- out=fopen(buffer,"w+");
- for(j=0;j<l->ln;j++)
- fprintf(out,"%g\n",lcurve[j]);
- fclose(out);
-
-#endif
-
- for(i=0;i<l->ln;i++)lcurve[i]/=l->barknorm[i];
- for(i=0;i<l->n;i++)curve[i]=lcurve[l->linearmap[i]];
-
-#ifdef ANALYSIS
-
- sprintf(buffer,"lpc%d.m",frameno-1);
- out=fopen(buffer,"w+");
- for(j=0;j<l->n;j++)
- fprintf(out,"%g\n",curve[j]);
- fclose(out);
- }
-#endif
-}
-
-/* subtract or add an lpc filter to data. Vorbis doesn't actually use this. */
-
-void vorbis_lpc_residue(double *coeff,double *prime,int m,
- double *data,long n){
-
- /* in: coeff[0...m-1] LPC coefficients
- prime[0...m-1] initial values
- data[0...n-1] data samples
- out: data[0...n-1] residuals from LPC prediction */
-
- long i,j;
- double *work=alloca(sizeof(double)*(m+n));
- double y;
-
- if(!prime)
- for(i=0;i<m;i++)
- work[i]=0;
- else
- for(i=0;i<m;i++)
- work[i]=prime[i];