"Fossies" - the Fresh Open Source Software Archive

Member "rawtherapee-5.7/rtengine/guidedfilter.cc" (10 Sep 2019, 6513 Bytes) of package /linux/misc/rawtherapee-5.7.tar.xz:


As a special service "Fossies" has tried to format the requested source page into HTML format using (guessed) C and C++ source code syntax highlighting (style: standard) with prefixed line numbers and code folding option. Alternatively you can here view or download the uninterpreted source code file. For more information about "guidedfilter.cc" see the Fossies "Dox" file reference documentation and the latest Fossies "Diffs" side-by-side code changes report: 5.6_vs_5.7.

    1 /* -*- C++ -*-
    2  *
    3  *  This file is part of RawTherapee.
    4  *
    5  *  Copyright (c) 2018 Alberto Griggio <alberto.griggio@gmail.com>
    6  *
    7  *  RawTherapee 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 3 of the License, or
   10  *  (at your option) any later version.
   11  *
   12  *  RawTherapee 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
   18  *  along with RawTherapee.  If not, see <https://www.gnu.org/licenses/>.
   19  */
   20 
   21 /**
   22  * This is a Fast Guided Filter implementation, derived directly from the
   23  * pseudo-code of the paper:
   24  *
   25  * Fast Guided Filter
   26  * by Kaiming He, Jian Sun
   27  *
   28  * available at https://arxiv.org/abs/1505.00996
   29  */
   30 
   31 #include "guidedfilter.h"
   32 #include "boxblur.h"
   33 #include "rescale.h"
   34 #include "imagefloat.h"
   35 
   36 namespace rtengine {
   37 
   38 #if 0
   39 #  define DEBUG_DUMP(arr)                                                 \
   40     do {                                                                \
   41         Imagefloat im(arr.width(), arr.height());                      \
   42         const char *out = "/tmp/" #arr ".tif";                     \
   43         for (int y = 0; y < im.getHeight(); ++y) {                      \
   44             for (int x = 0; x < im.getWidth(); ++x) {                   \
   45                 im.r(y, x) = im.g(y, x) = im.b(y, x) = arr[y][x] * 65535.f; \
   46             }                                                           \
   47         }                                                               \
   48         im.saveTIFF(out, 16);                                           \
   49     } while (false)
   50 #else
   51 #  define DEBUG_DUMP(arr)
   52 #endif
   53 
   54 
   55 namespace {
   56 
   57 int calculate_subsampling(int w, int h, int r)
   58 {
   59     if (r == 1) {
   60         return 1;
   61     }
   62     
   63     if (max(w, h) <= 600) {
   64         return 1;
   65     }
   66     
   67     for (int s = 5; s > 0; --s) {
   68         if (r % s == 0) {
   69             return s;
   70         }
   71     }
   72 
   73     return LIM(r / 2, 2, 4);
   74 }
   75 
   76 } // namespace
   77 
   78 
   79 void guidedFilter(const array2D<float> &guide, const array2D<float> &src, array2D<float> &dst, int r, float epsilon, bool multithread, int subsampling)
   80 {
   81 
   82     const int W = src.width();
   83     const int H = src.height();
   84 
   85     if (subsampling <= 0) {
   86         subsampling = calculate_subsampling(W, H, r);
   87     }
   88 
   89     enum Op { MUL, DIVEPSILON, ADD, SUB, ADDMUL, SUBMUL };
   90 
   91     const auto apply =
   92         [=](Op op, array2D<float> &res, const array2D<float> &a, const array2D<float> &b, const array2D<float> &c=array2D<float>()) -> void
   93         {
   94             const int w = res.width();
   95             const int h = res.height();
   96             
   97 #ifdef _OPENMP
   98             #pragma omp parallel for if (multithread)
   99 #endif
  100             for (int y = 0; y < h; ++y) {
  101                 for (int x = 0; x < w; ++x) {
  102                     float r;
  103                     float aa = a[y][x];
  104                     float bb = b[y][x];
  105                     switch (op) {
  106                     case MUL:
  107                         r = aa * bb;
  108                         break;
  109                     case DIVEPSILON:
  110                         r = aa / (bb + epsilon);
  111                         break;
  112                     case ADD:
  113                         r = aa + bb;
  114                         break;
  115                     case SUB:
  116                         r = aa - bb;
  117                         break;
  118                     case ADDMUL:
  119                         r = aa * bb + c[y][x];
  120                         break;
  121                     case SUBMUL:
  122                         r = c[y][x] - (aa * bb);
  123                         break;
  124                     default:
  125                         assert(false);
  126                         r = 0;
  127                         break;
  128                     }
  129                     res[y][x] = r;
  130                 }
  131             }
  132         };
  133 
  134     // use the terminology of the paper (Algorithm 2)
  135     const array2D<float> &I = guide;
  136     const array2D<float> &p = src;
  137     array2D<float> &q = dst;
  138 
  139     const auto f_subsample =
  140         [=](array2D<float> &d, const array2D<float> &s) -> void
  141         {
  142             rescaleBilinear(s, d, multithread);
  143         };
  144 
  145     const auto f_upsample = f_subsample;
  146     
  147     const size_t w = W / subsampling;
  148     const size_t h = H / subsampling;
  149 
  150     AlignedBuffer<float> blur_buf(w * h);
  151     const auto f_mean =
  152         [&](array2D<float> &d, array2D<float> &s, int rad) -> void
  153         {
  154             rad = LIM(rad, 0, (min(s.width(), s.height()) - 1) / 2 - 1);
  155             float **src = s;
  156             float **dst = d;
  157 #ifdef _OPENMP
  158             #pragma omp parallel if (multithread)
  159 #endif
  160             boxblur<float, float>(src, dst, blur_buf.data, rad, rad, s.width(), s.height());
  161         };
  162 
  163     array2D<float> I1(w, h);
  164     array2D<float> p1(w, h);
  165 
  166     f_subsample(I1, I);
  167     f_subsample(p1, p);
  168 
  169     DEBUG_DUMP(I);
  170     DEBUG_DUMP(p);
  171     DEBUG_DUMP(I1);
  172     DEBUG_DUMP(p1);
  173 
  174     float r1 = float(r) / subsampling;
  175 
  176     array2D<float> meanI(w, h);
  177     f_mean(meanI, I1, r1);
  178     DEBUG_DUMP(meanI);
  179 
  180     array2D<float> meanp(w, h);
  181     f_mean(meanp, p1, r1);
  182     DEBUG_DUMP(meanp);
  183 
  184     array2D<float> &corrIp = p1;
  185     apply(MUL, corrIp, I1, p1);
  186     f_mean(corrIp, corrIp, r1);
  187     DEBUG_DUMP(corrIp);
  188 
  189     array2D<float> &corrI = I1;
  190     apply(MUL, corrI, I1, I1);
  191     f_mean(corrI, corrI, r1);
  192     DEBUG_DUMP(corrI);
  193 
  194     array2D<float> &varI = corrI;
  195     apply(SUBMUL, varI, meanI, meanI, corrI);
  196     DEBUG_DUMP(varI);
  197 
  198     array2D<float> &covIp = corrIp;
  199     apply(SUBMUL, covIp, meanI, meanp, corrIp);
  200     DEBUG_DUMP(covIp);
  201 
  202     array2D<float> &a = varI;
  203     apply(DIVEPSILON, a, covIp, varI);
  204     DEBUG_DUMP(a);
  205 
  206     array2D<float> &b = covIp;
  207     apply(SUBMUL, b, a, meanI, meanp);
  208     DEBUG_DUMP(b);
  209 
  210     meanI.free(); // frees w * h * 4 byte
  211     meanp.free(); // frees w * h * 4 byte
  212 
  213     array2D<float> &meana = a;
  214     f_mean(meana, a, r1);
  215     DEBUG_DUMP(meana);
  216 
  217     array2D<float> &meanb = b;
  218     f_mean(meanb, b, r1);
  219     DEBUG_DUMP(meanb);
  220 
  221     blur_buf.resize(0); // frees w * h * 4 byte
  222 
  223     array2D<float> meanA(W, H);
  224     f_upsample(meanA, meana);
  225     DEBUG_DUMP(meanA);
  226 
  227     array2D<float> &meanB = q;
  228     f_upsample(meanB, meanb);
  229     DEBUG_DUMP(meanB);
  230 
  231     apply(ADDMUL, q, meanA, I, meanB);
  232     DEBUG_DUMP(q);
  233 }
  234 
  235 } // namespace rtengine