isodata.cpp

一个模式识别的源代码,新手可以看看

// ISODATA.cpp : Defines the entry point for the console application.
//

#include "stdafx.h"
#include "stdio.h"
#include "math.h"
#define  N 10
#define  eps   0.00001

struct Pointf
{
	int sequence;
	float x;
	float y;
};

struct PointZ
{
	float x;
	float y;
};

float CalDistancef(Pointf x1,Pointf x2)
{
	return sqrtf((x1.x-x2.x)*(x1.x-x2.x)+(x1.y-x2.y)*(x1.y-x2.y));
}

float CalDistanceZ(PointZ x1,PointZ x2)
{
	return sqrtf((x1.x-x2.x)*(x1.x-x2.x)+(x1.y-x2.y)*(x1.y-x2.y));
}

float CalDistancefZ(Pointf x1,PointZ x2)
{
	return sqrtf((x1.x-x2.x)*(x1.x-x2.x)+(x1.y-x2.y)*(x1.y-x2.y));
}




int main(int argc, char* argv[])
{
	Pointf pts[N]={
		{0,0.0,0.0},{1,3.0,8.0},{2,2.0,2.0},{3,1.0,1.0},{4,5.0,3.0},
		{5,4.0,8.0},{6,6.0,3.0},{7,5.0,4.0},{8,6.0,4.0},{9,7.0,5.0},
	};
	int i,j,m;

	printf("样本集为:\n");
	
	for(i=0;i<N;i++)
	{
		printf("X%d(%.1f,%.1f)  ",i,pts[i].x,pts[i].y);
		
		if((i+1)%5==0)
		{
			printf("\n");
		}
		
	}
	printf("\n");
	printf("\n");
	
	int Nc=0;
	printf("please input Nc(0-10): ");
	scanf("%d",&Nc);

	int Z[N];

	for(i=0;i<Nc;i++)
	{
		printf("输入初始第%d聚类中心的序号(0-9):",i);
		scanf("%d",&Z[i]);
	}
	
    int    Nj[N]; //记录每个类中元素的个数
	PointZ ZArray[N];
	Pointf SAArray[N][N];
	float  DjAv[N];
	float  Deltaj[N][2];
	float  Deltajmax[N];
	int  DeltajmaxCor[N];
	float  DAv;
	int    Nreal=N;
	int    count=0;

	float  Dij[N*N/2];
	int    Diji[N];
	int    Dijj[N];
	int    q=0;
	int    p=0;
    float ft;
	int   it;
	int   jt;
	int   flag;
	int ss=0;
	PointZ Ztp;
	PointZ ZArraytp[N];
	int    Nctp;

	char ch;
	
	int cur=0;

	for(i=0;i<N;i++)
	{
		Nj[i]=0; 
	}

	//聚类中心的特征值
	for(i=0;i<Nc;i++)
	{
		int ihere=Z[i];
		ZArray[i].x=pts[ihere].x;
		ZArray[i].y=pts[ihere].y;
	}

	int K,ThetaN;
	float ThetaS,ThetaC;
	int L,I;

Step1:
	printf("输入预期聚类中心数目 K :");
	scanf("%d",&K);
	printf("输入每个聚类域中最少的样本数ThetaN: ");
	scanf("%d",&ThetaN);

	

	printf("输入同一聚类域中样本标准差的最大值: ");
	scanf("%f",&ThetaS);

	printf("输入不同聚类域距离最小值: ");
	scanf("%f",&ThetaC);

	printf("输入一次可以合并的聚类中心的最多对数: ");
	scanf("%d",&L);

	printf("输入最大迭代次数: ");
	scanf("%d",&I);
Step2:
	for(i=0;i<Nc;i++)
	{
		Nj[i]=0;
	}

		printf("\n");
		printf("这是第%d次归类\n",count+1);
		
	for(i=0;i<N;i++)                          //将模式样本归类
	{
		if(pts[i].sequence==-1)continue;      //若该点的序号为－1则说明它是被剔除的
		float dis=1.0e+10;
		int xx=0;
		float ftemp;

		for(j=0;j<Nc;j++)                     
		{
			ftemp=CalDistancefZ(pts[i],ZArray[j]);
			if(ftemp<dis||fabs(dis-ftemp)<eps)
			{
				xx=j;
				dis=ftemp;
			}
		}
				
		SAArray[xx][Nj[xx]].x=pts[i].x;
		SAArray[xx][Nj[xx]].y=pts[i].y;
		SAArray[xx][Nj[xx]].sequence=pts[i].sequence;	
		Nj[xx]=Nj[xx]+1;	
	}

	   for(i=0;i<Nc;i++)
	   {
		   printf("第%d个聚类中心是:(%.2f,%.2f)   ",i,ZArray[i].x,ZArray[i].y);
		   printf("包含的元素有:",i);
		   for(j=0;j<Nj[i];j++)
		   {
			   printf(" X%d ",SAArray[i][j].sequence);
		   }
		   printf("\n");
		  // printf("Nj(%d) is %d\n",i,Nj[i]);
	   }

	count++;

Step3:

	for(j=0;j<Nc;j++)               //是否可以去掉一些数据
	{
		if(Nj[j]<ThetaN)
		{
			for(i=0;i<Nj[j];i++)
			{
				pts[SAArray[j][i].sequence].sequence=-1;
			}
			i=j;
			int tr=j;
			Nreal-=Nj[j];
			while(j<Nc-1)
			{
				
				for(m=0;m<Nj[j+1];m++)
				{
					SAArray[j][m].x=SAArray[j+1][m].x;
					SAArray[j][m].y=SAArray[j+1][m].y;
					SAArray[j][m].sequence=SAArray[j+1][m].sequence;
				} 
				j++;			
			}

			while(i<Nc-1)
			{
				Nj[i]=Nj[i+1];
				i++;
			}
			Nc--;
			j=tr;
		}		
	}

Step4:              //修正各聚类中心
	for(j=0;j<Nc;j++)
	{
		float temx=0,temy=0;
		for(i=0;i<Nj[j];i++)
		{
			temx+=SAArray[j][i].x;
			temy+=SAArray[j][i].y;		
		}
		ZArray[j].x=temx/Nj[j];
		ZArray[j].y=temy/Nj[j];
	}

Step5://计算各聚类域中诸样本与聚类中心的平均距离
	float temp=0.0;
	for(j=0;j<Nc;j++)
	{
		for(i=0;i<Nj[j];i++)
		{
			temp+=CalDistancefZ(SAArray[j][i],ZArray[j]);
		}

		DjAv[j]=temp/Nj[j];
		
		temp=0.0;
	}
Step6://计算全部模式样本对应聚类中心的总平均距离
	DAv=0;
	for(j=0;j<Nc;j++)
	{
		DAv+=Nj[j]*DjAv[j];
	}
	DAv/=Nreal;

Step7:

	if(count>=I) goto Step14;

	if(Nc<=K/2)goto Step8;

	if((count%2==0)||Nc>=2*K) goto Step11;
Step8://计算各聚类中样本距离标准差
	
	for(j=0;j<Nc;j++)
	{
		float temx=0.0,temy=0.0;
		for(i=0;i<Nj[j];i++)
		{
           temx+=(SAArray[j][i].x-ZArray[j].x)*(SAArray[j][i].x-ZArray[j].x);
		   temy+=(SAArray[j][i].y-ZArray[j].y)*(SAArray[j][i].y-ZArray[j].y);
		}
        Deltaj[j][0]=sqrtf(temx/Nj[j]);
		Deltaj[j][1]=sqrtf(temy/Nj[j]);	
		temx=0.0,temy=0.0;	
	}	
Step9://求每个标准差向量中的最大分量 
	for(j=0;j<Nc;j++)
	{
		Deltajmax[j]=Deltaj[j][0]>Deltaj[j][1]?Deltaj[j][0]:Deltaj[j][1];
		DeltajmaxCor[j]=Deltaj[j][0]>Deltaj[j][1]?0:1;
	}
Step10://分裂判断和计算
	for(j=0;j<Nc;j++)
	{
		if(Deltajmax[j]>ThetaS)
		{
			if((DjAv[j]>DAv&&Nj[j]>2*(ThetaN+1))||Nc<=K/2)
			{
				float Garma=0.5;
				PointZ Zj1,Zj2;

				if(DeltajmaxCor[j]==0)
				{
					Zj1.x=ZArray[j].x+Deltajmax[j]*Garma;
					Zj1.y=ZArray[j].y;
					Zj2.x=ZArray[j].x-Deltajmax[j]*Garma;
					Zj2.y=ZArray[j].y;
				}
				else if(DeltajmaxCor[j]==1)
				{
					Zj1.x=ZArray[j].x;
					Zj1.y=ZArray[j].y+Deltajmax[j]*Garma;
					Zj2.x=ZArray[j].x;
					Zj2.y=ZArray[j].y-Deltajmax[j]*Garma;
				}
				
                ZArray[j].x=Zj1.x;
				ZArray[j].y=Zj1.y;

				ZArray[Nc].x=Zj2.x;
				ZArray[Nc].y=Zj2.y;
				Nc++;
				goto Step2;
			}
		}
	}

Step11://计算全部聚类中心的距离
	ss=0;
	for(i=0;i<Nc-1;i++)
	{
		for(j=i+1;j<Nc;j++)
		{
			Dij[ss]=CalDistanceZ(ZArray[i],ZArray[j]);
			Diji[ss]=i;
			Dijj[ss]=j;
			ss++;		
		}
	}
		
Step12: //简单起见，只考虑一次只合并一对聚类中心的情况
	    //找出类间距离最小的
	ft=Dij[0];
	it=Diji[0];
	jt=Dijj[0];
	for(i=1;i<ss;i++)
	{
		if(Dij[i]<ft)
		{
			ft=Dij[i];
			it=Diji[i];
			jt=Dijj[i];
		}
	}
Step13:
	if(ft<ThetaC)
	{
		Ztp.x=(Nj[it]*ZArray[it].x+Nj[jt]*ZArray[jt].x)/(Nj[it]+Nj[jt]);
		Ztp.y=(Nj[it]*ZArray[it].y+Nj[jt]*ZArray[jt].y)/(Nj[it]+Nj[jt]);
		ZArray[it].x=Ztp.x;
		ZArray[jt].y=Ztp.y;

		j=jt;
		while(j<Nc-1)
		{
			ZArray[j].x=ZArray[j+1].x;
			ZArray[j].y=ZArray[j+1].y;
			j++;
		}
		j=jt;
		
		while(j<Nc-1)
		{
			Nj[j]=Nj[j+1];
			j++;
		}
		Nc--;
	}
Step14:
	if(count>=I) 
	{
		count=0;	
		printf("\n");
		printf("共分为%d类\n",Nc);
			return 0;
		
	}
	else
	{
		goto Step2;	
	}

}