Matlab实现HHT程序(源码-非常珍贵)

合集下载

matlab经典源程序带有注释(详细经典)

matlab经典源程序带有注释(详细经典)

2.1set 与get 函数 (1)2.2callback函数 (2)2.3元胞数组 (4)2.4结构数组 (6)2.5矩阵操作 (9)2.6字符串操作 (13)2.7判断函数使用大全 (16)2.11打开外部程序 (21)2.11程序运行时间 (22)2.14动画 (23)2.12动画 (24)2.23显示多行内容 (26)2.24 uitable 使用 (26)2.27鼠标操作 (27)2.28键盘操作 (27)2.32粘贴板 (28)2.1set 与get 函数set(edit_handle,'String','my value!'); %String为Edit控件的属性%%%2.1-1%创建figure对象hfig=figure(1);%创建坐标轴对象,指定其父对象为figure 1haxes1=axes('parent',hfig);prop.Color='b';prop.FontSize=12;set(haxes1,prop);%%%2.1-2hfig=figure(1);%查询其Units属性值get(hfig,'units')%其Units属性值为pixels(像素)% ans=% pixels%%%2.1-3%figure的Pointer属性标识了鼠标指针的形状set(gcf,'pointer');% 返回值为:[ crosshair | fullcrosshair | {arrow} | ibeam | watch | topl | topr | botl | botr | left | top | right | bottom | circle | cross | fleur | custom | hand ]%%%2.1-4%首先取得标识电脑屏幕大小的度量单位get(0,'units')% ans =% pixels%取得屏幕的尺寸get(0,'screensize')% ans =% 1 1 1280 8002.2callback函数%定义M文件的主函数名称为DefineCallback,不带输入和输出参数function DefineCallbackhFig= figure('units','normalize',...'position',[0.4 0.4 0.3 0.2]);%在窗口中创建按钮控件,并定义其Callback属性uicontrol('parent',hFig,...'style','pushbutton',...'String','Execute Callback',...'units','normalize',...'position',[0.4 0.4 0.3 0.2],...'callback',['figure;',...'x = 0:pi/20:2*pi;',...'y = sin(x);',...'plot(x,y);']);%定义M文件的主函数名称为DefineCallback,不带输入和输出参数function DefineCallback%创建界面窗口hFig= figure('units','normalize',...'position',[0.4 0.4 0.3 0.2]);%在窗口中创建按钮控件hpush=uicontrol('parent',hFig,...'style','pushbutton',...'String','Execute Callback',...'units','normalize',...'position',[0.4 0.4 0.3 0.2]);%设置按钮的Callback属性set(hpush,'callback',@mycallback);%定义回调函数为子函数function mycallback(hobj,event)figure;x = 0:pi/20:2*pi;y = sin(x);plot(x,y);2.3元胞数组a={'hello' [1 2 3;4 5 6];1 {'1''2'}}a ='hello' [2x3 double][ 1] {1x2 cell }%示例2:将元胞数组a中的元胞逐一赋值>> a{1,1}='hello';a{1,2}=[1 2 3;4 5 6];a{2,1}=1;a{2,2}={'1' '2'};>> aa ='hello' [2x3 double][ 1] {1x2 cell }%示例3:使用cell函数来创建元胞数组%生成2x3的元素为空数组的元胞数组>> a=cell(2,3)a =[] [] [][] [] []%示例4:判断数组A是否为元胞数组%定义一个元胞数组A>> A={1 2 3};%判断A是否为元胞数组,如果为元胞数组,则函数>> tf = iscell(A)tf =1%示例5:显示元胞数组C中的内容>> clear>> C={'Smith' [1 2;3 4] [12]};%直接显示元胞数组C中的内容>> celldisp(C)C{1} =SmithC{2} =1 23 4C{3} =12%显示元胞数组C中的内容,数组的名称用cellcontent代替>> celldisp(C,'cellcontent')cellcontent{1} =Smithcellcontent{2} =1 23 4cellcontent{3} =12%示例6:将字符数组转换为元胞数组>> S = ['abc '; 'defg'; 'hi m'];>> cellstr(S)ans ='abc'%原先abc后面的空格被清除'defg''hi m'%i和m之间的空格仍然保留%示例7:显示元胞数组S中的内容(包括空格和字符)>> S = {'abc ', 'defg','hi m'};>> cellplot(S)%示例8:将数字数组A按行或按列转换为元胞数组%A是4x3的数组>> A=[1 2 3;4 5 6;7 8 9;10 11 12];%把A的每一列转换为一个元胞,得到的C是1×3的元胞数组>> C=num2cell(A,1)C =[4x1 double] [4x1 double] [4x1 double] %把A的每一行转换为一个元胞,得到的C是4×1的元胞数组>> C=num2cell(A,2)C =[1x3 double][1x3 double][1x3 double][1x3 double]2.4结构数组%示例1:使用直接法来创建结构数组>> A(1).name = 'Pat';A(1).number = 176554;A(2).name = 'Tony';A(2).number = 901325;>> AA =1x2 struct array with fields:namenumber%示例2:利用struct函数来创建结构数组>> A(1)=struct('name','Pat','number',176554);A(2)=struct('name','Tony','number',901325);>> AA =1x2 struct array with fields:namenumber%示例3:使用deal函数来得到结构体中各结构域的值%定义结构数组A>> = 'Pat'; A.number = 176554;A(2).name = 'Tony';A(2).number = 901325;%得到结构数组中所有name结构域的数据>>[name1,name2] = deal(A(:).name)name1 =Patname2 =Tony%示例4:使用getfield函数来取得结构体中结构域的值%定义mystr结构数组>> mystr(1,1).name = 'alice';mystr(1,1).ID = 0;mystr(2,1).name = 'gertrude';mystr(2,1).ID = 1;%取得mystr(2,1)的结构域name的值>> f = getfield(mystr, {2,1}, 'name')f =gertrude%示例5:删除结构数组中的指定结构域%定义结构数组s>> s.field1=[1 2 3];s.field2='string';s.field3={1 2 3;4 5 6}; %删除结构域field1>> s=rmfield(s,'field1')s =field2: 'string'field3: {2x3 cell}%删除结构域'field2','field3'>> s=rmfield(s,{'field2','field3'})s =field1: [1 2 3]%示例6:%定义结构数组s>> s.field1=[1 2 3];s.field2='string';s.field3={1 2 3;4 5 6}; >> ss =field1: [1 2 3]field2: 'string'field3: {2x3 cell}%将结构数组转换为元胞数组>> c=struct2cell(s)c =[1x3 double]'string'{2x3 cell }%示例7:>>c = {'birch', 'betula', 65; 'maple', 'acer', 50}c ='birch''betula' [65]'maple''acer' [50]>>fields = {'name', 'genus', 'height'}; %fields包含struct中的结构域名>>s = cell2struct(c, fields, 2); %dim=2表示把c中的各行转换为struct数组s =2x1 struct array with fields:namegenusheight>> s(1)ans =name: 'birch'genus: 'betula'height: 65>> s(2)ans =name: 'maple'genus: 'acer'height: 50>> fields = {'field1', 'field2'};>> s = cell2struct(c, fields, 1); %dim=1表示把c中的各列转换为struct数组>> s(1)ans =field1: 'birch'field2: 'maple'>> s(2)ans =field1: 'betula'field2: 'acer'>> s(3)ans =field1:65field2:502.5矩阵操作%示例1: find函数的使用方法。

霍夫曼编码的MATLAB实现(完整版)

霍夫曼编码的MATLAB实现(完整版)

霍夫曼编码的MATLAB实现(完整版)%哈夫曼编码的 MATLAB 实现(基于 0、 1 编码):clc;clear;A=[0. 3, 0. 2, 0. 1, 0. 2, 0. 2]; 信源消息的概率序列A=fliplr(sort(A) ) ; %按降序排列T=A;[m, n]=size(A) ;B=zeros(n, n-1) ; %空的编码表(矩阵)for i=1: nB(i, 1) =T(i) ; %生成编码表的第一列endr=B(i, 1) +B(i-1, 1) ; %最后两个元素相加T(n-1) =r;T(n) =0;T=fliplr(sort(T) ) ;t=n-1;for j=2: n-1%生成编码表的其他各列for i=1: tB(i, j) =T(i) ;endK=find(T==r) ;B(n, j) =K(end) ; %从第二列开始,每列的最后一个元素记录特征元素在 %该列的位置r=(B(t-1, j) +B(t, j) ); %最后两个元素相加T(t-1) =r;T(t) =0;T=fliplr(sort(T) ) ;t=t-1;endB; %输出编码表END1=sym(' [0, 1] ' ) ; %给最后一列的元素编码END=END1;t=3;d=1;for j=n-2:-1:1%从倒数第二列开始依次对各列元素编码for i=1: t-2if i>1 & B(i, j) ==B(i-1, j)d=d+1;elsed=1;endB(B(n, j+1) , j+1) =-1;temp=B(: , j+1) ; x=find(temp==B(i, j) ) ; END(i) =END1(x(d) ) ; endy=B(n, j+1) ;END(t-1) =[char(END1(y) ) , ' 0' ]; END(t) =[char(END1(y) ) , ' 1' ]; t=t+1;END1=END;endA%排序后的原概率序列 END%编码结果for i=1: n[a, b]=size(char(END(i) ) ) ; L(i) =b;endavlen=sum(L. *A) %平均码长 H1=log2(A) ;H=-A*(H1' ) %熵P=H/avlen%编码效率。

matlab hht变换代码

matlab hht变换代码

matlab hht变换代码
在MATLAB中,可以使用以下代码实现希尔伯特-黄变换(HHT):
% 读取信号
Fs = 1000; % 采样频率
t = (0:1/Fs:length(x)-1/Fs); % 时间向量
x = sin(2*pi*50*t) + sin(2*pi*120*t); % 合成信号
% 希尔伯特-黄变换
[imf,t,A,f] = eemd(x,5); % 使用EEMD方法进行HHT变换
% 绘制原始信号和IMFs
figure;
subplot(2,1,1);
plot(t,x);
title('原始信号');
subplot(2,1,2);
for k = 1:length(imf)
plot(t,imf(k));
end
title('IMFs');
在上述代码中,x是一个合成信号,它由两个正弦波组成。

使用eemd函数对信号进行希尔伯特-黄变换,该函数使用经验模式分解(EMD)方法进行分解。

eemd函数的输出包括:
●imf:固有模式函数(IMFs)
●t:时间向量
●A:每个IMF的瞬时幅度
●f:每个IMF的瞬时频率
最后,使用subplot和plot函数绘制原始信号和IMFs。

程序1和2(HHT)

程序1和2(HHT)
function imf = emd(x)% Empiricial Mode Decomposition (Hilbert-Huang Transform)% imf = emd(x)% Func : findpeaksx= transpose(x(:));%转置为行矩阵imf = [];while ~ismonotonic(x) %当x不是单调函数,分解终止条件 x1 = x; sd = Inf;%均值%直到x1满足IMF条件,得c1 while (sd > 0.1) | ~isimf(x1) %当标准偏差系数sd大于0.1或x1不是固有模态函数时,分量终止条件 s1 = getspline(x1);%上包络线 s2 = -getspline(-x1);%下包络线 x2 = x1-(s1+s2)/2;%此处的x2为文章中的h sd = sum((x1-x2).^2)/sum(x1.^2); x1 = x2; end imf{end+1} = x1; x = x-x1;endimf{end+1} = x;% FUNCTIONSfunction u = ismonotonic(x)%u=0表示x不是单调函数,u=1表示x为单调的u1 = length(findpeaks(x))*length(findpeaks(-x));if u1 > 0, u = 0;else, u = 1; endfunction u = isimf(x)%u=0表示x不是固有模式函数,u=1表示x是固有模式函数N = length(x);u1 = sum(x(1:N-1).*x(2:N) < 0);u2 = length(findpeaks(x))+length(findpeaks(-x));if abs(u1-u2) > 1, u = 0;else, u = 1; endfunction s = getspline(x)%三次样条函数拟合成元数据包络线N = length(x);p = findpeaks(x);s = spline([0 p N+1],[0 x(p) 0],1:N);---------------------------------------------------------------------------------------------------------------------------------------------------------------function n = findpeaks(x)% Find peaks.找到极值 ,n为极值点所在位置% n = findpeaks(x)n = find(diff(diff(x) > 0) < 0);u = find(x(n+1) > x(n));n(u) = n(u)+1;----------------------------------------------------------------------------------------------------------------------------------------------------------------------------------function plot_hht00(x,Ts)% 双边带调幅信号的EMD分解% Plot the HHT.% plot_hht(x,Ts)%% :: Syntax% The array(列) x is the input signal and Ts is the sampling period(取样周期).% Example on use: [x,Fs] = wavread('Hum.wav');% plot_hht(x(1:6000),1/Fs);% Func : emd% Get HHT.clear all;close all;Ts=0.0005;t=0:Ts:10; % 采样率2000HZ% 调幅信号%x=sin(2*pi*t).*sin(40*pi*t);x=sin(2*pi*t);s1 = getspline(x);%上包络线s2 = -getspline(-x);%上包络线x1 = (s1+s2)/2;%此处的x2为文章中的hfigure;plot(t,x);xlabel('Time'), ylabel('Amplitude');title('双边带调幅信号');hold on;plot(t,s1,'-r');plot(t,s2,'-r');plot(t,x1,'g');imf = emd(x);for k = 1:length(imf) b(k) = sum(imf{k}.*imf{k}); th = angle(hilbert(imf{k})); d{k} = diff(th)/Ts/(2*pi);end[u,v] = sort(-b);b = 1-b/max(b);% Set time-frequency plots.N = length(x);c = linspace(0,(N-2)*Ts,N-1);%figure;for k = v(1:2) plot(c,d{k},'k.','Color',b([k k k]),'MarkerSize',3); hold on; set(gca,'FontSize',8,'XLim',[0 c(end)],'YLim',[0 50]);%设置x、y轴句柄 xlabel('Time'), ylabel('Frequency');title('原信号时频图');end% Set IMF plots.M = length(imf);N = length(x);c = linspace(0,(N-1)*Ts,N);for k1 = 0:4:M-1 figure for k2 = 1:min(4,M-k1), subplot(4,1,k2), plot(c,imf{k1+k2}); set(gca,'FontSize',8,'XLim',[0 c(end)]); title('EMD分解结果'); end xlabel('Time');end

EMD分解HHT变化matlab源代码

EMD分解HHT变化matlab源代码

EMD分解HHT变化matlab源代码第一篇:EMD分解HHT变化matlab源代码%function [] = UseEMD(x,t)function [] = UseEMD()%UNTITLED2 Summary of this function goes here %Detailed explanation goes hereN = 39;% # of data samplest = linspace(0,1,N);x=[20.3421.2520.6219.317.615.6815.4616.2717.9418.9718.7118.4719.1119.3118.6717.3216.2116.0916.0417.5419.320.1120.4221.3321.9422.1722.9723.0923.3923.9726.7129.3129.8230.3529.7827.9427.1727.8327.36];y=fft(x);plot(y,x);[imf,ort,nbits] = emd(x);emd_visu(x,t,imf,1);%¾ùÖµµÄƽ·½imfp2=mean(imf,2).^2%ƽ·½µÄ¾ùÖµimf2p=mean(imf.^2,2)%¸÷¸öIMFµÄ·½²îmse=imf2p-imfp2%·½²î°Ù·Ö±È£¬Ò²¾ÍÊÇ·½²î¹±Ï×ÂÊmseb=mse/sum(mse)*100%ÏÔʾ¸÷¸öIMFµÄ·½²îºÍ¹±Ï×ÂÊ[mse mseb]%¼ÆËãimfµÄÐÐÁÐάÊý[m,n] = size(imf)hx = hilbert(x)xf = real(hx);xi = imag(hx);%¼ÆËã˲ʱÕñ·ùA=sqrt(xf.^2 + xi.^2);figure(4)plot(t,A);title('˲ʱÕñ·ù')%¼ÆËã˲ʱÏàλp=angle(hx);figure(5);plot(t,p);title('˲ʱÏàλ')%¼ÆËã˲ʱƵÂÊ%dt = diff(t);%dx=diff(p);%sp = dx./dt;%figure(6);%plot(t(1:N-1),sp);title('˲ʱƵÂÊ') %imf1µÄhilbert±ä»»xn1 = hilbert(imf(1,:));xr1 = real(xn1);xi1 = imag(xn1);%imf1µÄ˲ʱÕñ·ùA1=sqrt(xr1.^2+xi1.^2);figure(7);subplot(2,1,1);plot(t,A1);xlabel('ʱ¼äª¡tª¢');ylabel('˲ʱÕñ·ù');title('imf1')%imf2µÄHilbert±ä»»xn2=hilbert(imf(2,:));xr2=real(xn2);xi2=imag(xn2);%imf2µÄ˲ʱÕñ·ùA2=sqrt(xr2.^2+xi2.^2); subplot(2,1,2);plot(t,A2);xlabel('ʱ¼äª¡tª¢');ylabel('˲ʱÕñ·ù');title('imf2')%imf1µÄ˲ʱÏàλP1=atan2(xi1,xr1);figure(8);subplot(2,1,1);plot(t,P1);xlabel('ʱ¼äª¡tª¢');ylabel('˲ʱÏàλ');title('imf1')%imf2µÄ˲ʱÏàλP2=atan2(xi2,xr2);subplot(2,1,2);plot(t,P2);xlabel('ʱ¼äª¡tª¢');ylabel('˲ʱÏàλ');title('imf2')%imf1˲ʱƵÂÊxh1=unwrap(P1);fs=1000;xhd1=fs*diff(xh1)/(2*pi);figure(9);subplot(2,1,1);plot(t(1:38),xhd1);xlabel('ʱ¼äª¡tª¢');ylabel('˲ʱƵÂÊ');title('imf1')%imf2˲ʱƵÂÊxh2=unwrap(P2);fs=1000;xhd2=fs*diff(xh2)/(2*pi);subplot(2,1,2);plot(t(1:38),xhd2);xlabel('ʱ¼äª¡tª¢');ylabel('˲ʱƵÂÊ');title('imf2')%¼ÆËãHHTµÄʱƵÆ×[A, fa, tt] = hhspectrum(imf);if size(imf,1)> 1[A,fa,tt] = hhspectrum(imf(1:end-1, :));else [A,fa,tt] = hhspectrum(imf);end[E,tt1]=toimage(A,fa,tt,length(tt));disp_hhs(E,[],fs)%¶þάͼÐÎÏÔʾHHTÊÓƵÆת£EÊÇÇóµÃµÄHHTÆ×for i=1:m-1faa=fa(i,:);[FA,TT1]=meshgrid(faa,tt1);%ÈýάͼÏÔʾHHTʱƵͼfigur e(11);surf(FA,TT1,E)title('HHTʱƵÆ×ÈýάÏÔʾ')end%»-±ß¼ÊÆ×E=flipud(E);for k=1:size(E,1)bjp(k)=sum(E(k,:))*1/fs;endf =(0:N-3)/N*(fs/2);figure(12)plot(f,bjp);end注意:matlab中需加载instfreq.m文件,从网上可下到第二篇:LU分解MatLab算法分析最近矩阵分析老师出了一道题目作为作业,是一道程序题,题目是对A矩阵做LU分解,要求能对A实现PA=LU的分解,并最终输出L,U,P矩阵。

FFT与HHT实例操作步骤1

FFT与HHT实例操作步骤1

FFT与HHT实例
1.打开Xstart,登陆到并行计算机。

2.在界面中输入pbsnodes,查看各节点运行情况。

3.第14节点loadsave=0,选择该节点,在界面中输入ssh –X node14 回车,输入密码:111111,再回车。

4.输入matlab ,运行软件。

5.下面构造一个信号,该信号由一个频率为5HZ 的正弦函数和一个频率随时间
变化的函数组成:210sin()sin(10)y t t ππ=+。

6.查看该函数的时域信号。

7.做FFT图。

先再本机上安装origin软件,然后打开。

8.双击matlab中我workspace的y,在array editor中选重整行右击copy。

9.在origin中选中B(y)右击粘贴。

10.选中B(y)列点击→分析→快速傅立叶变换。

11.调整坐标。

双击X轴刻度 刻度,把0.65改为0.5,应用。

12. 用除因子除0.001,确定。

13.双击曲线,颜色改为蓝色。

14.对y作EMD分解。

在matlab中输入imf=emd(y);
visuemd(imf);
15.HHT分解。

在matlab中输入:[A,f,tt]=hhspectrum(imf);
[im,tt]=toimage(A,f);
figure;
disp_hhs(im)。

地震波matlab hht变换代码

地震波matlab hht变换代码

地震波是指在地壳或地球内部产生的地震所携带的波动。

它是由地震破裂过程中的释放的能量引起的,具有复杂的波形和频谱特征。

地震波的研究对于地震学、地质学和工程地质学等领域有着重要的意义。

而Hilbert-Huang变换(HHT)则是一种能够有效处理非线性和非平稳信号的信号分析方法,它在地震波分析中具有广泛的应用。

在本文中,将介绍关于地震波以及HHT变换在Matlab中的代码实现,并对其进行深入探讨。

1. 地震波的特性地震波是地球内部或地表突然发生破裂时产生的振动波动,具有地球内部介质的特性。

地震波可分为P波、S波和面波等,它们在地震波传播过程中有着不同的特性表现。

地震波的振动特性呈现出复杂的波形和频谱结构,因此需要进行深入的分析和研究。

2. Hilbert-Huang变换(HHT)简介HHT是一种基于固有局部特征的自适应信号分析方法,由黄庭坚和黄淳在1996年提出。

它是一种多尺度、多分量和非线性的信号分析方法,适用于分析非线性和非平稳信号。

HHT主要包括经验模态分解(EMD)和 Hilbert 谱分析两部分,它可以将信号分解成若干本征模态函数(IMF),从而实现信号的时频特性分析。

HHT在地震波分析中可以有效地提取信号的时频特征,揭示信号内部的非线性和非平稳信息,具有很高的应用价值。

3. Matlab中的HHT变换代码实现在Matlab中,可以利用现有的工具箱或自行编写代码实现HHT变换。

需要对地震波信号进行预处理,包括去噪、滤波等操作,以保证HHT变换的准确性和稳定性。

接下来,可以利用Matlab内置的信号处理工具箱或自行编写代码实现HHT变换。

需要进行经验模态分解(EMD),将地震波信号分解成若干本征模态函数(IMF),然后对每个IMF进行Hilbert变换,得到相位和振幅信息,最终得到信号的时频分布。

通过Matlab的可视化工具,可以直观地展示地震波信号的时频特性,以及各个IMF的分布规律。

4. 个人观点和理解在地震波分析中,HHT变换具有独特的优势,能够有效地揭示地震波信号的非线性和非平稳特性。

(完整word版)用MATLAB实现计算器程序源代码(word文档良心出品)

(完整word版)用MATLAB实现计算器程序源代码(word文档良心出品)

function varargout = caculator(varargin)gui_Singleton = 1;gui_State = struct('gui_Name', mfilename, ...'gui_Singleton', gui_Singleton, ...'gui_OpeningFcn', @caculator_OpeningFcn, ...'gui_OutputFcn', @caculator_OutputFcn, ...'gui_LayoutFcn', [] , ...'gui_Callback', []);if nargin && ischar(varargin{1})gui_State.gui_Callback = str2func(varargin{1});endif nargout[varargout{1:nargout}] = gui_mainfcn(gui_State, varargin{:});elsegui_mainfcn(gui_State, varargin{:});endfunction caculator_OpeningFcn(hObject, eventdata, handles, varargin) handles.output = hObject;set(handles.edit1,'string','0');set(handles.edit5,'string','0');guidata(hObject, handles);function varargout = caculator_OutputFcn(hObject, eventdata, handles) varargout{1} = handles.output;function edit1_Callback(hObject, eventdata, handles)function edit1_CreateFcn(hObject, eventdata, handles)if ispc && isequal(get(hObject,'BackgroundColor'),get(0,'defaultUicontrolBackgroundColor'))set(hObject,'BackgroundColor','white');endfunction edit2_Callback(hObject, eventdata, handles)function edit2_CreateFcn(hObject, eventdata, handles)if ispc && isequal(get(hObject,'BackgroundColor'),get(0,'defaultUicontrolBackgroundColor'))set(hObject,'BackgroundColor','white');endfunction edit3_Callback(hObject, eventdata, handles)function edit3_CreateFcn(hObject, eventdata, handles)if ispc && isequal(get(hObject,'BackgroundColor'),get(0,'defaultUicontrolBackgroundColor'))set(hObject,'BackgroundColor','white');endfunction pushbutton1_Callback(hObject, eventdata, handles)handles.num1=strcat(get(handles.edit1,'string'),'+');set(handles.edit1,'string',handles.num1);guidata(hObject,handles);function pushbutton2_Callback(hObject, eventdata, handles) handles.num2=strcat(get(handles.edit1,'string'),'-');set(handles.edit1,'string',handles.num2);guidata(hObject,handles);function pushbutton3_Callback(hObject, eventdata, handles) handles.num3=strcat(get(handles.edit1,'string'),'*');set(handles.edit1,'string',handles.num3);guidata(hObject,handles);function pushbutton4_Callback(hObject, eventdata, handles) handles.num4=strcat(get(handles.edit1,'string'),'/');set(handles.edit1,'string',handles.num4);guidata(hObject,handles);function pushbutton_1_Callback(hObject, eventdata, handles) handles.shu1=get(hObject,'string');handles.yuanshu=get(handles.edit1,'string');if handles.yuanshu(end)=='N';handles.yuanshu(end)='';endif handles.yuanshu(1)=='0';handles.yuanshu=handles.yuanshu(2:end);endhandles.shu1=strcat(handles.yuanshu,handles.shu1);if length(handles.shu1)<2;elseif (length(handles.shu1)>=2)&&(handles.shu1(end-1)==')')&& (handles.shu1(1)=='l')temp=handles.shu1(end);handles.shu1(end)=handles.shu1(end-1);handles.shu1(end-1)=temp;endset(handles.edit1,'string',handles.shu1);guidata(hObject, handles);function pushbutton_2_Callback(hObject, eventdata, handles) handles.shu2=get(hObject,'string');handles.yuanshu=get(handles.edit1,'string');if handles.yuanshu(end)=='N';handles.yuanshu(end)='';endif handles.yuanshu(1)=='0';handles.yuanshu=handles.yuanshu(2:end);endhandles.shu2=strcat(handles.yuanshu,handles.shu2);if length(handles.shu2)<2;elseif (length(handles.shu2)>=2)&&(handles.shu2(end-1)==')')&& (handles.shu2(1)=='l')temp=handles.shu2(end);handles.shu2(end)=handles.shu2(end-1);handles.shu2(end-1)=temp;endset(handles.edit1,'string',handles.shu2);guidata(hObject, handles);function pushbutton_4_Callback(hObject, eventdata, handles) handles.shu4=get(hObject,'string');handles.yuanshu=get(handles.edit1,'string');if handles.yuanshu(end)=='N';handles.yuanshu(end)='';endif handles.yuanshu(1)=='0';handles.yuanshu=handles.yuanshu(2:end);endhandles.shu4=strcat(handles.yuanshu,handles.shu4);if length(handles.shu4)<2;elseif (length(handles.shu4)>=2)&&(handles.shu4(end-1)==')')&& (handles.shu4(1)=='l')temp=handles.shu4(end);handles.shu4(end)=handles.shu4(end-1);handles.shu4(end-1)=temp;endset(handles.edit1,'string',handles.shu4);guidata(hObject, handles);function pushbutton_3_Callback(hObject, eventdata, handles) handles.shu3=get(hObject,'string');handles.yuanshu=get(handles.edit1,'string');if handles.yuanshu(end)=='N';handles.yuanshu(end)='';endif handles.yuanshu(1)=='0';handles.yuanshu=handles.yuanshu(2:end);endhandles.shu3=strcat(handles.yuanshu,handles.shu3);if length(handles.shu3)<2;elseif (length(handles.shu3)>=2)&&(handles.shu3(end-1)==')')&& (handles.shu3(1)=='l')temp=handles.shu3(end);handles.shu3(end)=handles.shu3(end-1);handles.shu3(end-1)=temp;endset(handles.edit1,'string',handles.shu3);guidata(hObject, handles);function pushbutton_5_Callback(hObject, eventdata, handles) handles.shu5=get(hObject,'string');handles.yuanshu=get(handles.edit1,'string');if handles.yuanshu(end)=='N';handles.yuanshu(end)='';endif handles.yuanshu(1)=='0';handles.yuanshu=handles.yuanshu(2:end);endhandles.shu5=strcat(handles.yuanshu,handles.shu5);if length(handles.shu5)<2;elseif (length(handles.shu5)>=2)&&(handles.shu5(end-1)==')')&& (handles.shu5(1)=='l')temp=handles.shu5(end);handles.shu5(end)=handles.shu5(end-1);handles.shu5(end-1)=temp;endset(handles.edit1,'string',handles.shu5);guidata(hObject, handles);function pushbutton_9_Callback(hObject, eventdata, handles) handles.shu9=get(hObject,'string');handles.yuanshu=get(handles.edit1,'string');if handles.yuanshu(end)=='N';handles.yuanshu(end)='';endif handles.yuanshu(1)=='0';handles.yuanshu=handles.yuanshu(2:end);endhandles.shu9=strcat(handles.yuanshu,handles.shu9);if length(handles.shu9)<2;elseif (length(handles.shu9)>=2)&&(handles.shu9(end-1)==')')&& (handles.shu9(1)=='l')temp=handles.shu9(end);handles.shu9(end)=handles.shu9(end-1);handles.shu9(end-1)=temp;endset(handles.edit1,'string',handles.shu9);guidata(hObject, handles);function pushbutton_7_Callback(hObject, eventdata, handles) handles.shu7=get(hObject,'string');handles.yuanshu=get(handles.edit1,'string');if handles.yuanshu(end)=='N';handles.yuanshu(end)='';endif handles.yuanshu(1)=='0';handles.yuanshu=handles.yuanshu(2:end);endhandles.shu7=strcat(handles.yuanshu,handles.shu7);if length(handles.shu7)<2;elseif (length(handles.shu7)>=2)&&(handles.shu7(end-1)==')')&& (handles.shu7(1)=='l')temp=handles.shu7(end);handles.shu7(end)=handles.shu7(end-1);handles.shu7(end-1)=temp;endset(handles.edit1,'string',handles.shu7);guidata(hObject, handles);function pushbutton_8_Callback(hObject, eventdata, handles) handles.shu8=get(hObject,'string');handles.yuanshu=get(handles.edit1,'string');if handles.yuanshu(end)=='N';handles.yuanshu(end)='';endif handles.yuanshu(1)=='0';handles.yuanshu=handles.yuanshu(2:end);endhandles.shu8=strcat(handles.yuanshu,handles.shu8);if length(handles.shu8)<2;elseif (length(handles.shu8)>=2)&&(handles.shu8(end-1)==')')&& (handles.shu8(1)=='l')temp=handles.shu8(end);handles.shu8(end)=handles.shu8(end-1);handles.shu8(end-1)=temp;endset(handles.edit1,'string',handles.shu8);guidata(hObject, handles);function pushbutton_6_Callback(hObject, eventdata, handles) handles.shu6=get(hObject,'string');handles.yuanshu=get(handles.edit1,'string');if handles.yuanshu(end)=='N';handles.yuanshu(end)='';endif handles.yuanshu(1)=='0';handles.yuanshu=handles.yuanshu(2:end);endhandles.shu6=strcat(handles.yuanshu,handles.shu6);if length(handles.shu6)<2;elseif (length(handles.shu6)>=2)&&(handles.shu6(end-1)==')')&& (handles.shu6(1)=='l')temp=handles.shu6(end);handles.shu6(end)=handles.shu6(end-1);handles.shu6(end-1)=temp;endset(handles.edit1,'string',handles.shu6);guidata(hObject, handles);function pushbutton18_Callback(hObject, eventdata, handles) handles.jieguo=get(handles.edit1,'string');handles.jieguo=strcat('=',handles.jieguo);eval(['handles.result''1' handles.jieguo]);set(handles.edit5,'string',num2str(handles.result1));guidata(hObject,handles);function pushbutton_0_Callback(hObject, eventdata, handles) handles.shu0=get(hObject,'string');handles.yuanshu=get(handles.edit1,'string');if handles.yuanshu(end)=='N';handles.yuanshu(end)='';endif handles.yuanshu(1)=='0';handles.yuanshu=handles.yuanshu(2:end);endhandles.shu0=strcat(handles.yuanshu,handles.shu0);if length(handles.shu0)<2;elseif (length(handles.shu0)>=2)&&(handles.shu0(end-1)==')')&& (handles.shu0(1)=='l')temp=handles.shu0(end);handles.shu0(end)=handles.shu0(end-1);handles.shu0(end-1)=temp;endset(handles.edit1,'string',handles.shu0);guidata(hObject, handles);function pushbutton20_Callback(hObject, eventdata, handles) handles.shu10=get(hObject,'string');handles.yuanshu=get(handles.edit1,'string');handles.shu10=strcat(handles.yuanshu,handles.shu10);set(handles.edit1,'string',handles.shu10);guidata(hObject, handles);function pushbutton21_Callback(hObject, eventdata, handles) function edit5_Callback(hObject, eventdata, handles)function edit5_CreateFcn(hObject, eventdata, handles)if ispc && isequal(get(hObject,'BackgroundColor'),get(0,'defaultUicontrolBackgroundColor'))set(hObject,'BackgroundColor','white');endfunction pushbutton22_Callback(hObject, eventdata, handles) handles.yuanshu=get(handles.edit1,'string');handles.fanhao=strcat('-(',handles.yuanshu,')');set(handles.edit1,'string',handles.fanhao);guidata(hObject, handles);function pushbutton23_Callback(hObject, eventdata, handles)set(handles.edit1,'string','0');set(handles.edit5,'string','0');guidata(hObject, handles);function pushbutton24_Callback(hObject, eventdata, handles)result=questdlg('ÕæµÄÒªÍ˳ö£¿','Í˳öÈ·ÈÏ','È·¶¨','È¡Ïû','È¡Ïû'); if result=='È·¶¨', close(gcf); endfunction pushbutton25_Callback(hObject, eventdata, handles) handles.yuanshu=get(handles.edit1,'string');handles.sin=strcat('sin(',handles.yuanshu,')');set(handles.edit1,'string',handles.sin);guidata(hObject, handles);function pushbutton26_Callback(hObject, eventdata, handles) handles.yuanshu=get(handles.edit1,'string');handles.cos=strcat('cos(',handles.yuanshu,')');set(handles.edit1,'string',handles.cos);guidata(hObject, handles);function pushbutton27_Callback(hObject, eventdata, handles) handles.yuanshu=get(handles.edit1,'string');handles.tan=strcat('tan(',handles.yuanshu,')');set(handles.edit1,'string',handles.tan);guidata(hObject, handles);function pushbutton28_Callback(hObject, eventdata, handles) handles.yuanshu=get(handles.edit1,'string');handles.ln=strcat('reallog(',handles.yuanshu,')');set(handles.edit1,'string',handles.ln);guidata(hObject, handles);function pushbutton29_Callback(hObject, eventdata, handles) handles.yuanshu=get(handles.edit1,'string');handles.log=strcat('log',handles.yuanshu,'()');set(handles.edit1,'string',handles.log);guidata(hObject, handles);function pushbutton30_Callback(hObject, eventdata, handles) handles.yuanshu=get(handles.edit1,'string');handles.daoshu=strcat('1/(',handles.yuanshu,')');set(handles.edit1,'string',handles.daoshu);guidata(hObject, handles);function pushbutton31_Callback(hObject, eventdata, handles) handles.yuanshu=get(handles.edit1,'string');handles.sqrt=strcat('sqrt(',handles.yuanshu,')');set(handles.edit1,'string',handles.sqrt);guidata(hObject, handles);function pushbutton32_Callback(hObject, eventdata, handles) handles.yuanshu=get(handles.edit1,'string');handles.pingfang=strcat('(',handles.yuanshu,')^2');set(handles.edit1,'string',handles.pingfang);guidata(hObject, handles);function pushbutton33_Callback(hObject, eventdata, handles) handles.yuanshu=get(handles.edit1,'string');handles.ncifang=strcat('(',handles.yuanshu,')^N');set(handles.edit1,'string',handles.ncifang);guidata(hObject, handles);function pushbutton35_Callback(hObject, eventdata, handles) handles.yuanshu=get(handles.edit1,'string');handles.kuohao=strcat('(',handles.yuanshu,')');set(handles.edit1,'string',handles.kuohao);guidata(hObject, handles);function pushbutton36_Callback(hObject, eventdata, handles) handles.yuanshu=get(handles.edit1,'string');handles.yuanshu=handles.yuanshu(1:(end-1));if length(handles.yuanshu)<1handles.yuanshu='0';endset(handles.edit1,'string',handles.yuanshu);guidata(hObject, handles);% --- Executes on button press in pushbutton37.function pushbutton37_Callback(hObject, eventdata, handles) handles.yuanshu=get(handles.edit1,'string');handles.exp=strcat('exp(',handles.yuanshu,')');set(handles.edit1,'string',handles.exp);guidata(hObject, handles);% --- Executes when user attempts to close figure1.function figure1_CloseRequestFcn(hObject, eventdata, handles)% hObject handle to figure1 (see GCBO)% eventdata reserved - to be defined in a future version of MATLAB % handles structure with handles and user data (see GUIDATA)% Hint: delete(hObject) closes the figuredelete(hObject);%--------------------------------------------------------------------function Untitled_14_Callback(hObject, eventdata, handles)% hObject handle to Untitled_14 (see GCBO)% eventdata reserved - to be defined in a future version of MATLAB % handles structure with handles and user data (see GUIDATA)%--------------------------------------------------------------------function Untitled_15_Callback(hObject, eventdata, handles)set(gcf,'color','red')% hObject handle to Untitled_15 (see GCBO)% eventdata reserved - to be defined in a future version of MATLAB % handles structure with handles and user data (see GUIDATA)%--------------------------------------------------------------------function Untitled_16_Callback(hObject, eventdata, handles)set(gcf,'color','blue')% hObject handle to Untitled_16 (see GCBO)% eventdata reserved - to be defined in a future version of MATLAB % handles structure with handles and user data (see GUIDATA)%--------------------------------------------------------------------function Untitled_17_Callback(hObject, eventdata, handles)set(gcf,'color','green')% hObject handle to Untitled_17 (see GCBO)% eventdata reserved - to be defined in a future version of MATLAB % handles structure with handles and user data (see GUIDATA)%--------------------------------------------------------------------function Untitled_18_Callback(hObject, eventdata, handles)set(gcf,'color','black')% hObject handle to Untitled_18 (see GCBO)% eventdata reserved - to be defined in a future version of MATLAB % handles structure with handles and user data (see GUIDATA)%--------------------------------------------------------------------function Untitled_19_Callback(hObject, eventdata, handles)set(gcf,'color','yellow')% hObject handle to Untitled_19 (see GCBO)% eventdata reserved - to be defined in a future version of MATLAB % handles structure with handles and user data (see GUIDATA)%--------------------------------------------------------------------function Untitled_20_Callback(hObject, eventdata, handles)set(gcf,'color','m')% hObject handle to Untitled_20 (see GCBO)% eventdata reserved - to be defined in a future version of MATLAB % handles structure with handles and user data (see GUIDATA) function Untitled_1_Callback(hObject, eventdata, handles)function Untitled_2_Callback(hObject, eventdata, handles)function Untitled_3_Callback(hObject, eventdata, handles)function Untitled_4_Callback(hObject, eventdata, handles)function Untitled_4_CreateFcn(hObject, eventdata, handles)function Untitled_5_Callback(hObject, eventdata, handles)function Untitled_6_Callback(hObject, eventdata, handles)。

Matlab源程序代码

Matlab源程序代码

正弦波的源程序:(一),用到的函数1,f2t函数function x=f2t(X)global dt df t f T N%x=f2t(X)%x为时域的取样值矢量%X为x的傅氏变换%X与x长度相同并为2的整幂%本函数需要一个全局变量dt(时域取样间隔) X=[X(N/2+1:N),X(1:N/2)];x=ifft(X)/dt;end2,t2f函数。

function X=t2f(x)global dt df N t f T%X=t2f(x)%x为时域的取样值矢量%X为x的傅氏变换%X与x长度相同,并为2的整幂。

%本函数需要一个全局变量dt(时域取样间隔) H=fft(x);X=[H(N/2+1:N),H(1:N/2)]*dt;end(二),主程序。

1,%(1)绘出正弦信号波形及频谱global dt df t f Nclose allk=input('取样点数=2^k, k取10摆布');if isempty(k), k=10; endf0=input('f0=取1(kz)摆布');if isempty(f0), f0=1; endN=2^k;dt=0.01; %msdf=1/(N*dt); %KHzT=N*dt; %截短期Bs=N*df/2; %系统带宽f=[-Bs+df/2:df:Bs]; %频域横坐标t=[-T/2+dt/2:dt:T/2]; %时域横坐标s=sin(2*pi*f0*t); %输入的正弦信号S=t2f(s); %S是s的傅氏变换a=f2t(S); %a是S的傅氏反变换a=real(a);as=abs(S);subplot(2,1,1) %输出的频谱plot(f,as,'b');gridaxis([-2*f0,+2*f0,min(as),max(as)])xlabel('f (KHz)')ylabel('|S(f)| (V/KHz)') %figure(2)subplot(2,1,2)plot(t,a,'black') %输出信号波形画图gridaxis([-2/f0,+2/f0,-1.5,1.5])xlabel('t(ms)')ylabel('a(t)(V)')gtext('频谱图')最佳基带系统的源程序:(一),用到的函数f2t函数和t2f函数。

图论matlab程序源码

图论matlab程序源码
while m==0 %判断m的值,当m为0时进行循环执行,m不为0时终止
k=k+1;%k进行加1,k=2
B1=[B Y1(3:length(Y1))]; %取Y1中第三个到第六个元素,与B组合在一起赋给B1,B1=[6 5 6 8 9]
f=length(B1);%求B1的中元素的个数 f=5
B=Y(1)+Y(2);%取Y中第一个和第二个元素相加求和,Y(1)=2,Y(2)=4,B=6
W=[Y(1) Y(2) B];%给W赋值为[2 4 6],Y(1)=2,Y(2)=4,B=6, W=[2 4 6]
Y1=Y;%把Y的值赋给Y1,Y1=[2 4 5 6 7 8 9]
m=0;%对m进行初始化为0
W(j,k)=1; %给边的终点赋值为1W(2,1)=1
k=k+1;%对k进行加1
end
end
end
elseif f==1 %关联矩阵转换为邻接矩阵
m=size(F,2);
k=1;%初始化k=1;代表列的值
for i=1:n %i从1到4
for j=i:n %j从i到4 对F进行遍历
if F(i,j)~=0%判断F中(i,j)元素是否为0,不为0进行赋值,为0不进行操作,继续循环
W(i,k)=1; %给边的始点赋值为1 W(1,1)=1,
C(i+1)=t;%对C(2)赋值为2 C=[1 2 0]
W(C(1:(i+1)),C(1:(i+1)))=0;%对W([1 2],[1 2])赋值为0,W=[0 0 1;0 0 1;1 1 0]
end;
m=sum(sum(F))/2; %计算图的边数 m=6,sum(F)表示矩阵F中每一列求和,点度的总和是边数的2倍

Matlab实现HHT程序(源码)

Matlab实现HHT程序(源码)

clear all;x=load ('065140106、TXT');fs=;N=length(x);t=0:1/fs:(N-1)/fs;z=x;c=emd(z);%计算每个IMF分量及最后一个剩余分量residual与原始信号得相关性[m,n]=size(c);for i=1:m;a=corrcoef(c(i,:),z);xg(i)=a(1,2);endxg;for i=1:m-1%--------------------------------------------------------------------%计算各IMF得方差贡献率%定义:方差为平方得均值减去均值得平方%均值得平方%imfp2=mean(c(i,:),2)、^2%平方得均值%imf2p=mean(c(i,:)、^2,2)%各个IMF得方差mse(i)=mean(c(i,:)、^2,2)-mean(c(i,:),2)、^2;end;mmse=sum(mse);for i=1:m-1mse(i)=mean(c(i,:)、^2,2)-mean(c(i,:),2)、^2;%方差百分比,也就就是方差贡献率mseb(i)=mse(i)/mmse*100;%显示各个IMF得方差与贡献率end;%画出每个IMF分量及最后一个剩余分量residual得图形figure(1)for i=1:m-1disp(['imf',int2str(i)]) ;disp([mse(i) mseb(i)]);end;subplot(m+1,1,1)plot(t,z)set(gca,'fontname','times New Roman')set(gca,'fontsize',14、0)ylabel(['signal','Amplitude'])for i=1:m-1subplot(m+1,1,i+1);set(gcf,'color','w')plot(t,c(i,:),'k')set(gca,'fontname','times New Roman')set(gca,'fontsize',14、0)ylabel(['imf',int2str(i)])endsubplot(m+1,1,m+1);set(gcf,'color','w')plot(t,c(m,:),'k')set(gca,'fontname','times New Roman')set(gca,'fontsize',14、0)ylabel(['r',int2str(m-1)])%画出每个IMF分量及剩余分量residual得幅频曲线figure(2)subplot(m+1,1,1)set(gcf,'color','w')[f,z]=fft(t,z);plot(f,z,'k')set(gca,'fontname','times New Roman')set(gca,'fontsize',14、0)ylabel(['initial signal',int2str(m-1),'Amplitude'])for i=1:m-1subplot(m+1,1,i+1);set(gcf,'color','w')[f,z]=fft(t,c(i,:));plot(f,z,'k')set(gca,'fontname','times New Roman')set(gca,'fontsize',14、0)ylabel(['imf',int2str(i),'Amplitude'])endsubplot(m+1,1,m+1);set(gcf,'color','w')[f,z]=fft(t,c(m,:));plot(f,z,'k')set(gca,'fontname','times New Roman')set(gca,'fontsize',14、0)ylabel(['r',int2str(m-1),'Amplitude'])hx=hilbert(z);xr=real(hx);xi=imag(hx);%计算瞬时振幅sz=sqrt(xr、^2+xi、^2);%计算瞬时相位sx=angle(hx);%计算瞬时频率dt=diff(t);dx=diff(sx);sp=dx、/dt;figure(6)plot(t(1:N-1),sp)title('瞬时频率')%计算HHT时频谱与边际谱[A,fa,tt]=hhspectrum(c);[E,tt1]=toimage(A,fa,tt,length(tt));figure(3)disp_hhs(E,tt1) %二维图显示HHT时频谱,E就是求得得HHT谱pausefigure(4)for i=1:size(c,1)faa=fa(i,:);[FA,TT1]=meshgrid(faa,tt1);%三维图显示HHT时频图surf(FA,TT1,E)title('HHT时频谱三维显示')hold onendhold offE=flipud(E);for k=1:size(E,1)bjp(k)=sum(E(k,:))*1/fs;endf=(1:N-2)/N*(fs/2);figure(5)plot(f,bjp);xlabel('频率/ Hz');ylabel('信号幅值');title('信号边际谱')%要求边际谱必须先对信号进行EMD分解function [A,f,tt] = hhspectrum(x,t,l,aff)error(nargchk(1,4,nargin));if nargin < 2t=1:size(x,2);endif nargin < 3l=1;endif nargin < 4aff = 0;endif min(size(x)) == 1if size(x,2) == 1if nargin < 2t = 1:size(x,2);endendNmodes = 1;elseNmodes = size(x,1);endlt=length(t);tt=t((l+1):(lt-l));for i=1:Nmodesan(i,:)=hilbert(x(i,:)')';f(i,:)=instfreq(an(i,:)',tt,l)';A=abs(an(:,l+1:end-l));if affdisprog(i,Nmodes,max(Nmodes,100))endendfunction disp_hhs(im,t,inf)% DISP_HHS(im,t,inf)% displays in a new figure the spectrum contained in matrix "im" % (amplitudes in log)、%% inputs : - im : image matrix (e、g、, output of "toimage")% - t (optional) : time instants (e、g、, output of "toimage")% - inf (optional) : -dynamic range in dB (wrt max)% default : inf = -20%% utilisation : disp_hhs(im) ; disp_hhs(im,t) ; disp_hhs(im,inf) % disp_hhs(im,t,inf)figurecolormap(bone)colormap(1-colormap);if nargin==1inf=-20;t = 1:size(im,2);endif nargin == 2if length(t) == 1inf = t;t = 1:size(im,2);elseendendif inf >= 0error('inf doit etre < 0')endM=max(max(im));im = log10(im/M+1e-300);inf=inf/10;imagesc(t,fliplr((1:size(im,1))/(2*size(im,1))),im,[inf,0]);set(gca,'YDir','normal')xlabel(['time'])ylabel(['normalized frequency'])title('Hilbert-Huang spectrum')function [f,z]=fftfenxi(t,y)L=length(t);N=2^nextpow2(L);%fft默认计算得信号就是从0开始得t=linspace(t(1),t(L),N);deta=t(2)-t(1);m=0:N-1;f=1、/(N*deta)*m;%下面计算得Y就就是x(t)得傅里叶变换数值%Y=exp(i*4*pi*f)、*fft(y)%将计算出来得频谱乘以exp(i*4*pi*f)得到频移后[-2,2]之间得频谱值Y=fft(y);z=sqrt(Y、*conj(Y));。

(完整word版)matlab经典代码大全

(完整word版)matlab经典代码大全

(完整word版)matlab经典代码大全哈哈哈MATLAB显示正炫余炫图:plot(x,y1,'* r',x,y2,'o b')定义【0,2π】;t=0:pi/10:2*pi;定义函数文件:function [返回变量列表]=函数名(输入变量列表)顺序结构:选择结构1)if-else-end语句其格式为:if 逻辑表达式程序模块1;else程序模块2;End图片读取:%选择图片路径[filename, pathname] = ...uigetfile({'*.jpg';'*.bmp';'*.gif'},'选择图片');%合成路径+文件名str=[pathname,filename];%为什么pathname和filename要前面出现的位置相反才能运行呢%读取图片im=imread(str);%使用图片axes(handles.axes1);%显示图片imshow(im);边缘检测:global imstr=get(hObject,'string');axes (handles.axes1);switch strcase ' 原图'imshow(im);case 'sobel'BW = edge(rgb2gray(im),'sobel');imshow(BW);case 'prewitt'BW = edge(rgb2gray(im),'prewitt');imshow(BW);case 'canny'BW = edge(rgb2gray(im),'canny');imshow(BW);Canny算子边缘定位精确性和抗噪声能力效果较好,是一个折中方案end;开闭运算:se=[1,1,1;1,1,1;1,1,1;1,1,1]; %Structuring ElementI=rgb2gray(im);imshow(I,[]);title('Original Image');I=double(I);[im_height,im_width]=size(I);[se_height,se_width]=size(se);halfheight=floor(se_height/2);halfwidth=floor(se_width/2);[se_origin]=floor((size(se)+1)/2);image_dilation=padarray(I,se_origin,0,'both'); %Image to be used for dilationimage_erosion=padarray(I,se_origin,256,'both'); %Image to be used for erosion %%%%%%%%%%%%%%%%%% %%% Dilation %%%%%%%%%%%%%%%%%%%%%for k=se_origin(1)+1:im_height+se_origin(1)for kk=se_origin(2)+1:im_width+se_origin(2)dilated_image(k-se_origin(1),kk-se_origin(2))=max(max(se+image_dilation(k-se_origin(1):k+halfh eight-1,kk-se_origin(2):kk+halfwidth-1)));endendfigure;imshow(dilated_image,[]);title('Image after Dilation'); %%%%%%%%%%%%%%%%%%%% Erosion %%%%%%%%%%%%%%%%%%%%se=se';for k=se_origin(2)+1:im_height+se_origin(2)for kk=se_origin(1)+1:im_width+se_origin(1)eroded_image(k-se_origin(2),kk-se_origin(1))=min(min(image_erosion(k-se_origin(2):k+halfwidth -1,kk-se_origin(1):kk+halfheight-1)-se));endendfigure;imshow(eroded_image,[]);title('Image after Erosion'); %%%%%%%%%%%%%%%%%%%%%%%%%%%%% %%%%%%%%%%%%%%%%%%%%% Opening(Erosion first, then Dilation) %%% %%%%%%%%%%%%%%%%%%%%%%%%%% %%%%%%%%%%%%%%%%%%%%%se=se';image_dilation2=eroded_image; %Image to be used for dilationfor k=se_origin(1)+1:im_height-se_origin(1)for kk=se_origin(2)+1:im_width-se_origin(2)opening_image(k-se_origin(1),kk-se_origin(2))=max(max(se+image_dilation2(k-se_origin(1):k+hal fheight-1,kk-se_origin(2):kk+halfwidth-1)));endendfigure;imshow(opening_image,[]);title('OpeningImage'); %%%%%%%%%%%%%%%%%%%%%%%%%%%%%% %%%%%%%%%%%%%%%%%%%% Closing(Dilation first, then Erosion) %%% %%%%%%%%%%%%%%%%%%%%%%%%%% %%%%%%%%%%%%%%%%%%%%%se=se';image_erosion2=dilated_image; %Image to be used for erosionfor k=se_origin(2)+1:im_height-se_origin(2)for kk=se_origin(1)+1:im_width-se_origin(1)closing_image(k-se_origin(2),kk-se_origin(1))=min(min(image_erosion2(k-se_origin(2):k+halfwidt h-1,kk-se_origin(1):kk+halfheight-1)-se));endendfigure;imshow(closing_image,[]);title('Closing Image');Warning: Image is too big to fit on screen; displaying at 31% scale.> In truesize>Resize1 at 308In truesize at 44In imshow at 161图像的直方图归一化:I=imread(‘red.bmp’);%读入图像figure;%打开新窗口[M,N]=size(I);%计算图像大小[counts,x]=imhist(I,32);%计算有32个小区间的灰度直方图counts=counts/M/N;%计算归一化灰度直方图各区间的值stem(x,counts);%绘制归一化直方图图像平移:I=imread('shuichi.jpg');se=translate(strel(1),[180 190]);B=imdilate(I,se);figure;subplot(1,2,1),subimage(I);title('原图像');subplot(1,2,2),subimage(B);title('平移后图像');图像的转置;A=imread('nir.bmp');tform=maketform('affine',[0 1 0;1 0 0;0 0 1]);B=imtransform(A,tform,'nearest');figure;imshow(A);figure;imshow(B);imwrite(B,'nir转置后图像.bmp');图像滤波:B = imfilter(A,H,option1,option2,...)或写作g = imfilter(f, w, filtering_mode, boundary_options, size_options)其中,f为输入图像,w为滤波掩模,g为滤波后图像。

Hilbert-HuangTransform:matlab希尔伯特-黄变换:matlab实现

Hilbert-HuangTransform:matlab希尔伯特-黄变换:matlab实现

Hilbert-HuangTransform:matlab希尔伯特-黄变换:matlab实现关于Hilbert-Huang的matlab实现,材料汇总,⽐较杂...感谢所有⽹络上的贡献者们:)核⼼:以下代码计算HHT边际谱及其对应频率⼯具包要求:G-Rilling EMD Toolbox,TFTB Toolbox附:黄锷先⽣课题组开发的⼯具包(可以在找到),这⾥并未⽤到。

% Empirical mode decomposition, resulting in intrinc mode functions.% Without parameter 'MAXMODES' the process may be seriously delayed by% decompose original signals into too many IMFs (not necessary, 9 is% enough generally)imfs = emd(oriSig, 'MAXMODES', 9);% HHT spectrum: hhtS[A, f, t] = hhspectrum(imfs);[hhtS, ~, fCent] = toimage(A, f, t);% Marginal hilbert spectrum: hhtMS, xf: correspondig frequencyfor k = 1 : size(hhtS, 1)hhtMS(k) = sum(hhtS(k, : )) * 1 / fs;endxf = fCent(1, :) .* fs;简单来说,设置好路径之后输⼊install_emd即可。

hhspectrum 函数说明(8楼:⽼⽼的学⽣)% [A,f,tt] = HHSPECTRUM(imf,t,l,aff)% Input:%- imf : matrix with one IMF per row % 将emd分解得到的IMF代⼊就可以,就是你的程序中写的c变量,不⽤加最后⼀⾏的趋势项%- t : time instants % 瞬时时间或持续时间??(写[1:信号长度]就可以,真实的时间可以根据采样率转换%- l : estimation parameter for instfreq % 瞬时频率的估计参数??(写1就可以,决定瞬时频率估计时的边界从第⼏个点开始%- aff : if 1, displays the computation evolution % 显⽰计算进程选项,不想显⽰写0就可以%% Output:% - A : amplitudes of IMF rows% - f : instantaneous frequencies% - tt : truncated time instants % 截⽌时间??(截断时间,返回的是瞬时频率对应的时间,要⽐原来信号的时间按短,由前⾯的l值决定)%% calls:% - hilbert : computes the analytic signal% - instfreq : computes the instantaneous frequency % 瞬时频率%% Example:[A,f,tt] = hhspectrum(c(1:end-1,:),[1:n],1,0);[im,tt,ff] = toimage(A,f,tt,512);disp_hhs(im,tt,[],fs);9楼:⽼⽼的学⽣关于时频图的概念,我认为是与诸如⼩波,gabor等联合时频分析⽅法联系在⼀起的。

希尔伯特-黄变换说明及程序

希尔伯特-黄变换说明及程序

质模态函数(IMF)任何一个资料,满足下列两个条件即可称作本质模态函数。

⒈局部极大值(local maxima)以及局部极小值(local minima)的数目之和必须与零交越点(zero crossing)的数目相等或是最多只能差1,也就是说一个极值后面必需马上接一个零交越点。

⒉在任何时间点,局部最大值所定义的上包络线(upper envelope)与局部极小值所定义的下包络线,取平均要接近为零。

因此,一个函数若属于IMF,代表其波形局部对称于零平均值。

此类函数类似于弦波(sinusoid-like),但是这些类似于弦波的部分其周期与振幅可以不是固定。

因为,可以直接使用希尔伯特转换,求得有意义的瞬时频率。

经验模态分解(EMD)EMD算法流程图建立IMF是为了满足希尔伯特转换对于瞬时频率的限制条件之前置处理,也是一种转换的过程。

我们将IMF来做希尔伯特转换可以得到较良好的特性,不幸的是大部分的资料并不是IMF,而是由许多弦波所合成的一个组合。

如此一来,希尔伯特转换并不能得到正确的瞬时频率,我们便无法准确的分析资料。

为了解决非线性(non-linear)与非稳态(non-stationary)资料在分解成IMF时所遇到的困难,便发展出EMD。

经验模态分解是将讯号分解成IMF的组合。

经验模态分解是借着不断重复的筛选程序来逐步找出IMF。

以讯号为例,筛选程序的流程概述如下:步骤 1 : 找出中的所有局部极大值以及局部极小值,接着利用三次样条(cubic spline),分别将局部极大值串连成上包络线与局部极小值串连成下包络线。

步骤 2 : 求出上下包络线之平均,得到均值包络线。

步骤 3 : 原始信号与均值包络线相减,得到第一个分量。

步骤 4 : 检查是否符合IMF的条件。

如果不符合,则回到步骤1并且将当作原始讯号,进行第二次的筛选。

亦即重复筛选次直到符合IMF的条件,即得到第一个IMF分量,亦即步骤 5 : 原始讯号减去可得到剩余量,表示如下式步骤 6 : 将当作新的资料,重新执行步骤1至步骤5,得到新的剩余量。

遗传算法经典MATLAB代码

遗传算法经典MATLAB代码

遗传算法经典学习Matlab代码遗传算法实例:也是自己找来的,原代码有少许错误,本人都已更正了,调试运行都通过了的。

对于初学者,尤其是还没有编程经验的非常有用的一个文件遗传算法实例% 下面举例说明遗传算法%% 求下列函数的最大值%% f(x)=10*sin(5x)+7*cos(4x) x∈[0,10]%% 将x 的值用一个10位的二值形式表示为二值问题,一个10位的二值数提供的分辨率是每为(10-0)/(2^10-1)≈0.01。

%% 将变量域[0,10] 离散化为二值域[0,1023], x=0+10*b/1023, 其中 b 是[0,1023] 中的一个二值数。

%% %%--------------------------------------------------------------------------------------------------------------%%--------------------------------------------------------------------------------------------------------------%% 编程%-----------------------------------------------% 2.1初始化(编码)% initpop.m函数的功能是实现群体的初始化,popsize表示群体的大小,chromlength表示染色体的长度(二值数的长度),% 长度大小取决于变量的二进制编码的长度(在本例中取10位)。

%遗传算法子程序%Name: initpop.m%初始化function pop=initpop(popsize,chromlength)pop=round(rand(popsize,chromlength)); % rand随机产生每个单元为{0,1} 行数为popsize,列数为chromlength的矩阵,% roud对矩阵的每个单元进行圆整。

matlab源代码

matlab源代码

例错误!文档中没有指定样式的文字。

-1%周期信号(方波)的展开,fb_jinshi.mclose all;clear all;N=100; %取展开式的项数为2N+1项T=1;fs=1/T;N_sample=128; %为了画出波形,设置每个周期的采样点数dt = T/N_sample;t=0:dt:10*T-dt;n=-N:N;Fn = sinc(n/2).*exp(-j*n*pi/2);Fn(N+1)=0;ft = zeros(1,length(t));for m=-N:Nft = ft + Fn(m+N+1)*exp(j*2*pi*m*fs*t);endplot(t,ft)例错误!文档中没有指定样式的文字。

-4利用FFT计算信号的频谱并与信号的真实频谱的抽样比较。

脚本文件T2F.m定义了函数T2F,计算信号的傅立叶变换。

function [f,sf]= T2F(t,st)%This is a function using the FFT function to calculate a signal's Fourier %Translation%Input is the time and the signal vectors,the length of time must greater %than 2%Output is the frequency and the signal spectrumdt = t(2)-t(1);T=t(end);df = 1/T;N = length(st);f=-N/2*df:df:N/2*df-df;sf = fft(st);sf = T/N*fftshift(sf);脚本文件F2T.m定义了函数F2T,计算信号的反傅立叶变换。

function [t st]=F2T(f,sf)%This function calculate the time signal using ifft function for the input %signal's spectrumdf = f(2)-f(1);Fmx = ( f(end)-f(1) +df);dt = 1/Fmx;N = length(sf);T = dt*N;%t=-T/2:dt:T/2-dt;t = 0:dt:T-dt;sff = fftshift(sf);st = Fmx*ifft(sff);另写脚本文件fb_spec.m如下:%方波的傅氏变换, fb_spec.mclear all;close all;T=1;N_sample = 128;dt=T/N_sample;t=0:dt:T-dt;st=[ones(1,N_sample/2), -ones(1,N_sample/2)]; %方波一个周期subplot(211);plot(t,st);axis([0 1 -2 2]);xlabel('t'); ylabel('s(t)');subplot(212);[f sf]=T2F(t,st); %方波频谱plot(f,abs(sf)); hold on;axis([-10 10 0 1]);xlabel('f');ylabel('|S(f)|');%根据傅氏变换计算得到的信号频谱相应位置的抽样值sff= T^2*j*pi*f*0.5.*exp(-j*2*pi*f*T).*sinc(f*T*0.5).*sinc(f*T*0.5);plot(f,abs(sff),'r-')例错误!文档中没有指定样式的文字。

  1. 1、下载文档前请自行甄别文档内容的完整性,平台不提供额外的编辑、内容补充、找答案等附加服务。
  2. 2、"仅部分预览"的文档,不可在线预览部分如存在完整性等问题,可反馈申请退款(可完整预览的文档不适用该条件!)。
  3. 3、如文档侵犯您的权益,请联系客服反馈,我们会尽快为您处理(人工客服工作时间:9:00-18:30)。

Mat lab实现HHT程序(源码-非常珍贵)clear all;x=load (");fs=1000000;N=length(x);t=O:l/fs:(N-l)/fs;z=x;c=emd(z);%计•算每个IMF分量及最后一个剩余分量residual与原始信号的相关性[m,n]=size(c);for i=l:m;a=corrcoef(c(i,:),z);xg(i)=a(l,2);endxg;for i=l:m-l%计算各IMF的方差贡献率%定义:方差为平方的均值减去均值的平方%均值的平方%imfp2=mean(c(i,:),2).A2%平方的均值%imf2p=mean(c(i,:).A2/2)%各个IMF的方差mse(i)=mean(c(i/:).A2/2)-mean(c(i,:)z2).A2;end;mmse=sum(mse);for i=l:m-lmse(i)=mea n(c(i/:).A2/2)-mean(c(i/:)z2).A2;%方差百分比,也就是方差贡献率mseb(i)=mse(i)/mmse*100;%显示各个IMF的方差和贡献率end;%画出每个IMF分量及最后一个剩余分量residual的图形figure(l)for i=l:m-ldisp(['imf'/int2str(i)]) ;disp([mse(i) mseb(i)]);end;subplot(m+lU)plot(t’z)set(gca/fontname7times New RomarV)set(gca/f on tsize1,ylabelfl'signar/Amplitude'])for i=l:m-lsubplot(m+l z l,i+l);setfgcf/colorVw1)plot(t,c(i,:),k)set(gca/fontname1,'times New Roman‘)set(gca/f on tsize1,ylabel(['imf'/int2str(i)])endsubplot(m+l z l,m+l);setlgcf/color'/w')plot(tjC(m,:),'k‘)set(gca/fontname1,'times New Roman‘)set(gca/f on tsize1,ylabel(['r\int2str(m-l)])%画出每个IMF分量及剩余分量residual的幅频曲线figure(2)subplot(m+blj)setfgcf/color'/w1)[f,z]=fft(t,z);plot(f,z;k')set(gca/fontname1,'times New Roman‘)set(gca/f on tsize1,ylabel『initial signal',int2str(m-l),'Amplitude'])for i=l:m-lsubplot(m+l z l,i+l);settgcf/color'/w1)[f,z]=fft(t/c(i/:));plot(f,z,'k')set(gca/fontname7times New RomarV)set(gca/f on tsize1,ylabel(['imf'/int2str(i)/'Amplitude,])endsubplot(m+l,l,m+l);settgcf/color'/w1)[f/z]=fft(t,c(m/:));plot(f,z;k')set(gca/fontname7times New Roman‘)set(gca/f on tsize1,ylabel([,r,/int2str(m-l);Amplitude,])hx=hilbert(z);xr=real(hx);xi=imag(hx);%计算瞬时振幅sz=sqrt(xr.A2+xi.A2);%计算瞬时相位sx=a ng le(hx);%计算瞬时频率dt=diff(t);dx=diff(sx);sp=dx./dt;figure(6)plot(t(l:N-l)/sp)title/瞬时频率J%讣算HHT时频谱和边际谱[A,fa,tt]=hhspectrum(c);[E/ttl]=toimage(A,fa,tt,length(tt));figure(3)disp_hhs(E,ttl) %二维图显示HHT时频谱,E是求得的HHT谱pausefigure(4)for i=l:size(c,l)faa=fa(i,:);[FA/TTl]=meshgrid(faa,ttl);%三维图显示HHT 时频图surf(FA,TTl,E)title('HHT时频谱三维显示]hold onendhold offE=flipud(E);for k=l:size(E,l)bjp(k)=sum(E(k/))*l/fs;endf=(l:N-2)/N*(fs/2);figure(5)plot(f,bjp);xlabel('频率/Hz');ylabelf信号幅值,);title/信号边际谱J%要求边际谱必须先对信号进行EMD分解function [A,f,tt] = hhspectrum(x,t,I,aff)error(nargchk(l,4, nargin));讦nargin v 2t=l:size(x,2);end讦nargin v 31=1;end讦nargin < 4aff = 0;endif min(size(x)) == 1if size(x#2) == 1X = X ;讦nargin < 2t = l:size(x,2);endendNmodes = 1;elseNmodes = size(x,l);endlt=le ngth(t);tt=t((l+i):(lt-l));for i=l:Nmodesan(i/:)=hilbert(x(i/:)')';f(i/:)=instfreq(an(i/:),/tt/l)1;A=abs(an(:,l+l:end-l));if affdisprog(i,Nmodes,max(Nmodes/100))endendfunction disp_hhs(im,t,inf)%DISP_HHS(im/t,inf)% displays in a new figure the spectrum contained in matrix "im" % (amplitudes in log).%% inputs :・ im : image matrix output of "toimage")% -1 (optional): time instants output of '^oimage")% ・ inf (optional) : -dynamic range in dB (wrt max)% default: inf = -20%% utilisation : disp_hhs(im); disp_hhs(im,t); disp_hhs(im,inf)% disp_hhs(im,t,inf)figurecolormap(bone)colormap(l-colormap);讦nargiinf=-20;t = l:size(im,2);end讦nargin == 2if length(t) == 1inf = t;t = l:size(im,2);elseinf = -20;endendif inf >= 0errorfinf doit etre < 01)endM=max(max(im));im = logl0(im/M+le-300);inf=inf/10;imagesc(t,fliplr((l:size(im/l))/(2*size(im/l)))/im,[inf/0]);set(gca/YDir,/,normal1)xlabel(['time'])ylabel(「normalized frequency1])titlef^ilbert-Huang spectrum')function [f,z]=fftfenxi(t,y)L=length(t);N=2A nextpow2(L);%fft默认计算的信号是从0开始的t=linspace(t(l)/t(L),N);deta=t(2)-t(l);m=0:N-l;f=l./(N*deta)*m;%下面计算的Y就是x(t)的傅里叶变换数值%Y=exp(i*4*pi*f).*fft(y)%将计算出来的频谱乘以exp(i*4*pi*f)得到频移后卜2,2]之间的频谱值Y=fft(y);z=sqrt(Y.*conj(Y));。

相关文档
最新文档