ARTICLE DETAIL

资讯详情

深耕网站视觉设计与运营推广的一线实战洞察。

用Python手绘波特图:从传递函数到电路稳定性分析

用Python手绘波特图:从传递函数到电路稳定性分析 搞硬件这几年波特图一直是我评估电路稳定性的首选工具。以前在学校用Matlab画工作后想换Python又总觉得麻烦直到疫情期间宅家决定彻底把这套流程跑通。说实话真动手做完才发现用Python手绘波特图的难度被严重高估了但这套逻辑对理解电路稳定性分析的本质帮助却被严重低估了。这篇东西没有门槛只要你有一台装了Python的电脑哪怕只懂一点点基础语法也能跟着我一步步把从传递函数到波特图、再到相位裕度判稳的完整链路搭起来。我会把为什么这么画、为什么要关注这个点、以及我踩过的哪些坑全部写清楚。你最后拿到的不只是一段能出图的代码而是一套以后遇到任何电路都能自己上手分析稳定的方法论。1. 从零搭建环境Python绘图分析需要的东西1.1 安装Python解释器聊代码之前机器上得先有Python。很多新手卡在第一步说法、教程五花八门我直接给你最稳妥的答案去Python官网(python.org)下载对应系统的安装包。安装时有两点必须注意。一是勾选“Add Python to PATH”不勾的话后续在命令行用python命令大概率直接报python was not found这是新手日常遇到的第一道坎。二是安装路径尽量别带中文和空格有人装到C:\程序\python这种路径后面库安装大概率出各种诡异错误。装完在终端里敲一句python --version能显示版本号就说明环境通了。我是Windows 11 Python 3.10.x这套教程里的代码在3.8以上版本都跑得动。1.2 编辑器怎么选VS Code实操配置编辑器我用的是VS Code免费、轻量、插件生态好。去官网下载安装然后做两件关键配置安装官方Python插件在扩展商店搜“Python”认准微软出品那个。按CtrlShiftP输入“Python: Select Interpreter”选你刚装的解释器。这一步很容易忽略但很重要。不选解释器的话VS Code主界面不会推荐补全和调试功能就算代码能在终端跑编辑器也一直报错体验很磨人。还有个小技巧把终端默认改成PowerShell或CMD都行快捷键是Ctrl后续运行脚本直接用python 文件名.py就行。1.3 必须装的三个库我们这次核心用的是三个库numpy做数学运算matplotlib做绘图control用来做交叉验证。直接用pip安装pip install numpy matplotlib control三个库一起装慢的话加个国内镜像源具体怎么换镜像自己搜一下就行我不展开。装完可以跑个简单测试import numpy as np import matplotlib.pyplot as plt print(np.__version__) print(plt.__version__)打印出版本号说明导入成功。control库如果安装失败通常是依赖的scipy版本不匹配升级一下scipy再装就行。2. 波特图到底在画什么先看数学本质2.1 传递函数、幅频与相频一个类比说透很多新人对波特图的第一印象是一堆看不懂的曲线。我习惯用一个生活化类比来解释把电路当成一个音频均衡器你输入一个频率的声音它会告诉你“这个频率被放大了多少”和“这个频率被延迟了多少拍”。放大了多少体现在幅频特性上单位是dB计算公式是20 * log10(|H(jω)|)。延迟了多少拍体现在相频特性上单位是度或弧度本质是arg(H(jω))。波特图就是一条横轴以频率对数刻度、纵轴分别显示增益和相位的两张图。用对数横轴不是因为习惯而是因为电路响应覆盖的频率范围跨度极大可能从1Hz到100MHz只有对数轴才能把低频和高频行为放在同一张图里看清楚。2.2 极点和零点如何塑造曲线形状波特图的曲线形状不是随机跳变的它完全由传递函数的极点和零点决定。这个知识点是看懂图的关键。一个简单的一阶低通滤波器传递函数是H(s) 1 / (1 sRC)它有一个极点位置在s -1/(RC)对应频率f_p 1/(2πRC)。这个频率称为转折频率或极点频率。在它之前增益平坦在它之后增益以每十倍频程-20dB的速度滚降。相位在转折频率处正好是-45°远低于转折频率时趋近0°远高于时趋近-90°。二阶系统有两个极点情况复杂一些。如果没有阻尼或阻尼很弱在转折频率附近会出现增益尖峰相位下降也会更陡峭。这个尖峰往往就是不稳定或振荡隐患的前兆。所以看波特图首先就是找转折频率然后观察转折频率附近曲线的“斜率”和“尖峰形态”这能快速告诉你系统的响应模式。2.3 为什么只看增益还不够相位同样要命很多初学者只看幅频特性觉得增益曲线在系统带宽内保持平坦就行。这是一个会出大事的误区。负反馈系统的稳定性不只取决于增益还取决于环路在增益穿过0dB那一点的附加相移。如果附加相移接近-180°负反馈就会变成正反馈系统就会振荡。这就是所谓“180°相位瓶颈”的概念。所以稳定性的判断逻辑是这样的先找到增益曲线穿过0dB的频率叫穿越频率或增益交点频率然后看这个频率对应的相位和-180°差多少。这个差值就是相位裕度。工程上相位裕度大于45°才算稳定可靠典型目标在45°到60°之间。相位裕度过小系统的阶跃响应会有明显的振铃稍微有点外界扰动就可能导致振荡。这个逻辑是后面所有分析的主心骨。3. 用Python手绘波特图代码一步步拆解3.1 最简RC低通滤波器让图画出来先从最简单的RC开始。假设电阻R1kΩ电容C100nF转折频率约为1.59kHz。代码思路是这样的定义一个传递函数输入频率数组输出复数增益然后用np.logspace生成从10Hz到1MHz的对数分布频率点分别计算幅值和相位。import numpy as np import matplotlib.pyplot as plt def rc_lowpass(freq, R1e3, C100e-9): w 2 * np.pi * freq H 1 / (1 1j * w * R * C) return H f np.logspace(1, 6, 500) # 从10Hz到1MHz500个点 H rc_lowpass(f) mag 20 * np.log10(np.abs(H)) phase np.angle(H, degTrue) fig, (ax1, ax2) plt.subplots(2, 1, figsize(8, 6), sharexTrue) ax1.semilogx(f, mag, lw2, colorcrimson) ax1.axhline(y-3, colorgray, ls--, lw1) ax1.set_ylabel(Gain (dB)) ax1.grid(True, whichboth, ls:) ax2.semilogx(f, phase, lw2, colornavy) ax2.axhline(y-45, colorgray, ls--, lw1) ax2.set_xlabel(Frequency (Hz)) ax2.set_ylabel(Phase (deg)) ax2.grid(True, whichboth, ls:) plt.tight_layout() plt.show()跑起来你会看到典型的低通曲线。增益在低频接近0dB在1.59kHz附近正好到-3dB点之后以每十倍频程-20dB的斜率滚降。相位在转折频率处是-45°。这里有个细节容易让人困惑为什么相位在很低频率时不是正好0°曲线看着有些毛刺因为np.logspace生成的点是等对数间隔的低频段样本量相对少加上相位计算本身的浮点误差曲线上会出现微小锯齿。解决方法是增加点数比如500改到1000或者单独对低频段加密。3.2 画布美化与标注让图能放进报告里做工程分析图不能只给自己看还要能放进测试报告或技术评审里。所以必要的标注和美化不能省。我通常会在幅频图上标注0dB参考线在相频图上标注-180°参考线并把穿越频率、相位裕度、增益裕度直接标到图上。这个习惯在做稳定性分析时非常有用看图不用自己拿尺子量。def crossover_freq(freq, mag_db): idx np.where(mag_db 0)[0] if len(idx) 0: return None i idx[-1] f1, f2 freq[i], freq[i1] m1, m2 mag_db[i], mag_db[i1] return 10**(np.interp(0, [m1, m2], [np.log10(f1), np.log10(f2)]))这是一个线性插值找穿越频率的函数精度足够工程使用。原理很简单在增益大于0dB的最后一个点和小于0dB的第一个点之间做线性插值找到增益恰好等于0dB的频率。如果你觉得手写太麻烦也可以用scipy.interpolate做更精细的插值。不过对于伯德图这种平滑曲线线性插值已经非常精确了。3.3 从RC延伸到二阶系统更接近真实电路RC只是入门实际工程里电源环路、运放电路、锁相环大多是二阶或更高阶系统。我们用一个典型的二阶低通滤波电路举例H(s) ω₀² / (s² 2ζω₀·s ω₀²)假设自然频率f₀10kHzω₀2π×10⁴阻尼比ζ分别取0.3、0.5、1.0。同一张图画出三个阻尼比下的伯德图你能直观看到阻尼比对峰值和相位的影响。def second_order(freq, f010e3, zeta0.5): w0 2 * np.pi * f0 s 1j * 2 * np.pi * freq H w0**2 / (s**2 2*zeta*w0*s w0**2) return H f np.logspace(1, 6, 1000) fig, (ax1, ax2) plt.subplots(2, 1, figsize(9, 7), sharexTrue) for zeta, color in [(0.3, red), (0.5, blue), (1.0, green)]: H second_order(f, zetazeta) ax1.semilogx(f, 20*np.log10(np.abs(H)), lw2, colorcolor, labelfζ{zeta}) ax2.semilogx(f, np.angle(H, degTrue), lw2, colorcolor, labelfζ{zeta}) ax1.axhline(y0, colorblack, lw0.8, ls-) ax1.grid(True, whichboth, ls:) ax1.legend() ax1.set_ylabel(Gain (dB)) ax2.axhline(y-180, colorblack, lw0.8, ls-) ax2.grid(True, whichboth, ls:) ax2.legend() ax2.set_xlabel(Frequency (Hz)) ax2.set_ylabel(Phase (deg)) plt.tight_layout() plt.show()跑完你会发现阻尼比0.3时幅频曲线在10kHz附近明显凸起峰值超过0dB达好几个dB相位曲线在转折区也更陡。这个凸起就是谐振峰在实际电路中对应着Q值过高、容易振铃甚至自激的风险。而阻尼比1.0时曲线平坦得多相位变化也缓和系统更稳定但带宽相对保守。做电路稳定性分析时我通常会画出至少三组阻尼比或三组补偿参数下的曲线对比这样选参数时能直观看到趋势比单看一组数据可靠得多。4. 从波特图看电路稳定性实战案例推送4.1 增益裕度与相位裕度怎么读才稳判断稳定性的核心指标有两个相位裕度Phase Margin, PM和增益裕度Gain Margin, GM。相位裕度增益曲线穿越0dB时相位离-180°还有多远。记为PM args(H(jω_cross)) 180°其中ω_cross是穿越角频率。PM0是稳定必要条件工程上建议至少45°最好60°。增益裕度相位曲线穿越-180°时增益低于0dB的程度。记为GM 0dB - |H(jω_180)|。工程上GM至少6dB最好能到10dB以上。这两个指标一个看幅度快、一个看相位快箭在弦上。如果PM太小那么环路在带宽附近很容易还没衰减就先振荡起来如果GM太小意味着高频处存在一个潜在振荡条件噪声或者瞬态扰动就可能激发自激。4.2 实战案例一个带反馈的运放电路来看一个更有工程味的例子。假设有一个单位增益跟随器使用一个单极点运放开环增益A01000080dB单位增益带宽f_t1MHz。那么它的开环传递函数可以近似为A(s) A0 / (1 s/ω_p0)把闭环接成单位增益反馈系数β1环路增益T(s)A(s)·β。我们画出T(s)的波特图就是直接看开环特性。用上一节的方法计算出穿越频率然后找到对应相位算出PM。你会发现这个系统PM接近90°非常稳定。原因在于单极点系统的相位最多滞后90°永远不会到-180°。但如果运放本身在高频处还有一个次极点真实运放都有那么传递函数变成A(s) A0 / ((1 s/ω_p0)(1 s/ω_p1))这时PM就会下降。设f_p15MHz跑一次代码def opamp_loop(freq): A0 10000 fp0 f_t / A0 # 主极点: 100Hz fp1 5e6 s 1j * 2 * np.pi * freq A A0 / ((1 s/(2*np.pi*fp0)) * (1 s/(2*np.pi*fp1))) return A f np.logspace(1, 7, 10000) H opamp_loop(f) mag 20*np.log10(np.abs(H)) phase np.angle(H, degTrue) fc crossover_freq(f, mag) if fc: phase_fc np.interp(np.log10(fc), np.log10(f), phase) pm phase_fc 180 print(f穿越频率: {fc/1e3:.1f} kHz) print(f相位裕度: {pm:.1f}°)跑出来PM可能是50°到70°说明系统还能稳定工作但裕度已经下降。当你把次极点往低频推比如8MHz推到2MHzPM会进一步下降甚至低于30°这时候你去看时域阶跃响应就能看到明显的振铃。这个例子说明了一个关键道理波特图上的相位裕度不是凭空而来的它直接对应时域响应的阻尼程度。理解这个对应关系波特图就活了。4.3 用control库交叉验证避免自己算错的保险手写代码画图的好处是理解深但手写也有笔误风险。我通常会用control库做交叉验证两套结果一比对心里就有底了。import control as ct # 二阶系统 num [w0**2] den [1, 2*zeta*w0, w0**2] sys ct.TransferFunction(num, den) mag, phase, omega ct.bode(sys, PlotFalse, omeganp.logspace(1, 6, 1000))ct.bode(..., PlotFalse)返回频率响应数据画法跟matplotlib完全兼容。我用这个方法确认过很多次手写计算结果和库计算的结果基本完全重合差不了几个小数点。这里有个实用技巧omega参数不传的话control库会自动选择频率范围但往往不够密转折频率附近会显得粗糙。手动传入np.logspace生成的足够密集的频率点画出来的曲线更平滑更关键的是穿越频率的插值精度更高。5. 常见问题与排查技巧这些坑值得花时间记下来5.1 对数坐标的使用与频率范围选择很多新人在画波特图时容易犯一个错误直接用plot(f, mag)而不是semilogx画图。这样画出来的曲线在低频段挤成一团高频段拉得远远的完全看不出转折频率的位置。正确做法是横轴用对数刻度频率范围至少覆盖你关心的转折频率左右各两个十倍频程。比如转折频率1.59kHz那么频率范围至少要从16Hz到159kHz最稳妥的是10Hz到1MHz。还有一个小细节np.logspace(1, 6, 500)生成的是10^1到10^6之间的500个点在对数坐标下等距。如果你直接对频域数据做插值必须先把频率取对数再插值否则会出现严重误差。5.2 相位曲线的±180°跳跃问题这是up主们问得最多的一个问题。用np.angle(H, degTrue)算相位结果会限定在[-180°, 180°]之间。但很多系统的真实相位是连续下降的比如二阶系统从0°一路降到-360°在高频端这时np.angle会在-180°处突然跳到180°图形上看起来像一条“锯齿”或“断裂”线。解决办法是用np.unwrap。这个函数专门处理相位角度跳变它会自动把跳变的相位加上2π的整数倍恢复成连续曲线。phase np.unwrap(np.angle(H)) * 180 / np.pi经过unwrap之后相位曲线会从0°平滑连续地下降到-360°无论是计算相位裕度还是画图都舒服多了。还有一个细节unwrap默认是针对弧度制的所以要先把角度转成弧度unwrap之后再转回度顺序别搞反。5.3 增益曲线不平滑的元凶曲线如果感觉“毛刺”很多第一检查频率点数够不够。我习惯至少用500点复杂系统直接上1000点或更多。第二检查是不是有数值溢出特别是在计算1/(10j)这种复数运算时某些极端频率下浮点精度会下降增益出现小尖刺。还有一种情况是接近奈奎斯特频率或者频率范围过大时np.logspace在高频段的点数密度不足导致转折区曲线精度不够。解决办法是分两段绘制低频段用密一点的对数点高频段用另外一组点最后拼图。不过大多数情况下1000点已经足够。5.4 标注穿越频率的实用小工具稳定分析里最常用的一个操作就是找穿越频率和相位裕度。我写了一个很小的函数可以直接嵌入代码里复用def stability_metrics(freq, mag_db, phase_deg): idx np.where(mag_db 0)[0] if len(idx) 0: return None, None i idx[-1] f_cross 10**np.interp(0, [mag_db[i], mag_db[i1]], [np.log10(freq[i]), np.log10(freq[i1])]) pm np.interp(np.log10(f_cross), np.log10(freq), phase_deg) 180 idx180 np.where(phase_deg -180)[0] gm None if len(idx180) 0: j idx180[0] f_180 10**np.interp(-180, [phase_deg[j-1], phase_deg[j]], [np.log10(freq[j-1]), np.log10(freq[j])]) gm -np.interp(np.log10(f_180), np.log10(freq), mag_db) return f_cross, pm, f_180, gm return f_cross, pm, None, None这个函数返回四个值穿越频率、相位裕度、相位到达-180°时的频率、增益裕度。有了它每次分析数据只要调一行代码不用重复手算。不过要提醒一下这个函数假设数组是按频率单调递增排列的而且穿越频率发生的位置是你关心的高频段。如果遇到多穿越的情况比如复杂多环系统需要额外判断取哪一个穿越点当主穿越频率。5.5 别把噪声当振荡结合时域去验证最后分享一个实际的工程体会。只看波特图有时候会被某些现象迷惑比如曲线上有个小尖峰你以为是振荡前兆但实测时域却很干净。这时不要急着改电路先用示波器看下实际波形区分到底是真实的谐振峰还是测量噪声。我最近调试一个电源环路的案例波特图显示相位裕度只有35°理论上风险不低。但我的实测时域波形只有轻微振铃并没有持续振荡。后来深挖发现问题出在环路里一个电容的ESR在转折频率附近引入了一个额外的零点导致相位测量出现偏差实际相位裕度比计算值好很多。所以纯粹依赖理论计算和理想模型会在真实的电路里栽跟头。正确做法是先用Python手绘曲线得到理论预判再用网络分析仪或信号源示波器实测闭环波特图两者对比找差异归因到每一个寄生参数上。说到这我想起自己第一次画出的波特图歪歪扭扭的曲线里带着锯齿程序跑的几百行代码有近一半在处理格式问题但那之后我就不再怕稳定性分析。后来做电源环路补偿、运放电路调试、甚至给师弟讲怎么分析DCDC的环路用的都是这套思路。只是建议你答辩或出报告之前务必用真实的实测波形验证一版理论图那才是硬道理。
返回列表