66 float amp,
float ampoffset){
72 for(i=0;i<
m;i++)lsp[i]=vorbis_coslook(lsp[i]);
80 float w=vorbis_coslook(wdel*k);
103 q=vorbis_fromdBlook(amp*
105 vorbis_invsq2explook(qexp+m)-
122 static const int MLOOP_1[64]={
123 0,10,11,11, 12,12,12,12, 13,13,13,13, 13,13,13,13,
124 14,14,14,14, 14,14,14,14, 14,14,14,14, 14,14,14,14,
125 15,15,15,15, 15,15,15,15, 15,15,15,15, 15,15,15,15,
126 15,15,15,15, 15,15,15,15, 15,15,15,15, 15,15,15,15,
129 static const int MLOOP_2[64]={
130 0,4,5,5, 6,6,6,6, 7,7,7,7, 7,7,7,7,
131 8,8,8,8, 8,8,8,8, 8,8,8,8, 8,8,8,8,
132 9,9,9,9, 9,9,9,9, 9,9,9,9, 9,9,9,9,
133 9,9,9,9, 9,9,9,9, 9,9,9,9, 9,9,9,9,
136 static const int MLOOP_3[8]={0,1,2,2,3,3,3,3};
141 float amp,
float ampoffset){
147 int ampoffseti=rint(ampoffset*4096.
f);
148 int ampi=rint(amp*16.
f);
149 long *ilsp=alloca(m*
sizeof(*ilsp));
150 for(i=0;i<
m;i++)ilsp[i]=vorbis_coslook_i(lsp[i]/
M_PI*65536.
f+.5
f);
155 unsigned long pi=46341;
156 unsigned long qi=46341;
158 long wi=vorbis_coslook_i(k*65536/ln);
160 qi*=labs(ilsp[0]-wi);
161 pi*=labs(ilsp[1]-wi);
164 if(!(shift=MLOOP_1[(pi|qi)>>25]))
165 if(!(shift=MLOOP_2[(pi|qi)>>19]))
166 shift=MLOOP_3[(pi|qi)>>16];
167 qi=(qi>>shift)*labs(ilsp[j-1]-wi);
168 pi=(pi>>shift)*labs(ilsp[j]-wi);
171 if(!(shift=MLOOP_1[(pi|qi)>>25]))
172 if(!(shift=MLOOP_2[(pi|qi)>>19]))
173 shift=MLOOP_3[(pi|qi)>>16];
180 qi=(qi>>shift)*labs(ilsp[j-1]-wi);
184 if(!(shift=MLOOP_1[(pi|qi)>>25]))
185 if(!(shift=MLOOP_2[(pi|qi)>>19]))
186 shift=MLOOP_3[(pi|qi)>>16];
190 qexp+=shift-14*((m+1)>>1);
196 pi*=(1<<14)-((wi*wi)>>14);
227 while(qi && !(qi&0x8000)){
231 amp=vorbis_fromdBlook_i(ampi*
232 vorbis_invsqlook_i(qi,qexp)-
237 while(map[++i]==k)curve[
i]*=amp;
249 float amp,
float ampoffset){
252 for(i=0;i<
m;i++)lsp[i]=2.
f*
cos(lsp[i]);
259 float w=2.f*
cos(wdel*k);
276 q=
fromdB(amp/sqrt(p+q)-ampoffset);
279 while(map[++i]==k)curve[
i]*=
q;
290 for(i=2; i<= ord; i++) {
291 for(j=ord; j >=
i; j--) {
298 static int comp(
const void *
a,
const void *
b){
299 return (*(
float *)a<*(
float *)b)-(*(
float *)a>*(
float *)
b);
309 #define EPSILON 10e-7
312 double lastdelta=0.f;
313 double *defl=alloca(
sizeof(*defl)*(ord+1));
314 for(i=0;i<=ord;i++)defl[i]=a[i];
317 double new=0.f,delta;
321 double p=defl[
m],
pp=0.f,ppp=0.f,denom;
327 p =
new*p + defl[i-1];
331 denom=(m-1) * ((m-1)*pp*pp - m*p*ppp);
336 denom = pp + sqrt(denom);
339 denom = pp - sqrt(denom);
346 if(delta<0.
f)delta*=-1;
348 if(
fabs(delta/
new)<10
e-12)
break;
357 defl[i-1]+=
new*defl[i];
369 double *root=alloca(ord*
sizeof(*root));
371 for(i=0; i<ord;i++) root[i] = r[i];
376 for(i=0; i<ord; i++) {
378 double rooti=root[
i];
380 for(k=ord-1; k>= 0; k--) {
383 p = p * rooti + a[
k];
391 if(count>40)
return(-1);
399 for(i=0; i<ord;i++) r[i] = root[i];
407 int g1_order,g2_order;
408 float *g1=alloca(
sizeof(*g1)*(order2+1));
409 float *g2=alloca(
sizeof(*g2)*(order2+1));
410 float *g1r=alloca(
sizeof(*g1r)*(order2+1));
411 float *g2r=alloca(
sizeof(*g2r)*(order2+1));
424 for(i=1;i<=g1_order;i++) g1[g1_order-i] = lpc[i-1]+lpc[m-i];
426 for(i=1;i<=g2_order;i++) g2[g2_order-i] = lpc[i-1]-lpc[m-i];
428 if(g1_order>g2_order){
429 for(i=2; i<=g2_order;i++) g2[g2_order-i] += g2[g2_order-i+2];
431 for(i=1; i<=g1_order;i++) g1[g1_order-i] -= g1[g1_order-i+1];
432 for(i=1; i<=g2_order;i++) g2[g2_order-i] += g2[g2_order-i+1];
450 for(i=0;i<g1_order;i++)
451 lsp[i*2] = acos(g1r[i]);
453 for(i=0;i<g2_order;i++)
454 lsp[i*2+1] = acos(g2r[i]);
static int Newton_Raphson(float *a, int ord, float *r)
return Display return Display Bool Bool int int e
static void cheby(float *g, int ord)
GLboolean GLboolean GLboolean GLboolean a
static int comp(const void *a, const void *b)
#define vorbis_fpu_setround(vorbis_fpu_control)
int vorbis_lpc_to_lsp(float *lpc, float *lsp, int m)
void vorbis_lsp_to_curve(float *curve, int *map, int n, int ln, float *lsp, int m, float amp, float ampoffset)
#define vorbis_fpu_restore(vorbis_fpu_control)
GLdouble GLdouble GLdouble GLdouble q
static int Laguerre_With_Deflation(float *a, int ord, float *r)
GLdouble GLdouble GLdouble r
GLdouble GLdouble GLdouble b
GLint GLint GLint GLint GLint w
const GLubyte GLuint GLuint GLuint GLuint alpha GLboolean GLboolean GLboolean GLboolean alpha GLint GLint GLsizei GLsizei GLenum type GLenum GLint GLenum GLint GLint GLsizei GLsizei GLint border GLenum GLint GLint GLint GLint GLint GLsizei GLsizei height GLsizei GLsizei GLenum GLenum const GLvoid *pixels GLenum GLint GLint GLint GLint j2 GLdouble GLdouble GLdouble GLdouble GLdouble GLdouble zFar GLenum GLenum GLint *params GLenum GLenum GLint *params GLenum GLenum GLint *params GLenum GLenum GLfloat *params GLenum GLint GLenum GLenum GLvoid *pixels GLenum GLint GLenum GLint *params GLenum GLenum GLint *params GLenum GLsizei const GLvoid *pointer GLenum GLenum const GLint *params GLenum GLfloat GLfloat GLint GLint const GLfloat *points GLenum GLfloat GLfloat GLint GLint GLfloat GLfloat GLint GLint const GLfloat *points GLint GLfloat GLfloat GLint GLfloat GLfloat v2 GLenum GLenum const GLint *params GLdouble GLdouble GLdouble GLdouble GLdouble GLdouble zFar GLenum map
void qsort(void *base, size_t nmemb, size_t size, int(*compare)(const void *, const void *))