Coverage Report

Created: 2025-08-26 06:54

/src/speex/libspeex/ltp.c
Line
Count
Source (jump to first uncovered line)
1
/* Copyright (C) 2002-2006 Jean-Marc Valin
2
   File: ltp.c
3
   Long-Term Prediction functions
4
5
   Redistribution and use in source and binary forms, with or without
6
   modification, are permitted provided that the following conditions
7
   are met:
8
9
   - Redistributions of source code must retain the above copyright
10
   notice, this list of conditions and the following disclaimer.
11
12
   - Redistributions in binary form must reproduce the above copyright
13
   notice, this list of conditions and the following disclaimer in the
14
   documentation and/or other materials provided with the distribution.
15
16
   - Neither the name of the Xiph.org Foundation nor the names of its
17
   contributors may be used to endorse or promote products derived from
18
   this software without specific prior written permission.
19
20
   THIS SOFTWARE IS PROVIDED BY THE COPYRIGHT HOLDERS AND CONTRIBUTORS
21
   ``AS IS'' AND ANY EXPRESS OR IMPLIED WARRANTIES, INCLUDING, BUT NOT
22
   LIMITED TO, THE IMPLIED WARRANTIES OF MERCHANTABILITY AND FITNESS FOR
23
   A PARTICULAR PURPOSE ARE DISCLAIMED.  IN NO EVENT SHALL THE FOUNDATION OR
24
   CONTRIBUTORS BE LIABLE FOR ANY DIRECT, INDIRECT, INCIDENTAL, SPECIAL,
25
   EXEMPLARY, OR CONSEQUENTIAL DAMAGES (INCLUDING, BUT NOT LIMITED TO,
26
   PROCUREMENT OF SUBSTITUTE GOODS OR SERVICES; LOSS OF USE, DATA, OR
27
   PROFITS; OR BUSINESS INTERRUPTION) HOWEVER CAUSED AND ON ANY THEORY OF
28
   LIABILITY, WHETHER IN CONTRACT, STRICT LIABILITY, OR TORT (INCLUDING
29
   NEGLIGENCE OR OTHERWISE) ARISING IN ANY WAY OUT OF THE USE OF THIS
30
   SOFTWARE, EVEN IF ADVISED OF THE POSSIBILITY OF SUCH DAMAGE.
31
*/
32
33
#ifdef HAVE_CONFIG_H
34
#include "config.h"
35
#endif
36
37
#include <math.h>
38
#include "ltp.h"
39
#include "stack_alloc.h"
40
#include "filters.h"
41
#include "math_approx.h"
42
#include "os_support.h"
43
44
#ifndef NULL
45
#define NULL 0
46
#endif
47
48
49
#ifdef _USE_SSE
50
#include "ltp_sse.h"
51
#elif defined (ARM4_ASM) || defined(ARM5E_ASM)
52
#include "ltp_arm4.h"
53
#elif defined (BFIN_ASM)
54
#include "ltp_bfin.h"
55
#endif
56
57
#ifndef OVERRIDE_INNER_PROD
58
spx_word32_t inner_prod(const spx_word16_t *x, const spx_word16_t *y, int len)
59
6.39M
{
60
6.39M
   spx_word32_t sum=0;
61
6.39M
   len >>= 2;
62
141M
   while(len--)
63
135M
   {
64
135M
      spx_word32_t part=0;
65
135M
      part = MAC16_16(part,*x++,*y++);
66
135M
      part = MAC16_16(part,*x++,*y++);
67
135M
      part = MAC16_16(part,*x++,*y++);
68
135M
      part = MAC16_16(part,*x++,*y++);
69
      /* HINT: If you had a 40-bit accumulator, you could shift only at the end */
70
135M
      sum = ADD32(sum,SHR32(part,6));
71
135M
   }
72
6.39M
   return sum;
73
6.39M
}
74
#endif
75
76
#ifndef DISABLE_ENCODER
77
78
#ifndef OVERRIDE_PITCH_XCORR
79
#if 0 /* HINT: Enable this for machines with enough registers (i.e. not x86) */
80
static void pitch_xcorr(const spx_word16_t *_x, const spx_word16_t *_y, spx_word32_t *corr, int len, int nb_pitch, char *stack)
81
{
82
   int i,j;
83
   for (i=0;i<nb_pitch;i+=4)
84
   {
85
      /* Compute correlation*/
86
      /*corr[nb_pitch-1-i]=inner_prod(x, _y+i, len);*/
87
      spx_word32_t sum1=0;
88
      spx_word32_t sum2=0;
89
      spx_word32_t sum3=0;
90
      spx_word32_t sum4=0;
91
      const spx_word16_t *y = _y+i;
92
      const spx_word16_t *x = _x;
93
      spx_word16_t y0, y1, y2, y3;
94
      /*y0=y[0];y1=y[1];y2=y[2];y3=y[3];*/
95
      y0=*y++;
96
      y1=*y++;
97
      y2=*y++;
98
      y3=*y++;
99
      for (j=0;j<len;j+=4)
100
      {
101
         spx_word32_t part1;
102
         spx_word32_t part2;
103
         spx_word32_t part3;
104
         spx_word32_t part4;
105
         part1 = MULT16_16(*x,y0);
106
         part2 = MULT16_16(*x,y1);
107
         part3 = MULT16_16(*x,y2);
108
         part4 = MULT16_16(*x,y3);
109
         x++;
110
         y0=*y++;
111
         part1 = MAC16_16(part1,*x,y1);
112
         part2 = MAC16_16(part2,*x,y2);
113
         part3 = MAC16_16(part3,*x,y3);
114
         part4 = MAC16_16(part4,*x,y0);
115
         x++;
116
         y1=*y++;
117
         part1 = MAC16_16(part1,*x,y2);
118
         part2 = MAC16_16(part2,*x,y3);
119
         part3 = MAC16_16(part3,*x,y0);
120
         part4 = MAC16_16(part4,*x,y1);
121
         x++;
122
         y2=*y++;
123
         part1 = MAC16_16(part1,*x,y3);
124
         part2 = MAC16_16(part2,*x,y0);
125
         part3 = MAC16_16(part3,*x,y1);
126
         part4 = MAC16_16(part4,*x,y2);
127
         x++;
128
         y3=*y++;
129
130
         sum1 = ADD32(sum1,SHR32(part1,6));
131
         sum2 = ADD32(sum2,SHR32(part2,6));
132
         sum3 = ADD32(sum3,SHR32(part3,6));
133
         sum4 = ADD32(sum4,SHR32(part4,6));
134
      }
135
      corr[nb_pitch-1-i]=sum1;
136
      corr[nb_pitch-2-i]=sum2;
137
      corr[nb_pitch-3-i]=sum3;
138
      corr[nb_pitch-4-i]=sum4;
139
   }
140
141
}
142
#else
143
static void pitch_xcorr(const spx_word16_t *_x, const spx_word16_t *_y, spx_word32_t *corr, int len, int nb_pitch, char *stack)
144
42.0k
{
145
42.0k
   int i;
146
5.15M
   for (i=0;i<nb_pitch;i++)
147
5.11M
   {
148
      /* Compute correlation*/
149
5.11M
      corr[nb_pitch-1-i]=inner_prod(_x, _y+i, len);
150
5.11M
   }
151
152
42.0k
}
153
#endif
154
#endif
155
156
#ifndef OVERRIDE_COMPUTE_PITCH_ERROR
157
static inline spx_word32_t compute_pitch_error(spx_word16_t *C, spx_word16_t *g, spx_word16_t pitch_control)
158
10.0M
{
159
10.0M
   spx_word32_t sum = 0;
160
10.0M
   sum = ADD32(sum,MULT16_16(MULT16_16_16(g[0],pitch_control),C[0]));
161
10.0M
   sum = ADD32(sum,MULT16_16(MULT16_16_16(g[1],pitch_control),C[1]));
162
10.0M
   sum = ADD32(sum,MULT16_16(MULT16_16_16(g[2],pitch_control),C[2]));
163
10.0M
   sum = SUB32(sum,MULT16_16(MULT16_16_16(g[0],g[1]),C[3]));
164
10.0M
   sum = SUB32(sum,MULT16_16(MULT16_16_16(g[2],g[1]),C[4]));
165
10.0M
   sum = SUB32(sum,MULT16_16(MULT16_16_16(g[2],g[0]),C[5]));
166
10.0M
   sum = SUB32(sum,MULT16_16(MULT16_16_16(g[0],g[0]),C[6]));
167
10.0M
   sum = SUB32(sum,MULT16_16(MULT16_16_16(g[1],g[1]),C[7]));
168
10.0M
   sum = SUB32(sum,MULT16_16(MULT16_16_16(g[2],g[2]),C[8]));
169
10.0M
   return sum;
170
10.0M
}
171
#endif
172
173
#ifndef OVERRIDE_OPEN_LOOP_NBEST_PITCH
174
void open_loop_nbest_pitch(spx_word16_t *sw, int start, int end, int len, int *pitch, spx_word16_t *gain, int N, char *stack)
175
42.0k
{
176
42.0k
   int i,j,k;
177
42.0k
   VARDECL(spx_word32_t *best_score);
178
42.0k
   VARDECL(spx_word32_t *best_ener);
179
42.0k
   spx_word32_t e0;
180
42.0k
   VARDECL(spx_word32_t *corr);
181
42.0k
#ifdef FIXED_POINT
182
   /* In fixed-point, we need only one (temporary) array of 32-bit values and two (corr16, ener16)
183
      arrays for (normalized) 16-bit values */
184
42.0k
   VARDECL(spx_word16_t *corr16);
185
42.0k
   VARDECL(spx_word16_t *ener16);
186
42.0k
   spx_word32_t *energy;
187
42.0k
   int cshift=0, eshift=0;
188
42.0k
   int scaledown = 0;
189
42.0k
   ALLOC(corr16, end-start+1, spx_word16_t);
190
42.0k
   ALLOC(ener16, end-start+1, spx_word16_t);
191
42.0k
   ALLOC(corr, end-start+1, spx_word32_t);
192
42.0k
   energy = corr;
193
#else
194
   /* In floating-point, we need to float arrays and no normalized copies */
195
   VARDECL(spx_word32_t *energy);
196
   spx_word16_t *corr16;
197
   spx_word16_t *ener16;
198
   ALLOC(energy, end-start+2, spx_word32_t);
199
   ALLOC(corr, end-start+1, spx_word32_t);
200
   corr16 = corr;
201
   ener16 = energy;
202
#endif
203
204
42.0k
   ALLOC(best_score, N, spx_word32_t);
205
42.0k
   ALLOC(best_ener, N, spx_word32_t);
206
268k
   for (i=0;i<N;i++)
207
226k
   {
208
226k
        best_score[i]=-1;
209
226k
        best_ener[i]=0;
210
226k
        pitch[i]=start;
211
226k
   }
212
213
42.0k
#ifdef FIXED_POINT
214
7.90M
   for (i=-end;i<len;i++)
215
7.86M
   {
216
7.86M
      if (ABS16(sw[i])>16383)
217
11.3k
      {
218
11.3k
         scaledown=1;
219
11.3k
         break;
220
11.3k
      }
221
7.86M
   }
222
   /* If the weighted input is close to saturation, then we scale it down */
223
42.0k
   if (scaledown)
224
11.3k
   {
225
2.65M
      for (i=-end;i<len;i++)
226
2.64M
      {
227
2.64M
         sw[i]=SHR16(sw[i],1);
228
2.64M
      }
229
11.3k
   }
230
42.0k
#endif
231
42.0k
   energy[0]=inner_prod(sw-start, sw-start, len);
232
42.0k
   e0=inner_prod(sw, sw, len);
233
5.11M
   for (i=start;i<end;i++)
234
5.07M
   {
235
      /* Update energy for next pitch*/
236
5.07M
      energy[i-start+1] = SUB32(ADD32(energy[i-start],SHR32(MULT16_16(sw[-i-1],sw[-i-1]),6)), SHR32(MULT16_16(sw[-i+len-1],sw[-i+len-1]),6));
237
5.07M
      if (energy[i-start+1] < 0)
238
0
         energy[i-start+1] = 0;
239
5.07M
   }
240
241
42.0k
#ifdef FIXED_POINT
242
42.0k
   eshift = normalize16(energy, ener16, 32766, end-start+1);
243
42.0k
#endif
244
245
   /* In fixed-point, this actually overrites the energy array (aliased to corr) */
246
42.0k
   pitch_xcorr(sw, sw-end, corr, len, end-start+1, stack);
247
248
42.0k
#ifdef FIXED_POINT
249
   /* Normalize to 180 so we can square it and it still fits in 16 bits */
250
42.0k
   cshift = normalize16(corr, corr16, 180, end-start+1);
251
   /* If we scaled weighted input down, we need to scale it up again (OK, so we've just lost the LSB, who cares?) */
252
42.0k
   if (scaledown)
253
11.3k
   {
254
2.65M
      for (i=-end;i<len;i++)
255
2.64M
      {
256
2.64M
         sw[i]=SHL16(sw[i],1);
257
2.64M
      }
258
11.3k
   }
259
42.0k
#endif
260
261
   /* Search for the best pitch prediction gain */
262
5.15M
   for (i=start;i<=end;i++)
263
5.11M
   {
264
5.11M
      spx_word16_t tmp = MULT16_16_16(corr16[i-start],corr16[i-start]);
265
      /* Instead of dividing the tmp by the energy, we multiply on the other side */
266
5.11M
      if (MULT16_16(tmp,best_ener[N-1])>MULT16_16(best_score[N-1],ADD16(1,ener16[i-start])))
267
766k
      {
268
         /* We can safely put it last and then check */
269
766k
         best_score[N-1]=tmp;
270
766k
         best_ener[N-1]=ener16[i-start]+1;
271
766k
         pitch[N-1]=i;
272
         /* Check if it comes in front of others */
273
2.45M
         for (j=0;j<N-1;j++)
274
2.32M
         {
275
2.32M
            if (MULT16_16(tmp,best_ener[j])>MULT16_16(best_score[j],ADD16(1,ener16[i-start])))
276
638k
            {
277
2.92M
               for (k=N-1;k>j;k--)
278
2.28M
               {
279
2.28M
                  best_score[k]=best_score[k-1];
280
2.28M
                  best_ener[k]=best_ener[k-1];
281
2.28M
                  pitch[k]=pitch[k-1];
282
2.28M
               }
283
638k
               best_score[j]=tmp;
284
638k
               best_ener[j]=ener16[i-start]+1;
285
638k
               pitch[j]=i;
286
638k
               break;
287
638k
            }
288
2.32M
         }
289
766k
      }
290
5.11M
   }
291
292
   /* Compute open-loop gain if necessary */
293
42.0k
   if (gain)
294
18.3k
   {
295
128k
      for (j=0;j<N;j++)
296
110k
      {
297
110k
         spx_word16_t g;
298
110k
         i=pitch[j];
299
110k
         g = DIV32(SHL32(EXTEND32(corr16[i-start]),cshift), 10+SHR32(MULT16_16(spx_sqrt(e0),spx_sqrt(SHL32(EXTEND32(ener16[i-start]),eshift))),6));
300
         /* FIXME: g = max(g,corr/energy) */
301
110k
         if (g<0)
302
44.9k
            g = 0;
303
110k
         gain[j]=g;
304
110k
      }
305
18.3k
   }
306
307
308
42.0k
}
309
#endif
310
311
#ifndef OVERRIDE_PITCH_GAIN_SEARCH_3TAP_VQ
312
static int pitch_gain_search_3tap_vq(
313
  const signed char *gain_cdbk,
314
  int                gain_cdbk_size,
315
  spx_word16_t      *C16,
316
  spx_word16_t       max_gain
317
)
318
119k
{
319
119k
  const signed char *ptr=gain_cdbk;
320
119k
  int                best_cdbk=0;
321
119k
  spx_word32_t       best_sum=-VERY_LARGE32;
322
119k
  spx_word32_t       sum=0;
323
119k
  spx_word16_t       g[3];
324
119k
  spx_word16_t       pitch_control=64;
325
119k
  spx_word16_t       gain_sum;
326
119k
  int                i;
327
328
10.2M
  for (i=0;i<gain_cdbk_size;i++) {
329
330
10.0M
    ptr = gain_cdbk+4*i;
331
10.0M
    g[0]=ADD16((spx_word16_t)ptr[0],32);
332
10.0M
    g[1]=ADD16((spx_word16_t)ptr[1],32);
333
10.0M
    g[2]=ADD16((spx_word16_t)ptr[2],32);
334
10.0M
    gain_sum = (spx_word16_t)ptr[3];
335
336
10.0M
    sum = compute_pitch_error(C16, g, pitch_control);
337
338
10.0M
    if (sum>best_sum && gain_sum<=max_gain) {
339
410k
      best_sum=sum;
340
410k
      best_cdbk=i;
341
410k
    }
342
10.0M
  }
343
344
119k
  return best_cdbk;
345
119k
}
346
#endif
347
348
/** Finds the best quantized 3-tap pitch predictor by analysis by synthesis */
349
static spx_word32_t pitch_gain_search_3tap(
350
const spx_word16_t target[],       /* Target vector */
351
const spx_coef_t ak[],          /* LPCs for this subframe */
352
const spx_coef_t awk1[],        /* Weighted LPCs #1 for this subframe */
353
const spx_coef_t awk2[],        /* Weighted LPCs #2 for this subframe */
354
spx_sig_t exc[],                /* Excitation */
355
const signed char *gain_cdbk,
356
int gain_cdbk_size,
357
int   pitch,                    /* Pitch value */
358
int   p,                        /* Number of LPC coeffs */
359
int   nsf,                      /* Number of samples in subframe */
360
SpeexBits *bits,
361
char *stack,
362
const spx_word16_t *exc2,
363
const spx_word16_t *r,
364
spx_word16_t *new_target,
365
int  *cdbk_index,
366
int plc_tuning,
367
spx_word32_t cumul_gain,
368
int scaledown
369
)
370
119k
{
371
119k
   int i,j;
372
119k
   VARDECL(spx_word16_t *tmp1);
373
119k
   VARDECL(spx_word16_t *e);
374
119k
   spx_word16_t *x[3];
375
119k
   spx_word32_t corr[3];
376
119k
   spx_word32_t A[3][3];
377
119k
   spx_word16_t gain[3];
378
119k
   spx_word32_t err;
379
119k
   spx_word16_t max_gain=128;
380
119k
   int          best_cdbk=0;
381
382
119k
   ALLOC(tmp1, 3*nsf, spx_word16_t);
383
119k
   ALLOC(e, nsf, spx_word16_t);
384
385
119k
   if (cumul_gain > 262144)
386
895
      max_gain = 31;
387
388
119k
   x[0]=tmp1;
389
119k
   x[1]=tmp1+nsf;
390
119k
   x[2]=tmp1+2*nsf;
391
392
4.91M
   for (j=0;j<nsf;j++)
393
4.79M
      new_target[j] = target[j];
394
395
119k
   {
396
119k
      int bound;
397
119k
      VARDECL(spx_mem_t *mm);
398
119k
      int pp=pitch-1;
399
119k
      ALLOC(mm, p, spx_mem_t);
400
119k
      bound = nsf;
401
119k
      if (nsf-pp>0)
402
45.1k
         bound = pp;
403
4.21M
      for (j=0;j<bound;j++)
404
4.09M
         e[j]=exc2[j-pp];
405
119k
      bound = nsf;
406
119k
      if (nsf-pp-pitch>0)
407
15.6k
         bound = pp+pitch;
408
753k
      for (;j<bound;j++)
409
633k
         e[j]=exc2[j-pp-pitch];
410
186k
      for (;j<nsf;j++)
411
66.5k
         e[j]=0;
412
119k
#ifdef FIXED_POINT
413
      /* Scale target and excitation down if needed (avoiding overflow) */
414
119k
      if (scaledown)
415
27.2k
      {
416
1.11M
         for (j=0;j<nsf;j++)
417
1.09M
            e[j] = SHR16(e[j],1);
418
1.11M
         for (j=0;j<nsf;j++)
419
1.09M
            new_target[j] = SHR16(new_target[j],1);
420
27.2k
      }
421
119k
#endif
422
1.31M
      for (j=0;j<p;j++)
423
1.19M
         mm[j] = 0;
424
119k
      iir_mem16(e, ak, e, nsf, p, mm, stack);
425
1.31M
      for (j=0;j<p;j++)
426
1.19M
         mm[j] = 0;
427
119k
      filter10(e, awk1, awk2, e, nsf, mm, stack);
428
4.91M
      for (j=0;j<nsf;j++)
429
4.79M
         x[2][j] = e[j];
430
119k
   }
431
359k
   for (i=1;i>=0;i--)
432
239k
   {
433
239k
      spx_word16_t e0=exc2[-pitch-1+i];
434
239k
#ifdef FIXED_POINT
435
      /* Scale excitation down if needed (avoiding overflow) */
436
239k
      if (scaledown)
437
54.5k
         e0 = SHR16(e0,1);
438
239k
#endif
439
239k
      x[i][0]=MULT16_16_Q14(r[0], e0);
440
9.59M
      for (j=0;j<nsf-1;j++)
441
9.35M
         x[i][j+1]=ADD32(x[i+1][j],MULT16_16_P14(r[j+1], e0));
442
239k
   }
443
444
479k
   for (i=0;i<3;i++)
445
359k
      corr[i]=inner_prod(x[i],new_target,nsf);
446
479k
   for (i=0;i<3;i++)
447
1.07M
      for (j=0;j<=i;j++)
448
719k
         A[i][j]=A[j][i]=inner_prod(x[i],x[j],nsf);
449
450
119k
   {
451
119k
      spx_word32_t C[9];
452
119k
#ifdef FIXED_POINT
453
119k
      spx_word16_t C16[9];
454
#else
455
      spx_word16_t *C16=C;
456
#endif
457
119k
      C[0]=corr[2];
458
119k
      C[1]=corr[1];
459
119k
      C[2]=corr[0];
460
119k
      C[3]=A[1][2];
461
119k
      C[4]=A[0][1];
462
119k
      C[5]=A[0][2];
463
119k
      C[6]=A[2][2];
464
119k
      C[7]=A[1][1];
465
119k
      C[8]=A[0][0];
466
467
      /*plc_tuning *= 2;*/
468
119k
      if (plc_tuning<2)
469
0
         plc_tuning=2;
470
119k
      if (plc_tuning>30)
471
0
         plc_tuning=30;
472
119k
#ifdef FIXED_POINT
473
119k
      C[0] = SHL32(C[0],1);
474
119k
      C[1] = SHL32(C[1],1);
475
119k
      C[2] = SHL32(C[2],1);
476
119k
      C[3] = SHL32(C[3],1);
477
119k
      C[4] = SHL32(C[4],1);
478
119k
      C[5] = SHL32(C[5],1);
479
119k
      C[6] = MAC16_32_Q15(C[6],MULT16_16_16(plc_tuning,655),C[6]);
480
119k
      C[7] = MAC16_32_Q15(C[7],MULT16_16_16(plc_tuning,655),C[7]);
481
119k
      C[8] = MAC16_32_Q15(C[8],MULT16_16_16(plc_tuning,655),C[8]);
482
119k
      normalize16(C, C16, 32767, 9);
483
#else
484
      C[6]*=.5*(1+.02*plc_tuning);
485
      C[7]*=.5*(1+.02*plc_tuning);
486
      C[8]*=.5*(1+.02*plc_tuning);
487
#endif
488
489
119k
      best_cdbk = pitch_gain_search_3tap_vq(gain_cdbk, gain_cdbk_size, C16, max_gain);
490
491
119k
#ifdef FIXED_POINT
492
119k
      gain[0] = ADD16(32,(spx_word16_t)gain_cdbk[best_cdbk*4]);
493
119k
      gain[1] = ADD16(32,(spx_word16_t)gain_cdbk[best_cdbk*4+1]);
494
119k
      gain[2] = ADD16(32,(spx_word16_t)gain_cdbk[best_cdbk*4+2]);
495
      /*printf ("%d %d %d %d\n",gain[0],gain[1],gain[2], best_cdbk);*/
496
#else
497
      gain[0] = 0.015625*gain_cdbk[best_cdbk*4]  + .5;
498
      gain[1] = 0.015625*gain_cdbk[best_cdbk*4+1]+ .5;
499
      gain[2] = 0.015625*gain_cdbk[best_cdbk*4+2]+ .5;
500
#endif
501
119k
      *cdbk_index=best_cdbk;
502
119k
   }
503
504
119k
   SPEEX_MEMSET(exc, 0, nsf);
505
479k
   for (i=0;i<3;i++)
506
359k
   {
507
359k
      int j;
508
359k
      int tmp1, tmp3;
509
359k
      int pp=pitch+1-i;
510
359k
      tmp1=nsf;
511
359k
      if (tmp1>pp)
512
129k
         tmp1=pp;
513
12.7M
      for (j=0;j<tmp1;j++)
514
12.4M
         exc[j]=MAC16_16(exc[j],SHL16(gain[2-i],7),exc2[j-pp]);
515
359k
      tmp3=nsf;
516
359k
      if (tmp3>pp+pitch)
517
39.6k
         tmp3=pp+pitch;
518
2.17M
      for (j=tmp1;j<tmp3;j++)
519
1.81M
         exc[j]=MAC16_16(exc[j],SHL16(gain[2-i],7),exc2[j-pp-pitch]);
520
359k
   }
521
4.91M
   for (i=0;i<nsf;i++)
522
4.79M
   {
523
4.79M
      spx_word32_t tmp = ADD32(ADD32(MULT16_16(gain[0],x[2][i]),MULT16_16(gain[1],x[1][i])),
524
4.79M
                            MULT16_16(gain[2],x[0][i]));
525
4.79M
      new_target[i] = SUB16(new_target[i], EXTRACT16(PSHR32(tmp,6)));
526
4.79M
   }
527
119k
   err = inner_prod(new_target, new_target, nsf);
528
529
119k
   return err;
530
119k
}
531
532
/** Finds the best quantized 3-tap pitch predictor by analysis by synthesis */
533
int pitch_search_3tap(
534
spx_word16_t target[],                 /* Target vector */
535
spx_word16_t *sw,
536
spx_coef_t ak[],                     /* LPCs for this subframe */
537
spx_coef_t awk1[],                   /* Weighted LPCs #1 for this subframe */
538
spx_coef_t awk2[],                   /* Weighted LPCs #2 for this subframe */
539
spx_sig_t exc[],                    /* Excitation */
540
const void *par,
541
int   start,                    /* Smallest pitch value allowed */
542
int   end,                      /* Largest pitch value allowed */
543
spx_word16_t pitch_coef,               /* Voicing (pitch) coefficient */
544
int   p,                        /* Number of LPC coeffs */
545
int   nsf,                      /* Number of samples in subframe */
546
SpeexBits *bits,
547
char *stack,
548
spx_word16_t *exc2,
549
spx_word16_t *r,
550
int complexity,
551
int cdbk_offset,
552
int plc_tuning,
553
spx_word32_t *cumul_gain
554
)
555
29.1k
{
556
29.1k
   int i;
557
29.1k
   int cdbk_index, pitch=0, best_gain_index=0;
558
29.1k
   VARDECL(spx_sig_t *best_exc);
559
29.1k
   VARDECL(spx_word16_t *new_target);
560
29.1k
   VARDECL(spx_word16_t *best_target);
561
29.1k
   int best_pitch=0;
562
29.1k
   spx_word32_t err, best_err=-1;
563
29.1k
   int N;
564
29.1k
   const ltp_params *params;
565
29.1k
   const signed char *gain_cdbk;
566
29.1k
   int   gain_cdbk_size;
567
29.1k
   int scaledown=0;
568
569
29.1k
   VARDECL(int *nbest);
570
571
29.1k
   params = (const ltp_params*) par;
572
29.1k
   gain_cdbk_size = 1<<params->gain_bits;
573
29.1k
   gain_cdbk = params->gain_cdbk + 4*gain_cdbk_size*cdbk_offset;
574
575
29.1k
   N=complexity;
576
29.1k
   if (N>10)
577
0
      N=10;
578
29.1k
   if (N<1)
579
5.62k
      N=1;
580
581
29.1k
   ALLOC(nbest, N, int);
582
29.1k
   params = (const ltp_params*) par;
583
584
29.1k
   if (end<start)
585
1.63k
   {
586
1.63k
      speex_bits_pack(bits, 0, params->pitch_bits);
587
1.63k
      speex_bits_pack(bits, 0, params->gain_bits);
588
1.63k
      SPEEX_MEMSET(exc, 0, nsf);
589
1.63k
      return start;
590
1.63k
   }
591
592
27.4k
#ifdef FIXED_POINT
593
   /* Check if we need to scale everything down in the pitch search to avoid overflows */
594
981k
   for (i=0;i<nsf;i++)
595
959k
   {
596
959k
      if (ABS16(target[i])>16383)
597
4.69k
      {
598
4.69k
         scaledown=1;
599
4.69k
         break;
600
4.69k
      }
601
959k
   }
602
2.96M
   for (i=-end;i<0;i++)
603
2.94M
   {
604
2.94M
      if (ABS16(exc2[i])>16383)
605
4.31k
      {
606
4.31k
         scaledown=1;
607
4.31k
         break;
608
4.31k
      }
609
2.94M
   }
610
27.4k
#endif
611
27.4k
   if (N>end-start+1)
612
2.71k
      N=end-start+1;
613
27.4k
   if (end != start)
614
23.6k
      open_loop_nbest_pitch(sw, start, end, nsf, nbest, NULL, N, stack);
615
3.79k
   else
616
3.79k
      nbest[0] = start;
617
618
27.4k
   ALLOC(best_exc, nsf, spx_sig_t);
619
27.4k
   ALLOC(new_target, nsf, spx_word16_t);
620
27.4k
   ALLOC(best_target, nsf, spx_word16_t);
621
622
147k
   for (i=0;i<N;i++)
623
119k
   {
624
119k
      pitch=nbest[i];
625
119k
      SPEEX_MEMSET(exc, 0, nsf);
626
119k
      err=pitch_gain_search_3tap(target, ak, awk1, awk2, exc, gain_cdbk, gain_cdbk_size, pitch, p, nsf,
627
119k
                                 bits, stack, exc2, r, new_target, &cdbk_index, plc_tuning, *cumul_gain, scaledown);
628
119k
      if (err<best_err || best_err<0)
629
37.1k
      {
630
37.1k
         SPEEX_COPY(best_exc, exc, nsf);
631
37.1k
         SPEEX_COPY(best_target, new_target, nsf);
632
37.1k
         best_err=err;
633
37.1k
         best_pitch=pitch;
634
37.1k
         best_gain_index=cdbk_index;
635
37.1k
      }
636
119k
   }
637
   /*printf ("pitch: %d %d\n", best_pitch, best_gain_index);*/
638
27.4k
   speex_bits_pack(bits, best_pitch-start, params->pitch_bits);
639
27.4k
   speex_bits_pack(bits, best_gain_index, params->gain_bits);
640
27.4k
#ifdef FIXED_POINT
641
27.4k
   *cumul_gain = MULT16_32_Q13(SHL16(params->gain_cdbk[4*best_gain_index+3],8), MAX32(1024,*cumul_gain));
642
#else
643
   *cumul_gain = 0.03125*MAX32(1024,*cumul_gain)*params->gain_cdbk[4*best_gain_index+3];
644
#endif
645
   /*printf ("%f\n", cumul_gain);*/
646
   /*printf ("encode pitch: %d %d\n", best_pitch, best_gain_index);*/
647
27.4k
   SPEEX_COPY(exc, best_exc, nsf);
648
27.4k
   SPEEX_COPY(target, best_target, nsf);
649
27.4k
#ifdef FIXED_POINT
650
   /* Scale target back up if needed */
651
27.4k
   if (scaledown)
652
6.29k
   {
653
258k
      for (i=0;i<nsf;i++)
654
251k
         target[i]=SHL16(target[i],1);
655
6.29k
   }
656
27.4k
#endif
657
27.4k
   return pitch;
658
29.1k
}
659
#endif /* DISABLE_ENCODER */
660
661
#ifndef DISABLE_DECODER
662
void pitch_unquant_3tap(
663
spx_word16_t exc[],             /* Input excitation */
664
spx_word32_t exc_out[],         /* Output excitation */
665
int   start,                    /* Smallest pitch value allowed */
666
int   end,                      /* Largest pitch value allowed */
667
spx_word16_t pitch_coef,        /* Voicing (pitch) coefficient */
668
const void *par,
669
int   nsf,                      /* Number of samples in subframe */
670
int *pitch_val,
671
spx_word16_t *gain_val,
672
SpeexBits *bits,
673
char *stack,
674
int count_lost,
675
int subframe_offset,
676
spx_word16_t last_pitch_gain,
677
int cdbk_offset
678
)
679
0
{
680
0
   int i;
681
0
   int pitch;
682
0
   int gain_index;
683
0
   spx_word16_t gain[3];
684
0
   const signed char *gain_cdbk;
685
0
   int gain_cdbk_size;
686
0
   const ltp_params *params;
687
688
0
   params = (const ltp_params*) par;
689
0
   gain_cdbk_size = 1<<params->gain_bits;
690
0
   gain_cdbk = params->gain_cdbk + 4*gain_cdbk_size*cdbk_offset;
691
692
0
   pitch = speex_bits_unpack_unsigned(bits, params->pitch_bits);
693
0
   pitch += start;
694
0
   gain_index = speex_bits_unpack_unsigned(bits, params->gain_bits);
695
   /*printf ("decode pitch: %d %d\n", pitch, gain_index);*/
696
0
#ifdef FIXED_POINT
697
0
   gain[0] = ADD16(32,(spx_word16_t)gain_cdbk[gain_index*4]);
698
0
   gain[1] = ADD16(32,(spx_word16_t)gain_cdbk[gain_index*4+1]);
699
0
   gain[2] = ADD16(32,(spx_word16_t)gain_cdbk[gain_index*4+2]);
700
#else
701
   gain[0] = 0.015625*gain_cdbk[gain_index*4]+.5;
702
   gain[1] = 0.015625*gain_cdbk[gain_index*4+1]+.5;
703
   gain[2] = 0.015625*gain_cdbk[gain_index*4+2]+.5;
704
#endif
705
706
0
   if (count_lost && pitch > subframe_offset)
707
0
   {
708
0
      spx_word16_t gain_sum;
709
0
      if (1) {
710
0
#ifdef FIXED_POINT
711
0
         spx_word16_t tmp = count_lost < 4 ? last_pitch_gain : SHR16(last_pitch_gain,1);
712
0
         if (tmp>62)
713
0
            tmp=62;
714
#else
715
         spx_word16_t tmp = count_lost < 4 ? last_pitch_gain : 0.5 * last_pitch_gain;
716
         if (tmp>.95)
717
            tmp=.95;
718
#endif
719
0
         gain_sum = gain_3tap_to_1tap(gain);
720
721
0
         if (gain_sum > tmp)
722
0
         {
723
0
            spx_word16_t fact = DIV32_16(SHL32(EXTEND32(tmp),14),gain_sum);
724
0
            for (i=0;i<3;i++)
725
0
               gain[i]=MULT16_16_Q14(fact,gain[i]);
726
0
         }
727
728
0
      }
729
730
0
   }
731
732
0
   *pitch_val = pitch;
733
0
   gain_val[0]=gain[0];
734
0
   gain_val[1]=gain[1];
735
0
   gain_val[2]=gain[2];
736
0
   gain[0] = SHL16(gain[0],7);
737
0
   gain[1] = SHL16(gain[1],7);
738
0
   gain[2] = SHL16(gain[2],7);
739
0
   SPEEX_MEMSET(exc_out, 0, nsf);
740
0
   for (i=0;i<3;i++)
741
0
   {
742
0
      int j;
743
0
      int tmp1, tmp3;
744
0
      int pp=pitch+1-i;
745
0
      tmp1=nsf;
746
0
      if (tmp1>pp)
747
0
         tmp1=pp;
748
0
      for (j=0;j<tmp1;j++)
749
0
         exc_out[j]=MAC16_16(exc_out[j],gain[2-i],exc[j-pp]);
750
0
      tmp3=nsf;
751
0
      if (tmp3>pp+pitch)
752
0
         tmp3=pp+pitch;
753
0
      for (j=tmp1;j<tmp3;j++)
754
0
         exc_out[j]=MAC16_16(exc_out[j],gain[2-i],exc[j-pp-pitch]);
755
0
   }
756
   /*for (i=0;i<nsf;i++)
757
   exc[i]=PSHR32(exc32[i],13);*/
758
0
}
759
#endif /* DISABLE_DECODER */
760
761
#ifndef DISABLE_ENCODER
762
/** Forced pitch delay and gain */
763
int forced_pitch_quant(
764
spx_word16_t target[],                 /* Target vector */
765
spx_word16_t *sw,
766
spx_coef_t ak[],                     /* LPCs for this subframe */
767
spx_coef_t awk1[],                   /* Weighted LPCs #1 for this subframe */
768
spx_coef_t awk2[],                   /* Weighted LPCs #2 for this subframe */
769
spx_sig_t exc[],                    /* Excitation */
770
const void *par,
771
int   start,                    /* Smallest pitch value allowed */
772
int   end,                      /* Largest pitch value allowed */
773
spx_word16_t pitch_coef,               /* Voicing (pitch) coefficient */
774
int   p,                        /* Number of LPC coeffs */
775
int   nsf,                      /* Number of samples in subframe */
776
SpeexBits *bits,
777
char *stack,
778
spx_word16_t *exc2,
779
spx_word16_t *r,
780
int complexity,
781
int cdbk_offset,
782
int plc_tuning,
783
spx_word32_t *cumul_gain
784
)
785
47.0k
{
786
47.0k
   int i;
787
47.0k
   VARDECL(spx_word16_t *res);
788
47.0k
   ALLOC(res, nsf, spx_word16_t);
789
47.0k
#ifdef FIXED_POINT
790
47.0k
   if (pitch_coef>63)
791
544
      pitch_coef=63;
792
#else
793
   if (pitch_coef>.99)
794
      pitch_coef=.99;
795
#endif
796
1.55M
   for (i=0;i<nsf&&i<start;i++)
797
1.51M
   {
798
1.51M
      exc[i]=MULT16_16(SHL16(pitch_coef, 7),exc2[i-start]);
799
1.51M
   }
800
414k
   for (;i<nsf;i++)
801
367k
   {
802
367k
      exc[i]=MULT16_32_Q15(SHL16(pitch_coef, 9),exc[i-start]);
803
367k
   }
804
1.92M
   for (i=0;i<nsf;i++)
805
1.88M
      res[i] = EXTRACT16(PSHR32(exc[i], SIG_SHIFT-1));
806
47.0k
   syn_percep_zero16(res, ak, awk1, awk2, res, nsf, p, stack);
807
1.92M
   for (i=0;i<nsf;i++)
808
1.88M
      target[i]=EXTRACT16(SATURATE(SUB32(EXTEND32(target[i]),EXTEND32(res[i])),32700));
809
47.0k
   return start;
810
47.0k
}
811
#endif /* DISABLE_ENCODER */
812
813
#ifndef DISABLE_DECODER
814
/** Unquantize forced pitch delay and gain */
815
void forced_pitch_unquant(
816
spx_word16_t exc[],             /* Input excitation */
817
spx_word32_t exc_out[],         /* Output excitation */
818
int   start,                    /* Smallest pitch value allowed */
819
int   end,                      /* Largest pitch value allowed */
820
spx_word16_t pitch_coef,        /* Voicing (pitch) coefficient */
821
const void *par,
822
int   nsf,                      /* Number of samples in subframe */
823
int *pitch_val,
824
spx_word16_t *gain_val,
825
SpeexBits *bits,
826
char *stack,
827
int count_lost,
828
int subframe_offset,
829
spx_word16_t last_pitch_gain,
830
int cdbk_offset
831
)
832
0
{
833
0
   int i;
834
0
#ifdef FIXED_POINT
835
0
   if (pitch_coef>63)
836
0
      pitch_coef=63;
837
#else
838
   if (pitch_coef>.99)
839
      pitch_coef=.99;
840
#endif
841
0
   for (i=0;i<nsf;i++)
842
0
   {
843
0
      exc_out[i]=MULT16_16(exc[i-start],SHL16(pitch_coef,7));
844
0
      exc[i] = EXTRACT16(PSHR32(exc_out[i],13));
845
0
   }
846
0
   *pitch_val = start;
847
0
   gain_val[0]=gain_val[2]=0;
848
0
   gain_val[1] = pitch_coef;
849
0
}
850
#endif /* DISABLE_DECODER */