matlab 常见经典平差 程序

matlab 常见经典平差 程序
matlab 常见经典平差 程序

希望本人编写的经典平差可以对matlab初学者以及测绘专业的学生有用

By cowboy

function void()

sprintf('请选择平差类型')

sprintf('1: 参数平差\n2: 条件平差\n3:具有条件的参数平差\n4: 具有参数的条件平差\n') a=input('请输入对应号码\n');

switch (a)

case 1

% 参数平差V=AX-L

sprintf('1:参数平差V=AX-L')

n=input('请输入n=');

t=input('请输入t=');

A=input('请输入系数矩阵A=');

L=input('请输入常数项矩阵L=');

P=diag(input('请输入权阵P='))

%P=input('请输入权阵P=');

N=A'*P*A;

U=A'*P*L;

sprintf('开始计算')

X=inv(N)*(U);

V=A*X-L;

Qx=inv(N);

u=sqrt((V'*P*V)/(n-t));

fid=fopen('data.txt','wt'); %写入文件路径,文件输出

fprintf(fid,'=====================参数平差: V=AX-L =====================\n');

fprintf(fid,'输出n,t\n') ;

fprintf(fid,'%12.5f %12.5f\n',n ,t) ;%n,t

fprintf(fid,'输出矩阵A\n')

[m,n]=size(A); %A

for i=1:1:m

for j=1:1:n

if j==n

fprintf(fid,'%12.5f\n',A(i,j));

else

fprintf(fid,'%12.5f\t',A(i,j));

end

end

end

fprintf(fid,'\n');

fprintf(fid,'输出矩阵L\n') ;

[m,n]=size(L); %L

for i=1:1:m

for j=1:1:n

if j==n

fprintf(fid,'%12.5f\n',L(i,j));

else

fprintf(fid,'%12.5f\t',L(i,j));

end

end

end

fprintf(fid,'\n');

fprintf(fid,'输出矩阵P\n') ;

[m,n]=size(P); %P for i=1:1:m

for j=1:1:n

if j==n

fprintf(fid,'%12.5f\n',P(i,j));

else

fprintf(fid,'%12.5f\t',P(i,j));

end

end

end

fprintf(fid,'\n');

fprintf(fid,'输出矩阵X\n');

[m,n]=size(X); %X for i=1:1:m

for j=1:1:n

if j==n

fprintf(fid,'%12.5f\n',X(i,j));

else

fprintf(fid,'%12.5f\t',X(i,j));

end

end

end

fprintf(fid,'\n');

fprintf(fid,'输出矩阵V\n');

[m,n]=size(V); %V for i=1:1:m

for j=1:1:n

if j==n

fprintf(fid,'%12.5f\n',V(i,j));

else

fprintf(fid,'%12.5f\t',V(i,j));

end

end

end

fprintf(fid,'\n');

fprintf(fid,'输出矩阵Qx\n');

[m,n]=size(Qx); %Qx for i=1:1:m

for j=1:1:n

if j==n

fprintf(fid,'%12.5f\n',Qx(i,j));

else

fprintf(fid,'%12.5f\t',Qx(i,j));

end

end

end

fprintf(fid,'\n');

fprintf(fid,'输出单位权中误差u\n');

[m,n]=size(u); %u for i=1:1:m

for j=1:1:n

if j==n

fprintf(fid,'%12.5f\n',u(i,j));

else

fprintf(fid,'%12.5f\t',u(i,j));

end

end

end

fprintf(fid,'\n');

fclose(fid);

case 2

% 条件平差BV+W=0

sprintf('2:条件平差BV+W=0')

n=input('请输入n=');

t=input('请输入t=');

B=input('请输入系数矩阵B=');

W=input('请输入自由项向量W=');

P=input('请输入权阵P=');

sprintf('开始计算')

K=-inv((B*inv(P)*B'))*W;

V=inv(P)*B'*K;

u=sqrt((V'*P*V)/(n-t));

Ql=inv(P)-inv(P)*B'*inv(B*inv(P)*B')*B*inv(P);

fid=fopen('data.txt','wt'); %写入文件路径,文件输出

fprintf(fid,'=====================条件平差: BV+W=0 =====================\n'); fprintf(fid,'输出n,t\n') ;

fprintf(fid,'%12.5f %12.5f\n',n ,t) ;%n,t

fprintf(fid,'输出矩阵B\n');

[m,n]=size(B); %B

for i=1:1:m

for j=1:1:n

if j==n

fprintf(fid,'%12.5f\n',B(i,j));

else

fprintf(fid,'%12.5f\t',B(i,j));

end

end

end

fprintf(fid,'\n');

fprintf(fid,'输出矩阵W\n');

[m,n]=size(W); %W

for i=1:1:m

for j=1:1:n

if j==n

fprintf(fid,'%12.5f\n',W(i,j));

else

fprintf(fid,'%12.5f\t',W(i,j));

end

end

end

fprintf(fid,'\n');

fprintf(fid,'输出矩阵P\n');

[m,n]=size(P); %P

for i=1:1:m

for j=1:1:n

if j==n

fprintf(fid,'%12.5f\n',P(i,j));

else

fprintf(fid,'%12.5f\t',P(i,j));

end

end

end

fprintf(fid,'\n');

fprintf(fid,'输出矩阵V\n');

[m,n]=size(V); %V

for i=1:1:m

for j=1:1:n

if j==n

fprintf(fid,'%12.5f\n',V(i,j));

else

fprintf(fid,'%12.5f\t',V(i,j));

end

end

end

fprintf(fid,'\n');

fprintf(fid,'输出矩阵Qx\n');

[m,n]=size(Ql); %Ql for i=1:1:m

for j=1:1:n

if j==n

fprintf(fid,'%12.5f\n',Ql(i,j));

else

fprintf(fid,'%12.5f\t',Ql(i,j));

end

end

end

fprintf(fid,'\n');

fprintf(fid,'输出单位权中误差u\n');

[m,n]=size(u); %u for i=1:1:m

for j=1:1:n

if j==n

fprintf(fid,'%12.5f\n',u(i,j));

else

fprintf(fid,'%12.5f\t',u(i,j));

end

end

end

fprintf(fid,'\n');

fclose(fid);

case 3

% 具有条件的参数平差

% 参数平差V=AX-L

% 条件平差BX+W=0

sprintf('3:具有条件的参数平差V=AX-L;BX+W=0')

n=input('请输入n=');

t=input('请输入t=');

r=input('请输入r=');

A=input('请输入误差方程系数矩阵A=');

L=input('请输入常数项矩阵L=');

B=input('请输入条件方程系数矩阵B=');

W=input('请输入自由项向量W=');

P=input('请输入权阵P=');

sprintf('开始计算')

K=-(inv(B*inv(A'*P*A)*B'))*(W+B*inv(A'*P*A)*A'*P*L)

X=(inv(A'*P*A))*(B'*K+A'*P*L)

N=(A'*P*A)

M=B*inv(N)*B'

Qx=inv(N)-inv(N)*B'*inv(M)*B*inv(N)

u=sqrt((V'*P*V)/(n-t+r))

fid=fopen('data.txt','wt'); %写入文件路径,文件输出

fprintf(fid,'=====================具有条件的参数平差: V=AX-L;BX+W=0 ===========================\n');

fprintf(fid,'输出n,t,r\n') ;

fprintf(fid,'%12.5f %12.5f %12.5f\n',n ,t,r); %n,t,r

fprintf(fid,'输出矩阵A\n');

[m,n]=size(A); %A

for i=1:1:m

for j=1:1:n

if j==n

fprintf(fid,'%12.5f\n',A(i,j));

else

fprintf(fid,'%12.5f\t',A(i,j));

end

end

end

fprintf(fid,'\n');

fprintf(fid,'输出矩阵L\n');

[m,n]=size(L); %L

for i=1:1:m

for j=1:1:n

if j==n

fprintf(fid,'%12.5f\n',L(i,j));

else

fprintf(fid,'%12.5f\t',L(i,j));

end

end

end

fprintf(fid,'\n');

fprintf(fid,'输出矩阵B\n');

[m,n]=size(B); %B

for i=1:1:m

for j=1:1:n

if j==n

fprintf(fid,'%12.5f\n',B(i,j));

fprintf(fid,'%12.5f\t',B(i,j));

end

end

end

fprintf(fid,'\n');

fprintf(fid,'输出矩阵W\n');

[m,n]=size(W); %W for i=1:1:m

for j=1:1:n

if j==n

fprintf(fid,'%12.5f\n',W(i,j));

else

fprintf(fid,'%12.5f\t',W(i,j));

end

end

end

fprintf(fid,'\n');

fprintf(fid,'输出矩阵P\n');

[m,n]=size(P); %P for i=1:1:m

for j=1:1:n

if j==n

fprintf(fid,'%12.5f\n',P(i,j));

else

fprintf(fid,'%12.5f\t',P(i,j));

end

end

end

fprintf(fid,'\n');

fprintf(fid,'输出矩阵X\n');

[m,n]=size(X); %X for i=1:1:m

for j=1:1:n

if j==n

fprintf(fid,'%12.5f\n',X(i,j));

else

fprintf(fid,'%12.5f\t',X(i,j));

end

end

end

fprintf(fid,'\n');

fprintf(fid,'输出矩阵V\n');

[m,n]=size(V); %V

for j=1:1:n

if j==n

fprintf(fid,'%12.5f\n',V(i,j));

else

fprintf(fid,'%12.5f\t',V(i,j));

end

end

end

fprintf(fid,'\n');

fprintf(fid,'输出矩阵Qx\n');

[m,n]=size(Qx); %Qx for i=1:1:m

for j=1:1:n

if j==n

fprintf(fid,'%12.5f\n',Qx(i,j));

else

fprintf(fid,'%12.5f\t',Qx(i,j));

end

end

end

fprintf(fid,'\n');

fprintf(fid,'输出单位权中误差\n');

[m,n]=size(u); %u for i=1:1:m

for j=1:1:n

if j==n

fprintf(fid,'%12.5f\n',u(i,j));

else

fprintf(fid,'%12.5f\t',u(i,j));

end

end

end

fprintf(fid,'\n');

fclose(fid);

case 4

% 具有参数的条件平差

% 具有参数的条件方程式BV+AX+W=0

sprintf('4:具有参数的条件平差BV+AX+W=0 ')

n=input('请输入n=');

t=input('请输入t=');

r=input('请输入r=');

B=input('请输入条件方程中V的系数矩阵B=');

A=input('请输入条件方程中X的系数矩阵A=');

W=input('请输入条件方程自由向量W=');

P=input('请输入权阵P=');

sprintf('开始计算')

N=(B*inv(P)*B')

M=(A'*inv(N)*A)

X=-inv(M)*A'*inv(N)*W

K=-inv(N)*(A*X+W)

V=inv(P)*B'*K

Qx=inv(M)

u=sqrt((V'*P*V)/(r-t))

fid=fopen('data.txt','wt'); %写入文件路径,文件输出

fprintf(fid,'=====================具有参数的条件平差: BV+AX+W=0 ================================\n');

fprintf(fid,'输出n,t,r\n') ;

fprintf(fid,'%12.5f %12.5f %12.5f\n',n ,t,r) ; %n,t,r

fprintf(fid,'输出矩阵B\n');

[m,n]=size(B); %B

for i=1:1:m

for j=1:1:n

if j==n

fprintf(fid,'%12.5f\n',B(i,j));

else

fprintf(fid,'%12.5f\t',B(i,j));

end

end

end

fprintf(fid,'\n');

fprintf(fid,'输出矩阵A\n');

[m,n]=size(A); %A

for i=1:1:m

for j=1:1:n

if j==n

fprintf(fid,'%12.5f\n',A(i,j));

else

fprintf(fid,'%12.5f\t',A(i,j));

end

end

end

fprintf(fid,'\n');

fprintf(fid,'输出矩阵W\n');

[m,n]=size(W); %W

for i=1:1:m

for j=1:1:n

if j==n

fprintf(fid,'%12.5f\n',W(i,j));

else

fprintf(fid,'%12.5f\t',W(i,j));

end

end

end

fprintf(fid,'\n');

fprintf(fid,'输出矩阵P\n');

[m,n]=size(P); %P for i=1:1:m

for j=1:1:n

if j==n

fprintf(fid,'%12.5f\n',P(i,j));

else

fprintf(fid,'%12.5f\t',P(i,j));

end

end

end

fprintf(fid,'\n');

fprintf(fid,'输出矩阵X\n');

[m,n]=size(X); %X for i=1:1:m

for j=1:1:n

if j==n

fprintf(fid,'%12.5f\n',X(i,j));

else

fprintf(fid,'%12.5f\t',X(i,j));

end

end

end

fprintf(fid,'\n');

fprintf(fid,'输出矩阵V\n');

[m,n]=size(V); %V for i=1:1:m

for j=1:1:n

if j==n

fprintf(fid,'%12.5f\n',V(i,j));

else

fprintf(fid,'%12.5f\t',V(i,j));

end

end

end

fprintf(fid,'\n');

fprintf(fid,'输出矩阵Qx\n');

[m,n]=size(Qx); %Qx for i=1:1:m

for j=1:1:n

if j==n

fprintf(fid,'%12.5f\n',Qx(i,j));

else

fprintf(fid,'%12.5f\t',Qx(i,j));

end

end

end

fprintf(fid,'\n');

fprintf(fid,'输出单位权中误差\n');

[m,n]=size(u); %u for i=1:1:m

for j=1:1:n

if j==n

fprintf(fid,'%12.5f\n',u(i,j));

else

fprintf(fid,'%12.5f\t',u(i,j));

end

end

end

fprintf(fid,'\n');

fclose(fid);

otherwise

disp('代号输入有误')

end

几种常见窗函数及其MATLAB程序实现

几种常见窗函数及其MATLAB程序实现 2013-12-16 13:58 2296人阅读评论(0) 收藏举报 分类: Matlab(15) 数字信号处理中通常是取其有限的时间片段进行分析,而不是对无限长的信号进行测量和运算。具体做法是从信号中截取一个时间片段,然后对信号进行傅里叶变换、相关分析等数学处理。信号的截断产生了能量泄漏,而用FFT算法计算频谱又产生了栅栏效应,从原理上讲这两种误差都是不能消除的。在FFT分析中为了减少或消除频谱能量泄漏及栅栏效应,可采用不同的截取函数对信号进行截短,截短函数称为窗函数,简称为窗。 泄漏与窗函数频谱的两侧旁瓣有关,对于窗函数的选用总的原则是,要从保持最大信息和消除旁瓣的综合效果出发来考虑问题,尽可能使窗函数频谱中的主瓣宽度应尽量窄,以获得较陡的过渡带;旁瓣衰减应尽量大,以提高阻带的衰减,但通常都不能同时满足这两个要求。 频谱中的如果两侧瓣的高度趋于零,而使能量相对集中在主瓣,就可以较为接近于真实的频谱。不同的窗函数对信号频谱的影响是不一样的,这主要是因为不同的窗函数,产生泄漏的大小不一样,频率分辨能力也不一样。信号的加窗处理,重要的问题是在于根据信号的性质和研究目的来选用窗函数。图1是几种常用的窗函数的时域和频域波形,其中矩形窗主瓣窄,旁瓣大,频率识别精度最高,幅值识别精度最低,如果仅要求精确读出主瓣频率,而不考虑幅值精度,则可选用矩形窗,例如测量物体的自振频率等;布莱克曼窗主瓣宽,旁瓣小,频率识别精度最低,但幅值识别精度最高;如果分析窄带信号,且有较强的干扰噪声,则应选用旁瓣幅度小的窗函数,如汉宁窗、三角窗等;对于随时间按指数衰减的函数,可采用指数窗来提高信噪比。表1 是几种常用的窗函数的比较。 如果被测信号是随机或者未知的,或者是一般使用者对窗函数不大了解,要求也不是特别高时,可以选择汉宁窗,因为它的泄漏、波动都较小,并且选择性也较高。但在用于校准时选用平顶窗较好,因为它的通带波动非常小,幅度误差也较小。

matlab代码大全

MATLAB主要命令汇总 MATLAB函数参考 附录1.1 管理用命令 函数名功能描述函数名功能描述 addpath 增加一条搜索路径 rmpath 删除一条搜索路径 demo 运行Matlab演示程序 type 列出.M文件 doc 装入超文本文档 version 显示Matlab的版本号 help 启动联机帮助 what 列出当前目录下的有关文件 lasterr 显示最后一条信息 whatsnew 显示Matlab的新特性 lookfor 搜索关键词的帮助 which 造出函数与文件所在的目录 path 设置或查询Matlab路径 附录1.2管理变量与工作空间用命令 函数名功能描述函数名功能描述 clear 删除存中的变量与函数 pack 整理工作空间存 disp 显示矩阵与文本 save 将工作空间中的变量存盘 length 查询向量的维数 size 查询矩阵的维数 load 从文件中装入数据 who,whos 列出工作空间中的变量名 附录1.3文件与操作系统处理命令 函数名功能描述函数名功能描述 cd 改变当前工作目录 edit 编辑.M文件 delete 删除文件 matlabroot 获得Matlab的安装根目录 diary 将Matlab运行命令存盘 tempdir 获得系统的缓存目录 dir 列出当前目录的容 tempname 获得一个缓存(temp)文件 ! 执行操作系统命令 附录1.4窗口控制命令 函数名功能描述函数名功能描述 echo 显示文件中的Matlab中的命令 more 控制命令窗口的输出页面format 设置输出格式 附录1.5启动与退出命令 函数名功能描述函数名功能描述 matlabrc 启动主程序 quit 退出Matlab环境 startup Matlab自启动程序 附录2 运算符号与特殊字符附录 2.1运算符号与特殊字符 函数名功能描述函数名功能描述 + 加 ... 续行标志 - 减 , 分行符(该行结果不显示) * 矩阵乘 ; 分行符(该行结果显示) .* 向量乘 % 注释标志 ^ 矩阵乘方 ! 操作系统命令提示符 .^ 向量乘方矩阵转置 kron 矩阵kron积 . 向量转置 \ 矩阵左除 = 赋值运算 / 矩阵右除 == 关系运算之相等 .\ 向量左除 ~= 关系运算之不等 ./ 向量右除 < 关系运算之小于

一个简单的Matlab_GUI编程实例

Matlab GUI编程教程(适用于初学者) 1.首先我们新建一个GUI文件:如下图所示; 选择Blank GUI(Default) 2.进入GUI开发环境以后添加两个编辑文本框,6个静态文本框,和一个按钮,布置如下

图所示; 布置好各控件以后,我们就可以来为这些控件编写程序来实现两数相加的功能了。3.我们先为数据1文本框添加代码; 点击上图所示红色方框,选择edit1_Callback,光标便立刻移到下面这段代码的位置。 1. 2. 3.function edit1_Callback(hObject, eventdata, handles) 4.% hObject handle to edit1 (see GCBO) 5.% eventdata reserved - to be defined in a future version of MATLAB

6.% handles structure with handles and user data (see GUIDATA) 7.% Hints: get(hObject,'String') returns contents of edit1 as text 8.% str2double(get(hObject,'String')) returns contents of edit1 as a double 复制代码 然后在上面这段代码的下面插入如下代码: 1. 2.%以字符串的形式来存储数据文本框1的内容. 如果字符串不是数字,则现实空白内容input = str2num(get(hObject,'String')); %检查输入是否为空. 如果为空,则默认显示为0if (isempty(input)) set(hObject,'String','0')endguidata(hObject, handles); 复制代码 这段代码使得输入被严格限制,我们不能试图输入一个非数字。 4.为edit2_Callback添加同样一段代码 5 现在我们为计算按钮添加代码来实现把数据1和数据2相加的目的。 用3中同样的方法在m文件中找到pushbutton1_Callback代码段 如下; 1.function pushbutton1_Callback(hObject, eventdata, handles) 2.% hObject handle to pushbutton1 (see GCBO) 3.% eventdata reserved - to be defined in a future version of MATLAB 4.% handles structure with handles and user data (see GUIDATA) 复制代码

matlab编程实现求解最优解

《现代设计方法》课程 关于黄金分割法和二次插值法的Matlab语言实现在《现代设计方法》的第二章优化设计方法中有关一维搜索的最优化方法的 一节里,我们学习了黄金非分割法和二次插值法。它们都是建立在搜索区间的优先确定基础上实现的。 为了便于方便执行和比较,我将两种方法都写进了一个程序之内,以选择的方式实现执行其中一个。下面以《现代设计方法》课后习题为例。见课本70页,第2—7题。原题如下: 求函数f(x)=3*x^2+6*x+4的最优点,已知单谷区间为[-3,4],一维搜索精度为0.4。 1、先建立函数f(x),f(x)=3*x^2+6*x+4。函数文件保存为:lee.m 源代码为:function y=lee(x) y=3*x^2+6*x+4; 2、程序主代码如下,该函数文件保存为:ll.m clear; a=input('请输入初始点'); b=input('请输入初始步长'); Y1=lee(a);Y2=lee(a+b); if Y1>Y2 %Y1>Y2的情况 k=2; Y3=lee(a+2*b); while Y2>=Y3 %直到满足“大,小,大”为止 k=k+1; Y3=lee(a+k*b); end A=a+b;B=a+k*b; elseif Y1=Y3 %直到满足“大,小,大”为止 k=k+1; Y3=lee(a-k*b); end A=a-k*b;B=a; else A=a;B=a+b; %Y1=Y2的情况 end disp(['初始搜索区间为',num2str([A,B])])%输出符合的区间 xuanze=input('二次插值法输入0,黄金分割法输入1');%选择搜索方式 T=input('选定一维搜索精度'); if xuanze==1 while B-A>T %一维搜索法使精度符合要求 C=A+0.382*(B-A);D=A+0.618*(B-A); %黄金分割法选点

matlab实现中值滤波去除脉冲噪声matlab小程序

matlab实现中值滤波去除脉冲噪声matlab小程序(图像处理)2010-04-1612:58:44阅读8评论0字号:大中小 实验原理:中值滤波器是将领域内像素灰度的中值代替该像素的值,对处理脉冲噪声(椒盐噪声)非常有效。为了对一幅图像上的某个点进行中值滤波处理,必须先将掩模内欲求的像素及其领域的像素值排序,确定出中值,主要功能是使拥有不同灰度的点看起来更接近于它的邻近值。 程序说明:函数名为mid(pic_name,s)的函数,其中参数pic_name为读入的图像,s为掩模矩阵的边长,由用户自行决定。 实验说明:随着掩模矩阵的变大,我们可以看到脉冲噪声去除得更加理想,但同时图像会变得更模糊,因为各点像素与其邻域更为接近,因此,进行中值滤波时选择一个适合的掩模矩阵十分重要。另外,我们看到图像的边界处出现了黑色的斑点,这是由于我采用了0来直译边界,这种影响可用镜像反射方式对称地沿其边界扩展来减弱。 另附:其实本实验可以完全由matlab中的函数median或medfilt2简单实现,此处写出内部处理过程,主要是为了让大家理解中值滤波的具体处理过程。 程序源代码: function mid(pic_name,s) close all; s=double(s); X=imread(pic_name); Y1=imnoise(X,'salt&pepper',0.2);%对读入的图像加脉冲噪声 figure; imshow(uint8(Y1)); Y1=double(Y1); [m,n]=size(X); s2=round(s/2); s3=round(s*s/2);%中值像素点的位置

数字信号处理Matlab实现实例(推荐给学生)

数字信号处理Matlab 实现实例 第1章离散时间信号与系统 例1-1 用MATLAB计算序列{-2 0 1 –1 3}和序列{1 2 0 -1}的离散卷积。 解 MATLAB程序如下: a=[-2 0 1 -1 3]; b=[1 2 0 -1]; c=conv(a,b); M=length(c)-1; n=0:1:M; stem(n,c); xlabel('n'); ylabel('幅度'); 图1.1给出了卷积结果的图形,求得的结果存放在数组c中为:{-2 -4 1 3 1 5 1 -3}。 例1-2 用MATLAB计算差分方程 当输入序列为时的输出结果。 解 MATLAB程序如下: N=41; a=[0.8 -0.44 0.36 0.22]; b=[1 0.7 -0.45 -0.6]; x=[1 zeros(1,N-1)];

k=0:1:N-1; y=filter(a,b,x); stem(k,y) xlabel('n');ylabel('幅度') 图 1.2 给出了该差分方程的前41个样点的输出,即该系统的单位脉冲响应。 例1-3 用MATLAB 计算例1-2差分方程 所对应的系统函数的DTFT 。 解 例1-2差分方程所对应的系统函数为: 123 123 0.80.440.360.02()10.70.450.6z z z H z z z z -------++= +-- 其DTFT 为 23230.80.440.360.02()10.70.450.6j j j j j j j e e e H e e e e ωωωω ωωω--------++= +-- 用MATLAB 计算的程序如下: k=256; num=[0.8 -0.44 0.36 0.02]; den=[1 0.7 -0.45 -0.6]; w=0:pi/k:pi; h=freqz(num,den,w); subplot(2,2,1); plot(w/pi,real(h));grid title('实部') xlabel('\omega/\pi');ylabel('幅度')

CA码生成原理及matlab程序实现

作业:用Matlab写C/A码生成器程序,并画生成码的方波图。 C/A码生成原理 C/A 码是用m 序列优选对组合形成的Gold 码。Gold码是由两个长度相同而互相关极大值为最小的m 序列逐位模2 相加所得到的码序列。它是由两个10 级反馈移位寄存器组合产生的,其产生原理如图1 所示。 图1 C/A码生成原理 发生器的抽头号为3和10,发生器的抽头号为2、3、6、8、9、10;发生器的第10位输出的数字即为码,而码是由的两个抽头的输出结果进行模2相加得到。 卫星的PRN码与延时的量是相关联的,对C/A码来说,每颗卫星都有特别的延时,如第1颗GPS卫星的G2 抽为2、6,第2颗为3、7,第3 颗为4、8,第4 颗为5、9 等,如图2所示。通过G2 相位选择可以产生结构不同的伪随机码,从而可以实现不同卫星之间的码分多址技术与卫星识别。

图2 prn序号与G2抽头、时延对应关系 基于MATLAB的GPS信号实现 编写成“codegen”程序,输入[ca_used]=codegen(svnum),其中svnum为卫星号,ca_used 为得到的C/A码序列。程序具体实现流程如下: 在程序中定义一个数组,使得卫星号与G2的码片延时一一对应。 gs2=[5;6;7;8;17;18;139;140;141;251;252;254;255;256;257;258;469;470;471;472;473;474;509;512 ;513;514;515;516;859;860;861;862]; 定义两个1×1 023 的数组g1、g2 用来存放生成的Gold 码。定义一个全1 的10 位数组,作为移位寄存器,相当于G1、G2 生成模块的初值均置为全“1”。按原理式

matlab代码大全教学文案

m a t l a b代码大全

MATLAB主要命令汇总 MATLAB函数参考 附录1.1 管理用命令 函数名功能描述函数名功能描述 addpath 增加一条搜索路径 rmpath 删除一条搜索路径 demo 运行Matlab演示程序 type 列出.M文件 doc 装入超文本文档 version 显示Matlab的版本号 help 启动联机帮助 what 列出当前目录下的有关文件 lasterr 显示最后一条信息 whatsnew 显示Matlab的新特性 lookfor 搜索关键词的帮助 which 造出函数与文件所在的目录 path 设置或查询Matlab路径 附录1.2管理变量与工作空间用命令 函数名功能描述函数名功能描述 clear 删除内存中的变量与函数 pack 整理工作空间内存 disp 显示矩阵与文本 save 将工作空间中的变量存盘 length 查询向量的维数 size 查询矩阵的维数 load 从文件中装入数据 who,whos 列出工作空间中的变量名 附录1.3文件与操作系统处理命令 函数名功能描述函数名功能描述 cd 改变当前工作目录 edit 编辑.M文件 delete 删除文件 matlabroot 获得Matlab的安装根目录 diary 将Matlab运行命令存盘 tempdir 获得系统的缓存目录 dir 列出当前目录的内容 tempname 获得一个缓存(temp)文件 ! 执行操作系统命令 附录1.4窗口控制命令 函数名功能描述函数名功能描述 echo 显示文件中的Matlab中的命令 more 控制命令窗口的输出页面format 设置输出格式 附录1.5启动与退出命令 函数名功能描述函数名功能描述 matlabrc 启动主程序 quit 退出Matlab环境 startup Matlab自启动程序 附录2 运算符号与特殊字符附录

基于matlab程序实现人脸识别

基于m a t l a b程序实现 人脸识别 TYYGROUP system office room 【TYYUA16H-TYY-TYYYUA8Q8-

基于m a t l a b程序实现人脸识别 1.人脸识别流程 基于YCbCr颜色空间的肤色模型进行肤色分割。在YCbCr色彩空间内对肤色进行了建模发现,肤色聚类区域在Cb—Cr子平面上的投影将缩减,与中心区域显着不同。采用这种方法的图像分割已经能够较为精确的将人脸和非人脸分割开来。 人脸识别流程图 2.人脸识别程序 (1)人脸和非人脸区域分割程序 function result = skin(Y,Cb,Cr) %SKIN Summary of this function goes here % Detailed explanation goes here a=; b=; ecx=; ecy=; sita=; cx=; cy=; xishu=[cos(sita) sin(sita);-sin(sita) cos(sita)]; %如果亮度大于230,则将长短轴同时扩大为原来的倍 if(Y>230) a=*a; b=*b; end %根据公式进行计算 Cb=double(Cb); Cr=double(Cr);

t=[(Cb-cx);(Cr-cy)]; temp=xishu*t; value=(temp(1)-ecx)^2/a^2+(temp(2)-ecy)^2/b^2; %大于1则不是肤色,返回0;否则为肤色,返回1 if value>1 result=0; else result=1; end end (2)人脸的确认程序 function eye = findeye(bImage,x,y,w,h) %FINDEYE Summary of this function goes here % Detailed explanation goes here part=zeros(h,w); %二值化 for i=y:(y+h) for j=x:(x+w) if bImage(i,j)==0 part(i-y+1,j-x+1)=255; else part(i-y+1,j-x+1)=0; end end end [L,num]=bwlabel(part,8); %如果区域中有两个以上的矩形则认为有眼睛 if num<2 eye=0;

(仅供参考)Matlab编写与调用函数

MATLAB 学习指南 第六章.编写与调用函数 在这一章中,我们讨论如何用多源代码文件来构造一个程序。首先,解释代码文件在MATLAB中如何工作。在编译语言中,例如FORTRAN,C ,或C++,代码被存储在一个或多个源文件中,在进行编译的时候,这些源文件组合在一起 形成了一个单独的可执行文件。作为一种解释型语言,MATLAB以一种更广泛的方式来处理多个源文件。MATLAB代码被放入带有扩展名.m的ASCII文件(或称m-文件)中。MATLAB 6 有一个集成字处理与调试应用程序,尽管会用到其它编辑程序如vi或emacs,集成字处理与调试应用程序仍是编译m-文件的首选程序。 有两种不同的m-文件。一种是脚本文件,它是一种最简单的文件,仅仅将MATLAB中的指令收集在一起。当在交互提示符处输入文件名执行脚本文件时,MATLAB在m-文件内读取并执行指令,就好像指令是我们输入的。而且,似乎我们能够削减m-文件的内容并将削减过的内容传到MATLAB指令窗口中。这种m-文件的用法将在6.1节中给予概述。 在6.2节中要讨论的第二种m-文件包含一个单一函数,此函数名与此m-文件名相同。这种m-文件包含一段独立的代码,这段代码具有一个明确规定的输入/输出界面;那就是说,传给这段代码一列空变量arg1,arg2,…,这段独立代码就能够被调用,然后返回输出值out1,out2,…。一个函数m-文件的第一个非注释行包含函数标头,其形式如下: 此m-文件以返回指令结束,将执行程序返回到函数被调用的位置。或者在交互指令提示符处或者在另一个m-文件内,无论何时用下列指令调用函数代码,函数代码都将被执行。 输入映射到空变量:arg1=var1,arg2=var2,等等。在函数主体内,输出值被分配给了变量out1,out2,等等。当遇到返回值时,当前值out1,out2,…在函数被调用处被映射到变量outvar1,outvar2,…。在用可变长度自变量和输出变量列表编写函数时,MATLAB允许更多的自由。例如,也可以使用下列指令来调用函数。 在此情况下,仅返回一个单一输出变量,这个变量在出口处包含函数变量out1的值。输入和输出自变量可能是字符串,数值,向量,矩阵,或者更高级的数据结构。 为什么使用函数呢?因为从每门计算机科学课程中可知,把一个大的程序分割 成多个可以单独执行一个被明确规定的和被注释过的任务的小程序会使大程序 易读,易于修改,不易于出错。在MATLAB中,先为程序编写一个主文件,或者是一个脚本文件或者更好的话,是一个能够返回一个单一整数的函数m-文件(返回1表示程序执行成功,0表示不完全程序执行,负值表示出现运行误差),这个主文件是程序的进入点。通过把m-文件当作函数来调用,此程序文件可以

matlab的编码大全

附录Matlab源程序 附录A 信息熵 % 函数说明:% % H=entropy(P,r) 为信息熵函数% % P为信源的概率矢量, r为进制数% % H为信息熵% %****************************** % function H=entropy(P,r) if (length(find(P<=0))~=0) error('Not a prob.vector,negative component'); % 判断是否符合概率分布条件end if (abs(sum(P)-1)>10e-10) error('Not a prob.vector,component do not add up to 1'); end H=(sum(-P.*log2(P)))/(log2(r)+eps); 附录B 离散无记忆信道容量的迭代计算 % 信道容量C的迭代算法% % 函数说明:% % [CC,Paa]=ChannelCap(P,k) 为信道容量函数% % 变量说明:% % P:输入的正向转移概率矩阵,k:迭代计算精度% % CC:最佳信道容量,Paa:最佳输入概率矩阵% % Pa:初始输入概率矩阵,Pba:正向转移概率矩阵% % Pb:输出概率矩阵,Pab:反向转移概率矩阵% % C:初始信道容量,r:输入符号数,s:输出符号数% %************************************************** % function [CC,Paa]=ChannelCap(P,k) % 提示错误信息 if (length(find(P<0)) ~=0) error('Not a prob.vector,negative component'); % 判断是否符合概率分布条件end

matlab源代码实例

1.硬币模拟试验 源代码: clear; clc; head_count=0; p1_hist= [0]; p2_hist= [0]; n = 1000; p1 = 0.3; p2=0.03; head = figure(1); rand('seed',sum(100*clock)); fori = 1:n tmp = rand(1); if(tmp<= p1) head_count = head_count + 1; end p1_hist (i) = head_count /i; end figure(head); subplot(2,1,1); plot(p1_hist); grid on; hold on; xlabel('重复试验次数'); ylabel('正面向上的比率'); title('p=0.3试验次数N与正面向上比率的函数图'); head_count=0; fori = 1:n tmp = rand(1); if(tmp<= p2) head_count = head_count + 1; end p2_hist (i) = head_count /i; end figure(head); subplot(2,1,2); plot(p2_hist); grid on; hold on; xlabel('重复试验次数'); ylabel('正面向上的比率'); title('p=0.03试验次数N与正面向上比率的函数图'); 实验结果:

2.不同次数的随机试验均值方差比较 源代码: clear ; clc; close; rand('seed',sum(100*clock)); Titles = ['n=5时' 'n=20时' 'n=25时' 'n=50时' 'n=100时']; Titlestr = cellstr(Titles); X_n_bar=[0]; %the samples of the X_n_bar X_n=[0]; %the samples of X_n N=[5,10,25,50,100]; j=1; num_X_n = 100; num_X_n_bar = 100; h_X_n_bar = figure(1);

遗传算法的原理及MATLAB程序实现

遗传算法的原理及MATLAB程序实现 1 遗传算法的原理 1.1 遗传算法的基本思想 遗传算法(genetic algorithms,GA)是一种基于自然选择和基因遗传学原理,借鉴了生物进化优胜劣汰的自然选择机理和生物界繁衍进化的基因重组、突变的遗传机制的全局自适应概率搜索算法。 遗传算法是从一组随机产生的初始解(种群)开始,这个种群由经过基因编码的一定数量的个体组成,每个个体实际上是染色体带有特征的实体。染色体作为遗传物质的主要载体,其内部表现(即基因型)是某种基因组合,它决定了个体的外部表现。因此,从一开始就需要实现从表现型到基因型的映射,即编码工作。初始种群产生后,按照优胜劣汰的原理,逐代演化产生出越来越好的近似解。在每一代,根据问题域中个体的适应度大小选择个体,并借助于自然遗传学的遗传算子进行组合交叉和变异,产生出代表新的解集的种群。这个过程将导致种群像自然进化一样,后代种群比前代更加适应环境,末代种群中的最优个体经过解码,可以作为问题近似最优解。 计算开始时,将实际问题的变量进行编码形成染色体,随机产生一定数目的个体,即种群,并计算每个个体的适应度值,然后通过终止条件判断该初始解是否是最优解,若是则停止计算输出结果,若不是则通过遗传算子操作产生新的一代种群,回到计算群体中每个个体的适应度值的部分,然后转到终止条件判断。这一过程循环执行,直到满足优化准则,最终产生问题的最优解。图1-1给出了遗传算法的基本过程。 1.2 遗传算法的特点 1.2.1 遗传算法的优点

遗传算法具有十分强的鲁棒性,比起传统优化方法,遗传算法有如下优点: 1. 遗传算法以控制变量的编码作为运算对象。传统的优化算法往往直接利用控制变量的实际值的本身来进行优化运算,但遗传算法不是直接以控制变量的值,而是以控制变量的特定形式的编码为运算对象。这种对控制变量的编码处理方式,可以模仿自然界中生物的遗传和进化等机理,也使得我们可以方便地处理各种变量和应用遗传操作算子。 2. 遗传算法具有内在的本质并行性。它的并行性表现在两个方面,一是遗传 开始 初始化,输入原始参 数及给定参数,gen=1 染色体编码,产生初始群体 计算种群中每个个体的适应值 终止条件的判断, N gen=gen+1 选择 交叉 Y 变异 新种群 输出结果 结束 图1-1 简单遗传算法的基本过程

MATLAB小程序:将TXT中十六进制数转为十进制输出

matlab小程序:将txt中十六进制数转为十进制输出function htod(filename) clc [n]=textread(filename,'%2c'); [a b]=size(n) m=zeros(a,b); mm=zeros(a,1); for i=1:a for j=1:b switch n(i,j) case{'0'}m(i,j)=0; case{'1'}m(i,j)=1; case{'2'}m(i,j)=2; case{'3'}m(i,j)=3; case{'4'}m(i,j)=4; case{'5'}m(i,j)=5; case{'6'}m(i,j)=6; case{'7'}m(i,j)=7; case{'8'}m(i,j)=8; case{'9'}m(i,j)=9; case{'A'}m(i,j)=10; case{'B'}m(i,j)=11; case{'C'}m(i,j)=12; case{'D'}m(i,j)=13; case{'E'}m(i,j)=14; case{'F'}m(i,j)=15; otherwise m(i,j)=nan; end end end %m for i=1:a for j=1:b mm(i)=mm(i)+m(i,j)*16^(j-1); end end %mm [a b]=size(mm); size_mm=a mmm=mm'; savefile='C:\Documents and Settings\Administrator\桌面\test.txt'; fid=fopen(savefile,'w');

fprintf(fid,'%4d',mmm) fclose(fid); matlab如何读取二进制、十六进制txt文档 发现matlab如何读取十六进制的和二进制的txt文章不多。今天刚想了一种方法,所以在这里小结一下,所以matlab中文论坛共享一下,没有参考其他的文章哦,觉得好用就帮顶,不好用提意见。 原帖地址https://www.360docs.net/doc/db14130687.html,/thread-23226-1-1.html 本方法同样适合读取十六进制和二进制以外的其他进制文件, txt使用一个最简单的命令就可以读取textread这是一个十分有用,简便的函数(对于fopen fscanf而言) 读取二进制txt文件: 假如txt文档中内容为00010010001101001000,保存在pin.txt文档中 使用a=textread('pin.txt','%s')' a= '0001''0010''0011''0100''1000' 可以看到数据保存为了char格式。 使用bin2dec b=bin2dec(a)' b= 12348 可以看到成功地转换成了十进制文件。 十六进制文件: 00010010001101001000A B C AA a=textread('pin.txt','%s')' a= '0001''0010''0011''0100''1000''A''B''C''AA' 可以看到成功读取了文件。 b=hex2dec(a)' b= 11617256409610 1112170 读取完毕。 小结:本方法以简单使用方便的方法读取二进制、十六进制的txt文档,欢迎大家提出意见

matlab 指令大全

分享 我的分享 当前分享 返回分享首页? 分享 matlab命令,应该很全了!来源:李家叶的日志 matlab命令 一、常用对象操作:除了一般windows窗口的常用功能键外。 1、!dir 可以查看当前工作目录的文件。!dir& 可以在dos状态下查看。 2、who 可以查看当前工作空间变量名,whos 可以查看变量名细节。 3、功能键: 功能键快捷键说明 方向上键Ctrl+P 返回前一行输入 方向下键Ctrl+N 返回下一行输入 方向左键Ctrl+B 光标向后移一个字符 方向右键Ctrl+F 光标向前移一个字符 Ctrl+方向右键Ctrl+R 光标向右移一个字符 Ctrl+方向左键Ctrl+L 光标向左移一个字符 home Ctrl+A 光标移到行首 End Ctrl+E 光标移到行尾 Esc Ctrl+U 清除一行 Del Ctrl+D 清除光标所在的字符 Backspace Ctrl+H 删除光标前一个字符 Ctrl+K 删除到行尾 Ctrl+C 中断正在执行的命令 4、clc可以命令窗口显示的内容,但并不清除工作空间。 二、函数及运算 1、运算符: +:加,-:减,*:乘,/:除,\:左除^:幂,‘:复数的共轭转置,():制定运算顺序。 2、常用函数表: sin( ) 正弦(变量为弧度) Cot( ) 余切(变量为弧度) sind( ) 正弦(变量为度数) Cotd( ) 余切(变量为度数) asin( ) 反正弦(返回弧度) acot( ) 反余切(返回弧度) Asind( ) 反正弦(返回度数) acotd( ) 反余切(返回度数) cos( ) 余弦(变量为弧度) exp( ) 指数 cosd( ) 余弦(变量为度数)

基于Matlab的动态规划程序实现

动态规划方法的Matlab 实现与应用 动态规划(Dynamic Programming)是求解决策过程最优化的有效数学方法,它是根据“最优决策的任何截断仍是最优的”这最优性原理,通过将多阶段决策过程转化为一系列单段决策问题,然后从最后一段状态开始逆向递推到初始状态为止的一套最优化求解方法。 1.动态规划基本组成 (1) 阶段 整个问题的解决可分为若干个阶段依次进行,描述阶段的变量称为阶段变量,记为k (2) 状态 状态表示每个阶段开始所处的自然状况或客观条件,它描述了研究问题过程的状况。各阶段状态通常用状态变量描述,用k x 表示第k 阶段状态变量,n 个阶段决策过程有n+ 1个状态。 (3) 决策 从一确定的状态作出各种选择从而演变到下一阶段某一状态,这种选择手段称为决策。描述决策的变量称为决策变量,决策变量限制的取值范围称为允许决策集合。用()k k u x 表示第k 阶段处于状态k x 时的决策变量,它是k x 的函数。用()k k D x Dk(xk)表示k x 的允许决策的集合。 (4) 策略 每个阶段的决策按顺序组成的集合称为策略。由第k 阶段的状态k x 开始到终止状态的后部子过程的策略记为{}11(),(),,()k k k k n n u x u x u x ++ 。可供选择的策略的范围称为允许策略集合,允许策略集合中达到最优效果的策略称为最优策略。从初始状态* 11()x x =出发,过程按照最优策略和状态转移方程演变所经历的状态序列{ } **** 121,,,,n n x x x x + 称为最优轨线。 (5) 状态转移方程 如果第k 个阶段状态变量为k x ,作出的决策为k u ,那么第k+ 1阶段的状态变量1k x +也被完全确定。用状态转移方程表示这种演变规律,记为1(,)k k k x T x u +=。 (6) 指标函数 指标函数是系统执行某一策略所产生结果的数量表示,是衡量策略优劣的数量指标,它定义在全过程和所有后部子过程上,用()k k f x 表示。过程在某阶段j 的阶段指标函数是衡量该阶段决策优劣数量指标,取决于状态j x 和决策j u ,用(,)j j j v x u 表示。 2.动态规划基本方程 (){} 11()min ,,(),()k k k k k k k k k k f x g v x u f x u D x ++=∈???? Matlab 实现 (dynprog.m 文件) function [p_opt,fval]=dynprog (x,DecisFun,SubObjFun,TransFun,ObjFun) % x 是状态变量,一列代表一个阶段的所有状态; % M-函数DecisFun(k,x) 由阶段k 的状态变量x 求出相应的允许决策变量; % M-函数SubObjFun(k,x,u) 是阶段指标函数, % M-函数ObjFun(v,f) 是第k 阶段至最后阶段的总指标函数 % M-函数TransFun(k,x,u) 是状态转移函数, 其中x 是阶段k 的某状态变量, u 是相应的决策变量; %输出 p_opt 由4列构成,p_opt=[序号组;最优策略组;最优轨线组;指标函数值组]; %输出 fval 是一个列向量,各元素分别表示p_opt 各最优策略组对应始端状态x 的最优函数值。

数学建模基础入门小程序文件

自己整理MATLAB知识 1入门 例1-1 绘制正弦曲线和余弦曲线。 x=[0:0.5:360]*pi/180; plot(x,sin(x),x,cos(x)); 例1-2 求方程3x4+7x3+9x2-23=0的全部根。 p=[3,7,9,0,-23]; %建立多项式系数向量 x=roots(p) %求根 例1-3 求积分 quad('x.*log(1+x)',0,1) %‘里是被积函数’0,1分 别是积分上下限 例1-4 求解线性方程组。 a=[2,-3,1;8,3,2;45,1,-9]; %方程左面系数 b=[4;2;17]; %方程右面系数 x=inv(a)*b %也可是x=a\b的形式 例1-5 水仙花 for m=100:999 m1=fix(m/100); %求m的百位数字 m2=rem(fix(m/10),10); %求m的十位数字 m3=rem(m,10); %求m的个位数字 if m==m1*m1*m1+m2*m2*m2+m3*m3*m3 disp(m)

end end 例1-6 已知,当n=100时,求y的值。程序如下: y=0; n=100; for i=1:n y=y+1/(2*i-1); end y 例1-7 求[100,200]之间第一个能被21整除的整数 for n=100:200 if rem(n,21)~=0 continue end break end n 例1-8 若一个数等于它的各个真因子之和,则称该数为完数,如6=1+2+3,所以6是完数。求[1,500]之间的全部完数。for m=1:500 s=0; for k=1:m/2

实验一 典型环节的MATLAB仿真汇总

实验一 典型环节的MATLAB 仿真 一、实验目的 1.熟悉MATLAB 桌面和命令窗口,初步了解SIMULINK 功能模块的使用方法。 2.通过观察典型环节在单位阶跃信号作用下的动态特性,加深对各典型环节响应曲线的理解。 3.定性了解各参数变化对典型环节动态特性的影响。 二、SIMULINK 的使用 MATLAB 中SIMULINK 是一个用来对动态系统进行建模、仿真和分析的软件包。利用SIMULINK 功能模块可以快速的建立控制系统的模型,进行仿真和调试。 1.运行MATLAB 软件,在命令窗口栏“>>”提示符下键入simulink 命令,按Enter 键或在工具栏单击按钮,即可进入如图1-1所示的SIMULINK 仿真 环境下。 2.选择File 菜单下New 下的Model 命令,新建一个simulink 仿真环境常规模板。 3.在simulink 仿真环境下,创建所需要的系统 三、实验内容 按下列各典型环节的传递函数,建立相应的SIMULINK 仿真模型,观察并记录其单位阶跃响应波形。 ① 比例环节1)(1=s G 和2)(1=s G 实验处理:1)(1=s G SIMULINK 仿真模型

波形图为: 实验处理:2)(1=s G SIMULINK 仿真模型 波形图为: 实验结果分析:增加比例函数环节以后,系统的输出型号将输入信号成倍数放大. ② 惯性环节11)(1+= s s G 和15.01)(2+=s s G 实验处理:1 1 )(1+=s s G SIMULINK 仿真模型

波形图为: 实验处理:1 5.01 )(2+= s s G SIMULINK 仿真模型 波形图为: 实验结果分析:当1 1 )(1+= s s G 时,系统达到稳定需要时间接近5s,当

如何编写MATLAB程序才能实现对

关闭文件用fclose函数,调用格式为:sta=fclose(fid)说明:该函数关闭fid所表示的文件。其调用格式为:[A,COUNT]=fscanf(fid,format,size)说明:其中A用来存放读取的数据,COUNT返回所读取的数据元素个数,fid为文件句柄,format用来控制读取的数据格式,由%加上格式符组成,常见的格式符有:d(整型)、f(浮点型)、s(字符串型)、c(字符型)等,在%与格式符之间还可以插入附加格式说明符,如数据宽度说明等。 matlab fprintf.数据的格式化输出:fprintf(fid, format, variables)fprintf(fid,format,A)说明:fid为文件句柄,指定要写入数据的文件,format是用来控制所写数据格式的格式符,与fscanf函数相同,A是用来存放数据的矩阵。>> fid=fopen(""d:\char1.txt"",""w"");>> fid1=fopen(""d:\char1.txt"",""rt"");matlab读txt文件fid=fopen(""fx.txt"",""r"");%得到文件号[f,count]=fscanf(fid,""%f %f"",[12,90]);%把文件号1的数据读到f中。 matlab函数fgetl和fgets:按行读取格式文本函数Matlab提供了两个函数fgetl和fgets来从格式文本文件读取行,并存储到字符串向量中。这两个函数集几乎相同;不同之处是,fgets拷贝新行字符到字符向量,而fgetl则不。下面的M-file函数说明了fgetl的一个可能用法。此函数使用fgetl一次读取一整行。while f eof(fid) == 0 tline = fgetl(fid); %用Fourier变换求取信号的功率谱---周期图法 clf; Fs=1000; N=256;Nfft=256;%数据的长度和FFT所用的数据长度 n=0:N-1;t=n/Fs;%采用的时间序列 xn=sin(2*pi*50*t)+2*sin(2*pi*120*t)+randn(1,N); Pxx=10*log10(abs(fft(xn,Nfft).^2)/N);%Fourier振幅谱平方的平均值,并转化为dB f=(0:length(Pxx)-1)*Fs/length(Pxx);%给出频率序列 subplot(2,1,1),plot(f,Pxx);%绘制功率谱曲线 xlabel('频率/Hz');ylabel('功率谱/dB'); title('周期图N=256');grid on; Fs=1000; N=1024;Nfft=1024;%数据的长度和FFT所用的数据长度 n=0:N-1;t=n/Fs;%采用的时间序列 xn=sin(2*pi*50*t)+2*sin(2*pi*120*t)+randn(1,N); Pxx=10*log10(abs(fft(xn,Nfft).^2)/N);%Fourier振幅谱平方的平均值,并转化为dB f=(0:length(Pxx)-1)*Fs/length(Pxx);%给出频率序列 subplot(2,1,2),plot(f,Pxx);%绘制功率谱曲线 xlabel('频率/Hz');ylabel('功率谱/dB'); title('周期图N=256');grid on; %用Fourier变换求取信号的功率谱---分段周期图法 %思想:把信号分为重叠或不重叠的小段,对每小段信号序列进行功率谱估计,然后取平均值作为整个序列的功率谱 clf;

相关文档
最新文档