clear
clc

fid = fopen('10.prn', 'r');
fgetl(fid);  % Пропускаем первую строку
y = fscanf(fid, '%g', [3 Inf]);
fclose(fid);

y = y(3, :) - mean(y(3, :));

N = size(y, 2);

S = abs(fftshift(fft(y .* parzenwin(N)')));
S = S / max(S);
S = S(fix(N/2+1):N);

Thr = 0.25;            % Порог
NB = 5;                % Количество соседних спектральных составляющих

ind1 = find(S > Thr);  % Всё, что выше порога, считаем сигналом

a = zeros(size(S));    % Добавляем сюда соседние отсчёты - NB штук с каждой стороны
a(ind1) = ones(size(ind1));
b = zeros(size(S));
for i=-NB:NB
    b = b | circshift(a, [0 i]);
end

ind1 = find(b == 1);   % Отсчёты, соответствующие сигналу
ind2 = find(b == 0);   % Отсчёты, соответствующие шуму

if size(ind1, 2) > 6*NB
    fprintf('Результат сомнителен\n');
end

sinad = 10*log10(sum(S(ind1).^2) / sum(S(ind2).^2)) + 2.7;

N = size(y, 2);
f = 125/2*((fix(N/2):N-1)-N/2)/(N/2);
plot(f, 20*log10(S), 'r')
hold on
plot(f(ind2), 20*log10(S(ind2)), 'b-')
grid on
hold off
set(gcf, 'color', 'white');
xlabel('f, MHz');
ylabel('S, dB');

fprintf('SINAD: %g dB\n', sinad);