Remove trailing whitespace.
[opus.git] / silk / float / burg_modified_FLP.c
index b5862a5..2f6edfa 100644 (file)
@@ -30,29 +30,29 @@ OF THIS SOFTWARE, EVEN IF ADVISED OF THE POSSIBILITY OF SUCH DAMAGE.
 #endif
 
 #include "SigProc_FLP.h"
+#include "tuning_parameters.h"
+#include "define.h"
 
 #define MAX_FRAME_SIZE              384 /* subfr_length * nb_subfr = ( 0.005 * 16000 + 16 ) * 4 = 384*/
-#define MAX_NB_SUBFR                4
 
 /* Compute reflection coefficients from input signal */
 silk_float silk_burg_modified_FLP(          /* O    returns residual energy                                     */
     silk_float          A[],                /* O    prediction coefficients (length order)                      */
     const silk_float    x[],                /* I    input signal, length: nb_subfr*(D+L_sub)                    */
+    const silk_float    minInvGain,         /* I    minimum inverse prediction gain                             */
     const opus_int      subfr_length,       /* I    input signal subframe length (incl. D preceeding samples)   */
     const opus_int      nb_subfr,           /* I    number of subframes stacked in x                            */
-    const silk_float    WhiteNoiseFrac,     /* I    fraction added to zero-lag autocorrelation                  */
     const opus_int      D                   /* I    order                                                       */
 )
 {
-    opus_int         k, n, s;
-    double          C0, num, nrg_f, nrg_b, rc, Atmp, tmp1, tmp2;
+    opus_int         k, n, s, reached_max_gain;
+    double           C0, invGain, num, nrg_f, nrg_b, rc, Atmp, tmp1, tmp2;
     const silk_float *x_ptr;
-    double          C_first_row[ SILK_MAX_ORDER_LPC ], C_last_row[ SILK_MAX_ORDER_LPC ];
-    double          CAf[ SILK_MAX_ORDER_LPC + 1 ], CAb[ SILK_MAX_ORDER_LPC + 1 ];
-    double          Af[ SILK_MAX_ORDER_LPC ];
+    double           C_first_row[ SILK_MAX_ORDER_LPC ], C_last_row[ SILK_MAX_ORDER_LPC ];
+    double           CAf[ SILK_MAX_ORDER_LPC + 1 ], CAb[ SILK_MAX_ORDER_LPC + 1 ];
+    double           Af[ SILK_MAX_ORDER_LPC ];
 
     silk_assert( subfr_length * nb_subfr <= MAX_FRAME_SIZE );
-    silk_assert( nb_subfr <= MAX_NB_SUBFR );
 
     /* Compute autocorrelations, added over subframes */
     C0 = silk_energy_FLP( x, nb_subfr * subfr_length );
@@ -66,8 +66,9 @@ silk_float silk_burg_modified_FLP(          /* O    returns residual energy
     silk_memcpy( C_last_row, C_first_row, SILK_MAX_ORDER_LPC * sizeof( double ) );
 
     /* Initialize */
-    CAb[ 0 ] = CAf[ 0 ] = C0 + WhiteNoiseFrac * C0 + 1e-9f;
-
+    CAb[ 0 ] = CAf[ 0 ] = C0 + FIND_LPC_COND_FAC * C0 + 1e-9f;
+    invGain = 1.0f;
+    reached_max_gain = 0;
     for( n = 0; n < D; n++ ) {
         /* Update first row of correlation matrix (without first element) */
         /* Update last row of correlation matrix (without last element, stored in reversed order) */
@@ -93,7 +94,7 @@ silk_float silk_burg_modified_FLP(          /* O    returns residual energy
         tmp2 = C_last_row[ n ];
         for( k = 0; k < n; k++ ) {
             Atmp = Af[ k ];
-            tmp1 += C_last_row[ n - k - 1 ]  * Atmp;
+            tmp1 += C_last_row[  n - k - 1 ] * Atmp;
             tmp2 += C_first_row[ n - k - 1 ] * Atmp;
         }
         CAf[ n + 1 ] = tmp1;
@@ -116,6 +117,21 @@ silk_float silk_burg_modified_FLP(          /* O    returns residual energy
         rc = -2.0 * num / ( nrg_f + nrg_b );
         silk_assert( rc > -1.0 && rc < 1.0 );
 
+        /* Update inverse prediction gain */
+        tmp1 = invGain * ( 1.0 - rc * rc );
+        if( tmp1 <= minInvGain ) {
+            /* Max prediction gain exceeded; set reflection coefficient such that max prediction gain is exactly hit */
+            rc = sqrt( 1.0 - minInvGain / invGain );
+            if( num > 0 ) {
+                /* Ensure adjusted reflection coefficients has the original sign */
+                rc = -rc;
+            }
+            invGain = minInvGain;
+            reached_max_gain = 1;
+        } else {
+            invGain = tmp1;
+        }
+
         /* Update the AR coefficients */
         for( k = 0; k < (n + 1) >> 1; k++ ) {
             tmp1 = Af[ k ];
@@ -125,6 +141,14 @@ silk_float silk_burg_modified_FLP(          /* O    returns residual energy
         }
         Af[ n ] = rc;
 
+        if( reached_max_gain ) {
+            /* Reached max prediction gain; set remaining coefficients to zero and exit loop */
+            for( k = n + 1; k < D; k++ ) {
+                Af[ k ] = 0.0;
+            }
+            break;
+        }
+
         /* Update C * Af and C * Ab */
         for( k = 0; k <= n + 1; k++ ) {
             tmp1 = CAf[ k ];
@@ -133,16 +157,30 @@ silk_float silk_burg_modified_FLP(          /* O    returns residual energy
         }
     }
 
-    /* Return residual energy */
-    nrg_f = CAf[ 0 ];
-    tmp1 = 1.0;
-    for( k = 0; k < D; k++ ) {
-        Atmp = Af[ k ];
-        nrg_f += CAf[ k + 1 ] * Atmp;
-        tmp1  += Atmp * Atmp;
-        A[ k ] = (silk_float)(-Atmp);
+    if( reached_max_gain ) {
+        /* Convert to silk_float */
+        for( k = 0; k < D; k++ ) {
+            A[ k ] = (silk_float)( -Af[ k ] );
+        }
+        /* Subtract energy of preceeding samples from C0 */
+        for( s = 0; s < nb_subfr; s++ ) {
+            C0 -= silk_energy_FLP( x + s * subfr_length, D );
+        }
+        /* Approximate residual energy */
+        nrg_f = C0 * invGain;
+    } else {
+        /* Compute residual energy and store coefficients as silk_float */
+        nrg_f = CAf[ 0 ];
+        tmp1 = 1.0;
+        for( k = 0; k < D; k++ ) {
+            Atmp = Af[ k ];
+            nrg_f += CAf[ k + 1 ] * Atmp;
+            tmp1  += Atmp * Atmp;
+            A[ k ] = (silk_float)(-Atmp);
+        }
+        nrg_f -= FIND_LPC_COND_FAC * C0 * tmp1;
     }
-    nrg_f -= WhiteNoiseFrac * C0 * tmp1;
 
+    /* Return residual energy */
     return (silk_float)nrg_f;
 }