Частотная передаточная функция (FRD)

Частотная передаточная функция (частотный отклик) задается с помощью наборов частот и комплексных чисел, выражающих реакцию системы на этих частотах.
Это удобный способ хранить данные частотной характеристики системы для визуализации или использования в расчетах.

neededLibs = [
    "EngeeControlSystems", # Для использования функции frd
    "ControlSystemsBase",  # Для построения ЛАФЧХ
    "CSV"                  # Для загрузки экспериментальных данных из таблицы
]

for lib in neededLibs
    try
        eval(Meta.parse("using $lib"))
    catch ex
        Pkg.add(lib)
        eval(Meta.parse("using $lib"))
    end
end

Загрузка данных из таблицы

df = CSV.read("frd_data.csv", DataFrame)
freq = df.freq
resp = [parse(ComplexF64, r) for r in df.resp]; # Преобразование из String в ComplexF64

Создание Модели частотной характеристики

sys = frd(resp, freq)

frd(ComplexF64[0.23841984895937737 - 0.0019508505988724586im;;; 0.25400219150798586 - 0.0035114829943812854im;;; 0.24295100787892412 + 0.0018065135835565785im;;; … ;;; -7.792719980301814e-5 + 0.0008045520700888559im;;; -0.0001896945323066172 + 0.0006147494316645394im;;; 1.116473000890296e-5 + 0.0006156621012700552im], [0.1, 0.10274594854461801, 0.10556729942333293, 0.10846612314544048, 0.11144454707535631, 0.11450475699382817, 0.11764899870201859, 0.12087957966963436, 0.12419887072831304, 0.12760930781150914  …  78.36418966217519, 80.51602998770538, 82.72695874133697, 84.99859846090146, 87.33261623828437, 89.73072494285637, 92.19468447849992, 94.72630307515246, 97.32743861581504, 100.0], nothing, "seconds", "rad/timeunit")

Построение ЛАФЧХ

magnitude, phase, _ = bode(sys, freq)

([0.2384278301622347;;; 0.25402646280196745;;; 0.24295772414293115;;; … ;;; 0.0008083171914251253;;; 0.0006433512876484359;;; 0.0006157633263977576], [-0.46880749163629404;;; -0.7920417471113614;;; 0.42602703045454404;;; … ;;; -264.4677100629009;;; -252.85123290043435;;; -271.038916961609], [0.1, 0.10274594854461801, 0.10556729942333293, 0.10846612314544048, 0.11144454707535631, 0.11450475699382817, 0.11764899870201859, 0.12087957966963436, 0.12419887072831304, 0.12760930781150914  …  78.36418966217519, 80.51602998770538, 82.72695874133697, 84.99859846090146, 87.33261623828437, 89.73072494285637, 92.19468447849992, 94.72630307515246, 97.32743861581504, 100.0])

magnitude = reshape(magnitude, size(magnitude, 3),)
phase = reshape(phase, size(phase, 3),)

256-element Vector{Float64}:
   -0.46880749163629404
   -0.7920417471113614
    0.42602703045454404
   -1.0400851545567138
   -1.757285611061378
   -0.4572532835944419
    2.0688061063360683
   -0.5731712493801628
    0.22512626521159854
   -0.9097604824609797
   -3.0282154711474245
   -0.0793722602111739
    2.327323942243346
    ⋮
 -254.55952891164728
 -259.66789284233107
 -253.52620714314904
 -265.78529244104374
 -256.87843361975564
 -271.8782200891723
 -262.35584773156495
 -256.0261319236719
 -270.1132862026374
 -264.4677100629009
 -252.85123290043435
 -271.038916961609

p1 = Plots.plot(freq, 20*log10.(magnitude), ylabel = "L, дБ", legend = :none)
p2 = Plots.plot(freq, phase, xlabel = "ω, рад/с", ylabel = "Ψ, град", legend = :none)
Plots.plot(p1, p2, layout = (2,1), grid = false, framestyle = :box, xscale=:log10, xlims = (0.1,100))