简介：MATLAB环境下基于粒子群优化（PSO）与部分传输序列（PTS）的OFDM峰均功率比（PAPR）抑制仿真源码，面向无线通信、信号处理方向的工程师与研究者，可用于学习OFDM系统中降低PAPR的优化实现。压缩包内共4个文件，均为.m脚本，包含主运行程序Runme.m、子函数binary2factor.m、get80216map.m以及核心算法PSO_PTS.m，整体仅4KB，结构精简，便于快速阅读和二次开发。目前已有232人学习下载。通过源码可完整跟踪OFDM符号生成、PTS分割相位旋转、PAPR计算及PSO迭代寻优的全流程，直观对比优化前后的PAPR曲线；同时还能掌握粒子群算法在通信信号优化问题中的参数设置与收敛判断方法，为后续改进其他PAPR抑制技术（如SLM、限幅等）提供基础。适合具备一定MATLAB和OFDM基础的学习者用来实验验证和算法扩展。
1. 为什么OFDM系统的PAPR抑制要选PSO优化PTS
做OFDM物理层仿真的人，早晚要撞上PAPR这个指标。子载波数一到256，IFFT之后时域信号的峰值功率普遍比平均功率高10dB以上，射频功放稍微进入压缩区，EVM和邻道泄漏一起恶化。PTS是经典的抑制手段，可它的相位组合搜索量随分块数指数膨胀，V=12时穷举2^11种组合，单个符号就要算两千多次合成信号，整套CCDF仿真跑下来非常吃力。把PSO粒子群优化拿来做相位因子搜索，不保证全局最优，但能顺着群体搜索的思路在几十次迭代里逼近穷举效果，这是工程仿真里最常见的折中方案。下面把OFDM原理、PTS分块、PSO寻优串成一条可复现的MATLAB链路，代码、参数和坑位都摆出来，适合正在做PAPR抑制仿真或者需要交付仿真结果的读者。
2. PTS分块与PSO寻优：OFDM系统中PAPR抑制的机理拆解
2.1 OFDM的PAPR从哪来：IFFT叠加的统计本质
OFDM发射端把频域符号X_k经N点IFFT变换到时域，第n个采样点是x_n = (1/sqrt(N)) Σ X_k·e^{j2πkn/N}。N较大且各子载波符号独立同分布时，中心极限定理使得x_n的实部和虚部逼近高斯分布，幅度服从瑞利分布。这意味着某个采样点上所有子载波恰好同相叠加的概率并不低，峰值功率可以远超平均功率。PAPR定义为10·log10(max|x_n|² / E[|x_n|²])，理论极限是10·log10(N) dB，N=256时上限约24dB，但即便在99.9%概率下，实际信号也会出现10dB以上的瞬时峰值。功放若按平均功率设计，这些峰值会被削掉，带内失真和带外辐射一起冒出来。
有个细节容易被忽略：真实DAC和ADC在连续时间上看到峰值，按符号率采样会漏掉峰。所以仿真里普遍对IFFT输出做L倍过采样，L=4是工程默认值。过采样不改变PAPR的定义，但改变统计结果，后面CCDF对比必须统一过采样倍数，否则两条曲线对应的是两套不同口径的系统。
2.2 PTS分块策略与相位因子搜索空间
PTS的核心思路是把N个子载波分成V个互不重叠的子块，各子块分别IFFT，再用相位因子b_v加权叠加。IFFT是线性运算，合成信号可以写成x = Σ b_v·ifft(X_v)，搜索目标就是找一组b_v让合成信号的峰值最小。分块方式有三种常见选择，对仿真结论影响不小：
| 分块方式 | 实现复杂度 | 频域相关性 | PAPR抑制上限 | PSO收敛难度 |
|---|---|---|---|---|
| 相邻分块 | 最低 | 高 | 偏低 | 容易 |
| 交织分块 | 低 | 低 | 高 | 中等 |
| 伪随机分块 | 中 | 低 | 高 | 中等 |
b_v一般从W个等间隔相位里取，比如{1,-1}对应W=2，{1,-1,j,-j}对应W=4。搜索空间是W^(V-1)，减一是因为整体乘一个常数相位不影响PAPR，第一个子块的相位可以固定为1。V=8、W=2时共128种组合，穷举没有问题；V=16、W=2时是32768种，穷举还能忍；V=8、W=4时是4^7=16384种，V=16、W=4直接到4^15约10亿量级。搜索空间一旦过千，逐符号穷举就跑不动了，这正是PSO介入的时机。
2.3 PSO算法在相位空间里的搜索方式
PSO算法把每个候选解抽象成空间里的一个粒子，位置向量就是一组合法的相位因子，适应度就是合成信号的PAPR。每次迭代，粒子向自己的历史最优pbest和群体最优gbest靠拢，核心更新可以写成下面的MATLAB形式：
% PSO单步更新，r1、r2是[0,1]均匀随机数矩阵 r1 = rand(size(pos)); r2 = rand(size(pos)); vel = w*vel + c1*r1.*(pbest - pos) + c2*r2.*(gbest - pos); pos = pos + vel;w是惯性权重，c1、c2是加速常数。r1、r2逐维独立随机，避免所有粒子沿着同一条路径冲向当前最优，保住群体多样性。PTS问题有个特性让PSO特别划算：V个子块的时域信号在搜索前一次性算好并缓存，评估适应度时只需做V次复数加权求和再算一次峰值，复杂度是O(N·L·V)。粒子数30、迭代30次的搜索总量在千次量级，远小于穷举上万次组合的代价。
实际实现里通常做连续松弛：粒子位置先按连续值在[-1,1]里搜索，最后量化到最近的合法相位。连续搜索比直接在离散集合上跳变收敛更平滑，但量化会引入微小偏差，第4章讲处理办法。
3. MATLAB实现：从标准PTS到PSO优化的完整仿真脚本
3.1 仿真环境与参数表
下面的仿真只用MATLAB基础函数，不依赖通信工具箱、优化工具箱或全局优化工具箱，手写IFFT链路和PSO主循环，方便移植到Octave或改写成C风格代码。参数统一放在脚本头部，改配置时只动这一块。表里是一组常用起点，对应256子载波、QPSK、4倍过采样的OFDM基带系统。
| 参数 | 取值 | 说明 |
|---|---|---|
| N | 256 | 子载波总数 |
| V | 4 | PTS分块数 |
| W | 2 | 每块可选相位数，即{1,-1} |
| L | 4 | 过采样倍数 |
| 调制方式 | QPSK | 频域符号能量归一化 |
| M | 2000 | CCDF统计的OFDM符号数 |
| nP | 30 | PSO粒子数 |
| maxIter | 30 | PSO最大迭代数 |
3.2 OFDM发射链路与PTS分块的MATLAB代码
N = 256; V = 4; L = 4; M = 2000; % CCDF统计样本数 % QPSK符号映射，能量归一化到1 data = randi([0 3], N, 1); X = (1/sqrt(2)) * ((2*(bitget(data,2)-0.5)) + 1j*(2*(bitget(data,1)-0.5))); % 相邻分块：把N个子载波均分为V段，分别放入V个子块 block_len = N / V; Xv = zeros(N, V); for v = 1:V idx = (v-1)*block_len + 1 : v*block_len; Xv(idx, v) = X(idx); end % 每个子块做L倍过采样IFFT，结果缓存为xv矩阵 xv = zeros(N*L, V); for v = 1:V Xup = oversample(Xv(:, v), L); % 频谱映射到N*L点 xv(:, v) = ifft(Xup, N*L) * sqrt(L); % sqrt(L)保持时域功率一致 endoversample辅助函数把N点频域符号映射到N*L点频谱：前半段子载波放低频区，后半段放高频区，中间补零完成过采样插值。IFFT输出乘以sqrt(L)是因为补零后Parseval能量摊到了更长的序列上，不乘这个系数，时域平均功率会被压到1/L，PAPR计算跟着失真。xv矩阵是整个PTS搜索的公共缓存，V列波形只算一次，后面所有候选解都在这个基上做加权叠加，这是PTS能高效跑起来的根基。
3.3 PSO寻优主循环与适应度函数
nP = 30; maxIter = 30; w = 0.9; c1 = 1.8; c2 = 1.8; % 粒子位置是V-1维连续向量，第v维对应第v+1个子块的相位 pos = 2*rand(nP, V-1) - 1; vel = zeros(nP, V-1); pbest = pos; gbest = pos(1, :); % 初始化个体最优和群体最优 fit = zeros(nP, 1); for i = 1:nP fit(i) = deal_papr(pos(i, :), xv); end pbest_fit = fit; [gbest_fit, gi] = min(fit); gbest = pos(gi, :); for t = 1:maxIter for i = 1:nP r1 = rand(1, V-1); r2 = rand(1, V-1); vel(i, :) = w * vel(i, :) + c1 * r1 .* (pbest(i, :) - pos(i, :)) ... + c2 * r2 .* (gbest - pos(i, :)); pos(i, :) = pos(i, :) + vel(i, :); pos(i, :) = max(-1, min(1, pos(i, :))); % 位置钳位到[-1,1] f = deal_papr(pos(i, :), xv); if f < pbest_fit(i) pbest_fit(i) = f; pbest(i, :) = pos(i, :); end if f < gbest_fit gbest_fit = f; gbest = pos(i, :); end end % 惯性权重从0.9线性衰减到0.4 w = max(0.4, w - (0.9 - 0.4) / maxIter); end % 连续位置量化到{+1,-1}，得到最终相位因子 b_opt = [1, sign(gbest)]; x_pts = xv * b_opt.'; % PTS合成后的时域信号function papr = deal_papr(b, xv) bv = [1, b]; % 补上固定的第一个子块相位 x = xv * bv.'; % V个时域子块加权叠加 papr = 10 * log10(max(abs(x).^2) / mean(abs(x).^2)); end搜索空间是V-1维不是V维，因为第一个子块相位固定为1省略了全局相位模糊。位置钳位到[-1,1]是因为最终量化目标是{+1,-1}，连续区间越贴近离散集合，量化损失越小。速度更新里r1、r2逐维独立，防止粒子群过早抱团；惯性权重前期大保持探索能力，后期小做局部精修，这是PSO参数设计里最经典的节奏。
3.4 CCDF统计与曲线对比
papr_orig = zeros(M, 1); papr_pts = zeros(M, 1); for m = 1:M data = randi([0 3], N, 1); X = (1/sqrt(2)) * ((2*(bitget(data,2)-0.5)) + 1j*(2*(bitget(data,1)-0.5))); x_orig = ifft(oversample(X, L), N*L) * sqrt(L); papr_orig(m) = 10*log10(max(abs(x_orig).^2) / mean(abs(x_orig).^2)); % run_pso_pts是3.2和3.3代码封装成的函数，对单个符号跑完整PSO papr_pts(m) = run_pso_pts(X, N, V, L, nP, maxIter); end th = 0:0.2:12; ccdf_orig = arrayfun(@(p0) mean(papr_orig > p0), th); ccdf_pts = arrayfun(@(p0) mean(papr_pts > p0), th); figure; semilogy(th, ccdf_orig, 'k-', 'LineWidth', 1.5); hold on; semilogy(th, ccdf_pts, 'r-', 'LineWidth', 1.5); grid on; xlabel('PAPR_0 (dB)'); ylabel('Pr(PAPR > PAPR_0)'); legend('原始OFDM', 'PSO-PTS', 'Location', 'southwest');CCDF的统计口径是：每个OFDM符号独立跑一次完整PSO，M个符号的PAPR值构成经验分布，超过门限p0的比例就是该点的CCDF估计。M=2000时曲线已经平滑，再往上到1e4需要多等几分钟。两个分支必须复用同一个oversample实现，原始OFDM和PTS分支如果用了不同的过采样方式，对比的就是两套系统，结论没有意义。
4. 参数怎么设：PSO-PTS在OFDM仿真里的收敛性和复杂度权衡
4.1 PSO参数取值表与收敛性判断
| 参数 | 常见区间 | 默认值 | 调参说明 |
|---|---|---|---|
| w | 0.4~0.9线性衰减 | 0.9→0.4 | 固定小w容易早熟，大w收敛慢 |
| c1, c2 | 1.5~2.5 | 1.8, 1.8 | c1大偏向个体经验，c2大偏向群体导向 |
| nP | 20~60 | 30 | V增大时按比例加大到40~60 |
| maxIter | 15~50 | 30 | 看gbest_fit末段是否还在下降 |
| 位置范围 | [-1,1] | [-1,1] | 对应二值相位{+1,-1}的连续松弛 |
收敛性判断不要只看最终gbest_fit的值，要看每一轮gbest_fit的下降曲线。V=8、W=2时通常10~20次迭代就把PAPR压到接近穷举值，后面继续迭代收益很小；V=16时搜索空间大了两个数量级，30次迭代后曲线仍缓慢下行，这时候把maxIter加到50~80，或者把粒子数翻倍，效果比单方面增大迭代数更明显。记录每轮的gbest_fit画出来，是判断参数是否合理最直接的手段。
4.2 分块数V和相位备选集W怎么搭配
V和W直接决定搜索空间大小，也决定PTS需要发给接收端的边带信息量。每多一个子块，接收端就要多知道一份相位信息，边带是真实系统里要付出的代价。V=4、W=2时总共只有2^3=8种组合，穷举一次只要算8次合成信号，PSO没有任何优势，老老实实穷举更稳。PSO真正划算的区间是搜索空间大于1000的配置，比如V=8、W=4对应16384种，或V=16、W=2对应32768种。
分块方式也影响PSO的搜索难度。相邻分块的子块频域相关性高，PAPR抑制上限偏低，但搜索曲面平滑，PSO容易收敛；交织分块和伪随机分块的上限更高，但曲面更崎岖，容易陷进局部最优。实测中V=8、W=2时交织分块配合PSO能比相邻分块多0.3~0.5dB的增益，代价是收敛需要的迭代数多一半左右。
4.3 连续松弛量化偏差和早熟收敛的处理
连续PSO搜到的是[-1,1]里的实数值，最后sign()一量化，最终解的PAPR通常比搜索过程中看到的gbest_fit高0.1dB左右。V较小时偏差不明显，V大了粒子维度高、变量间相关性复杂，偏差会放大。我一般做第二轮离散局部搜索补偿量化损失：
% 围绕gbest做离散局部搜索：逐维翻转相位，保留改善的解 b_try = sign([1, gbest]); best = deal_papr(b_try(2:end), xv); for v = 2:V b_cand = b_try; b_cand(v) = -b_cand(v); % 翻转第v个子块相位 f_cand = deal_papr(b_cand(2:end), xv); if f_cand < best best = f_cand; b_try = b_cand; end end b_opt = b_try;这轮搜索只做V-1次适应度评估，开销可以忽略，但能把sign量化丢掉的那部分增益找补回来。早熟收敛是PSO-PTS最常见的坑，现象是gbest_fit前5次迭代就定型，后面不再下降，但PAPR明显高于穷举基准。常见原因有两个：粒子数太少导致群体多样性不足，或者w衰减太快把探索能力提前掐断。处理办法是把w衰减推迟到后三分之一迭代再进行，前三分之二保持w=0.9；还跳不出就在每次迭代随机重置10%粒子的位置，让群体重新获得多样性。
提示：对比不同PSO配置时，始终保持V、L、调制方式和CCDF符号数相同。PAPR曲线相差0.1dB以内基本是随机波动，判断配置好坏要看多次独立实验的均值，而不是单次结果。
5. 用穷举基准校验PSO-PTS的优化效果上限
5.1 小规模穷举gap测试
建议养成一个习惯：在小参数下用穷举算出理论最优PAPR，再跑同样配置的PSO，两者之差就是PSO的优化gap。以V=8、W=2为例，穷举128种组合，对每个OFDM符号记录穷举最优PAPR和PSO结果，统计两者的平均差。这个gap通常小于0.1dB，说明参数设置合理；大于0.3dB就需要回到第4章调参。gap测试的关键价值在大V场景：V=16时无法穷举，但可以拿小V下的gap大致预估PSO在大V下的优化余量，心里有底。习惯用Simulink搭OFDM调制解调模型的话，把run_pso_pts封装成MATLAB Function模块，就能在模型里直接观测PAPR波形。
5.2 仿真时间预算与迭代早停
仿真耗时的瓶颈在每符号的迭代搜索。M=2000、nP=30、maxIter=30时单符号评估900次适应度，V=8、L=4、N=256的配置跑完整个CCDF大约几分钟，机器不同差异很大。批量调参时先用M=200粗跑确认趋势，再上完整统计，能省大量时间。用tic/toc包住主循环实测运行时间，同时记录gbest_fit的变化轨迹：如果最后10次迭代完全没有变化，说明maxIter给多了，砍半再跑结果一致；如果结果变差，说明之前刚好在临界迭代数附近，往回调到20~25。这种时间换精度的边界测量，比照搬经验值可靠得多。
本文还有配套的精品资源，点击获取
