most of leakage estimation converted
authorjm <jm@0101bb08-14d6-0310-b084-bc0e0c8e3800>
Sat, 17 Dec 2005 02:26:28 +0000 (02:26 +0000)
committerjm <jm@0101bb08-14d6-0310-b084-bc0e0c8e3800>
Sat, 17 Dec 2005 02:26:28 +0000 (02:26 +0000)
git-svn-id: http://svn.xiph.org/trunk/speex@10609 0101bb08-14d6-0310-b084-bc0e0c8e3800

libspeex/mdf.c
libspeex/pseudofloat.h

index efaebcd..7607308 100644 (file)
 #endif
 
 #ifdef FIXED_POINT
-static const spx_float_t MAX_ALPHA = {16777, -21};
-static const spx_float_t ALPHA0 = {26214, -19};
-static const spx_float_t MIN_LEAK = {16777, -24};
-static const spx_float_t SPEC_AVERAGE = {20972, -20};
-static const spx_float_t SPEC_AVERAGE_1 = {32113,-15};
+static const spx_float_t MAX_ALPHA = ((spx_float_t){16777, -21});
+static const spx_float_t ALPHA0 = ((spx_float_t){26214, -19});
+static const spx_float_t MIN_LEAK = ((spx_float_t){16777, -24});
+static const spx_float_t SPEC_AVERAGE = ((spx_float_t){20972, -20});
+static const spx_float_t SPEC_AVERAGE_1 = ((spx_float_t){32113,-15});
 #else
-static const spx_float_t MAX_ALPHA = .008;
-static const spx_float_t ALPHA0 = .05;
-static const spx_float_t MIN_LEAK = .001;
-static const spx_float_t SPEC_AVERAGE = .02;
-static const spx_float_t SPEC_AVERAGE_1 = .98;
+static const spx_float_t MAX_ALPHA = .008f;
+static const spx_float_t ALPHA0 = .05f;
+static const spx_float_t MIN_LEAK = .001f;
+static const spx_float_t SPEC_AVERAGE = .02f;
+static const spx_float_t SPEC_AVERAGE_1 = .98f;
 #endif
 
 /** Speex echo cancellation state. */
@@ -101,8 +101,8 @@ struct SpeexEchoState_ {
    spx_word32_t *Xf;
    spx_word32_t *Eh;
    spx_word32_t *Yh;
-   float Pey;
-   float Pyy;
+   spx_float_t Pey;
+   spx_float_t Pyy;
    float *window;
    /*struct drft_lookup *fft_lookup;*/
    void *fft_table;
@@ -248,7 +248,8 @@ SpeexEchoState *speex_echo_state_init(int frame_size, int filter_length)
    st->memX=st->memD=st->memE=0;
    st->preemph = QCONST16(.9,15);
    st->adapted = 0;
-   st->Pey = st->Pyy = 0;
+   st->Pey = FLOAT_ZERO;
+   st->Pyy = FLOAT_ZERO;
    return st;
 }
 
@@ -269,7 +270,7 @@ void speex_echo_state_reset(SpeexEchoState *st)
    
    st->adapted = 0;
    st->sum_adapt = 0;
-   st->Pey = st->Pyy = 0;
+   st->Pey = st->Pyy = FLOAT_ZERO;
 
 }
 
@@ -313,7 +314,7 @@ void speex_echo_cancel(SpeexEchoState *st, short *ref, short *echo, short *out,
    spx_word16_t ss, ss_1;
    float adapt_rate=0;
    spx_float_t Pey = FLOAT_ZERO, Pyy=FLOAT_ZERO;
-   float alpha;
+   spx_float_t alpha;
    float RER;
    
    N = st->window_size;
@@ -419,16 +420,33 @@ void speex_echo_cancel(SpeexEchoState *st, short *ref, short *echo, short *out,
       st->Yh[j] = .98*st->Yh[j] + .02*st->Yf[j];
 #endif
    }
+   
+#ifdef FIXED_POINT
+   spx_word32_t tmp;
+   tmp = MULT16_32_Q15(QCONST16(.05,15),Syy);
+   if (tmp > MULT16_32_Q15(QCONST16(.008,15),See))
+      tmp = MULT16_32_Q15(QCONST16(.008,15),See);
+   alpha = FLOAT_DIV32(tmp, SHR(10000,6)+See);
+   spx_float_t alpha_1 = FLOAT_SUB(FLOAT_ONE, alpha);
+   /*printf ("%f %f\n", REALFLOAT(alpha), REALFLOAT(alpha_1));*/
+   st->Pey = FLOAT_ADD(FLOAT_MULT(alpha_1,st->Pey) , FLOAT_MULT(alpha,Pey));
+   st->Pyy = FLOAT_ADD(FLOAT_MULT(alpha_1,st->Pyy) , FLOAT_MULT(alpha,Pyy));
+   /* We don't really hope to get better than 33 dB attenuation anyway */
+   if (FLOAT_LT(st->Pey, FLOAT_MULT(MIN_LEAK,st->Pyy)))
+      st->Pey = FLOAT_MULT(MIN_LEAK,st->Pyy);
+   leak_estimate = REALFLOAT(st->Pey) / (1+REALFLOAT(st->Pyy));
+#else
    alpha = .05*Syy / (SHR(10000,6)+See);
    if (alpha > .008)
       alpha = .008;
    st->Pey = (1-alpha)*st->Pey + alpha*REALFLOAT(Pey);
    st->Pyy = (1-alpha)*st->Pyy + alpha*REALFLOAT(Pyy);
-   
    /* We don't really hope to get better than 33 dB attenuation anyway */
    if (st->Pey< .001*st->Pyy)
       st->Pey = .001*st->Pyy;
    leak_estimate = st->Pey / (1+st->Pyy);
+#endif
+
    if (leak_estimate > 1)
       leak_estimate = 1;
    /*printf ("%f\n", leak_estimate);*/
index 9ab0b49..e0c5d84 100644 (file)
@@ -45,7 +45,9 @@ typedef struct {
    spx_int16_t e;
 } spx_float_t;
 
-#define FLOAT_ZERO {0,0}
+#define FLOAT_ZERO ((spx_float_t){0,0})
+#define FLOAT_ONE ((spx_float_t){16384,-14})
+#define FLOAT_HALF ((spx_float_t){16384,-15})
 
 #define MIN(a,b) ((a)<(b)?(a):(b))
 static inline spx_float_t PSEUDOFLOAT(float x)
@@ -107,6 +109,46 @@ static inline spx_float_t FLOAT_ADD(spx_float_t a, spx_float_t b)
    return r;
 }
 
+static inline spx_float_t FLOAT_SUB(spx_float_t a, spx_float_t b)
+{
+   if (a.m==0)
+      return b;
+   else if (b.m==0)
+      return a;
+   spx_float_t r = (a).e > (b).e ? (spx_float_t) {((a).m>>1) - ((b).m>>MIN(15,(a).e-(b).e+1)),(a).e+1} : (spx_float_t) {((a).m>>MIN(15,(b).e-(a).e+1)) - ((b).m>>1) ,(b).e+1};
+   if (r.m>0)
+   {
+      if (r.m<16384)
+      {
+         r.m<<=1;
+         r.e-=1;
+      }
+   } else {
+      if (r.m>-16384)
+      {
+         r.m<<=1;
+         r.e-=1;
+      }
+   }
+   /*printf ("%f + %f = %f\n", REALFLOAT(a), REALFLOAT(b), REALFLOAT(r));*/
+   return r;
+}
+
+static inline int FLOAT_LT(spx_float_t a, spx_float_t b)
+{
+   if (a.m==0)
+      return b.m<0;
+   else if (b.m==0)
+      return a.m>0;
+   spx_float_t r = (a).e > (b).e ? (spx_float_t) {((a).m>>1) - ((b).m>>MIN(15,(a).e-(b).e+1)),(a).e+1} : (spx_float_t) {((b).m>>1) - ((a).m>>MIN(15,(b).e-(a).e+1)),(b).e+1};
+   
+   if ((a).e > (b).e)
+      return ((a).m>>1) < ((b).m>>MIN(15,(a).e-(b).e+1));
+   else 
+      return ((b).m>>1) > ((a).m>>MIN(15,(b).e-(a).e+1));
+
+}
+
 static inline spx_float_t FLOAT_MULT(spx_float_t a, spx_float_t b)
 {
    spx_float_t r = (spx_float_t) {(spx_int16_t)((spx_int32_t)(a).m*(b).m>>15), (a).e+(b).e+15};