首先实现P(x|y)=高斯函数(正态分布),即y是正态分布,是条件,我们把它还原到图像界面I(i,j)每一个像素,是结果。
这个实际就是生成高斯模板,我们前面有3*3的高斯生成算法,但是根本没有办法显示,3*3=9个像素,针尖一样,我们在这个基础上,生成11*11的高斯模板,并放大16倍,生成图像,在界面显示出来,看是否是同心圆,二维高斯像一个大钟,我们用平面平行xy平面,切割大钟很多次,感到很像核磁共振人体切片,把这些等高线,让z=0,在xy平面画出来,验证一下,发现图像的高斯离散型和高斯的连续函数还是有差别(颗粒感很强,我们把模板数据变成了灰度),请看图:
c#完整代码如下:
namespace 画高斯同心圆
{
public partial class Form1 : Form
{
public Form1()
{
InitializeComponent();
}
private void Form1_Load(object sender, EventArgs e)
{
//生成二维高斯滤波系数//
// int nWidowSize = 3;
int nWidowSize = 11;
// double nSigma = 0.4;
double nSigma = 3.3;
// double nSigma = 1.4;
// double* pdKernal_2 = new double[nWidowSize*nWidowSize]; //定义二维高斯核数组
double[] pdKernal_2 = new double[11 * 11];
//double[] pdKernal_2 = new double[3 * 3];
///二维高斯函数公式
// x*x+y*y ///
// -1*-------------- ///
// 1 2*nSigma*nSigma ///
// ---------------- e ///
// 2*pi*nSigma*nSigma ///
///
int nCenter = (nWidowSize) / 2;
for (int i = 0; i < nWidowSize; i++)
{
for (int j = 0; j < nWidowSize; j++)
{
int nDis_x = i - nCenter;
int nDis_y = j - nCenter;
double temp = (nDis_x * nDis_x + nDis_y * nDis_y);
temp = temp / (nSigma * nSigma);
double r = Math.Sqrt(temp);//马氏距离
double 指数 = Math.Exp(-(1 / 2.0) * temp);
double temp2 = 指数 / (2 * 3.1415926 * nSigma * nSigma);
pdKernal_2[i + j * nWidowSize] = (int)(temp2 * 1000);
}
}
//因为是11*11,不妥,X12,所以改为44*44*3*3,最终放大了16倍
byte[] gaos = new byte[11 * 11 ];
for (int j = 0; j < 11; j++)
{
for (int i = 0; i < 11; i = i + 1)
{
int nn = j * 11 + i;
int temp3 = (int)(pdKernal_2[nn] * 10);
gaos[nn] = temp3 > 255 ? (byte)255 : (byte)temp3;//
}
}
int ww = 132 + 44;
int hh = 132 + 44;
//int ww = 132 ;
//int hh = 132 ;
byte[] glob_buffer8 = new byte[ww * hh];
//ZoomNormal(gaos, ref glob_buffer8, 11, 11,
// ref ww, ref hh, 12, 12);
ZoomNormal(gaos, ref glob_buffer8, 11, 11,
ref ww, ref hh, 16, 16);
///测试8位显示
Bitmap cutPic88 = new Bitmap(ww, hh, PixelFormat.Format8bppIndexed);
System.Drawing.Imaging.ColorPalette cp = cutPic88.Palette;
for (int i = 0; i < 256; i++)
{
cp.Entries[i] = Color.FromArgb(255, i, i, i);
}
cutPic88.Palette = cp;
BitmapData _cutPic8 = cutPic88.LockBits(new Rectangle(0, 0, ww, hh), ImageLockMode.ReadWrite,
cutPic88.PixelFormat);
IntPtr ptr00 = _cutPic8.Scan0;//得到首地址
//把cutvalues数组给ptr
System.Runtime.InteropServices.Marshal.Copy(glob_buffer8, 0, ptr00, ww * hh);
cutPic88.UnlockBits(_cutPic8);
pictureBox显示gaos.Image = cutPic88;//pictureBox显示gaos
}
void ZoomNormal(byte[] image0, ref byte[] image1, int w, int h,
ref int outwidth, ref int outheight, double ZoomX, double ZoomY)
{
outwidth = (int)(w * ZoomX + 0.5);
outheight = (int)(h * ZoomY + 0.5);
int newbuffersize = outwidth * outheight ;
image1 = new byte[newbuffersize];
int x = 0;
int y = 0;
int tempY;
int tempJ;
for (int j = 0; j < outheight; j++)
{
y = (int)(j / ZoomY + 0.5);
if (y >= h) y--;
tempY = y * w ;
tempJ = j * outwidth ;
for (int i = 0; i < outwidth; i++)
{
x = (int)(i / ZoomX + 0.5);
if (x >= w)
x--;
image1[tempJ + i ] = image0[tempY + x ];
}
}
}
private void pictureBox显示gaos_Paint(object sender, PaintEventArgs e)
{
Graphics g = e.Graphics;//马氏距离r=1.21,70.10*10E-4(概率),位置(1,5),放大16倍后如图
g.DrawEllipse(new Pen(Color.Red), new Rectangle( 16, 16, 4 * 16 * 2, 4 * 16 * 2));
//马氏距离r=0.9,96.67*10E-4(概率),位置(2,5),放大16倍后如图
g.DrawEllipse(new Pen(Color.Red), new Rectangle(16*2, 16*2, 3 * 16 * 2, 3 * 16 * 2));
}
}
}
说明:当 double nSigma = 3.3;更改为:double nSigma = 1;图像变成如下:
这是什么意思呢?也就是说3*nSigma+1=3*1+1=4单元,模板矩阵4*4没有中心,所以模板变成5*5,对比前边,5*5之外,高斯模板就不起作用了,为什么我们nSigma=3.3?原因在此,3.3*3+1=11,所以我们使用11*11的模板。3*3的模板,nSigma=0.4.
另外,当nSigma=1时,重点是:协方差矩阵=协方差矩阵的逆==I=单位矩阵,马氏距离的平方=(x y),正态分布变成了Gaos(x,y)=0.5*(1/3.1416)*exp((-0.5)(x y)),这里有一个非凡的发现,就是中的1变成了特征值,什么是特征值,就是椭圆的长短轴a,b,为什么这样说?因为(x y)=+=f(x,y),这个曲线是什么?同心圆,椭圆长短轴为1的特殊现象,即/1+/1=f(x,y),若f(x,y)为任意实数Z,则/Z+/Z=1,像极了椭圆。
那么怎样让同心圆变成同心椭圆?改变协方差矩阵,比如,改为,这里一个小小的改变,简直是天大的变化,不信请看协方差矩阵的定义,P(协方差系数)不等于零了,即变量x和y不独立了,相关了,但这仍然是高斯分布,不再是同心圆,而是同心椭圆了。
回到我们的程序中的协方差矩阵=,特征值a+b=A+C=矩阵的迹,a*b=矩阵的行列式值=A*C-B*B,可以解得a=3.3,b=3.3,长轴等于短轴,圆一个,变化为=,可以解得a=3.3+0.5=3.8,b=3.3-0.5=2.8;长短轴不等,圆回归常态椭圆。
下一节,把这个同心椭圆画出来。