📄 hh.cpp
字号:
#include "math.h"
#include "stdio.h"
#include "string.h"
#include"stdlib.h"
#include"conio.h"
#include"iostream.h"
#include "fstream.h"
#define G 3 /*类别*/
#define P 2 /*指标*/
#define ql 1/G /*先验概率*/
#define N 11 /*总数*/
#define M 100
struct layer
{ int num;
char depth[20];
int category;
float curves[20];
} layer[N];
struct C
{
float aa;
float bb;
}C[M];
double s2[P][P]={0.0,0.0};
double s1[P][P]={0.0,0.0};
double s3[P][P]={0.0,0.0};
double s4[P][P]={0.0,0.0};
double s[P][P]={0.0,0.0};
double u1[P]={0.0,0.0};
double u2[P]={0.0,0.0};
double u3[P]={0.0,0.0};
double u4[P]={0.0,0.0};
double max=0.0;
void read()
{
int i,j,flag;
char ch1[80],ch2[80];
FILE *fp,*fp1,*fp2,*fp3,*fp4;
if((fp=fopen("井段数据1.txt","r+"))==NULL)
{
printf("cannot open this file\n");
return;
}
fgets(ch1,80,fp);
printf("%s\n\n",ch1);
fgets(ch2,80,fp);
printf("%s",ch2);
if((fp1=fopen("井段数据2.txt","w+"))==NULL)
{
printf("cannot open this file\n");
return;
}
fprintf(fp1,"%s",ch1);
fprintf(fp1,"%s",ch2);
if ((fp2=fopen("井段数据3.txt","w+"))==NULL)
{
printf("cannot open this file\n");
return;
}
fprintf(fp2,"%s",ch1);
fprintf(fp2,"%s",ch2);
if ((fp3=fopen("井段数据4.txt","w+"))==NULL)
{
printf("cannot open this file\n");
return;
}
fprintf(fp3,"%s",ch1);
fprintf(fp3,"%s",ch2);
if ((fp4=fopen("井段数据5.txt","w+"))==NULL)
{
printf("cannot open this file\n");
return;
}
fprintf(fp4,"%s",ch1);
fprintf(fp4,"%s",ch2);
for(i=0;i<N; i++)
{
fscanf(fp,"%d%s%d",&layer[i].num,&layer[i].depth,&layer[i].category);
printf("%d %s %d ",layer[i].num,layer[i].depth,layer[i].category);
for(j=0;j<P;j++)
{
fscanf(fp,"%f",&layer[i].curves[j]);
printf(" %4.1f ",layer[i].curves[j]);
}
printf("\n");
}
for(i=0;i<N; i++)
{
flag=layer[i].category;
if (flag==1)
{
fprintf(fp1,"%d %s %d ",layer[i].num,layer[i].depth,layer[i].category);
for(j=0;j<P;j++)fprintf(fp1," %4.1f",layer[i].curves[j]);
fprintf(fp1,"\n");
}
else if (flag==2)
{
fprintf(fp2,"%d %s %d ",layer[i].num,layer[i].depth,layer[i].category);
for(j=0;j<P;j++)fprintf(fp2," %4.1f",layer[i].curves[j]);
fprintf(fp2,"\n");
}
else if (flag==3)
{
fprintf(fp3,"%d %s %d ",layer[i].num,layer[i].depth,layer[i].category);
for(j=0;j<P;j++)fprintf(fp3," %4.1f",layer[i].curves[j]);
fprintf(fp3,"\n");
}
else if (flag==4)
{
fprintf(fp4,"%d %s %d ",layer[i].num,layer[i].depth,layer[i].category);
for(j=0;j<P;j++)fprintf(fp4," %4.1f",layer[i].curves[j]);
fprintf(fp4,"\n");
}
}
fclose(fp);
fclose(fp1);
fclose(fp2);
fclose(fp3);
fclose(fp4);
}
void average()
{
int i,j,flag;
int g=0,h=0,k=0,f=0;
for(i=0;i<N;i++)
{
flag=layer[i].category;
if(flag==1)
{
for(j=0;j<P;j++)
u1[j]=u1[j]+layer[i].curves[j];
g++;
}
if (flag==2)
{
for(j=0;j<P;j++)
u2[j]=u2[j]+layer[i].curves[j];
h++;
}
if (flag==3)
{
for(j=0;j<P;j++)
u3[j]=u3[j]+layer[i].curves[j];
k++;
}
if (flag==4)
{
for(j=0;j<P;j++)
u4[j]=u4[j]+layer[i].curves[j];
f++;
}
}
printf("--------------------------第一类均值向量-------------------------------------\n");
for(i=0;i<P;i++)
{
u1[i]=u1[i]/g;
printf("%10.2f",u1[i]);
}
printf("\n");
printf("--------------------------第二类均值向量----------------------------------------\n");
for(i=0;i<P;i++)
{
u2[i]=u2[i]/h;
printf("%10.2f",u2[i]);
}
printf("\n");
printf("--------------------------第三类均值向量---------------------------------------\n");
for(i=0;i<P;i++)
{
u3[i]=u3[i]/k;
printf("%10.2f",u3[i]);
}
printf("\n");
printf("--------------------------第四类均值向量---------------------------------------\n");
for(i=0;i<P;i++)
{
u4[i]=u4[i]/f;
printf("%10.2f",u4[i]);
}
printf("\n");
}
void array()
{
// read();
//average();
int i,j,flag=0,t;
for(i=0;i<N;i++)
{
flag=layer[i].category;
if(flag==1)
{
for(j=0;j<P;j++)
{
for(t=0;t<P;t++)
s1[j][t]=s1[j][t]+(layer[i].curves[j]-u1[j])*(layer[i].curves[t]-u1[t]);
}
}
if(flag==2)
{
for(j=0;j<P;j++)
{
for(t=0;t<P;t++)
s2[j][t]=s2[j][t]+(layer[i].curves[j]-u2[j])*(layer[i].curves[t]-u2[t]);
}
}
if(flag==3)
{
for(j=0;j<P;j++)
{
for(t=0;t<P;t++)
s3[j][t]=s3[j][t]+(layer[i].curves[j]-u3[j])*(layer[i].curves[t]-u3[t]);
}
}
if(flag==4)
{
for(j=0;j<P;j++)
{
for(t=0;t<P;t++)
s4[j][t]=s4[j][t]+(layer[i].curves[j]-u4[j])*(layer[i].curves[t]-u4[t]);
}
}
}
for(i=0;i<P;i++)
{
for(j=0;j<P;j++)
{
s1[i][j]= s1[i][j]/(N-G);
s2[i][j]=s2[i][j]/(N-G);
s3[i][j]=s3[i][j]/(N-G);
s4[i][j]=s4[i][j]/(N-G);
s[i][j]=s[i][j]+s1[i][j]+s2[i][j]+s3[i][j]+ s4[i][j];
}
}
printf("-------------------------协方差矩阵*********-------------------------------------\n");
for(i=0;i<P;i++)
{
for(j=0;j<P;j++)
{
printf("%10.4f",s[i][j]);
}
printf("\n");
}
}
double * MatrixInver(double A[],int m,int n) /*矩阵转置*/
{
int i,j;
double *B=NULL;
B=(double *)malloc(m*n*sizeof(double));
for(i=0;i<n;i++)
for(j=0;j<m;j++)
B[i*m+j]=A[j*n+i];
return B;
}
double Surplus(double A[],int m,int n) /*求矩阵行列式*/
{
int i,j,k,p,r;
double X,temp=1,temp1=1,s=0,s1=0;
if(n==2)
{for(i=0;i<m;i++)
for(j=0;j<n;j++)
if((i+j)%2) temp1*=A[i*n+j];
else temp*=A[i*n+j];
X=temp-temp1;}
else{
for(k=0;k<n;k++)
{for(i=0,j=k;i<m,j<n;i++,j++)
temp*=A[i*n+j];
if(m-i)
{for(p=m-i,r=m-1;p>0;p--,r--)
temp*=A[r*n+p-1];}
s+=temp;
temp=1;
}
for(k=n-1;k>=0;k--)
{for(i=0,j=k;i<m,j>=0;i++,j--)
temp1*=A[i*n+j];
if(m-i)
{for(p=m-1,r=i;r<m;p--,r++)
temp1*=A[r*n+p];}
s1+=temp1;
temp1=1;
}
X=s-s1;}
return X;
}
double * MatrixOpp(double A[],int m,int n) /*矩阵求逆*/
{
int i,j,x,y,k;
double *SP=NULL,*AB=NULL,*B=NULL,X,*C;
SP=(double *)malloc(m*n*sizeof(double));
AB=(double *)malloc(m*n*sizeof(double));
B=(double *)malloc(m*n*sizeof(double));
X=Surplus(A,m,n);
X=1/X;
for(i=0;i<m;i++)
for(j=0;j<n;j++)
{for(k=0;k<m*n;k++)
B[k]=A[k];
{for(x=0;x<n;x++)
B[i*n+x]=0;
for(y=0;y<m;y++)
B[m*y+j]=0;
B[i*n+j]=1;
SP[i*n+j]=Surplus(B,m,n);
AB[i*n+j]=X*SP[i*n+j];
}
}
C=MatrixInver(AB,m,n);
return C;
}
void matrix(double X[] ,double *p ,int m,int n ,int k,double c[P]) /*向量与矩阵相乘*/
{
double *w,*b;
int i,j,l,t;
w=X;
b=p;
for(i=0;i<m;i++)
{
for(j=0;j<k;j++)
{
t=i*k+j;
c[t]=0.0;
for(l=0;l<n;l++)
c[t]+=w[i*n+l]*b[l*k+j];
}
}
}
void maxvalue(double *array,int n)
{
double *p,*array_end;
array_end=array+n;
max=*array;
for(p=array+1;p<array_end;p++)
if(*p>max)max=*p;
return;
}
void main()
{
read();
average();
array();
double * MatrixOpp(double A[],int m,int n);
double F1[P]={0.0},F2[P]={0.0},F3[P]={0.0},F4[P]={0.0};
double c1[1]={0.0},c2[1]={0.0},c3[1]={0.0},c4[1]={0.0};
double w1[1]={0.0},w2[1]={0.0},w3[1]={0.0},w4[1]={0.0};
double sum[20]={0.0,0.0,0.0,0.0};
double *p1,*p2,*p3,*p4,*p5,*p6,*p7,*p8,*p9;
FILE *fp1,*fp2;
int i,jj=0;
char dd[1][20];
char ff[1][20];
double X[20]={0.0,0.0};/*输入待判的值*/
printf("--------------------------s协方差矩阵的可逆矩阵-----------------------------\n");
p1=MatrixOpp(s[0],P,P);
for(i=0;i<P*P;i++)
{ printf("%10.4f ",p1[i]);
if(i==1)
printf("\n");
}
printf("\n");
matrix(u1,p1,1,P,P,F1);
printf("%10.4f ",F1[1]);
/******************************************/
p1=MatrixOpp(s[0],P,P);
matrix(u2,p1,1,P,P,F2);
printf("%10.4f",F2[1]);
/******************************************/
p1=MatrixOpp(s[0],P,P);
matrix(u3,p1,1,P,P,F3);
printf("%10.4f ",F3[1]);
/*******************************************/
p1=MatrixOpp(s[0],P,P);
matrix(u4,p1,1,P,P,F4);
printf("%10.4f",F4[1]);
printf("\n");
p2=MatrixInver(u1,P,1);
p3=MatrixInver(u2,P,1);
p4=MatrixInver(u3,P,1);
p5=MatrixInver(u4,P,1);
double h1=0.0,h2=0.0;
matrix(F1,u1,P,P,1,c1);
printf("%10.4f",c1[0]);
h1=c1[0];
matrix(F2,u2,P,P,1,c2);
printf("%10.4f",c2[0]);
printf("\n");
h2=c2[0];
matrix(F3,u3,P,P,1,c3);
printf("%10.4f",c3[0]);
double h3=0.0,h4=0.0;
h3=c3[0];
matrix(F4,u4,P,P,1,c4);
printf("%10.4f",c4[0]);
h4=c4[0];
printf("\n");
sum[0]=sum[0]+log(0.5)-0.5*h1;
sum[1]=sum[1]+log(0.5)-0.5*h2;
sum[2]=sum[2]+log(0.5)-0.5*h3;
sum[3]=sum[3]+log(0.5)-0.5*h4;
printf("\n");
printf("-------------------------**********-------------------------------\n");
printf(" %10.4f ",sum[0]);
printf(" %10.4f ",sum[1]);
printf(" %10.4f",sum[2]);
printf(" %10.4f ",sum[3]);
/**************************/
p6=MatrixInver(F1,P,1);
p7=MatrixInver(F2,P,1);
p8=MatrixInver(F3,P,1);
p9=MatrixInver(F4,P,1);
/********************************************************************************************************************/
printf("\n");
printf("请输入你要处理的文件名：");
scanf("%s",dd[1]);
if((fp1=fopen(dd[1],"r+"))==NULL)
{
printf("Can't find file\n");
exit(0);
}
printf("请输入要写入的文件名：");
scanf("%s",ff[1]);
if((fp2=fopen(ff[1],"w+"))==NULL)
{
perror("Cannot open data file");
exit(0);
}
for(jj=0;feof(fp1)==0;jj++)
{
double kk=0.0,v=0.0, tt=0.0;
double ss1=0.0, ss2=0.0, ss3=0.0, ss4=0.0;
fscanf(fp1,"%10f %10f",&C[jj].aa,&C[jj].bb);
X[0]=C[jj].aa;
X[1]=C[jj].bb;
double sum[9]={log(0.5)-0.5*h1,log(0.5)-0.5*h2,log(0.5)-0.5*h3};
matrix(X,p6,1,P,P,w1);
sum[0]=sum[0]+w1[0];
ss1=sum[0];
printf(" %10.4f ",ss1);
matrix(X,p7,1,P,P,w2);
sum[1]=sum[1]+w2[0];
ss2=sum[1];
printf(" %10.4f ",ss2);
matrix(X,p8,1,P,P,w3);
sum[2]=sum[2]+w3[0];
ss3=sum[2];
printf(" %10.4f ",ss3);
matrix(X,p9,1,P,P,w4);
sum[3]=sum[3]+w4[0];
ss4=sum[3];
printf(" %10.4f ",ss4);
double *w ;
w=sum;
maxvalue(w,4);
if((ss1/max)==1)
{ printf("第%d数据属于第1大类\n",jj+1);
if(feof(fp1)==0)
fprintf(fp2,"第%d个数据属于第1大类",jj+1);
}
else if((ss2/max)==1)
{ printf("第%d数据属于第2大类\n",jj+1);
if(feof(fp1)==0)
fprintf(fp2,"第%d个数据属于第2大类",jj+1);
}
else if((ss3/max)==1)
{
printf("第%d数据属于第3大类\n",jj+1);
if(feof(fp1)==0)
fprintf(fp2,"第%d个数据属于第3大类",jj+1);
}
else if((ss4/max)==1)
{
printf("第%d数据属于第四大类\n",jj+1);
if(feof(fp1)==0)
fprintf(fp2,"第%d个数据属于第四大类\n",jj+1);
}
printf("最大值=%10.4f\n",max);
kk=max*max;
tt=max*sum[0]+max*sum[1]+max*sum[2];
if(feof(fp1)==0) fprintf(fp2," 后验概率为:");
v=(kk/tt)*100;
if(v>100 ) v=100;
if(v<0) v=0;
printf("%10.4f\n",v);
if(feof(fp1)==0) fprintf(fp2,"%10.4f\n",v);
}
}
⌨️ 快捷键说明
复制代码
Ctrl + C
搜索代码
Ctrl + F
全屏模式
F11
切换主题
Ctrl + Shift + D
显示快捷键
?
增大字号
Ctrl + =
减小字号
Ctrl + -