FFmpeg
afir_template.c
Go to the documentation of this file.
1 /*
2  * Copyright (c) 2017 Paul B Mahol
3  *
4  * This file is part of FFmpeg.
5  *
6  * FFmpeg is free software; you can redistribute it and/or
7  * modify it under the terms of the GNU Lesser General Public
8  * License as published by the Free Software Foundation; either
9  * version 2.1 of the License, or (at your option) any later version.
10  *
11  * FFmpeg is distributed in the hope that it will be useful,
12  * but WITHOUT ANY WARRANTY; without even the implied warranty of
13  * MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the GNU
14  * Lesser General Public License for more details.
15  *
16  * You should have received a copy of the GNU Lesser General Public
17  * License along with FFmpeg; if not, write to the Free Software
18  * Foundation, Inc., 51 Franklin Street, Fifth Floor, Boston, MA 02110-1301 USA
19  */
20 
21 #include "libavutil/tx.h"
22 #include "avfilter.h"
23 #include "formats.h"
24 #include "internal.h"
25 #include "audio.h"
26 
27 #undef ctype
28 #undef ftype
29 #undef SQRT
30 #undef HYPOT
31 #undef SAMPLE_FORMAT
32 #undef TX_TYPE
33 #if DEPTH == 32
34 #define SAMPLE_FORMAT float
35 #define SQRT sqrtf
36 #define HYPOT hypotf
37 #define ctype AVComplexFloat
38 #define ftype float
39 #define TX_TYPE AV_TX_FLOAT_RDFT
40 #else
41 #define SAMPLE_FORMAT double
42 #define SQRT sqrt
43 #define HYPOT hypot
44 #define ctype AVComplexDouble
45 #define ftype double
46 #define TX_TYPE AV_TX_DOUBLE_RDFT
47 #endif
48 
49 #define fn3(a,b) a##_##b
50 #define fn2(a,b) fn3(a,b)
51 #define fn(a) fn2(a, SAMPLE_FORMAT)
52 
53 static void fn(draw_response)(AVFilterContext *ctx, AVFrame *out)
54 {
55  AudioFIRContext *s = ctx->priv;
56  ftype *mag, *phase, *delay, min = FLT_MAX, max = FLT_MIN;
57  ftype min_delay = FLT_MAX, max_delay = FLT_MIN;
58  int prev_ymag = -1, prev_yphase = -1, prev_ydelay = -1;
59  char text[32];
60  int channel, i, x;
61 
62  for (int y = 0; y < s->h; y++)
63  memset(out->data[0] + y * out->linesize[0], 0, s->w * 4);
64 
65  phase = av_malloc_array(s->w, sizeof(*phase));
66  mag = av_malloc_array(s->w, sizeof(*mag));
67  delay = av_malloc_array(s->w, sizeof(*delay));
68  if (!mag || !phase || !delay)
69  goto end;
70 
71  channel = av_clip(s->ir_channel, 0, s->ir[s->selir]->ch_layout.nb_channels - 1);
72  for (i = 0; i < s->w; i++) {
73  const ftype *src = (const ftype *)s->ir[s->selir]->extended_data[channel];
74  double w = i * M_PI / (s->w - 1);
75  double div, real_num = 0., imag_num = 0., real = 0., imag = 0.;
76 
77  for (x = 0; x < s->nb_taps[s->selir]; x++) {
78  real += cos(-x * w) * src[x];
79  imag += sin(-x * w) * src[x];
80  real_num += cos(-x * w) * src[x] * x;
81  imag_num += sin(-x * w) * src[x] * x;
82  }
83 
84  mag[i] = hypot(real, imag);
85  phase[i] = atan2(imag, real);
86  div = real * real + imag * imag;
87  delay[i] = (real_num * real + imag_num * imag) / div;
88  min = fminf(min, mag[i]);
89  max = fmaxf(max, mag[i]);
90  min_delay = fminf(min_delay, delay[i]);
91  max_delay = fmaxf(max_delay, delay[i]);
92  }
93 
94  for (i = 0; i < s->w; i++) {
95  int ymag = mag[i] / max * (s->h - 1);
96  int ydelay = (delay[i] - min_delay) / (max_delay - min_delay) * (s->h - 1);
97  int yphase = (0.5 * (1. + phase[i] / M_PI)) * (s->h - 1);
98 
99  ymag = s->h - 1 - av_clip(ymag, 0, s->h - 1);
100  yphase = s->h - 1 - av_clip(yphase, 0, s->h - 1);
101  ydelay = s->h - 1 - av_clip(ydelay, 0, s->h - 1);
102 
103  if (prev_ymag < 0)
104  prev_ymag = ymag;
105  if (prev_yphase < 0)
106  prev_yphase = yphase;
107  if (prev_ydelay < 0)
108  prev_ydelay = ydelay;
109 
110  draw_line(out, i, ymag, FFMAX(i - 1, 0), prev_ymag, 0xFFFF00FF);
111  draw_line(out, i, yphase, FFMAX(i - 1, 0), prev_yphase, 0xFF00FF00);
112  draw_line(out, i, ydelay, FFMAX(i - 1, 0), prev_ydelay, 0xFF00FFFF);
113 
114  prev_ymag = ymag;
115  prev_yphase = yphase;
116  prev_ydelay = ydelay;
117  }
118 
119  if (s->w > 400 && s->h > 100) {
120  drawtext(out, 2, 2, "Max Magnitude:", 0xDDDDDDDD);
121  snprintf(text, sizeof(text), "%.2f", max);
122  drawtext(out, 15 * 8 + 2, 2, text, 0xDDDDDDDD);
123 
124  drawtext(out, 2, 12, "Min Magnitude:", 0xDDDDDDDD);
125  snprintf(text, sizeof(text), "%.2f", min);
126  drawtext(out, 15 * 8 + 2, 12, text, 0xDDDDDDDD);
127 
128  drawtext(out, 2, 22, "Max Delay:", 0xDDDDDDDD);
129  snprintf(text, sizeof(text), "%.2f", max_delay);
130  drawtext(out, 11 * 8 + 2, 22, text, 0xDDDDDDDD);
131 
132  drawtext(out, 2, 32, "Min Delay:", 0xDDDDDDDD);
133  snprintf(text, sizeof(text), "%.2f", min_delay);
134  drawtext(out, 11 * 8 + 2, 32, text, 0xDDDDDDDD);
135  }
136 
137 end:
138  av_free(delay);
139  av_free(phase);
140  av_free(mag);
141 }
142 
143 static int fn(get_power)(AVFilterContext *ctx, AudioFIRContext *s,
144  int cur_nb_taps, int ch,
145  ftype *time)
146 {
147  ftype ch_gain = 1;
148 
149  switch (s->gtype) {
150  case -1:
151  ch_gain = 1;
152  break;
153  case 0:
154  {
155  ftype sum = 0;
156 
157  for (int i = 0; i < cur_nb_taps; i++)
158  sum += FFABS(time[i]);
159  ch_gain = 1. / sum;
160  }
161  break;
162  case 1:
163  {
164  ftype sum = 0;
165 
166  for (int i = 0; i < cur_nb_taps; i++)
167  sum += time[i];
168  ch_gain = 1. / sum;
169  }
170  break;
171  case 2:
172  {
173  ftype sum = 0;
174 
175  for (int i = 0; i < cur_nb_taps; i++)
176  sum += time[i] * time[i];
177  ch_gain = 1. / SQRT(sum);
178  }
179  break;
180  case 3:
181  case 4:
182  {
183  ftype *inc, *outc, scale, power;
184  AVTXContext *tx;
185  av_tx_fn tx_fn;
186  int ret, size;
187 
188  size = 1 << av_ceil_log2_c(cur_nb_taps);
189  inc = av_calloc(size + 2, sizeof(SAMPLE_FORMAT));
190  outc = av_calloc(size + 2, sizeof(SAMPLE_FORMAT));
191  if (!inc || !outc) {
192  av_free(outc);
193  av_free(inc);
194  break;
195  }
196 
197  scale = 1.;
198  ret = av_tx_init(&tx, &tx_fn, TX_TYPE, 0, size, &scale, 0);
199  if (ret < 0) {
200  av_free(outc);
201  av_free(inc);
202  break;
203  }
204 
205  {
206  memcpy(inc, time, cur_nb_taps * sizeof(SAMPLE_FORMAT));
207  tx_fn(tx, outc, inc, sizeof(SAMPLE_FORMAT));
208 
209  power = 0;
210  if (s->gtype == 3) {
211  for (int i = 0; i < size / 2 + 1; i++)
212  power = FFMAX(power, HYPOT(outc[i * 2], outc[i * 2 + 1]));
213  } else {
214  ftype sum = 0;
215  for (int i = 0; i < size / 2 + 1; i++)
216  sum += HYPOT(outc[i * 2], outc[i * 2 + 1]);
217  power = SQRT(sum / (size / 2 + 1));
218  }
219 
220  ch_gain = 1. / power;
221  }
222 
223  av_tx_uninit(&tx);
224  av_free(outc);
225  av_free(inc);
226  }
227  break;
228  default:
229  return AVERROR_BUG;
230  }
231 
232  if (ch_gain != 1. || s->ir_gain != 1.) {
233  ftype gain = ch_gain * s->ir_gain;
234 
235  av_log(ctx, AV_LOG_DEBUG, "ch%d gain %f\n", ch, gain);
236 #if DEPTH == 32
237  s->fdsp->vector_fmul_scalar(time, time, gain, FFALIGN(cur_nb_taps, 4));
238 #else
239  s->fdsp->vector_dmul_scalar(time, time, gain, FFALIGN(cur_nb_taps, 8));
240 #endif
241  }
242 
243  return 0;
244 }
245 
246 static void fn(convert_channel)(AVFilterContext *ctx, AudioFIRContext *s, int ch,
247  AudioFIRSegment *seg, int coeff_partition, int selir)
248 {
249  const int coffset = coeff_partition * seg->coeff_size;
250  const int nb_taps = s->nb_taps[selir];
251  ftype *time = (ftype *)s->norm_ir[selir]->extended_data[ch];
252  ftype *tempin = (ftype *)seg->tempin->extended_data[ch];
253  ftype *tempout = (ftype *)seg->tempout->extended_data[ch];
254  ctype *coeff = (ctype *)seg->coeff->extended_data[ch];
255  const int remaining = nb_taps - (seg->input_offset + coeff_partition * seg->part_size);
256  const int size = remaining >= seg->part_size ? seg->part_size : remaining;
257 
258  memset(tempin + size, 0, sizeof(*tempin) * (seg->block_size - size));
259  memcpy(tempin, time + seg->input_offset + coeff_partition * seg->part_size,
260  size * sizeof(*tempin));
261  seg->ctx_fn(seg->ctx[ch], tempout, tempin, sizeof(*tempin));
262  memcpy(coeff + coffset, tempout, seg->coeff_size * sizeof(*coeff));
263 
264  av_log(ctx, AV_LOG_DEBUG, "channel: %d\n", ch);
265  av_log(ctx, AV_LOG_DEBUG, "nb_partitions: %d\n", seg->nb_partitions);
266  av_log(ctx, AV_LOG_DEBUG, "partition size: %d\n", seg->part_size);
267  av_log(ctx, AV_LOG_DEBUG, "block size: %d\n", seg->block_size);
268  av_log(ctx, AV_LOG_DEBUG, "fft_length: %d\n", seg->fft_length);
269  av_log(ctx, AV_LOG_DEBUG, "coeff_size: %d\n", seg->coeff_size);
270  av_log(ctx, AV_LOG_DEBUG, "input_size: %d\n", seg->input_size);
271  av_log(ctx, AV_LOG_DEBUG, "input_offset: %d\n", seg->input_offset);
272 }
273 
274 static void fn(fir_fadd)(AudioFIRContext *s, ftype *dst, const ftype *src, int nb_samples)
275 {
276  if ((nb_samples & 15) == 0 && nb_samples >= 8) {
277 #if DEPTH == 32
278  s->fdsp->vector_fmac_scalar(dst, src, 1.f, nb_samples);
279 #else
280  s->fdsp->vector_dmac_scalar(dst, src, 1.0, nb_samples);
281 #endif
282  } else {
283  for (int n = 0; n < nb_samples; n++)
284  dst[n] += src[n];
285  }
286 }
287 
288 static int fn(fir_quantum)(AVFilterContext *ctx, AVFrame *out, int ch, int ioffset, int offset, int selir)
289 {
290  AudioFIRContext *s = ctx->priv;
291  const ftype *in = (const ftype *)s->in->extended_data[ch] + ioffset;
292  ftype *blockout, *ptr = (ftype *)out->extended_data[ch] + offset;
293  const int min_part_size = s->min_part_size;
294  const int nb_samples = FFMIN(min_part_size, out->nb_samples - offset);
295  const int nb_segments = s->nb_segments[selir];
296  const float dry_gain = s->dry_gain;
297  const float wet_gain = s->wet_gain;
298 
299  for (int segment = 0; segment < nb_segments; segment++) {
300  AudioFIRSegment *seg = &s->seg[selir][segment];
301  ftype *src = (ftype *)seg->input->extended_data[ch];
302  ftype *dst = (ftype *)seg->output->extended_data[ch];
303  ftype *sumin = (ftype *)seg->sumin->extended_data[ch];
304  ftype *sumout = (ftype *)seg->sumout->extended_data[ch];
305  ftype *tempin = (ftype *)seg->tempin->extended_data[ch];
306  ftype *buf = (ftype *)seg->buffer->extended_data[ch];
307  int *output_offset = &seg->output_offset[ch];
308  const int nb_partitions = seg->nb_partitions;
309  const int input_offset = seg->input_offset;
310  const int part_size = seg->part_size;
311  int j;
312 
313  seg->part_index[ch] = seg->part_index[ch] % nb_partitions;
314  if (dry_gain == 1.f) {
315  memcpy(src + input_offset, in, nb_samples * sizeof(*src));
316  } else if (min_part_size >= 8) {
317 #if DEPTH == 32
318  s->fdsp->vector_fmul_scalar(src + input_offset, in, dry_gain, FFALIGN(nb_samples, 4));
319 #else
320  s->fdsp->vector_dmul_scalar(src + input_offset, in, dry_gain, FFALIGN(nb_samples, 8));
321 #endif
322  emms_c();
323  } else {
324  ftype *src2 = src + input_offset;
325  for (int n = 0; n < nb_samples; n++)
326  src2[n] = in[n] * dry_gain;
327  }
328 
329  output_offset[0] += min_part_size;
330  if (output_offset[0] >= part_size) {
331  output_offset[0] = 0;
332  } else {
333  memmove(src, src + min_part_size, (seg->input_size - min_part_size) * sizeof(*src));
334 
335  dst += output_offset[0];
336  fn(fir_fadd)(s, ptr, dst, nb_samples);
337  continue;
338  }
339 
340  memset(sumin, 0, sizeof(*sumin) * seg->fft_length);
341 
342  blockout = (ftype *)seg->blockout->extended_data[ch] + seg->part_index[ch] * seg->block_size;
343  memset(tempin + part_size, 0, sizeof(*tempin) * (seg->block_size - part_size));
344  memcpy(tempin, src, sizeof(*src) * part_size);
345  seg->tx_fn(seg->tx[ch], blockout, tempin, sizeof(ftype));
346 
347  j = seg->part_index[ch];
348  for (int i = 0; i < nb_partitions; i++) {
349  const int input_partition = j;
350  const int coeff_partition = i;
351  const int coffset = coeff_partition * seg->coeff_size;
352  const ftype *blockout = (const ftype *)seg->blockout->extended_data[ch] + input_partition * seg->block_size;
353  const ctype *coeff = ((const ctype *)seg->coeff->extended_data[ch]) + coffset;
354 
355  if (j == 0)
356  j = nb_partitions;
357  j--;
358 
359 #if DEPTH == 32
360  s->afirdsp.fcmul_add(sumin, blockout, (const ftype *)coeff, part_size);
361 #else
362  s->afirdsp.dcmul_add(sumin, blockout, (const ftype *)coeff, part_size);
363 #endif
364  }
365 
366  seg->itx_fn(seg->itx[ch], sumout, sumin, sizeof(ctype));
367 
368  fn(fir_fadd)(s, buf, sumout, part_size);
369  memcpy(dst, buf, part_size * sizeof(*dst));
370  memcpy(buf, sumout + part_size, part_size * sizeof(*buf));
371 
372  fn(fir_fadd)(s, ptr, dst, nb_samples);
373 
374  if (part_size != min_part_size)
375  memmove(src, src + min_part_size, (seg->input_size - min_part_size) * sizeof(*src));
376 
377  seg->part_index[ch] = (seg->part_index[ch] + 1) % nb_partitions;
378  }
379 
380  if (wet_gain == 1.f)
381  return 0;
382 
383  if (min_part_size >= 8) {
384 #if DEPTH == 32
385  s->fdsp->vector_fmul_scalar(ptr, ptr, wet_gain, FFALIGN(nb_samples, 4));
386 #else
387  s->fdsp->vector_dmul_scalar(ptr, ptr, wet_gain, FFALIGN(nb_samples, 8));
388 #endif
389  emms_c();
390  } else {
391  for (int n = 0; n < nb_samples; n++)
392  ptr[n] *= wet_gain;
393  }
394 
395  return 0;
396 }
av_clip
#define av_clip
Definition: common.h:95
draw_response
static void fn() draw_response(AVFilterContext *ctx, AVFrame *out)
Definition: afir_template.c:53
convert_channel
static void fn() convert_channel(AVFilterContext *ctx, AudioFIRContext *s, int ch, AudioFIRSegment *seg, int coeff_partition, int selir)
Definition: afir_template.c:246
AudioFIRSegment::block_size
int block_size
Definition: af_afir.h:36
out
FILE * out
Definition: movenc.c:54
ctype
#define ctype
Definition: afir_template.c:44
TX_TYPE
#define TX_TYPE
Definition: afir_template.c:46
AVTXContext
Definition: tx_priv.h:234
SQRT
#define SQRT
Definition: afir_template.c:42
AVFrame
This structure describes decoded (raw) audio or video data.
Definition: frame.h:330
AudioFIRSegment::buffer
AVFrame * buffer
Definition: af_afir.h:50
w
uint8_t w
Definition: llviddspenc.c:38
HYPOT
#define HYPOT
Definition: afir_template.c:43
AudioFIRSegment::input_offset
int input_offset
Definition: af_afir.h:40
AudioFIRSegment::tx_fn
av_tx_fn tx_fn
Definition: af_afir.h:56
max
#define max(a, b)
Definition: cuda_runtime.h:33
FFMAX
#define FFMAX(a, b)
Definition: macros.h:47
AudioFIRSegment::part_size
int part_size
Definition: af_afir.h:35
AudioFIRSegment::input_size
int input_size
Definition: af_afir.h:39
av_tx_init
av_cold int av_tx_init(AVTXContext **ctx, av_tx_fn *tx, enum AVTXType type, int inv, int len, const void *scale, uint64_t flags)
Initialize a transform context with the given configuration (i)MDCTs with an odd length are currently...
Definition: tx.c:883
formats.h
AudioFIRSegment::coeff
AVFrame * coeff
Definition: af_afir.h:51
fn
#define fn(a)
Definition: afir_template.c:51
scale
static av_always_inline float scale(float x, float s)
Definition: vf_v360.c:1389
AudioFIRSegment::blockout
AVFrame * blockout
Definition: af_afir.h:47
AudioFIRSegment
Definition: af_afir.h:33
AudioFIRSegment::tx
AVTXContext ** tx
Definition: af_afir.h:55
ftype
#define ftype
Definition: afir_template.c:45
av_tx_fn
void(* av_tx_fn)(AVTXContext *s, void *out, void *in, ptrdiff_t stride)
Function pointer to a function to perform the transform.
Definition: tx.h:127
emms_c
#define emms_c()
Definition: internal.h:51
s
#define s(width, name)
Definition: cbs_vp9.c:256
fminf
float fminf(float, float)
AV_LOG_DEBUG
#define AV_LOG_DEBUG
Stuff which is only useful for libav* developers.
Definition: log.h:201
ctx
AVFormatContext * ctx
Definition: movenc.c:48
FFABS
#define FFABS(a)
Absolute value, Note, INT_MIN / INT64_MIN result in undefined behavior as they are not representable ...
Definition: common.h:64
draw_line
static void draw_line(AVFrame *out, int x0, int y0, int x1, int y1, uint32_t color)
Definition: af_afir.c:72
AudioFIRSegment::itx_fn
av_tx_fn itx_fn
Definition: af_afir.h:56
AudioFIRSegment::output
AVFrame * output
Definition: af_afir.h:53
f
f
Definition: af_crystalizer.c:122
fmaxf
float fmaxf(float, float)
hypot
static av_const double hypot(double x, double y)
Definition: libm.h:366
size
int size
Definition: twinvq_data.h:10344
AudioFIRSegment::sumin
AVFrame * sumin
Definition: af_afir.h:45
fir_fadd
static void fn() fir_fadd(AudioFIRContext *s, ftype *dst, const ftype *src, int nb_samples)
Definition: afir_template.c:274
offset
it s the only field you need to keep assuming you have a context There is some magic you don t need to care about around this just let it vf offset
Definition: writing_filters.txt:86
AudioFIRSegment::tempin
AVFrame * tempin
Definition: af_afir.h:48
M_PI
#define M_PI
Definition: mathematics.h:67
av_tx_uninit
av_cold void av_tx_uninit(AVTXContext **ctx)
Frees a context and sets *ctx to NULL, does nothing when *ctx == NULL.
Definition: tx.c:294
internal.h
i
#define i(width, name, range_min, range_max)
Definition: cbs_h2645.c:269
AudioFIRSegment::input
AVFrame * input
Definition: af_afir.h:52
get_power
static int fn() get_power(AVFilterContext *ctx, AudioFIRContext *s, int cur_nb_taps, int ch, ftype *time)
Definition: afir_template.c:143
AudioFIRSegment::coeff_size
int coeff_size
Definition: af_afir.h:38
AVFrame::extended_data
uint8_t ** extended_data
pointers to the data planes/channels.
Definition: frame.h:391
src2
const pixel * src2
Definition: h264pred_template.c:422
av_malloc_array
#define av_malloc_array(a, b)
Definition: tableprint_vlc.h:31
AudioFIRSegment::nb_partitions
int nb_partitions
Definition: af_afir.h:34
FFMIN
#define FFMIN(a, b)
Definition: macros.h:49
av_calloc
void * av_calloc(size_t nmemb, size_t size)
Definition: mem.c:262
ret
ret
Definition: filter_design.txt:187
AudioFIRSegment::itx
AVTXContext ** itx
Definition: af_afir.h:55
AudioFIRSegment::fft_length
int fft_length
Definition: af_afir.h:37
power
static float power(float r, float g, float b, float max)
Definition: preserve_color.h:45
AudioFIRSegment::sumout
AVFrame * sumout
Definition: af_afir.h:46
AudioFIRContext
Definition: af_afir.h:59
avfilter.h
segment
Definition: hls.c:75
AVFilterContext
An instance of a filter.
Definition: avfilter.h:397
audio.h
av_free
#define av_free(p)
Definition: tableprint_vlc.h:33
FFALIGN
#define FFALIGN(x, a)
Definition: macros.h:78
src
INIT_CLIP pixel * src
Definition: h264pred_template.c:418
AVERROR_BUG
#define AVERROR_BUG
Internal bug, also see AVERROR_BUG2.
Definition: error.h:52
coeff
static const double coeff[2][5]
Definition: vf_owdenoise.c:78
av_log
#define av_log(a,...)
Definition: tableprint_vlc.h:27
snprintf
#define snprintf
Definition: snprintf.h:34
AudioFIRSegment::output_offset
int * output_offset
Definition: af_afir.h:42
channel
channel
Definition: ebur128.h:39
av_ceil_log2_c
static av_always_inline av_const int av_ceil_log2_c(int x)
Compute ceil(log2(x)).
Definition: common.h:417
drawtext
static void drawtext(AVFrame *pic, int x, int y, const char *txt, uint32_t color)
Definition: af_afir.c:49
tx.h
SAMPLE_FORMAT
#define SAMPLE_FORMAT
Definition: afir_template.c:41
fir_quantum
static int fn() fir_quantum(AVFilterContext *ctx, AVFrame *out, int ch, int ioffset, int offset, int selir)
Definition: afir_template.c:288
min
float min
Definition: vorbis_enc_data.h:429
AudioFIRSegment::part_index
int * part_index
Definition: af_afir.h:43
Generated on Fri Jun 9 2023 19:22:05 for FFmpeg by   doxygen 1.8.17