Line data Source code
1 : /******************************************************************************************************
2 :
3 : (C) 2022-2025 IVAS codec Public Collaboration with portions copyright Dolby International AB, Ericsson AB,
4 : Fraunhofer-Gesellschaft zur Foerderung der angewandten Forschung e.V., Huawei Technologies Co. LTD.,
5 : Koninklijke Philips N.V., Nippon Telegraph and Telephone Corporation, Nokia Technologies Oy, Orange,
6 : Panasonic Holdings Corporation, Qualcomm Technologies, Inc., VoiceAge Corporation, and other
7 : contributors to this repository. All Rights Reserved.
8 :
9 : This software is protected by copyright law and by international treaties.
10 : The IVAS codec Public Collaboration consisting of Dolby International AB, Ericsson AB,
11 : Fraunhofer-Gesellschaft zur Foerderung der angewandten Forschung e.V., Huawei Technologies Co. LTD.,
12 : Koninklijke Philips N.V., Nippon Telegraph and Telephone Corporation, Nokia Technologies Oy, Orange,
13 : Panasonic Holdings Corporation, Qualcomm Technologies, Inc., VoiceAge Corporation, and other
14 : contributors to this repository retain full ownership rights in their respective contributions in
15 : the software. This notice grants no license of any kind, including but not limited to patent
16 : license, nor is any license granted by implication, estoppel or otherwise.
17 :
18 : Contributors are required to enter into the IVAS codec Public Collaboration agreement before making
19 : contributions.
20 :
21 : This software is provided "AS IS", without any express or implied warranties. The software is in the
22 : development stage. It is intended exclusively for experts who have experience with such software and
23 : solely for the purpose of inspection. All implied warranties of non-infringement, merchantability
24 : and fitness for a particular purpose are hereby disclaimed and excluded.
25 :
26 : Any dispute, controversy or claim arising under or in relation to providing this software shall be
27 : submitted to and settled by the final, binding jurisdiction of the courts of Munich, Germany in
28 : accordance with the laws of the Federal Republic of Germany excluding its conflict of law rules and
29 : the United Nations Convention on Contracts on the International Sales of Goods.
30 :
31 : *******************************************************************************************************/
32 :
33 : /*====================================================================================
34 : EVS Codec 3GPP TS26.443 Nov 04, 2021. Version 12.14.0 / 13.10.0 / 14.6.0 / 15.4.0 / 16.3.0
35 : ====================================================================================*/
36 :
37 : /*! @file jbm_pcmdsp_similarityestimation.c Algorithms for correlation and similarity estimation. */
38 :
39 : /* system headers */
40 : #include <stdlib.h>
41 : #include <math.h>
42 : #include <stdint.h>
43 : #include <assert.h>
44 : #include "options.h"
45 : #include "wmc_auto.h"
46 : #include "basop_util.h"
47 :
48 : /* local headers */
49 : #include "jbm_pcmdsp_similarityestimation.h"
50 :
51 :
52 : /* Calculates cross correlation coefficient for template segment. */
53 566 : void scaleSignal16( const Word16 *src, Word16 *dst, Word16 n, Word16 rightShift )
54 : {
55 : Word16 i;
56 :
57 1826166 : FOR( i = 0; i < n; i++ )
58 : {
59 1825600 : dst[i] = shr_r( src[i], rightShift );
60 1825600 : move16();
61 : }
62 566 : }
63 : /* Calculates cross correlation coefficient for template segment. */
64 97200 : Word32 cross_correlation_subsampled_self_fx( const Word16 *signal,
65 : Word16 x,
66 : Word16 y,
67 : Word16 corr_len,
68 : Word16 subsampling )
69 : {
70 : Word32 sum;
71 : Word16 i;
72 :
73 97200 : sum = 0;
74 97200 : move32();
75 7873200 : FOR( i = 0; i < corr_len; i += subsampling )
76 : {
77 7776000 : sum = L_mac0( sum, signal[x + i], signal[y + i] );
78 : }
79 :
80 97200 : return sum;
81 : }
82 :
83 : /* Calculates cross correlation coefficient for template segment. */
84 :
85 3616 : Word16 normalized_cross_correlation_self_fx( const Word16 *signal,
86 : Word16 x,
87 : Word16 y,
88 : Word16 corr_len,
89 : Word16 subsampling,
90 : Word32 *energy )
91 : {
92 : const Word16 *signalX, *signalY;
93 : Word32 sumXY, sumXX, sumYY, product;
94 : Word16 sqrtXY, cc;
95 : Word16 i, normX, normY, normXY, normCC;
96 :
97 3616 : signalX = &signal[x];
98 3616 : signalY = &signal[y];
99 3616 : sumXY = 0;
100 3616 : sumXX = 0;
101 3616 : sumYY = 0;
102 3616 : move32();
103 3616 : move32();
104 3616 : move32();
105 760096 : FOR( i = 0; i < corr_len; i += subsampling )
106 : {
107 756480 : sumXY = L_mac0( sumXY, signalX[i], signalY[i] );
108 756480 : sumXX = L_mac0( sumXX, signalX[i], signalX[i] );
109 756480 : sumYY = L_mac0( sumYY, signalY[i], signalY[i] );
110 : }
111 :
112 3616 : normX = norm_l( sumXX );
113 3616 : sumXX = L_shl( sumXX, normX );
114 3616 : normY = norm_l( sumYY );
115 3616 : sumYY = L_shl( sumYY, normY );
116 3616 : product = L_mult0( extract_h( sumXX ), extract_h( sumYY ) );
117 3616 : normXY = add( normX, normY );
118 3616 : normXY = sub( normXY, 32 );
119 :
120 : /* change norm to factor of 2 */
121 3616 : IF( s_and( normXY, 0x1 ) != 0 )
122 : {
123 1577 : product = L_shr( product, 1 );
124 1577 : normXY = sub( normXY, 1 );
125 : }
126 3616 : sqrtXY = getSqrtWord32( product );
127 3616 : normXY = shr( normXY, 1 );
128 :
129 3616 : IF( sqrtXY != 0 )
130 : {
131 3552 : move32();
132 3552 : normCC = 0;
133 3552 : move16();
134 3552 : cc = BASOP_Util_Divide3216_Scale( sumXY, sqrtXY, &normCC );
135 3552 : normCC = add( normCC, 16 );
136 : /* scale to Q15 with saturation */
137 : BASOP_SATURATE_WARNING_OFF
138 : #ifdef ISSUE_1751_replace_shl_ro
139 3552 : cc = shr_r_sat( cc, negate( add( normXY, normCC ) ) );
140 : #else
141 : Flag Overflow;
142 : cc = shl_ro( cc, add( normXY, normCC ), &Overflow );
143 : #endif
144 : BASOP_SATURATE_WARNING_ON
145 3552 : *energy = L_shr_r( L_deposit_l( sqrtXY ), normXY );
146 : }
147 : ELSE /* conceal silent frames */
148 : {
149 64 : cc = 0;
150 64 : move16();
151 64 : *energy = L_deposit_l( 1 );
152 64 : move32();
153 : }
154 :
155 3616 : return cc; /* Q15 */
156 : }
157 :
158 27 : Word16 getSignalScaleForCorrelation( Word32 sampleRate )
159 : {
160 : Word16 ret;
161 :
162 27 : IF( LT_32( sampleRate, 16000 ) )
163 : {
164 0 : ret = 2;
165 0 : move16();
166 : }
167 27 : ELSE IF( GE_32( sampleRate, 32000 ) )
168 : {
169 24 : ret = 4;
170 24 : move16();
171 : }
172 : ELSE
173 : {
174 3 : ret = 3;
175 3 : move16();
176 : }
177 :
178 27 : return ret;
179 : }
180 :
181 0 : Word32 cross_correlation_self_fx( const Word16 *signal,
182 : Word16 x,
183 : Word16 y,
184 : Word16 corr_len )
185 : {
186 : Word32 sum;
187 : Word16 i;
188 :
189 0 : sum = 0;
190 0 : move32();
191 0 : FOR( i = 0; i < corr_len; i++ )
192 : {
193 0 : sum = L_mac0( sum, signal[x + i], signal[y + i] );
194 : }
195 :
196 0 : return sum;
197 : }
198 :
199 0 : Word8 isSilence_fx( const Word16 *signal, Word16 len, Word16 segments )
200 : {
201 : Word16 i, j, samplesPerSegment;
202 : Word32 energy, maxEnergy;
203 : Word8 ret;
204 :
205 0 : assert( len > 0 );
206 0 : assert( segments > 0 );
207 :
208 : /* Every segment is checked using the following formula:
209 : * 10 * log10(sum_i(signal[i]*signal[i]))) > -65
210 : * For simplification/complexity, this is replaced by:
211 : * 20 * log10(sum_i(abs(signal[i]))) > -65
212 : */
213 :
214 0 : ret = 1;
215 0 : move16();
216 0 : energy = 0;
217 0 : move32();
218 0 : samplesPerSegment = idiv1616U( len, segments );
219 : /* calculate maxEnergy with factor 2 to reduce rounding error */
220 0 : maxEnergy = L_mult0( samplesPerSegment, 37 ); /* 37 = 2 * exp10(-65.0 / 20) * 32768 */
221 0 : maxEnergy = L_shr( maxEnergy, 1 );
222 0 : j = samplesPerSegment;
223 0 : move16();
224 : /* check all but last segment */
225 0 : FOR( i = 0; i < len; i++ )
226 : {
227 : /* division by 32768 is done later */
228 0 : energy = L_add( energy, L_abs( L_deposit_l( signal[i] ) ) );
229 0 : IF( EQ_16( i, j ) )
230 : {
231 : /* check energy of current segment */
232 : /* 20 * log10(energy / 32768 / samplesPerSegment) > -65
233 : * => energy > samplesPerSegment * 10 ^ (-65 / 20) * 32768 */
234 0 : IF( GT_32( energy, maxEnergy ) )
235 : {
236 0 : ret = 0;
237 0 : move16();
238 0 : BREAK;
239 : }
240 0 : energy = 0;
241 0 : move32();
242 0 : j = add( j, samplesPerSegment );
243 : }
244 : }
245 : /* check last segment */
246 0 : if ( GT_32( energy, maxEnergy ) )
247 : {
248 0 : ret = 0;
249 0 : move16();
250 : }
251 0 : return ret;
252 : }
253 :
254 566 : Word8 isSilence_ivas_fx( const Word16 *signal, Word16 q_sig, Word16 len, Word16 segments )
255 : {
256 : Word16 i, j, samplesPerSegment;
257 : Word64 energy, maxEnergy;
258 566 : Word16 sig_shift = shl( q_sig, 1 );
259 :
260 566 : assert( len > 0 );
261 566 : assert( segments > 0 );
262 :
263 : /* Every segment is checked using the following formula:
264 : * 10 * log10(sum_i(signal[i] / 32768.0 *signal[i] / 32768.0)) / samples_per_segment) > -65
265 : * For simplification/complexity, this is replaced by:
266 : * sum_i(signal[i] / 32768.0 * signal[i] / 32768.0) / samples_per_segment > 10^(-65 / 10)
267 : * sum_i(signal[i] * signal[i]) > 10^(-65 / 10) * 32768.0 * 32768.0 * samples_per_segment
268 : */
269 :
270 566 : energy = 0;
271 566 : move64();
272 566 : samplesPerSegment = idiv1616U( len, segments );
273 566 : maxEnergy = W_shl( W_mac0_16_16( 0, samplesPerSegment, 340 ), sig_shift ); /* 340 = 10^(-65.0 / 10) * 32768 * 32768.0 */
274 566 : j = samplesPerSegment;
275 566 : move16();
276 : /* check all but last segment */
277 1145 : FOR( i = 0; i < len; i++ )
278 : {
279 : /* division by 32768 is done later */
280 1145 : energy = W_mac0_16_16( energy, signal[i], signal[i] );
281 : /* check energy of current segment
282 : * => energy > samplesPerSegment * 10 ^ (-65 / 10) * 32768 * 32768.0*/
283 1145 : IF( GT_64( energy, maxEnergy ) )
284 : {
285 566 : return 0;
286 : }
287 579 : IF( EQ_16( i, j ) )
288 : {
289 3 : energy = 0;
290 3 : move64();
291 3 : j = add( j, samplesPerSegment );
292 : }
293 : }
294 :
295 0 : return 1;
296 : }
|