When building circle blur profile evaluate kernel vertically once per column

src/effects/GrCircleBlurFragmentProcessor.cpp

 /*
  * Copyright 2015 Google Inc.
  *
  * Use of this source code is governed by a BSD-style license that can be
  * found in the LICENSE file.
  */
 
 #include "GrCircleBlurFragmentProcessor.h"
 
 #if SK_SUPPORT_GPU
@@ skipping 118 matching lines @@
     }
     float sum = 0.f;
     // The half kernel should sum to 0.5 not 1.0.
     tot *= 2.f;
     for (int i = 0; i < halfKernelSize; ++i) {
         halfKernel[i] /= tot;
         sum += halfKernel[i];
         summedHalfKernel[i] = sum;
     }
 }
 

[robertphillips, 2016/05/20 14:44:18]

at a points ?

Maybe something more like:

Create a table to store the ...

?

Should the name of this method be changed?

[bsalomon, 2016/05/20 15:17:14]

Done.

-// Applies the 1D half kernel vertically at a point (x, 0) to a circle centered at the origin with
+// Applies the 1D half kernel vertically at a points (x, 0) to a circle centered at the origin with
 // radius circleR.
-static float eval_vertically(float x, float circleR, const float* summedHalfKernelTable,
-                             int halfKernelSize) {
-    // Given x find the positive y that is on the edge of the circle.
-    float y = sqrtf(fabs(circleR * circleR - x * x));
-    // In the column at x we exit the circle at +y and -y
-    // table entry j is actually the kernel evaluated at j + 0.5.
-    y -= 0.5f;
-    int yInt = SkScalarFloorToInt(y);
-    SkASSERT(yInt >= -1);
-    if (y < 0) {
-        return (y + 0.5f) * summedHalfKernelTable[0];
-    } else if (yInt >= halfKernelSize - 1) {
-        return 0.5f;
-    } else {
-        float yFrac = y - yInt;
-        return (1.f - yFrac) * summedHalfKernelTable[yInt] +
-                       yFrac * summedHalfKernelTable[yInt + 1];
+void apply_kernel_in_y(float* results, int numSteps, float firstX, float circleR,
+                       int halfKernelSize, const float* summedHalfKernelTable) {
+    float x = firstX;
+    for (int i = 0; i < numSteps; ++i, x += 1.f) {
+        if (x < -circleR || x > circleR) {

[bsalomon, 2016/05/20 00:24:45]

I tried breaking this loop up into five parts:

1) x fully out of the circle to the left
2) x inside the circle, half kernel extends vertically outside the circle
3) x inside the circle, half kernel contained in y
4) 2 again
5) x fully outside the circle to the right

However, I wasn't able to measure a speedup even with nullgpu.

+            results[i] = 0;
+            continue;
+        }
+        float y = sqrtf(circleR * circleR - x * x);
+        // In the column at x we exit the circle at +y and -y
+        // The summed table entry j is actually reflects an offset of j + 0.5.
+        y -= 0.5f;
+        int yInt = SkScalarFloorToInt(y);
+        SkASSERT(yInt >= -1);
+        if (y < 0) {
+            results[i] = (y + 0.5f) * summedHalfKernelTable[0];
+        } else if (yInt >= halfKernelSize - 1) {
+            results[i] = 0.5f;
+        } else {
+            float yFrac = y - yInt;
+            results[i] = (1.f - yFrac) * summedHalfKernelTable[yInt] +
+                         yFrac * summedHalfKernelTable[yInt + 1];
+        }
     }
 }
 

[robertphillips, 2016/05/20 14:44:19]

// Conceptually combine the 2D Gaussian kernel positioned at (evalX, 0) with a
circle centered at the origin with radius circleR. In reality, combine the look
ups in the halfKernel table (for X) with look ups from the yKernelEvaluations
(for Y).

[bsalomon, 2016/05/20 15:17:14]

Done.

-// Apply the kernel at point (t, 0) to a circle centered at the origin with radius circleR.
-static uint8_t eval_at(float t, float circleR, const float* halfKernel,
-                       const float* summedHalfKernelTable, int halfKernelSize) {
+// Apply the kernel at point (evalX, 0) to a circle centered at the origin with radius circleR.
+static uint8_t eval_at(float evalX, float circleR, const float* halfKernel, int halfKernelSize,
+                       const float* yKernelEvaluations) {
     float acc = 0;
 
-    for (int i = 0; i < halfKernelSize; ++i) {
-        float x = t - i - 0.5f;
+    float x = evalX - halfKernelSize;
+    for (int i = 0; i < halfKernelSize; ++i, x += 1.f) {
         if (x < -circleR || x > circleR) {
             continue;
         }
-        float verticalEval = eval_vertically(x, circleR, summedHalfKernelTable, halfKernelSize);
-        acc += verticalEval * halfKernel[i];
+        float verticalEval = yKernelEvaluations[i];
+        acc += verticalEval * halfKernel[halfKernelSize - i - 1];
     }
-    for (int i = 0; i < halfKernelSize; ++i) {
-        float x = t + i + 0.5f;
+    for (int i = 0; i < halfKernelSize; ++i, x += 1.f) {
         if (x < -circleR || x > circleR) {
             continue;
         }
-        float verticalEval = eval_vertically(x, circleR, summedHalfKernelTable, halfKernelSize);
+        float verticalEval = yKernelEvaluations[i + halfKernelSize];
         acc += verticalEval * halfKernel[i];
     }
     // Since we applied a half kernel in y we multiply acc by 2 (the circle is symmetric about the
     // x axis).
     return SkUnitScalarClampToByte(2.f * acc);
 }
 
 static inline void compute_profile_offset_and_size(float circleR, float sigma,
                                                    float* offset, int* size) {
     if (3*sigma <= circleR) {
         // The circle is bigger than the Gaussian. In this case we know the interior of the
         // blurred circle is solid.
         *offset = circleR - 3 * sigma; // This location maps to 0.5f in the weights texture.
                                        // It should always be 255.
         *size = SkScalarCeilToInt(6*sigma);
     } else {
         // The Gaussian is bigger than the circle.
         *offset = 0.0f;
         *size = SkScalarCeilToInt(circleR + 3*sigma);
     }
 }
 

[robertphillips, 2016/05/20 14:44:19]

Seems like this might need updating

[bsalomon, 2016/05/20 15:17:14]

Done.

 // This function creates a profile of a blurred circle. It does this by computing a kernel for
 // half the Gaussian and a matching summed area table. To compute a profile value at x = r it steps
 // outward in x from (r, 0) in both directions. There is a step for each direction for each entry
 // in the half kernel. The y contribution at each step is computed from the summed area table using
 // the height of the circle above the step point. Each y contribution is multiplied by the half
 // kernel value corresponding to the step in x.
 static uint8_t* create_profile(float circleR, float sigma) {
     float offset;
     int numSteps;
     compute_profile_offset_and_size(circleR, sigma, &offset, &numSteps);
 
     uint8_t* weights = new uint8_t[numSteps];
 
     // The full kernel is 6 sigmas wide.
     int halfKernelSize = SkScalarCeilToInt(6.0f*sigma);
     // round up to next multiple of 2 and then divide by 2
     halfKernelSize = ((halfKernelSize + 1) & ~1) >> 1;
-    SkAutoTArray<float> halfKernel(halfKernelSize);
-    SkAutoTArray<float> summedKernel(halfKernelSize);
-    make_half_kernel_and_summed_table(halfKernel.get(), summedKernel.get(), halfKernelSize,
-                                      sigma);
+
+    // Number of x steps at which to apply kernel in y to cover all the profile samples in x.
+    int numYSteps = numSteps + 2 * halfKernelSize;
+
+    SkAutoTArray<float> bulkAlloc(halfKernelSize + halfKernelSize + numYSteps);
+    float* halfKernel = bulkAlloc.get();
+    float* summedKernel = bulkAlloc.get() + halfKernelSize;
+    float* yEvals = bulkAlloc.get() + 2 * halfKernelSize;
+    make_half_kernel_and_summed_table(halfKernel, summedKernel, halfKernelSize, sigma);
+
+    float firstX = offset - halfKernelSize + 0.5;
+    apply_kernel_in_y(yEvals, numYSteps, firstX, circleR, halfKernelSize, summedKernel);
+
     for (int i = 0; i < numSteps - 1; ++i) {
-        weights[i] = eval_at(offset+i, circleR, halfKernel.get(), summedKernel.get(),
-                             halfKernelSize);
+        float evalX = offset + i + 0.5f;
+        weights[i] = eval_at(evalX, circleR, halfKernel, halfKernelSize, yEvals + i);
     }
     // Ensure the tail of the Gaussian goes to zero.
     weights[numSteps - 1] = 0;
     return weights;
 }
 
 GrTexture* GrCircleBlurFragmentProcessor::CreateCircleBlurProfileTexture(
                                                                 GrTextureProvider* textureProvider,
                                                                 const SkRect& circle,
                                                                 float sigma,
@@ skipping 35 matching lines @@
 GR_DEFINE_FRAGMENT_PROCESSOR_TEST(GrCircleBlurFragmentProcessor);
 
 const GrFragmentProcessor* GrCircleBlurFragmentProcessor::TestCreate(GrProcessorTestData* d) {
     SkScalar wh = d->fRandom->nextRangeScalar(100.f, 1000.f);
     SkScalar sigma = d->fRandom->nextRangeF(1.f,10.f);
     SkRect circle = SkRect::MakeWH(wh, wh);
     return GrCircleBlurFragmentProcessor::Create(d->fContext->textureProvider(), circle, sigma);
 }
 
 #endif