[GMavr 0]

Вроде бы, теоретически, комплексный коэффициент передачи K(w) = jwG(w), где w - круговая частота, G(w) - преобразование Фурье от переходной функции. Так?

Но вот рассмотрим в Matlab-е переходную функцию RC цепочки - y(t). tau=RC.

fs = 10;
N = 1024;
tau = 10;
t = 0:1/fs:(N-1)/fs;
f = 0:fs/N:fs*(N-1)/N;
y= 1-exp(-t/tau);

Коэффициент передачи всем известен: Kr = 1./(1+1i*2*pi*f*tau);

Как я ни пытался прикладывать ДПФ к массиву y, никак нужная зависимость не получается. Лучшая попытка:

Y = fft(y);

K=2i*pi*(0:N-1)/N.*Y+ y(N-1);

Здесь предполагалось, что до прихода ступеньки на выходе 0, а последний элемент массива y(N-1) уже вышел на постоянное значение 1.

В начале abs(K) и abs(Kr) практически одинаковые, но довольно быстро расходятся. Почему?

Как правильно получить из y массив K с отсчетами комплексного коэффициента передачи?

[Alex-lab 4]

По моему, заменить преобразование Лапласа на Фурье нельзя, замена p на jw делается уже после преобразования для нахождения афчх.

[V_G 11]

ЧХ - преобразование Фурье не от переходной, а от импульсной характеристики. Импульсная - производная от переходной.

[GMavr 0]

 

10 hours ago, Alex-lab said:

По моему, заменить преобразование Лапласа на Фурье нельзя, замена p на jw делается уже после преобразования для нахождения афчх.

Можно подробнее? Формула из книжки и там именно Фурье. Лаплас вроде здесь не нужен.

 

8 hours ago, V_G said:

ЧХ - преобразование Фурье не от переходной, а от импульсной характеристики. Импульсная - производная от переходной.

Так вроде в формуле K(w) = jwG(w) умножение на jw это учитывает? Нет?

Изменено 29 декабря, 2019 пользователем GMavr

[Grizzly 0]

@GMavr

Для аналогового случая у вас переходная характеристика g(t) = 1 - exp(-t/(RC)) -> импульсная характеристика h(t) = 1/(RC)* exp(-t/(RC)), t >= 0. Применяя преобразование Фурье (свойство затухающей экспоненты), получаем K(f) = 1/(1+j*2pi*f*RC).

В дискретной системе: https://www.dummies.com/education/science/science-engineering/working-across-domains-example-take-the-rc-low-pass-filter-to-the-z-domain/

В примере АЧХ аналогового и дискретного фильтров совпадают в соответствующей полосе.

Поскольку у вас коэффициенты импульсной характеристики в дискретной форме известны (2 штуки в числителе, знаменателя нет), то немного изменив код, сможете получить вашу АЧХ с помощью fft: https://ccrma.stanford.edu/~jos/filters/Practical_Frequency_Response_Analysis.html

[GMavr 0]

@Grizzly

Спасибо за ответ и за ссылки. Но вопрос состоял немного не в этом.

Есть массив значений переходной функции системы. Полученный экспериментально. Нужно по этим значениям оценить АЧХ системы. С экспериментальными данными ничего не выходит. Попробовал заполнить массив значениями переходной функции RC цепочки, для которой все хорошо известно. И тоже ничего не вышло.

Можно, конечно, по известным входным и выходным данным построить фильтр (Виннеровское оценивание). И потом вычислить АЧХ этого фильтра. Но нет ли пути попроще? Да и результат такого подхода меня как-то не впечатлил, хотя этот способ более или менее рабочий.

6 hours ago, Grizzly said:

Поскольку у вас коэффициенты импульсной характеристики в дискретной форме известны (2 штуки в числителе, знаменателя нет)

Вроде 1 коэффициент в числителе и 2 в знаменателе для цифрового аналога RC фильтра.

[Grizzly 0]

1 час назад, GMavr сказал:

Вроде 1 коэффициент в числителе и 2 в знаменателе для цифрового аналога RC фильтра.

Прошу прощения. Конечно же, у вас БИХ фильтр.

По виду переходной характеристики вы можете оценить tau = RC, а затем рассчитать частотную характеристику. Или у вас произвольная система?

Как вариант - найти производную от переходной характеристики (diff), а затем посчитать FFT. Пример на стр. 14-16: https://ethz.ch/content/dam/ethz/special-interest/mavt/dynamic-systems-n-control/idsc-dam/Lectures/Signals-and-Systems/Lectures/Fall2018/Lecture11_sigsys.pdf

Но хорошие результаты всё равно не получатся в таком случае.

[thermit 1]

Импульсная х-ка аналогового фильтра - реакция на воздействие дельта-функции дирака. Импульсная х-ка дискретной системы - реакция на воздействие дискретной дельта-функции. Разница в воздействиях, в принципе самих систем порождают когнитивный диссонанс. Ну и умножение на jw этот диссонанс только усугубляет. Не проще почитать про дискретизацию аналоговых фильтров? Многое прояснится.

[V_G 11]

10 часов назад, GMavr сказал:

Есть массив значений переходной функции системы. Полученный экспериментально. Нужно по этим значениям оценить АЧХ системы.

Импульсная характеристика получается из переходной дифференцированием, т.е. подчеркивает (усиливает) верхние частоты. Если есть проблемы с адекватностью АЧХ, вычисленной по экспериментальной ПХ , значит, есть проблемы с точностью дифференцирования и/или с частотой Найквиста.

Попробуйте брать точки на переходной характеристике почаще, эти проблемы постепенно уйдут за точность вычислений

[GMavr 0]

Спасибо.

On 12/30/2019 at 1:10 AM, Grizzly said:

Как вариант - найти производную от переходной характеристики (diff), а затем посчитать FFT.

Спасибо. Пока нет шума, полученный данным способом график АЧХ почти совпадает с тем, что должно быть в идеале. Пробовал также численное интегрирование. Результат чуть лучше, но не принципиально. А вот добавление даже небольшого шума заметно портит картину.

fs = 10;
N = 1024;
tau = 10;
t = 0:1/fs:(N-1)/fs;
f = 0:fs/N:fs*(N-1)/N;
y = 1-exp(-t/tau); % Реакция RC цепи на ступеньку (переходная функция)
yn = y + (rand(size(t))-0.5)*0.01; % + шум
% Комплексный коэффициент передачи, вычисленный по переходной функции
K = fft(diff(y));
% Тоже + шум
Kn =  fft(diff(yn));
% Такой коэффициент передачи должен быть
Kr = 1./(1+1i*2*pi*f*tau);
plot(f(1:N/2), abs(K(1:N/2)), 'b',f(1:N/2), abs(Kn(1:N/2)), 'g', f(1:N/2), abs(Kr(1:N/2)), 'r');
legend('K', 'Kn', 'Kr');

Тут, наверно, ничего сделать нельзя. Придется ограничится максимальной частотой 1/20 от частоты дискретизации.

[Grizzly 0]

14 часов назад, GMavr сказал:

А вот добавление даже небольшого шума заметно портит картину.

Изначально у вас про шум ничего не говорилось. Это уже совершенно иная задача. Необходимо использовать различные методы идентификации систем. В принципе вы можете попытаться как-то сгладить полученную нынешнюю оценку.

Как вам посоветовал @V_G, остаётся только увеличивать частоту дискретизации, чтобы получать больше точек (чаще) переходной характеристики.

[Tanya 4]

18 часов назад, GMavr сказал:

Спасибо.

Спасибо. Пока нет шума, полученный данным способом график АЧХ почти совпадает с тем, что должно быть в идеале. Пробовал также численное интегрирование. Результат чуть лучше, но не принципиально. А вот добавление даже небольшого шума заметно портит картину.

 

Тут, наверно, ничего сделать нельзя. Придется ограничится максимальной частотой 1/20 от частоты дискретизации.

Преобразование одной функции в другую относится только к линейным системам. Добавление шума убивает линейность.

[Tano 0]

On 1/6/2020 at 6:14 PM, Tanya said:

Преобразование одной функции в другую относится только к линейным системам. Добавление шума убивает линейность.

to TC.
Если шум аддитивный, (а у Вас "yn = y + (rand(size(t))...."), то все с ним боряться усреднением/накоплением. 
Запустите Вашу функцию: "yn = y + (rand(size(t))-0.5)*0.01; % + шум"  в цикле и будет Вам точность...
Кстати, а почему у Вас  "rand" , а не "randn" ?                                                                            
to Tanya: "Добавление шума убивает линейность." -- ну если ОН ну оооочень большой...