python画psd图

Question

已知输出电压和输入加速度的Laplace变换之前的比值，以及count和加速度之间的关系，怎么构造response的gain和sensitivity

Answer 1

上午好☀️☀️☀️️
本答案参考通义千问

在Python中绘制PSD（Power Spectral Density）图时，我们需要根据系统的输入和输出之间的传递函数（即Laplace变换的比值），来计算系统的响应增益（Gain）和灵敏度（Sensitivity）。以下是一个完整的流程说明和代码示例。

---

一、问题理解

你提到“已知输出电压和输入加速度的Laplace变换之前的比值为上图”，这表示我们有系统的传递函数：

$$ H(s) = \frac{V_{out}(s)}{A_{in}(s)} $$

其中：

• $ V_{out}(s) $ 是输出电压的Laplace变换
• $ A_{in}(s) $ 是输入加速度的Laplace变换

我们要基于这个传递函数构造 response的gain 和 sensitivity。

---

二、关键概念解释

1. Gain（增益）

系统对特定频率输入的放大能力，通常在频域中表示为：

$$ |H(j\omega)|^2 $$

它代表了系统在不同频率下对输入信号的放大程度。

2. Sensitivity（灵敏度）

灵敏度表示系统对输入扰动或参数变化的敏感程度。在控制系统中，灵敏度函数通常定义为：

$$ S(\omega) = \frac{1}{1 + H(j\omega)G(j\omega)} $$

但如果你只是需要计算系统的响应对输入的敏感性（比如对加速度输入的响应），可以简单地使用：

$$ S(\omega) = |H(j\omega)| $$

---

三、解决方案步骤

1. 导入必要的库

import numpy as np
import matplotlib.pyplot as plt
from scipy.signal import TransferFunction, lsim2

1. 定义传递函数 $ H(s) $

假设你的传递函数是如下形式（例如一个简单的二阶系统）：

$$ H(s) = \frac{1}{s^2 + 2s + 1} $$

你可以将其表示为分子和分母多项式系数：

num = [1]         # 分子系数
den = [1, 2, 1]   # 分母系数
H = TransferFunction(num, den)

1. 计算频率响应

w, mag = H.bode(n=1000)

1. 计算 Gain 和 Sensitivity

gain = mag ** 2     # 增益（平方）
sensitivity = mag   # 灵敏度（幅度）

1. 绘制PSD图（功率谱密度）

注意：PSD通常是对信号的功率随频率分布的描述。如果你没有实际数据，可以模拟一个加速度输入信号，并通过系统得到输出电压，然后进行PSD分析。

# 模拟输入加速度信号（例如白噪声）
t = np.linspace(0, 10, 1000)
A_in = np.random.normal(size=len(t))  # 白噪声作为输入

# 仿真系统响应
t_out, y_out, _ = lsim2(H, A_in, t)

# 计算PSD
frequencies, psd = plt.psd(y_out, Fs=100, sides='twosided')

1. 绘图展示

plt.figure(figsize=(12, 6))

plt.subplot(2, 1, 1)
plt.semilogx(w, gain)
plt.title('Gain (Magnitude Squared)')
plt.xlabel('Frequency [rad/s]')
plt.ylabel('Gain')

plt.subplot(2, 1, 2)
plt.semilogx(w, sensitivity)
plt.title('Sensitivity (Magnitude)')
plt.xlabel('Frequency [rad/s]')
plt.ylabel('Sensitivity')

plt.tight_layout()
plt.show()

---

四、总结

重点部分：

• Gain 是系统对输入信号的放大能力，计算方式为 $ |H(j\omega)|^2 $
• Sensitivity 表示系统对输入的响应敏感程度，计算方式为 $ |H(j\omega)| $
• 使用 scipy.signal.TransferFunction 可以方便地定义和分析传递函数
• 使用 plt.psd() 或 scipy.signal.welch() 可以计算并绘制PSD图

---

五、修改后的代码（完整）

import numpy as np
import matplotlib.pyplot as plt
from scipy.signal import TransferFunction, lsim2, welch

# 定义传递函数 H(s) = 1/(s^2 + 2s + 1)
num = [1]
den = [1, 2, 1]
H = TransferFunction(num, den)

# 计算频率响应
w, mag = H.bode(n=1000)

# 计算 Gain 和 Sensitivity
gain = mag ** 2
sensitivity = mag

# 模拟输入加速度信号（白噪声）
t = np.linspace(0, 10, 1000)
A_in = np.random.normal(size=len(t))

# 仿真系统响应
t_out, y_out, _ = lsim2(H, A_in, t)

# 计算PSD
frequencies, psd = welch(y_out, fs=100, nperseg=256)

# 绘制图形
plt.figure(figsize=(12, 6))

plt.subplot(2, 1, 1)
plt.semilogx(w, gain)
plt.title('Gain (Magnitude Squared)')
plt.xlabel('Frequency [rad/s]')
plt.ylabel('Gain')

plt.subplot(2, 1, 2)
plt.semilogx(w, sensitivity)
plt.title('Sensitivity (Magnitude)')
plt.xlabel('Frequency [rad/s]')
plt.ylabel('Sensitivity')

plt.tight_layout()
plt.show()

---

如需进一步定制传递函数或输入信号，请提供具体的系统模型或数据格式。