使用伪随机二进制序列信号发生器块确定频率响应
要构建系统的频率特性，需要获得模型的频率响应：将测试信号应用于系统的输入并测量其输出。 为了减少实验的时间，而不损失结果的质量，可以发送伪随机二进制序列（PRBS）信号而不是谐波信号。 PRBS是一个周期信号,只能取两个特定的值.
在此示例中，在model_analysis_with_prbs中。engee模型，使用特殊块发送PRBS信号并记录系统响应。 在此脚本中，处理数据并构建LFC。
Pkg.add("ControlSystems")
Pkg.add("DSP")
In [ ]:
using ControlSystems
using DSP
using EngeeDSP
算法的采样周期
In [ ]:
Ts = 0.001 # @param {type:"number"}
;
将以对数标度收集频率响应的频率阵列。
wmin-建立LFCH的最小频率
wmax是建立的最大频率
LFCH n_ponts是图形上的点数。
In [ ]:
wmin = 0.01 # @param {type:"number"}
wmax = 100 # @param {type:"number"}
n_points = 20 # @param {type:"number"}
w = collect(logrange(wmin, wmax, length=n_points));
运行模拟
In [ ]:
modelName = "model_analysis_with_prbs" # @param {type:"string"}
;
In [ ]:
if modelName ∉ getfield.(engee.get_all_models(), :name)
engee.load("$(@__DIR__)/$(modelName).engee");
end
engee.run(modelName);
将仿真数据写入变量
In [ ]:
model_data = collect(PRBSSigGenData);
t = model_data[:,1];
data = model_data[:,2]
nd = size(data, 1)
nw = length(w)
ready = zeros(nd)
perturbation = zeros(nd)
input = zeros(nd)
output = zeros(nd)
for i = 1:nd
ready[i] = data[i][1]
perturbation[i] = data[i][2]
input[i] = data[i][3]
output[i] = data[i][4]
end
我们在计算系统的频率响应之前对信号进行处理
In [ ]:
ts = t[2] - t[1]
dupidx = findall(diff(t) .< ts/2)
deleteat!(ready, dupidx)
deleteat!(perturbation, dupidx)
deleteat!(input, dupidx)
idx = findall(ready .== 1.)
d = perturbation[idx]
u = input[idx]
y = output[idx]
nfft = sum(idx .> 0.)
if !iszero(rem(nfft, 2))
freq = (2*pi/ts) .* collect(range(0., 1., nfft + 1))
lastel = Int64((nfft + 1) / 2)
freq = freq[1:lastel]
else
freq = (pi/ts) .* collect(range(0., 1., Int64(nfft/2) + 1))
end
pd = 2*pi / (nfft - 1)
filterwindow = [0.5 .* ones(nfft-1) .- 0.5 .* cos.(pd .* collect(0:nfft-2)); 0.]
fd = EngeeDSP.Functions.fft(d .* filterwindow);
fu = EngeeDSP.Functions.fft(u .* filterwindow);
fy = EngeeDSP.Functions.fft(y .* filterwindow);
d2u = fu ./ fd;
d2y = fy ./ fd;
u2y = d2y ./ d2u;
构造频率响应
我们得到系统的频率响应
请注意，文件夹必须包含具有函数的文件 interpfreqresp.jl
In [ ]:
# 将文件与服务功能连接
include("$(@__DIR__)/interpfreqresp.jl")
plant_response = interpfreqresp(freq, u2y[1:length(freq)], w);
我们根据频率响应计算系统的频率响应
In [ ]:
magnitude_estimated = 20 * log10.(abs.(plant_response));
phase_estimated = rad2deg.(DSP.unwrap((atan.(imag.(plant_response), real.(plant_response)))));
我们正在建立一个低频响应系统
In [ ]:
gr()
p1 = plot(w, magnitude_estimated, linealpha = 0.0, markershape = :circle, markersize = 4, ylabel = "L，dB", legend = :none)
p2 = plot(w, phase_estimated, linealpha = 0.0, markershape = :circle, markersize = 4, label = "模型", xlabel = "哦，我很高兴", ylabel = ",冰雹", legend = :bottomleft)
plot(p1, p2, layout = (2,1), grid = false, framestyle = :box, xscale=:log10)
Out[0]:
.png)
我们计算系统传递函数的频率响应。
让我们将实验LCF与传递函数构建的特性进行比较。
In [ ]:
# 管理对象
L = 0.006 # 电感Gn
k = 0.10 # 反电动势系数（单位为△s）/rad
R = 0.25 # 欧姆电阻
J = 0.015 # 转子的转动惯量kg〇m^2
K = 0.10 # 扭矩常数n∶m/A
tfnum = [K];
tfden = [L*J, R*J, K*k];
sys = tf(tfnum, tfden)
W = sys
Out[0]:
In [ ]:
W = c2d(sys, Ts) # 我们用等于算法采样周期的量化周期离散化传递函数。
magnitude_ideal, phase_ideal, w_ideal = bode(W);
magnitude_ideal = reshape(magnitude_ideal, size(magnitude_ideal, 3),)
phase_ideal = reshape(phase_ideal, size(phase_ideal, 3),);
LFCH的比较
In [ ]:
p1 = plot(w, magnitude_estimated, linealpha = 0.0, markershape = :circle, markersize = 4, ylabel = "L，dB", legend = :none)
plot!(w_ideal, 20 * log10.(magnitude_ideal), xlims = (wmin, wmax))
p2 = plot(w, phase_estimated, linealpha = 0.0, markershape = :circle, markersize = 4, label = "模型", xlabel = "哦，我很高兴", ylabel = ",冰雹", legend = :bottomleft)
plot!(w_ideal, phase_ideal, xlims = (wmin, wmax), label = "传递函数")
plot(p1, p2, layout = (2,1), grid = false, framestyle = :box, xscale=:log10)
Out[0]:
.png)