matlab中fspecial中生成高斯模板的解释以及c语言实现

下面这个函数是模拟fspecial中的生成高斯模板的函数。

function h=makefilter(size,sigma)
%size为模板大小
%sigma为标准差



%下面的代码其实是从fspecial中摘录出来的，我做了一些更改放到自己写的函数里面便于解释
%计算高斯模板的中心位置
siz    = ([size size]-1)/2;
sig    = sigma;


%用meshgrid是为了加速，不用for循环
[x y]  = meshgrid(-siz(2):siz(2),-siz(1):siz(1));


%计算exp(-(x^2+y^2)/(2*sig^2))
%我想你肯定有一个疑问，那就是为什么不除以2*pi*sig^2
%因为不除也没有关系，因为最后归一化之后会约掉
arg    = -(x.*x+y.*y)./(2*sig*sig);
h      = exp(arg);




%不知道它为什么要这样，忘懂得人解释一下
h(h<eps*max(h(:))) = 0;


%求和，用来归一化
sumh   = sum(h(:));


%防止求和之后出现为0的情况，然后再归一化一下使高斯，模板为小数
if sumh ~=0
    h=h/sumh;
end






end

运行结果：

>> makefilter(3,0.5)

ans =

    0.0113    0.0838    0.0113
    0.0838    0.6193    0.0838
    0.0113    0.0838    0.0113

下面给出C的代码

#include <conio.h>
#include <stdio.h>
#include <math.h>
 
  double** makeGaussianFilter(int iSize, double sigma);
  double** malloc2DArray(int iRow, int iCol);
  void free2DArray(double** p,int iRow);
  int _tmain(int argc, _TCHAR* argv[])
  {
	  double** filter = makeGaussianFilter(3,0.5);
	  for (int i=0; i<3; i++)
	  {
		  for (int j=0; j<3; j++)
		  {
			  printf("%lf\t",filter[i][j]);
		  }
		  printf("\n");
	  }
	  _getch();
	  free2DArray(filter,3);
	  return 0;
  }

  double** makeGaussianFilter(int iSize, double sigma)
  {
	  double** filter = malloc2DArray(iSize,iSize);
	  int center = (iSize-1)/2;
	  double sum = 0;
	  double x2 = 0;
	  double y2 = 0;
	  for (int i=0; i<iSize; i++)
	  {
		  x2 = pow(double(i-center),2);
		  for (int j=0; j<iSize; j++)
		  {
			  y2 = pow(double(j-center),2);
			  sum += filter[i][j] = exp(-(x2+y2)/(2*sigma*sigma));
		  }
	  }


	  if (sum!=0)
	  {
		  //归一化
		  for (int i=0; i<iSize; i++)
		  {
			  for (int j=0; j<iSize; j++)
			  {
				  filter[i][j] /= sum;
			  }
		  }
	  }

	  return filter;
  }


  
  double** malloc2DArray(int iRow, int iCol)
  {
	  double **filter = new double*[iRow];
	  for (int i=0; i<iRow; i++)
	  {
		  filter[i] = new double[iCol];
	  }
	  return filter;
  }
  void free2DArray(double** p,int iRow)
  {
	  for (int i=0; i<iRow; i++)
	  {
		  delete []p[i];
	  }
	  delete []p;
  }

结果：

总结：有时候不一定按照公式的算法去做，可以简化计算。