0

我试图找到一个函数的零点。请参阅下面的代码。

因为fft需要一个数值数组，所以我没有定义要使用的符号函数fzero。

但是，这种方法并不准确，并且依赖于step. 你有更好的主意吗？

step=2000;
x=0:pi/step:2*pi;
y= 4+5*cos(10*x)+20*cos(40*x)+cos(100*x);
fy = fft(y');
fy(1:8) =0;
fy(12:end-10)=0; 
fy(end-6:end)=0;
ffy = ifft(fy);
t=diff(ffy);
x=0:pi/step:2*pi-pi/step;
plot(x,t)
indices= find(t<5e-4 & t>-5e-4);
4

1 回答 1

2

您可以沿着数组继续t并查找值更改符号的点。这将表明零的存在。

实际上，MATLAB 的fzero函数使用了类似的方法。您说您没有使用它是因为您需要一个数组，而不是匿名函数，但是您可以使用简单的线性插值将数组转换为匿名函数，如下所示：

func = @(k) interp1(x,t,k);  % value of t(x) interpolated at x=k
fzero(func,initial_value);


编辑：只是为了澄清我的意思。如果你有一个数组t并且你想找到它的零点......

f = 5;                            % frequency of wave in Hz
x = 0:0.01:1;                     % time index
t = cos( 2*pi*f*x );              % cosine wave of frequency f

zeroList = [];                    % start with an empty list of zeros
for i = 2:length(t)               % loop through the array t

    current_t = t(i);             % current value of t
    previous_t = t(i-1);          % previous value of t

    if current_t == 0                  % the case where the zero is exact
        newZero = x(i);
        zeroList = [zeroList,newZero];
    elseif current_t*previous_t < 0    % a and b have opposite sign if a*b is -ve
        % do a linear interpolation to find the zero (solve y=mx+b)
        slope = (current_t-previous_t)/(x(i)-x(i-1));
        newZero = x(i) - current_t/slope;
        zeroList = [zeroList,newZero];
    end

end

figure(1); hold on;
axis([ min(x) max(x) -(max(abs(t))) (max(abs(t))) ]); 
plot(x,t,'b');
plot(x,zeros(length(x)),'k-.');
scatter(zeroList,zeros(size(zeroList)),'ro');

我得到的零是正确的：

在此处输入图像描述

于 2014-12-11T19:52:41.893 回答