media/libopus/silk/float/find_LTP_FLP.c

Thu, 22 Jan 2015 13:21:57 +0100

author
Michael Schloh von Bennewitz <michael@schloh.com>
date
Thu, 22 Jan 2015 13:21:57 +0100
branch
TOR_BUG_9701
changeset 15
b8a032363ba2
permissions
-rw-r--r--

Incorporate requested changes from Mozilla in review:
https://bugzilla.mozilla.org/show_bug.cgi?id=1123480#c6

michael@0 1 /***********************************************************************
michael@0 2 Copyright (c) 2006-2011, Skype Limited. All rights reserved.
michael@0 3 Redistribution and use in source and binary forms, with or without
michael@0 4 modification, are permitted provided that the following conditions
michael@0 5 are met:
michael@0 6 - Redistributions of source code must retain the above copyright notice,
michael@0 7 this list of conditions and the following disclaimer.
michael@0 8 - Redistributions in binary form must reproduce the above copyright
michael@0 9 notice, this list of conditions and the following disclaimer in the
michael@0 10 documentation and/or other materials provided with the distribution.
michael@0 11 - Neither the name of Internet Society, IETF or IETF Trust, nor the
michael@0 12 names of specific contributors, may be used to endorse or promote
michael@0 13 products derived from this software without specific prior written
michael@0 14 permission.
michael@0 15 THIS SOFTWARE IS PROVIDED BY THE COPYRIGHT HOLDERS AND CONTRIBUTORS "AS IS"
michael@0 16 AND ANY EXPRESS OR IMPLIED WARRANTIES, INCLUDING, BUT NOT LIMITED TO, THE
michael@0 17 IMPLIED WARRANTIES OF MERCHANTABILITY AND FITNESS FOR A PARTICULAR PURPOSE
michael@0 18 ARE DISCLAIMED. IN NO EVENT SHALL THE COPYRIGHT OWNER OR CONTRIBUTORS BE
michael@0 19 LIABLE FOR ANY DIRECT, INDIRECT, INCIDENTAL, SPECIAL, EXEMPLARY, OR
michael@0 20 CONSEQUENTIAL DAMAGES (INCLUDING, BUT NOT LIMITED TO, PROCUREMENT OF
michael@0 21 SUBSTITUTE GOODS OR SERVICES; LOSS OF USE, DATA, OR PROFITS; OR BUSINESS
michael@0 22 INTERRUPTION) HOWEVER CAUSED AND ON ANY THEORY OF LIABILITY, WHETHER IN
michael@0 23 CONTRACT, STRICT LIABILITY, OR TORT (INCLUDING NEGLIGENCE OR OTHERWISE)
michael@0 24 ARISING IN ANY WAY OUT OF THE USE OF THIS SOFTWARE, EVEN IF ADVISED OF THE
michael@0 25 POSSIBILITY OF SUCH DAMAGE.
michael@0 26 ***********************************************************************/
michael@0 27
michael@0 28 #ifdef HAVE_CONFIG_H
michael@0 29 #include "config.h"
michael@0 30 #endif
michael@0 31
michael@0 32 #include "main_FLP.h"
michael@0 33 #include "tuning_parameters.h"
michael@0 34
michael@0 35 void silk_find_LTP_FLP(
michael@0 36 silk_float b[ MAX_NB_SUBFR * LTP_ORDER ], /* O LTP coefs */
michael@0 37 silk_float WLTP[ MAX_NB_SUBFR * LTP_ORDER * LTP_ORDER ], /* O Weight for LTP quantization */
michael@0 38 silk_float *LTPredCodGain, /* O LTP coding gain */
michael@0 39 const silk_float r_lpc[], /* I LPC residual */
michael@0 40 const opus_int lag[ MAX_NB_SUBFR ], /* I LTP lags */
michael@0 41 const silk_float Wght[ MAX_NB_SUBFR ], /* I Weights */
michael@0 42 const opus_int subfr_length, /* I Subframe length */
michael@0 43 const opus_int nb_subfr, /* I number of subframes */
michael@0 44 const opus_int mem_offset /* I Number of samples in LTP memory */
michael@0 45 )
michael@0 46 {
michael@0 47 opus_int i, k;
michael@0 48 silk_float *b_ptr, temp, *WLTP_ptr;
michael@0 49 silk_float LPC_res_nrg, LPC_LTP_res_nrg;
michael@0 50 silk_float d[ MAX_NB_SUBFR ], m, g, delta_b[ LTP_ORDER ];
michael@0 51 silk_float w[ MAX_NB_SUBFR ], nrg[ MAX_NB_SUBFR ], regu;
michael@0 52 silk_float Rr[ LTP_ORDER ], rr[ MAX_NB_SUBFR ];
michael@0 53 const silk_float *r_ptr, *lag_ptr;
michael@0 54
michael@0 55 b_ptr = b;
michael@0 56 WLTP_ptr = WLTP;
michael@0 57 r_ptr = &r_lpc[ mem_offset ];
michael@0 58 for( k = 0; k < nb_subfr; k++ ) {
michael@0 59 lag_ptr = r_ptr - ( lag[ k ] + LTP_ORDER / 2 );
michael@0 60
michael@0 61 silk_corrMatrix_FLP( lag_ptr, subfr_length, LTP_ORDER, WLTP_ptr );
michael@0 62 silk_corrVector_FLP( lag_ptr, r_ptr, subfr_length, LTP_ORDER, Rr );
michael@0 63
michael@0 64 rr[ k ] = ( silk_float )silk_energy_FLP( r_ptr, subfr_length );
michael@0 65 regu = 1.0f + rr[ k ] +
michael@0 66 matrix_ptr( WLTP_ptr, 0, 0, LTP_ORDER ) +
michael@0 67 matrix_ptr( WLTP_ptr, LTP_ORDER-1, LTP_ORDER-1, LTP_ORDER );
michael@0 68 regu *= LTP_DAMPING / 3;
michael@0 69 silk_regularize_correlations_FLP( WLTP_ptr, &rr[ k ], regu, LTP_ORDER );
michael@0 70 silk_solve_LDL_FLP( WLTP_ptr, LTP_ORDER, Rr, b_ptr );
michael@0 71
michael@0 72 /* Calculate residual energy */
michael@0 73 nrg[ k ] = silk_residual_energy_covar_FLP( b_ptr, WLTP_ptr, Rr, rr[ k ], LTP_ORDER );
michael@0 74
michael@0 75 temp = Wght[ k ] / ( nrg[ k ] * Wght[ k ] + 0.01f * subfr_length );
michael@0 76 silk_scale_vector_FLP( WLTP_ptr, temp, LTP_ORDER * LTP_ORDER );
michael@0 77 w[ k ] = matrix_ptr( WLTP_ptr, LTP_ORDER / 2, LTP_ORDER / 2, LTP_ORDER );
michael@0 78
michael@0 79 r_ptr += subfr_length;
michael@0 80 b_ptr += LTP_ORDER;
michael@0 81 WLTP_ptr += LTP_ORDER * LTP_ORDER;
michael@0 82 }
michael@0 83
michael@0 84 /* Compute LTP coding gain */
michael@0 85 if( LTPredCodGain != NULL ) {
michael@0 86 LPC_LTP_res_nrg = 1e-6f;
michael@0 87 LPC_res_nrg = 0.0f;
michael@0 88 for( k = 0; k < nb_subfr; k++ ) {
michael@0 89 LPC_res_nrg += rr[ k ] * Wght[ k ];
michael@0 90 LPC_LTP_res_nrg += nrg[ k ] * Wght[ k ];
michael@0 91 }
michael@0 92
michael@0 93 silk_assert( LPC_LTP_res_nrg > 0 );
michael@0 94 *LTPredCodGain = 3.0f * silk_log2( LPC_res_nrg / LPC_LTP_res_nrg );
michael@0 95 }
michael@0 96
michael@0 97 /* Smoothing */
michael@0 98 /* d = sum( B, 1 ); */
michael@0 99 b_ptr = b;
michael@0 100 for( k = 0; k < nb_subfr; k++ ) {
michael@0 101 d[ k ] = 0;
michael@0 102 for( i = 0; i < LTP_ORDER; i++ ) {
michael@0 103 d[ k ] += b_ptr[ i ];
michael@0 104 }
michael@0 105 b_ptr += LTP_ORDER;
michael@0 106 }
michael@0 107 /* m = ( w * d' ) / ( sum( w ) + 1e-3 ); */
michael@0 108 temp = 1e-3f;
michael@0 109 for( k = 0; k < nb_subfr; k++ ) {
michael@0 110 temp += w[ k ];
michael@0 111 }
michael@0 112 m = 0;
michael@0 113 for( k = 0; k < nb_subfr; k++ ) {
michael@0 114 m += d[ k ] * w[ k ];
michael@0 115 }
michael@0 116 m = m / temp;
michael@0 117
michael@0 118 b_ptr = b;
michael@0 119 for( k = 0; k < nb_subfr; k++ ) {
michael@0 120 g = LTP_SMOOTHING / ( LTP_SMOOTHING + w[ k ] ) * ( m - d[ k ] );
michael@0 121 temp = 0;
michael@0 122 for( i = 0; i < LTP_ORDER; i++ ) {
michael@0 123 delta_b[ i ] = silk_max_float( b_ptr[ i ], 0.1f );
michael@0 124 temp += delta_b[ i ];
michael@0 125 }
michael@0 126 temp = g / temp;
michael@0 127 for( i = 0; i < LTP_ORDER; i++ ) {
michael@0 128 b_ptr[ i ] = b_ptr[ i ] + delta_b[ i ] * temp;
michael@0 129 }
michael@0 130 b_ptr += LTP_ORDER;
michael@0 131 }
michael@0 132 }

mercurial