FFmpeg coverage

Directory: ../../../ffmpeg/
File: src/libswscale/filters.c
Date: 2026-09-24 20:08:24
Exec Total Coverage
Lines: 154 229 67.2%
Functions: 10 24 41.7%
Branches: 69 99 69.7%
Line Branch Exec Source
1 /*
2 * Copyright (C) 2026 Niklas Haas
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 <math.h>
22 #include <stdbool.h>
23
24 #include <libavutil/attributes.h>
25 #include <libavutil/avassert.h>
26 #include <libavutil/mem.h>
27
28 #include "filters.h"
29
30 #ifdef _WIN32
31 # define j1 _j1
32 #endif
33
34 /* Maximum (pre-stretching) radius (for tunable filters) */
35 #define RADIUS_MAX 10.0
36
37 /* Defined only on [0, radius]. */
38 typedef double (*SwsFilterKernel)(double x, const double *params);
39
40 typedef struct SwsFilterFunction {
41 char name[16];
42 double radius; /* negative means resizable */
43 SwsFilterKernel kernel;
44 SwsFilterKernel window; /* optional */
45 double params[SWS_NUM_SCALER_PARAMS]; /* default params */
46 } SwsFilterFunction;
47
48 static const SwsFilterFunction filter_functions[SWS_SCALE_NB];
49
50 73579320 static double scaler_sample(const SwsFilterFunction *f, double x)
51 {
52 73579320 x = fabs(x);
53
2/2
✓ Branch 0 taken 1930179 times.
✓ Branch 1 taken 71649141 times.
 73579320 if (x > f->radius)
54 1930179 return 0.0;
55
56 71649141 double w = f->kernel(x, f->params);
57
1/2
✗ Branch 0 not taken.
✓ Branch 1 taken 71649141 times.
 71649141 if (f->window)
58 ✗ w *= f->window(x / f->radius, f->params);
59 71649141 return w;
60 }
61
62 21797760 static void compute_row(SwsFilterWeights *f, const SwsFilterFunction *fun,
63 double radius, double ratio_inv, double stretch_inv,
64 int dst_pos, double *tmp)
65 {
66 21797760 int *out = &f->weights[dst_pos * f->filter_size];
67 21797760 int *pos = &f->offsets[dst_pos];
68
69 /**
70 * Explanation of the 0.5 offsets: Normally, pixel samples are assumed
71 * to be representative of the center of their containing area; e.g. for
72 * a 2x2 image, the samples are located at {0.5, 1.5}^2. However, with
73 * integer indexing, we round sample positions down (0-based indexing).
74 * So the (0, 0) sample is actually located at (0.5, 0.5) and represents
75 * the entire square from (0,0) to (1,1). When normalizing between different
76 * image sizes, we therefore need to add/subtract off these 0.5 offsets.
77 */
78 21797760 const double src_pos = (dst_pos + 0.5) * ratio_inv - 0.5 + f->offset;
79
2/2
✓ Branch 0 taken 524032 times.
✓ Branch 1 taken 21273728 times.
 21797760 if (f->filter_size == 1) {
80 524032 *pos = fmin(fmax(round(src_pos), 0.0), f->src_size - 1);
81 524032 *out = SWS_FILTER_SCALE;
82 524032 return;
83 }
84
85 /* First pixel that is actually within the filter envelope */
86 21273728 const double start_pos = src_pos - radius;
87 21273728 int64_t start_idx = ceil(start_pos);
88 21273728 start_idx = FFMAX(start_idx, 0); /* edge clamping */
89 21273728 start_idx = FFMIN(start_idx, f->src_size - f->filter_size);
90 21273728 const double offset = start_idx - src_pos;
91 21273728 *pos = start_idx;
92
93 /**
94 * Generate raw filter weights with maximum precision. Sum the positive
95 * and negative weights separately to avoid catastrophic cancellation. This
96 * summation order should already give the best precision because abs(w)
97 * is monotonically decreasing
98 */
99 21273728 const double base = stretch_inv * offset;
100 21273728 double wsum_pos = 0.0, wsum_neg = 0.0;
101
2/2
✓ Branch 0 taken 53372704 times.
✓ Branch 1 taken 21273728 times.
 74646432 for (int i = 0; i < f->filter_size; i++) {
102 53372704 tmp[i] = scaler_sample(fun, base + stretch_inv * i);
103
2/2
✓ Branch 0 taken 44404140 times.
✓ Branch 1 taken 8968564 times.
 53372704 if (tmp[i] >= 0)
104 44404140 wsum_pos += tmp[i];
105 else
106 8968564 wsum_neg += tmp[i];
107 }
108
109 21273728 const double wsum = wsum_pos + wsum_neg;
110
1/2
✗ Branch 0 not taken.
✓ Branch 1 taken 21273728 times.
 21273728 av_assert0(wsum > 0);
111
112 /* Generate correctly rounded filter weights with error diffusion */
113 21273728 double error = 0.0;
114 21273728 int sum_pos = 0, sum_neg = 0;
115
2/2
✓ Branch 0 taken 53372704 times.
✓ Branch 1 taken 21273728 times.
 74646432 for (int i = 0; i < f->filter_size; i++) {
116
2/2
✓ Branch 0 taken 21273728 times.
✓ Branch 1 taken 32098976 times.
 53372704 if (i == f->filter_size - 1) {
117 /* Ensure weights sum to exactly SWS_FILTER_SCALE */
118 21273728 out[i] = SWS_FILTER_SCALE - sum_pos - sum_neg;
119 } else {
120 32098976 const double w = tmp[i] / wsum + error;
121 32098976 out[i] = round(w * SWS_FILTER_SCALE);
122 32098976 error = w - (double) out[i] / SWS_FILTER_SCALE;
123 }
124
2/2
✓ Branch 0 taken 44510988 times.
✓ Branch 1 taken 8861716 times.
 53372704 if (out[i] >= 0)
125 44510988 sum_pos += out[i];
126 else
127 8861716 sum_neg += out[i];
128 }
129
130
2/2
✓ Branch 0 taken 711548 times.
✓ Branch 1 taken 20562180 times.
 21273728 if (sum_pos > f->sum_positive)
131 711548 f->sum_positive = sum_pos;
132
2/2
✓ Branch 0 taken 206805 times.
✓ Branch 1 taken 21066923 times.
 21273728 if (sum_neg < f->sum_negative)
133 206805 f->sum_negative = sum_neg;
134 }
135
136 603762 static void sws_filter_free(AVRefStructOpaque opaque, void *obj)
137 {
138 603762 SwsFilterWeights *filter = obj;
139 603762 av_refstruct_unref(&filter->weights);
140 603762 av_refstruct_unref(&filter->offsets);
141 603762 }
142
143 603762 static bool validate_params(const SwsFilterFunction *fun, SwsScaler scaler)
144 {
145
2/4
✗ Branch 0 not taken.
✗ Branch 1 not taken.
✓ Branch 2 taken 90943 times.
✓ Branch 3 taken 512819 times.
 603762 switch (scaler) {
146 ✗ case SWS_SCALE_GAUSSIAN:
147 ✗ return fun->params[0] >= 0.0; /* sigma */
148 ✗ case SWS_SCALE_LANCZOS:
149 ✗ return fun->params[0] >= 1.0 && fun->params[0] <= RADIUS_MAX; /* radius */
150 90943 case SWS_SCALE_BICUBIC:
151 90943 return fun->params[0] < 3.0; /* B param (division by zero) */
152 512819 default:
153 512819 return true;
154 }
155 }
156
157 /**
158 * Numerically estimate the last intersection between the function value
159 * and the cutoff domain [-SWS_MAX_REDUCE_CUTOFF, SWS_MAX_REDUCE_CUTOFF].
160 */
161 92959 static double est_filter_radius(const SwsFilterFunction *fun)
162 {
163 92959 const double bound = fun->radius;
164 92959 const double step = 1e-2;
165
166 92959 double radius = bound;
167 92959 double prev = 0.0, fprev = 1.0; /* f(0) is always 1.0 */
168 92959 double integral = 0.0;
169
2/2
✓ Branch 0 taken 20206616 times.
✓ Branch 1 taken 92959 times.
 20299575 for (double x = step; x < bound + step; x += step) {
170 20206616 const double fx = scaler_sample(fun, x);
171 20206616 integral += (fprev + fx) * step; /* trapezoidal rule (mirrored) */
172 20206616 double cutoff = SWS_MAX_REDUCE_CUTOFF * integral;
173
8/8
✓ Branch 0 taken 10084156 times.
✓ Branch 1 taken 10122460 times.
✓ Branch 2 taken 9983133 times.
✓ Branch 3 taken 101023 times.
✓ Branch 4 taken 9534466 times.
✓ Branch 5 taken 10571127 times.
✓ Branch 6 taken 101023 times.
✓ Branch 7 taken 9433443 times.
 20206616 if ((fprev > cutoff && fx <= cutoff) || (fprev < -cutoff && fx >= -cutoff)) {
174 /* estimate crossing with secant method; note that we have to
175 * bias by the cutoff to find the actual cutoff radius */
176
2/2
✓ Branch 0 taken 101023 times.
✓ Branch 1 taken 101023 times.
 202046 double estimate = fx + (fx > fprev ? cutoff : -cutoff);
177 202046 double root = x - estimate * (x - prev) / (fx - fprev);
178 202046 radius = fmin(root, bound);
179 }
180 20206616 prev = x;
181 20206616 fprev = fx;
182 }
183
184 92959 return radius;
185 }
186
187 603762 int ff_sws_filter_generate(void *log, const SwsFilterParams *params,
188 SwsFilterWeights **out)
189 {
190 603762 SwsScaler scaler = params->scaler;
191
1/2
✗ Branch 0 not taken.
✓ Branch 1 taken 603762 times.
 603762 if (scaler >= SWS_SCALE_NB)
192 ✗ return AVERROR(EINVAL);
193
194
2/2
✓ Branch 0 taken 12340 times.
✓ Branch 1 taken 591422 times.
 603762 if (scaler == SWS_SCALE_AUTO)
195 12340 scaler = SWS_SCALE_BICUBIC;
196
197 603762 double virtual_size = params->virtual_size;
198
1/2
✓ Branch 0 taken 603762 times.
✗ Branch 1 not taken.
 603762 if (!virtual_size)
199 603762 virtual_size = params->dst_size;
200
201 603762 const double ratio = virtual_size / params->src_size;
202 603762 double stretch = 1.0;
203
4/4
✓ Branch 0 taken 18931 times.
✓ Branch 1 taken 584831 times.
✓ Branch 2 taken 16763 times.
✓ Branch 3 taken 2168 times.
 603762 if (ratio < 1.0 && scaler != SWS_SCALE_POINT) {
204 /* Widen filter for downscaling (anti-aliasing) */
205 16763 stretch = 1.0 / ratio;
206 }
207
208
2/2
✓ Branch 0 taken 4226 times.
✓ Branch 1 taken 599536 times.
 603762 if (scaler == SWS_SCALE_AREA) {
209 /**
210 * SWS_SCALE_AREA is a pseudo-filter that is equivalent to bilinear
211 * filtering for upscaling (since bilinear just evenly mixes samples
212 * according to the relative distance), and equivalent to (anti-aliased)
213 * point sampling for downscaling.
214 */
215
2/2
✓ Branch 0 taken 2176 times.
✓ Branch 1 taken 2050 times.
 4226 scaler = ratio >= 1.0 ? SWS_SCALE_BILINEAR : SWS_SCALE_POINT;
216 }
217
218 603762 SwsFilterFunction fun = filter_functions[scaler];
219
1/2
✗ Branch 0 not taken.
✓ Branch 1 taken 603762 times.
 603762 if (!fun.kernel)
220 ✗ return AVERROR(EINVAL);
221
222
2/2
✓ Branch 0 taken 1207524 times.
✓ Branch 1 taken 603762 times.
 1811286 for (int i = 0; i < SWS_NUM_SCALER_PARAMS; i++) {
223
1/2
✗ Branch 0 not taken.
✓ Branch 1 taken 1207524 times.
 1207524 if (params->scaler_params[i] != SWS_PARAM_DEFAULT)
224 ✗ fun.params[i] = params->scaler_params[i];
225 }
226
227
1/2
✗ Branch 1 not taken.
✓ Branch 2 taken 603762 times.
 603762 if (!validate_params(&fun, scaler)) {
228 ✗ av_log(log, AV_LOG_ERROR, "Invalid parameters for scaler %s: {%f, %f}\n",
229 fun.name, fun.params[0], fun.params[1]);
230 ✗ return AVERROR(EINVAL);
231 }
232
233
1/2
✗ Branch 0 not taken.
✓ Branch 1 taken 603762 times.
 603762 if (fun.radius < 0.0) /* tunable width kernels like lanczos */
234 ✗ fun.radius = fun.params[0];
235
236 double radius;
237
3/3
✓ Branch 0 taken 8446 times.
✓ Branch 1 taken 502357 times.
✓ Branch 2 taken 92959 times.
 603762 switch (scaler) {
238 8446 case SWS_SCALE_POINT:
239 8446 radius = 0.5;
240 8446 break;
241 502357 case SWS_SCALE_BILINEAR:
242 502357 radius = 1.0 - SWS_MAX_REDUCE_CUTOFF;
243 502357 break;
244 92959 default:
245 /* Numerically estimate radius of nontrivial or parametric kernels */
246 92959 radius = est_filter_radius(&fun);
247 92959 break;
248 }
249 603762 radius *= stretch;
250
251 603762 int filter_size = ceil(radius * 2.0);
252 603762 filter_size = FFMIN(filter_size, params->src_size);
253
1/2
✗ Branch 0 not taken.
✓ Branch 1 taken 603762 times.
 603762 av_assert0(filter_size >= 1);
254
1/2
✗ Branch 0 not taken.
✓ Branch 1 taken 603762 times.
 603762 if (filter_size > SWS_FILTER_SIZE_MAX)
255 ✗ return AVERROR(ENOTSUP);
256
257 SwsFilterWeights *filter;
258 603762 filter = av_refstruct_alloc_ext(sizeof(*filter), 0, NULL, sws_filter_free);
259
1/2
✗ Branch 0 not taken.
✓ Branch 1 taken 603762 times.
 603762 if (!filter)
260 ✗ return AVERROR(ENOMEM);
261 603762 memcpy(filter->name, fun.name, sizeof(filter->name));
262 603762 filter->src_size = params->src_size;
263 603762 filter->dst_size = params->dst_size;
264 603762 filter->virtual_size = virtual_size;
265 603762 filter->offset = params->offset;
266 603762 filter->filter_size = filter_size;
267
2/2
✓ Branch 0 taken 6732 times.
✓ Branch 1 taken 597030 times.
 603762 if (filter->filter_size == 1)
268 6732 filter->sum_positive = SWS_FILTER_SCALE;
269
270 603762 av_log(log, AV_LOG_DEBUG, "Generating %s filter with %d taps (radius = %f)\n",
271 603762 filter->name, filter->filter_size, radius);
272
273 603762 filter->num_weights = (size_t) params->dst_size * filter->filter_size;
274 603762 filter->weights = av_refstruct_allocz(filter->num_weights * sizeof(*filter->weights));
275
1/2
✗ Branch 0 not taken.
✓ Branch 1 taken 603762 times.
 603762 if (!filter->weights) {
276 ✗ av_refstruct_unref(&filter);
277 ✗ return AVERROR(ENOMEM);
278 }
279
280 603762 filter->offsets = av_refstruct_allocz(params->dst_size * sizeof(*filter->offsets));
281
1/2
✗ Branch 0 not taken.
✓ Branch 1 taken 603762 times.
 603762 if (!filter->offsets) {
282 ✗ av_refstruct_unref(&filter);
283 ✗ return AVERROR(ENOMEM);
284 }
285
286 603762 double *tmp = av_malloc(filter->filter_size * sizeof(*tmp));
287
1/2
✗ Branch 0 not taken.
✓ Branch 1 taken 603762 times.
 603762 if (!tmp) {
288 ✗ av_refstruct_unref(&filter);
289 ✗ return AVERROR(ENOMEM);
290 }
291
292 603762 const double ratio_inv = 1.0 / ratio, stretch_inv = 1.0 / stretch;
293
2/2
✓ Branch 0 taken 21797760 times.
✓ Branch 1 taken 603762 times.
 22401522 for (int i = 0; i < params->dst_size; i++)
294 21797760 compute_row(filter, &fun, radius, ratio_inv, stretch_inv, i, tmp);
295 603762 av_free(tmp);
296
297 603762 *out = filter;
298 603762 return 0;
299 }
300
301 /*
302 * Some of the filter code originally derives (via libplacebo/mpv) from Glumpy:
303 * # Copyright (c) 2009-2016 Nicolas P. Rougier. All rights reserved.
304 * # Distributed under the (new) BSD License.
305 * (https://github.com/glumpy/glumpy/blob/master/glumpy/library/build-spatial-filters.py)
306 *
307 * The math underlying each filter function was written from scratch, with
308 * some algorithms coming from a number of different sources, including:
309 * - https://en.wikipedia.org/wiki/Window_function
310 * - https://en.wikipedia.org/wiki/Jinc
311 * - http://vector-agg.cvs.sourceforge.net/viewvc/vector-agg/agg-2.5/include/agg_image_filters.h
312 * - Vapoursynth plugin fmtconv (WTFPL Licensed), which is based on
313 * dither plugin for avisynth from the same author:
314 * https://github.com/vapoursynth/fmtconv/tree/master/src/fmtc
315 * - Paul Heckbert's "zoom"
316 * - XBMC: ConvolutionKernels.cpp etc.
317 * - https://github.com/AviSynth/jinc-resize (only used to verify the math)
318 */
319
320 262400 av_unused static double box(double x, const double *params)
321 {
322 262400 return 1.0;
323 }
324
325 32939030 av_unused static double triangle(double x, const double *params)
326 {
327 32939030 return 1.0 - x;
328 }
329
330 ✗ av_unused static double cosine(double x, const double *params)
331 {
332 ✗ return cos(x);
333 }
334
335 ✗ av_unused static double hann(double x, const double *params)
336 {
337 ✗ return 0.5 + 0.5 * cos(M_PI * x);
338 }
339
340 ✗ av_unused static double hamming(double x, const double *params)
341 {
342 ✗ return 0.54 + 0.46 * cos(M_PI * x);
343 }
344
345 ✗ av_unused static double welch(double x, const double *params)
346 {
347 ✗ return 1.0 - x * x;
348 }
349
350 ✗ av_unused static double bessel_i0(double x)
351 {
352 ✗ double s = 1.0;
353 ✗ double y = x * x / 4.0;
354 ✗ double t = y;
355 ✗ int i = 2;
356 ✗ while (t > 1e-12) {
357 ✗ s += t;
358 ✗ t *= y / (i * i);
359 ✗ i += 1;
360 }
361 ✗ return s;
362 }
363
364 ✗ av_unused static double kaiser(double x, const double *params)
365 {
366 ✗ double alpha = fmax(params[0], 0.0);
367 ✗ double scale = bessel_i0(alpha);
368 ✗ return bessel_i0(alpha * sqrt(1.0 - x * x)) / scale;
369 }
370
371 ✗ av_unused static double blackman(double x, const double *params)
372 {
373 ✗ double a = params[0];
374 ✗ double a0 = (1 - a) / 2.0, a1 = 1 / 2.0, a2 = a / 2.0;
375 ✗ x *= M_PI;
376 ✗ return a0 + a1 * cos(x) + a2 * cos(2 * x);
377 }
378
379 ✗ av_unused static double bohman(double x, const double *params)
380 {
381 ✗ double pix = M_PI * x;
382 ✗ return (1.0 - x) * cos(pix) + sin(pix) / M_PI;
383 }
384
385 ✗ av_unused static double gaussian(double x, const double *params)
386 {
387 ✗ return exp(-params[0] * x * x);
388 }
389
390 ✗ av_unused static double quadratic(double x, const double *params)
391 {
392 ✗ if (x < 0.5) {
393 ✗ return 1.0 - 4.0/3.0 * (x * x);
394 } else {
395 ✗ return 2.0 / 3.0 * (x - 1.5) * (x - 1.5);
396 }
397 }
398
399 2771328 av_unused static double sinc(double x, const double *params)
400 {
401
2/2
✓ Branch 0 taken 13440 times.
✓ Branch 1 taken 2757888 times.
 2771328 if (x < 1e-8)
402 13440 return 1.0;
403 2757888 x *= M_PI;
404 2757888 return sin(x) / x;
405 }
406
407 ✗ av_unused static double jinc(double x, const double *params)
408 {
409 ✗ if (x < 1e-8)
410 ✗ return 1.0;
411 ✗ x *= M_PI;
412 ✗ return 2.0 * j1(x) / x;
413 }
414
415 ✗ av_unused static double sphinx(double x, const double *params)
416 {
417 ✗ if (x < 1e-8)
418 ✗ return 1.0;
419 ✗ x *= M_PI;
420 ✗ return 3.0 * (sin(x) - x * cos(x)) / (x * x * x);
421 }
422
423 35676383 av_unused static double cubic(double x, const double *params)
424 {
425 35676383 const double b = params[0], c = params[1];
426 35676383 double p0 = 6.0 - 2.0 * b,
427 35676383 p2 = -18.0 + 12.0 * b + 6.0 * c,
428 35676383 p3 = 12.0 - 9.0 * b - 6.0 * c,
429 35676383 q0 = 8.0 * b + 24.0 * c,
430 35676383 q1 = -12.0 * b - 48.0 * c,
431 35676383 q2 = 6.0 * b + 30.0 * c,
432 35676383 q3 = -b - 6.0 * c;
433
434
2/2
✓ Branch 0 taken 17953887 times.
✓ Branch 1 taken 17722496 times.
 35676383 if (x < 1.0) {
435 17953887 return (p0 + x * x * (p2 + x * p3)) / p0;
436 } else {
437 17722496 return (q0 + x * (q1 + x * (q2 + x * q3))) / p0;
438 }
439 }
440
441 ✗ static double spline_coeff(double a, double b, double c, double d, double x)
442 {
443 ✗ if (x <= 1.0) {
444 ✗ return ((d * x + c) * x + b) * x + a;
445 } else {
446 ✗ return spline_coeff(0.0,
447 ✗ b + 2.0 * c + 3.0 * d,
448 ✗ c + 3.0 * d,
449 ✗ -b - 3.0 * c - 6.0 * d,
450 x - 1.0);
451 }
452 }
453
454 ✗ av_unused static double spline(double x, const double *params)
455 {
456 ✗ const double p = -2.196152422706632;
457 ✗ return spline_coeff(1.0, 0.0, p, -p - 1.0, x);
458 }
459
460 static const SwsFilterFunction filter_functions[SWS_SCALE_NB] = {
461 [SWS_SCALE_BILINEAR] = { "bilinear", 1.0, triangle },
462 [SWS_SCALE_BICUBIC] = { "bicubic", 2.0, cubic, .params = { 0.0, 0.6 } },
463 [SWS_SCALE_POINT] = { "point", 0.5, box },
464 [SWS_SCALE_GAUSSIAN] = { "gaussian", 4.0, gaussian, .params = { 3.0 } },
465 [SWS_SCALE_SINC] = { "sinc", RADIUS_MAX, sinc },
466 [SWS_SCALE_LANCZOS] = { "lanczos", -1.0, sinc, sinc, .params = { 3.0 } },
467 [SWS_SCALE_SPLINE] = { "spline", RADIUS_MAX, spline },
468 /* SWS_SCALE_AREA is a pseudo-filter, see code above */
469 };
470