Revert "Gauss filter calculation"
This reverts commit 53ec7dc7cb523f220a9f5cd713b241c706779c81.
Reason for revert: Segv on very specific machines.
Original change's description:
> Gauss filter calculation
>
> Change-Id: I921ef815d4f788c312aa729f353b6ea154140555
> Reviewed-on: https://skia-review.googlesource.com/67723
> Commit-Queue: Herb Derby <herb@google.com>
> Reviewed-by: Robert Phillips <robertphillips@google.com>
TBR=herb@google.com,robertphillips@google.com
Change-Id: I15164809d081dee0076e815b40fbfdbc6374cfba
No-Presubmit: true
No-Tree-Checks: true
No-Try: true
Reviewed-on: https://skia-review.googlesource.com/69641
Reviewed-by: Herb Derby <herb@google.com>
Commit-Queue: Herb Derby <herb@google.com>
diff --git a/gn/core.gni b/gn/core.gni
index 358c387..21dcff4 100644
--- a/gn/core.gni
+++ b/gn/core.gni
@@ -134,8 +134,6 @@
"$_src/core/SkFindAndPlaceGlyph.h",
"$_src/core/SkArenaAlloc.cpp",
"$_src/core/SkArenaAlloc.h",
- "$_src/core/SkGaussFilter.cpp",
- "$_src/core/SkGaussFilter.h",
"$_src/core/SkFlattenable.cpp",
"$_src/core/SkFlattenableSerialization.cpp",
"$_src/core/SkFont.cpp",
diff --git a/gn/tests.gni b/gn/tests.gni
index 3f8283d..b051326 100644
--- a/gn/tests.gni
+++ b/gn/tests.gni
@@ -212,7 +212,6 @@
"$_tests/SkColor4fTest.cpp",
"$_tests/SkDOMTest.cpp",
"$_tests/SkFixed15Test.cpp",
- "$_tests/SkGaussFilterTest.cpp",
"$_tests/SkImageTest.cpp",
"$_tests/SkLiteDLTest.cpp",
"$_tests/SkNxTest.cpp",
diff --git a/src/core/SkGaussFilter.cpp b/src/core/SkGaussFilter.cpp
deleted file mode 100644
index 548ff43..0000000
--- a/src/core/SkGaussFilter.cpp
+++ /dev/null
@@ -1,152 +0,0 @@
-/*
- * Copyright 2017 Google Inc.
- *
- * Use of this source code is governed by a BSD-style license that can be
- * found in the LICENSE file.
- */
-
-
-#include "SkGaussFilter.h"
-
-#include <cmath>
-#include "SkTypes.h"
-
-static constexpr double kPi = 3.14159265358979323846264338327950288;
-
-// The value when we can stop expanding the filter. The spec implies that 3% is acceptable, but
-// we just use 1%.
-static constexpr double kGoodEnough = 1.0 / 100.0;
-
-// Normalize the values of gauss to 1.0, and make sure they add to one.
-// NB if n == 1, then this will force gauss[0] == 1.
-static void normalize(int n, double* gauss) {
- // Carefully add from smallest to largest to calculate the normalizing sum.
- double sum = 0;
- for (int i = n-1; i >= 1; i--) {
- sum += 2 * gauss[i];
- }
- sum += gauss[0];
-
- // Normalize gauss.
- for (int i = 0; i < n; i++) {
- gauss[i] /= sum;
- }
-
- // The factors should sum to 1. Take any remaining slop, and add it to gauss[0]. Add the
- // values in such a way to maintain the most accuracy.
- sum = 0;
- for (int i = n - 1; i >= 1; i--) {
- sum += 2 * gauss[i];
- }
-
- gauss[0] = 1 - sum;
-}
-
-static int calculate_bessel_factors(double sigma, double *gauss) {
- auto var = sigma * sigma;
-
- // The two functions below come from the equations in "Handbook of Mathematical Functions"
- // by Abramowitz and Stegun. Specifically, equation 9.6.10 on page 375. Bessel0 is given
- // explicitly as 9.6.12
- // BesselI_0 for 0 <= sigma < 2.
- // NB the k = 0 factor is just sum = 1.0.
- auto besselI_0 = [](double t) -> double {
- auto tSquaredOver4 = t * t / 4.0;
- auto sum = 1.0;
- auto factor = 1.0;
- auto k = 1;
- // Use a variable number of loops. When sigma is small, this only requires 3-4 loops, but
- // when sigma is near 2, it could require 10 loops. The same holds for BesselI_1.
- while(factor > 1.0/1000000.0) {
- factor *= tSquaredOver4 / (k * k);
- sum += factor;
- k += 1;
- }
- return sum;
- };
- // BesselI_1 for 0 <= sigma < 2.
- auto besselI_1 = [](double t) -> double {
- auto tSquaredOver4 = t * t / 4.0;
- auto sum = t / 2.0;
- auto factor = sum;
- auto k = 1;
- while (factor > 1.0/1000000.0) {
- factor *= tSquaredOver4 / (k * (k + 1));
- sum += factor;
- k += 1;
- }
- return sum;
- };
-
- // The following formula for calculating the Gaussian kernel is from
- // "Scale-Space for Discrete Signals" by Tony Lindeberg.
- // gauss(n; var) = besselI_n(var) / (e^var)
- auto d = std::exp(var);
- double b[6] = {besselI_0(var), besselI_1(var)};
- gauss[0] = b[0]/d;
- gauss[1] = b[1]/d;
-
- int n = 1;
- // The recurrence relation below is from "Numerical Recipes" 3rd Edition.
- // Equation 6.5.16 p.282
- while (gauss[n] > kGoodEnough) {
- b[n+1] = -(2*n/var) * b[n] + b[n-1];
- gauss[n+1] = b[n+1] / d;
- n += 1;
- }
-
- normalize(n, gauss);
-
- return n;
-}
-
-static int calculate_gauss_factors(double sigma, double* gauss) {
- SkASSERT(0 <= sigma && sigma < 2);
-
- // From the SVG blur spec: 8.13 Filter primitive <feGaussianBlur>.
- // H(x) = exp(-x^2/ (2s^2)) / sqrt(2π * s^2)
- auto var = sigma * sigma;
- auto expGaussDenom = -2 * var;
- auto normalizeDenom = std::sqrt(2 * kPi) * sigma;
-
- // Use the recursion relation from "Incremental Computation of the Gaussian" by Ken
- // Turkowski in GPUGems 3. Page 877.
- double g0 = 1.0 / normalizeDenom;
- double g1 = std::exp(1.0 / expGaussDenom);
- double g2 = g1 * g1;
-
- gauss[0] = g0;
- g0 *= g1;
- g1 *= g2;
- gauss[1] = g0;
-
- int n = 1;
- while (gauss[n] > kGoodEnough) {
- g0 *= g1;
- g1 *= g2;
- gauss[n+1] = g0;
- n += 1;
- }
-
- normalize(n, gauss);
-
- return n;
-}
-
-SkGaussFilter::SkGaussFilter(double sigma, Type type) {
- SkASSERT(0 <= sigma && sigma < 2);
-
- if (type == Type::Bessel) {
- fN = calculate_bessel_factors(sigma, fBasis);
- } else {
- fN = calculate_gauss_factors(sigma, fBasis);
- }
-}
-
-int SkGaussFilter::filterDouble(double* values) const {
- for (int i = 0; i < fN; i++) {
- values[i] = fBasis[i];
- }
- return fN;
-}
-
diff --git a/src/core/SkGaussFilter.h b/src/core/SkGaussFilter.h
deleted file mode 100644
index 9af45c8..0000000
--- a/src/core/SkGaussFilter.h
+++ /dev/null
@@ -1,41 +0,0 @@
-/*
- * Copyright 2017 Google Inc.
- *
- * Use of this source code is governed by a BSD-style license that can be
- * found in the LICENSE file.
- */
-
-#ifndef SkGaussFilter_DEFINED
-#define SkGaussFilter_DEFINED
-
-#include <cstdint>
-
-// Define gaussian filters for values of sigma < 2. Produce values good to 1 part in 1,000,000.
-// Gaussian produces values as defined in the SVG 1.1 spec:
-// https://www.w3.org/TR/SVG/filters.html#feGaussianBlurElement
-// Bessel produces values as defined in "Scale-Space for Discrete Signals" by Tony Lindeberg
-class SkGaussFilter {
-public:
- enum class Type : bool {
- Gaussian,
- Bessel
- };
-
- // Type selects which method is used to calculate the gaussian factors.
- SkGaussFilter(double sigma, Type type);
-
- int radius() const { return fN - 1; }
- int width() const { return 2 * this->radius() + 1; }
-
- // Take an array of values where the gaussian factors will be placed. Return the number of
- // values filled.
- int filterDouble(double values[5]) const;
-
-private:
- double fBasis[5];
- int fN;
-};
-
-#endif // SkGaussFilter_DEFINED
-
-
diff --git a/tests/SkGaussFilterTest.cpp b/tests/SkGaussFilterTest.cpp
deleted file mode 100644
index 958e811..0000000
--- a/tests/SkGaussFilterTest.cpp
+++ /dev/null
@@ -1,88 +0,0 @@
-/*
- * Copyright 2017 Google Inc.
- *
- * Use of this source code is governed by a BSD-style license that can be
- * found in the LICENSE file.
- */
-
-#include "SkGaussFilter.h"
-
-#include <cmath>
-#include <tuple>
-#include <vector>
-#include "Test.h"
-
-// one part in a million
-static constexpr double kEpsilon = 0.000001;
-
-static double careful_add(int n, double* gauss) {
- // Sum smallest to largest to retain precision.
- double sum = 0;
- for (int i = n - 1; i >= 1; i--) {
- sum += 2.0 * gauss[i];
- }
- sum += gauss[0];
- return sum;
-}
-
-DEF_TEST(SkGaussFilterCommon, r) {
- using Test = std::tuple<double, SkGaussFilter::Type, std::vector<double>>;
-
- auto golden_check = [&](const Test& test) {
- double sigma; SkGaussFilter::Type type; std::vector<double> golden;
- std::tie(sigma, type, golden) = test;
- SkGaussFilter filter{sigma, type};
- double result[5];
- size_t n = filter.filterDouble(result);
- REPORTER_ASSERT(r, n == golden.size());
- double sum = careful_add(n, result);
- REPORTER_ASSERT(r, sum == 1.0);
- for (size_t i = 0; i < golden.size(); i++) {
- REPORTER_ASSERT(r, std::abs(golden[i] - result[i]) < kEpsilon);
- }
- };
-
- // The following two sigmas account for about 85% of all sigmas used for masks.
- // Golden values generated using Mathematica.
- auto tests = {
- // 0.788675 - most common mask sigma.
- // GaussianMatrix[{{Automatic}, {.788675}}, Method -> "Gaussian"]
- Test{0.788675, SkGaussFilter::Type::Gaussian, {0.506205, 0.226579, 0.0203189}},
-
- // GaussianMatrix[{{Automatic}, {.788675}}]
- Test{0.788675, SkGaussFilter::Type::Bessel, {0.593605, 0.176225, 0.0269721}},
-
- // 1.07735 - second most common mask sigma.
- // GaussianMatrix[{{Automatic}, {1.07735}}, Method -> "Gaussian"]
- Test{1.07735, SkGaussFilter::Type::Gaussian, {0.376362, 0.244636, 0.0671835}},
-
- // GaussianMatrix[{{4}, {1.07735}}, Method -> "Bessel"]
- Test{1.07735, SkGaussFilter::Type::Bessel, {0.429537, 0.214955, 0.059143, 0.0111337}},
- };
-
- for (auto& test : tests) {
- golden_check(test);
- }
-}
-
-DEF_TEST(SkGaussFilterSweep, r) {
- // The double just before 2.0.
- const double maxSigma = nextafter(2.0, 0.0);
- for (auto type : {SkGaussFilter::Type::Gaussian, SkGaussFilter::Type::Bessel}) {
-
- auto check = [&](double sigma) {
- SkGaussFilter filter{sigma, type};
- double result[5];
- int n = filter.filterDouble(result);
- REPORTER_ASSERT(r, n <= 5);
- double sum = careful_add(n, result);
- REPORTER_ASSERT(r, sum == 1.0);
- };
-
- for (double sigma = 0.0; sigma < 2.0; sigma += 0.1) {
- check(sigma);
- }
-
- check(maxSigma);
- }
-}