Mercurial > hg > audiostuff

comparison spandsp-0.0.6pre17/src/g726.c @ 4:26cd8f1ef0b1

Find changesets by keywords (author, files, the commit message), revision number or hash, or revset expression.
import spandsp-0.0.6pre17
author Peter Meerwald <pmeerw@cosy.sbg.ac.at>
date Fri, 25 Jun 2010 15:50:58 +0200
parents
children
comparison
equal deleted inserted replaced
3:c6c5a16ce2f2 4:26cd8f1ef0b1
1 /*
2 * SpanDSP - a series of DSP components for telephony
3 *
4 * g726.c - The ITU G.726 codec.
5 *
6 * Written by Steve Underwood <steveu@coppice.org>
7 *
8 * Copyright (C) 2006 Steve Underwood
9 *
10 * All rights reserved.
11 *
12 * This program is free software; you can redistribute it and/or modify
13 * it under the terms of the GNU Lesser General Public License version 2.1,
14 * as published by the Free Software Foundation.
15 *
16 * This program is distributed in the hope that it will be useful,
17 * but WITHOUT ANY WARRANTY; without even the implied warranty of
18 * MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the
19 * GNU Lesser General Public License for more details.
20 *
21 * You should have received a copy of the GNU Lesser General Public
22 * License along with this program; if not, write to the Free Software
23 * Foundation, Inc., 675 Mass Ave, Cambridge, MA 02139, USA.
24 *
25 * Based on G.721/G.723 code which is:
26 *
27 * This source code is a product of Sun Microsystems, Inc. and is provided
28 * for unrestricted use. Users may copy or modify this source code without
29 * charge.
30 *
31 * SUN SOURCE CODE IS PROVIDED AS IS WITH NO WARRANTIES OF ANY KIND INCLUDING
32 * THE WARRANTIES OF DESIGN, MERCHANTIBILITY AND FITNESS FOR A PARTICULAR
33 * PURPOSE, OR ARISING FROM A COURSE OF DEALING, USAGE OR TRADE PRACTICE.
34 *
35 * Sun source code is provided with no support and without any obligation on
36 * the part of Sun Microsystems, Inc. to assist in its use, correction,
37 * modification or enhancement.
38 *
39 * SUN MICROSYSTEMS, INC. SHALL HAVE NO LIABILITY WITH RESPECT TO THE
40 * INFRINGEMENT OF COPYRIGHTS, TRADE SECRETS OR ANY PATENTS BY THIS SOFTWARE
41 * OR ANY PART THEREOF.
42 *
43 * In no event will Sun Microsystems, Inc. be liable for any lost revenue
44 * or profits or other special, indirect and consequential damages, even if
45 * Sun has been advised of the possibility of such damages.
46 *
47 * Sun Microsystems, Inc.
48 * 2550 Garcia Avenue
49 * Mountain View, California 94043
50 *
51 * $Id: g726.c,v 1.28.4.1 2009/12/28 12:20:46 steveu Exp $
52 */
53
54 /*! \file */
55
56 #if defined(HAVE_CONFIG_H)
57 #include "config.h"
58 #endif
59
60 #include <inttypes.h>
61 #include <memory.h>
62 #include <stdlib.h>
63 #if defined(HAVE_TGMATH_H)
64 #include <tgmath.h>
65 #endif
66 #if defined(HAVE_MATH_H)
67 #include <math.h>
68 #endif
69 #include "floating_fudge.h"
70
71 #include "spandsp/telephony.h"
72 #include "spandsp/dc_restore.h"
73 #include "spandsp/bitstream.h"
74 #include "spandsp/bit_operations.h"
75 #include "spandsp/g711.h"
76 #include "spandsp/g726.h"
77
78 #include "spandsp/private/bitstream.h"
79 #include "spandsp/private/g726.h"
80
81 /*
82 * Maps G.726_16 code word to reconstructed scale factor normalized log
83 * magnitude values.
84 */
85 static const int g726_16_dqlntab[4] =
86 {
87 116, 365, 365, 116
88 };
89
90 /* Maps G.726_16 code word to log of scale factor multiplier. */
91 static const int g726_16_witab[4] =
92 {
93 -704, 14048, 14048, -704
94 };
95
96 /*
97 * Maps G.726_16 code words to a set of values whose long and short
98 * term averages are computed and then compared to give an indication
99 * how stationary (steady state) the signal is.
100 */
101 static const int g726_16_fitab[4] =
102 {
103 0x000, 0xE00, 0xE00, 0x000
104 };
105
106 static const int qtab_726_16[1] =
107 {
108 261
109 };
110
111 /*
112 * Maps G.726_24 code word to reconstructed scale factor normalized log
113 * magnitude values.
114 */
115 static const int g726_24_dqlntab[8] =
116 {
117 -2048, 135, 273, 373, 373, 273, 135, -2048
118 };
119
120 /* Maps G.726_24 code word to log of scale factor multiplier. */
121 static const int g726_24_witab[8] =
122 {
123 -128, 960, 4384, 18624, 18624, 4384, 960, -128
124 };
125
126 /*
127 * Maps G.726_24 code words to a set of values whose long and short
128 * term averages are computed and then compared to give an indication
129 * how stationary (steady state) the signal is.
130 */
131 static const int g726_24_fitab[8] =
132 {
133 0x000, 0x200, 0x400, 0xE00, 0xE00, 0x400, 0x200, 0x000
134 };
135
136 static const int qtab_726_24[3] =
137 {
138 8, 218, 331
139 };
140
141 /*
142 * Maps G.726_32 code word to reconstructed scale factor normalized log
143 * magnitude values.
144 */
145 static const int g726_32_dqlntab[16] =
146 {
147 -2048, 4, 135, 213, 273, 323, 373, 425,
148 425, 373, 323, 273, 213, 135, 4, -2048
149 };
150
151 /* Maps G.726_32 code word to log of scale factor multiplier. */
152 static const int g726_32_witab[16] =
153 {
154 -384, 576, 1312, 2048, 3584, 6336, 11360, 35904,
155 35904, 11360, 6336, 3584, 2048, 1312, 576, -384
156 };
157
158 /*
159 * Maps G.726_32 code words to a set of values whose long and short
160 * term averages are computed and then compared to give an indication
161 * how stationary (steady state) the signal is.
162 */
163 static const int g726_32_fitab[16] =
164 {
165 0x000, 0x000, 0x000, 0x200, 0x200, 0x200, 0x600, 0xE00,
166 0xE00, 0x600, 0x200, 0x200, 0x200, 0x000, 0x000, 0x000
167 };
168
169 static const int qtab_726_32[7] =
170 {
171 -124, 80, 178, 246, 300, 349, 400
172 };
173
174 /*
175 * Maps G.726_40 code word to ructeconstructed scale factor normalized log
176 * magnitude values.
177 */
178 static const int g726_40_dqlntab[32] =
179 {
180 -2048, -66, 28, 104, 169, 224, 274, 318,
181 358, 395, 429, 459, 488, 514, 539, 566,
182 566, 539, 514, 488, 459, 429, 395, 358,
183 318, 274, 224, 169, 104, 28, -66, -2048
184 };
185
186 /* Maps G.726_40 code word to log of scale factor multiplier. */
187 static const int g726_40_witab[32] =
188 {
189 448, 448, 768, 1248, 1280, 1312, 1856, 3200,
190 4512, 5728, 7008, 8960, 11456, 14080, 16928, 22272,
191 22272, 16928, 14080, 11456, 8960, 7008, 5728, 4512,
192 3200, 1856, 1312, 1280, 1248, 768, 448, 448
193 };
194
195 /*
196 * Maps G.726_40 code words to a set of values whose long and short
197 * term averages are computed and then compared to give an indication
198 * how stationary (steady state) the signal is.
199 */
200 static const int g726_40_fitab[32] =
201 {
202 0x000, 0x000, 0x000, 0x000, 0x000, 0x200, 0x200, 0x200,
203 0x200, 0x200, 0x400, 0x600, 0x800, 0xA00, 0xC00, 0xC00,
204 0xC00, 0xC00, 0xA00, 0x800, 0x600, 0x400, 0x200, 0x200,
205 0x200, 0x200, 0x200, 0x000, 0x000, 0x000, 0x000, 0x000
206 };
207
208 static const int qtab_726_40[15] =
209 {
210 -122, -16, 68, 139, 198, 250, 298, 339,
211 378, 413, 445, 475, 502, 528, 553
212 };
213
214 /*
215 * returns the integer product of the 14-bit integer "an" and
216 * "floating point" representation (4-bit exponent, 6-bit mantessa) "srn".
217 */
218 static int16_t fmult(int16_t an, int16_t srn)
219 {
220 int16_t anmag;
221 int16_t anexp;
222 int16_t anmant;
223 int16_t wanexp;
224 int16_t wanmant;
225 int16_t retval;
226
227 anmag = (an > 0) ? an : ((-an) & 0x1FFF);
228 anexp = (int16_t) (top_bit(anmag) - 5);
229 anmant = (anmag == 0) ? 32 : (anexp >= 0) ? (anmag >> anexp) : (anmag << -anexp);
230 wanexp = anexp + ((srn >> 6) & 0xF) - 13;
231
232 wanmant = (anmant*(srn & 0x3F) + 0x30) >> 4;
233 retval = (wanexp >= 0) ? ((wanmant << wanexp) & 0x7FFF) : (wanmant >> -wanexp);
234
235 return (((an ^ srn) < 0) ? -retval : retval);
236 }
237 /*- End of function --------------------------------------------------------*/
238
239 /*
240 * Compute the estimated signal from the 6-zero predictor.
241 */
242 static __inline__ int16_t predictor_zero(g726_state_t *s)
243 {
244 int i;
245 int sezi;
246
247 sezi = fmult(s->b[0] >> 2, s->dq[0]);
248 /* ACCUM */
249 for (i = 1; i < 6; i++)
250 sezi += fmult(s->b[i] >> 2, s->dq[i]);
251 return (int16_t) sezi;
252 }
253 /*- End of function --------------------------------------------------------*/
254
255 /*
256 * Computes the estimated signal from the 2-pole predictor.
257 */
258 static __inline__ int16_t predictor_pole(g726_state_t *s)
259 {
260 return (fmult(s->a[1] >> 2, s->sr[1]) + fmult(s->a[0] >> 2, s->sr[0]));
261 }
262 /*- End of function --------------------------------------------------------*/
263
264 /*
265 * Computes the quantization step size of the adaptive quantizer.
266 */
267 static int step_size(g726_state_t *s)
268 {
269 int y;
270 int dif;
271 int al;
272
273 if (s->ap >= 256)
274 return s->yu;
275 y = s->yl >> 6;
276 dif = s->yu - y;
277 al = s->ap >> 2;
278 if (dif > 0)
279 y += (dif*al) >> 6;
280 else if (dif < 0)
281 y += (dif*al + 0x3F) >> 6;
282 return y;
283 }
284 /*- End of function --------------------------------------------------------*/
285
286 /*
287 * Given a raw sample, 'd', of the difference signal and a
288 * quantization step size scale factor, 'y', this routine returns the
289 * ADPCM codeword to which that sample gets quantized. The step
290 * size scale factor division operation is done in the log base 2 domain
291 * as a subtraction.
292 */
293 static int16_t quantize(int d, /* Raw difference signal sample */
294 int y, /* Step size multiplier */
295 const int table[], /* quantization table */
296 int quantizer_states) /* table size of int16_t integers */
297 {
298 int16_t dqm; /* Magnitude of 'd' */
299 int16_t exp; /* Integer part of base 2 log of 'd' */
300 int16_t mant; /* Fractional part of base 2 log */
301 int16_t dl; /* Log of magnitude of 'd' */
302 int16_t dln; /* Step size scale factor normalized log */
303 int i;
304 int size;
305
306 /*
307 * LOG
308 *
309 * Compute base 2 log of 'd', and store in 'dl'.
310 */
311 dqm = (int16_t) abs(d);
312 exp = (int16_t) (top_bit(dqm >> 1) + 1);
313 /* Fractional portion. */
314 mant = ((dqm << 7) >> exp) & 0x7F;
315 dl = (exp << 7) + mant;
316
317 /*
318 * SUBTB
319 *
320 * "Divide" by step size multiplier.
321 */
322 dln = dl - (int16_t) (y >> 2);
323
324 /*
325 * QUAN
326 *
327 * Search for codword i for 'dln'.
328 */
329 size = (quantizer_states - 1) >> 1;
330 for (i = 0; i < size; i++)
331 {
332 if (dln < table[i])
333 break;
334 }
335 if (d < 0)
336 {
337 /* Take 1's complement of i */
338 return (int16_t) ((size << 1) + 1 - i);
339 }
340 if (i == 0 && (quantizer_states & 1))
341 {
342 /* Zero is only valid if there are an even number of states, so
343 take the 1's complement if the code is zero. */
344 return (int16_t) quantizer_states;
345 }
346 return (int16_t) i;
347 }
348 /*- End of function --------------------------------------------------------*/
349
350 /*
351 * Returns reconstructed difference signal 'dq' obtained from
352 * codeword 'i' and quantization step size scale factor 'y'.
353 * Multiplication is performed in log base 2 domain as addition.
354 */
355 static int16_t reconstruct(int sign, /* 0 for non-negative value */
356 int dqln, /* G.72x codeword */
357 int y) /* Step size multiplier */
358 {
359 int16_t dql; /* Log of 'dq' magnitude */
360 int16_t dex; /* Integer part of log */
361 int16_t dqt;
362 int16_t dq; /* Reconstructed difference signal sample */
363
364 dql = (int16_t) (dqln + (y >> 2)); /* ADDA */
365
366 if (dql < 0)
367 return ((sign) ? -0x8000 : 0);
368 /* ANTILOG */
369 dex = (dql >> 7) & 15;
370 dqt = 128 + (dql & 127);
371 dq = (dqt << 7) >> (14 - dex);
372 return ((sign) ? (dq - 0x8000) : dq);
373 }
374 /*- End of function --------------------------------------------------------*/
375
376 /*
377 * updates the state variables for each output code
378 */
379 static void update(g726_state_t *s,
380 int y, /* quantizer step size */
381 int wi, /* scale factor multiplier */
382 int fi, /* for long/short term energies */
383 int dq, /* quantized prediction difference */
384 int sr, /* reconstructed signal */
385 int dqsez) /* difference from 2-pole predictor */
386 {
387 int16_t mag;
388 int16_t exp;
389 int16_t a2p; /* LIMC */
390 int16_t a1ul; /* UPA1 */
391 int16_t pks1; /* UPA2 */
392 int16_t fa1;
393 int16_t ylint;
394 int16_t dqthr;
395 int16_t ylfrac;
396 int16_t thr;
397 int16_t pk0;
398 int i;
399 int tr;
400
401 a2p = 0;
402 /* Needed in updating predictor poles */
403 pk0 = (dqsez < 0) ? 1 : 0;
404
405 /* prediction difference magnitude */
406 mag = (int16_t) (dq & 0x7FFF);
407 /* TRANS */
408 ylint = (int16_t) (s->yl >> 15); /* exponent part of yl */
409 ylfrac = (int16_t) ((s->yl >> 10) & 0x1F); /* fractional part of yl */
410 /* Limit threshold to 31 << 10 */
411 thr = (ylint > 9) ? (31 << 10) : ((32 + ylfrac) << ylint);
412 dqthr = (thr + (thr >> 1)) >> 1; /* dqthr = 0.75 * thr */
413 if (!s->td) /* signal supposed voice */
414 tr = FALSE;
415 else if (mag <= dqthr) /* supposed data, but small mag */
416 tr = FALSE; /* treated as voice */
417 else /* signal is data (modem) */
418 tr = TRUE;
419
420 /*
421 * Quantizer scale factor adaptation.
422 */
423
424 /* FUNCTW & FILTD & DELAY */
425 /* update non-steady state step size multiplier */
426 s->yu = (int16_t) (y + ((wi - y) >> 5));
427
428 /* LIMB */
429 if (s->yu < 544)
430 s->yu = 544;
431 else if (s->yu > 5120)
432 s->yu = 5120;
433
434 /* FILTE & DELAY */
435 /* update steady state step size multiplier */
436 s->yl += s->yu + ((-s->yl) >> 6);
437
438 /*
439 * Adaptive predictor coefficients.
440 */
441 if (tr)
442 {
443 /* Reset the a's and b's for a modem signal */
444 s->a[0] = 0;
445 s->a[1] = 0;
446 s->b[0] = 0;
447 s->b[1] = 0;
448 s->b[2] = 0;
449 s->b[3] = 0;
450 s->b[4] = 0;
451 s->b[5] = 0;
452 }
453 else
454 {
455 /* Update the a's and b's */
456 /* UPA2 */
457 pks1 = pk0 ^ s->pk[0];
458
459 /* Update predictor pole a[1] */
460 a2p = s->a[1] - (s->a[1] >> 7);
461 if (dqsez != 0)
462 {
463 fa1 = (pks1) ? s->a[0] : -s->a[0];
464 /* a2p = function of fa1 */
465 if (fa1 < -8191)
466 a2p -= 0x100;
467 else if (fa1 > 8191)
468 a2p += 0xFF;
469 else
470 a2p += fa1 >> 5;
471
472 if (pk0 ^ s->pk[1])
473 {
474 /* LIMC */
475 if (a2p <= -12160)
476 a2p = -12288;
477 else if (a2p >= 12416)
478 a2p = 12288;
479 else
480 a2p -= 0x80;
481 }
482 else if (a2p <= -12416)
483 a2p = -12288;
484 else if (a2p >= 12160)
485 a2p = 12288;
486 else
487 a2p += 0x80;
488 }
489
490 /* TRIGB & DELAY */
491 s->a[1] = a2p;
492
493 /* UPA1 */
494 /* Update predictor pole a[0] */
495 s->a[0] -= s->a[0] >> 8;
496 if (dqsez != 0)
497 {
498 if (pks1 == 0)
499 s->a[0] += 192;
500 else
501 s->a[0] -= 192;
502 }
503 /* LIMD */
504 a1ul = 15360 - a2p;
505 if (s->a[0] < -a1ul)
506 s->a[0] = -a1ul;
507 else if (s->a[0] > a1ul)
508 s->a[0] = a1ul;
509
510 /* UPB : update predictor zeros b[6] */
511 for (i = 0; i < 6; i++)
512 {
513 /* Distinguish 40Kbps mode from the others */
514 s->b[i] -= s->b[i] >> ((s->bits_per_sample == 5) ? 9 : 8);
515 if (dq & 0x7FFF)
516 {
517 /* XOR */
518 if ((dq ^ s->dq[i]) >= 0)
519 s->b[i] += 128;
520 else
521 s->b[i] -= 128;
522 }
523 }
524 }
525
526 for (i = 5; i > 0; i--)
527 s->dq[i] = s->dq[i - 1];
528 /* FLOAT A : convert dq[0] to 4-bit exp, 6-bit mantissa f.p. */
529 if (mag == 0)
530 {
531 s->dq[0] = (dq >= 0) ? 0x20 : 0xFC20;
532 }
533 else
534 {
535 exp = (int16_t) (top_bit(mag) + 1);
536 s->dq[0] = (dq >= 0)
537 ? ((exp << 6) + ((mag << 6) >> exp))
538 : ((exp << 6) + ((mag << 6) >> exp) - 0x400);
539 }
540
541 s->sr[1] = s->sr[0];
542 /* FLOAT B : convert sr to 4-bit exp., 6-bit mantissa f.p. */
543 if (sr == 0)
544 {
545 s->sr[0] = 0x20;
546 }
547 else if (sr > 0)
548 {
549 exp = (int16_t) (top_bit(sr) + 1);
550 s->sr[0] = (int16_t) ((exp << 6) + ((sr << 6) >> exp));
551 }
552 else if (sr > -32768)
553 {
554 mag = (int16_t) -sr;
555 exp = (int16_t) (top_bit(mag) + 1);
556 s->sr[0] = (exp << 6) + ((mag << 6) >> exp) - 0x400;
557 }
558 else
559 {
560 s->sr[0] = (uint16_t) 0xFC20;
561 }
562
563 /* DELAY A */
564 s->pk[1] = s->pk[0];
565 s->pk[0] = pk0;
566
567 /* TONE */
568 if (tr) /* this sample has been treated as data */
569 s->td = FALSE; /* next one will be treated as voice */
570 else if (a2p < -11776) /* small sample-to-sample correlation */
571 s->td = TRUE; /* signal may be data */
572 else /* signal is voice */
573 s->td = FALSE;
574
575 /* Adaptation speed control. */
576 /* FILTA */
577 s->dms += ((int16_t) fi - s->dms) >> 5;
578 /* FILTB */
579 s->dml += (((int16_t) (fi << 2) - s->dml) >> 7);
580
581 if (tr)
582 s->ap = 256;
583 else if (y < 1536) /* SUBTC */
584 s->ap += (0x200 - s->ap) >> 4;
585 else if (s->td)
586 s->ap += (0x200 - s->ap) >> 4;
587 else if (abs((s->dms << 2) - s->dml) >= (s->dml >> 3))
588 s->ap += (0x200 - s->ap) >> 4;
589 else
590 s->ap += (-s->ap) >> 4;
591 }
592 /*- End of function --------------------------------------------------------*/
593
594 static int16_t tandem_adjust_alaw(int16_t sr, /* decoder output linear PCM sample */
595 int se, /* predictor estimate sample */
596 int y, /* quantizer step size */
597 int i, /* decoder input code */
598 int sign,
599 const int qtab[],
600 int quantizer_states)
601 {
602 uint8_t sp; /* A-law compressed 8-bit code */
603 int16_t dx; /* prediction error */
604 int id; /* quantized prediction error */
605 int sd; /* adjusted A-law decoded sample value */
606
607 if (sr <= -32768)
608 sr = -1;
609 sp = linear_to_alaw((sr >> 1) << 3);
610 /* 16-bit prediction error */
611 dx = (int16_t) ((alaw_to_linear(sp) >> 2) - se);
612 id = quantize(dx, y, qtab, quantizer_states);
613 if (id == i)
614 {
615 /* No adjustment of sp required */
616 return (int16_t) sp;
617 }
618 /* sp adjustment needed */
619 /* ADPCM codes : 8, 9, ... F, 0, 1, ... , 6, 7 */
620 /* 2's complement to biased unsigned */
621 if ((id ^ sign) > (i ^ sign))
622 {
623 /* sp adjusted to next lower value */
624 if (sp & 0x80)
625 sd = (sp == 0xD5) ? 0x55 : (((sp ^ 0x55) - 1) ^ 0x55);
626 else
627 sd = (sp == 0x2A) ? 0x2A : (((sp ^ 0x55) + 1) ^ 0x55);
628 }
629 else
630 {
631 /* sp adjusted to next higher value */
632 if (sp & 0x80)
633 sd = (sp == 0xAA) ? 0xAA : (((sp ^ 0x55) + 1) ^ 0x55);
634 else
635 sd = (sp == 0x55) ? 0xD5 : (((sp ^ 0x55) - 1) ^ 0x55);
636 }
637 return (int16_t) sd;
638 }
639 /*- End of function --------------------------------------------------------*/
640
641 static int16_t tandem_adjust_ulaw(int16_t sr, /* decoder output linear PCM sample */
642 int se, /* predictor estimate sample */
643 int y, /* quantizer step size */
644 int i, /* decoder input code */
645 int sign,
646 const int qtab[],
647 int quantizer_states)
648 {
649 uint8_t sp; /* u-law compressed 8-bit code */
650 int16_t dx; /* prediction error */
651 int id; /* quantized prediction error */
652 int sd; /* adjusted u-law decoded sample value */
653
654 if (sr <= -32768)
655 sr = 0;
656 sp = linear_to_ulaw(sr << 2);
657 /* 16-bit prediction error */
658 dx = (int16_t) ((ulaw_to_linear(sp) >> 2) - se);
659 id = quantize(dx, y, qtab, quantizer_states);
660 if (id == i)
661 {
662 /* No adjustment of sp required. */
663 return (int16_t) sp;
664 }
665 /* ADPCM codes : 8, 9, ... F, 0, 1, ... , 6, 7 */
666 /* 2's complement to biased unsigned */
667 if ((id ^ sign) > (i ^ sign))
668 {
669 /* sp adjusted to next lower value */
670 if (sp & 0x80)
671 sd = (sp == 0xFF) ? 0x7E : (sp + 1);
672 else
673 sd = (sp == 0x00) ? 0x00 : (sp - 1);
674 }
675 else
676 {
677 /* sp adjusted to next higher value */
678 if (sp & 0x80)
679 sd = (sp == 0x80) ? 0x80 : (sp - 1);
680 else
681 sd = (sp == 0x7F) ? 0xFE : (sp + 1);
682 }
683 return (int16_t) sd;
684 }
685 /*- End of function --------------------------------------------------------*/
686
687 /*
688 * Encodes a linear PCM, A-law or u-law input sample and returns its 3-bit code.
689 */
690 static uint8_t g726_16_encoder(g726_state_t *s, int16_t amp)
691 {
692 int y;
693 int16_t sei;
694 int16_t sezi;
695 int16_t se;
696 int16_t d;
697 int16_t sr;
698 int16_t dqsez;
699 int16_t dq;
700 int16_t i;
701
702 sezi = predictor_zero(s);
703 sei = sezi + predictor_pole(s);
704 se = sei >> 1;
705 d = amp - se;
706
707 /* Quantize prediction difference */
708 y = step_size(s);
709 i = quantize(d, y, qtab_726_16, 4);
710 dq = reconstruct(i & 2, g726_16_dqlntab[i], y);
711
712 /* Reconstruct the signal */
713 sr = (dq < 0) ? (se - (dq & 0x3FFF)) : (se + dq);
714
715 /* Pole prediction difference */
716 dqsez = sr + (sezi >> 1) - se;
717
718 update(s, y, g726_16_witab[i], g726_16_fitab[i], dq, sr, dqsez);
719 return (uint8_t) i;
720 }
721 /*- End of function --------------------------------------------------------*/
722
723 /*
724 * Decodes a 2-bit CCITT G.726_16 ADPCM code and returns
725 * the resulting 16-bit linear PCM, A-law or u-law sample value.
726 */
727 static int16_t g726_16_decoder(g726_state_t *s, uint8_t code)
728 {
729 int16_t sezi;
730 int16_t sei;
731 int16_t se;
732 int16_t sr;
733 int16_t dq;
734 int16_t dqsez;
735 int y;
736
737 /* Mask to get proper bits */
738 code &= 0x03;
739 sezi = predictor_zero(s);
740 sei = sezi + predictor_pole(s);
741
742 y = step_size(s);
743 dq = reconstruct(code & 2, g726_16_dqlntab[code], y);
744
745 /* Reconstruct the signal */
746 se = sei >> 1;
747 sr = (dq < 0) ? (se - (dq & 0x3FFF)) : (se + dq);
748
749 /* Pole prediction difference */
750 dqsez = sr + (sezi >> 1) - se;
751
752 update(s, y, g726_16_witab[code], g726_16_fitab[code], dq, sr, dqsez);
753
754 switch (s->ext_coding)
755 {
756 case G726_ENCODING_ALAW:
757 return tandem_adjust_alaw(sr, se, y, code, 2, qtab_726_16, 4);
758 case G726_ENCODING_ULAW:
759 return tandem_adjust_ulaw(sr, se, y, code, 2, qtab_726_16, 4);
760 }
761 return (sr << 2);
762 }
763 /*- End of function --------------------------------------------------------*/
764
765 /*
766 * Encodes a linear PCM, A-law or u-law input sample and returns its 3-bit code.
767 */
768 static uint8_t g726_24_encoder(g726_state_t *s, int16_t amp)
769 {
770 int16_t sei;
771 int16_t sezi;
772 int16_t se;
773 int16_t d;
774 int16_t sr;
775 int16_t dqsez;
776 int16_t dq;
777 int16_t i;
778 int y;
779
780 sezi = predictor_zero(s);
781 sei = sezi + predictor_pole(s);
782 se = sei >> 1;
783 d = amp - se;
784
785 /* Quantize prediction difference */
786 y = step_size(s);
787 i = quantize(d, y, qtab_726_24, 7);
788 dq = reconstruct(i & 4, g726_24_dqlntab[i], y);
789
790 /* Reconstruct the signal */
791 sr = (dq < 0) ? (se - (dq & 0x3FFF)) : (se + dq);
792
793 /* Pole prediction difference */
794 dqsez = sr + (sezi >> 1) - se;
795
796 update(s, y, g726_24_witab[i], g726_24_fitab[i], dq, sr, dqsez);
797 return (uint8_t) i;
798 }
799 /*- End of function --------------------------------------------------------*/
800
801 /*
802 * Decodes a 3-bit CCITT G.726_24 ADPCM code and returns
803 * the resulting 16-bit linear PCM, A-law or u-law sample value.
804 */
805 static int16_t g726_24_decoder(g726_state_t *s, uint8_t code)
806 {
807 int16_t sezi;
808 int16_t sei;
809 int16_t se;
810 int16_t sr;
811 int16_t dq;
812 int16_t dqsez;
813 int y;
814
815 /* Mask to get proper bits */
816 code &= 0x07;
817 sezi = predictor_zero(s);
818 sei = sezi + predictor_pole(s);
819
820 y = step_size(s);
821 dq = reconstruct(code & 4, g726_24_dqlntab[code], y);
822
823 /* Reconstruct the signal */
824 se = sei >> 1;
825 sr = (dq < 0) ? (se - (dq & 0x3FFF)) : (se + dq);
826
827 /* Pole prediction difference */
828 dqsez = sr + (sezi >> 1) - se;
829
830 update(s, y, g726_24_witab[code], g726_24_fitab[code], dq, sr, dqsez);
831
832 switch (s->ext_coding)
833 {
834 case G726_ENCODING_ALAW:
835 return tandem_adjust_alaw(sr, se, y, code, 4, qtab_726_24, 7);
836 case G726_ENCODING_ULAW:
837 return tandem_adjust_ulaw(sr, se, y, code, 4, qtab_726_24, 7);
838 }
839 return (sr << 2);
840 }
841 /*- End of function --------------------------------------------------------*/
842
843 /*
844 * Encodes a linear input sample and returns its 4-bit code.
845 */
846 static uint8_t g726_32_encoder(g726_state_t *s, int16_t amp)
847 {
848 int16_t sei;
849 int16_t sezi;
850 int16_t se;
851 int16_t d;
852 int16_t sr;
853 int16_t dqsez;
854 int16_t dq;
855 int16_t i;
856 int y;
857
858 sezi = predictor_zero(s);
859 sei = sezi + predictor_pole(s);
860 se = sei >> 1;
861 d = amp - se;
862
863 /* Quantize the prediction difference */
864 y = step_size(s);
865 i = quantize(d, y, qtab_726_32, 15);
866 dq = reconstruct(i & 8, g726_32_dqlntab[i], y);
867
868 /* Reconstruct the signal */
869 sr = (dq < 0) ? (se - (dq & 0x3FFF)) : (se + dq);
870
871 /* Pole prediction difference */
872 dqsez = sr + (sezi >> 1) - se;
873
874 update(s, y, g726_32_witab[i], g726_32_fitab[i], dq, sr, dqsez);
875 return (uint8_t) i;
876 }
877 /*- End of function --------------------------------------------------------*/
878
879 /*
880 * Decodes a 4-bit CCITT G.726_32 ADPCM code and returns
881 * the resulting 16-bit linear PCM, A-law or u-law sample value.
882 */
883 static int16_t g726_32_decoder(g726_state_t *s, uint8_t code)
884 {
885 int16_t sezi;
886 int16_t sei;
887 int16_t se;
888 int16_t sr;
889 int16_t dq;
890 int16_t dqsez;
891 int y;
892
893 /* Mask to get proper bits */
894 code &= 0x0F;
895 sezi = predictor_zero(s);
896 sei = sezi + predictor_pole(s);
897
898 y = step_size(s);
899 dq = reconstruct(code & 8, g726_32_dqlntab[code], y);
900
901 /* Reconstruct the signal */
902 se = sei >> 1;
903 sr = (dq < 0) ? (se - (dq & 0x3FFF)) : (se + dq);
904
905 /* Pole prediction difference */
906 dqsez = sr + (sezi >> 1) - se;
907
908 update(s, y, g726_32_witab[code], g726_32_fitab[code], dq, sr, dqsez);
909
910 switch (s->ext_coding)
911 {
912 case G726_ENCODING_ALAW:
913 return tandem_adjust_alaw(sr, se, y, code, 8, qtab_726_32, 15);
914 case G726_ENCODING_ULAW:
915 return tandem_adjust_ulaw(sr, se, y, code, 8, qtab_726_32, 15);
916 }
917 return (sr << 2);
918 }
919 /*- End of function --------------------------------------------------------*/
920
921 /*
922 * Encodes a 16-bit linear PCM, A-law or u-law input sample and retuens
923 * the resulting 5-bit CCITT G.726 40Kbps code.
924 */
925 static uint8_t g726_40_encoder(g726_state_t *s, int16_t amp)
926 {
927 int16_t sei;
928 int16_t sezi;
929 int16_t se;
930 int16_t d;
931 int16_t sr;
932 int16_t dqsez;
933 int16_t dq;
934 int16_t i;
935 int y;
936
937 sezi = predictor_zero(s);
938 sei = sezi + predictor_pole(s);
939 se = sei >> 1;
940 d = amp - se;
941
942 /* Quantize prediction difference */
943 y = step_size(s);
944 i = quantize(d, y, qtab_726_40, 31);
945 dq = reconstruct(i & 0x10, g726_40_dqlntab[i], y);
946
947 /* Reconstruct the signal */
948 sr = (dq < 0) ? (se - (dq & 0x7FFF)) : (se + dq);
949
950 /* Pole prediction difference */
951 dqsez = sr + (sezi >> 1) - se;
952
953 update(s, y, g726_40_witab[i], g726_40_fitab[i], dq, sr, dqsez);
954 return (uint8_t) i;
955 }
956 /*- End of function --------------------------------------------------------*/
957
958 /*
959 * Decodes a 5-bit CCITT G.726 40Kbps code and returns
960 * the resulting 16-bit linear PCM, A-law or u-law sample value.
961 */
962 static int16_t g726_40_decoder(g726_state_t *s, uint8_t code)
963 {
964 int16_t sezi;
965 int16_t sei;
966 int16_t se;
967 int16_t sr;
968 int16_t dq;
969 int16_t dqsez;
970 int y;
971
972 /* Mask to get proper bits */
973 code &= 0x1F;
974 sezi = predictor_zero(s);
975 sei = sezi + predictor_pole(s);
976
977 y = step_size(s);
978 dq = reconstruct(code & 0x10, g726_40_dqlntab[code], y);
979
980 /* Reconstruct the signal */
981 se = sei >> 1;
982 sr = (dq < 0) ? (se - (dq & 0x7FFF)) : (se + dq);
983
984 /* Pole prediction difference */
985 dqsez = sr + (sezi >> 1) - se;
986
987 update(s, y, g726_40_witab[code], g726_40_fitab[code], dq, sr, dqsez);
988
989 switch (s->ext_coding)
990 {
991 case G726_ENCODING_ALAW:
992 return tandem_adjust_alaw(sr, se, y, code, 0x10, qtab_726_40, 31);
993 case G726_ENCODING_ULAW:
994 return tandem_adjust_ulaw(sr, se, y, code, 0x10, qtab_726_40, 31);
995 }
996 return (sr << 2);
997 }
998 /*- End of function --------------------------------------------------------*/
999
1000 SPAN_DECLARE(g726_state_t *) g726_init(g726_state_t *s, int bit_rate, int ext_coding, int packing)
1001 {
1002 int i;
1003
1004 if (bit_rate != 16000 && bit_rate != 24000 && bit_rate != 32000 && bit_rate != 40000)
1005 return NULL;
1006 if (s == NULL)
1007 {
1008 if ((s = (g726_state_t *) malloc(sizeof(*s))) == NULL)
1009 return NULL;
1010 }
1011 s->yl = 34816;
1012 s->yu = 544;
1013 s->dms = 0;
1014 s->dml = 0;
1015 s->ap = 0;
1016 s->rate = bit_rate;
1017 s->ext_coding = ext_coding;
1018 s->packing = packing;
1019 for (i = 0; i < 2; i++)
1020 {
1021 s->a[i] = 0;
1022 s->pk[i] = 0;
1023 s->sr[i] = 32;
1024 }
1025 for (i = 0; i < 6; i++)
1026 {
1027 s->b[i] = 0;
1028 s->dq[i] = 32;
1029 }
1030 s->td = FALSE;
1031 switch (bit_rate)
1032 {
1033 case 16000:
1034 s->enc_func = g726_16_encoder;
1035 s->dec_func = g726_16_decoder;
1036 s->bits_per_sample = 2;
1037 break;
1038 case 24000:
1039 s->enc_func = g726_24_encoder;
1040 s->dec_func = g726_24_decoder;
1041 s->bits_per_sample = 3;
1042 break;
1043 case 32000:
1044 default:
1045 s->enc_func = g726_32_encoder;
1046 s->dec_func = g726_32_decoder;
1047 s->bits_per_sample = 4;
1048 break;
1049 case 40000:
1050 s->enc_func = g726_40_encoder;
1051 s->dec_func = g726_40_decoder;
1052 s->bits_per_sample = 5;
1053 break;
1054 }
1055 bitstream_init(&s->bs, (s->packing != G726_PACKING_LEFT));
1056 return s;
1057 }
1058 /*- End of function --------------------------------------------------------*/
1059
1060 SPAN_DECLARE(int) g726_release(g726_state_t *s)
1061 {
1062 return 0;
1063 }
1064 /*- End of function --------------------------------------------------------*/
1065
1066 SPAN_DECLARE(int) g726_free(g726_state_t *s)
1067 {
1068 free(s);
1069 return 0;
1070 }
1071 /*- End of function --------------------------------------------------------*/
1072
1073 SPAN_DECLARE(int) g726_decode(g726_state_t *s,
1074 int16_t amp[],
1075 const uint8_t g726_data[],
1076 int g726_bytes)
1077 {
1078 int i;
1079 int samples;
1080 uint8_t code;
1081 int sl;
1082
1083 for (samples = i = 0; ; )
1084 {
1085 if (s->packing != G726_PACKING_NONE)
1086 {
1087 /* Unpack the code bits */
1088 if (s->packing != G726_PACKING_LEFT)
1089 {
1090 if (s->bs.residue < s->bits_per_sample)
1091 {
1092 if (i >= g726_bytes)
1093 break;
1094 s->bs.bitstream |= (g726_data[i++] << s->bs.residue);
1095 s->bs.residue += 8;
1096 }
1097 code = (uint8_t) (s->bs.bitstream & ((1 << s->bits_per_sample) - 1));
1098 s->bs.bitstream >>= s->bits_per_sample;
1099 }
1100 else
1101 {
1102 if (s->bs.residue < s->bits_per_sample)
1103 {
1104 if (i >= g726_bytes)
1105 break;
1106 s->bs.bitstream = (s->bs.bitstream << 8) | g726_data[i++];
1107 s->bs.residue += 8;
1108 }
1109 code = (uint8_t) ((s->bs.bitstream >> (s->bs.residue - s->bits_per_sample)) & ((1 << s->bits_per_sample) - 1));
1110 }
1111 s->bs.residue -= s->bits_per_sample;
1112 }
1113 else
1114 {
1115 if (i >= g726_bytes)
1116 break;
1117 code = g726_data[i++];
1118 }
1119 sl = s->dec_func(s, code);
1120 if (s->ext_coding != G726_ENCODING_LINEAR)
1121 ((uint8_t *) amp)[samples++] = (uint8_t) sl;
1122 else
1123 amp[samples++] = (int16_t) sl;
1124 }
1125 return samples;
1126 }
1127 /*- End of function --------------------------------------------------------*/
1128
1129 SPAN_DECLARE(int) g726_encode(g726_state_t *s,
1130 uint8_t g726_data[],
1131 const int16_t amp[],
1132 int len)
1133 {
1134 int i;
1135 int g726_bytes;
1136 int16_t sl;
1137 uint8_t code;
1138
1139 for (g726_bytes = i = 0; i < len; i++)
1140 {
1141 /* Linearize the input sample to 14-bit PCM */
1142 switch (s->ext_coding)
1143 {
1144 case G726_ENCODING_ALAW:
1145 sl = alaw_to_linear(((const uint8_t *) amp)[i]) >> 2;
1146 break;
1147 case G726_ENCODING_ULAW:
1148 sl = ulaw_to_linear(((const uint8_t *) amp)[i]) >> 2;
1149 break;
1150 default:
1151 sl = amp[i] >> 2;
1152 break;
1153 }
1154 code = s->enc_func(s, sl);
1155 if (s->packing != G726_PACKING_NONE)
1156 {
1157 /* Pack the code bits */
1158 if (s->packing != G726_PACKING_LEFT)
1159 {
1160 s->bs.bitstream |= (code << s->bs.residue);
1161 s->bs.residue += s->bits_per_sample;
1162 if (s->bs.residue >= 8)
1163 {
1164 g726_data[g726_bytes++] = (uint8_t) (s->bs.bitstream & 0xFF);
1165 s->bs.bitstream >>= 8;
1166 s->bs.residue -= 8;
1167 }
1168 }
1169 else
1170 {
1171 s->bs.bitstream = (s->bs.bitstream << s->bits_per_sample) | code;
1172 s->bs.residue += s->bits_per_sample;
1173 if (s->bs.residue >= 8)
1174 {
1175 g726_data[g726_bytes++] = (uint8_t) ((s->bs.bitstream >> (s->bs.residue - 8)) & 0xFF);
1176 s->bs.residue -= 8;
1177 }
1178 }
1179 }
1180 else
1181 {
1182 g726_data[g726_bytes++] = (uint8_t) code;
1183 }
1184 }
1185 return g726_bytes;
1186 }
1187 /*- End of function --------------------------------------------------------*/
1188 /*- End of file ------------------------------------------------------------*/

Repositories maintained by Peter Meerwald, pmeerw@pmeerw.net.