www.pudn.com > bsiftC_.rar > GaussianConvolution.cs

/* GaussianConvolution.cs
 *
 * Gaussian convolution filter, separated passes.
 *
 * (C) Copyright 2004 -- Sebastian Nowozin (nowozin@cs.tu-berlin.de)
 */

// TODO: There might be a way to do this faster by creating a different mask
// and using just one pass from left to right per row or up to down from
// column, instead of (filterSize * dim) passes.

using System;

public class GaussianConvolution
{
	ConvLinearMask mask;

	public GaussianConvolution (double sigma)
		// From "Image Processing, Analysis and Machine Vision", pp. 84:
		// 'Pixels more distant from the center of the operator have smaller
		// influence, and pixels farther than 3 \sigma from the center have
		// neglible influence.'
		//
		// So, find the kernel size by rounding twice 3 \sigma up.
		: this (sigma, 1 + 2 * ((int) (3.0 * sigma)))
	{
	}

	// Like Generate, but manually specifying the kernel size.
	public GaussianConvolution (double sigma, int dim)
	{
		// Assert the kernel size is odd, so we have a clear center pixel.
		dim |= 1;
		mask = new ConvLinearMask (dim);

		double sigma2sq = 2 * sigma * sigma;
		double normalizeFactor = 1.0 / (Math.Sqrt (2.0 * Math.PI) * sigma);

		for (int n = 0 ; n < dim ; ++n) {
			int relPos = n - mask.Middle;

			double G = (relPos * relPos) / sigma2sq;
			G = Math.Exp (-G);
			G *= normalizeFactor;
			mask[n] = G;
			mask.MaskSum += G;
		}
	}

	// Apply the gaussian filter.
	public ImageMap Convolve (ImageMap img)
	{
		return (ConvolutionFilter.Convolve (img, mask));
	}

}

// Static utility functions for convolution
//
public class ConvolutionFilter
{
	public static ImageMap Convolve (ImageMap img, ConvLinearMask mask)
	{
		ImageMap res = new ImageMap (img.XDim, img.YDim);
		ImageMap res2 = new ImageMap (img.XDim, img.YDim);

		Convolve1D (res, mask, img, Direction.Vertical);
		Convolve1D (res2, mask, res, Direction.Horizontal);

		return (res2);
	}

	public enum Direction {
		Vertical,
		Horizontal,
	};

	public static void Convolve1D (ImageMap dest, ConvLinearMask mask,
		ImageMap src, Direction dir)
	{
		int maxN;	// outer loop maximum index
		int maxP;	// inner loop maximum index

		if (dir == Direction.Vertical) {
			maxN = src.XDim;
			maxP = src.YDim;
		} else if (dir == Direction.Horizontal) {
			maxN = src.YDim;
			maxP = src.XDim;
		} else
			throw (new Exception ("TODO: invalid direction"));

		for (int n = 0 ; n < maxN ; ++n) {
			for (int p = 0 ; p < maxP ; ++p) {
				double val = ConvolutionFilter.CalculateConvolutionValue1D (src, mask,
					n, p, maxN, maxP, dir);

				if (dir == Direction.Vertical)
					dest[n, p] = val;
				else
					dest[p, n] = val;
			}
		}
	}

	internal static double CalculateConvolutionValue1D (ImageMap src,
		ConvLinearMask mask, int n, int p, int maxN, int maxP, Direction dir)
	{
		double sum = 0.0;

		bool isOut = false;
		double outBound = 0.0;	// values that are out of bound

		for (int xw = 0 ; xw < mask.Count ; ++xw) {
			int curAbsP = xw - mask.Middle + p;

			if (curAbsP < 0 || curAbsP >= maxP) {
				isOut = true;
				outBound += mask[xw];

				continue;
			}

			if (dir == Direction.Vertical)
				sum += mask[xw] * src[n, curAbsP];
			else
				sum += mask[xw] * src[curAbsP, n];
		}

		// if part of the mask was outside, correct the resulting value by the
		// in/out ratio.
		if (isOut)
			sum *= 1.0 / (1.0 - outBound);

		return (sum);
	}
}

public class ConvLinearMask
{
	int dim;
	public int Count {
		get {
			return (dim);
		}
	}
	int middle;
	public int Middle {
		get {
			return (middle);
		}
	}
	double[] mask;
	public double this[int idx] {
		get {
			return (mask[idx]);
		}
		set {
			mask[idx] = value;
		}
	}
	double maskSum;
	public double MaskSum {
		get {
			return (maskSum);
		}
		set {
			maskSum = value;
		}
	}

	private ConvLinearMask ()
	{
	}

	public ConvLinearMask (int dim)
	{
		mask = new double[dim];
		this.dim = dim;
		this.middle = dim / 2;
	}
}