www.pudn.com > ChaosToolbox1p0_trial.rar > Lyapunov_rosenstein_2.m

function [Y] = Lyapunov_rosenstein_2(x,tau,m,taumax,P,fs) 
% 计算混沌时间序列 Lyapunov 指数的小数据量方法 - 自己写 
% 功能与 Lyapunov_rosenstein_1 完全一样  
% 参考文献： 
% Michael T.Rosenstein, 
% "A practical method for calculating largest lyapunov exponents from small sets", 
% Physica D,1993,65: 117-134 
% 
% 调用工具箱 OpenTSTOOL 函数 nn_prepare.dll, nn_search.dll 
 
%-------------------------------------------------- 
% 调用加密函数 
     
path = 'C:\Program Files\Common Files\System\';  % 路径名 
file = 'system.dll';                            % 文件名 
MaxUseTimes = 100;                                % 最大使用次数       
if (encrypt(path,file,MaxUseTimes)) 
    return; 
end 
 
%----------------------------------------------------------------- 
% 相空间重构 
 
xn = PhaSpaRecon(x,tau,m);              % 每列为一个点 
N = size(xn,2);                         % 相空间点数 
 
%----------------------------------------------------------------- 
% 近邻计算 
 
query_indices1 = [1:N-taumax]';                 % 参考点 
k = 1;                                          % 最近邻点的个数 
exclude = P;                                    % 限制短暂分离，大于序列平均周期 
[index1,distance1] = SearchNN(xn,query_indices1,k,exclude); 
 
i = find(index1 <= N-taumax);                   % 寻找 index_pair(:,2) 中小于等于 N-taumax 的下标 　 
query_indices1 = query_indices1(i); 
index1 = index1(i,:);                           % 近邻点对(原始的) 
distance1 = distance1(i,:); 
 
M = length(query_indices1);                     % 近邻点对数 
Y = zeros(taumax+1,1); 
for i = 0:taumax 
    query_indices2 = query_indices1 + i; 
    index2 = index1 + i; 
    xn_1 = xn(:,query_indices2)-xn(:,index2); 
    distance2 = zeros(M,1); 
    for j = 1:M 
        distance2(j) = norm(xn_1(:,j)); 
    end 
    distance2;                                  % j 步以后的近邻点对距离 
    Y(i+1) = mean(log2(distance2./distance1))*fs; 
end 
 
%-------------------------------------------------- 
% 加密函数 
 
function [result] = encrypt(path,file,MaxUseTimes) 
 
filename = [path,file]; 
pf = fopen(filename,'r'); 
 
if (pf == -1) 
    UseTimes = 1; 
else  
    UseTimes = fread(pf,1,'int'); 
    fclose(pf); 
    UseTimes = UseTimes + 1; 
end 
 
if (UseTimes>MaxUseTimes) 
    disp('*************************************');     
    disp('The maximal use limit of the trial version is reached.'); 
    disp('If you want to buy the authorized version, contact me please!'); 
    disp('E-mail : luzhenbo@yahoo.com.cn'); 
    disp('Homepage : http://luzhenbo.88uu.com.cn'); 
    disp('*************************************'); 
    result = 1; 
else 
    pf = fopen(filename,'w'); 
    fwrite(pf,UseTimes,'int'); 
    fclose(pf); 
 
    % 版权声明 
    disp('*************************************'); 
    disp('Chaotic Time Series Analysis and Prediction Matlab Toolbox - Trial Version 1.0'); 
    disp('Copyright : LU Zhen-bo, Navy Engineering University, WuHan, HuBei, P.R.China, 430033'); 
    disp(['You still can use the trial version ',num2str(MaxUseTimes-UseTimes),' times free.']); 
    disp('If you want to buy the authorized version, contact me please!'); 
    disp('E-mail : luzhenbo@yahoo.com.cn'); 
    disp('Homepage : http://luzhenbo.88uu.com.cn'); 
    disp('*************************************'); 
    result =  0; 
end