Дифференциация от FFT нахождения экстремумов

Question

Я пытаюсь найти нули функции. См. Мой код ниже.

Потому что 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); 

Answer 1

Вы можете продолжить вдоль массива t и искать точки, где значения меняющих знак. Это указывает на наличие нуля.

Actaully, функция fzero MATLAB использует аналогичный метод. Вы сказали, что вы не использовали, потому что вам требуется массив, а не анонимная функции, но вы могли бы преобразовать массив в анонимную функцию с помощью простой линейной интерполяции так:

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

EDIT: Просто, чтобы уточнить, что я имею в виду. Если у вас есть массив 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'); 

нули я получаю правильно: