FreeCalypso > hg > gsm-codec-lib

comparison libtwamr/pitch_fr.c @ 391:be8edf9e6bc1

Find changesets by keywords (author, files, the commit message), revision number or hash, or revset expression.
libtwamr: integrate pitch_fr.c
author Mychaela Falconia <falcon@freecalypso.org>
date Mon, 06 May 2024 18:30:00 +0000
parents
children
comparison
equal deleted inserted replaced
390:bde9f5804670 391:be8edf9e6bc1
1 /*
2 ********************************************************************************
3 *
4 * GSM AMR-NB speech codec R98 Version 7.6.0 December 12, 2001
5 * R99 Version 3.3.0
6 * REL-4 Version 4.1.0
7 *
8 ********************************************************************************
9 *
10 * File : pitch_fr.c
11 * Purpose : Find the pitch period with 1/3 or 1/6 subsample
12 * : resolution (closed loop).
13 *
14 ********************************************************************************
15 */
16 /*
17 ********************************************************************************
18 * MODULE INCLUDE FILE AND VERSION ID
19 ********************************************************************************
20 */
21 #include "namespace.h"
22 #include "pitch_fr.h"
23 /*
24 ********************************************************************************
25 * INCLUDE FILES
26 ********************************************************************************
27 */
28 #include "typedef.h"
29 #include "basic_op.h"
30 #include "oper_32b.h"
31 #include "no_count.h"
32 #include "cnst.h"
33 #include "enc_lag3.h"
34 #include "enc_lag6.h"
35 #include "inter_36.h"
36 #include "inv_sqrt.h"
37 #include "convolve.h"
38
39 /*
40 ********************************************************************************
41 * LOCAL VARIABLES AND TABLES
42 ********************************************************************************
43 */
44
45 /*
46 * mode dependent parameters used in Pitch_fr()
47 * Note: order of MRxx in 'enum Mode' is important!
48 */
49 static const struct {
50 Word16 max_frac_lag; /* lag up to which fractional lags are used */
51 Word16 flag3; /* enable 1/3 instead of 1/6 fract. resolution */
52 Word16 first_frac; /* first fractional to check */
53 Word16 last_frac; /* last fractional to check */
54 Word16 delta_int_low; /* integer lag below TO to start search from */
55 Word16 delta_int_range; /* integer range around T0 */
56 Word16 delta_frc_low; /* fractional below T0 */
57 Word16 delta_frc_range; /* fractional range around T0 */
58 Word16 pit_min; /* minimum pitch */
59 } mode_dep_parm[8] = {
60 /* MR475 */ { 84, 1, -2, 2, 5, 10, 5, 9, PIT_MIN },
61 /* MR515 */ { 84, 1, -2, 2, 5, 10, 5, 9, PIT_MIN },
62 /* MR59 */ { 84, 1, -2, 2, 3, 6, 5, 9, PIT_MIN },
63 /* MR67 */ { 84, 1, -2, 2, 3, 6, 5, 9, PIT_MIN },
64 /* MR74 */ { 84, 1, -2, 2, 3, 6, 5, 9, PIT_MIN },
65 /* MR795 */ { 84, 1, -2, 2, 3, 6, 10, 19, PIT_MIN },
66 /* MR102 */ { 84, 1, -2, 2, 3, 6, 5, 9, PIT_MIN },
67 /* MR122 */ { 94, 0, -3, 3, 3, 6, 5, 9, PIT_MIN_MR122 }
68 };
69
70 /*
71 ********************************************************************************
72 * LOCAL PROGRAM CODE
73 ********************************************************************************
74 */
75 /*************************************************************************
76 *
77 * FUNCTION: Norm_Corr()
78 *
79 * PURPOSE: Find the normalized correlation between the target vector
80 * and the filtered past excitation.
81 *
82 * DESCRIPTION:
83 * The normalized correlation is given by the correlation between the
84 * target and filtered past excitation divided by the square root of
85 * the energy of filtered excitation.
86 * corr[k] = <x[], y_k[]>/sqrt(y_k[],y_k[])
87 * where x[] is the target vector and y_k[] is the filtered past
88 * excitation at delay k.
89 *
90 *************************************************************************/
91 static void Norm_Corr (Word16 exc[], Word16 xn[], Word16 h[], Word16 L_subfr,
92 Word16 t_min, Word16 t_max, Word16 corr_norm[])
93 {
94 Word16 i, j, k;
95 Word16 corr_h, corr_l, norm_h, norm_l;
96 Word32 s;
97
98 /* Usally dynamic allocation of (L_subfr) */
99 Word16 excf[L_SUBFR];
100 Word16 scaling, h_fac, *s_excf, scaled_excf[L_SUBFR];
101
102 k = -t_min; move16 ();
103
104 /* compute the filtered excitation for the first delay t_min */
105
106 Convolve (&exc[k], h, excf, L_subfr);
107
108 /* scale "excf[]" to avoid overflow */
109
110 for (j = 0; j < L_subfr; j++) {
111 scaled_excf[j] = shr (excf[j], 2); move16 ();
112 }
113
114 /* Compute 1/sqrt(energy of excf[]) */
115
116 s = 0; move32 ();
117 for (j = 0; j < L_subfr; j++) {
118 s = L_mac (s, excf[j], excf[j]);
119 }
120 test ();
121 if (L_sub (s, 67108864L) <= 0) { /* if (s <= 2^26) */
122 s_excf = excf; move16 ();
123 h_fac = 15 - 12; move16 ();
124 scaling = 0; move16 ();
125 }
126 else {
127 /* "excf[]" is divided by 2 */
128 s_excf = scaled_excf; move16 ();
129 h_fac = 15 - 12 - 2; move16 ();
130 scaling = 2; move16 ();
131 }
132
133 /* loop for every possible period */
134
135 for (i = t_min; i <= t_max; i++) {
136 /* Compute 1/sqrt(energy of excf[]) */
137
138 s = 0; move32 ();
139 for (j = 0; j < L_subfr; j++) {
140 s = L_mac (s, s_excf[j], s_excf[j]);
141 }
142
143 s = Inv_sqrt (s);
144 L_Extract (s, &norm_h, &norm_l);
145
146 /* Compute correlation between xn[] and excf[] */
147
148 s = 0; move32 ();
149 for (j = 0; j < L_subfr; j++) {
150 s = L_mac (s, xn[j], s_excf[j]);
151 }
152 L_Extract (s, &corr_h, &corr_l);
153
154 /* Normalize correlation = correlation * (1/sqrt(energy)) */
155
156 s = Mpy_32 (corr_h, corr_l, norm_h, norm_l);
157
158 corr_norm[i] = extract_h (L_shl (s, 16));
159 move16 ();
160
161 /* modify the filtered excitation excf[] for the next iteration */
162
163 test ();
164 if (sub (i, t_max) != 0) {
165 k--;
166 for (j = L_subfr - 1; j > 0; j--) {
167 s = L_mult (exc[k], h[j]);
168 s = L_shl (s, h_fac);
169 s_excf[j] = add (extract_h (s), s_excf[j - 1]); move16 ();
170 }
171 s_excf[0] = shr (exc[k], scaling); move16 ();
172 }
173 }
174 return;
175 }
176
177 /*************************************************************************
178 *
179 * FUNCTION: searchFrac()
180 *
181 * PURPOSE: Find fractional pitch
182 *
183 * DESCRIPTION:
184 * The function interpolates the normalized correlation at the
185 * fractional positions around lag T0. The position at which the
186 * interpolation function reaches its maximum is the fractional pitch.
187 * Starting point of the search is frac, end point is last_frac.
188 * frac is overwritten with the fractional pitch.
189 *
190 *************************************************************************/
191 static void searchFrac (
192 Word16 *lag, /* i/o : integer pitch */
193 Word16 *frac, /* i/o : start point of search -
194 fractional pitch */
195 Word16 last_frac, /* i : endpoint of search */
196 Word16 corr[], /* i : normalized correlation */
197 Word16 flag3 /* i : subsample resolution
198 (3: =1 / 6: =0) */
199 )
200 {
201 Word16 i;
202 Word16 max;
203 Word16 corr_int;
204
205 /* Test the fractions around T0 and choose the one which maximizes */
206 /* the interpolated normalized correlation. */
207
208 max = Interpol_3or6 (&corr[*lag], *frac, flag3); move16 (); /* function result */
209
210 for (i = add (*frac, 1); i <= last_frac; i++) {
211 corr_int = Interpol_3or6 (&corr[*lag], i, flag3);
212 move16 ();
213 test ();
214 if (sub (corr_int, max) > 0) {
215 max = corr_int; move16 ();
216 *frac = i; move16 ();
217 }
218 }
219
220 test();
221 if (flag3 == 0) {
222 /* Limit the fraction value in the interval [-2,-1,0,1,2,3] */
223
224 test ();
225 if (sub (*frac, -3) == 0) {
226 *frac = 3; move16 ();
227 *lag = sub (*lag, 1);
228 }
229 }
230 else {
231 /* limit the fraction value between -1 and 1 */
232
233 test ();
234 if (sub (*frac, -2) == 0) {
235 *frac = 1; move16 ();
236 *lag = sub (*lag, 1);
237 }
238 test ();
239 if (sub (*frac, 2) == 0) {
240 *frac = -1; move16 ();
241 *lag = add (*lag, 1);
242 }
243 }
244 }
245
246 /*************************************************************************
247 *
248 * FUNCTION: getRange()
249 *
250 * PURPOSE: Sets range around open-loop pitch or integer pitch of last subframe
251 *
252 * DESCRIPTION:
253 * Takes integer pitch T0 and calculates a range around it with
254 * t0_min = T0-delta_low and t0_max = (T0-delta_low) + delta_range
255 * t0_min and t0_max are bounded by pitmin and pitmax
256 *
257 *************************************************************************/
258 static void getRange (
259 Word16 T0, /* i : integer pitch */
260 Word16 delta_low, /* i : search start offset */
261 Word16 delta_range, /* i : search range */
262 Word16 pitmin, /* i : minimum pitch */
263 Word16 pitmax, /* i : maximum pitch */
264 Word16 *t0_min, /* o : search range minimum */
265 Word16 *t0_max) /* o : search range maximum */
266 {
267 *t0_min = sub(T0, delta_low);
268 test ();
269 if (sub(*t0_min, pitmin) < 0) {
270 *t0_min = pitmin; move16();
271 }
272 *t0_max = add(*t0_min, delta_range);
273 test ();
274 if (sub(*t0_max, pitmax) > 0) {
275 *t0_max = pitmax; move16();
276 *t0_min = sub(*t0_max, delta_range);
277 }
278 }
279
280 /*
281 ********************************************************************************
282 * PUBLIC PROGRAM CODE
283 ********************************************************************************
284 */
285
286 /*************************************************************************
287 *
288 * Function: Pitch_fr_reset
289 * Purpose: Initializes state memory to zero
290 *
291 **************************************************************************
292 */
293 void Pitch_fr_reset (Pitch_frState *state)
294 {
295 state->T0_prev_subframe = 0;
296 }
297
298 /*************************************************************************
299 *
300 * FUNCTION: Pitch_fr()
301 *
302 * PURPOSE: Find the pitch period with 1/3 or 1/6 subsample resolution
303 * (closed loop).
304 *
305 * DESCRIPTION:
306 * - find the normalized correlation between the target and filtered
307 * past excitation in the search range.
308 * - select the delay with maximum normalized correlation.
309 * - interpolate the normalized correlation at fractions -3/6 to 3/6
310 * with step 1/6 around the chosen delay.
311 * - The fraction which gives the maximum interpolated value is chosen.
312 *
313 *************************************************************************/
314 Word16 Pitch_fr ( /* o : pitch period (integer) */
315 Pitch_frState *st, /* i/o : State struct */
316 enum Mode mode, /* i : codec mode */
317 Word16 T_op[], /* i : open loop pitch lags */
318 Word16 exc[], /* i : excitation buffer Q0 */
319 Word16 xn[], /* i : target vector Q0 */
320 Word16 h[], /* i : impulse response of synthesis and
321 weighting filters Q12 */
322 Word16 L_subfr, /* i : Length of subframe */
323 Word16 i_subfr, /* i : subframe offset */
324 Word16 *pit_frac, /* o : pitch period (fractional) */
325 Word16 *resu3, /* o : subsample resolution 1/3 (=1) or 1/6 (=0) */
326 Word16 *ana_index /* o : index of encoding */
327 )
328 {
329 Word16 i;
330 Word16 t_min, t_max;
331 Word16 t0_min, t0_max;
332 Word16 max, lag, frac;
333 Word16 tmp_lag;
334 Word16 *corr;
335 Word16 corr_v[40]; /* Total length = t0_max-t0_min+1+2*L_INTER_SRCH */
336
337 Word16 max_frac_lag;
338 Word16 flag3, flag4;
339 Word16 last_frac;
340 Word16 delta_int_low, delta_int_range;
341 Word16 delta_frc_low, delta_frc_range;
342 Word16 pit_min;
343 Word16 frame_offset;
344 Word16 delta_search;
345
346 /*-----------------------------------------------------------------------*
347 * set mode specific variables *
348 *-----------------------------------------------------------------------*/
349
350 max_frac_lag = mode_dep_parm[mode].max_frac_lag; move16 ();
351 flag3 = mode_dep_parm[mode].flag3; move16 ();
352 frac = mode_dep_parm[mode].first_frac; move16 ();
353 last_frac = mode_dep_parm[mode].last_frac; move16 ();
354 delta_int_low = mode_dep_parm[mode].delta_int_low; move16 ();
355 delta_int_range = mode_dep_parm[mode].delta_int_range; move16 ();
356
357 delta_frc_low = mode_dep_parm[mode].delta_frc_low; move16 ();
358 delta_frc_range = mode_dep_parm[mode].delta_frc_range; move16 ();
359 pit_min = mode_dep_parm[mode].pit_min; move16 ();
360
361 /*-----------------------------------------------------------------------*
362 * decide upon full or differential search *
363 *-----------------------------------------------------------------------*/
364
365 delta_search = 1; move16 ();
366
367 test (); test ();
368 if ((i_subfr == 0) || (sub(i_subfr,L_FRAME_BY2) == 0)) {
369
370 /* Subframe 1 and 3 */
371
372 test (); test (); test ();
373 if (((sub(mode, MR475) != 0) && (sub(mode, MR515) != 0)) ||
374 (sub(i_subfr,L_FRAME_BY2) != 0)) {
375
376 /* set t0_min, t0_max for full search */
377 /* this is *not* done for mode MR475, MR515 in subframe 3 */
378
379 delta_search = 0; /* no differential search */ move16 ();
380
381 /* calculate index into T_op which contains the open-loop */
382 /* pitch estimations for the 2 big subframes */
383
384 frame_offset = 1; move16 ();
385 test ();
386 if (i_subfr == 0)
387 frame_offset = 0; move16 ();
388
389 /* get T_op from the corresponding half frame and */
390 /* set t0_min, t0_max */
391
392 getRange (T_op[frame_offset], delta_int_low, delta_int_range,
393 pit_min, PIT_MAX, &t0_min, &t0_max);
394 }
395 else {
396
397 /* mode MR475, MR515 and 3. Subframe: delta search as well */
398 getRange (st->T0_prev_subframe, delta_frc_low, delta_frc_range,
399 pit_min, PIT_MAX, &t0_min, &t0_max);
400 }
401 }
402 else {
403
404 /* for Subframe 2 and 4 */
405 /* get range around T0 of previous subframe for delta search */
406
407 getRange (st->T0_prev_subframe, delta_frc_low, delta_frc_range,
408 pit_min, PIT_MAX, &t0_min, &t0_max);
409 }
410
411 /*-----------------------------------------------------------------------*
412 * Find interval to compute normalized correlation *
413 *-----------------------------------------------------------------------*/
414
415 t_min = sub (t0_min, L_INTER_SRCH);
416 t_max = add (t0_max, L_INTER_SRCH);
417
418 corr = &corr_v[-t_min]; move16 ();
419
420 /*-----------------------------------------------------------------------*
421 * Compute normalized correlation between target and filtered excitation *
422 *-----------------------------------------------------------------------*/
423
424 Norm_Corr (exc, xn, h, L_subfr, t_min, t_max, corr);
425
426 /*-----------------------------------------------------------------------*
427 * Find integer pitch *
428 *-----------------------------------------------------------------------*/
429
430 max = corr[t0_min]; move16 ();
431 lag = t0_min; move16 ();
432
433 for (i = t0_min + 1; i <= t0_max; i++) {
434 test ();
435 if (sub (corr[i], max) >= 0) {
436 max = corr[i]; move16 ();
437 lag = i; move16 ();
438 }
439 }
440
441 /*-----------------------------------------------------------------------*
442 * Find fractional pitch *
443 *-----------------------------------------------------------------------*/
444 test (); test ();
445 if ((delta_search == 0) && (sub (lag, max_frac_lag) > 0)) {
446
447 /* full search and integer pitch greater than max_frac_lag */
448 /* fractional search is not needed, set fractional to zero */
449
450 frac = 0; move16 ();
451 }
452 else {
453
454 /* if differential search AND mode MR475 OR MR515 OR MR59 OR MR67 */
455 /* then search fractional with 4 bits resolution */
456
457 test (); test (); test (); test (); test ();
458 if ((delta_search != 0) &&
459 ((sub (mode, MR475) == 0) ||
460 (sub (mode, MR515) == 0) ||
461 (sub (mode, MR59) == 0) ||
462 (sub (mode, MR67) == 0))) {
463
464 /* modify frac or last_frac according to position of last */
465 /* integer pitch: either search around integer pitch, */
466 /* or only on left or right side */
467
468 tmp_lag = st->T0_prev_subframe; move16 ();
469 test ();
470 if ( sub( sub(tmp_lag, t0_min), 5) > 0)
471 tmp_lag = add (t0_min, 5);
472 test ();
473 if ( sub( sub(t0_max, tmp_lag), 4) > 0)
474 tmp_lag = sub (t0_max, 4);
475
476 test (); test ();
477 if ((sub (lag, tmp_lag) == 0) ||
478 (sub (lag, sub(tmp_lag, 1)) == 0)) {
479
480 /* normal search in fractions around T0 */
481
482 searchFrac (&lag, &frac, last_frac, corr, flag3);
483
484 }
485 else if (sub (lag, sub (tmp_lag, 2)) == 0) {
486 test ();
487 /* limit search around T0 to the right side */
488 frac = 0; move16 ();
489 searchFrac (&lag, &frac, last_frac, corr, flag3);
490 }
491 else if (sub (lag, add(tmp_lag, 1)) == 0) {
492 test (); test ();
493 /* limit search around T0 to the left side */
494 last_frac = 0; move16 ();
495 searchFrac (&lag, &frac, last_frac, corr, flag3);
496 }
497 else {
498 test (); test ();
499 /* no fractional search */
500 frac = 0; move16 ();
501 }
502 }
503 else
504 /* test the fractions around T0 */
505 searchFrac (&lag, &frac, last_frac, corr, flag3);
506 }
507
508 /*-----------------------------------------------------------------------*
509 * encode pitch *
510 *-----------------------------------------------------------------------*/
511
512 test ();
513 if (flag3 != 0) {
514 /* flag4 indicates encoding with 4 bit resolution; */
515 /* this is needed for mode MR475, MR515 and MR59 */
516
517 flag4 = 0; move16 ();
518 test (); test (); test (); test ();
519 if ( (sub (mode, MR475) == 0) ||
520 (sub (mode, MR515) == 0) ||
521 (sub (mode, MR59) == 0) ||
522 (sub (mode, MR67) == 0) ) {
523 flag4 = 1; move16 ();
524 }
525
526 /* encode with 1/3 subsample resolution */
527
528 *ana_index = Enc_lag3(lag, frac, st->T0_prev_subframe,
529 t0_min, t0_max, delta_search, flag4); move16 (); /* function result */
530
531 }
532 else
533 {
534 /* encode with 1/6 subsample resolution */
535
536 *ana_index = Enc_lag6(lag, frac, t0_min, delta_search); move16 (); /* function result */
537 }
538
539 /*-----------------------------------------------------------------------*
540 * update state variables *
541 *-----------------------------------------------------------------------*/
542
543 st->T0_prev_subframe = lag; move16 ();
544
545 /*-----------------------------------------------------------------------*
546 * update output variables *
547 *-----------------------------------------------------------------------*/
548
549 *resu3 = flag3; move16 ();
550
551 *pit_frac = frac; move16 ();
552
553 return (lag);
554 }