#include<iostream.h>
#include<math.h>
#include<iomanip.h>
#define n 19//樣品數
#define m 9//樣品屬性個數
void main()
{
	int i,j,k;
	double a[n][m]={{86.76,12.25,0.00,0.40,46.20,5.36,35.55,3.58,1.30},
{92.52,7.32,0.00,0.16,47.52,3.68,37.70,1.95,0.59},
{96.76,3.07,0.17,0.00,47.86,7.82,37.77,1.38,0.63},
{85.56,13.03,0.94,0.47,45.04,7.44,36.06,3.43,0.65},
{87.75,10.26,1.23,0.76,36.22,24.29,29.25,3.23,1.15},
{75.69,24.06,0.00,2.50,47.40,3.76,40.00,0.97,1.15},
{98.21,0.63,0.16,0.00,86.82,18.87,31.75,2.82,0.96},
{84.50,14.46,1.04,0.00,46.94,4.59,37.35,2.20,0.89},
{94.14,5.86,0.00,0.00,48.66,8.41,38.42,0.51,0.67},
{90.50,6.72,2.78,0.00,35.18,31.58,30.11,0.51,0.59},
{72.05,26.49,1.46,0.00,49.04,5.19,39.28,0.05,0.74},
{98.10,1.71,0.00,0.00,54.22,5.79,32.04,1.33,0.63},
{96.35,3.13,0.26,0.00,44.76,6.54,36.20,0.76,0.70},
{98.48,0.00,1.08,0.44,28.60,17.77,24.37,3.37,5.56},
{88.51,2.54,4.33,0.29,49.12,25.27,18.49,2.46,1.07},
{79.75,19.463,0.00,0.00,53.98,8.27,24.66,7.25,2.37},
{94.54,0.76,2.05,0.31,35.64,34.81,11.18,6.50,0.78},
{91.51,7.88,0.61,0.00,46.06,8.42,39.71,0.82,0.74},
{90.89,3.04,0.23,1.17,62.48,5.11,25.73,0.72,0.70}};//原始數據表
	double d[m][m];//距離系數陣
	double R[m][m];//R型聚類
	double max_min[2];//在變量正規化中，存儲a[][]中的最大與最小值
	double t;//在計算距離系數陣中作臨時變量
	double result[m-1][3];//做m-1次聚類，第一、二個存儲聚類對象的號碼，第三個存儲與兩者相關的距離
	int sp=0;//大循環,共聚類m-1次
//	char c[m][50]={"56","83","80","58","79","98","102"};//樣本名稱
	int temp=0;//在制作圖表時作臨時變量
	double sum1,sum2,x_t,x_s;//R型聚類
	for(i=0;i<m;i++)//賦初值
	{
		for(j=0;j<m;j++)
		{
			d[i][j]=0;
		}
	}
	for(i=0;i<m;i++)
	{
		for(j=0;j<m;j++)
		{
			R[i][j]=0;
		}
	}
	//變量正規化
    for(i=0;i<n;i++)//t
	{
		max_min[0]=0;////////////////////////給一個最小值
		max_min[1]=1000;////////////////////給一個最大值
		for(j=0;j<m;j++)//t
		{
			for(k=0;k<n;k++)
			{
				if(a[k][j]>max_min[0])
				{
					max_min[0]=a[k][j];
				}
				if(a[k][j]<max_min[1])
				{
					max_min[1]=a[k][j];
				}
			}
			a[i][j]=(a[i][j]-max_min[1])/(max_min[0]-max_min[1]);
		}	
	}
	cout<<"對樣本數據標準化處理得："<<endl;
	for(i=0;i<n;i++)
	{
		for(j=0;j<m;j++)
		{
			cout<<setw(8)<<setprecision(3)<<a[i][j];
		}
		cout<<endl;
	}
	cout<<endl;
	//R型聚類？？？？？？？？？？？？？？？？？
	for(i=0;i<m;i++)
	{
		for(j=0;j<m;j++)
		{
			sum1=sum2=0;
			for(k=0;k<n;k++)
			{
				sum1+=a[i][k];
				sum2+=a[j][k];
			}
			x_t=sum1/n;
			x_s=sum2/n;
			sum1=sum2=0;
			for(int f=0;f<n;f++)
			{
				R[i][j]+=(a[f][i]-x_t)*(a[f][j]-x_s);
				sum1+=pow((a[f][i]-x_t),2);
				sum2+=pow((a[f][j]-x_s),2);
			}
			R[i][j]=R[i][j]/pow((sum1*sum2),0.5);
		}
	}
	cout<<endl<<"R聚類："<<endl;
	for(i=0;i<m;i++)//輸出R聚類結果
	{
		for(j=0;j<m;j++)
		{
			cout<<setw(8)<<setprecision(3)<<R[i][j];
		}
		cout<<endl;
	}
	for(;sp<m-1;)
	{
		//計算距離系數陣（樣品之間的距離）
		for(i=0;i<n-sp;i++)
		{
			for(j=0;j<n-sp;j++)
			{
				t=0;
				for(k=0;k<m;k++)
				{
					t+=pow((a[k][i]-a[k][j]),2);
				}
				d[i][j]=pow((t/m),0.5);
			}
		}
		cout<<"距離系數陣:"<<endl;
		for(i=0;i<n-sp;i++)
		{
			for(j=0;j<n-sp;j++)
			{
				cout<<setw(8)<<setprecision(3)<<d[i][j];
			}
			cout<<endl;
		}
		// 開始聚類
		//（1）找最小值
		t=1000;///////////////////////////給一個最大值
		for(i=0;i<n-sp;i++)
		{
			for(j=i;j<n-sp;j++)
			{
				if(d[i][j]<t&&d[i][j]!=0)
				{
					t=d[i][j];
					result[sp][0]=i;
					result[sp][1]=j;
				}
			}
		}
		result[sp][2]=t;
		//cout<<endl<<result[0][2]<<endl<<result[0][0]<<result[0][1];
		for(i=0;i<n;i++)//
		{
			a[i][int(result[sp][0])]=(a[i][int(result[sp][0])]+a[i][int(result[sp][1])])/2;
		}
		for(i=0;i<n;i++)
		{
			for(j=int(result[sp][1])+1;j<m-sp;j++)
			{
				a[i][j-1]=a[i][j];
			}
		}
		cout<<endl;
		sp++;
		cout<<endl<<"第"<<sp<<"次聚類得到矩陣："<<endl;
		for(i=0;i<n;i++)
		{
			for(j=0;j<m-sp;j++)
			{
				cout<<setw(8)<<setprecision(3)<<a[i][j];
			}
			cout<<endl;
		}
	}
	//結論
	/*
	cout<<endl<<setw(11)<<"連接順序"<<setw(30)<<"-----------連接巖體---------"<<setw(15)<<"距離函數"<<endl;	
	for(i=0;i<m-1;i++)//做m-1次聚類
	{
		cout<<setw(8)<<i+1<<setw(20)<<c[int(result[i][0])]<<setw(15)<<c[int(result[i][1])]<<setw(13)<<setprecision(3)<<result[i][2]<<endl;
		for(j=0;j<50;j++)
		{
			if(c[int(result[i][0])][j]!='\0')
			{
				temp++;
			}
		}
		c[int(result[i][0])][temp]=',';
		temp++;
		for(j=temp;j<50;j++)
		{
			c[int(result[i][0])][j]=c[int(result[i][1])][j-temp];
		}
		for(j=int(result[i][1]);j<m-1;j++)
		{
			for(k=0;k<50;k++)
			{
				c[j][k]=c[j+1][k];
			}
		}
		temp=0;
	}*/
}