粗大-系统-随机误差处理
- 1、下载文档前请自行甄别文档内容的完整性,平台不提供额外的编辑、内容补充、找答案等附加服务。
- 2、"仅部分预览"的文档,不可在线预览部分如存在完整性等问题,可反馈申请退款(可完整预览的文档不适用该条件!)。
- 3、如文档侵犯您的权益,请联系客服反馈,我们会尽快为您处理(人工客服工作时间:9:00-18:30)。
课程设计用仪器设备名称
此次课程设计用到的仪器设备和软件包括: (1) 个人计算机; (2) Matlab 软件。
课程设计过程
1、课程设计处理原理:
此次课程开展的数据处理包:(1)粗大误差处理;(2)系统误差处理;(3)随机误差处理。
他们的原理分别分析如下:
(1)粗大误差处理
对于粗大误差,采用莱以特准则和罗曼诺夫斯基准则。
莱以特准则:求出数据的算数平均值x 和标准差σ,将残差的绝对值i x v 和
3σ进行比较,大于3σ的值都认为是粗大误差。
罗曼诺夫斯基准则:首先剔除该数据中的最大值,然后再按照t 分布检验,
求出该项与剔除后平均值的差,即d x x −,再与()2,K n a σ−进行比较,如果前者大于等于后者,那么该数据有系统误差。
(2)系统误差处理
对于系统误差,我们采用了残差总和判断法,阿贝-赫梅判别法,标准差比较法,他们的原理如下:
残差总和判断法:对于等精度的系统测量数据12,,...n x x x ,设相对的残差
分别是12,,...n v v v ,若有12n
i i v =>∑,则怀疑测量数据有系统误差
阿贝-赫梅判别法:对于等精度的系统测量数据12,,...n x x x ,设相对的残,
分别是12,,...n v v v ,122311
1
...n
n n i i i u v v v v v v v v
−+==+++=
∑,如果2u >,
则判定该组数据含有系统误差。
标准差比较法:对于等精度的系统测量数据12,,...n x x x ,设相对的残差分
别是12,
,...n v v v ,用不同的公式计算标准差,通过比较可以发现存在的系
统误差。
用贝塞尔公式计算,1s
=
,用别捷尔斯公式计算,
1s
=211s s ≥,则怀疑测量中存在系统误差。
(3)随机误差处理
我们考虑了正态分布和t 分布两种情况,通过置信概率和自由度分别在正态分布积分表和t 分布表中找到对应的t 值,再求出极限误差lim x t ςσ=+。
2
课程设计的整体流程图如图(一)所示。
在图(一)中,粗大误差分析,系统误差分析,随机误差分析都作为子程序存在。
首先我们是将存储在txt 文件中的测量数据导入到matlab 中,然后进行在子程序中用两种方式进行粗大误差分析,并返回剔除异常值以及剔除异常值后的测量数据。
接着进行系统误差分析,用了三种方法检测是否具有系统误差,并返回测量结果。
之后进行随机误差分析,返回两种分布的极限误差。
最后将本次测量结果都写入到txt 文件中。
斯基准则进行分析。
首先看莱以特准则,计算测量数据第i项的残差,之后比较残差与三倍标准差的的值大小,如果残差大于三倍标准差就剔除该项,否则进行下一项的比较,比较完成后返回剔除异常值后剩余的数据。
另外是罗曼诺夫斯基准则,首先将最大值作为怀疑对象,从测量数据中剔除出来。
然后求该值与剔除后数据的平均值之差,与t分布下的标准差进行比较,如果前者大,则进行下一次循环,否则结束循环,并把该次循环中的最大值重新放入测量数据中,最后返回剔除异常值后剩余的数据。
(3)系统误差处理
计算判断方法所需要的数值,然后对各个方法的判断条件进行比较,如果满足就返回有系统误差,否则返回没有系统误差。
如(图四)所示,首先接收已经剔除异常值的数据,然后计算平均值和标准差,进行正态分布和t 分
布的查表,我将正态分布和t 分布的表格数据存储到matlab 的数组中,并将它们封装成独立的
子m 程序,只要输入对应的自由度和置信概率就
可以得到相应的t 值。
然后将两个t 值乘上标准
差得到对应的极限误差,并返回两个极限误差。
3、实验规划
如何输入数据:将测量结果以行存储在txt文件中,,在matlab中利用textread函数读取txt数据,并将其存储在一维数组中。
如何使用matlab:首先新建一个主文件并保存,接着建立了5个函数,如下所示: Thick_Error.m 粗大误差处理
Systematic_Error.m 系统误差处理
Random_Error.m 随机误差处理
Normal_distribution.m正态分布t值查找
T_FenBu.m t分布t值查找
如何编写程序:依照误差处理的原理一步一步进行数据分析,将所得结果进行打印在命令行窗口,并将所有数据打印在txt文件中
如何程序运行、调试程序:运行上边三角运行按钮,还可以选中某些行代码右键运行。
运行错误时会设置断点,一步一步运行,并看着工作区的显示的数值是否是我想得到的那个值。
通常我还在程序中加一些打印字符
串,用于发现错误之处。
4、代码分析
(1)数据输入:
path="C:\Users\dingshuai\Desktop\误差\课程设计\误差处理.txt";
date=textread(path,'%f',15)';
首先将txt文件路径存储在path中,利用textread函数读取存储在txt文件中的15个测量数据,并存储在数组date中
(2)主函数
粗大误差处理
[result_1,rid_value_1,result_2,rid_value_2]=Thick_Error(date);
%莱以特准则:rid_value_1表示剔除值 result_1表示剔除后剩余数据
%罗曼诺夫斯基准则:rid_value_2表示剔除值 result_2表示剔除后剩余数据
disp("莱以特准则剔除数据");disp(rid_value_1);
disp("罗曼诺夫斯基准则");disp(rid_value_2);
判断系统误差
[answer_1,answer_2,answer_3]=Systematic_Error(result_2); %answer 是否具
有系统误差
% answer_1,answer_2,answer_3分别对应残差总和判断法,阿贝-赫梅判别法,标准差
比较法的判断
disp("残差总和判断法是否有系统误差");disp(answer_1);
disp("阿贝-赫梅判别法是否有系统误差");disp(answer_2);
disp("标准差比较法是否有系统误差");disp(answer_3);
随机误差
[limit_error_zt,limit_error_t]=Random_Error(result_2);
%分别表示正态分布和t分布的极限误差
disp("正态分布极限误差")
disp(limit_error_zt);
disp("t分布的极限误差");
disp(limit_error_t);
写入文件
fileID=fopen(path,'w');
fprintf(fileID,'%.2f ',date); fprintf(fileID,'\n');
fprintf(fileID,'算数平均值: %.2f\n',mean(result_2));
fprintf(fileID,'莱以特准则粗大误差剔除数据: %.2f\n',rid_value_1);
fprintf(fileID,'%s',"罗曼诺夫斯基准则粗大误差剔除数据: ");
fprintf(fileID,'%.2f ',rid_value_2); fprintf(fileID,'\n');
fprintf(fileID,'残差总和判断法是否有系统误差是否有系统误差: %s\n',answer_1); fprintf(fileID,'阿贝-赫梅判别法是否有系统误差是否有系统误
差: %s\n',answer_1);
fprintf(fileID,'标准差比较法是否有系统误差是否有系统误差: %s\n',answer_1); fprintf(fileID,'按正态分布结
果: %.2f+-%.2f\n',mean(result_2),limit_error_zt);
fprintf(fileID,'按t分布结果: %.2f+-%.2f',mean(result_2),limit_error_t); fclose(fileID);
利用子函数求出各种数值结果,将结果打印在命令行窗口,并将结果写入txt
文件
(3)莱以特准则
a=1; %用于剔除数据的索引
tem_date=date;
date_avg=mean(tem_date); %求平均值
date_std=std(tem_date); %求标准差
residual=tem_date-date_avg; %求残差
for i=1:length(tem_date) %大于三倍标准差元素进行剔除
if(residual(i)>(3*date_std))
rid_value_1(a)=tem_date(i); %求剔除值
a=a+1; %索引加一,用于存储下次异常值
tem_date(i)=[]; %剔除粗大误差
end
end
result_1=tem_date;
首先定义a值,作为异常值数组的索引,将测量数据存储到tem_date中,然后求其平均值date_avg,标准差date_std,再求其残差residual,然后进行for循环,如果满足莱以特准则的条件,就让该索引对应的值为空,即剔除该数据,最后将剔除异常值的数据存储到result_1中,并返回到主函数。
(4)罗曼诺夫斯基准则
b=1; %用于剔除数据的索引
while(1)
tem_date=date; %临时保存date原本数据
[max_value,i]=max(tem_date); %求其最大值及其索引,作为怀疑对象
rid_value_2(b)=max_value; %存储怀疑对象
tem_date(i)=[]; %剔除怀疑对象
avg=mean(date); %求平均值
std_=std(date); %求标准差
if(abs(max_value-avg)<(std_*T_FenBu((length(tem_date)-1),0.05)))
rid_value_2(b)=[]; %如果不满足粗大误差条件就将怀疑对象清除break;
else
b=b+1; %满足粗大误差条件就将怀疑对象增一
date=tem_date;
end
end
result_2=date;
首先定义b值,作为异常值数组的索引,接着进入while循环,其中有break 语句,所以不是死循环。
提取其最大值及其索引,将怀疑对象存储在
rid_value_2中,计算平均值和标准差,再进行罗曼诺夫斯基准则的判断条件,如果该值满足,则将异常值索引加1,方便下次异常值存储,如果不满足,将该次怀疑的异常值对象去掉,并退出循环,最后再将测量数据返回主函数
(5)判断系统误差
输入: date 表示要进行系统误差检验的数据
返回值: answer_1,answer_2,answer_3分别对应残差总和判断法,阿贝-赫梅判别法,标准差比较法的判断
残差总和判断法
function [answer_1,answer_2,answer_3]=Systematic_Error(date)
date_avg=mean(date);%求平均值
residual=abs(date-date_avg);%求残差绝对值
sum_residual=sum(residual);%求残差绝对值和
s=sqrt((sum(residual.*residual))/(length(date)));%贝塞尔公式求标准差
if(sum_residual>(2*s*sqrt(length(date))))
answer_1="数据有系统误差";
else
answer_1="数据没有系统误差";
end
阿贝-赫梅判别法
date_1=date(1:length(date)-1);
date_2=date(2:length(date));
u=abs(sum(date_1.*date_2));
if(u>(sqrt(length(date)-1)*std(date)*std(date)))
answer_2="数据有系统误差";
else
answer_2="数据没有系统误差";
end
标准差比较法
s1=std(date);
s2=1.253*sum_residual/sqrt(length(date)*(length(date)-1));
if((s2/s1)<(1+2/sqrt(length(date)-1)))
answer_3="数据没有系统误差";
else
answer_3="数据有系统误差";
end
从主函数中获取已经剔除异常值的数据,利用三种方法进行判断,并且返回三种方法各自分析的结果。
(6)随机误差处理
输入: date 表示要进行系统误差检验的数据
返回值: 正态分布和t分布的极限误差
function [limit_error_zt,limit_error_t]=Random_Error(date)
date_avg=mean(date);%求平均值
residual=abs(date-date_avg);%求残差绝对值
s=sqrt((sum(residual.*residual))/(length(date)));
s_x=s/sqrt(length(date));%求单次重复标准差
正态分布
t_zt=Normal_distribution(0.95/2); % 置信概率为%95,自由度为length(date)-1 limit_error_zt=t_zt*s_x;
t分布
t_t=T_FenBu(length(date)-1,0.05);
limit_error_t=t_t*s_x;
接收剔除异常值后的数据,计算平均值和标准差,再通过
Normal_distribution和T_FenBu函数获取对应的t值,计算两种情况的极限误差,返回主函数中。
两种分布获取t值代码如下:
输入: Fiducial_Probability 测量数据的置信概率
返回值: t 正态分布的t值
正态分布
function t=Normal_distribution(Fiducial_Probability)
a=[0.00 0.05 0.10 0.15 0.20 0.25 0.30 0.35 1.60 1.65 1.70 1.75 1.80 1.85
1.90 0.40 0.45 0.50 0.55 0.60 0.65 0.70 0.75 1.95
2.00 2.10 2.20 2.30 2.40
2.50 0.80 0.85 0.90 0.95 1.00 1.05 1.10 1.15 2.60 2.70 2.80 2.90
3.5 3.20
3.40 1.20 1.25 1.30 1.35 1.40 1.45 1.50 1.55 3.60 3.80
4.00 4.50
5.00 ]; b=[0.0000 0.0199 0.0398 0.0596 0.0793 0.0987 0.1179 0.1368 0.4452 0.4505 0.4554 0.4599 0.4641 0.4678 0.4713 0.1554 0.1736 0.1915 0.2088 0.2257 0.2422 0.2580 0.2734 0.4744 0.4772 0.4821 0.4861 0.4893 0.4918 0.4938 0.2881 0.3023 0.3159 0.3289 0.3413 0.3531 0.3643 0.3740 0.4953 0.4965 0.4974 0.4981 0.49865 0.49931 0.49966 0.3949 0.3944 0.4032 0.4115 0.4192 0.4265 0.4332 0.4394 0.499841 0.49928 0.49968 0.499997 0.49999997];
[~,i]=min(abs(b-Fiducial_Probability));%i表示最小值在向量中的索引
t=a(i);
输入: dof 代表自由度a 代表1-置信概率
返回值: t t分布的t值
t分布
function t=T_FenBu(dof,a)
v=[1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 40 50 60 70 80 90 100];
a1=[12.71 4.30 3.18 2.78 2.57 2.45 2.36 2.31 3.26 .23 2.20 2.18 2.16 2.14 2.13 2.12 2.11 2.10 2.09 2.09 2.08 2.07 2.07 2.06 2.06 2.06 2.05 2.05 2.05 2.04 2.02 2.01 2.00 1.99 1.99 1.99 1.98];
a2=[63.66 9.92 5.84 4.60 4.03 3.71 3.50 3.36 3.25 3.17 3.11 3.05 3.01 2.98 2.95 2.92 2.90 2.88 2.86 2.85 2.83 2.82 2.81 2.80 2.79 2.78 2.77 2.76 2.76 2.75 2.70 2.68 2.66 2.65 2.64 2.63 2.63];
a3=[235.0 19.21 9.21 5.51 4.90 4.53 4.28 4.09 3.96 3.85 3.76 3.69 3.64 3.59 3.54 3.51 3.48 3.45 3.42 3.40 3.38 3.36 3.34 3.33 3.32 3.30 3.29 3.28 3.27 3.20 3.18 3.13 3.11 3.10 3.09 3.08 3.00];
for i=1:length(v)
if(v(i)==dof) %首先判断自由度一样
if(a==0.05) %再判断a一样
t=a1(i);
elseif(a==0.01)
t=a2(i);
elseif(a==0.0027)
t=a3(i);
else
t=0;
end
end
end
(7)正态分布输出结果
按正态分布结果: 20.40+-0.03
(8)t分布输出结果
按t分布结果: 20.40+-0.03
5、实验结果及分析
实验结果如左所示,将测量数据的平均值、
标准差、粗大误差,系统误差以及极限误差显示
在命令行窗口,txt文件如下图所示,如有其它
数据需要,可对写入代码稍作修改将其它数据写
入。
输出的数据均有较多的位数保留,如不需
要,可适当进行取舍。
另外,程序对相关的查表值利用数组进行写
入,可满足大部分情况下的数据处理,进行数据
修改或者置信概率修改仍可得到正确的处理结
果。
课程设计心得
通过这次学习,对MATLAB的使用更加熟练了,虽然以前数学建模的时候用过这个软件,但是很多操作依旧不熟悉,通过在百度上查找,现在我也掌握的了如何从txt文件中读写数据。
但是也有很多不足的地方,比如对于正态分布和t 分布的表格的存储我是用数组存储的,感觉这个方法效率不高,但是上网查也没查到很有效的方法。
另外,我也用到了visio这个流程图软件,以前从来没有用过,现在对它的使用虽然不熟练,但是也不陌生了。
这次的难点我觉得有两个,一个是在编写程序过程中部分函数不是很了解,再加上对于函数的使用有些在网上的解释的不是很清楚。
第二个就是测量误差的方法记得不是很牢,依然需要去重新翻阅课件。
这次课设学到了不少,也发现自己对于知识的理解不够深刻,像是依旧停留在记公式的状态。
我需要进一步的加深理解,把这些知识争取能够成为自己内化的东西
附录
说明:代码是以实时函数方式打开进行复制
主函数:Error.m
clear;clc;
path="C:\Users\dingshuai\Desktop\误差\课程设计\误差处理.txt";
date=textread(path,'%f',15)';
粗大误差处理
[result_1,rid_value_1,result_2,rid_value_2]=Thick_Error(date);
%莱以特准则:rid_value_1表示剔除值 result_1表示剔除后剩余数据
%罗曼诺夫斯基准则:rid_value_2表示剔除值 result_2表示剔除后剩余数据
disp("莱以特准则剔除数据");disp(rid_value_1);
disp("罗曼诺夫斯基准则");disp(rid_value_2);
判断系统误差
[answer_1,answer_2,answer_3]=Systematic_Error(result_2); %answer 是否具有系统误差
% answer_1,answer_2,answer_3分别对应残差总和判断法,阿贝-赫梅判别法,标准差
比较法的判断
disp("残差总和判断法是否有系统误差");disp(answer_1);
disp("阿贝-赫梅判别法是否有系统误差");disp(answer_2);
disp("标准差比较法是否有系统误差");disp(answer_3);
随机误差
[limit_error_zt,limit_error_t]=Random_Error(result_2);
%分别表示正态分布和t分布的极限误差
disp("剔除异常值后的平均值");
disp(mean(result_2));
disp("正态分布极限误差")
disp(limit_error_zt);
disp("t分布的极限误差");
disp(limit_error_t);
写入文件
fileID=fopen(path,'w');
fprintf(fileID,'%.2f ',date); fprintf(fileID,'\n');
fprintf(fileID,'算数平均值: %.2f\n',mean(result_2));
fprintf(fileID,'莱以特准则粗大误差剔除数据: %.2f\n',rid_value_1);
fprintf(fileID,'%s',"罗曼诺夫斯基准则粗大误差剔除数据: ");
fprintf(fileID,'%.2f ',rid_value_2); fprintf(fileID,'\n');
fprintf(fileID,'残差总和判断法是否有系统误差是否有系统误差: %s\n',answer_1); fprintf(fileID,'阿贝-赫梅判别法是否有系统误差是否有系统误
差: %s\n',answer_1);
fprintf(fileID,'标准差比较法是否有系统误差是否有系统误差: %s\n',answer_1); fprintf(fileID,'按正态分布结
果: %.2f+-%.2f\n',mean(result_2),limit_error_zt);
fprintf(fileID,'按t分布结果: %.2f+-%.2f',mean(result_2),limit_error_t); fclose(fileID);
粗大误差处理:
输入: date 表示要进行粗大误差检验的数据
返回值:莱以特准则:rid_value_1表示剔除值result_1表示剔除后剩余数据返回值:罗曼诺夫斯基准则:rid_value_2表示剔除值result_2表示剔除后剩余数据
莱以特准则
function [result_1,rid_value_1,result_2,rid_value_2]=Thick_Error(date)
a=1; %用于剔除数据的索引
tem_date=date;
date_avg=mean(tem_date); %求平均值
date_std=std(tem_date); %求标准差
residual=tem_date-date_avg; %求残差
for i=1:length(tem_date) %大于三倍标准差元素进行剔除
if(residual(i)>(3*date_std))
rid_value_1(a)=tem_date(i); %求剔除值
a=a+1; %索引加一,用于存储下次异常值
tem_date(i)=[]; %剔除粗大误差
end
end
result_1=tem_date;
罗曼诺夫斯基准则
b=1; %用于剔除数据的索引
while(1)
tem_date=date; %临时保存date原本数据
[max_value,i]=max(tem_date); %求其最大值及其索引,作为怀疑对象
rid_value_2(b)=max_value; %存储怀疑对象
tem_date(i)=[]; %剔除怀疑对象
avg=mean(date); %求平均值
std_=std(date); %求标准差
if(abs(max_value-avg)<(std_*T_FenBu((length(tem_date)-1),0.05)))
rid_value_2(b)=[]; %如果不满足粗大误差条件就将怀疑对象清除break;
else
b=b+1; %满足粗大误差条件就将怀疑对象增一
date=tem_date;
end
end
result_2=date;
系统误差处理:
输入: date 表示要进行系统误差检验的数据
返回值: answer_1,answer_2,answer_3分别对应残差总和判断法,阿贝-赫梅判别法,标准差比较法的判断
残差总和判断法
function [answer_1,answer_2,answer_3]=Systematic_Error(date)
date_avg=mean(date);%求平均值
residual=abs(date-date_avg);%求残差绝对值
sum_residual=sum(residual);%求残差绝对值和
s=sqrt((sum(residual.*residual))/(length(date)));%贝塞尔公式求标准差
if(sum_residual>(2*s*sqrt(length(date))))
answer_1="数据有系统误差";
else
answer_1="数据没有系统误差";
end
阿贝-赫梅判别法
date_1=date(1:length(date)-1);
date_2=date(2:length(date));
u=abs(sum(date_1.*date_2));
if(u>(sqrt(length(date)-1)*std(date)*std(date)))
answer_2="数据有系统误差";
else
answer_2="数据没有系统误差";
end
标准差比较法
s1=std(date);
s2=1.253*sum_residual/sqrt(length(date)*(length(date)-1));
if((s2/s1)<(1+2/sqrt(length(date)-1)))
answer_3="数据没有系统误差";
else
answer_3="数据有系统误差";
end
随机误差处理:
输入: date 表示要进行系统误差检验的数据
返回值: 正态分布和t分布的极限误差
function [limit_error_zt,limit_error_t]=Random_Error(date)
date_avg=mean(date);%求平均值
residual=abs(date-date_avg);%求残差绝对值
s=sqrt((sum(residual.*residual))/(length(date)));
s_x=s/sqrt(length(date));%求单次重复标准差
正态分布
t_zt=Normal_distribution(0.95/2); % 置信概率为%95,自由度为length(date)-1 limit_error_zt=t_zt*s_x;
t分布
t_t=T_FenBu(length(date)-1,0.05);
limit_error_t=t_t*s_x;
正态分布查表:
输入: Fiducial_Probability 测量数据的置信概率
返回值: t 正态分布的t值
正态分布
function t=Normal_distribution(Fiducial_Probability)
a=[0.00 0.05 0.10 0.15 0.20 0.25 0.30 0.35 1.60 1.65 1.70 1.75 1.80 1.85
1.90 0.40 0.45 0.50 0.55 0.60 0.65 0.70 0.75 1.95
2.00 2.10 2.20 2.30 2.40
2.50 0.80 0.85 0.90 0.95 1.00 1.05 1.10 1.15 2.60 2.70 2.80 2.90
3.5 3.20
3.40 1.20 1.25 1.30 1.35 1.40 1.45 1.50 1.55 3.60 3.80
4.00 4.50
5.00 ];
b=[0.0000 0.0199 0.0398 0.0596 0.0793 0.0987 0.1179 0.1368 0.4452 0.4505 0.4554 0.4599 0.4641 0.4678 0.4713 0.1554 0.1736 0.1915 0.2088 0.2257 0.2422 0.2580 0.2734 0.4744 0.4772 0.4821 0.4861 0.4893 0.4918 0.4938 0.2881 0.3023 0.3159 0.3289 0.3413 0.3531 0.3643 0.3740 0.4953 0.4965 0.4974 0.4981 0.49865 0.49931 0.49966 0.3949 0.3944 0.4032 0.4115 0.4192 0.4265 0.4332 0.4394 0.499841 0.49928 0.49968 0.499997 0.49999997];
[~,i]=min(abs(b-Fiducial_Probability));%i表示最小值在向量中的索引
t=a(i);
t分布查表:
输入: dof 代表自由度a 代表1-置信概率
返回值: t t分布的t值
t分布
function t=T_FenBu(dof,a)
v=[1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 40 50 60 70 80 90 100];
a1=[12.71 4.30 3.18 2.78 2.57 2.45 2.36 2.31 3.26 .23 2.20 2.18 2.16 2.14 2.13 2.12 2.11 2.10 2.09 2.09 2.08 2.07 2.07 2.06 2.06 2.06 2.05 2.05 2.05 2.04 2.02 2.01 2.00 1.99 1.99 1.99 1.98];
a2=[63.66 9.92 5.84 4.60 4.03 3.71 3.50 3.36 3.25 3.17 3.11 3.05 3.01 2.98 2.95 2.92 2.90 2.88 2.86 2.85 2.83 2.82 2.81 2.80 2.79 2.78 2.77 2.76 2.76 2.75 2.70 2.68 2.66 2.65 2.64 2.63 2.63];
a3=[235.0 19.21 9.21 5.51 4.90 4.53 4.28 4.09 3.96 3.85 3.76 3.69 3.64 3.59 3.54 3.51 3.48 3.45 3.42 3.40 3.38 3.36 3.34 3.33 3.32 3.30 3.29 3.28 3.27 3.20 3.18 3.13 3.11 3.10 3.09 3.08 3.00];
for i=1:length(v)
if(v(i)==dof) %首先判断自由度一样
if(a==0.05) %再判断a一样
t=a1(i);
elseif(a==0.01)
t=a2(i);
elseif(a==0.0027)
t=a3(i);
else
t=0;
end
end
end。