- int j;
- FILE *out;
- char buffer[80];
- static int frameno=0;
-
- sprintf(buffer,"preloglpc%d.m",frameno++);
- out=fopen(buffer,"w+");
- for(j=0;j<l->ln;j++)
- fprintf(out,"%g\n",work[j]);
- fclose(out);
- }
-#endif
-
- return vorbis_lpc_from_spectrum(work,lpc,l);
-}
-
-
-/* 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 on 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;
-
- _vlpc_de_helper(lcurve,lpc,amp,l);
-
-#ifdef ANALYSIS
- {
- int j;
- FILE *out;
- char buffer[80];
- static int frameno=0;
-
- 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
-
- if(amp==0)return;
-
- for(i=0;i<l->n;i++){
- int ii=l->iscale[i];
- curve[i]=((1.-l->ifrac[i])*lcurve[ii]+
- l->ifrac[i]*lcurve[ii+1])*l->norm[i];
- }
-
-}
-
-void vorbis_lpc_apply(double *residue,double *lpc,double amp,lpc_lookup *l){
- double *lcurve=alloca(sizeof(double)*((l->ln+l->n)*2));
- int i;
-
- if(amp==0){
- memset(residue,0,l->n*sizeof(double));
- }else{
-
- _vlpc_de_helper(lcurve,lpc,amp,l);
-
- for(i=0;i<l->n;i++){
- if(residue[i]!=0){
- int ii=l->iscale[i];
- residue[i]*=((1.-l->ifrac[i])*lcurve[ii]+
- l->ifrac[i]*lcurve[ii+1])*l->norm[i];
- }