sift算法—图像的高斯金字塔的生成（1）

代码并没有优化，只是草草生成，

#include "stdafx.h"

#include<iostream>
#include<Windows.h>
#include"highgui.h"
#include"cv.h"
#include"math.h" 
#include<vector>
#include<string>
using namespace std;
//使用floor函数。floor(x)返回的是小于或等于x的最大整数。
#define pi 3.1415
#define S   3 //高斯差分金字塔的组内层数
//#define dem_gauss 3 
//高斯模板的维数
//高斯模板
vector<CvMat*> Image_Hash;
double** gauss_model(double sigma)
{
//int dem_gauss = ceil(6 * sigma + 1);//使用ceil函数。ceil(x)返回的是大于x的最小整数。
int dem_gauss = int(6 * sigma + 1);//使用ceil函数。ceil(x)返回的是大于x的最小整数。
if (dem_gauss % 2 == 0)
{
dem_gauss++;
}
int mid = dem_gauss / 2;
double guiyi_total=0;
double** template_su = new double*[dem_gauss];
for (int i = 0; i < dem_gauss; i++)
{
template_su[i] = new double[dem_gauss];
}
for (int m = 0; m < dem_gauss; m++)
{
for (int n = 0; n < dem_gauss; n++)
{
template_su[m][n] = 1 / (2 * pi*sigma*sigma);
template_su[m][n] = template_su[m][n] * exp((-pow((mid - m), 2) - pow((mid - n), 2)) / (2 * pow(sigma, 2)));
guiyi_total = guiyi_total + template_su[m][n];
}
}
cout<<endl;
//归一化
for (int i = 0; i < dem_gauss; i++)
{
for (int j = 0; j < dem_gauss; j++)
{
template_su[i][j] = template_su[i][j] / guiyi_total;
}
}
//for (int k = 0; k < dem_gauss; k++)
//{ 
// for (int v = 0; v < dem_gauss; v++)
// {
// cout<<'\t' << &template_su[k][v] <<",";
// }
// cout << endl;
//}
return template_su;
}
//图像的灰度化
CvMat* grey_img(IplImage* image)
{
IplImage  *GrayImage;
CvMat* pGrayMat = NULL;
GrayImage = cvCreateImage(cvGetSize(image), 8, 1);
pGrayMat = cvCreateMat(image->height, image->width, CV_32FC1);
BYTE data_b;   //BYTE类型的头文件为windows.h，相当于unsigned char
BYTE data_g;
BYTE data_r;
BYTE data_gray;
for (int j = 0; j<image->height; j++)
{
for (int i = 0; i<image->width; i++)
{
data_b = (BYTE)image->imageData[j*image->widthStep + i * 3];     //B分量  
data_g = (BYTE)image->imageData[j*image->widthStep + i * 3 + 1]; //G分量  
data_r = (BYTE)image->imageData[j*image->widthStep + i * 3 + 2]; //R分量
//Gray= 0.072169B+ 0.715160G+ 0.212671R
data_gray = (BYTE)(0.072169*data_b + 0.715160*data_g + 0.212671*data_r);
cvmSet(pGrayMat, j, i, data_gray);
}
}
cvConvert(pGrayMat, GrayImage);
cvNamedWindow("Image", CV_WINDOW_AUTOSIZE);
cvShowImage("Image", GrayImage);
//cvWaitKey(0);
//cvReleaseImage(&GrayImage);
//cvDestroyWindow("Image_1");
return pGrayMat;
}
 //图像高斯化
CvMat* gauss_image(CvMat *Gray_array, double sigma_float)
{
IplImage  *GrayImage;
CvMat* pGrayMat = NULL;
GrayImage = cvCreateImage(cvGetSize(Gray_array), 8, 1);
pGrayMat = cvCreateMat(Gray_array->height, Gray_array->width, CV_32FC1);


int dem_gauss = int(6 * sigma_float+1);//使用ceil函数。ceil(x)返回的是大于x的最小整数。
if (dem_gauss%2==0)
{
dem_gauss++;
}
int mid_num = dem_gauss / 2;
double** template_gauss = NULL;
//double** template_gauss = new double*[dem_gauss];
//for (int i = 0; i < dem_gauss; i++)
//{
// template_gauss[i] = new double[dem_gauss];
//}/


template_gauss=gauss_model(sigma_float);
//cout << template_gauss[0][0] << endl;
//cout << &template_gauss;
// for (int k = 0; k < dem_gauss; k++)
// {
// for (int v = 0; v < dem_gauss; v++)
// {
// cout << '\t' << template_gauss[k][v] << ",";
// }
// cout << endl;
// }
int gray_height = Gray_array->height;
int gray_width = Gray_array->width;
//double filter_num = 0.0;
//double filter_sum = 0.0;
for (int i = 0; i < gray_height; i++)
{
for (int j = 0; j < gray_width; j++)
{
double dFilter = 0.0;
double dSum = 0.0;
int a1, a2;
for (int x = -mid_num; x <= mid_num; x++)
{
for (int y = -mid_num; y <= mid_num; y++)
{
if ((j + x) >= 0 && (j + x)<gray_width && (i + y) >= 0 && (i + y)<gray_height)  //判断边缘  
{
double Imgaedata = cvmGet(Gray_array, i + y, j + x);//访问矩阵元素:
//cout << Imgaedata << endl;
///int y1 = (y + mid_num)*dem_gauss;
//int x1=x + mid_num;
//cout << template_gauss[y1][x1];
//cout << Imgaedata << endl;           dem_gauss
//a1=*(*(template_gauss+(y + mid_num))+(x + mid_num));
//a2 = Imgaedata;
//dFilter += a1*a2;
dFilter += Imgaedata *template_gauss[(y + mid_num) ][x + mid_num];
dSum += template_gauss[(y + mid_num)][x + mid_num];
}
}
}
cvmSet(pGrayMat, i, j, dFilter / dSum);
}

}
//cvConvert(pGrayMat, GrayImage);
//cvNamedWindow("Image_2", CV_WINDOW_AUTOSIZE);
//cvShowImage("Image_2", GrayImage);
//cvWaitKey(0);
//cvReleaseImage(&GrayImage);
//cvDestroyWindow("Image_2");
//释放内存
for (int b = 0; b < dem_gauss; b++)
{
delete[] template_gauss[b];
}
delete[] template_gauss;
return pGrayMat;
}
//构建高斯金字塔
CvMat* convert(CvMat*  Gray_array)
{
CvMat* pGrayMat = NULL;
pGrayMat = cvCreateMat(Gray_array->height, Gray_array->width, CV_32FC1);
int gray_height = Gray_array->height;
int gray_width = Gray_array->width;

for (int i = 0; i < gray_height; i++)
{
for (int j = 0; j < gray_width; j++)
{
int  v = cvmGet(Gray_array, i, j);
cvmSet(pGrayMat, i, j, v);
}
}


return pGrayMat;
}
//图像降采样函数
CvMat* downsample(CvMat* pMat)
{
int dowidth = (int)(pMat->width/2);
int doheight = (int)(pMat->height/2);
  CvMat* cpMat = NULL;
cpMat = cvCreateMat(doheight, dowidth, CV_32FC1);
for (int i = 0; i < doheight; i++)
{
for (int j = 0; j < dowidth; j++)
{
int  y = 2 * i;
int x = 2 * j;
cvmSet(cpMat, i, j, cvmGet(pMat, y, x));
}
}
return cpMat;
}
void Hashgausss(CvMat*  Gray_array, int minSize, double sigma_float)
{

int nheight = Gray_array->height;
int nwidth = Gray_array->width;
int nsmall = nheight < nwidth ? nheight : nwidth;
int octave = log(nsmall) / log(2) - log(minSize) / log(2);
double* sigma_array = new double[S + 3];
sigma_array[0] = sigma_float;
//int v = pow(2.0, 1.0 / S);
for (int k = 1; k < S + 3; k++)
{
double sigma_pre = pow(2,(double)(k - 1)*1.0/S)*sigma_float;
double sigma_tai = sigma_pre*pow(2,(double)(1.0/S));
sigma_array[k] = (double)sqrt((double)abs(sigma_pre *sigma_pre - sigma_tai*sigma_tai));
cout << sigma_array[k] << endl;
}
CvMat* pMat = NULL;
pMat = convert(Gray_array);
for (int m = 0; m < octave; m++)//组数，金字塔层数
{
for (int n = 0; n < S + 3; n++)
{
pMat = gauss_image(pMat, sigma_array[n]);
Image_Hash.push_back(pMat);
}
vector<CvMat*>::iterator iter = Image_Hash.begin();
pMat = *(iter+(Image_Hash.size() - 3));
pMat=convert(pMat);
//降采样
pMat = downsample(pMat);
}
vector<CvMat*>::iterator ster = Image_Hash.begin();

char um = 0;
for (; ster != Image_Hash.end(); ster++)
{

char str[] = "Image_";
char srr[1];
itoa(um, srr, 10);
strcat(str, srr);
um++;
//CvMat* bb = *ster;
//CvSize size = cvSize(bb->width, bb->height);
IplImage  *GrayImage=NULL;
GrayImage = cvCreateImage(cvGetSize(*ster), 8, 1);
cvConvert(*ster, GrayImage);
cvNamedWindow(str, CV_WINDOW_AUTOSIZE);

cvShowImage(str, GrayImage);
cvWaitKey(0);
}
}
int _tmain(int argc, _TCHAR* argv[])
{
IplImage* Img_1;
CvMat*  Gray_array;
vector<CvMat*> Image_Hash;
Img_1 = cvLoadImage("d://psb.jpg");
Gray_array = cvCreateMat(Img_1->height, Img_1->width, CV_32FC1);
if (Img_1 == NULL)
{
cout << "载入图片错误" << endl;
return 0;
}
//图像的灰度化
Gray_array=grey_img(Img_1);
//图像的简单高斯化
//double sigma_float=0.9;
//gauss_image(Gray_array, sigma_float);
Hashgausss(Gray_array,64,0.9);


return 0;
}