FFmpeg  4.0
dcaadpcm.c
Go to the documentation of this file.
1 /*
2  * DCA ADPCM engine
3  * Copyright (C) 2017 Daniil Cherednik
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 
23 #include "dcaadpcm.h"
24 #include "dcaenc.h"
25 #include "dca_core.h"
26 #include "mathops.h"
27 
29 
30 //assume we have DCA_ADPCM_COEFFS values before x
31 static inline int64_t calc_corr(const int32_t *x, int len, int j, int k)
32 {
33  int n;
34  int64_t s = 0;
35  for (n = 0; n < len; n++)
36  s += MUL64(x[n-j], x[n-k]);
37  return s;
38 }
39 
40 static inline int64_t apply_filter(const int16_t a[DCA_ADPCM_COEFFS], const int64_t corr[15], const int32_t aa[10])
41 {
42  int64_t err = 0;
43  int64_t tmp = 0;
44 
45  err = corr[0];
46 
47  tmp += MUL64(a[0], corr[1]);
48  tmp += MUL64(a[1], corr[2]);
49  tmp += MUL64(a[2], corr[3]);
50  tmp += MUL64(a[3], corr[4]);
51 
52  tmp = norm__(tmp, 13);
53  tmp += tmp;
54 
55  err -= tmp;
56  tmp = 0;
57 
58  tmp += MUL64(corr[5], aa[0]);
59  tmp += MUL64(corr[6], aa[1]);
60  tmp += MUL64(corr[7], aa[2]);
61  tmp += MUL64(corr[8], aa[3]);
62 
63  tmp += MUL64(corr[9], aa[4]);
64  tmp += MUL64(corr[10], aa[5]);
65  tmp += MUL64(corr[11], aa[6]);
66 
67  tmp += MUL64(corr[12], aa[7]);
68  tmp += MUL64(corr[13], aa[8]);
69 
70  tmp += MUL64(corr[14], aa[9]);
71 
72  tmp = norm__(tmp, 26);
73 
74  err += tmp;
75 
76  return llabs(err);
77 }
78 
79 static int64_t find_best_filter(const DCAADPCMEncContext *s, const int32_t *in, int len)
80 {
81  const premultiplied_coeffs *precalc_data = s->private_data;
82  int i, j, k = 0;
83  int vq = -1;
84  int64_t err;
85  int64_t min_err = 1ll << 62;
86  int64_t corr[15];
87 
88  for (i = 0; i <= DCA_ADPCM_COEFFS; i++)
89  for (j = i; j <= DCA_ADPCM_COEFFS; j++)
90  corr[k++] = calc_corr(in+4, len, i, j);
91 
92  for (i = 0; i < DCA_ADPCM_VQCODEBOOK_SZ; i++) {
93  err = apply_filter(ff_dca_adpcm_vb[i], corr, *precalc_data);
94  if (err < min_err) {
95  min_err = err;
96  vq = i;
97  }
98  precalc_data++;
99  }
100 
101  return vq;
102 }
103 
104 static inline int64_t calc_prediction_gain(int pred_vq, const int32_t *in, int32_t *out, int len)
105 {
106  int i;
107  int32_t error;
108 
109  int64_t signal_energy = 0;
110  int64_t error_energy = 0;
111 
112  for (i = 0; i < len; i++) {
113  error = in[DCA_ADPCM_COEFFS + i] - ff_dcaadpcm_predict(pred_vq, in + i);
114  out[i] = error;
115  signal_energy += MUL64(in[DCA_ADPCM_COEFFS + i], in[DCA_ADPCM_COEFFS + i]);
116  error_energy += MUL64(error, error);
117  }
118 
119  if (!error_energy)
120  return -1;
121 
122  return signal_energy / error_energy;
123 }
124 
126 {
127  int pred_vq, i;
128  int32_t input_buffer[16 + DCA_ADPCM_COEFFS];
129  int32_t input_buffer2[16 + DCA_ADPCM_COEFFS];
130 
131  int32_t max = 0;
132  int shift_bits;
133  uint64_t pg = 0;
134 
135  for (i = 0; i < len + DCA_ADPCM_COEFFS; i++)
136  max |= FFABS(in[i]);
137 
138  // normalize input to simplify apply_filter
139  shift_bits = av_log2(max) - 11;
140 
141  for (i = 0; i < len + DCA_ADPCM_COEFFS; i++) {
142  input_buffer[i] = norm__(in[i], 7);
143  input_buffer2[i] = norm__(in[i], shift_bits);
144  }
145 
146  pred_vq = find_best_filter(s, input_buffer2, len);
147 
148  if (pred_vq < 0)
149  return -1;
150 
151  pg = calc_prediction_gain(pred_vq, input_buffer, diff, len);
152 
153  // Greater than 10db (10*log(10)) prediction gain to use ADPCM.
154  // TODO: Tune it.
155  if (pg < 10)
156  return -1;
157 
158  for (i = 0; i < len; i++)
159  diff[i] <<= 7;
160 
161  return pred_vq;
162 }
163 
165 {
166  int i, j, k;
167 
168  for (i = 0; i < DCA_ADPCM_VQCODEBOOK_SZ; i++) {
169  int id = 0;
170  int32_t t = 0;
171  for (j = 0; j < DCA_ADPCM_COEFFS; j++) {
172  for (k = j; k < DCA_ADPCM_COEFFS; k++) {
173  t = (int32_t)ff_dca_adpcm_vb[i][j] * (int32_t)ff_dca_adpcm_vb[i][k];
174  if (j != k)
175  t *= 2;
176  (*data)[id++] = t;
177  }
178  }
179  data++;
180  }
181 }
182 
183 int ff_dcaadpcm_do_real(int pred_vq_index,
184  softfloat quant, int32_t scale_factor, int32_t step_size,
185  const int32_t *prev_hist, const int32_t *in, int32_t *next_hist, int32_t *out,
186  int len, int32_t peak)
187 {
188  int i;
189  int64_t delta;
190  int32_t dequant_delta;
191  int32_t work_bufer[16 + DCA_ADPCM_COEFFS];
192 
193  memcpy(work_bufer, prev_hist, sizeof(int32_t) * DCA_ADPCM_COEFFS);
194 
195  for (i = 0; i < len; i++) {
196  work_bufer[DCA_ADPCM_COEFFS + i] = ff_dcaadpcm_predict(pred_vq_index, &work_bufer[i]);
197 
198  delta = (int64_t)in[i] - ((int64_t)work_bufer[DCA_ADPCM_COEFFS + i] << 7);
199 
200  out[i] = quantize_value(av_clip64(delta, -peak, peak), quant);
201 
202  ff_dca_core_dequantize(&dequant_delta, &out[i], step_size, scale_factor, 0, 1);
203 
204  work_bufer[DCA_ADPCM_COEFFS+i] += dequant_delta;
205  }
206 
207  memcpy(next_hist, &work_bufer[len], sizeof(int32_t) * DCA_ADPCM_COEFFS);
208 
209  return 0;
210 }
211 
213 {
214  if (!s)
215  return -1;
216 
218  if (!s->private_data)
219  return AVERROR(ENOMEM);
220 
221  precalc(s->private_data);
222  return 0;
223 }
224 
226 {
227  if (!s)
228  return;
229 
230  av_freep(&s->private_data);
231 }
#define MUL64(a, b)
Definition: mathops.h:54
av_cold int ff_dcaadpcm_init(DCAADPCMEncContext *s)
Definition: dcaadpcm.c:212
const char * s
Definition: avisynth_c.h:768
int32_t premultiplied_coeffs[10]
Definition: dcaadpcm.c:28
static int64_t ff_dcaadpcm_predict(int pred_vq_index, const int32_t *input)
Definition: dcaadpcm.h:33
static int64_t calc_corr(const int32_t *x, int len, int j, int k)
Definition: dcaadpcm.c:31
int ff_dcaadpcm_do_real(int pred_vq_index, softfloat quant, int32_t scale_factor, int32_t step_size, const int32_t *prev_hist, const int32_t *in, int32_t *next_hist, int32_t *out, int len, int32_t peak)
Definition: dcaadpcm.c:183
static int32_t quantize_value(int32_t value, softfloat quant)
Definition: dcaenc.h:149
#define av_cold
Definition: attributes.h:82
#define av_malloc(s)
float delta
const char data[16]
Definition: mxf.c:90
static void ff_dca_core_dequantize(int32_t *output, const int32_t *input, int32_t step_size, int32_t scale, int residual, int len)
Definition: dca_core.h:227
#define DCA_ADPCM_VQCODEBOOK_SZ
Definition: dcadata.h:29
#define AVERROR(e)
Definition: error.h:43
#define DCA_ADPCM_COEFFS
Definition: dcadata.h:28
void * private_data
Definition: dcaadpcm.h:30
int32_t
#define FFABS(a)
Absolute value, Note, INT_MIN / INT64_MIN result in undefined behavior as they are not representable ...
Definition: common.h:72
int n
Definition: avisynth_c.h:684
static void error(const char *err)
#define av_log2
Definition: intmath.h:83
static int64_t apply_filter(const int16_t a[DCA_ADPCM_COEFFS], const int64_t corr[15], const int32_t aa[10])
Definition: dcaadpcm.c:40
uint8_t pi<< 24) CONV_FUNC_GROUP(AV_SAMPLE_FMT_FLT, float, AV_SAMPLE_FMT_U8, uint8_t,(*(const uint8_t *) pi - 0x80) *(1.0f/(1<< 7))) CONV_FUNC_GROUP(AV_SAMPLE_FMT_DBL, double, AV_SAMPLE_FMT_U8, uint8_t,(*(const uint8_t *) pi - 0x80) *(1.0/(1<< 7))) CONV_FUNC_GROUP(AV_SAMPLE_FMT_U8, uint8_t, AV_SAMPLE_FMT_S16, int16_t,(*(const int16_t *) pi >> 8)+0x80) CONV_FUNC_GROUP(AV_SAMPLE_FMT_FLT, float, AV_SAMPLE_FMT_S16, int16_t, *(const int16_t *) pi *(1.0f/(1<< 15))) CONV_FUNC_GROUP(AV_SAMPLE_FMT_DBL, double, AV_SAMPLE_FMT_S16, int16_t, *(const int16_t *) pi *(1.0/(1<< 15))) CONV_FUNC_GROUP(AV_SAMPLE_FMT_U8, uint8_t, AV_SAMPLE_FMT_S32, int32_t,(*(const int32_t *) pi >> 24)+0x80) CONV_FUNC_GROUP(AV_SAMPLE_FMT_FLT, float, AV_SAMPLE_FMT_S32, int32_t, *(const int32_t *) pi *(1.0f/(1U<< 31))) CONV_FUNC_GROUP(AV_SAMPLE_FMT_DBL, double, AV_SAMPLE_FMT_S32, int32_t, *(const int32_t *) pi *(1.0/(1U<< 31))) CONV_FUNC_GROUP(AV_SAMPLE_FMT_U8, uint8_t, AV_SAMPLE_FMT_FLT, float, av_clip_uint8(lrintf(*(const float *) pi *(1<< 7))+0x80)) CONV_FUNC_GROUP(AV_SAMPLE_FMT_S16, int16_t, AV_SAMPLE_FMT_FLT, float, av_clip_int16(lrintf(*(const float *) pi *(1<< 15)))) CONV_FUNC_GROUP(AV_SAMPLE_FMT_S32, int32_t, AV_SAMPLE_FMT_FLT, float, av_clipl_int32(llrintf(*(const float *) pi *(1U<< 31)))) CONV_FUNC_GROUP(AV_SAMPLE_FMT_U8, uint8_t, AV_SAMPLE_FMT_DBL, double, av_clip_uint8(lrint(*(const double *) pi *(1<< 7))+0x80)) CONV_FUNC_GROUP(AV_SAMPLE_FMT_S16, int16_t, AV_SAMPLE_FMT_DBL, double, av_clip_int16(lrint(*(const double *) pi *(1<< 15)))) CONV_FUNC_GROUP(AV_SAMPLE_FMT_S32, int32_t, AV_SAMPLE_FMT_DBL, double, av_clipl_int32(llrint(*(const double *) pi *(1U<< 31)))) #define SET_CONV_FUNC_GROUP(ofmt, ifmt) static void set_generic_function(AudioConvert *ac) { } void ff_audio_convert_free(AudioConvert **ac) { if(! *ac) return;ff_dither_free(&(*ac) ->dc);av_freep(ac);} AudioConvert *ff_audio_convert_alloc(AVAudioResampleContext *avr, enum AVSampleFormat out_fmt, enum AVSampleFormat in_fmt, int channels, int sample_rate, int apply_map) { AudioConvert *ac;int in_planar, out_planar;ac=av_mallocz(sizeof(*ac));if(!ac) return NULL;ac->avr=avr;ac->out_fmt=out_fmt;ac->in_fmt=in_fmt;ac->channels=channels;ac->apply_map=apply_map;if(avr->dither_method !=AV_RESAMPLE_DITHER_NONE &&av_get_packed_sample_fmt(out_fmt)==AV_SAMPLE_FMT_S16 &&av_get_bytes_per_sample(in_fmt) > 2) { ac->dc=ff_dither_alloc(avr, out_fmt, in_fmt, channels, sample_rate, apply_map);if(!ac->dc) { av_free(ac);return NULL;} return ac;} in_planar=ff_sample_fmt_is_planar(in_fmt, channels);out_planar=ff_sample_fmt_is_planar(out_fmt, channels);if(in_planar==out_planar) { ac->func_type=CONV_FUNC_TYPE_FLAT;ac->planes=in_planar ? ac->channels :1;} else if(in_planar) ac->func_type=CONV_FUNC_TYPE_INTERLEAVE;else ac->func_type=CONV_FUNC_TYPE_DEINTERLEAVE;set_generic_function(ac);if(ARCH_AARCH64) ff_audio_convert_init_aarch64(ac);if(ARCH_ARM) ff_audio_convert_init_arm(ac);if(ARCH_X86) ff_audio_convert_init_x86(ac);return ac;} int ff_audio_convert(AudioConvert *ac, AudioData *out, AudioData *in) { int use_generic=1;int len=in->nb_samples;int p;if(ac->dc) { av_log(ac->avr, AV_LOG_TRACE, "%d samples - audio_convert: %s to %s (dithered)\", len, av_get_sample_fmt_name(ac->in_fmt), av_get_sample_fmt_name(ac->out_fmt));return ff_convert_dither(ac-> in
static void precalc(premultiplied_coeffs *data)
Definition: dcaadpcm.c:164
const int16_t ff_dca_adpcm_vb[DCA_ADPCM_VQCODEBOOK_SZ][DCA_ADPCM_COEFFS]
Definition: dcadata.c:60
int ff_dcaadpcm_subband_analysis(const DCAADPCMEncContext *s, const int32_t *in, int len, int *diff)
Definition: dcaadpcm.c:125
const uint8_t * quant
static int64_t calc_prediction_gain(int pred_vq, const int32_t *in, int32_t *out, int len)
Definition: dcaadpcm.c:104
static av_always_inline int diff(const uint32_t a, const uint32_t b)
static int64_t find_best_filter(const DCAADPCMEncContext *s, const int32_t *in, int len)
Definition: dcaadpcm.c:79
int len
static int32_t norm__(int64_t a, int bits)
Definition: dcamath.h:27
FILE * out
Definition: movenc.c:54
#define av_freep(p)
av_cold void ff_dcaadpcm_free(DCAADPCMEncContext *s)
Definition: dcaadpcm.c:225
static uint8_t tmp[11]
Definition: aes_ctr.c:26