Subversion Repositories Kolibri OS

Rev

Go to most recent revision | Blame | Last modification | View Log | RSS feed

  1. /*
  2.  * Copyright (c) 2013 Clément Bœsch
  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. /**
  22.  * A simple, relatively efficient and extremely slow DCT image denoiser.
  23.  * @see http://www.ipol.im/pub/art/2011/ys-dct/
  24.  */
  25.  
  26. #include "libavcodec/avfft.h"
  27. #include "libavutil/eval.h"
  28. #include "libavutil/opt.h"
  29. #include "drawutils.h"
  30. #include "internal.h"
  31.  
  32. #define NBITS 4
  33. #define BSIZE (1<<(NBITS))
  34.  
  35. static const char *const var_names[] = { "c", NULL };
  36. enum { VAR_C, VAR_VARS_NB };
  37.  
  38. typedef struct {
  39.     const AVClass *class;
  40.  
  41.     /* coefficient factor expression */
  42.     char *expr_str;
  43.     AVExpr *expr;
  44.     double var_values[VAR_VARS_NB];
  45.  
  46.     int pr_width, pr_height;    // width and height to process
  47.     float sigma;                // used when no expression are st
  48.     float th;                   // threshold (3*sigma)
  49.     float color_dct[3][3];      // 3x3 DCT for color decorrelation
  50.     float *cbuf[2][3];          // two planar rgb color buffers
  51.     float *weights;             // dct coeff are cumulated with overlapping; these values are used for averaging
  52.     int p_linesize;             // line sizes for color and weights
  53.     int overlap;                // number of block overlapping pixels
  54.     int step;                   // block step increment (BSIZE - overlap)
  55.     DCTContext *dct, *idct;     // DCT and inverse DCT contexts
  56.     float *block, *tmp_block;   // two BSIZE x BSIZE block buffers
  57. } DCTdnoizContext;
  58.  
  59. #define OFFSET(x) offsetof(DCTdnoizContext, x)
  60. #define FLAGS AV_OPT_FLAG_FILTERING_PARAM|AV_OPT_FLAG_VIDEO_PARAM
  61. static const AVOption dctdnoiz_options[] = {
  62.     { "sigma",   "set noise sigma constant",               OFFSET(sigma),    AV_OPT_TYPE_FLOAT,  {.dbl=0},            0, 999,          .flags = FLAGS },
  63.     { "s",       "set noise sigma constant",               OFFSET(sigma),    AV_OPT_TYPE_FLOAT,  {.dbl=0},            0, 999,          .flags = FLAGS },
  64.     { "overlap", "set number of block overlapping pixels", OFFSET(overlap),  AV_OPT_TYPE_INT,    {.i64=(1<<NBITS)-1}, 0, (1<<NBITS)-1, .flags = FLAGS },
  65.     { "expr",    "set coefficient factor expression",      OFFSET(expr_str), AV_OPT_TYPE_STRING, {.str=NULL},                          .flags = FLAGS },
  66.     { "e",       "set coefficient factor expression",      OFFSET(expr_str), AV_OPT_TYPE_STRING, {.str=NULL},                          .flags = FLAGS },
  67.     { NULL }
  68. };
  69.  
  70. AVFILTER_DEFINE_CLASS(dctdnoiz);
  71.  
  72. static float *dct_block(DCTdnoizContext *ctx, const float *src, int src_linesize)
  73. {
  74.     int x, y;
  75.     float *column;
  76.  
  77.     for (y = 0; y < BSIZE; y++) {
  78.         float *line = ctx->block;
  79.  
  80.         memcpy(line, src, BSIZE * sizeof(*line));
  81.         src += src_linesize;
  82.         av_dct_calc(ctx->dct, line);
  83.  
  84.         column = ctx->tmp_block + y;
  85.         column[0] = line[0] * (1. / sqrt(BSIZE));
  86.         column += BSIZE;
  87.         for (x = 1; x < BSIZE; x++) {
  88.             *column = line[x] * sqrt(2. / BSIZE);
  89.             column += BSIZE;
  90.         }
  91.     }
  92.  
  93.     column = ctx->tmp_block;
  94.     for (x = 0; x < BSIZE; x++) {
  95.         av_dct_calc(ctx->dct, column);
  96.         column[0] *= 1. / sqrt(BSIZE);
  97.         for (y = 1; y < BSIZE; y++)
  98.             column[y] *= sqrt(2. / BSIZE);
  99.         column += BSIZE;
  100.     }
  101.  
  102.     for (y = 0; y < BSIZE; y++)
  103.         for (x = 0; x < BSIZE; x++)
  104.             ctx->block[y*BSIZE + x] = ctx->tmp_block[x*BSIZE + y];
  105.  
  106.     return ctx->block;
  107. }
  108.  
  109. static void idct_block(DCTdnoizContext *ctx, float *dst, int dst_linesize)
  110. {
  111.     int x, y;
  112.     float *block = ctx->block;
  113.     float *tmp = ctx->tmp_block;
  114.  
  115.     for (y = 0; y < BSIZE; y++) {
  116.         block[0] *= sqrt(BSIZE);
  117.         for (x = 1; x < BSIZE; x++)
  118.             block[x] *= 1./sqrt(2. / BSIZE);
  119.         av_dct_calc(ctx->idct, block);
  120.         block += BSIZE;
  121.     }
  122.  
  123.     block = ctx->block;
  124.     for (y = 0; y < BSIZE; y++) {
  125.         tmp[0] = block[y] * sqrt(BSIZE);
  126.         for (x = 1; x < BSIZE; x++)
  127.             tmp[x] = block[x*BSIZE + y] * (1./sqrt(2. / BSIZE));
  128.         av_dct_calc(ctx->idct, tmp);
  129.         for (x = 0; x < BSIZE; x++)
  130.             dst[x*dst_linesize + y] += tmp[x];
  131.     }
  132. }
  133.  
  134. static int config_input(AVFilterLink *inlink)
  135. {
  136.     AVFilterContext *ctx = inlink->dst;
  137.     DCTdnoizContext *s = ctx->priv;
  138.     int i, x, y, bx, by, linesize, *iweights;
  139.     const float dct_3x3[3][3] = {
  140.         { 1./sqrt(3),  1./sqrt(3),  1./sqrt(3) },
  141.         { 1./sqrt(2),           0, -1./sqrt(2) },
  142.         { 1./sqrt(6), -2./sqrt(6),  1./sqrt(6) },
  143.     };
  144.     uint8_t rgba_map[4];
  145.  
  146.     ff_fill_rgba_map(rgba_map, inlink->format);
  147.     for (y = 0; y < 3; y++)
  148.         for (x = 0; x < 3; x++)
  149.             s->color_dct[y][x] = dct_3x3[rgba_map[y]][rgba_map[x]];
  150.  
  151.     s->pr_width  = inlink->w - (inlink->w - BSIZE) % s->step;
  152.     s->pr_height = inlink->h - (inlink->h - BSIZE) % s->step;
  153.     if (s->pr_width != inlink->w)
  154.         av_log(ctx, AV_LOG_WARNING, "The last %d horizontal pixels won't be denoised\n",
  155.                inlink->w - s->pr_width);
  156.     if (s->pr_height != inlink->h)
  157.         av_log(ctx, AV_LOG_WARNING, "The last %d vertical pixels won't be denoised\n",
  158.                inlink->h - s->pr_height);
  159.  
  160.     s->p_linesize = linesize = FFALIGN(s->pr_width, 32);
  161.     for (i = 0; i < 2; i++) {
  162.         s->cbuf[i][0] = av_malloc(linesize * s->pr_height * sizeof(*s->cbuf[i][0]));
  163.         s->cbuf[i][1] = av_malloc(linesize * s->pr_height * sizeof(*s->cbuf[i][1]));
  164.         s->cbuf[i][2] = av_malloc(linesize * s->pr_height * sizeof(*s->cbuf[i][2]));
  165.         if (!s->cbuf[i][0] || !s->cbuf[i][1] || !s->cbuf[i][2])
  166.             return AVERROR(ENOMEM);
  167.     }
  168.  
  169.     s->weights = av_malloc(s->pr_height * linesize * sizeof(*s->weights));
  170.     if (!s->weights)
  171.         return AVERROR(ENOMEM);
  172.     iweights = av_calloc(s->pr_height, linesize * sizeof(*iweights));
  173.     if (!iweights)
  174.         return AVERROR(ENOMEM);
  175.     for (y = 0; y < s->pr_height - BSIZE + 1; y += s->step)
  176.         for (x = 0; x < s->pr_width - BSIZE + 1; x += s->step)
  177.             for (by = 0; by < BSIZE; by++)
  178.                 for (bx = 0; bx < BSIZE; bx++)
  179.                     iweights[(y + by)*linesize + x + bx]++;
  180.     for (y = 0; y < s->pr_height; y++)
  181.         for (x = 0; x < s->pr_width; x++)
  182.             s->weights[y*linesize + x] = 1. / iweights[y*linesize + x];
  183.     av_free(iweights);
  184.  
  185.     return 0;
  186. }
  187.  
  188. static av_cold int init(AVFilterContext *ctx)
  189. {
  190.     DCTdnoizContext *s = ctx->priv;
  191.  
  192.     if (s->expr_str) {
  193.         int ret = av_expr_parse(&s->expr, s->expr_str, var_names,
  194.                                 NULL, NULL, NULL, NULL, 0, ctx);
  195.         if (ret < 0)
  196.             return ret;
  197.     }
  198.  
  199.     s->th   = s->sigma * 3.;
  200.     s->step = BSIZE - s->overlap;
  201.     s->dct  = av_dct_init(NBITS, DCT_II);
  202.     s->idct = av_dct_init(NBITS, DCT_III);
  203.     s->block     = av_malloc(BSIZE * BSIZE * sizeof(*s->block));
  204.     s->tmp_block = av_malloc(BSIZE * BSIZE * sizeof(*s->tmp_block));
  205.  
  206.     if (!s->dct || !s->idct || !s->tmp_block || !s->block)
  207.         return AVERROR(ENOMEM);
  208.  
  209.     return 0;
  210. }
  211.  
  212. static int query_formats(AVFilterContext *ctx)
  213. {
  214.     static const enum AVPixelFormat pix_fmts[] = {
  215.         AV_PIX_FMT_BGR24, AV_PIX_FMT_RGB24,
  216.         AV_PIX_FMT_NONE
  217.     };
  218.     ff_set_common_formats(ctx, ff_make_format_list(pix_fmts));
  219.     return 0;
  220. }
  221.  
  222. static void color_decorrelation(float dct3ch[3][3], float **dst, int dst_linesize,
  223.                                 const uint8_t *src, int src_linesize, int w, int h)
  224. {
  225.     int x, y;
  226.     float *dstp_r = dst[0];
  227.     float *dstp_g = dst[1];
  228.     float *dstp_b = dst[2];
  229.  
  230.     for (y = 0; y < h; y++) {
  231.         const uint8_t *srcp = src;
  232.  
  233.         for (x = 0; x < w; x++) {
  234.             dstp_r[x] = srcp[0] * dct3ch[0][0] + srcp[1] * dct3ch[0][1] + srcp[2] * dct3ch[0][2];
  235.             dstp_g[x] = srcp[0] * dct3ch[1][0] + srcp[1] * dct3ch[1][1] + srcp[2] * dct3ch[1][2];
  236.             dstp_b[x] = srcp[0] * dct3ch[2][0] + srcp[1] * dct3ch[2][1] + srcp[2] * dct3ch[2][2];
  237.             srcp += 3;
  238.         }
  239.         src += src_linesize;
  240.         dstp_r += dst_linesize;
  241.         dstp_g += dst_linesize;
  242.         dstp_b += dst_linesize;
  243.     }
  244. }
  245.  
  246. static void color_correlation(float dct3ch[3][3], uint8_t *dst, int dst_linesize,
  247.                               float **src, int src_linesize, int w, int h)
  248. {
  249.     int x, y;
  250.     const float *src_r = src[0];
  251.     const float *src_g = src[1];
  252.     const float *src_b = src[2];
  253.  
  254.     for (y = 0; y < h; y++) {
  255.         uint8_t *dstp = dst;
  256.  
  257.         for (x = 0; x < w; x++) {
  258.             dstp[0] = av_clip_uint8(src_r[x] * dct3ch[0][0] + src_g[x] * dct3ch[1][0] + src_b[x] * dct3ch[2][0]);
  259.             dstp[1] = av_clip_uint8(src_r[x] * dct3ch[0][1] + src_g[x] * dct3ch[1][1] + src_b[x] * dct3ch[2][1]);
  260.             dstp[2] = av_clip_uint8(src_r[x] * dct3ch[0][2] + src_g[x] * dct3ch[1][2] + src_b[x] * dct3ch[2][2]);
  261.             dstp += 3;
  262.         }
  263.         dst += dst_linesize;
  264.         src_r += src_linesize;
  265.         src_g += src_linesize;
  266.         src_b += src_linesize;
  267.     }
  268. }
  269.  
  270. static void filter_plane(AVFilterContext *ctx,
  271.                          float *dst, int dst_linesize,
  272.                          const float *src, int src_linesize,
  273.                          int w, int h)
  274. {
  275.     int x, y, bx, by;
  276.     DCTdnoizContext *s = ctx->priv;
  277.     float *dst0 = dst;
  278.     const float *weights = s->weights;
  279.  
  280.     // reset block sums
  281.     memset(dst, 0, h * dst_linesize * sizeof(*dst));
  282.  
  283.     // block dct sums
  284.     for (y = 0; y < h - BSIZE + 1; y += s->step) {
  285.         for (x = 0; x < w - BSIZE + 1; x += s->step) {
  286.             float *ftb = dct_block(s, src + x, src_linesize);
  287.  
  288.             if (s->expr) {
  289.                 for (by = 0; by < BSIZE; by++) {
  290.                     for (bx = 0; bx < BSIZE; bx++) {
  291.                         s->var_values[VAR_C] = FFABS(*ftb);
  292.                         *ftb++ *= av_expr_eval(s->expr, s->var_values, s);
  293.                     }
  294.                 }
  295.             } else {
  296.                 for (by = 0; by < BSIZE; by++) {
  297.                     for (bx = 0; bx < BSIZE; bx++) {
  298.                         if (FFABS(*ftb) < s->th)
  299.                             *ftb = 0;
  300.                         ftb++;
  301.                     }
  302.                 }
  303.             }
  304.             idct_block(s, dst + x, dst_linesize);
  305.         }
  306.         src += s->step * src_linesize;
  307.         dst += s->step * dst_linesize;
  308.     }
  309.  
  310.     // average blocks
  311.     dst = dst0;
  312.     for (y = 0; y < h; y++) {
  313.         for (x = 0; x < w; x++)
  314.             dst[x] *= weights[x];
  315.         dst += dst_linesize;
  316.         weights += dst_linesize;
  317.     }
  318. }
  319.  
  320. static int filter_frame(AVFilterLink *inlink, AVFrame *in)
  321. {
  322.     AVFilterContext *ctx = inlink->dst;
  323.     DCTdnoizContext *s = ctx->priv;
  324.     AVFilterLink *outlink = inlink->dst->outputs[0];
  325.     int direct, plane;
  326.     AVFrame *out;
  327.  
  328.     if (av_frame_is_writable(in)) {
  329.         direct = 1;
  330.         out = in;
  331.     } else {
  332.         direct = 0;
  333.         out = ff_get_video_buffer(outlink, outlink->w, outlink->h);
  334.         if (!out) {
  335.             av_frame_free(&in);
  336.             return AVERROR(ENOMEM);
  337.         }
  338.         av_frame_copy_props(out, in);
  339.     }
  340.  
  341.     color_decorrelation(s->color_dct, s->cbuf[0], s->p_linesize,
  342.                         in->data[0], in->linesize[0], s->pr_width, s->pr_height);
  343.     for (plane = 0; plane < 3; plane++)
  344.         filter_plane(ctx, s->cbuf[1][plane], s->p_linesize,
  345.                           s->cbuf[0][plane], s->p_linesize,
  346.                           s->pr_width, s->pr_height);
  347.     color_correlation(s->color_dct, out->data[0], out->linesize[0],
  348.                       s->cbuf[1], s->p_linesize, s->pr_width, s->pr_height);
  349.  
  350.     if (!direct) {
  351.         int y;
  352.         uint8_t *dst = out->data[0];
  353.         const uint8_t *src = in->data[0];
  354.         const int dst_linesize = out->linesize[0];
  355.         const int src_linesize = in->linesize[0];
  356.         const int hpad = (inlink->w - s->pr_width) * 3;
  357.         const int vpad = (inlink->h - s->pr_height);
  358.  
  359.         if (hpad) {
  360.             uint8_t       *dstp = dst + s->pr_width * 3;
  361.             const uint8_t *srcp = src + s->pr_width * 3;
  362.  
  363.             for (y = 0; y < s->pr_height; y++) {
  364.                 memcpy(dstp, srcp, hpad);
  365.                 dstp += dst_linesize;
  366.                 srcp += src_linesize;
  367.             }
  368.         }
  369.         if (vpad) {
  370.             uint8_t       *dstp = dst + s->pr_height * dst_linesize;
  371.             const uint8_t *srcp = src + s->pr_height * src_linesize;
  372.  
  373.             for (y = 0; y < vpad; y++) {
  374.                 memcpy(dstp, srcp, inlink->w * 3);
  375.                 dstp += dst_linesize;
  376.                 srcp += src_linesize;
  377.             }
  378.         }
  379.  
  380.         av_frame_free(&in);
  381.     }
  382.  
  383.     return ff_filter_frame(outlink, out);
  384. }
  385.  
  386. static av_cold void uninit(AVFilterContext *ctx)
  387. {
  388.     int i;
  389.     DCTdnoizContext *s = ctx->priv;
  390.  
  391.     av_dct_end(s->dct);
  392.     av_dct_end(s->idct);
  393.     av_free(s->block);
  394.     av_free(s->tmp_block);
  395.     av_free(s->weights);
  396.     for (i = 0; i < 2; i++) {
  397.         av_free(s->cbuf[i][0]);
  398.         av_free(s->cbuf[i][1]);
  399.         av_free(s->cbuf[i][2]);
  400.     }
  401.     av_expr_free(s->expr);
  402. }
  403.  
  404. static const AVFilterPad dctdnoiz_inputs[] = {
  405.     {
  406.         .name         = "default",
  407.         .type         = AVMEDIA_TYPE_VIDEO,
  408.         .filter_frame = filter_frame,
  409.         .config_props = config_input,
  410.     },
  411.     { NULL }
  412. };
  413.  
  414. static const AVFilterPad dctdnoiz_outputs[] = {
  415.     {
  416.         .name = "default",
  417.         .type = AVMEDIA_TYPE_VIDEO,
  418.     },
  419.     { NULL }
  420. };
  421.  
  422. AVFilter avfilter_vf_dctdnoiz = {
  423.     .name          = "dctdnoiz",
  424.     .description   = NULL_IF_CONFIG_SMALL("Denoise frames using 2D DCT."),
  425.     .priv_size     = sizeof(DCTdnoizContext),
  426.     .init          = init,
  427.     .uninit        = uninit,
  428.     .query_formats = query_formats,
  429.     .inputs        = dctdnoiz_inputs,
  430.     .outputs       = dctdnoiz_outputs,
  431.     .priv_class    = &dctdnoiz_class,
  432.     .flags         = AVFILTER_FLAG_SUPPORT_TIMELINE_GENERIC,
  433. };
  434.