Itti算法c代码
//Itti視覺注意模型的opencv實現
//代碼參照Matlab的代碼進行編寫，結構上進行了一定得優化
//Matlab代碼的下載地址見本人CSDN資源
//http://download.csdn.net/detail/megaarthur/7402441#include<opencv2/opencv.hpp>
#include<opencv/highgui.h>
#include<opencv/cv.h>
#include<stdio.h>
#include<iostream>
#include<Windows.h>
#define PI 3.1415926
using namespace std;//高斯金字塔結構體
typedef struct GaussPyr{IplImage *lev[9];}GaussPyr;//將金字塔尺寸定義為全局變量
CvSize PyrSize[9] = {NULL};//初始化金字塔結構體
void initPyr(GaussPyr *p)
{for(int i = 0; i < 9; i++)p->lev[i] = cvCreateImage(PyrSize[i], IPL_DEPTH_64F, 1);
}//根據層數，求第i層的尺寸
int downsample(int x, int level)
{if(level-- > 0){if(x%2 == 0)x = x/2;elsex = (x+1)/2;downsample(x,level);}if(level == -1)return x;
}//計算并產生一幅圖的高斯金字塔 每層的圖像
void Gscale(GaussPyr *p, IplImage *data, int level, double sigma)
{for(int i = 0; i<level; i++){if(i == 0)cvSmooth(data, p->lev[0], CV_GAUSSIAN, 5, 5, sigma, 0);else{IplImage *tem = cvCreateImage(PyrSize[i-1], IPL_DEPTH_64F, 1);cvSmooth(p->lev[i-1], tem, CV_GAUSSIAN, 5, 5, sigma, 0);for(int a = 0; a < PyrSize[i].height; a++)for(int b = 0; b < PyrSize[i].width; b++)((double *)(p->lev[i]->imageData + a*p->lev[i]->widthStep))[b] = ((double *)(tem->imageData + 2*a*tem->widthStep))[2*b];}}
}//c-s過程中用到的跨尺度相減
void overScaleSub(IplImage *s1, IplImage *s2, IplImage *dst)
{cvResize(s2, dst, CV_INTER_LINEAR);cvAbsDiff(s1, dst, dst);
}//求圖像的局部最大值
void getLocalMaxima(IplImage *scr, double thresh, double *lm_sum, int *lm_num, double*lm_avg)
{*lm_sum = 0.0;*lm_num = 0;*lm_avg = 0.0;int count = 0;//查找局部最大值 for(int a = 1; a < ((scr->height) - 1); a++)for(int b = 1; b< ((scr->width) - 1); b++){double val = ((double *)(scr->imageData + a*scr->widthStep))[b];if((val >= thresh) &&(val >= ((double *)(scr->imageData + (a-1)*scr->widthStep))[b]) &&(val >= ((double *)(scr->imageData + (a+1)*scr->widthStep))[b]) &&(val >= ((double *)(scr->imageData + a*scr->widthStep))[b-1]) &&(val >= ((double *)(scr->imageData + a*scr->widthStep))[b+1])){if (val == 10) count++;//因為設定圖像的最大值是10 *lm_sum += val;//局部最大值加上全局最大值(*lm_num) ++;//總的個數 包含局部和全局}}if(*lm_num > count){*lm_sum = *lm_sum - 10*count;//局部最大值的總和*lm_num = *lm_num - count;//局部最大值的個數if(*lm_num > 0)*lm_avg = *lm_sum / *lm_num;//局部最大值的平均值else*lm_avg = 0;}else*lm_avg = 0;
}//N操作，包括下采樣到第5層的尺寸
void N_operation(IplImage *scr, IplImage *dst)
{cvNormalize(scr, scr, 1.0, 0.0, CV_MINMAX, NULL);cvConvertScale(scr, scr, 10, 0);int lm_num;double lm_sum;double lm_avg;getLocalMaxima(scr ,1, &lm_sum, &lm_num, &lm_avg);if(lm_num > 0)cvConvertScale(scr, scr, (10-lm_avg)*(10-lm_avg), 0);//下采樣，先判斷目前的層數if( scr->height == PyrSize[2].height )//scr在第三層{for(int a = 0; a < PyrSize[4].height; a++)for(int b = 0; b < PyrSize[4].width; b++)((double *)(dst->imageData + a*dst->widthStep))[b] = ((double *)(scr->imageData + 4*a*scr->widthStep))[4*b];}else if( scr->height == PyrSize[3].height )//scr在第四層{ for(int a = 0; a < PyrSize[4].height; a++)for(int b = 0; b < PyrSize[4].width; b++)((double *)(dst->imageData + a*dst->widthStep))[b] = ((double *)(scr->imageData + 2*a*scr->widthStep))[2*b];}else //scr在第五層{for(int a = 0; a < PyrSize[4].height; a++)for(int b = 0; b < PyrSize[4].width; b++)((double *)(dst->imageData + a*dst->widthStep))[b] = ((double *)(scr->imageData + a*scr->widthStep))[b];}
}//C_S過程，金字塔最底層為第0層(實際上（c,s）=(3,6),(3,7),(4,7),(4,8),(5,8),(5,9),)
void CS_operation(GaussPyr *p1, GaussPyr *p2, IplImage *dst)
{for(int c = 2; c < 5; c++)for(int delta = 3, s = c + delta; delta < 5; delta++, s = c + delta){IplImage *tem_c = cvCreateImage(PyrSize[c], IPL_DEPTH_64F, 1);overScaleSub(p1->lev[c], p2->lev[s], tem_c);IplImage *tem_5lev = cvCreateImage(PyrSize[4], IPL_DEPTH_64F, 1);N_operation(tem_c, tem_5lev);cvAdd(tem_5lev, dst, dst, NULL);}
}//金字塔相減，在求顏色特征時使用
void PyrSub(GaussPyr *s1, GaussPyr *s2, GaussPyr *dst)
{for(int i = 0; i<9; i++)cvSub(s1->lev[i], s2->lev[i], dst->lev[i], NULL);
}//計算Gabor濾波
void cvGabor(IplImage *scr, IplImage *dst, int width, double lamda, double theta)
{CvMat *gabor_kernel = cvCreateMat(width, width, CV_32FC1);double tmp1,tmp2,xtmp,ytmp,re;int i,j,x,y;for(i = 0; i < width; i++)for(j = 0; j < width; j++){x = (i*16/(width - 1)) - 8;y = (j*16/(width - 1)) - 8;xtmp = (float)x*cos(theta) + (float)y*sin(theta);ytmp = (float)(-x)*sin(theta) + (float)y*cos(theta);tmp1 = (1/(PI*pow(lamda,2)))*exp(-((pow(xtmp,2) + pow(ytmp,2))/pow(lamda,2)));tmp2 = cos(2*PI*xtmp/lamda);re = tmp1*tmp2;cvSetReal2D((CvMat*)gabor_kernel,i,j,re);}cvFilter2D(scr, dst, gabor_kernel, cvPoint(-1, -1));cvAbs(dst, dst);double max = 0;for(int i = 0; i<dst->height; i++)for(int j = 0; j<dst->width; j++)if( ((double *)(dst->imageData + i*dst->widthStep))[j] >= max )max = ((double *)(dst->imageData + i*dst->widthStep))[j];cvConvertScale(dst, dst, 1/max, 0);
}void runItti(const char* filename)
{DWORD time_start, time_end;time_start = GetTickCount();//讀入原圖IplImage *origin = cvLoadImage(filename, CV_LOAD_IMAGE_UNCHANGED);//將原圖插值為兩倍CvSize newsize;newsize.height = 2*origin->height;newsize.width = 2*origin->width;IplImage *originx2 = cvCreateImage(newsize, origin->depth, origin->nChannels);cvResize(origin, originx2 , CV_INTER_LINEAR);//分離bgr分量IplImage *rgb[4] = {NULL};for(int i=0; i < 3; i++)rgb[i] = cvCreateImage(newsize, IPL_DEPTH_8U, 1);cvSplit(originx2, rgb[2], rgb[1], rgb[0], rgb[3]);//計算亮度I和R、G、B、Y分量IplImage *Intensity = cvCreateImage(newsize, IPL_DEPTH_64F, 1);IplImage *R = cvCreateImage(newsize, IPL_DEPTH_64F, 1);IplImage *G = cvCreateImage(newsize, IPL_DEPTH_64F, 1);IplImage *B = cvCreateImage(newsize, IPL_DEPTH_64F, 1);IplImage *Y = cvCreateImage(newsize, IPL_DEPTH_64F, 1);uchar *data_r = (uchar *)rgb[0]->imageData;uchar *data_g = (uchar *)rgb[1]->imageData;uchar *data_b = (uchar *)rgb[2]->imageData;for(int i=0; i < newsize.height; i++)for(int j=0; j < newsize.width; j++){int index_8 = i*(rgb[0]->widthStep) + j;((double *)(Intensity->imageData + i*Intensity->widthStep))[j] = (double(data_r[index_8]) + double(data_g[index_8]) + double(data_b[index_8]))/3;((double *)(R->imageData + i*Intensity->widthStep))[j] = double(data_r[index_8]) - (double(data_g[index_8]) + double(data_b[index_8]))/2;((double *)(G->imageData + i*Intensity->widthStep))[j] = double(data_g[index_8]) - (double(data_r[index_8]) + double(data_b[index_8]))/2;((double *)(B->imageData + i*Intensity->widthStep))[j] = double(data_b[index_8]) - (double(data_r[index_8]) + double(data_g[index_8]))/2;((double *)(Y->imageData + i*Intensity->widthStep))[j] = (double(data_r[index_8]) + double(data_g[index_8]))/2 - (abs(double(data_r[index_8]) - double(data_g[index_8])))/2 - double(data_b[index_8]);}//releasecvReleaseImage(&originx2); cvReleaseImage(&rgb[0]); cvReleaseImage(&rgb[1]);cvReleaseImage(&rgb[2]);cvReleaseImage(&rgb[3]);//初始化高斯金字塔尺寸for(int i = 0; i < 9; i++){PyrSize[i].height = downsample(newsize.height, i);PyrSize[i].width = downsample(newsize.width, i);}亮度特征提取//計算I金字塔GaussPyr Pyr_I;initPyr(&Pyr_I);Gscale(&Pyr_I, Intensity, 9, 1.6);//I的CS過程IplImage *I_mean = cvCreateImage(PyrSize[4], IPL_DEPTH_64F, 1);cvSetZero(I_mean);CS_operation(&Pyr_I, &Pyr_I, I_mean);N_operation(I_mean, I_mean);//releasefor(int i = 0; i < 9; i++){cvReleaseImage(&(Pyr_I.lev[i]));}方向特征提取//計算O金字塔IplImage *orient[4];for(int i = 0; i < 4; i++)orient[i] = cvCreateImage(cvGetSize(Intensity), IPL_DEPTH_64F,1);cvGabor(Intensity, orient[0], 11, 5, 0);cvGabor(Intensity, orient[1], 11, 5, PI/4);cvGabor(Intensity, orient[2], 11, 5, PI/2);cvGabor(Intensity, orient[3], 11, 5, 3*PI/4);//releasecvReleaseImage(&Intensity);GaussPyr Pyr_0, Pyr_45, Pyr_90, Pyr_135;initPyr(&Pyr_0);initPyr(&Pyr_45);initPyr(&Pyr_90);initPyr(&Pyr_135);Gscale(&Pyr_0, orient[0], 9, 0.5);Gscale(&Pyr_45, orient[1], 9, 0.5);Gscale(&Pyr_90, orient[2], 9, 0.5);Gscale(&Pyr_135, orient[3], 9, 0.5);//releasecvReleaseImage(&orient[0]); cvReleaseImage(&orient[1]); cvReleaseImage(&orient[2]); cvReleaseImage(&orient[3]);//O的CS過程IplImage *O_mean = cvCreateImage(PyrSize[4], IPL_DEPTH_64F, 1);IplImage *O_1 = cvCreateImage(PyrSize[4], IPL_DEPTH_64F, 1);IplImage *O_2 = cvCreateImage(PyrSize[4], IPL_DEPTH_64F, 1);IplImage *O_3 = cvCreateImage(PyrSize[4], IPL_DEPTH_64F, 1);cvSetZero(O_mean);cvSetZero(O_1);cvSetZero(O_2);cvSetZero(O_3);CS_operation(&Pyr_0, &Pyr_0, O_mean);CS_operation(&Pyr_45, &Pyr_45, O_1);CS_operation(&Pyr_90, &Pyr_90, O_2);CS_operation(&Pyr_135, &Pyr_135, O_3);//releasefor(int i = 0; i < 9; i++){cvReleaseImage(&(Pyr_0.lev[i]));cvReleaseImage(&(Pyr_45.lev[i]));cvReleaseImage(&(Pyr_90.lev[i]));cvReleaseImage(&(Pyr_135.lev[i]));}N_operation(O_mean, O_mean);N_operation(O_1, O_1);N_operation(O_2, O_2);N_operation(O_3, O_3);cvAdd(O_mean, O_1, O_mean, NULL);cvAdd(O_mean, O_2, O_mean, NULL);cvAdd(O_mean, O_3, O_mean, NULL);N_operation(O_mean, O_mean);//releasecvReleaseImage(&O_1); cvReleaseImage(&O_2); cvReleaseImage(&O_3); 顏色特征提取//計算RGBY金字塔GaussPyr Pyr_R, Pyr_G, Pyr_B, Pyr_Y;initPyr(&Pyr_R);initPyr(&Pyr_G);initPyr(&Pyr_B);initPyr(&Pyr_Y);Gscale(&Pyr_R, R, 9, 0.5);Gscale(&Pyr_G, G, 9, 0.5);Gscale(&Pyr_B, B, 9, 0.5);Gscale(&Pyr_Y, Y, 9, 0.5);//releasecvReleaseImage(&R); cvReleaseImage(&G); cvReleaseImage(&B); cvReleaseImage(&Y); //計算RG—BY金字塔GaussPyr Pyr_RG, Pyr_GR, Pyr_BY, Pyr_YB;initPyr(&Pyr_RG);initPyr(&Pyr_GR);initPyr(&Pyr_BY);initPyr(&Pyr_YB);PyrSub(&Pyr_R, &Pyr_G, &Pyr_RG);PyrSub(&Pyr_G, &Pyr_R, &Pyr_GR);PyrSub(&Pyr_B, &Pyr_Y, &Pyr_BY);PyrSub(&Pyr_Y, &Pyr_B, &Pyr_YB);//releasefor(int i = 0; i < 9; i++){cvReleaseImage(&(Pyr_R.lev[i]));cvReleaseImage(&(Pyr_G.lev[i]));cvReleaseImage(&(Pyr_B.lev[i]));cvReleaseImage(&(Pyr_Y.lev[i]));}//C的CS過程IplImage *C_mean = cvCreateImage(PyrSize[4], IPL_DEPTH_64F, 1);IplImage *C_1 = cvCreateImage(PyrSize[4], IPL_DEPTH_64F, 1);cvSetZero(C_mean);cvSetZero(C_1);CS_operation(&Pyr_RG, &Pyr_GR, C_mean);CS_operation(&Pyr_BY, &Pyr_YB, C_1);cvAdd(C_mean, C_1, C_mean, NULL);//releasefor(int i = 0; i < 9; i++){cvReleaseImage(&(Pyr_RG.lev[i]));cvReleaseImage(&(Pyr_GR.lev[i]));cvReleaseImage(&(Pyr_BY.lev[i]));cvReleaseImage(&(Pyr_YB.lev[i]));}cvReleaseImage(&C_1);N_operation(C_mean, C_mean);//整合所有特征，生成顯著性圖IplImage *all = cvCreateImage(PyrSize[4], IPL_DEPTH_64F, 1);//在第四層生成的顯著圖 融合后的顯著圖cvSetZero(all);cvAdd(I_mean, C_mean, all, NULL);cvAdd(all, O_mean, all, NULL);//releasecvReleaseImage(&I_mean); cvReleaseImage(&C_mean); cvReleaseImage(&O_mean); cvConvertScale(all, all, 0.33333, 0);cvNormalize(all, all, 1.0, 0.0, CV_MINMAX, NULL);IplImage *all_8U = cvCreateImage(cvGetSize(all), IPL_DEPTH_8U, 1);cvConvertScale(all, all_8U, 255, 0);//releasecvReleaseImage(&all);IplImage *saliency = cvCreateImage(cvSize(origin->width, origin->height), IPL_DEPTH_8U, 1);//最后調整到原圖像大小cvResize(all_8U, saliency, CV_INTER_LINEAR);//releasecvReleaseImage(&all_8U);time_end = GetTickCount();cvNamedWindow("原圖", CV_WINDOW_AUTOSIZE);cvNamedWindow("顯著性圖", CV_WINDOW_AUTOSIZE);cvShowImage("原圖", origin);cvShowImage("顯著性圖", saliency);cout << "總用時：" << time_end - time_start << "ms" << endl;cvWaitKey(0);cvReleaseImage(&origin);cvReleaseImage(&saliency);cvDestroyWindow("原圖"); cvDestroyWindow("顯著性圖"); }int main()
{string loop = "y";while(loop == "y"){cout << "輸入圖片名稱：" << endl;string originFILENAME = "NULL";cin >> originFILENAME;const char *filename = originFILENAME.data();runItti(filename);cout << "繼續？ y/n" << endl;cin >> loop;}return 0;}
總結
