// Copyright (c) 2013 The Chromium Authors. All rights reserved. // Use of this source code is governed by a BSD-style license that can be // found in the LICENSE file. #include #include #include #include "base/check_op.h" #include "base/notreached.h" #include "skia/ext/recursive_gaussian_convolution.h" namespace skia { namespace { // Takes the value produced by accumulating element-wise product of image with // a kernel and brings it back into range. // All of the filter scaling factors are in fixed point with kShiftBits bits of // fractional part. template inline unsigned char FloatTo8(float f) { int a = static_cast(f + 0.5f); if (take_absolute) a = std::abs(a); else if (a < 0) return 0; if (a < 256) return a; return 255; } template inline float ForwardFilter(float in_n_1, float in_n, float in_n1, const std::vector& w, int n, const float* b) { switch (order) { case RecursiveFilter::FUNCTION: return b[0] * in_n + b[1] * w[n-1] + b[2] * w[n-2] + b[3] * w[n-3]; case RecursiveFilter::FIRST_DERIVATIVE: return b[0] * 0.5f * (in_n1 - in_n_1) + b[1] * w[n-1] + b[2] * w[n-2] + b[3] * w[n-3]; case RecursiveFilter::SECOND_DERIVATIVE: return b[0] * (in_n - in_n_1) + b[1] * w[n-1] + b[2] * w[n-2] + b[3] * w[n-3]; } NOTREACHED(); return 0.0f; } template inline float BackwardFilter(const std::vector& out, int n, float w_n, float w_n1, const float* b) { switch (order) { case RecursiveFilter::FUNCTION: case RecursiveFilter::FIRST_DERIVATIVE: return b[0] * w_n + b[1] * out[n + 1] + b[2] * out[n + 2] + b[3] * out[n + 3]; case RecursiveFilter::SECOND_DERIVATIVE: return b[0] * (w_n1 - w_n) + b[1] * out[n + 1] + b[2] * out[n + 2] + b[3] * out[n + 3]; } NOTREACHED(); return 0.0f; } template unsigned char SingleChannelRecursiveFilter( const unsigned char* const source_data, int source_pixel_stride, int source_row_stride, int row_width, int row_count, unsigned char* const output, int output_pixel_stride, int output_row_stride, const float* b) { const int intermediate_buffer_size = row_width + 6; std::vector w(intermediate_buffer_size); const unsigned char* in = source_data; unsigned char* out = output; unsigned char max_output = 0; for (int r = 0; r < row_count; ++r, in += source_row_stride, out += output_row_stride) { // Compute forward filter. // First initialize start of the w (temporary) vector. if (order == RecursiveFilter::FUNCTION) w[0] = w[1] = w[2] = in[0]; else w[0] = w[1] = w[2] = 0.0f; // Note that special-casing of w[3] is needed because of derivatives. w[3] = ForwardFilter( in[0], in[0], in[source_pixel_stride], w, 3, b); int n = 4; int c = 1; int byte_index = source_pixel_stride; for (; c < row_width - 1; ++c, ++n, byte_index += source_pixel_stride) { w[n] = ForwardFilter(in[byte_index - source_pixel_stride], in[byte_index], in[byte_index + source_pixel_stride], w, n, b); } // The value of w corresponding to the last image pixel needs to be computed // separately, again because of derivatives. w[n] = ForwardFilter(in[byte_index - source_pixel_stride], in[byte_index], in[byte_index], w, n, b); // Now three trailing bytes set to the same value as current w[n]. w[n + 1] = w[n]; w[n + 2] = w[n]; w[n + 3] = w[n]; // Now apply the back filter. float w_n1 = w[n + 1]; int output_index = (row_width - 1) * output_pixel_stride; for (; c >= 0; output_index -= output_pixel_stride, --c, --n) { float w_n = BackwardFilter(w, n, w[n], w_n1, b); w_n1 = w[n]; w[n] = w_n; out[output_index] = FloatTo8(w_n); max_output = std::max(max_output, out[output_index]); } } return max_output; } unsigned char SingleChannelRecursiveFilter( const unsigned char* const source_data, int source_pixel_stride, int source_row_stride, int row_width, int row_count, unsigned char* const output, int output_pixel_stride, int output_row_stride, const float* b, RecursiveFilter::Order order, bool absolute_values) { if (absolute_values) { switch (order) { case RecursiveFilter::FUNCTION: return SingleChannelRecursiveFilter( source_data, source_pixel_stride, source_row_stride, row_width, row_count, output, output_pixel_stride, output_row_stride, b); case RecursiveFilter::FIRST_DERIVATIVE: return SingleChannelRecursiveFilter< RecursiveFilter::FIRST_DERIVATIVE, true>( source_data, source_pixel_stride, source_row_stride, row_width, row_count, output, output_pixel_stride, output_row_stride, b); case RecursiveFilter::SECOND_DERIVATIVE: return SingleChannelRecursiveFilter< RecursiveFilter::SECOND_DERIVATIVE, true>( source_data, source_pixel_stride, source_row_stride, row_width, row_count, output, output_pixel_stride, output_row_stride, b); } } else { switch (order) { case RecursiveFilter::FUNCTION: return SingleChannelRecursiveFilter( source_data, source_pixel_stride, source_row_stride, row_width, row_count, output, output_pixel_stride, output_row_stride, b); case RecursiveFilter::FIRST_DERIVATIVE: return SingleChannelRecursiveFilter< RecursiveFilter::FIRST_DERIVATIVE, false>( source_data, source_pixel_stride, source_row_stride, row_width, row_count, output, output_pixel_stride, output_row_stride, b); case RecursiveFilter::SECOND_DERIVATIVE: return SingleChannelRecursiveFilter< RecursiveFilter::SECOND_DERIVATIVE, false>( source_data, source_pixel_stride, source_row_stride, row_width, row_count, output, output_pixel_stride, output_row_stride, b); } } NOTREACHED(); return 0; } } float RecursiveFilter::qFromSigma(float sigma) { DCHECK_GE(sigma, 0.5f); if (sigma <= 2.5f) return 3.97156f - 4.14554f * std::sqrt(1.0f - 0.26891f * sigma); return 0.98711f * sigma - 0.96330f; } void RecursiveFilter::computeCoefficients(float q, float b[4]) { b[0] = 1.57825f + 2.44413f * q + 1.4281f * q * q + 0.422205f * q * q * q; b[1] = 2.4413f * q + 2.85619f * q * q + 1.26661f * q * q * q; b[2] = - 1.4281f * q * q - 1.26661f * q * q * q; b[3] = 0.422205f * q * q * q; // The above is exactly like in the paper. To cut down on computations, // we can fix up these numbers a bit now. float b_norm = 1.0f - (b[1] + b[2] + b[3]) / b[0]; b[1] /= b[0]; b[2] /= b[0]; b[3] /= b[0]; b[0] = b_norm; } RecursiveFilter::RecursiveFilter(float sigma, Order order) : order_(order), q_(qFromSigma(sigma)) { computeCoefficients(q_, b_); } unsigned char SingleChannelRecursiveGaussianX(const unsigned char* source_data, int source_byte_row_stride, int input_channel_index, int input_channel_count, const RecursiveFilter& filter, const SkISize& image_size, unsigned char* output, int output_byte_row_stride, int output_channel_index, int output_channel_count, bool absolute_values) { return SingleChannelRecursiveFilter(source_data + input_channel_index, input_channel_count, source_byte_row_stride, image_size.width(), image_size.height(), output + output_channel_index, output_channel_count, output_byte_row_stride, filter.b(), filter.order(), absolute_values); } unsigned char SingleChannelRecursiveGaussianY(const unsigned char* source_data, int source_byte_row_stride, int input_channel_index, int input_channel_count, const RecursiveFilter& filter, const SkISize& image_size, unsigned char* output, int output_byte_row_stride, int output_channel_index, int output_channel_count, bool absolute_values) { return SingleChannelRecursiveFilter(source_data + input_channel_index, source_byte_row_stride, input_channel_count, image_size.height(), image_size.width(), output + output_channel_index, output_byte_row_stride, output_channel_count, filter.b(), filter.order(), absolute_values); } } // namespace skia