图像处理 高斯模板生成器

项目中用到不同大小的高斯模板，为了方便自己写了个高斯模板生成器，以下代码实现高斯模板的生成：

/**************************************************
***A program with function to generate gauss kernel
***author: zhang jie
***Date:2014.04.15
**************************************************/
#include<stdio.h>
#include<stdlib.h>
#include<assert.h>
#include<math.h>
#define LENGTH(radius) 2*radius+1
typedef struct
{
	int nrows;
	int ncols;
	int radius; //the length of templates
	float variance; //variance of 2D gauss function
	float ratio; //template ratio
	int** data;
}gausskernel2D,*Gauss_Kernel_2D;
enum ImageType
{
	GAUSS_IMAGE_GAUSS,
	GAUSS_IMAGE_HORIZONTAL,
	GAUSS_IMAGE_VERTICAL
};
void** Create2DArray(int nrows,int ncols,int nBytes)
{
	char** tt;
	int i;
	tt=(char**)malloc(nrows*sizeof(void*)+ncols*nrows*nBytes);
	for(int i=0;i<nrows;i++)
	{
		tt[i]=((char*)tt)+(nrows*sizeof(void*)+i*ncols*nBytes);
	}
	return ((void**)tt);
}
Gauss_Kernel_2D InitGaussKernel(int radius,float variance)
{
	Gauss_Kernel_2D g_kernel_2d;
	/* Allocate memory */
    g_kernel_2d = (Gauss_Kernel_2D)malloc(sizeof(gausskernel2D));
	g_kernel_2d->radius=radius;
	g_kernel_2d->variance=variance;
	g_kernel_2d->nrows=LENGTH(radius);
	g_kernel_2d->ncols=LENGTH(radius);
	g_kernel_2d->ratio=0;
	int nBytes=g_kernel_2d->nrows*g_kernel_2d->ncols*sizeof(int);
	g_kernel_2d->data=(int**)Create2DArray(g_kernel_2d->nrows,g_kernel_2d->ncols,nBytes);
	for(int i=0;i<g_kernel_2d->nrows;i++)
	 for(int j=0;j<g_kernel_2d->ncols;j++)
		 g_kernel_2d->data[i][j]=0;	
	return g_kernel_2d;
}
bool Generate_Gauss_Kernel(Gauss_Kernel_2D kernel,int radius)
{
	  int i,j;
      float t,*Ftemp;
	  Ftemp=(float*)malloc((radius+1)*(radius+1)*sizeof(float));
	  if(!Ftemp) return false;
	  for(i=0;i<radius+1;i++) 
	  for(j=0;j<radius+1;j++)
	  {
		t=(float)(i*i+j*j)/(float)kernel->variance; // smooth scale
		Ftemp[i*(radius+1)+j]=(float)(1.0/exp(t/2)); 
      }
      int C=(int)ceil(1/Ftemp[(radius+1)*(radius+1)-1]+0.5); //compute normalized coefficient
	  /*printf("C=%d\n",C);*/
    
	  for(i=0;i<radius+1;i++) //assign valute to lt part
	  for(j=0;j<radius+1;j++)
      kernel->data[radius+i][radius+j]=(int)(Ftemp[i*(radius+1)+j]*C+0.5); 
	

      for(i=radius;i<2*radius+1;i++) //assign valute to ld part
	  for(j=0;j<radius+1;j++)
	  kernel->data[i][j]=kernel->data[i][2*radius-j];


      for(i=0;i<radius;i++) //assign valute to up  
      for(j=0;j<2*radius+1;j++)
	  kernel->data[i][j]=kernel->data[2*radius-i][j]; 
	
      for(i=0;i<2*radius+1;i++) //compute total coefficient
	  for(j=0;j<2*radius+1;j++)
	  kernel->ratio +=(float )kernel->data[i][j]; 
	  printf("%f\n",kernel->ratio);
      kernel->ratio=(float)( 1.0/kernel->ratio);
	  return true;
}

bool Generate_Gauss_Kernel_Horizontal(Gauss_Kernel_2D kernel,int radius)
{
	  int i,j;
      float t,*Ftemp;
	  Ftemp=(float*)malloc((2*radius+1)*(2*radius+1)*sizeof(float));
	  float den=0.0;
	  for(int i=-radius;i<=radius;i++)
		for(int j=-radius;j<=radius;j++)
		{
		    t=(float)(i*i+j*j)/(float)kernel->variance; // smooth scale
			Ftemp[(j+radius)*(2*radius+1)+(i+radius)]=-i*exp((-t)/2);
			den -=i*Ftemp[(j+radius)*(2*radius+1)+(i+radius)];
		}
	   int C=(int)ceil(1.0/Ftemp[(2*radius+1)*(2*radius+1)-1]+0.5); //compute normalized coefficient
	   printf("C=%d\n",C);
       printf("Ftemp=%f\n",Ftemp[(2*radius+1)*(2*radius+1)-1]);
		/*for(int j=-radius;j<=radius;j++)
		  for(int i=-radius;i<=radius;i++)
		{
			if((i+radius+1)%(2*radius+1)==0)
			printf("%f\n",Ftemp[(j+radius)*(2*radius+1)+(i+radius)]);
			else
            printf("%f ",Ftemp[(j+radius)*(2*radius+1)+(i+radius)]);
			kernel->data[i+radius][j+radius]=(int)(Ftemp[(j+radius)*(2*radius+1)+(i+radius)]/den);
		}*/
	   for(int j=-radius;j<=radius;j++)
		  for(int i=-radius;i<=radius;i++)
		{
			kernel->data[i+radius][j+radius]=(int)(Ftemp[(j+radius)*(2*radius+1)+i+radius]*C+0.5); 
			printf("%d\n",kernel->data[i+radius][j+radius]);
		}
      for(i=0;i<2*radius+1;i++) //compute total coefficient
	  for(j=0;j<2*radius+1;j++)
	  kernel->ratio +=(float)kernel->data[i][j]; 

      kernel->ratio=(float)( 1.0/kernel->ratio);
	  return true;
}

bool Generate_Gauss_Kernel_Vertical(Gauss_Kernel_2D kernel,int radius)
{
	  int i,j;
      float t,*Ftemp;
	  Ftemp=(float*)malloc((2*radius+1)*(2*radius+1)*sizeof(float));
	  float den=0.0;
	  for(int i=-radius;i<=radius;i++)
		for(int j=-radius;j<=radius;j++)
		{
		    t=(float)(i*i+j*j)/(float)kernel->variance; // smooth scale
			Ftemp[(j+radius)*(2*radius+1)+(i+radius)]=-j*exp((-t)/2);
			den -=j*Ftemp[(j+radius)*(2*radius+1)+(i+radius)];
		}
        printf("%f\n",den);
		for(int j=-radius;j<=radius;j++)
		  for(int i=-radius;i<=radius;i++)
		{
			if((i+radius+1)%(2*radius+1)==0)
			printf("%f\n",Ftemp[(j+radius)*(2*radius+1)+(i+radius)]);
			else
            printf("%f ",Ftemp[(j+radius)*(2*radius+1)+(i+radius)]);
			kernel->data[i+radius][j+radius]=(int)(Ftemp[(j+radius)*(2*radius+1)+(i+radius)]/den);
		}
	   for(int i=-radius;i<=radius;i++)
		for(int j=-radius;j<=radius;j++)
		{
			kernel->data[i+radius][j+radius]=(int)(Ftemp[(j+radius)*(2*radius+1)+(i+radius)]/den);
		}
      for(i=0;i<2*radius+1;i++) //compute total coefficient
	  for(j=0;j<2*radius+1;j++)
	  kernel->ratio +=(float)kernel->data[i][j]; 

      kernel->ratio=(float)( 1.0/kernel->ratio);
	  return true;
}

bool Store_Gauss_kernel(Gauss_Kernel_2D kernel,const char* ch)
{
	FILE* fp;
	if((fp=fopen(ch,"w"))==NULL) return false;
	fprintf(fp,"the length of kernel is:%d\n",kernel->nrows);
	fprintf(fp,"the coefficient of kernel is:%f\n",kernel->ratio);
	fprintf(fp,"the variance of kernel is:%f\n",kernel->variance);
	fprintf(fp,"the template data of kernel is:\n");
	for(int i=0;i<kernel->nrows;i++)
	 for(int j=0;j<kernel->ncols;j++)
	 {
		 if((j+1)%kernel->ncols==0) {fprintf(fp,"%d\n",kernel->data[i][j]); printf("%d\n",kernel->data[i][j]);}
		 else                        { fprintf(fp,"%d ",kernel->data[i][j]); printf("%d ",kernel->data[i][j]);}
	 }
	 fclose(fp);
	 return true;
}
int main(int argc,char* argv[])
{
	//Input variance,radius
	float variance;
	int   radius;
	ImageType imageType;
	imageType=GAUSS_IMAGE_VERTICAL;
	printf("Input the variance and radius of gauss kernel:\n");
	scanf("%f,%d",&variance,&radius);
	assert(variance>0);
	assert(radius>0);
	//
	Gauss_Kernel_2D g_kernel_2d=InitGaussKernel(radius,variance);
	bool success;
	switch(imageType)
	{
	case GAUSS_IMAGE_GAUSS: 
	success=Generate_Gauss_Kernel(g_kernel_2d,radius);
	break;
    case GAUSS_IMAGE_HORIZONTAL:
	success=Generate_Gauss_Kernel_Horizontal(g_kernel_2d,radius);
	break;
	case GAUSS_IMAGE_VERTICAL:
	success=Generate_Gauss_Kernel_Vertical(g_kernel_2d,radius);
	break;
	}
	if(!success) printf("faile to generate\n");
	//Store the gauss kernel into file
	bool ret=Store_Gauss_kernel(g_kernel_2d,"zhj.txt");
	//if(!ret) printf("fail to store\n");
	printf("finish \n");
	system("pause");
	//exit(0);
	return 0;
}