📄 fusionalgrithm.cpp
字号:
// FusionAlgrithm.cpp: implementation of the CFusionAlgrithm class.
//
//////////////////////////////////////////////////////////////////////
#include "stdafx.h"
#include "FusionAlgrithm.h"
#include "FeatureAlgrithm.h"
#include "dibapi.h"
#include "math.h"
#ifdef _DEBUG
#undef THIS_FILE
static char THIS_FILE[]=__FILE__;
#define new DEBUG_NEW
#endif
//////////////////////////////////////////////////////////////////////
// Construction/Destruction
//////////////////////////////////////////////////////////////////////
CFusionAlgrithm::CFusionAlgrithm()
{
}
CFusionAlgrithm::~CFusionAlgrithm()
{
}
HDIB CFusionAlgrithm::HISFUSION(HDIB HCOLOR, HDIB HPAN)
{
/////////////////进行HIS变换时当前图像为空图像（即输出图像）/////////////////////////
LPSTR muldib,pandib,outdib; //指向图像的指针
BYTE* colorimage;//指向彩色图像的数据区位置
BYTE* panimage;//指向全色图像的数据区位置
muldib= (LPSTR) ::GlobalLock((HGLOBAL) HCOLOR);
colorimage=(BYTE*)FindDIBBits(muldib);
pandib= (LPSTR) ::GlobalLock((HGLOBAL) HPAN);
panimage =(BYTE*)FindDIBBits(pandib);
int wide,height,DibWidth,WIDE;
wide=::DIBWidth(muldib);//取宽的像素个数
height=::DIBHeight(muldib);
DibWidth=3*wide+( sizeof(DWORD)- 3*wide %sizeof(DWORD) )%sizeof(DWORD);//取得原图的每行总位数
WIDE=((wide*8)+31)/32*4;//取宽为4的整数倍
BYTE * outimage; //输出图像的数据区位置
HDIB HOUT;
HOUT=(HDIB)::GlobalAlloc(GHND,sizeof(BITMAPINFOHEADER)+DibWidth*height);//生成一个新的图像句柄
if(!HOUT){ return NULL; }
outdib = (LPSTR)::GlobalLock((HGLOBAL)HOUT);
memcpy(outdib,muldib,sizeof(BITMAPINFOHEADER));//sizeof(BITMAPINFOHEADER)=40
outimage= (BYTE*)FindDIBBits (outdib);
if(muldib==NULL || pandib==NULL)
{
AfxMessageBox("必须选择两幅影像!",MB_ICONSTOP);
return NULL;
}
if(muldib==pandib)
{
AfxMessageBox("不可以出现两幅相同的影像!",MB_ICONSTOP);
return NULL;
}
if(::DIBWidth(muldib)!=::DIBWidth(pandib) || ::DIBHeight(muldib)!= ::DIBHeight(pandib))
{
AfxMessageBox("两幅影像的宽度高度必须相等!",MB_ICONSTOP);
return NULL;
}
float * HIS;//用于HIS空间
float * Temp;
float * Im;
float * R;
float * G;
float * B;
HIS = new float[DibWidth*height] ;//开辟一个新内存
Temp = new float[WIDE*height] ;//处理全色图像
Im = new float[WIDE*height] ;//用于存放I分量
B = new float[WIDE*height] ;
R = new float[WIDE*height] ;
G = new float[WIDE*height] ;
//////////////// HIS TRANSFORM ///////////
for(int j=0;j<height;j++)
for (int i=0;i<DibWidth;i+=3)
{
B[j*WIDE+i/3] =colorimage[j*DibWidth+i ]; //得到蓝值
G[j*WIDE+i/3] =colorimage[j*DibWidth+i+1]; //得到绿值
R[j*WIDE+i/3] =colorimage[j*DibWidth+i+2]; //得到红值
float H,S,I;
I=(R[j*WIDE+i/3]+G[j*WIDE+i/3]+B[j*WIDE+i/3])/3;
if (B[j*WIDE+i/3]<=G[j*WIDE+i/3]&&B[j*WIDE+i/3]<=R[j*WIDE+i/3])
{
H=(G[j*WIDE+i/3]-B[j*WIDE+i/3])/((I-B[j*WIDE+i/3])*3);
S=1-(B[j*WIDE+i/3]/I);
}
else if (R[j*WIDE+i/3]<=G[j*WIDE+i/3]&&R[j*WIDE+i/3]<=B[j*WIDE+i/3])
{
H=(B[j*WIDE+i/3]-R[j*WIDE+i/3])/((I-R[j*WIDE+i/3])*3);
S=1-(R[j*WIDE+i/3]/I);
}
else if (G[j*WIDE+i/3]<=B[j*WIDE+i/3]&&G[j*WIDE+i/3]<=R[j*WIDE+i/3])
{
H=(R[j*WIDE+i/3]-G[j*WIDE+i/3])/((I-G[j*WIDE+i/3])*3);
S=1-(G[j*WIDE+i/3]/I);
}
/**/
if (H<0) H=(-1)*H;
if (S<0) S=(-1)*S;
if (I<0) I=(-1)*I;
HIS[i+j*DibWidth] = I;
HIS[i+1+j*DibWidth] = H;
HIS[i+2+j*DibWidth] = S;
Im[j*WIDE+i/3] = I;
Temp[j*WIDE+i/3] = (float)panimage[j*WIDE+i/3];
}
/////////////////// 融合 //////////////////
float *imh;
imh=new float[WIDE*height];
HistogramMatch(Temp,Im,WIDE,height,imh);//用I分量图像对全色图像进行直方图均衡
if(imh == NULL)
{
AfxMessageBox("直方图均衡失败!",MB_ICONSTOP);
return NULL;
}
free(Temp);
free(Im) ; //自动分配内存，不用释放
///////////////////// HIS逆变换 /////////////
for(int m=0;m<height;m++)
for (int n=0;n<DibWidth;n+=3)
{
float I,H,S,r,g,b;
I = imh[m*WIDE+n/3] ; //HIS内存赋给I
H = HIS[m*DibWidth+n+1] ; //HIS内存赋给H
S = HIS[m*DibWidth+n+2] ; //HIS内存赋给S
if (B[m*WIDE+n/3]<=G[m*WIDE+n/3]&&B[m*WIDE+n/3]<=R[m*WIDE+n/3])
{
b=I*(1-S);
g=3*H*(I-b)+b;
r=3*I-g-b;
}
else if (R[m*WIDE+n/3]<=G[m*WIDE+n/3]&&R[m*WIDE+n/3]<=B[m*WIDE+n/3])
{
r=I*(1-S);
b=3*H*(I-r)+r;
g=3*I-r-b;
}
else if(G[m*WIDE+n/3]<=B[m*WIDE+n/3]&&G[m*WIDE+n/3]<=R[m*WIDE+n/3])
{
g=I*(1-S);
r=3*H*(I-g)+g;
b=3*I-r-g;
}
if (b>=255) b=255;
else if (b<=0) b=0;
else b=b;
if (r>=255) r=255;
else if (r<=0) r=0;
else r=r;
if (g>=255) g=255;
else if (g<=0) g=0;
else g=g;
HIS[m*DibWidth+n ] = b ;//得到的b赋给蓝值
HIS[m*DibWidth+n+1] = g ;//得到的g赋给绿值
HIS[m*DibWidth+n+2] = r ;//得到的r赋给红值
}
///////////////// OUTIMAGE中输出融合图像/////
for(int u=0;u<height;u++)
for (int v=0;v<DibWidth;v+=3)
{
outimage[u*DibWidth+v ] = (BYTE)HIS[u*DibWidth+v ];
outimage[u*DibWidth+v+1] = (BYTE)HIS[u*DibWidth+v+1];
outimage[u*DibWidth+v+2] = (BYTE)HIS[u*DibWidth+v+2];
}
free(HIS);
free(R);
free(B);
free(G);
delete []imh;
::GlobalUnlock((HGLOBAL) HCOLOR);
::GlobalUnlock((HGLOBAL) HPAN );
return (HDIB)HOUT;
}
float * CFusionAlgrithm::HistogramMatch(float *img1, float *img2, int width, int height,float *imgResult)
{
//////////////////// 直方图匹配分析 /////////////////////////
//img1全色波段,img2为I分量
if( img1 == NULL || img2 == NULL)
{
AfxMessageBox("获取影像指针失败!",MB_ICONSTOP);
return NULL;
}
static int K = 256;
int i,j,k,l;
double min; BYTE gray;
int *histo1 = new int[K]; int *histo2 = new int[K];
double *p1 = new double[K]; double *p2 = new double[K];
double *f = new double[K]; double *g = new double[K];
int *cpd = new int[K];//
for(i=0;i<K;i++)
{ histo1[i] = 0; histo2[i] = 0; }
//遍历图象的每一个像素* /
for(i=0;i<height;i++)
for(j=0;j<width;j++)
{
//计算图象的直方图
histo1[(BYTE)img1[i*width+j]]++;
histo2[(BYTE)img2[i*width+j]]++;
}
//计算图象的直方图
for(i=0;i<K;i++)
{
p1[i] = (double)histo1[i] / (double)(width*height);//p1[i]灰度值i的概率
p2[i] = (double)histo2[i] / (double)(width*height);
//直方图均衡化* /
f[i] = g[i] = 0.0;
for(l=0;l<=i;l++)
{
f[i] += p1[l];
g[i] += p2[l];
}
}
//建立直方图1和直方图2的对应关系(使用单映射规则)* /
for(k=0;k<K;k++)
{
min = fabs(f[k] - g[0]); cpd[k] = 0;
for(l=1;l<=k;l++)
if(fabs(f[k] - g[l]) < min)
{ min = fabs(f[k] - g[l]); cpd[k] = l; }
}
//遍历图象的每一个像素* /
//计算直方图变换后的图象* /
for(i=0;i<height;i++)
for(j=0;j<width;j++)
{
//将灰度值赋给新数组的该像素* /
gray = (BYTE)img1[i*width+j];
imgResult[i*width+j] = (float)cpd[gray];
}
delete []histo1;
delete []histo2;
delete []cpd;
delete [] p1;
delete [] p2;
delete [] f;
delete [] g;
return imgResult;
}
HDIB CFusionAlgrithm::YIQFUSION(HDIB HMUL, HDIB HPAN)
{
LPSTR muldib,pandib,outdib; //指向图像的指针
BYTE * colorimage;//指向彩色图像的数据区位置
BYTE * panimage;//指向全色图像的数据区位置
muldib= (LPSTR) ::GlobalLock((HGLOBAL) HMUL);
colorimage=(BYTE *)FindDIBBits(muldib);
pandib= (LPSTR) ::GlobalLock((HGLOBAL) HPAN);
panimage =(BYTE *)FindDIBBits(pandib);
int wide,height,DibWidth,WIDE;
wide=::DIBWidth(muldib);//取宽的像素个数
height=::DIBHeight(muldib);
DibWidth=3*wide+( sizeof(DWORD)- 3*wide %sizeof(DWORD) )%sizeof(DWORD);//取得原图的每行总位数
WIDE=((wide*8)+31)/32*4;//取宽为4的整数倍
BYTE * outimage; //输出图像的数据区位置
HDIB HOUT;
HOUT=(HDIB)::GlobalAlloc(GHND,sizeof(BITMAPINFOHEADER)+DibWidth*height);//生成一个新的图像句柄
if(!HOUT){ return NULL; }
outdib = (LPSTR)::GlobalLock((HGLOBAL)HOUT);
memcpy(outdib,muldib,sizeof(BITMAPINFOHEADER));
outimage= (BYTE *) FindDIBBits (outdib);
if(muldib==NULL || pandib==NULL)
{
AfxMessageBox("必须选择两幅影像!",MB_ICONSTOP);
return NULL;
}
if(muldib==pandib)
{
AfxMessageBox("不可以出现两幅相同的影像!",MB_ICONSTOP);
return NULL;
}
if(::DIBWidth(muldib)!=::DIBWidth(pandib) || ::DIBHeight(muldib)!= ::DIBHeight(pandib))
{
AfxMessageBox("两幅影像的宽度高度必须相等!",MB_ICONSTOP);
return NULL;
}
double * Y;
double * Q;
double * I;
double * Temp;
Y = new double[WIDE*height] ;//用于YIQ空间
Q = new double[WIDE*height] ;//用于YIQ空间
I = new double[WIDE*height] ;//用于存放I分量
Temp = new double[WIDE*height] ;//处理全色图像
for(int i=0;i<height;i++)
{
for(int j=0;j<WIDE;j++)
{
double R,G,B,Ym,Im,Qm;
R = colorimage[i*DibWidth+3*j +2];//2,1,0用来控制RGB分量
G = colorimage[i*DibWidth+3*j +1];
B = colorimage[i*DibWidth+3*j +0];
Ym= 0.299*R+0.587*G+0.114*B;
Im= 0.596*R-0.274*G-0.322*B;
Qm= 0.211*R-0.522*G+0.311*B;
Y[i*WIDE+j] = Ym;
I[i*WIDE+j] = Im;
Q[i*WIDE+j] = Qm;
Temp[(DWORD)j*WIDE+i] = (float)panimage[(DWORD)j*WIDE+i];
}
}
/////////////////// 融合 //////////////////
double *imh;
imh=new double[WIDE*height];
HistogramMatch(Temp,Y,WIDE,height,imh);//对全色图像进行直方图均衡
if(imh == NULL)
{
AfxMessageBox("直方图均衡失败!",MB_ICONSTOP);
return NULL;
}
free(Temp);
free(Y);//自动分配内存，不用释放
///////////////////// YIQ逆变换 /////////////
//用全色图像来代替YIQ空间的Y亮度分量
for(int m=0;m<height;m++)
for (int n=0;n<WIDE;n++)
{
double Yn,In,Qn,r,g,b;
Yn=imh[m*WIDE+n];
In=I[m*WIDE+n];
Qn=Q[m*WIDE+n];
r=Yn+0.956*In+0.623*Qn;
g=Yn-0.272*In-0.648*Qn;
b=Yn-1.105*In+0.705*Qn;
if (b>=256) b=256;
if (b<=0) b=0;
if (r>=256) r=256;
if (r<=0) r=0;
if (g>=256) g=256;
if (g<=0) g=0;
outimage[m*DibWidth+n*3 ]=(BYTE)b;
outimage[m*DibWidth+n*3+1]=(BYTE)g;
outimage[m*DibWidth+n*3+2]=(BYTE)r;
}
free (I);
free (Q);
delete []imh;
return (HDIB)HOUT;
}
double * CFusionAlgrithm::HistogramMatch(double *img1, double *img2, int width, int height,double *imgResult)
{
//////////////////// 直方图匹配分析 /////////////////////////
//img1全色波段,img2为I分量
if( img1 == NULL || img2 == NULL)
{
AfxMessageBox("获取影像指针失败!",MB_ICONSTOP);
return NULL;
}
static int K = 256;
int i,j,k,l;
double min; BYTE gray;
int *histo1 = new int[K]; int *histo2 = new int[K];
double *p1 = new double[K]; double *p2 = new double[K];
double *f = new double[K]; double *g = new double[K];
int *cpd = new int[K];//
for(i=0;i<K;i++)
{ histo1[i] = 0; histo2[i] = 0; }
//遍历图象的每一个像素* /
for(i=0;i<height;i++)
for(j=0;j<width;j++)
{
//计算图象的直方图
histo1[(BYTE)img1[i*width+j]]++;
histo2[(BYTE)img2[i*width+j]]++;
}
//计算图象的直方图
for(i=0;i<K;i++)
{
p1[i] = (double)histo1[i] / (double)(width*height);//p1[i]灰度值i的概率
p2[i] = (double)histo2[i] / (double)(width*height);
//直方图均衡化* /
f[i] = g[i] = 0.0;
for(l=0;l<=i;l++)
{
f[i] += p1[l];
g[i] += p2[l];
}
}
//建立直方图1和直方图2的对应关系(使用单映射规则)* /
for(k=0;k<K;k++)
{
min = fabs(f[k] - g[0]); cpd[k] = 0;
for(l=1;l<=k;l++)
if(fabs(f[k] - g[l]) < min)
{ min = fabs(f[k] - g[l]); cpd[k] = l; }
}
//遍历图象的每一个像素* /
//计算直方图变换后的图象* /
for(i=0;i<height;i++)
for(j=0;j<width;j++)
{
⌨️ 快捷键说明
复制代码
Ctrl + C
搜索代码
Ctrl + F
全屏模式
F11
切换主题
Ctrl + Shift + D
显示快捷键
?
增大字号
Ctrl + =
减小字号
Ctrl + -