简介面向无线通信与认知无线电方向的学习者和研究人员这份MATLAB仿真资源聚焦动态频谱接入中的频谱感知与接入策略设计针对频谱资源紧张与利用不足的矛盾提供了一套可直接运行的仿真基线可用于理解能量检测、匹配滤波等常见感知方法并验证不同参数条件下的系统接入性能。压缩包共2个文件main.m为仿真主程序README.md对运行逻辑与模块划分做了说明整体仅6KB轻量易用适合快速修改与调试。目前已有26人学习浏览。借助该代码读者可快速复现一个包含主用户与认知用户的动态频谱接入场景观察感知周期、检测阈值、信号调制方式等参数对频谱感知准确性和通信性能的影响并以此为基线扩展算法或增加信道干扰模型是快速上手认知无线电仿真实验的实用参考。1. 认知无线电动态频谱接入这套 MATLAB 仿真到底能复现什么认知无线电Cognitive Radio, CR里动态频谱接入Dynamic Spectrum Access, DSA是最常被拿来开题、做课程设计和毕设的方向但多数人卡在同一个地方论文里的检测概率、虚警概率、吞吐量曲线用 MATLAB 到底怎么一步步跑出来。这套仿真系统把整条链路拆成四个可独立运行的模块——频谱感知、主用户活动建模、接入决策、性能统计跑一轮蒙特卡洛能同时输出检测概率、虚警概率、次用户有效接入率和碰撞率。它给的不是摆好的成品曲线而是可改的感知算法、可调的主用户模型参数和完整事件流程适合两类人赶毕设或课程设计、需要快速拿到可复现结果的本科生以及刚接触 CR、想弄懂能量检测阈值推导和 ON-OFF 主用户模型关系的入门研究者。下面按模块逐层拆每个文件管什么、参数怎么改、坑在哪里一次讲透。2. 频谱感知层能量检测的 MATLAB 实现与阈值推导频谱感知是整个 DSA 闭环的第一环感知结果错了后面的接入决策全是空中楼阁。这一章先把能量检测的数学模型讲清楚再给可运行的蒙特卡洛仿真最后用参数表说明阈值怎么选。2.1 能量检测模型先从判决统计量说起认知无线电频谱感知里能量检测是复杂度最低、也最常被当作基准的算法。思路很直接主用户存在H1时接收信号是信号叠加噪声平均功率比纯噪声H0大拿接收信号能量和门限比一比就能判决。假设接收采样序列为 x[n]n 1, 2, ..., N能量检测器的判决统计量定义为T (1/N) * sum(|x[n]|^2)N 足够大时由中心极限定理T 在 H0 下近似服从高斯分布 N(σw², 2σw⁴/N)在 H1 下近似服从 N(σw² σs², 2(σw² σs²)²/N)。σw² 是噪声方差σs² 是主用户信号在接收端的功率。判决规则就是 T 大于门限 λ 判 H1占用否则判 H0空闲。这套系统里能量检测核心函数如下function [decision, T] energy_detector(x, N, lambda) % 能量检测判决函数 % x : 接收采样序列长度 N % N : 观测窗口内采样点数 % lambda : 判决门限 T sum(abs(x).^2) / N; % 判决统计量归一化平均能量 decision T lambda; % 1 判 H1主用户占用0 判 H0 end逻辑说明函数先算平均能量再和门限比较。T 用平均能量而不是总能量目的是把门限和采样点数解耦改变 N 时门限不用跟着换。参数说明两点容易错。第一abs(x).^2 对实信号、复基带信号都成立如果接收机给的是复信号漏掉 abs 直接 x.^2 会把虚部丢掉统计量分布整个错掉。第二N 的取值一般 32 到 256太小中心极限近似不成立太大单次仿真变慢。2.2 阈值推导从目标虚警概率反解门限门限 λ 不能拍脑袋定。虚警概率 Pf 的定义是 H0 下被误判成 H1 的概率按上面的高斯近似Pf Q((λ - σw²) / sqrt(2*σw⁴/N))反解出门限λ σw² sqrt(2*σw⁴/N) * Q⁻¹(Pf)工程上通常先定一个可容忍的 Pf 目标再反解 λ。MATLAB 没有直接的 Q⁻¹但 Q(x) 0.5·erfc(x/√2)所以 Q⁻¹(Pf) √2·erfinv(1 − 2Pf)。写成函数function lambda compute_threshold(Pf, sigma2, N) % 由目标虚警概率反解能量检测门限 % Pf : 目标虚警概率取值 (0,1) % sigma2 : 噪声方差必须与仿真中的噪声功率对齐 % N : 采样点数 qinv sqrt(2) * erfinv(1 - 2 * Pf); % Q^{-1}(Pf) lambda sigma2 sqrt(2 * sigma2^4 / N) * qinv; end逻辑说明erfinv 是误差函数的逆MATLAB 自带不需要额外工具箱。公式推导里 sigma2^4 是标量运算写成 sigma2^4 没问题如果 sigma2 是向量就改成 sigma2.^4。注意门限计算、噪声生成、SNR 换算三处必须使用同一个噪声方差变量分开写死是这类仿真最常见的翻车源。参数说明Pf 目标一般取 0.01 到 0.1太小门限抬太高Pd 会被压下去N 越小门限对 Pf 越敏感因为 sqrt(2σw⁴/N) 这项变大同样 Pf 对应的门限离噪声均值越远。2.3 蒙特卡洛仿真Pd 和 Pf 的估计与 ROC 曲线有了检测函数和门限下一步用蒙特卡洛估计检测概率 Pd 和虚警概率 Pf。核心是重复生成 H0、H1 两组样本统计判对次数function [Pd, Pf] mc_energy_detection(SNR_dB, N, lambda, N_MC) % 蒙特卡洛估计能量检测的检测概率与虚警概率 % SNR_dB : 信噪比dB % N : 采样点数 % lambda : 判决门限 % N_MC : 蒙特卡洛重复次数 sigma2 1; % 噪声方差归一化 Ps sigma2 * 10^(SNR_dB / 10); % 由 SNR 反推信号功率 hits_det 0; hits_fa 0; for k 1:N_MC x0 sqrt(sigma2) * randn(N, 1); % H0: 纯噪声 T0 sum(abs(x0).^2) / N; if T0 lambda hits_fa hits_fa 1; % 虚警计数 end x1 sqrt(Ps sigma2) * randn(N, 1); % H1: 信号噪声 T1 sum(abs(x1).^2) / N; if T1 lambda hits_det hits_det 1; % 正确检测计数 end end Pd hits_det / N_MC; Pf hits_fa / N_MC; end逻辑说明每轮循环生成一对样本。H0 样本是零均值、方差 σw² 的高斯噪声H1 样本是零均值、方差 Ps σw² 的高斯序列代表信号加噪声的总功率。Pd、Pf 就是两类事件发生的频率N_MC 越大估计越稳。参数说明SNR_dB 扫 -20 到 0 dB 最有区分度太低 Pd 贴着 Pf_target 看不出算法差异N_MC 不低于 1000020000 更稳。跑完一组 SNR 就能画 Pd-SNR 曲线SNR_list -20:2:0; N 64; N_MC 20000; Pf_target 0.1; sigma2 1; lambda compute_threshold(Pf_target, sigma2, N); Pd zeros(size(SNR_list)); Pf_est zeros(size(SNR_list)); for i 1:length(SNR_list) [Pd(i), Pf_est(i)] mc_energy_detection(SNR_list(i), N, lambda, N_MC); end plot(SNR_list, Pd, -o, LineWidth, 1.5); xlabel(SNR (dB)); ylabel(Detection Probability P_d); grid on;逻辑说明先按目标 Pf 算好门限再对每个 SNR 跑一次蒙特卡洛。画出来的曲线通常是 S 型SNR 从 -20 dB 抬到 0 dBPd 从接近 Pf_target 爬升到 1。如果 0 dB 附近 Pd 还没到 0.9先查门限是不是按同一个噪声方差算的再查 N 是否太小。参数建议取值作用N32 ~ 256采样点数决定高斯近似精度N_MC10000 ~ 20000蒙特卡洛次数决定估计方差Pf_target0.01 ~ 0.1目标虚警概率决定门限高低SNR_dB-20 ~ 0感知信噪比扫描范围3. 主用户活动建模ON-OFF 过程与频谱空洞生成DSA 仿真的第二块地基是主用户活动模型。感知回答的是“此刻信道是否被占”接入决策需要知道“未来一段时间信道可不可用”这就要靠主用户活动的统计模型来支撑。3.1 ON-OFF 指数模型主用户活动的标准假设文献里最常见的主用户活动模型是两状态 ON-OFF 过程ON 代表主用户正在发射OFF 代表信道空闲。ON、OFF 持续时长都假设服从指数分布均值记为 T_on、T_off。主用户占空比就是 T_on / (T_on T_off)这个值直接决定频谱空洞总量。指数分布的无记忆特性让这个模型特别好用不管当前时刻处在什么状态剩余持续时间都服从同一参数分布。MATLAB 里生成指数随机数可以直接用 exprnd也可以用 -log(rand)*mean两者等价。生成整段活动序列的代码如下function [state, t] pu_activity_sequence(T_on, T_off, T_total, Ts) % 生成主用户 ON/OFF 活动序列 % T_on : ON 状态平均持续时间秒 % T_off : OFF 状态平均持续时间秒 % T_total: 仿真总时长秒 % Ts : 时间轴采样间隔秒 n round(T_total / Ts); state zeros(1, n); % 0 空闲, 1 占用 t (0:n-1) * Ts; pos 1; is_on false; % 初始从空闲开始 while pos n dur exprnd(is_on * T_on (~is_on) * T_off); steps round(dur / Ts); % 时长换算成采样点数 if is_on seg_len min(steps, n - pos 1); state(pos:pos seg_len - 1) 1; end pos pos steps; is_on ~is_on; % 状态翻转 end end逻辑说明循环每轮生成一段持续时间ON 段写入 1OFF 段跳过然后翻转状态直到填满整个时间轴。exprnd(is_on * T_on (~is_on) * T_off)用 0/1 开关在两档均值间切换比 if-else 分支紧凑。参数说明Ts 决定时间分辨率一般取感知周期的十分之一左右比如感知周期 100 msTs 取 10 ms。T_on、T_off 典型值 1~5 秒比值决定占空比T_on2、T_off3 时占空比 40%。初始状态选 OFF 是惯例也可以加参数控制因为初始状态会影响前几个空洞的时长分布。序列生成后第一件事是验证占空比统计 state 里 1 的比例应该接近 T_on/(T_onT_off)。偏差来源主要是有限时长下的首尾截断这部分在避坑章展开。3.2 从占用序列提取频谱空洞占用序列生成后下一步把空闲区间提取成频谱空洞列表。空洞是连续的 OFF 区间每个空洞有起止时刻这就是动态频谱接入能用的资源function holes find_spectrum_holes(state, Ts) % 从占用序列提取频谱空洞 % state : 0/1 占用序列 % Ts : 采样间隔秒 % 返回 : holes 为两列矩阵每行 [起始时刻, 结束时刻] holes []; n length(state); idx 1; while idx n if state(idx) 0 start_t (idx - 1) * Ts; % 空洞起始时刻 while idx n state(idx) 0 idx idx 1; end end_t (idx - 2) * Ts; % 空洞结束时刻 holes [holes; start_t, end_t]; else idx idx 1; end end end逻辑说明外层循环扫描遇到 0 就用内层循环吃掉整个连续空闲段记录起止时刻后继续遇到 1 直接跳过。边界上 idx 走到 n1 时(idx-2)*Ts 正好是最后一个空洞的末尾不会越界。拿到空洞列表可以做两个有用的统计一是空洞数量和平均空洞时长理论上平均空洞时长应接近 T_off二是把空洞时长画成直方图和指数分布概率密度叠在一起看拟合程度。这两步能快速暴露状态翻转类的 bug——漏翻转会让空洞时长翻倍这类问题从数值上看很容易被当成参数设错而漏掉。3.3 多信道场景独立 ON-OFF 过程的组合单信道模型只是第一步DSA 的多信道场景是 M 个独立的 ON-OFF 过程并行。这里“独立”是关键假设各信道的主用户活动互不影响。实现上对每个信道调用一次活动序列生成堆成 M×n 矩阵function [state_mat, t] gen_multi_channel(M, T_on, T_off, T_total, Ts) % 生成 M 个独立信道的 ON/OFF 活动矩阵 % state_mat : M×n 矩阵每行对应一个信道 state_mat zeros(M, round(T_total / Ts)); for ch 1:M [state_mat(ch, :), t] pu_activity_sequence(T_on, T_off, T_total, Ts); end end逻辑说明循环里每次调用都从 MATLAB 全局随机流取数各信道的 ON/OFF 序列统计独立这正是多信道场景需要的性质。多信道场景有三个留意点。第一各信道参数可以不同比如有的信道主用户活跃T_on 大有的信道常年空闲T_off 大这会让后面的信道排序决策有意义。第二M 取 4~16 比较合适太少体现不出排序价值太多矩阵占内存、决策循环变慢。第三MATLAB 随机流是全局的多信道生成时各信道消耗的随机数个数不同要复现结果必须固定 rng 种子。场景MT_on秒T_off秒占空比单信道基准12340%多信道均匀82340%多信道异构81~52~630%~60%4. 动态频谱接入决策从感知结果到信道分配感知和主用户模型就位后核心的 DSA 决策逻辑上场。这一章给出决策函数的输入输出设计、信道排序实现以及吞吐量和碰撞率的统计口径。4.1 决策输入感知结果向量与历史空闲概率DSA 决策的输入分三类本次感知结果向量、各信道历史空闲概率、当前占用状态。输出是一个动作——接入哪个信道或者等待。历史空闲概率用指数滑动平均更新每完成一次感知就刷新一次function idle_prob update_idle_prob(idle_prob, sense_vec, alpha) % 更新各信道历史空闲概率指数滑动平均 % idle_prob : 上一轮的空闲概率向量 % sense_vec : 本轮感知结果1 占用 / 0 空闲 % alpha : 平滑因子0~1越大越依赖近期观测 idle_prob (1 - alpha) * idle_prob alpha * (1 - sense_vec); end逻辑说明指数滑动平均是新旧观测的加权组合alpha 取 0.1~0.3 在响应速度和稳定性之间比较平衡。注意 sense_vec 为 1 时 (1-sense_vec) 为 0表示观测到占用就拉低空闲概率方向不能写反写反了等于“越占用越觉得空闲”整个排序逻辑都会颠倒。参数说明alpha 是这轮唯一的超参数。alpha 太小空闲概率更新迟钝信道状态已经变了排序还停留在旧印象alpha 太大一次虚警就把排序彻底打乱。多信道异构场景里alpha 取 0.2 左右通常表现最稳。4.2 信道排序与接入判决函数有了空闲概率决策就变成排序问题在感知为空闲、且当前未被占用的信道里选历史空闲概率最高的接入。直觉是空闲概率高的信道次用户接入后撞上主用户回退的概率更低function [action, ch_idx] dsa_decision(sense_vec, idle_prob, ch_busy) % 动态频谱接入决策 % sense_vec : 本轮感知结果向量1 占用 / 0 空闲 % idle_prob : 各信道历史空闲概率向量 % ch_busy : 各信道当前占用状态接入前检查 % action : 返回 access 或 wait % ch_idx : 接入的信道编号wait 时为空 free find(sense_vec 0 ch_busy 0); if isempty(free) action wait; ch_idx []; return; end [~, best] max(idle_prob(free)); % 空闲概率最高的信道 ch_idx free(best); action access; end逻辑说明find 里两个条件是交集——感知为空闲且当前未被占用。max 返回值里 best 是 free 数组内的下标必须映射回 free(best) 才是真实信道编号这个映射漏掉是常见笔误。wait 状态下次用户不发射进入下一轮感知。主循环的结构是每 Tsense 时隙调用一次感知更新 idle_prob调用 dsa_decision 决定接入或等待然后推进一个时隙。感知周期、决策周期、统计周期三者要一致别让决策函数访问到还没更新的状态。接入后要处理主用户回退次用户正在传输时主用户回来就发生碰撞。碰撞检测要借助主用户活动序列的 ground truthif strcmp(action, access) su_tx_time su_tx_time 1; if state_mat(ch_idx, t_idx) 1 % 主用户实际占用 coll_time coll_time 1; ch_busy(ch_idx) 1; % 标记被挤占触发频谱切换 end end逻辑说明这段代码放在仿真主循环里每个时隙执行一次。碰撞判断依据是 ground truth 而不是感知结果——感知可能虚警也可能漏检用 ground truth 评估碰撞率才公平。ch_busy 置 1 让次用户立即停止传输这是 DSA 里“频谱移动性”的最小实现。4.3 吞吐量与碰撞率的统计口径DSA 仿真最后输出两个核心指标但统计口径不同结果差异很大。这套系统里吞吐量定义为次用户成功传输的时隙数除以总时隙数成功传输指次用户发射且主用户不在场。更严谨的说法是“有效接入率”它不是信息论意义的 bps。碰撞率定义碰撞时隙数除以次用户发射时隙数。分母很重要——如果拿总时隙数当分母碰撞率会被大量空闲时隙稀释数值好看但没意义。文献里两种口径都有人用写论文时必须在图表说明里写清楚。仿真结束时的汇总输出throughput su_success / total_slots; % 有效接入率 collision coll_time / su_tx_time; % 碰撞率按发射时隙归一 fprintf(Throughput %.3f, Collision rate %.3f\n, throughput, collision);逻辑说明total_slots T_total / Ts 是总时隙数su_success 每时隙累计一次条件是次用户发射且 ground truth 为空闲。两个指标要一起看只调低感知门限更激进接入能提升吞吐量碰撞率必然上升DSA 的本质就是在吞吐量和碰撞率之间做权衡。举个例子感知门限调低后虚警变多次用户把空闲信道误判成占用选择等待吞吐量下降门限调高后漏检变多次用户大胆接入本来就占用的信道碰撞率上升。所以门限不能只按感知性能来定还要结合 DSA 层的代价函数一起看这是把感知和接入串成闭环后才会意识到的问题。参数说明碰撞率分母用 su_tx_time 而不是 total_slots这是仿真统计里最容易口径不一的地方换人复现时经常对不上数建议在代码注释里写明归一化方式。5. 认知无线电仿真避坑五个高频问题的现象与根因这块是实践里最容易卡人的部分我把五个高频问题按“现象 → 原因 → 解决”整理出来都是实际跑仿真攒下来的血泪经验。5.1 虚警概率实测值与目标值对不上现象目标 Pf 设 0.1蒙特卡洛跑出来 Pf_est 是 0.2 甚至更高而且门限调了也没改善。原因门限计算用的噪声方差和仿真生成噪声用的方差不一致。最常见的是门限函数里写死 sigma2 1而仿真里噪声用了别的方差或者信号先做了归一化导致噪声方差不再是 1。比如生成噪声时用 randn(N,1) 再乘 sqrt(2)实际噪声方差是 2门限却按 1 算虚警概率直接翻倍不止。解决把噪声方差做成全局统一变量计算门限、生成噪声、SNR 换算都引用同一个 sigma2。我在仿真入口加一行断言校验噪声生成的实际方差与门限计算值一致不一致直接报错。排查时可以在蒙特卡洛循环外统计 H0 样本的平均能量它应该收敛到噪声方差收敛值不对就是生成环节的问题。5.2 Pd 曲线抖动剧烈、趋势不稳定现象同样的参数跑两次Pd 在 -10 dB 附近一次 0.6、一次 0.45曲线像锯齿看不出单调趋势。原因N_MC 太小统计涨落主导了结果。Pd 落在 0.5 附近时估计标准差约为 sqrt(Pd(1-Pd)/N_MC)N_MC 1000 时单点标准差约 1.6%叠加到整条曲线上视觉上就很明显。解决N_MC 提到 20000 以上。如果嫌慢先用 2 dB 的 SNR 步长快速看趋势确认 S 型形状后加密扫描点。想进一步压抖动可以用多条不同种子的曲线做平均或对曲线做滑动平均后处理但论文主图建议用原始曲线后处理过的曲线要在说明里交代清楚。每次跑之前固定 rng 种子论文里的曲线才能复现。5.3 频谱空洞数量与平均时长偏离理论值现象T_off 3 秒统计出的平均空洞时长只有 1.8 秒空洞数量却比理论值多一倍。原因时间轴采样 Ts 太粗短空洞被采样吃掉导致时长偏小更隐蔽的原因是状态翻转逻辑写错比如循环末尾漏掉 is_on ~is_on导致 ON 段之后还跟着 ON 段空洞被合并。此外有限时长下首尾截断也有影响序列从 OFF 开始、可能在 ON 段被截断第一个空洞和最后一个 ON 段的期望都偏离理论值。解决先用极端参数自测——T_on 设 0序列应全 0T_off 设 0序列应全 1。两个用例过了再检查 TsTs 至少比 T_off 小一个数量级。仿真时长取 T_on T_off 的 50 倍以上截断影响就小到可以忽略。5.4 SNR 口径混乱导致感知曲线异常现象SNR 扫到 0 dB 时 Pd 还不到 0.5或者 -20 dB 时就接近 1曲线完全违背直觉。原因dB 和线性功率的换算搞混。Ps 10^(SNR_dB/10) 是功率口径有人误写成 10^(SNR_dB/20)那是电压幅值口径用在功率里信号功率直接差了一倍。如果仿真里用了复基带信号噪声方差是实数基带的 2 倍SNR 换算公式也要跟着调整这类口径问题在复信号模型里更隐蔽。解决在代码注释里写明 SNR_dB 是功率比并在生成 H1 样本时打印一次 Ps 和 sigma2 的比值回验。Ps/sigma2 10^(SNR_dB/10) 说明口径正确否则就是 20 那个系数写错了。5.5 多信道仿真每次运行结果不同现象同样的脚本、同样的参数两次运行的空洞直方图、吞吐量都对不上想复现论文数据很痛苦。原因MATLAB 全局随机数流没有复位。多信道生成、噪声生成、蒙特卡洛都从同一个流里取数调用顺序一变结果就全变。解决仿真主函数第一行强制 rng 固定种子且必须放在所有随机数生成之前包括主用户序列生成、噪声生成、蒙特卡洛循环放在中间等于没固定。我一般把种子写进参数表每次实验保留“参数 种子”的完整记录回头查结果时不用靠猜。6. 进阶把固定感知周期改成自适应感知间隔6.1 固定感知周期的浪费在哪里前面所有仿真里次用户每 Tsense 固定感知一次。这在主用户活动忽快忽慢的场景里很浪费信道长时间空闲时高频感知白耗能量信道频繁切换时固定周期又容易错过主用户回退。自适应的思路是让感知间隔跟着状态走——状态越稳感知越稀变化越频感知越密。6.2 自适应感知间隔的实现与对比实现是在主循环里维护动态间隔 T_sense每次感知后按结果与预测是否一致来调整一致说明模型准拉长间隔不符说明状态在变收紧间隔。surprise (sense_vec(ch) ~ pred_occupy(ch)); % 与预测不符视为意外 if surprise T_sense max(T_min, T_sense / shrink); % 收紧更快感知 else T_sense min(T_max, T_sense * expand); % 放松减少感知 end逻辑说明surprise 为 1 时收紧为 0 时放松shrink 取 2、expand 取 1.5 是常用组合T_min、T_max 负责钳位防止间隔失去意义。pred_occupy 由上一轮感知结合空闲概率给出它本身就是一个简单的状态预测器。对比实验跑同一组主用户序列固定基准 T_sense 100 ms自适应从 100 ms 起步T_min 20 ms、T_max 500 ms。现象有两个主用户频繁切换时段自适应感知更密碰撞率下降长空闲时段感知间隔拉长感知开销下降吞吐量几乎不受影响。收益在碰撞率和能耗上对比时别只盯吞吐量。注意 T_min 不能小于感知周期的物理下限采样、判决、切换都耗时我会在参数表里加一行 min_sense_time 校验超了直接报错。从那以后我每次跑 DSA 仿真都强制走一遍固定流程先验主用户占空比再查 SNR 口径固定 rng 种子最后才看吞吐量和碰撞率这套习惯帮我省掉大量重复调试。上面这些模块都是这套仿真系统里实际跑通的部分照着搭就能复现核心曲线。希望帮到你。本文还有配套的精品资源点击获取