From c47e2adf4ece851a8e4bf88ecf251e6cea9cb7de Mon Sep 17 00:00:00 2001 From: Paul B Mahol Date: Wed, 24 Aug 2016 22:13:03 +0200 Subject: [PATCH] avfilter: add vaguedenoiser filter --- libavfilter/Makefile | 1 + libavfilter/allfilters.c | 1 + libavfilter/vf_vaguedenoiser.c | 543 +++++++++++++++++++++++++++++++++++++++++ 3 files changed, 545 insertions(+) create mode 100644 libavfilter/vf_vaguedenoiser.c diff --git a/libavfilter/Makefile b/libavfilter/Makefile index d00b40c..4ec7d8a 100644 --- a/libavfilter/Makefile +++ b/libavfilter/Makefile @@ -284,6 +284,7 @@ OBJS-$(CONFIG_TRANSPOSE_FILTER) += vf_transpose.o OBJS-$(CONFIG_TRIM_FILTER) += trim.o OBJS-$(CONFIG_UNSHARP_FILTER) += vf_unsharp.o OBJS-$(CONFIG_USPP_FILTER) += vf_uspp.o +OBJS-$(CONFIG_VAGUEDENOISER_FILTER) += vf_vaguedenoiser.o OBJS-$(CONFIG_VECTORSCOPE_FILTER) += vf_vectorscope.o OBJS-$(CONFIG_VFLIP_FILTER) += vf_vflip.o OBJS-$(CONFIG_VIDSTABDETECT_FILTER) += vidstabutils.o vf_vidstabdetect.o diff --git a/libavfilter/allfilters.c b/libavfilter/allfilters.c index ca1bfb7..1ca2679 100644 --- a/libavfilter/allfilters.c +++ b/libavfilter/allfilters.c @@ -300,6 +300,7 @@ void avfilter_register_all(void) REGISTER_FILTER(TRIM, trim, vf); REGISTER_FILTER(UNSHARP, unsharp, vf); REGISTER_FILTER(USPP, uspp, vf); + REGISTER_FILTER(VAGUEDENOISER, vaguedenoiser, vf); REGISTER_FILTER(VECTORSCOPE, vectorscope, vf); REGISTER_FILTER(VFLIP, vflip, vf); REGISTER_FILTER(VIDSTABDETECT, vidstabdetect, vf); diff --git a/libavfilter/vf_vaguedenoiser.c b/libavfilter/vf_vaguedenoiser.c new file mode 100644 index 0000000..b8bc6e8 --- /dev/null +++ b/libavfilter/vf_vaguedenoiser.c @@ -0,0 +1,543 @@ +/* + * Copyright (c) 2003 LeFunGus, lefungus@altern.org + * + * This file is part of FFmpeg + * + * FFmpeg is free software; you can redistribute it and/or modify + * it under the terms of the GNU General Public License as published by + * the Free Software Foundation; either version 2 of the License, or + * (at your option) any later version. + * + * FFmpeg is distributed in the hope that it will be useful, + * but WITHOUT ANY WARRANTY; without even the implied warranty of + * MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the + * GNU General Public License for more details. + * + * You should have received a copy of the GNU General Public License along + * with FFmpeg; if not, write to the Free Software Foundation, Inc., + * 51 Franklin Street, Fifth Floor, Boston, MA 02110-1301 USA. + */ + +#include + +#include "config.h" +#include "libavutil/attributes.h" +#include "libavutil/common.h" +#include "libavutil/pixdesc.h" +#include "libavutil/intreadwrite.h" +#include "libavutil/opt.h" + +#include "avfilter.h" +#include "formats.h" +#include "internal.h" +#include "video.h" + +typedef struct VagueContext { + const AVClass *class; + + float threshold; + float percent; + int method; + int nsteps; + int planes; + + int nb_planes; + int planeheight[4]; + int planewidth[4]; + + float *block; + float *in; + float *out; + float *tmp; + + float lower[4]; + float upper[4]; + + int peak; + int hlowsize[32]; + int hhighsize[32]; + int vlowsize[32]; + int vhighsize[32]; + + float analysis_low[32]; + float analysis_high[32]; + float synthesis_low[32]; + float synthesis_high[32]; +} VagueContext; + +#define OFFSET(x) offsetof(VagueContext, x) +#define FLAGS AV_OPT_FLAG_VIDEO_PARAM | AV_OPT_FLAG_FILTERING_PARAM +static const AVOption vaguedenoiser_options[] = { + { "threshold", "set filtering strength", OFFSET(threshold), AV_OPT_TYPE_FLOAT, {.dbl=2.}, 0, DBL_MAX, FLAGS }, + { "method", "set filtering method", OFFSET(method), AV_OPT_TYPE_INT, {.i64=0 }, 0, 2, FLAGS, "method" }, + { "hard", "hard thresholding", 0, AV_OPT_TYPE_CONST, {.i64=0}, 0, 0, FLAGS, "method" }, + { "soft", "soft thresholding", 0, AV_OPT_TYPE_CONST, {.i64=1}, 0, 0, FLAGS, "method" }, + { "garrote", "garotte thresholding", 0, AV_OPT_TYPE_CONST, {.i64=2}, 0, 0, FLAGS, "method" }, + { "nsteps", "set number of steps", OFFSET(nsteps), AV_OPT_TYPE_INT, {.i64=6 }, 1, 32, FLAGS }, + { "percent", "set percent of full denoising", OFFSET(percent), AV_OPT_TYPE_FLOAT, {.dbl=85}, 0,100, FLAGS }, + { "planes", "set planes to filter", OFFSET(planes), AV_OPT_TYPE_INT, {.i64=1 }, 0, 15, FLAGS }, + { NULL } +}; + +AVFILTER_DEFINE_CLASS(vaguedenoiser); + +#define NPAD 10 + +static const float analysis_low[9] = { + 0.037828455506995f, -0.023849465019380f, -0.110624404418423f, 0.377402855612654f, + 0.852698679009403f, 0.377402855612654f, -0.110624404418423f, -0.023849465019380f, 0.037828455506995f +}; + +static const float analysis_high[7] = { + -0.064538882628938f, 0.040689417609558f, 0.418092273222212f, -0.788485616405664f, + 0.418092273222212f, 0.040689417609558f, -0.064538882628938f +}; + +static const float synthesis_low[7] = { + -0.064538882628938f, -0.040689417609558f, 0.418092273222212f, 0.788485616405664f, + 0.418092273222212f, -0.040689417609558f, -0.064538882628938f +}; + +static const float synthesis_high[9] = { + -0.037828455506995f, -0.023849465019380f, 0.110624404418423f, 0.377402855612654f, + -0.852698679009403f, 0.377402855612654f, 0.110624404418423f, -0.023849465019380f, -0.037828455506995f +}; + +static int query_formats(AVFilterContext *ctx) +{ + static const enum AVPixelFormat pix_fmts[] = { + AV_PIX_FMT_GRAY8, + AV_PIX_FMT_GRAY16, + AV_PIX_FMT_YUV410P, AV_PIX_FMT_YUV411P, + AV_PIX_FMT_YUV420P, AV_PIX_FMT_YUV422P, + AV_PIX_FMT_YUV440P, AV_PIX_FMT_YUV444P, + AV_PIX_FMT_YUVJ420P, AV_PIX_FMT_YUVJ422P, + AV_PIX_FMT_YUVJ440P, AV_PIX_FMT_YUVJ444P, + AV_PIX_FMT_YUVJ411P, + AV_PIX_FMT_YUV420P9, AV_PIX_FMT_YUV422P9, AV_PIX_FMT_YUV444P9, + AV_PIX_FMT_YUV420P10, AV_PIX_FMT_YUV422P10, AV_PIX_FMT_YUV444P10, + AV_PIX_FMT_YUV440P10, + AV_PIX_FMT_YUV444P12, AV_PIX_FMT_YUV422P12, AV_PIX_FMT_YUV420P12, + AV_PIX_FMT_YUV440P12, + AV_PIX_FMT_YUV444P14, AV_PIX_FMT_YUV422P14, AV_PIX_FMT_YUV420P14, + AV_PIX_FMT_YUV420P16, AV_PIX_FMT_YUV422P16, AV_PIX_FMT_YUV444P16, + AV_PIX_FMT_GBRP, AV_PIX_FMT_GBRP9, AV_PIX_FMT_GBRP10, + AV_PIX_FMT_GBRP12, AV_PIX_FMT_GBRP14, AV_PIX_FMT_GBRP16, + AV_PIX_FMT_NONE + }; + AVFilterFormats *fmts_list = ff_make_format_list(pix_fmts); + if (!fmts_list) + return AVERROR(ENOMEM); + return ff_set_common_formats(ctx, fmts_list); +} + +static int config_input(AVFilterLink *inlink) +{ + VagueContext *s = inlink->dst->priv; + const AVPixFmtDescriptor *desc = av_pix_fmt_desc_get(inlink->format); + int depth; + + depth = desc->comp[0].depth; + s->nb_planes = desc->nb_components; + + s->planeheight[1] = s->planeheight[2] = AV_CEIL_RSHIFT(inlink->h, desc->log2_chroma_h); + s->planeheight[0] = s->planeheight[3] = inlink->h; + s->planewidth[1] = s->planewidth[2] = AV_CEIL_RSHIFT(inlink->w, desc->log2_chroma_w); + s->planewidth[0] = s->planewidth[3] = inlink->w; + + s->block = av_malloc_array(inlink->w * inlink->h, sizeof(*s->block)); + s->in = av_malloc_array(inlink->w * inlink->h, sizeof(*s->in)); + s->out = av_malloc_array(inlink->w * inlink->h, sizeof(*s->out)); + s->tmp = av_malloc_array(inlink->w * inlink->h, sizeof(*s->tmp)); + + if (!s->block || !s->in || !s->out || !s->tmp) + return AVERROR(ENOMEM); + + s->threshold *= 1 << (depth - 8); + s->peak = (1 << depth) - 1; + + s->lower[0] = 0; + s->upper[0] = 1; + s->lower[1] = -0.5; + s->upper[1] = 0.5; + s->lower[2] = -0.5; + s->upper[2] = 0.5; + s->lower[3] = 1; + s->upper[3] = 1; + + memcpy(s->analysis_low, analysis_low, sizeof(analysis_low)); + memcpy(s->analysis_high, analysis_high, sizeof(analysis_high)); + memcpy(s->synthesis_low, synthesis_low, sizeof(synthesis_low)); + memcpy(s->synthesis_high, synthesis_high, sizeof(synthesis_high)); + + return 0; +} + +static inline void copy(const float *p1, float *p2, const int length) +{ + memcpy(p2, p1, length * sizeof(float)); +} + +static inline void copyv(const float *p1, const int stride1, float *p2, const int length) +{ + int i; + + for (i = 0; i < length; i++) { + p2[i] = *p1; + p1 += stride1; + } +} + +static inline void copyh(const float *p1, float *p2, const int stride2, const int length) +{ + int i; + + for (i = 0; i < length; i++) { + *p2 = p1[i]; + p2 += stride2; + } +} + +// Do symmetric extension of data using prescribed symmetries +// Original values are in output[npad] through output[npad+size-1] +// New values will be placed in output[0] through output[npad] and in output[npad+size] through output[2*npad+size-1] (note: end values may not be filled in) +// extension at left bdry is ... 3 2 1 0 | 0 1 2 3 ... +// same for right boundary +// if rightExt=1 then ... 3 2 1 0 | 1 2 3 +static void symmetricExtension(float *output, const int size, const int leftExt, const int rightExt) +{ + int first = NPAD; + int last = NPAD - 1 + size; + + const int originalLast = last; + + if (leftExt == 2) + output[--first] = output[NPAD]; + if (rightExt == 2) + output[++last] = output[originalLast]; + + // extend left end + int nextend = first; + for (int i = 0; i < nextend; i++) + output[--first] = output[NPAD + 1 + i]; + + const int idx = NPAD + NPAD - 1 + size; + + // extend right end + nextend = idx - last; + for (int i = 0; i < nextend; i++) + output[++last] = output[originalLast - 1 - i]; +} + +static void transformStep(float *input, float *output, const int size, const int lowSize, VagueContext *s) +{ + symmetricExtension(input, size, 1, 1); + + for (int i = NPAD; i < NPAD + lowSize; i++) { + const float a = input[2 * i - 14] * s->analysis_low[0]; + const float b = input[2 * i - 13] * s->analysis_low[1]; + const float c = input[2 * i - 12] * s->analysis_low[2]; + const float d = input[2 * i - 11] * s->analysis_low[3]; + const float e = input[2 * i - 10] * s->analysis_low[4]; + const float f = input[2 * i - 9] * s->analysis_low[3]; + const float g = input[2 * i - 8] * s->analysis_low[2]; + const float h = input[2 * i - 7] * s->analysis_low[1]; + const float k = input[2 * i - 6] * s->analysis_low[0]; + + output[i] = a + b + c + d + e + f + g + h + k; + } + + for (int i = NPAD; i < NPAD + lowSize; i++) { + const float a = input[2 * i - 12] * s->analysis_high[0]; + const float b = input[2 * i - 11] * s->analysis_high[1]; + const float c = input[2 * i - 10] * s->analysis_high[2]; + const float d = input[2 * i - 9] * s->analysis_high[3]; + const float e = input[2 * i - 8] * s->analysis_high[2]; + const float f = input[2 * i - 7] * s->analysis_high[1]; + const float g = input[2 * i - 6] * s->analysis_high[0]; + + output[i + lowSize] = a + b + c + d + e + f + g; + } +} + +static void invertStep(const float *input, float *output, float *temp, const int size, VagueContext *s) +{ + const int lowSize = (size + 1) >> 1; + const int highSize = size >> 1; + int findex; + + memcpy(temp + NPAD, input + NPAD, lowSize * sizeof(float)); + + int leftExt = 1; + int rightExt = (size % 2 == 0) ? 2 : 1; + symmetricExtension(temp, lowSize, leftExt, rightExt); + + memset(output, 0, (NPAD + NPAD + size) * sizeof(float)); + findex = (size + 2) >> 1; + + for (int i = 9; i < findex + 11; i++) { + const float a = temp[i] * s->synthesis_low[0]; + const float b = temp[i] * s->synthesis_low[1]; + const float c = temp[i] * s->synthesis_low[2]; + const float d = temp[i] * s->synthesis_low[3]; + + output[2 * i - 13] += a; + output[2 * i - 12] += b; + output[2 * i - 11] += c; + output[2 * i - 10] += d; + output[2 * i - 9] += c; + output[2 * i - 8] += b; + output[2 * i - 7] += a; + } + + memcpy(temp + NPAD, input + NPAD + lowSize, highSize * sizeof(float)); + + leftExt = 2; + rightExt = (size % 2 == 0) ? 1 : 2; + symmetricExtension(temp, highSize, leftExt, rightExt); + + for (int i = 8; i < findex + 11; i++) { + const float a = temp[i] * s->synthesis_high[0]; + const float b = temp[i] * s->synthesis_high[1]; + const float c = temp[i] * s->synthesis_high[2]; + const float d = temp[i] * s->synthesis_high[3]; + const float e = temp[i] * s->synthesis_high[4]; + + output[2 * i - 13] += a; + output[2 * i - 12] += b; + output[2 * i - 11] += c; + output[2 * i - 10] += d; + output[2 * i - 9] += e; + output[2 * i - 8] += d; + output[2 * i - 7] += c; + output[2 * i - 6] += b; + output[2 * i - 5] += a; + } +} + +static void hardThresholding(float *block, const int width, const int height, const int stride, const float threshold, const float percent) +{ + const float frac = 1.f - percent * 0.01f; + + for (int y = 0; y < height; y++) { + for (int x = 0; x < width; x++) { + if (FFABS(block[x]) <= threshold) + block[x] *= frac; + } + block += stride; + } +} + +static void softThresholding(float *block, const int width, const int height, const int stride, const float threshold, const float percent, const int nsteps) +{ + const float frac = 1.f - percent * 0.01f; + const float shift = threshold * 0.01f * percent; + + int w = width; + int h = height; + for (int l = 0; l < nsteps; l++) { + w = (w + 1) >> 1; + h = (h + 1) >> 1; + } + + for (int y = 0; y < height; y++) { + const int x0 = (y < h) ? w : 0; + for (int x = x0; x < width; x++) { + const float temp = FFABS(block[x]); + if (temp <= threshold) + block[x] *= frac; + else + block[x] = (block[x] < 0.f ? -1.f : (block[x] > 0.f ? 1.f : 0.f)) * (temp - shift); + } + block += stride; + } +} + +static void qianThresholding(float *block, const int width, const int height, + const int stride, const float threshold, + const float percent) +{ + const float percent01 = percent * 0.01f; + const float tr2 = threshold * threshold * percent01; + const float frac = 1.f - percent01; + + for (int y = 0; y < height; y++) { + for (int x = 0; x < width; x++) { + const float temp = FFABS(block[x]); + if (temp <= threshold) { + block[x] *= frac; + } else { + const float tp2 = temp * temp; + block[x] *= (tp2 - tr2) / tp2; + } + } + block += stride; + } +} + +static void filter(VagueContext *s, AVFrame *in, AVFrame *out) +{ + int p, y, x, i, j; + + for (p = 0; p < s->nb_planes; p++) { + const int height = s->planeheight[p]; + const int width = s->planewidth[p]; + const uint8_t *srcp = in->data[p]; + uint8_t *dstp = out->data[p]; + float *output = s->block; + + for (y = 0; y < height; y++) { + for (x = 0; x < width; x++) + output[x] = srcp[x]; + srcp += in->linesize[p]; + output += width; + } + + int hLowSize0 = width; + int vLowSize0 = height; + int nstepsTransform = s->nsteps; + while (nstepsTransform--) { + int lowSize = (hLowSize0 + 1) >> 1; + float *input = s->block; + for (j = 0; j < vLowSize0; j++) { + copy(input, s->in + NPAD, hLowSize0); + transformStep(s->in, s->out, hLowSize0, lowSize, s); + copy(s->out + NPAD, input, hLowSize0); + input += width; + } + + lowSize = (vLowSize0 + 1) >> 1; + input = s->block; + for (j = 0; j < hLowSize0; j++) { + copyv(input, width, s->in + NPAD, vLowSize0); + transformStep(s->in, s->out, vLowSize0, lowSize, s); + copyh(s->out + NPAD, input, width, vLowSize0); + input++; + } + + hLowSize0 = (hLowSize0 + 1) >> 1; + vLowSize0 = (vLowSize0 + 1) >> 1; + } + + if (s->method == 0) + hardThresholding(s->block, width, height, width, s->threshold, s->percent); + else if (s->method == 1) + softThresholding(s->block, width, height, width, s->threshold, s->percent, s->nsteps); + else + qianThresholding(s->block, width, height, width, s->threshold, s->percent); + + s->hlowsize[0] = (width + 1) >> 1; + s->hhighsize[0] = width >> 1; + s->vlowsize[0] = (height + 1) >> 1; + s->vhighsize[0] = height >> 1; + for (i = 1; i < s->nsteps; i++) { + s->hlowsize[i] = (s->hlowsize[i - 1] + 1) >> 1; + s->hhighsize[i] = s->hlowsize[i - 1] >> 1; + s->vlowsize[i] = (s->vlowsize[i - 1] + 1) >> 1; + s->vhighsize[i] = s->vlowsize[i - 1] >> 1; + } + + int nstepsInvert = s->nsteps; + while (nstepsInvert--) { + const int idx = s->vlowsize[nstepsInvert] + s->vhighsize[nstepsInvert]; + const int idx2 = s->hlowsize[nstepsInvert] + s->hhighsize[nstepsInvert]; + float * idx3 = s->block; + for (i = 0; i < idx2; i++) { + copyv(idx3, width, s->in + NPAD, idx); + invertStep(s->in, s->out, s->tmp, idx, s); + copyh(s->out + NPAD, idx3, width, idx); + idx3++; + } + + idx3 = s->block; + for (i = 0; i < idx; i++) { + copy(idx3, s->in + NPAD, idx2); + invertStep(s->in, s->out, s->tmp, idx2, s); + copy(s->out + NPAD, idx3, idx2); + idx3 += width; + } + } + + const float *input = s->block; + for (y = 0; y < height; y++) { + for (x = 0; x < width; x++) + dstp[x] = av_clip_uint8(input[x] + 0.5f); + input += width; + dstp += out->linesize[p]; + } + } +} + +static int filter_frame(AVFilterLink *inlink, AVFrame *in) +{ + AVFilterContext *ctx = inlink->dst; + VagueContext *s = ctx->priv; + AVFilterLink *outlink = ctx->outputs[0]; + AVFrame *out; + int direct = av_frame_is_writable(in) && !ctx->is_disabled; + + if (direct) { + out = in; + } else { + out = ff_get_video_buffer(outlink, outlink->w, outlink->h); + if (!out) { + av_frame_free(&in); + return AVERROR(ENOMEM); + } + + av_frame_copy_props(out, in); + } + + filter(s, in, out); + + if (ctx->is_disabled) { + av_frame_free(&out); + return ff_filter_frame(outlink, in); + } + + if (!direct) + av_frame_free(&in); + + return ff_filter_frame(outlink, out); +} + +static av_cold void uninit(AVFilterContext *ctx) +{ + VagueContext *s = ctx->priv; + + av_freep(&s->block); + av_freep(&s->in); + av_freep(&s->out); + av_freep(&s->tmp); +} + +static const AVFilterPad vaguedenoiser_inputs[] = { + { + .name = "default", + .type = AVMEDIA_TYPE_VIDEO, + .config_props = config_input, + .filter_frame = filter_frame, + }, + { NULL } +}; + + +static const AVFilterPad vaguedenoiser_outputs[] = { + { + .name = "default", + .type = AVMEDIA_TYPE_VIDEO + }, + { NULL } +}; + +AVFilter ff_vf_vaguedenoiser = { + .name = "vaguedenoiser", + .description = NULL_IF_CONFIG_SMALL("Apply a Wavelet Denoiser."), + .priv_size = sizeof(VagueContext), + .priv_class = &vaguedenoiser_class, + .uninit = uninit, + .query_formats = query_formats, + .inputs = vaguedenoiser_inputs, + .outputs = vaguedenoiser_outputs, + .flags = AVFILTER_FLAG_SUPPORT_TIMELINE_INTERNAL, +}; -- 2.9.1