1 /*
2 * Copyright (c) 2007 Michael Niedermayer <michaelni@gmx.at>
3 * Copyright (c) 2013 Clément Bœsch <u pkh me>
4 *
5 * This file is part of FFmpeg.
6 *
7 * FFmpeg is free software; you can redistribute it and/or modify
8 * it under the terms of the GNU General Public License as published by
9 * the Free Software Foundation; either version 2 of the License, or
10 * (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
15 * GNU General Public License for more details.
16 *
17 * You should have received a copy of the GNU General Public License along
18 * with FFmpeg; if not, write to the Free Software Foundation, Inc.,
19 * 51 Franklin Street, Fifth Floor, Boston, MA 02110-1301 USA.
20 */
21
22 // The code written by Michael Niedermayer in 70024b6b47b9eacfe01e8f92349ca9bf1ccd7d5a:libavfilter/vf_owdenoise.c
23 // can also be used under the LGPL due to:
24 // <michaelni> durandal_1707, if you do all the "todo" points from vf_owdenoise.c that are in that file since 2013 then sure i would be more than happy to relicense my part of it to LGPL
25 // <durandal_1707> michaelni: first relicense than work
26
27 /**
28 * @todo try to change to int
29 * @todo try lifting based implementation
30 * @todo optimize optimize optimize
31 * @todo hard thresholding
32 * @todo use QP to decide filter strength
33 * @todo wavelet normalization / least squares optimal signal vs. noise thresholds
34 */
35
42
53
54 #define OFFSET(x) offsetof(OWDenoiseContext, x)
55 #define FLAGS AV_OPT_FLAG_FILTERING_PARAM|AV_OPT_FLAG_VIDEO_PARAM
63 };
64
66
68 { 0, 48, 12, 60, 3, 51, 15, 63 },
69 { 32, 16, 44, 28, 35, 19, 47, 31 },
70 { 8, 56, 4, 52, 11, 59, 7, 55 },
71 { 40, 24, 36, 20, 43, 27, 39, 23 },
72 { 2, 50, 14, 62, 1, 49, 13, 61 },
73 { 34, 18, 46, 30, 33, 17, 45, 29 },
74 { 10, 58, 6, 54, 9, 57, 5, 53 },
75 { 42, 26, 38, 22, 41, 25, 37, 21 },
76 };
77
78 static const double coeff[2][5] = {
79 {
85 },{
90 }
91 };
92
94 {
99 },{
102 -0.07822326652898785 *
M_SQRT2,
105 }
106 };
107
108
109 static inline void decompose(
float *dst_l,
float *dst_h,
const float *
src,
111 {
113 for (x = 0; x <
w; x++) {
114 double sum_l =
src[x * linesize] *
coeff[0][0];
115 double sum_h =
src[x * linesize] *
coeff[1][0];
116 for (
i = 1;
i <= 4;
i++) {
119
122 }
123 dst_l[x * linesize] = sum_l;
124 dst_h[x * linesize] = sum_h;
125 }
126 }
127
128 static inline void compose(
float *dst,
const float *src_l,
const float *src_h,
130 {
132 for (x = 0; x <
w; x++) {
133 double sum_l = src_l[x * linesize] *
icoeff[0][0];
134 double sum_h = src_h[x * linesize] *
icoeff[1][0];
135 for (
i = 1;
i <= 4;
i++) {
138
139 sum_l +=
icoeff[0][
i] * (src_l[x0] + src_l[x1]);
140 sum_h +=
icoeff[1][
i] * (src_h[x0] + src_h[x1]);
141 }
142 dst[x * linesize] = (sum_l + sum_h) * 0.5;
143 }
144 }
145
147 int xlinesize, int ylinesize,
149 {
150 int y, x;
151 for (y = 0; y <
h; y++)
152 for (x = 0; x <
step; x++)
153 decompose(dst_l + ylinesize*y + xlinesize*x,
154 dst_h + ylinesize*y + xlinesize*x,
155 src + ylinesize*y + xlinesize*x,
157 }
158
159 static inline void compose2D(
float *dst,
const float *src_l,
const float *src_h,
160 int xlinesize, int ylinesize,
162 {
163 int y, x;
164 for (y = 0; y <
h; y++)
165 for (x = 0; x <
step; x++)
166 compose(dst + ylinesize*y + xlinesize*x,
167 src_l + ylinesize*y + xlinesize*x,
168 src_h + ylinesize*y + xlinesize*x,
170 }
171
173 int linesize,
int step,
int w,
int h)
174 {
178 }
179
181 int linesize,
int step,
int w,
int h)
182 {
186 }
187
189 uint8_t *dst, int dst_linesize,
190 const uint8_t *
src,
int src_linesize,
192 {
193 int x, y,
i, j, depth =
s->depth;
194
196 depth--;
197
198 if (
s->pixel_depth <= 8) {
199 for (y = 0; y <
height; y++)
200 for(x = 0; x <
width; x++)
201 s->plane[0][0][y*
s->linesize + x] =
src[y*src_linesize + x];
202 } else {
203 const uint16_t *src16 = (
const uint16_t *)
src;
204
205 src_linesize /= 2;
206 for (y = 0; y <
height; y++)
207 for(x = 0; x <
width; x++)
208 s->plane[0][0][y*
s->linesize + x] = src16[y*src_linesize + x];
209 }
210
211 for (
i = 0;
i < depth;
i++)
213
214 for (
i = 0;
i < depth;
i++) {
215 for (j = 1; j < 4; j++) {
216 for (y = 0; y <
height; y++) {
217 for (x = 0; x <
width; x++) {
218 double v =
s->plane[
i + 1][j][y*
s->linesize + x];
219 if (v > strength) v -= strength;
220 else if (v < -strength) v += strength;
221 else v = 0;
222 s->plane[
i + 1][j][x + y*
s->linesize] = v;
223 }
224 }
225 }
226 }
227 for (
i = depth-1;
i >= 0;
i--)
229
230 if (
s->pixel_depth <= 8) {
231 for (y = 0; y <
height; y++) {
232 for (x = 0; x <
width; x++) {
233 i =
s->plane[0][0][y*
s->linesize + x] +
dither[x&7][y&7]*(1.0/64) + 1.0/128;
// yes the rounding is insane but optimal :)
234 if ((
unsigned)
i > 255
U)
i = ~(
i >> 31);
235 dst[y*dst_linesize + x] =
i;
236 }
237 }
238 } else {
239 uint16_t *dst16 = (uint16_t *)dst;
240
241 dst_linesize /= 2;
242 for (y = 0; y <
height; y++) {
243 for (x = 0; x <
width; x++) {
244 i =
s->plane[0][0][y*
s->linesize + x];
245 dst16[y*dst_linesize + x] =
i;
246 }
247 }
248 }
249 }
250
252 {
259
262
263 if (
s->luma_strength > 0)
265 if (
s->chroma_strength > 0) {
268 }
269 } else {
274 }
276
277 if (
s->luma_strength > 0) {
279 } else {
281 }
282 if (
s->chroma_strength > 0) {
285 } else {
288 }
289
295 }
296
298 }
299
314 };
315
317 {
322
323 s->hsub =
desc->log2_chroma_w;
324 s->vsub =
desc->log2_chroma_h;
325 s->pixel_depth =
desc->comp[0].depth;
326
328 for (j = 0; j < 4; j++) {
329 for (
i = 0;
i <=
s->depth;
i++) {
333 }
334 }
335 return 0;
336 }
337
339 {
342
343 for (j = 0; j < 4; j++)
344 for (
i = 0;
i <=
s->depth;
i++)
346 }
347
349 {
354 },
355 };
356
358 {
361 },
362 };
363
372 .priv_class = &owdenoise_class,
374 };