ARTICLE DETAIL

资讯详情

深耕编程入门与网站建设的一线实战洞察。

基于决策树的数字调制识别MATLAB实现与实战

基于决策树的数字调制识别MATLAB实现与实战 简介本资源是一套面向通信工程专业本科生及信号处理初学者的MATLAB实践代码包聚焦数字调制样式自动识别这一典型通信信号分析任务。程序完整覆盖2ASK、4ASK、2PSK、4PSK、2FSK、4FSK和16QAM七类常见调制信号的仿真生成、加性高斯白噪声信道建模、瞬时幅度/频率/相位特征提取与分类识别全流程特别适合课程设计、毕设验证及算法原理理解。压缩包共14个文件10个.url为配套学习链接4个.m为主程序模块Digit_Modul.m为总控入口Feature.m负责特征计算channel.m模拟信道recognition.m执行判决识别总容量仅11KB轻量易部署。已有1034人学习下载代码结构清晰、模块职责分明附带可调门限机制便于用户基于实际数据统计优化识别阈值是深入掌握调制识别核心思想与MATLAB工程实现的理想参考范例。 有一次我在实验室里调一套SDR接收链路信号进来了界面上的星座图却在慢慢打转。旁边的师弟问我这到底是QPSK还是16QAM我盯着屏幕愣了几秒最后只能说“不知道”。从那时候起我就意识到真正的麻烦不是解调而是在不知道对方用什么格式发信号的情况下怎么把这个“不知道”变成“知道”。后来我用MATLAB把数字调制解调样式识别这套流程完整写成了程序源代码从信号生成、特征提取到分类决策全部打通。这篇内容就是把这套实现思路完整拆开信号怎么造、特征怎么算、阈值怎么定、代码长什么样以及实际接入信号后最容易坑人的几个地方。适合两类人看一类是刚接触软件无线电、需要快速给接收链路补上盲识别能力的学生另一类是工作中要处理不明信源、需要评估和落地自动调制识别AMC模块的工程师。1. 为什么接收机需要“不认识信号的自动识别”1.1 样式识别在整套通信链路里的位置传统的数字接收机工作模式是“先知道再解调”。发送端用什么调制方式、符号速率是多少、载波频率在哪、滚降系数取多少这些参数在接收机开始工作之前就已经写死在配置里了。接收机做的事情无非是下变频、匹配滤波、同步、判决然后把比特流吐出来。但实际场景往往不给你这个“先知道”。频谱监测设备扫到一个未知频点信号就在那里没有任何信令告诉你它是什么格式认知无线电要动态接入频谱也必须在短时间内判断当前频段里正在跑的是哪种调制方式甚至一些故障排查场景里你面对的是自家产线上下来的设备但配置丢了接收端完全不知道发送端设成了什么参数。在这些场景里调制样式识别就成了整套链路里绕不开的一环。它处在什么位置呢简单说在同步完成之后、正式解调判决之前。接收机先通过盲估计把载波频偏、符号定时这些基础参数抓出来然后对这个“已经初步同步好的信号”做样式识别判断出调制方式最后再用这个判断结果去配置真正的解调器。换句话说样式识别是一个“给解调器做确认”的模块。1.2 决策树特征工程和深度学习工程上到底怎么选现在做调制样式识别的路子大概分成两派。一派是经典特征参数决策树我这次的源码就是这条路。另一派是深度学习端到端识别把IQ采样直接喂给卷积网络或者循环网络让模型自己学特征。两派各有各的适用场景。我直接说结论如果你在实验室做原型验证、要写一版看得见摸得着、方便调试的代码经典决策树是最好上手的如果你手里有大量真实采集数据而且信道的恶劣程度超出常规模型假设那深度学习的上限确实更高。但深度学习有个很现实的问题——需要标注数据。调制样式识别里的标注本身就是重活你要先确认每条数据的真实调制格式这在很多场景下恰恰是最难解决的问题。没有干净标注深度模型就是空中楼阁。经典决策树的优势在于不需要训练数据、计算量小、单次判决只需要计算十几个统计量、而且每一步都可解释。这个“可解释”在工程上太重要了。识别错了你能顺着决策树看是哪一层判断出了问题是同步没做好还是特征阈值标定有问题。深度学习给一个概率分布出来你很难定位是哪个环节坏了。所以在工程落地时我的习惯是先用决策树跑通链路再按需引入更复杂的特征或模型。2. 底层设计预处理、特征提取与阈值逻辑2.1 预处理三板斧载波同步、符号定时、匹配滤波特征提取的前提是信号已经完成了基本的同步这一步没做好后面所有特征值都会失真。我在这套代码里做了三层预处理。第一层是载波同步。SDR采集下来的复基带信号往往带有残余频偏这个频偏会让星座图缓慢旋转也会让瞬时相位特征完全崩掉。我在预处理函数里先做一次粗频偏估计然后用Costas环或判决导向环做细同步。粗估计用周期图法对信号做FFT找频谱峰值细同步在仿真里用理想参数在真实采集数据上会切到判决导向锁相环。第二层是符号定时同步。采样点的位置不能落在两个符号的跳变沿上否则同一个符号采出来的幅度和相位都是错的。仿真里我直接用成型滤波器和过采样来保证符号中心对齐真实场景则用Gardner定时环。第三层是匹配滤波。发送端的成型滤波器通常用根升余弦接收端也要配备同样参数匹配滤波器才能把信噪比拉回最优。这个滤波器的滚降系数、抽头数都会影响后面特征值的数值分布尤其是包络类特征所以参数要固定下来不能每轮实验随便改。2.2 核心统计特征包络、相位、频率三个维度抓住格式差异特征提取是整个识别器的灵魂。我用的特征都是从瞬时幅度、瞬时相位、瞬时频率这三个维度推导出来的经典文献里叫Nandi-Azzouz特征集。下面逐个说清楚它们到底在捕捉什么。第一个特征是零中心归一化瞬时幅度谱密度最大值记作γmax。计算时先求信号的Hilbert包络a(n)然后做归一化a_cn(n) a(n)/m_a - 1其中m_a是整个包络的均值。对a_cn做FFT取幅度谱平方的最大值除以均值得到γmax。为什么这个特征能区分调制方式因为PSK信号的包络理论上恒定归一化之后a_cn非常接近零它的频谱就没什么像样的谱峰而ASK信号的包络随码元内容变化会存在明显的调制频率分量频谱上会出现突出的尖峰。所以γmax大的一侧是“包络有起伏”的ASK和FSK小的一侧是“包络恒定型”的PSK。第二个特征是零中心归一化瞬时幅度绝对值的标准差记作σaa。它用来在ASK家族内部定阶2ASK只存在两种幅度值4ASK存在四种归一化之后取绝对值两者的离散程度不一样。4ASK的幅度层次更多σaa会更大。第三个特征是零中心非线性相位标准差记作σdp。它的作用是区分BPSK和QPSK。做法是先对信号的瞬时相位做解卷绕去掉由载频产生的线性相位项剩下的就是调制相位变化。BPSK的相位跳变只有0和π两个值非线性相位在零中心后分布非常集中QPSK有0、±π/2、π四个相位值散布范围更大σdp自然更大。这里有个关键细节计算σdp前一定要剔除包络幅度过小的采样点。因为幅度接近零时噪声会把瞬时相位打得乱飞这些点的相位完全是噪声会直接污染统计量。第四个特征是零中心归一化瞬时频率标准差记作σaf。它用来区分FSK和其他调制格式。FSK的瞬时频率在不同码元之间跳变频率标准差天然很大而ASK和PSK的频率都集中在载频附近σaf较小。同理2FSK和4FSK之间也可以用σaf的大小继续细分。这四个特征单独看都有一定区分能力组合进决策树之后能把六种常见调制方式切干净。2.3 决策树结构先分大类再定阶数我用的决策树分两步走。第一步用γmax把信号分成“包络起伏类”和“包络恒定类”。“包络起伏类”包含ASK和FSK“包络恒定类”包含BPSK和QPSK。这一步是鲁棒性最强的分流因为γmax的计算不需要精确的相位信息对频偏和相位噪声的敏感度最低。第二步在大类内部继续细分。“包络恒定类”里用σdp区分BPSK和QPSK“包络起伏类”里用σaf先区分FSK和ASK然后用σaa给ASK定阶再用σaf的频率散布程度给FSK定阶。这个“先大类后细类”的设计不是随意的。如果一上来就用相位特征去分所有调制方式频偏稍微没消干净BPSK和QPSK的区分就全乱了。而γmax在频偏存在时依然稳定所以让最稳的特征打头阵把大类切对后面细分类的压力就小了。这也是工程上特别重要的一点决策树的排列顺序本质上就是按特征对信道损伤的鲁棒性排序。3. MATLAB源码实现从信号产生到分类决策全链路3.1 文件结构与主脚本整个工程我拆成了五个文件结构如下modulation_classifier/ ├── main_demo.m % 主脚本蒙特卡洛仿真 ├── generate_signal.m % 信号发生模块 ├── preprocess_rx.m % 接收预处理 ├── extract_features.m % 特征提取 └── classify_modulation.m % 分类决策主脚本的作用是循环仿真六种调制方式2ASK、4ASK、2FSK、4FSK、BPSK、QPSK在设定的信噪比下各跑几百次蒙特卡洛统计识别正确率。它同时负责调用其他四个函数把整个识别链路串起来。%% main_demo.m modTypes {BPSK, QPSK, 2ASK, 4ASK, 2FSK, 4FSK}; snrVec 0:2:20; Ntrials 500; N_symbols 1024; sps 8; % 每符号采样点数 fs 400e3; % 采样率 accMat zeros(length(modTypes), length(snrVec)); for m 1:length(modTypes) for s 1:length(snrVec) okCount 0; for trial 1:Ntrials x generate_signal(modTypes{m}, N_symbols, sps, fs); rx awgn(x, snrVec(s), measured); rxSync preprocess_rx(rx, fs, sps); feats extract_features(rxSync, fs); predType classify_modulation(feats, thresholds); if strcmp(predType, modTypes{m}) okCount okCount 1; end end accMat(m, s) okCount / Ntrials; end end这段代码里有一个隐藏点thresholds结构体要在主脚本里预先定义。我建议把阈值集中放在一个地方不要散落在各个函数里后面标定时改起来方便。3.2 信号发生模块确保仿真数据覆盖真实通信条件generate_signal函数负责生成六种调制信号。它的核心逻辑是按调制类型生成符号序列然后过成型滤波器并做过采样。这里有一个容易忽略的问题FSK信号的生成方式跟PSK/ASK不一样不能用upfirdn直接处理复数符号而是要在频率域上累加相位。%% generate_signal.m 片段 function x generate_signal(modType, N_symbols, sps, fs) switch modType case BPSK data randi([0 1], N_symbols, 1); symb pskmod(data, 2); case QPSK data randi([0 3], N_symbols, 1); symb pskmod(data, 4, pi/4); case 2ASK data randi([0 1], N_symbols, 1); symb data.; case 4ASK data randi([0 3], N_symbols, 1); symb (2*data - 3) / 3; case 2FSK data randi([0 1], N_symbols, 1); fDev 20e3; freq (2*data - 1) * fDev; freqUp upsample(freq, sps); freqUp conv(freqUp, ones(sps,1), same); phase 2*pi * cumsum(freqUp) / fs; symb exp(1j * phase); case 4FSK data randi([0 3], N_symbols, 1); fDev 10e3; freq (2*data - 3) * fDev; freqUp upsample(freq, sps); freqUp conv(freqUp, ones(sps,1), same); phase 2*pi * cumsum(freqUp) / fs; symb exp(1j * phase); end % 成型滤波 rrc rcosdesign(0.35, 6, sps); x upfirdn(symb, rrc, sps); x x(1:N_symbols*sps); end注意这里的upsample加上conv的做法等于是对频率序列做了一个零阶保持保证FSK在一个符号周期内频率恒定这样相位累积是线性的。滚降因子我取了0.35这个数值要跟接收端的匹配滤波器保持一致。3.3 特征提取核心函数代码和公式一一对应extract_features是整个程序里最核心的函数。它的实现跟公式一一对应最好逐行读。%% extract_features.m function feats extract_features(x, fs) N length(x); analytic hilbert(x); % 解析信号 env abs(analytic); % 瞬时包络 m_a mean(env); a_cn env ./ m_a - 1; % 零中心归一化瞬时幅度 % --- 特征1: gamma_max --- spec abs(fft(a_cn)).^2; spec spec(2:floor(N/2)1); % 去掉直流 gamma_max max(spec) / (mean(spec) eps); % --- 相位特征 --- phase unwrap(angle(analytic)); valid env 0.9 * m_a; % 剔除低幅度点 phi_nl phase(valid); phi_nl phi_nl - mean(phi_nl); % 零中心化 sigma_dp std(phi_nl); % --- 瞬时频率特征 --- f_inst diff(phi_nl) * fs / (2*pi); m_f mean(f_inst); sigma_af std(f_inst ./ (m_f eps)); % --- sigma_aa --- a_v a_cn(valid); sigma_aa sqrt(mean(a_v.^2) - (mean(abs(a_v)))^2); feats [gamma_max, sigma_dp, sigma_af, sigma_aa]; end三个关键细节值得展开。第一计算γmax时FFT前要去掉直流。因为a_cn的直流分量直接对应包络均值如果不去掉频谱里会出现一个巨大的零频峰这个峰没有任何区分度还会把整个频谱的平均值抬高导致γmax被压低。第二相位解卷绕用unwrap但并不是所有采样点都能保留。幅度小于0.9倍平均包络的点相位完全被噪声主导必须剔掉。这个门限我试过0.7、0.8、0.9最后0.9在低信噪比下表现最稳。门限太低会把噪声相位放进来门限太高会丢掉大量有效符号导致统计量方差变大。第三σaf计算的是瞬时频率的归一化标准差。这里用diff求差分等于做了频率解调。分母上的m_f是平均频率这个平均值在理想情况下应该接近零频偏。如果残余频偏没消干净m_f会变大σaf会被压小这就会影响FSK和ASK的区分。这就是为什么预处理里频偏消除必须做得足够干净。3.4 分类决策与阈值标定classify_modulation函数按决策树结构逐层判断。%% classify_modulation.m function modType classify_modulation(feats, thr) gamma_max feats(1); sigma_dp feats(2); sigma_af feats(3); sigma_aa feats(4); if gamma_max thr.gamma_ps if sigma_af thr.sigma_af_askfsk if sigma_af thr.sigma_af_fsk_order modType 4FSK; else modType 2FSK; end else if sigma_aa thr.sigma_aa_ask_order modType 4ASK; else modType 2ASK; end end else if sigma_dp thr.sigma_dp_psk_order modType QPSK; else modType BPSK; end end end阈值不是拍脑袋定的。我在这套仿真参数下符号数1024、每符号8个采样点、SNR15dB实测标定出来的参考阈值如下特征阈值字段参考值作用γmaxgamma_ps2.5区分PSK与ASK/FSKσdpsigma_dp_psk_order0.5区分BPSK与QPSKσafsigma_af_askfsk0.35区分FSK与ASKσafsigma_af_fsk_order0.6区分2FSK与4FSKσaasigma_aa_ask_order0.3区分2ASK与4ASK必须强调这些阈值是相对值不是绝对值。换一套采样率、符号数、滚降系数甚至换一个信噪比工作区间阈值都要重新标定。标定的方法很简单跑一遍已知标签的仿真数据把每个特征值和真实标签打印出来画出箱线图在两类分布的中间位置取阈值。我在工程里就是这么干的不依赖任何理论公式硬算因为理论推导很难覆盖滤波器实现、定时误差、噪声模型等实际损失。4. 实测踩坑频偏、定时与低信噪比下的特征失真4.1 频偏没消干净QPSK被识别成了2FSK我第一次跑完整个程序发现识别率在SNR高于10dB的时候也不是100%QPSK偶尔会被判成2FSK。这个结果非常反直觉因为这两个调制方式在星座图上差异巨大怎么会被混淆我把出错的样本拉出来分析先看星座图发现星座点不是聚成四个点而是在画圆——这是典型残余频偏的现象。再看瞬时相位曲线整条曲线带着一个线性上升的斜坡。问题一下就清楚了残余频偏让瞬时相位多了一项线性项σdp被撑大了同时这个线性相位映射到瞬时频率上又产生了非零的频率分量σaf也偏大。两个特征同时失真决策树就把QPSK推进了FSK的分类支路。排查链路是这样的先画出错样本的星座图和相位曲线确认是频偏问题再在仿真里人为把频偏从0慢慢增加到几个Hz观察特征值变化趋势最后回到预处理函数把细同步环路的收敛精度提升。修完之后QPSK的σaf明显回落分类就恢复了。这个坑给我最大的教训是特征提取是否可靠完全取决于预处理是否够狠。频偏残余哪怕只有一个符号周期的百分之几对相位类特征都是致命的。4.2 符号定时偏差让4ASK变成了FSK另一个坑出在定时同步上。仿真里我本来用upfirdn做的过采样每个符号固定采8个点按理说不该有定时问题。但后来我为了让仿真更接近真实情况在发射端加了一个随机时延结果4ASK的识别率暴跌。看了一眼误判矩阵4ASK大量被判成了2FSK。为什么随机时延导致采样点落在符号跳变沿包络在跳变沿出现突然的凹坑和尖峰。这些包络瞬变会被Hilbert变换捕捉到在频谱上产生额外的谱峰γmax变大这本身没问题因为ASK本来就属于γmax大的那一类。但问题是这个随机扰动也被带进了瞬时频率的计算里σaf被异常增大于是4ASK被推给了FSK支路。解决办法是在预处理函数里加上定时同步。仿真里最简单可靠的方法是在每个符号周期内搜索包络最大值的位置以该位置作为符号中心重新采样。这个办法对PSK和ASK都有效FSK因为包络恒定需要改用Gardner定时环。把定时同步加进去之后这个误判基本消失了。4.3 低信噪比下的阈值漂移固定阈值会失效最后一个坑是特征阈值本身会随信噪比漂移。我一开始把所有阈值固定在15dB下标定的值然后把SNR从15dB降到0dB跑完整套仿真识别率掉得很难看。画了各特征值在不同SNR下的箱线图之后发现γmax在低SNR下会整体抬升。原因是噪声让所有信号的包络都产生了随机起伏原本包络恒定的PSK信号在低幅度采样点上也被噪声“掰出”了幅度变化。这样一来PSK类的γmax分布和ASK/FSK类的γmax分布开始重叠固定阈值自然就切不开了。解决思路有三个层次。最低成本的方案是按SNR分段设阈值把工作区间分成若干个SNR档位每档用一组标定好的阈值。这个方案简单粗暴但要保证系统能比较准确地估计当前SNR否则阈值切换时机反而是新问题。第二个方案是对包络做平滑降低噪声对瞬时幅度的扰动但平滑会造成符号间串扰需要权衡。第三个方案是换更鲁棒的特征比如高阶累积量或者循环谱特征这些特征对加性白噪声的敏感性要低得多。我在这套代码里最终采用了“分段阈值平滑后置”的组合方案低信噪比识别率提升了十几个百分点。5. 往下走高阶调制、深度学习与设备落地5.1 识别范围扩充高阶累积量的思路如果需要识别的调制方式有16QAM、64QAM、8PSK这些前面四个特征就有点不够用了。16QAM和64QAM的包络都有起伏但靠σaa很难稳定区分因为两者的幅度层次很多归一化后统计特征重叠严重。这个场景我建议引入高阶累积量特征。高阶累积量对高斯噪声天然免疫而且不同调制格式的累积量理论值是分离的。比如对复基带信号C40、C41、C42的组合可以比较好地区分BPSK、QPSK、8PSK、16QAM、64QAM。实际计算也不复杂就是对信号的几个矩做组合运算MATLAB里用moment函数加上自定义公式就可以实现。加入高阶累积量之后决策树的前两级可以保持不变只是在细分类支路里增加判断节点。5.2 深度学习端到端路线什么时候值得换深度学习的优势我之前说过要在数据充足、信道复杂的情况下才值得换。具体到调制样式识别我见过效果不错的做法是把IQ两路采样拼成一个2×N的矩阵当作“图像”输入CNN输出层用softmax做分类。这种方法的识别精度通常比经典决策树高尤其是在多径衰落信道下。但工程上有个很实际的成本训练数据需要覆盖足够的信道场景否则模型在真实数据上的泛化能力很差。仿真数据训练出来的模型换到真实SDR环境通常会掉点。我的建议是如果场景相对固定、需要快速上线经典决策树足够如果要做长期的平台能力建设可以两条腿走路决策树做保底深度模型做增强。5.3 从MATLAB到SDR设备实时落地的几个关键点仿真代码跑通之后真正的问题才刚开始怎么搬到实时的SDR设备上跑。MATLAB提供了Coder工具箱可以把特征提取和分类函数转成C代码但有几个地方要提前优化。第一hilbert函数在Coder里支持有限最好自己用FFT实现解析信号构造。第二FFT长度不需要做完整长度加窗截断到1024点或2048点足够能省不少计算时间。第三浮点转定点时γmax和σdp这些统计量的动态范围比较大要对每一级运算单独做定点范围标定否则截断误差会在低SNR时放大。第四如果目标平台是嵌入式处理器可以考虑把分类器换成查表法把特征空间离散化之后用查找表直接输出结果这样分类部分零运算量。最后说一个我自己的调试习惯这个习惯帮我省了非常多时间不管仿真还是实测我都会把每个样本的特征值、真实标签、当前SNR一起打印出来或者画成散点图。特征分布一旦可视化阈值怎么定、哪个特征在什么SNR下失效全都一目了然。很多同学跑来问我“为什么识别率上不去”我第一反应永远是先把特征分布图发我看看。大多数时候问题根本不在分类器而在特征本身已经糊成一团了。本文还有配套的精品资源点击获取
返回列表