| 1 | /* |
| 2 | * LPC utility code |
| 3 | * Copyright (c) 2006 Justin Ruggles <justin.ruggles@gmail.com> |
| 4 | * |
| 5 | * This file is part of FFmpeg. |
| 6 | * |
| 7 | * FFmpeg is free software; you can redistribute it and/or |
| 8 | * modify it under the terms of the GNU Lesser General Public |
| 9 | * License as published by the Free Software Foundation; either |
| 10 | * version 2.1 of the License, or (at your option) any later version. |
| 11 | * |
| 12 | * FFmpeg is distributed in the hope that it will be useful, |
| 13 | * but WITHOUT ANY WARRANTY; without even the implied warranty of |
| 14 | * MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the GNU |
| 15 | * Lesser General Public License for more details. |
| 16 | * |
| 17 | * You should have received a copy of the GNU Lesser General Public |
| 18 | * License along with FFmpeg; if not, write to the Free Software |
| 19 | * Foundation, Inc., 51 Franklin Street, Fifth Floor, Boston, MA 02110-1301 USA |
| 20 | */ |
| 21 | |
| 22 | #ifndef AVCODEC_LPC_H |
| 23 | #define AVCODEC_LPC_H |
| 24 | |
| 25 | #include <stdint.h> |
| 26 | #include "libavutil/avassert.h" |
| 27 | |
| 28 | #define ORDER_METHOD_EST 0 |
| 29 | #define ORDER_METHOD_2LEVEL 1 |
| 30 | #define ORDER_METHOD_4LEVEL 2 |
| 31 | #define ORDER_METHOD_8LEVEL 3 |
| 32 | #define ORDER_METHOD_SEARCH 4 |
| 33 | #define ORDER_METHOD_LOG 5 |
| 34 | |
| 35 | #define MIN_LPC_ORDER 1 |
| 36 | #define MAX_LPC_ORDER 32 |
| 37 | |
| 38 | /** |
| 39 | * LPC analysis type |
| 40 | */ |
| 41 | enum FFLPCType { |
| 42 | FF_LPC_TYPE_DEFAULT = -1, ///< use the codec default LPC type |
| 43 | FF_LPC_TYPE_NONE = 0, ///< do not use LPC prediction or use all zero coefficients |
| 44 | FF_LPC_TYPE_FIXED = 1, ///< fixed LPC coefficients |
| 45 | FF_LPC_TYPE_LEVINSON = 2, ///< Levinson-Durbin recursion |
| 46 | FF_LPC_TYPE_CHOLESKY = 3, ///< Cholesky factorization |
| 47 | FF_LPC_TYPE_NB , ///< Not part of ABI |
| 48 | }; |
| 49 | |
| 50 | typedef struct LPCContext { |
| 51 | int blocksize; |
| 52 | int max_order; |
| 53 | enum FFLPCType lpc_type; |
| 54 | double *windowed_buffer; |
| 55 | double *windowed_samples; |
| 56 | |
| 57 | /** |
| 58 | * Apply a Welch window to an array of input samples. |
| 59 | * The output samples have the same scale as the input, but are in double |
| 60 | * sample format. |
| 61 | * @param data input samples |
| 62 | * @param len number of input samples |
| 63 | * @param w_data output samples |
| 64 | */ |
| 65 | void (*lpc_apply_welch_window)(const int32_t *data, int len, |
| 66 | double *w_data); |
| 67 | /** |
| 68 | * Perform autocorrelation on input samples with delay of 0 to lag. |
| 69 | * @param data input samples. |
| 70 | * constraints: no alignment needed, but must have at |
| 71 | * least lag*sizeof(double) valid bytes preceding it, and |
| 72 | * size must be at least (len+1)*sizeof(double) if data is |
| 73 | * 16-byte aligned or (len+2)*sizeof(double) if data is |
| 74 | * unaligned. |
| 75 | * @param len number of input samples to process |
| 76 | * @param lag maximum delay to calculate |
| 77 | * @param autoc output autocorrelation coefficients. |
| 78 | * constraints: array size must be at least lag+1. |
| 79 | */ |
| 80 | void (*lpc_compute_autocorr)(const double *data, int len, int lag, |
| 81 | double *autoc); |
| 82 | } LPCContext; |
| 83 | |
| 84 | |
| 85 | /** |
| 86 | * Calculate LPC coefficients for multiple orders |
| 87 | */ |
| 88 | int ff_lpc_calc_coefs(LPCContext *s, |
| 89 | const int32_t *samples, int blocksize, int min_order, |
| 90 | int max_order, int precision, |
| 91 | int32_t coefs[][MAX_LPC_ORDER], int *shift, |
| 92 | enum FFLPCType lpc_type, int lpc_passes, |
| 93 | int omethod, int max_shift, int zero_shift); |
| 94 | |
| 95 | int ff_lpc_calc_ref_coefs(LPCContext *s, |
| 96 | const int32_t *samples, int order, double *ref); |
| 97 | |
| 98 | /** |
| 99 | * Initialize LPCContext. |
| 100 | */ |
| 101 | int ff_lpc_init(LPCContext *s, int blocksize, int max_order, |
| 102 | enum FFLPCType lpc_type); |
| 103 | void ff_lpc_init_x86(LPCContext *s); |
| 104 | |
| 105 | /** |
| 106 | * Uninitialize LPCContext. |
| 107 | */ |
| 108 | void ff_lpc_end(LPCContext *s); |
| 109 | |
| 110 | #ifdef LPC_USE_DOUBLE |
| 111 | #define LPC_TYPE double |
| 112 | #else |
| 113 | #define LPC_TYPE float |
| 114 | #endif |
| 115 | |
| 116 | /** |
| 117 | * Schur recursion. |
| 118 | * Produces reflection coefficients from autocorrelation data. |
| 119 | */ |
| 120 | static inline void compute_ref_coefs(const LPC_TYPE *autoc, int max_order, |
| 121 | LPC_TYPE *ref, LPC_TYPE *error) |
| 122 | { |
| 123 | int i, j; |
| 124 | LPC_TYPE err; |
| 125 | LPC_TYPE gen0[MAX_LPC_ORDER], gen1[MAX_LPC_ORDER]; |
| 126 | |
| 127 | for (i = 0; i < max_order; i++) |
| 128 | gen0[i] = gen1[i] = autoc[i + 1]; |
| 129 | |
| 130 | err = autoc[0]; |
| 131 | ref[0] = -gen1[0] / err; |
| 132 | err += gen1[0] * ref[0]; |
| 133 | if (error) |
| 134 | error[0] = err; |
| 135 | for (i = 1; i < max_order; i++) { |
| 136 | for (j = 0; j < max_order - i; j++) { |
| 137 | gen1[j] = gen1[j + 1] + ref[i - 1] * gen0[j]; |
| 138 | gen0[j] = gen1[j + 1] * ref[i - 1] + gen0[j]; |
| 139 | } |
| 140 | ref[i] = -gen1[0] / err; |
| 141 | err += gen1[0] * ref[i]; |
| 142 | if (error) |
| 143 | error[i] = err; |
| 144 | } |
| 145 | } |
| 146 | |
| 147 | /** |
| 148 | * Levinson-Durbin recursion. |
| 149 | * Produce LPC coefficients from autocorrelation data. |
| 150 | */ |
| 151 | static inline int compute_lpc_coefs(const LPC_TYPE *autoc, int max_order, |
| 152 | LPC_TYPE *lpc, int lpc_stride, int fail, |
| 153 | int normalize) |
| 154 | { |
| 155 | int i, j; |
| 156 | LPC_TYPE err; |
| 157 | LPC_TYPE *lpc_last = lpc; |
| 158 | |
| 159 | av_assert2(normalize || !fail); |
| 160 | |
| 161 | if (normalize) |
| 162 | err = *autoc++; |
| 163 | |
| 164 | if (fail && (autoc[max_order - 1] == 0 || err <= 0)) |
| 165 | return -1; |
| 166 | |
| 167 | for(i=0; i<max_order; i++) { |
| 168 | LPC_TYPE r = -autoc[i]; |
| 169 | |
| 170 | if (normalize) { |
| 171 | for(j=0; j<i; j++) |
| 172 | r -= lpc_last[j] * autoc[i-j-1]; |
| 173 | |
| 174 | r /= err; |
| 175 | err *= 1.0 - (r * r); |
| 176 | } |
| 177 | |
| 178 | lpc[i] = r; |
| 179 | |
| 180 | for(j=0; j < (i+1)>>1; j++) { |
| 181 | LPC_TYPE f = lpc_last[ j]; |
| 182 | LPC_TYPE b = lpc_last[i-1-j]; |
| 183 | lpc[ j] = f + r * b; |
| 184 | lpc[i-1-j] = b + r * f; |
| 185 | } |
| 186 | |
| 187 | if (fail && err < 0) |
| 188 | return -1; |
| 189 | |
| 190 | lpc_last = lpc; |
| 191 | lpc += lpc_stride; |
| 192 | } |
| 193 | |
| 194 | return 0; |
| 195 | } |
| 196 | |
| 197 | #endif /* AVCODEC_LPC_H */ |