简介这份资源面向通信工程、电子信息类专业学生及数字通信初学者聚焦连续相位调制CPM在MATLAB环境下的原理验证与仿真实现帮助读者理解调制指数、符号速率与信息速率之间的关系并掌握MSK、GMSK等典型CPM方案的建模思路。压缩包共7个文件全部为m脚本整体约2KB涵盖参数配置、信号生成、性能分析与结果可视化等模块结构紧凑便于逐文件阅读与二次修改。已有296人学习下载适合作为课程实验或自学练手的参考素材。读者可借助这些脚本完成CPM信号的相位累加、滤波处理与星座图绘制并结合高斯白噪声信道模拟观察误码表现从而把抽象的调制理论落到可运行的代码上快速建立对CPM通信链路的直观认识。1. 从一份 CPM.rar 说起连续相位调制到底解决了什么问题如果你手头正好有一个叫CPM.rar的压缩包里面大概率躺着一套 MATLAB 的 CPM 调制仿真脚本——可能是某个通信原理课程设计也可能是水声通信、卫星通信项目里抽出来的调制模块。CPM全称 Continuous Phase Modulation连续相位调制它要解决的核心问题很朴素如何在功率受限、带宽受限的信道上把数据送得更远、更省频谱同时不让功放因为包络跳变而失真。和 PSK、QAM 这类相位可以突变的调制方式不同CPM 的相位是连续变化的这个「连续」二字直接决定了它的频谱旁瓣衰减极快带外泄漏小对非线性功放友好。它适合谁做水声通信、深空通信、物联网低功耗链路、软件无线电的工程师以及正在啃《现代水声通信原理与 MATLAB 应用》这类教材、需要把公式落成波形的研究生。这一章先把 CPM 的数学骨架和 MATLAB 里的落地路径讲清楚后面几章再一步步把CPM.rar里那套东西复现出来。CPM 的一般表达式可以写成s(t) sqrt(2E/T) * cos(2πf_c t φ(t, I) φ0)其中 φ(t, I) 是携带信息的相位它由调制指数 h、脉冲成形函数 g(t) 和符号序列 I 共同决定。关键点在于 φ(t, I) 是连续的哪怕符号在切换相位也不会跳变。这个连续性带来两个直接好处一是频谱效率二是可以用限幅器或非线性功放而不产生严重的频谱再生。代价是接收端要做序列检测因为当前符号的相位状态和前面若干符号有关这就是 CPM 的「记忆性」。理解这一点后面写 MATLAB 代码时才知道为什么要用 Viterbi 或者 BCJR 做解调而不是简单的一符号一判决。在 MATLAB 里做 CPM常见做法有三条路一是用 Communications Toolbox 自带的comm.CPMModulator和comm.CPMDemodulator系统对象适合快速验证二是自己从相位递推公式写起适合理解原理和改参数三是用cpmmod、cpmdemod这类老函数兼容旧代码但灵活性差。CPM.rar里那套脚本我猜多半是第二条路——自己写相位累加因为只有这样才能把 h、L、g(t) 这些参数完全握在手里。接下来就从参数选型开始把这条路走通。2. CPM 参数选型h、L、M 和脉冲成形怎么定2.1 调制指数 h 和关联长度 L 的取舍调制指数 h 决定了相位每符号增量的大小。h 越大相位变化越快频谱越宽但符号间距离越大误码性能越好。h 越小频谱越紧凑但接收端区分状态越困难。常见取值有 h0.5MSK 类、h1/3、h1/4以及各种有理数 hp/q因为有理数 h 能让相位状态有限接收端可以用有限状态机处理。如果 h 是无理数相位状态无限多实际系统没法做网格搜索所以工程上几乎只用有理数 h。关联长度 L 表示一个符号的脉冲成形影响多少个符号周期。L1 叫全响应 CPM相位只和当前符号有关接收端简单但频谱一般L1 叫部分响应 CPM频谱更紧凑但接收端要记忆 L 个符号复杂度指数上升。我一般建议新手从 L1、h0.5 开始跑通后再试 L2 或 L3观察频谱和误码率的变化。CPM.rar里如果同时有L1和L3的脚本那作者大概率是在对比部分响应的收益。2.2 脉冲成形函数 g(t) 的三种常见形式g(t) 是相位脉冲它的积分就是频率脉冲。常见形式有类型表达式特点频谱特性实现难度REC矩形持续 L 个符号旁瓣衰减慢最低RC升余弦平滑旁瓣衰减快中等GMSK高斯滤波后的矩形极紧凑但引入 ISI较高在 MATLAB 里REC 脉冲最容易写就是一个长度为 L 的矩形窗RC 需要按升余弦公式生成GMSK 则要先做高斯滤波再积分。CPM.rar里如果出现gmsk字样那多半是 GMSK 的变体。选哪种取决于你的频谱模板有多严。水声通信里常用 RC 或 GMSK因为带外泄漏必须压得很低。2.3 用 MATLAB 生成相位脉冲和相位轨迹下面这段代码演示如何生成 REC 和 RC 的相位脉冲并画出相位轨迹。你可以直接抄进 MATLAB 跑。% 参数设置 M 2; % 进制数2 表示二进制 CPM h 0.5; % 调制指数 L 3; % 关联长度 sps 8; % 每符号采样点数 N 20; % 符号数 T 1; % 符号周期归一化 % 生成 REC 相位脉冲 t_pulse linspace(0, L*T, L*sps); g_rec ones(1, L*sps) / (L*sps); % 矩形面积归一化为 1/2? 注意归一化 % 生成 RC 相位脉冲 beta 0.3; % 滚降系数 g_rc zeros(1, L*sps); for i 1:length(t_pulse) ti t_pulse(i); if ti 0 || ti L*T g_rc(i) 0; elseif ti (1-beta)*T/2 g_rc(i) 1; elseif ti (1beta)*T/2 g_rc(i) 0.5 * (1 cos(pi/(beta*T) * (ti - (1-beta)*T/2))); else g_rc(i) 0; end end g_rc g_rc / sum(g_rc); % 归一化使总面积为 1 % 生成随机符号序列 data randi([0 M-1], 1, N); % 映射为 ±1二进制 CPM symbols 2*data - 1; % 相位累加 phase zeros(1, N*sps); current_phase 0; for k 1:N % 当前符号对应的相位增量 for j 1:sps idx (k-1)*sps j; % 简化用 g_rec 的累加近似 current_phase current_phase 2*pi*h*symbols(k)*g_rec(min(j, length(g_rec))); phase(idx) current_phase; end end % 画相位轨迹 figure; plot(phase); xlabel(采样点); ylabel(相位 (rad)); title(CPM 相位轨迹 (REC, h0.5, L3)); grid on;这段代码的关键在相位累加那几行。2*pi*h*symbols(k)*g_rec(...)是相位增量的离散近似实际连续时间公式是积分离散化后用求和代替。注意g_rec的归一化相位脉冲的积分应该是 1/2对于二进制 CPM但这里为了画图方便做了简单归一化真正做误码率仿真时要严格按sum(g)*dt 1/2来。参数sps决定采样率一般取 8 或 16太小会引入离散化误差太大浪费计算。L越大g_rec越长相位累加时要注意索引不要越界。提示如果你用comm.CPMModulator它内部已经处理好了归一化但自己写代码时归一化错了误码率曲线会整体偏移这是最常见的翻车点之一。3. 从零写一个 CPM 调制器相位递推与波形生成3.1 相位状态网格的构建CPM 的相位状态不是连续的而是离散的。对于有理数 hp/q相位状态在 2π 周期内有 q 个可能值如果考虑符号记忆还要乘以 M^(L-1)。构建网格是写解调器的前提。在 MATLAB 里可以用一个结构体数组表示每个状态字段包括当前相位、历史符号。下面是一个简化的网格构建示例。% 构建 CPM 相位状态网格 p 1; q 2; % h p/q 0.5 M 2; % 二进制 L 3; % 关联长度 num_states q * M^(L-1); % 状态数 states struct(phase, {}, history, {}); for i 0:q-1 for j 0:M^(L-1)-1 idx i * M^(L-1) j 1; states(idx).phase 2*pi*i/q; % 将 j 转换为 M 进制历史符号 hist zeros(1, L-1); temp j; for k 1:L-1 hist(k) mod(temp, M); temp floor(temp / M); end states(idx).history hist; end end % 显示前几个状态 for i 1:min(5, num_states) fprintf(状态 %d: 相位 %.2f rad, 历史 %s\n, ... i, states(i).phase, mat2str(states(i).history)); end这段代码构建了所有可能的相位状态。q是 h 的分母M^(L-1)是历史符号的组合数。状态总数随 L 指数增长L3、M2、q2 时是 8 个状态还能接受L5 时就是 32 个L7 时 128 个再大就不适合用网格搜索了。实际工程中如果 L 很大会用简化算法如 reduced-state sequence detection。CPM.rar里如果 L 不超过 3那网格法完全够用。3.2 用相位递推生成调制波形有了状态网格调制就是根据输入符号在网格上转移并生成对应的波形。下面这段代码把符号序列映射成 CPM 波形。% 生成 CPM 调制波形 sps 8; % 每符号采样点数 N 100; % 符号数 data randi([0 M-1], 1, N); symbols 2*data - 1; % 二进制映射为 ±1 % 预生成相位脉冲RECL3 g ones(1, L*sps) / (L*sps); % 注意这里面积归一化为 1实际应为 1/2 % 相位累加 phase zeros(1, N*sps); current_phase 0; for k 1:N for j 1:sps idx (k-1)*sps j; % 当前符号对相位的贡献 current_phase current_phase 2*pi*h*symbols(k)*g(min(j, length(g))); phase(idx) current_phase; end end % 生成复基带波形 fc 10; % 载波频率归一化 t (0:N*sps-1) / sps; waveform exp(1j * phase); % 画频谱 figure; pwelch(waveform, [], [], [], sps, centered); title(CPM 信号功率谱密度);这里waveform exp(1j * phase)是复基带信号实际发射时要乘上exp(1j*2*pi*fc*t)再取实部。pwelch用来估计功率谱你可以看到 CPM 的旁瓣衰减比 QPSK 快很多。参数sps影响频谱估计的精度一般取 8 以上。g的长度是L*sps相位累加时用min(j, length(g))防止越界但更严谨的做法是把g对齐到每个符号的起始位置而不是简单截断。这个细节在CPM.rar里如果处理得不好频谱会出现异常毛刺。注意相位累加时每个符号的贡献应该从该符号周期的起始点开始持续 L 个符号周期。上面代码简化成每个采样点都加实际实现要用卷积或状态机。如果你发现频谱旁瓣降不下去先检查这里。4. CPM 解调Viterbi 算法在 MATLAB 里怎么落地4.1 分支度量的计算CPM 解调的核心是 Viterbi 算法它在状态网格上搜索最优路径。分支度量是接收信号与预期信号之间的欧氏距离。对于每个状态转移预期信号由当前状态相位、历史符号和当前输入符号共同决定。下面代码计算分支度量。% 计算分支度量 % 假设接收信号为 rx长度 N*sps已同步 % 对于每个时刻 k每个状态 s每个输入符号 m计算度量 num_states q * M^(L-1); metrics zeros(num_states, N); % 累积度量 survivor zeros(num_states, N); % 幸存路径 % 初始化 metrics(:, 1) 0; for k 2:N for s 1:num_states best_metric inf; best_prev 0; for m 0:M-1 % 根据状态 s 和历史符号计算预期相位 % 这里简化假设状态 s 的相位已知输入 m 产生新相位 expected_phase states(s).phase 2*pi*h*(2*m-1); expected_symbol exp(1j * expected_phase); % 取接收信号对应段 rx_seg rx((k-1)*sps1 : k*sps); % 计算欧氏距离 metric sum(abs(rx_seg - expected_symbol).^2); % 更新累积度量 total_metric metrics(s, k-1) metric; if total_metric best_metric best_metric total_metric; best_prev s; end end metrics(s, k) best_metric; survivor(s, k) best_prev; end end % 回溯 [~, final_state] min(metrics(:, N)); decoded zeros(1, N); current_state final_state; for k N:-1:2 decoded(k) current_state; current_state survivor(current_state, k); end decoded(1) current_state;这段代码是 Viterbi 的骨架但做了大量简化。实际实现中状态转移不是简单加一个相位而是要根据历史符号和输入符号重新计算相位增量并且要考虑脉冲成形的影响。expected_phase的计算需要查表或实时计算不能只加一个常数。rx_seg的长度是sps但 CPM 的符号间干扰会跨越多个符号周期所以预期信号应该用完整的相位轨迹生成而不是单个复指数。这些细节决定了误码率性能CPM.rar里如果解调部分写得比较粗糙误码率曲线可能比理论值差 2-3 dB。4.2 用 MATLAB 自带对象做交叉验证自己写的解调器对不对可以用comm.CPMDemodulator交叉验证。下面代码演示如何用系统对象做同样的解调。% 使用 MATLAB 自带 CPM 调制解调器 mod comm.CPMModulator(ModulationOrder, M, ... BitInput, true, ... FrequencyPulse, Rectangular, ... ModulationIndex, h, ... PulseLength, L, ... SymbolMapping, Binary); demod comm.CPMDemodulator(ModulationOrder, M, ... BitOutput, true, ... FrequencyPulse, Rectangular, ... ModulationIndex, h, ... PulseLength, L, ... SymbolMapping, Binary, ... TracebackDepth, 16); % 生成数据 data randi([0 1], 1000, 1); modSignal mod(data); % 加噪声 rxSignal awgn(modSignal, 10, measured); % 解调 demodData demod(rxSignal); % 计算误码率 [~, ber] biterr(data(1:length(demodData)), demodData); fprintf(误码率: %.4f\n, ber);comm.CPMModulator和comm.CPMDemodulator是 MATLAB 官方实现参数设置和你的自定义代码应该一致。如果两者误码率差很多说明你的自定义代码有问题。TracebackDepth是 Viterbi 的回溯深度一般取 5L 到 10L太小会损失性能太大增加延迟。FrequencyPulse可以选Rectangular、Raised Cosine、Gaussian等和你的g(t)对应。SymbolMapping选Binary或Gray二进制 CPM 一般用Binary。提示用awgn加噪声时measured选项会根据信号功率自动计算噪声功率但 CPM 是恒包络信号功率恒定所以也可以直接指定 SNR。如果误码率曲线在高 SNR 下出现地板检查相位同步和定时同步这两个是 CPM 解调的玄学问题。5. 避坑与排查CPM 仿真里最容易翻车的 5 个地方5.1 相位脉冲归一化错误导致误码率整体偏移现象误码率曲线形状正常但整体比理论值差 3-6 dB或者在高 SNR 下无法降到零。 原因g(t)的积分没有归一化到 1/2对于二进制 CPM导致实际调制指数偏离设定值。调制指数偏大或偏小都会让接收端的状态网格和实际信号不匹配。 解决在生成g(t)后强制g g / sum(g) * 0.5确保sum(g) 0.5。如果是多进制归一化到(M-1)/2或其他正确值。用trapz做数值积分验证。5.2 采样率不足导致频谱混叠现象频谱在高频端出现异常抬升或者误码率随sps增加而改善。 原因sps太小比如取 2 或 4相位累加的离散化误差大频谱旁瓣被混叠掩盖。 解决sps至少取 8推荐 16。对于 GMSK 这种带外衰减极快的sps要取 16 以上。可以用resample做上采样但最好在生成阶段就设够。5.3 Viterbi 回溯深度不够导致误码率地板现象低 SNR 时误码率正常高 SNR 时误码率不再下降出现地板。 原因TracebackDepth太小Viterbi 没有足够深度回溯到正确路径。 解决TracebackDepth设为5*L到10*L。对于 L3取 15 到 30。如果还不行检查状态网格是否完整有没有漏掉某些历史符号组合。5.4 定时同步偏差导致星座图旋转现象解调输出的星座图复基带出现旋转误码率对定时偏差敏感。 原因CPM 是恒包络但相位对定时非常敏感。采样点偏离最佳时刻会引入额外相位偏移。 解决在解调前做定时同步可以用早迟门或 Mueller-Muller 算法。MATLAB 里可以用comm.SymbolSynchronizer。如果只是仿真确保接收端和发送端的采样时钟完全一致或者用interp1做插值对齐。5.5 用错脉冲成形类型导致频谱不达标现象频谱旁瓣衰减比预期慢或者带外泄漏超标。 原因FrequencyPulse选错比如该用Gaussian却用了Rectangular或者 RC 的滚降系数beta设得太小。 解决根据频谱模板选脉冲。水声通信常用 RCbeta取 0.3-0.5GSM 用 GMSKBT取 0.3。在 MATLAB 里用pwelch看频谱和模板对比不达标就换脉冲或调参数。6. 进阶技巧用相位轨迹验证和加速仿真6.1 用相位轨迹图快速定位调制问题相位轨迹是 CPM 最直观的调试工具。把发送端和接收端解调后重建的相位轨迹画在一起如果两条线分叉说明解调出错。下面代码演示如何画对比图。% 画发送和接收相位轨迹对比 figure; subplot(2,1,1); plot(phase_tx, b); hold on; plot(phase_rx, r--); legend(发送相位, 接收相位); title(相位轨迹对比); xlabel(采样点); ylabel(相位 (rad)); grid on; subplot(2,1,2); plot(phase_tx - phase_rx); title(相位误差); xlabel(采样点); ylabel(误差 (rad)); grid on;如果相位误差在某个符号后突然变大检查那个符号对应的状态转移。常见问题是历史符号索引错位或者g(t)的截断位置不对。这个技巧比看误码率曲线快得多尤其适合调试部分响应 CPM。6.2 用查表法加速 ViterbiViterbi 的瓶颈是分支度量计算。如果每次都要重新生成预期波形仿真会非常慢。我一般会预计算所有状态转移的预期波形存成矩阵运行时直接查表。下面是一个示例。% 预计算分支度量表 num_trans num_states * M; branch_waveforms zeros(num_trans, sps); for s 1:num_states for m 0:M-1 idx (s-1)*M m 1; % 根据状态 s 和输入 m 生成预期波形 expected_phase states(s).phase 2*pi*h*(2*m-1); branch_waveforms(idx, :) exp(1j * expected_phase * ones(1, sps)); end end % 运行时直接查表 for k 2:N rx_seg rx((k-1)*sps1 : k*sps); for s 1:num_states for m 0:M-1 idx (s-1)*M m 1; metric sum(abs(rx_seg - branch_waveforms(idx, :)).^2); % ... 更新累积度量 end end end查表法能把仿真速度提高 5-10 倍尤其当sps较大时。注意branch_waveforms的生成要严格按相位递推公式不能简化。如果内存不够可以只存相位值运行时再算复指数但那样会慢一些。6.3 用 MATLAB 的parfor并行化误码率扫描误码率仿真通常要扫多个 SNR 点每个点跑很多帧。用parfor可以并行化。下面是一个框架。snr_range 0:2:12; ber zeros(size(snr_range)); parfor i 1:length(snr_range) snr snr_range(i); error_count 0; total_bits 0; for frame 1:100 data randi([0 1], 1000, 1); modSignal mod(data); rxSignal awgn(modSignal, snr, measured); demodData demod(rxSignal); [~, e] biterr(data(1:length(demodData)), demodData); error_count error_count e; total_bits total_bits length(demodData); end ber(i) error_count / total_bits; end semilogy(snr_range, ber, o-); xlabel(SNR (dB)); ylabel(BER); grid on;parfor要求循环体独立这里每个 SNR 点独立所以没问题。注意mod和demod是系统对象在parfor里每个 worker 会复制一份内存够就行。如果帧数很多可以把frame循环也并行化但那样通信开销大一般只并行 SNR 点。6.4 一个我常犯的错误我刚开始做 CPM 时总以为相位连续就是「相位不能跳」结果在写代码时把每个符号的相位增量直接累加忘了脉冲成形会跨符号。后来发现频谱怎么调都不对才回去检查g(t)的卷积。现在我的习惯是每写一个 CPM 脚本先画相位轨迹再画频谱最后跑误码率。三步都对了才敢说这个脚本能用。希望帮到你。本文还有配套的精品资源点击获取