简介：一份MATLAB源码包，用于GSM移动通信系统的全链路仿真，适合通信工程专业学生、科研人员及数字蜂窝系统初学者。程序覆盖GSM核心物理层流程：信道卷积编码、GMSK调制、FDMA/TDMA多址接入、AWGN与瑞利衰落信道模拟、维特比解调与解码，以及误码率和符号错误率评估；运行demo脚本还可看到星座图、频谱图等可视化结果，帮助把理论公式转化为可操作的仿真实验。其中卷积编码、交织、突发脉冲生成、解交织等模块均可独立复用，配合调制解调函数灵活组合不同传输条件。包内共31个文件，以28个.m功能脚本为主，另有3个.sim仿真模型文件，资源包仅30KB，代码量精炼、结构清晰，便于按模块阅读和二次修改，可用于参数调节与算法验证。目前已有470人学习参考，适合课程设计、毕业设计或自学入门。
1. 一个能跑通的GSM物理层链路：从编解码到误码率曲线
用MATLAB写GSM仿真程序，最尴尬的不是不会写GMSK调制，而是写了半天发现整条链路对不上：发射端的突发脉冲格式错了，接收端解调的参考相位没对准，交织器解交织的顺序反了。你盯着误码率曲线发呆，根本不知道是哪一环出了问题。这份资源解决的正是这个问题——它把GSM物理层的一条完整收发链路拆成可运行的MATLAB模块，从信道编码、交织、突发脉冲格式化、GMSK调制，到多径信道建模和接收端解调，每一步都有独立脚本和可验证的中间结果。
这套程序特别适合两类人：一是通信方向正在做课程设计或毕业设计的学生，需要把GSM的物理层流程真正跑起来出BER曲线；二是刚接触无线通信仿真、想理解GMSK调制和均衡器在实际代码里长什么样的开发者。程序用MATLAB脚本实现，不需要Simulink，也不需要额外的通信工具箱之外的复杂依赖，改参数就能换场景，比如把静态信道换成典型城区多径信道，直接观察误码率恶化情况。下文按“调制发射—成帧编码—信道传输—接收解调”的顺序拆开讲，最后给出一份能直接改的二次开发路线。
2. 核心模块拆解：GMSK调制与突发脉冲的生成链路
2.1 为什么GSM选GMSK而不是QPSK
GSM选GMSK的根本原因不是频谱效率高，而是它包络恒定，对功率放大器的线性度要求低。早期的手机功放基本都工作在饱和区，如果信号包络有起伏，过零点的幅度会被削掉，误码率直接恶化。GMSK的本质是MSK（最小频移键控）前面加一个高斯低通滤波器，把相位路径的突变磨平，让频谱旁瓣快速衰减，从而满足GSM的邻道干扰指标。
在MATLAB里实现GMSK，常见做法是用comm.GMSKModulator对象，但更清晰的算法教学实现是手动完成高斯滤波加调频。资源里的实现路径如下：先把二进制比特映射为+1/-1的NRZ序列，过高斯滤波器得到平滑后的基带波形，再对波形积分得到相位，最后用相位生成I/Q两路信号。这个过程让你能随时抽看中间信号，比黑匣子调制器更直观。
% gmsk_modulate.m —— 手动实现GMSK调制 function [I, Q, phase] = gmsk_modulate(bits, Bt, sps, over) % bits : 输入比特流（0/1） % Bt : 高斯滤波器带宽与符号速率之比，GSM标准为0.3 % sps : 每个符号的采样点数，典型值8 % over : 过采样倍数，用于成形滤波器 nrz = 2*bits - 1; % 0映射为-1，1映射为+1 up = upsample(nrz, sps); % 插值为采样域 up = conv(up, ones(1, sps)/sps, 'same'); % 简易矩形保持 % 高斯滤波器时域响应 sigma = sqrt(log(2)) / (2*pi*Bt); t = -3:1/sps:3; h = (sqrt(2*pi)/sigma) * exp(-2*pi^2*sigma^2*t.^2); h = h / sum(h); % 归一化 filtered = conv(up, h, 'same'); % 高斯滤波 phase = cumsum(filtered) * pi/2 / sps; % 相位累加，MSK频偏 I = cos(phase); Q = sin(phase); end这段代码的核心在于phase = cumsum(filtered) * pi/2 / sps这一行。GMSK每个符号的相位变化量是±π/2，但经过高斯滤波后相位路径是平滑过渡的，累积相位就是调制信号的瞬时相位。BT=0.3是GSM标准参数，代表高斯滤波器3dB带宽和符号速率的比值，这个值越小频谱越窄，但符号间干扰越严重，0.3是工程上的折中值。
实际仿真中你会发现，如果在接收端直接用comm.GMSKDemodulator解调，误码率不容易做低，原因是GMSK引入的码间串扰在高斯滤波后已经存在。所以完整的GSM接收机通常要加均衡器，最常用的是基于Viterbi的MLSE均衡，这一点在第4章展开讲。
2.2 突发脉冲格式化：把数据装进GSM的时隙容器
GSM的物理层不是裸比特流直接上空中，而是以突发（burst）为单位组织数据。一个普通突发（Normal Burst，NB）包含3个起始尾比特、57+57个数据比特、26个训练序列比特、3个结束尾比特和8.25个保护间隔比特，共148比特加保护段。训练序列用于接收端的信道估计和均衡器参考，这就是GSM能对抗多径的关键设计之一。
资源里的build_burst.m把编码后的数据比特按照GSM的时隙格式组装，过程包括：插入训练序列、加尾比特和保护间隔，最后按突发长度输出。这里有个关键点：训练序列在发射端是已知的，接收端利用收到的训练序列受信道影响后的失真来估计信道冲激响应，再反向补偿数据段。
% build_burst.m —— 构造GSM普通突发 function burst = build_burst(data_bits, train_seq) % data_bits : 114位编码后数据（两段各57位） % train_seq : 26位训练序列，GSM标准定义8种可选 tail = [0 0 0]; guard = [0 0 0 0 0 0 0 0.25]; % 8.25位保护间隔，末尾0.25位是虚拟 burst = [tail, data_bits(1:57), train_seq, data_bits(58:114), tail, guard]; end这个函数是纯拼接逻辑，但容易踩坑的地方在于：数据比特的排列顺序和交织器的映射关系必须一致。GSM的语音业务里，57+57的两段数据分别来自两个不同交织块，接收端解交织时必须知道哪段数据对应哪个交织块。资源里在interleave.m里专门维护了一个置换表，这个表与GSM规范中的交织矩阵严格对应，改错一个索引后面全部错位。
需要注意的是，保护间隔（guard period）的8.25位在发射端并不是真的发射信号，而是留空——它用来吸收不同手机到达基站的时延差，避免相邻时隙碰撞。仿真里如果不加这个空闲位，接收端同步时会把相邻burst的数据误当作本burst延长段，导致Viterbi均衡器的路径度量计算错误。
3. 信道编码与交织：语音数据如何被保护
3.1 卷积码加CRC：GSM语音信道的级联编码结构
GSM的全速率语音信道采用级联编码：先对语音帧做CRC校验，再做卷积编码。CRC用于检测严重错误，卷积码用于纠错。仿真的完整流程是：语音帧（50个重要比特+132个次要比特）经过CRC、卷积码（码率1/2，约束长度5）、重排和交织后，形成456比特的编码块，再拆进4个突发里发送。
资源里conv_encoder.m实现的卷积码生成多项式为G0=1011和G1=1101（八进制为133和171），约束长度K=5。这是GSM规范里定义的经典多项式，删余（puncturing）后会得到不同码率，但基础仿真里常用未删余的全速率编码。
% conv_encoder.m —— GSM语音信道1/2码率卷积编码 function out = conv_encoder(bits, poly0, poly1) % poly0, poly1 : 生成多项式，以二进制向量表示 % 例如 poly0 = [1 0 1 1], poly1 = [1 1 0 1] 对应GSM标准 reg = zeros(1, 4); % K=5，寄存器长度4 out = zeros(1, 2*length(bits)); for n = 1:length(bits) reg = [bits(n), reg(1:end-1)]; x = mod(sum(reg .* poly0(2:end)), 2); % 监督位1 y = mod(sum(reg .* poly1(2:end)), 2); % 监督位2 out(2*n-1) = x; out(2*n) = y; end end这段实现里最容易被忽视的是reg = [bits(n), reg(1:end-1)]这一行——它是把新比特移入寄存器时，寄存器从高位到低位的顺序必须和生成多项式从左到右的系数一致，否则编码出的监督位完全错位，且看起来很像实现正确。拿标准的GSM测试序列跑一遍，用接收端译码复核就能暴露问题。
卷积编码后是交织。交织的目的不是纠错，而是把连续比特错误打散到不同的突发中。GSM的语音帧交织深度为8个突发，一个456比特的编码块分散在8个连续burst里，使得一个突发被深度衰落打掉时，单块内的错误比特数不超过16个，落在卷积码纠错能力范围内。interleave.m用一个456×8的映射矩阵实现吞吐顺序和交织顺序的变换，这个矩阵来自GSM规范中的公式，直接查表映射即可。
3.2 Viterbi译码：硬判决与软判决的差异
接收端译码用的是Viterbi算法，这是卷积码最大似然译码的标准做法。资源里viterbi_decoder.m的函数实现包括：从接收符号计算分支度量、按蝶形结构做ACS（加-比-选）运算、最终回溯得到译码比特。
% viterbi_decoder.m —— 约束长度5的Viterbi译码 function decoded = viterbi_decoder(recv, poly0, poly1) % recv : 接收端符号序列（如果做软判决，传入软信息） states = 16; % 2^(K-1) = 16种状态 path_metric = Inf(1, states); path_metric(1) = 0; % 从全零状态开始 survivor = zeros(length(recv)/2, states); for n = 1:length(recv)/2 % 提取当前两个监督比特 r0 = recv(2*n-1); r1 = recv(2*n); new_metric = Inf(1, states); for s = 0:states-1 % 每个状态有两条入边，分别对应输入0和1 for inp = 0:1 prev = floor(s/2) + inp * 8; % 回溯前序状态 % 计算理论输出并比较距离 exp0 = mod(sum(dec2bin(prev,4) .* poly0(2:end)), 2); exp1 = mod(sum(dec2bin(prev,4) .* poly1(2:end)), 2); dist = (r0 - (2*exp0-1))^2 + (r1 - (2*exp1-1))^2; if path_metric(prev+1) + dist < new_metric(s+1) new_metric(s+1) = path_metric(prev+1) + dist; survivor(n, s+1) = prev; end end end path_metric = new_metric; end % 回溯省略：取最终最小度量状态并逆推 end这里需要提醒一点：硬判决在输入r0、r1时通常用±1的实数值代表比特1和0，距离度量用欧氏距离平方；如果信道是AWGN且做了软解调，可以把接收信号直接当软信息输入，增益在0.5～2dB之间。这段代码性能不会特别好，原因是遍历状态的方式没有用蝶形结构，但可读性很好，适合理解原理。工程里常见的做法是用comm.ConvolutionalDecoder替代手写实现，但手写代码的价值在于你可以加上traceback length、tail biting这些细节，这是调误码率曲线时的关键参数。
4. 信道模型与接收端同步：从高斯噪声到典型城区多径
4.1 AWGN与瑞利衰落信道的仿真差异
GSM系统在实际部署中面对的主要挑战是多径衰落。仿真的第一步通常从AWGN信道开始，用于验证编码和调制的正确性；第二步换成频率选择性衰落信道，观察均衡器带来的增益。资源里channel_model.m同时实现了这两种信道模式，用参数切换。
% channel_model.m —— AWGN与多径信道二选一 function rx = channel_model(tx, snr_dB, channel_type) % tx : 发射端I/Q信号 % snr_dB : 信噪比（符号级） % channel_type : 'awgn' 或 'tu50' if strcmp(channel_type, 'awgn') noise = sqrt(10^(-snr_dB/10)) * randn(size(tx)); rx = tx + noise; elseif strcmp(channel_type, 'tu50') % 典型城区信道，50km/h移动速度 % 抽头延迟线：6条径，相对时延和平均功率按GSM标准 delays = [0 0.2 0.4 0.6 0.8 1.2]*1e-6; % 微秒级时延 powers = [0 -3 -6 -8 -10 -14]; % 相对功率dB % 每个抽头加独立瑞利衰落 rx = zeros(size(tx)); for k = 1:length(delays) h = rayleighchan(1/1e6, 50, delays(k), powers(k)); rx = rx + filter(h, tx); end % 再叠加AWGN noise = sqrt(10^(-snr_dB/10)) * randn(size(tx)); rx = rx + noise; end end这里rayleighchan是MATLAB通信工具箱的对象，在较新版本中推荐使用comm.RayleighChannel。但资源里保留老接口是为了兼容性——很多教学环境中MATLAB版本较旧，换成新接口反而会报错。注意TU50信道的多径时延扩展接近5微秒，而GSM符号周期约3.69微秒，这意味着延迟径的能量会跨越接近1.5个符号周期。符号间干扰是必然的，这正好引出了均衡器的必要性。
4.2 训练序列做信道估计：MLSE均衡的前置步骤
GSM接收端的均衡器基于最大似然序列估计（MLSE），本质是Viterbi算法在信道响应上的应用。它把信道的码间串扰当成一个有限状态机，用训练序列估计信道冲激响应后，对所有可能的发送序列做最大似然搜索。这个搜索的复杂度是星座点数量乘以信道记忆长度的指数，但GMSK的相位状态加上ISI状态，实际实现会做简化。
% mlse_equalizer.m —— 基于Viterbi的MLSE均衡 function est_bits = mlse_equalizer(rx, train_seq, chan_est) % rx : 接收数据段 % train_seq : 本地训练序列 % chan_est : 估计出的信道冲激响应 % 构建卷积状态：状态数 = 2^(chan_len-1) chan_len = length(chan_est); states = 2^(chan_len-1); path_metric = inf(1, states); path_metric(1) = 0; for n = 1:length(rx) new_metric = inf(1, states); for s = 0:states-1 for inp = 0:1 % 构造该状态对应的发送符号序列 seq = [inp, dec2bin(floor(s/2), chan_len-1) == '1']; % 计算理想接收值 ideal = seq * chan_est.'; dist = abs(rx(n) - ideal)^2; prev = ...; % 回溯前一状态 if path_metric(prev+1) + dist < new_metric(s+1) new_metric(s+1) = path_metric(prev+1) + dist; end end end path_metric = new_metric; end % 回溯得到估计比特 endMLSE均衡的代价是复杂度随信道长度指数上升。GSM规范中训练序列长度为26位，信道估计窗长通常取5～7个抽头，状态数为16～64，这个规模在MATLAB里仿真可接受。工程手机里为了省电，通常用DFSE（判决反馈序列估计）做近似，牺牲少量性能换取功耗降低。仿真改成DFSE的做法是：保留前两个抽头的Viterbi搜索，后面的抽头用硬判决反馈抵消。
4.3 频偏估计与补偿：仿真的隐藏难点
GSM仿真里最容易被忽略的环节是频偏。发射端晶振偏差、接收端本振偏移，在真实系统里会有几百赫兹的载波频偏。GMSK信号对频偏非常敏感，1kHz的频偏就能让误码率从1e-4恶化到0.1以上。资源里的freq_offset_est.m利用训练序列的前后对称性计算相位旋转，实现简单的频偏估计。
% freq_offset_est.m —— 基于训练序列的频偏估计 function f_est = freq_offset_est(rx_train, train_seq, sps, fs) % rx_train : 接收端提取的训练序列段 % train_seq : 本地训练序列 % sps : 每符号采样数 % fs : 采样率 % 方法：用本地序列的调制结果与接收信号共轭相乘，提取相位斜率 local = gmsk_modulate(train_seq, 0.3, sps, 4); z = rx_train .* conj(local); phase_diff = unwrap(angle(z)); p = polyfit((0:length(phase_diff)-1)/fs, phase_diff, 1); f_est = p(1) / (2*pi); end这段代码的原理是：如果接收信号带频偏，那么共轭相乘后的相位会随时间线性增长，线性拟合的斜率就是频偏的角频率。unwrap是必须的，因为相位有2π跳变，直接diff会得到错误的频偏值。用此估计值对接收数据做乘上exp(-j*2*pi*f_est*t)即可补偿。仿真时不加这一步，BER曲线在高信噪比区域会掉不下去，读者拿到代码后第一个要确认的就是有没有做频偏补偿。
5. 常见问题与排查：跑仿真时最常踩的四个坑
5.1 误码率曲线下不去，平台期在1e-2附近
这是仿真中最常见的现象：BER曲线在低信噪比时正常下降，但在高信噪比区域卡住，不再下降。原因几乎都是残留频偏或残留符号定时误差，而不是信道编码的问题。解决方法是检查频偏估计后的残余量，或者直接在接收端加入一个简单的相位跟踪环。我一般会在频偏补偿后打印残余相位的标准差，如果大于0.05弧度，说明估计精度不够，需要提高训练序列段的长度或者改用二阶环路。
5.2 训练序列发对了，信道估计结果却是错的
训练序列的位置和接收端提取的位置必须严格对齐。GSM突发里的训练序列前后各有一个57位的数据段，如果需要做信道估计，必须先用滑动相关找到训练序列的准确起点，再做LS信道估计。如果直接用已知位置切割，一旦定时偏差超过半个符号，估计出的信道冲激响应会混入相邻符号的能量。解决方法是先做定时同步，用本地训练序列和接收信号做互相关峰值定位。
5.3 交织器和解交织器的索引对不上
这个问题非常隐蔽：单独跑交织器时，输出看起来完全随机，无法验证正确性。如果发送端用随机比特做交织，接收端解交织后能恢复原序列，说明索引一致。但很多人的实现是照抄公式，公式里矩阵下标从0开始还是从1开始弄反了，导致奇偶位置互换。解决方法是构造一个固定序列如1:456，过一遍交织再解交织，输出与输入完全一致才算通过。这个自检过程应该放在任何调制解调之前。
5.4 采样率设置错误导致频谱折叠
GMSK信号的带宽约为200kHz，仿真中的采样率至少要400kHz，工程里通常取1.6MHz或更高。新手常犯的错误是采样率刚好等于符号速率的2倍，结果频谱混叠严重，误码率莫名其妙地高。资源里的代码建议sps设为8或以上，这样既能看清调制波形，又能保证卷积滤波器的时域分辨率。改采样率时要注意，所有滤波器的截止频率、频偏估计都隐含依赖fs，不能只改一个地方。
6. 二次开发思路：把静态仿真改成可配置的实验平台
6.1 参数化入口让不同实验一键切换
这份资源的核心改进空间在工程化易用性上。当前每个函数都是独立脚本，把它们统一收敛到一个顶层入口文件run_gsm_sim.m里，用结构体配置参数，这样换信道类型、改信噪比区间、切换编解码方式都只需要改一行。实用做法是定义params结构体，包含channel_type、snr_range、enable_freq_offset、eq_type等字段，主循环按参数分派到不同的处理函数。
% run_gsm_sim.m —— 统一仿真入口 params.channel_type = 'tu50'; % 'awgn' 或 'tu50' params.snr_db = 0:2:20; % 扫描信噪比 params.freq_offset_hz = 100; % 注入频偏 params.eq_mode = 'mlse'; % 'mlse' 或 'dfse' results.ber = zeros(size(params.snr_db)); for k = 1:length(params.snr_db) % 完整链路调用 [tx_bits, tx_signal] = gsm_transmit(data_bits, params); rx_signal = channel_model(tx_signal, params.snr_db(k), params.channel_type); [rx_bits, metric] = gsm_receive(rx_signal, params); results.ber(k) = sum(rx_bits ~= data_bits) / length(data_bits); results.metric(k) = metric; end把链路封装成gsm_transmit和gsm_receive两个函数后，你可以用同一套流程跑不同配置，比如对比硬判决与软判决Viterbi的增益、对比MLSE与DFSE的复杂度和性能差异、查看频偏估计前后的星座图变化。这样改造之后，这套仿真程序就不再是死代码，而是一个能支撑多种实验的框架。
6.2 扩展方向：跳频、设备指纹与上行链路
GSM物理层仿真还可以继续往上叠功能。跳频是GSM的重要抗干扰手段，实现思路是多个突发在不同频点上发射，但频谱结构相同，仿真里只需在突发间切换载波频率，无需改动调制函数。设备指纹识别则是把发射端的相位误差和幅度误差当作“指纹”，通过比较接收信号和理想信号的偏差来区分设备，这在GSM安全研究里很有价值。上行链路的仿真与下行存在差异：上行需要考虑不同用户的定时提前和功率控制，但核心的GMSK调制、均衡和译码模块完全复用。如果你把这份资源和实际的GSM空口抓包工具配合，还能离线验证信道编码的比特级正确性，这是改代码时最扎实的参考依据。
6.3 一处值得替换的性能瓶颈
我个人的建议是优先把mlse_equalizer.m里的全网格Viterbi换成蝶形结构实现。当前版本函数可读性好，但复杂度是O(N*2^L)且每一层遍历全部状态时有很多冗余乘加操作。改成蝶形结构后，状态转移从“每个状态遍历所有入边”变成“每两个状态配对更新”，复杂度降一半，且代码性能更接近工程实现。完成后用同一组信道参数对比BER曲线，如果改动正确，曲线应该完全重合。这个替换本身也是一个很好的算法练习题，能帮你彻底理解Viterbi的几何结构。从那以后，我每次在仿真平台里测试新改动，都会先跑一遍1:456序列的自检交织映射，强制确认索引一致，再谈误码率。希望帮到你。
本文还有配套的精品资源，点击获取
