bd.c

在定时器中断中做LED的PWM输出 AT89C2051实现A/D转换的C51程序 单片机开发系统 指令系统 程序设计 定时与中断 系统扩展 接口技术 串行口

//仪器系数自动标定算法。
//标定的结果:
//系数: COEFF[3][3]={{4.287,-2.185,-0.1945},
//                {-0.3181,10.56,-0.2134},
//                {-0.3802,-0.1864,0.3893}} 
//
//本底: cps0[3]={0.2782,0.1101,0.6737}

float former[5][3]={{137.4,9.61,0.59},
         {5.05,284.3,0.37},
         {2.43,8.42,5.71},
         {118.9,270.9,3.05},
         {2.89,7.3,0.9}}; //五个模型的含量数据（已知）

float cps[5][3]={{35.396,2.828,37.787},
	{17.154,28.164,31.537},
	{2.236,1.313,17.828},
	{45.424,28.446,66.166},
	{1.543,0.919,4.608}}; //五个模型的测量数据（已知）

float coef[3][3];  //9个小写字母系数（中间数据）
float COEFF[3][3];  //9个大写字母系数（待标定）
float cps0[3];    //3个本底计数率（待标定）

float A[3][3] ;	//存放三元一次方程组的系数矩阵
float B[3];	//存放三元一次方程组的常数项
float X[3];	//存放三元一次方程组的解:X1=1.537067,X2=-0.896159,X3=3.592657

float DETA ( ) //行列式计算
{
	return (A[0][0]*A[1][1]*A[2][2]+A[1][0]*A[2][1]*A[0][2]+A[2][0]*A[0][1]*A[1][2]
	 -A[2][0]*A[1][1]*A[0][2]-A[0][0]*A[2][1]*A[1][2]-A[1][0]*A[0][1]*A[2][2]);	
}

void swap (int k)//常数项与k列系数交换
{
	int i;
	float t;
	for (i=0;i<3;i++) {
		t=B[i];
		B[i]=A[i][k];
		A[i][k]=t;
		}
}	

void XYZ ( )	//求解三元一次方程组
{
	float d;
	d=DETA( );	//求系数行列式的值
	swap (0) ;	//将x1的系数与常数项交换
	X[0]=DETA()/d;	//解出x1
	swap (0) ;	//恢复原来的系数与常数项
	swap (1) ;	//将x2的系数与常数项交换
	X[1]=DETA()/d;	//解出x2
	swap (1) ;	//恢复原来的系数与常数项
	swap (2) ;	//将x3的系数与常数项交换
	X[2]=DETA()/d;	//解出x3
	swap (2) ;	//恢复原来的系数与常数项
}

void RTK (int k)	//求解第k行的3个小写字母系数
{
	int i,j;
	for (i=0;i<3;i++)	//准备方程组的系数
		for (j=0;j<3;j++) A[i][j]=former[i][j];
	for (i=0;i<3;i++) B[i]=cps[i][k]-cps0[k];//准备方程组的常数项
	XYZ ();	//解出3个小写字母系数
	for (i=0;i<3;i++) coef[k][i]=X[i];//保存3个小写字母系数
}

void conv ( )	//将9个小写字母系数转换成9个大写字母系数
{
	float d;
	int i,j;
	for (i=0;i<3;i++) //取9个小写字母系数	
		for (j=0;j<3;j++) A[i][j]=coef[i][j];
	d=DETA( );	//求系数行列式的值
	COEFF[0][0]=(A[1][1]*A[2][2]-A[2][1]*A[1][2])/d;//计算9个大写字母系数
	COEFF[0][1]=(A[2][1]*A[0][2]-A[0][1]*A[2][2])/d;
	COEFF[0][2]=(A[0][1]*A[1][2]-A[0][2]*A[1][1])/d;
	COEFF[1][0]=(A[1][2]*A[2][0]-A[1][0]*A[2][2])/d;
	COEFF[1][1]=(A[0][0]*A[2][2]-A[2][0]*A[0][2])/d;
	COEFF[1][2]=(A[1][0]*A[0][2]-A[0][0]*A[1][2])/d;
	COEFF[2][0]=(A[1][0]*A[2][1]-A[1][1]*A[2][0])/d;
	COEFF[2][1]=(A[0][1]*A[2][0]-A[0][0]*A[2][1])/d;
	COEFF[2][2]=(A[0][0]*A[1][1]-A[0][1]*A[1][0])/d;
}
	
 main ( )
{
	int i,j,k=1;
	float c;
	for (i=0;i<3;i++) cps0[i]=0; //3个本底计数率初始化为0。
	for (i=0;i<5;i++) {  //进行5次迭代运算。
		for (j=0;j<3;j++) RTK (j);//计算9个小写字母系数
		for (j=0;j<3;j++) 	//计算零值模型各道的净计数率
			B[j]=coef[j][0]*former[4][0]+
			coef[j][1]*former[4][1]+coef[j][2]*former[4][2];
		for (j=0;j<3;j++) cps0[j]=cps[4][j]-B[j];//计算各道的本底计数率
		}
	conv ( );	//将9个小写字母系数转换成9个大写字母系数。
	for (i=0;i<3;i++) 	//计算混合模型各元素的含量
		B[i]=COEFF[i][0]*(cps[3][0]-cps0[0])+
		COEFF[i][1]*(cps[3][1]-cps0[1])+COEFF[i][2]*(cps[3][2]-cps0[2]);
	c=former[3][0]/B[0];
	if ( c<0.96 || c>1.04 ) k=0;//铀含量误差超过百分之四,标定失败
	c=former[3][1]/B[1];
	if ( c<0.94 || c>1.06 ) k=0;//钍含量误差超过百分之六,标定失败
	c=former[3][2]/B[2];
	if ( c<0.88 || c>1.12 ) k=0;//钾含量误差超过百分之十二,标定失败
	while (1) ; //在这一行设置断点,中止程序运行,以便观察程序运行的结果 
}