简介面向MATLAB毕业设计的Voigt模型拟合项目以简洁的代码与文档展示如何通过自定义Voigt函数结合最小二乘思路对光谱类数据做谱形拟合。Voigt模型融合Lorentzian与Gaussian分布适用于光谱学、核磁共振和声学等宽峰信号分析场景项目重点给出拟合主程序与模型函数并配套演示脚本、说明文档及示例效果图便于从算法到应用快速理解。压缩包共8个文件以3个m源码文件为主辅以md说明、txt文本、png样图及Git配置信息整体仅31KB结构轻量目前已有98人学习适合需要完成相关课题或补强MATLAB数值拟合能力的本科生作为参考。通过阅读主程序与演示脚本可掌握Voigt模型建模、参数拟合与结果可视化整个流程也可迁移至其他线型拟合任务。1. 从光谱峰形到Voigt拟合这个MATLAB毕业设计在解决什么实际光谱分析中我们测到的峰形往往既不是纯高斯也不是纯洛伦兹而是两者卷积后的形态。Voigt模型能同时描述多普勒展宽带来的高斯成分和碰撞展宽带来的洛伦兹成分在红外光谱、核磁共振、X射线衍射等领域非常常用。这个zip包里装的正是一套用MATLAB实现Voigt模型拟合的毕业设计源码myvoigt.m负责计算Voigt分布曲线voigtfit.m是拟合主函数voigtfit_demo.m告诉你如何把两者串起来最后还会输出一张spectrumfit.png拟合结果图。如果你正在准备“matlab 毕业设计”题目或者手头有需要分离高斯-洛伦兹贡献的光谱数据这套代码值得逐行拆解。2. Voigt模型的数学本质与myvoigt.m的两种实现路径2.1 Voigt分布是Gaussian与Lorentzian的卷积不是简单加和Voigt分布的定义是高斯分布与洛伦兹分布的卷积积分数学上写成V(x) ∫ G(x) L(x - x) dx其中G是高斯函数L是洛伦兹函数。这个积分没有封闭的解析表达式只能通过数值积分或近似公式来计算。很多初学者会直接把高斯和洛伦兹曲线做线性叠加但那只是pseudo-Voigt不是真正的Voigt。在工程上两种做法都有应用但原理和误差边界完全不同需要根据拟合的物理场景来选择。2.2 数值积分实现基于quadgk的精确卷积myvoigt.m中一种直接实现方式是利用MATLAB的integral函数做无穷积分。下面这段代码展示了最朴素的写法function y myvoigt(x, sigma, gamma) % 数值积分实现的Voigt函数 % x: 频率或波长偏移sigma: 高斯展宽参数gamma: 洛伦兹展宽参数 y zeros(size(x)); for i 1:length(x) f (t) exp(-t.^2 / (2 * sigma^2)) ./ ... (1 ((x(i) - t) ./ gamma).^2); y(i) (1 / (sqrt(2 * pi) * sigma)) * ... (gamma / pi) * integral(f, -Inf, Inf); end end这个实现的好处是概念上严格对应数学定义适合用来做理论校验。但实际运行时会发现它非常慢原因是integral需要多次迭代求积分而且for循环遍历每个x点整体复杂度很高。如果光谱数据有上千个点拟合一次可能要数十秒甚至更久这在批量拟合场景中不可接受。另外当sigma或gamma趋近于0时被积函数变化过于尖锐数值积分可能不收敛。2.3 工程近似pseudo-Voigt加权组合更常用的做法是pseudo-Voigt模型它用高斯和洛伦兹的加权和来近似Voigt。权重参数eta控制两者比例eta0时是纯高斯eta1时是纯洛伦兹。这个近似在拟合精度要求不高时误差通常在2%以内计算速度却快三个数量级。function y myvoigt(x, sigma, gamma, eta) % 伪Voigt近似Gaussian与Lorentzian线性加权 % eta 为0时是纯Gaussian为1时是纯Lorentzian % 注意这里未包含幅值系数调用时可在外部乘amplitude g exp(-(x.^2) / (2 * sigma^2)); l 1 ./ (1 (x ./ gamma).^2); y (1 - eta) .* g eta .* l; end这段代码中的sigma对应高斯宽度gamma对应洛伦兹半宽eta做插值。实际拟合时为了让曲线峰值与数据匹配通常会在外层乘以一个幅值参数。值得注意的是pseudo-Voigt的eta与真实Voigt的sigma/gamma之间存在一个经验关系更严谨的做法是使用Olivero和Longbothum提出的近似公式但多数光谱拟合场景直接让eta自由拟合也能得到可接受的曲线形状。2.4 两种实现的精度与耗时对比表实现方式理论精度单次计算相对耗时适用场景数值积分高约千倍理论验证、生成仿真数据pseudo-Voigt近似较高极快光谱拟合、批量处理、实时分析从voigtfit_demo.m的实际运行效果看默认使用的应该是pseudo-Voigt近似因为spectrumfit.png中拟合曲线平滑且计算很快。若你在自己的数据上发现拟合残差呈现系统性的S形偏差可以尝试切回数值积分版本对比两种实现的结果来检查是否近似误差导致的。3. voigtfit.m 核心拟合算法与MATLAB优化工具箱对接3.1 函数接口设计与调用约定voigtfit.m是整套代码的拟合入口。按照MATLAB函数文件惯例它的调用方式通常类似[param, resnorm, residual] voigtfit(x, y, initParams, lb, ub)参数含义分别为x是自变量向量y是观测光谱强度initParams是模型参数的初始猜测lb和ub是参数下界和上界。返回的param是拟合得到的最优参数向量resnorm是残差平方和。如果调用时没有提供lb和ub则内部默认使用足够宽松的边界例如所有参数的下界为0上界为数据最大值的10倍。这样的设计让初级用户可以直接运行demo但做正式实验时建议显式传入边界。3.2 用lsqcurvefit做非线性最小二乘拟合MATLAB的Optimization Toolbox中lsqcurvefit是求解这类非线性曲线拟合的标准函数。voigtfit.m内部的常见实现方式是把myvoigt包成一个模型函数句柄然后交给lsqcurvefit迭代优化。下面是一段核心代码框架function [param, resnorm] voigtfit(x, y, init, lb, ub) % 定义带幅值的Voigt模型函数句柄 model (p, xdata) myvoigt(xdata, p(1), p(2), p(3)) .* p(4); % 配置优化选项显示迭代信息设置高精度和最大迭代次数 opts optimoptions(lsqcurvefit, ... Display, final, ... TolFun, 1e-10, ... TolX, 1e-10, ... MaxIter, 2000); % 执行拟合 [param, resnorm] lsqcurvefit(model, init, x, y, lb, ub, opts); end这里模型句柄中p(1)是sigmap(2)是gammap(3)是etap(4)是幅值amplitude。之所以把幅值单独放在myvoigt外面是因为myvoigt本身体现的是归一化形状而实际光谱的强度是任意单位的必须有一个幅值缩放项。注意lsqcurvefit要求x和y都是列向量如果传入行向量MATLAB不会报错但会得到错误结果这是常见坑点之一。3.3 参数初值、上下界与归一化策略初始值设置得不好lsqcurvefit很容易收敛到局部极小值。对于Voigt拟合合理的初值估计方法是先找到数据中的最大峰值位置和峰值强度用峰值强度作为幅值的初值然后测量半高宽FWHM由于FWHM与sigma、gamma的关系不是线性的可以设置sigmaFWHM/2.354gammaFWHM/2eta0.5。如果数据中噪声很强可以先用smooth函数做一次平滑再计算FWHM。在myvoigt的x不包含中心位置时需要先将数据以峰位为中心进行平移也就是把x减去峰位后再传入拟合函数。若希望拟合过程中峰位也能自由调整就需要在模型函数里加入center参数例如model (p, xdata) myvoigt(xdata - p(5), p(1), p(2), p(3)) .* p(4);其中p(5)为峰位。此时初值p(5)设为最大峰对应的x坐标上下界则根据数据范围设定比如峰位只允许在数据最左侧和最右侧之间变化。边界设置的原则是既不能太紧导致正确解在边界外也不能太松导致算法搜索过大区间而浪费迭代。3.4 拟合质量指标SSE、R²与残差分析lsqcurvefit返回的resnorm就是SSE但单看SSE不能判断拟合是否良好。更常用的指标是决定系数R²计算方式如下SSE resnorm; SST sum((y - mean(y)).^2); R2 1 - SSE / SST;R²越接近1表示拟合越好但要注意R²高不代表模型正确。如果残差residual y - model(param, x)存在明显的周期性或形状说明模型缺少基线项。很多实际光谱数据会有线性漂移需要在模型里增加一个一次项model (p, xdata) myvoigt(xdata - p(5), p(1), p(2), p(3)) .* p(4) p(6) p(7) * xdata;这样扩展后p(6)是常数基线p(7)是斜率。添加基线项后SSE通常会显著下降但要防止过拟合可以用赤池信息准则AIC来权衡。4. voigtfit_demo.m 实战演示从数据到拟合并保存spectrumfit.png4.1 生成模拟Voigt谱线数据voigtfit_demo.m的第一步通常是用已知参数模拟一组光谱数据用来验证拟合代码的正确性。由于真实数据的真值是未知的模拟数据可以让我们对比拟合结果与设定值的偏差。常见做法如下x linspace(-8, 8, 600); sigma 0.8; gamma 0.6; eta 0.4; amp 3.0; y_clean myvoigt(x, sigma, gamma, eta) * amp; rng(2024); % 固定随机种子使结果可重现 y_noisy y_clean 0.1 * randn(size(x));这里x取[-8,8]共600点保证峰两侧有足够长的尾部。噪声标准差设为0.1模拟真实测量中的随机误差。固定随机种子后每次运行结果一致方便在答辩中演示时保证可复现性。4.2 调用voigtfit.m完成拟合有了模拟数据后直接调用voigtfit.minit [1.0, 0.5, 0.5, max(y_noisy)]; lb [0.01, 0.01, 0, 0]; ub [5, 5, 1, 20]; p voigtfit(x, y_noisy, init, lb, ub);这里初值选择依据是sigma和gamma分别取1.0和0.5eta取0.5为中性值幅值用数据最大值。经过几次迭代后p应当收敛到接近[0.8, 0.6, 0.4, 3.0]的值。如果发现p与设定值相差超过20%说明算法可能陷入了局部极小或者数据点数不够。此时可以尝试将初值改为sigma0.7, gamma0.7, eta0.3看是否收敛到不同结果。4.3 绘制并保存结果图拟合完成后demo脚本会绘制原始数据点和拟合曲线并保存为spectrumfit.png。绘图代码通常为y_fit myvoigt(x, p(1), p(2), p(3)) * p(4); figure(Color, w); plot(x, y_noisy, ob, MarkerSize, 4, LineWidth, 1); hold on; plot(x, y_fit, -r, LineWidth, 2); legend(原始数据, Voigt拟合, Location, north); xlabel(相对位置 (cm^{-1})); ylabel(强度 (a.u.)); grid on; saveas(gcf, spectrumfit.png);这一段代码涉及MATLAB绘图的基本操作。plot的第一个参数是x第二个是y第三个是样式字符串ob表示蓝色圆圈标记-r表示红色实线。legend的Location设为north意思是图例放在顶部中间。saveas支持png、pdf、eps等格式写论文时建议输出矢量图saveas(gcf, spectrumfit, epsc)。绘图是毕业设计中的加分项清晰美观的图比参数表格更直观。4.4 对demo脚本做参数敏感性测试为了让答辩更有说服力可以修改demo脚本做一个简单的敏感性分析把初值在真实值附近上下浮动50%循环100次调用voigtfit记录每次的收敛结果。如果100次中绝大多数都收敛到同一组参数说明拟合的稳定性好如果发散或收敛到两种不同结果说明模型在该数据条件下过参数化。这个测试用一个小循环就能完成也是评审老师比较关心的部分。5. 毕业设计答辩前必看的边界条件与排错指南5.1 数据量、采样率对拟合稳定性的影响如果数据点数太少比如不足50个那么Voigt函数的尾部信息不足以约束gamma和eta拟合结果会非常不稳定。尤其是gamma参数它对峰尾的拖尾形状非常敏感采样范围不够时gamma会倾向于变大或变小以吸收噪声。一般建议采样范围至少覆盖峰中心两侧各3个半高宽数据点不少于200。如果原始数据不满足这个条件可以通过插值增密但插值会引入额外误差更好的做法是采集时就保证范围足够。5.2 初始化不当导致的局部极小值与规避lsqcurvefit是局部优化算法初值离全局最优太远时容易陷入局部极小。规避方法有三个第一先用smooth函数做滑动平均找到平滑数据中的峰值和半高宽用平滑结果计算初值第二先固定eta为0和1分别用两参数模型sigmaamp和gammaamp拟合取残差较小的初始eta第三使用全局优化工具箱中的MultiStart或GlobalSearch配合lsqcurvefit虽然计算量更大但能显著提高找到全局解的概率。对于毕业设计而言使用MultiStart可以作为一个展示点代码量也不大。5.3 工具箱版本差异与setparam注意事项不同MATLAB版本中优化选项的设置方式有差异。R2013a之前的版本使用optimset之后推荐optimoptions。如果你打开老代码看到optimset(Display,final)在R2024b上运行通常没问题但如果是optimset(Algorithm,trust-region-reflective)则有可能因算法选择受限而报错。另一个容易踩的坑是变量名冲突如果工作区中存在名为param或model的变量而函数内部也使用相同变量名会导致MATLAB提示“Subscript indices must either be real positive integers or logicals”之类的错误。强烈建议运行demo之前先执行clear; clc;清理工作区。5.4 批量处理多个峰时的策略如果你的光谱数据有多个峰不能直接使用现在的voigtfit因为模型函数只支持单个Voigt。常见做法是用findpeaks函数找到所有局部极大值然后对每个峰截取一个区间单独拟合最后合并结果。区间范围可以用相邻峰谷作为分界。这样做的问题在于重叠峰的尾部会相互干扰严格的处理方式是把所有峰放在一个模型函数中同时拟合此时参数数量成倍增加需要给每个峰分配5个参数sigma、gamma、eta、amplitude、center。批量拟合时不同峰的参数初值可以复用同一个先验值但边界可能不同需要动态生成上下界向量。6. 让这个项目往后延伸的3个实用技巧6.1 用fittype自定义模型让fit函数更灵活如果你的MATLAB安装了Curve Fitting Toolbox可以用fittype把myvoigt包装成可识别的模型直接使用fit函数。这个方法适合快速批处理因为fit会自动计算参数的95%置信区间ft fittype(amp * myvoigt(x - center, sigma, gamma, eta), ... independent, x, ... problem, center); % 如果不希望center自由拟合设为problem opt fitoptions(ft); opt.StartPoint [1, 0.5, 0.3, 3]; [fresult, gof] fit(x, y_noisy, ft, opt);这里gof结构体包含sse、rsquare、adjrsquare等统计量比手算更方便。注意fit要求自变量单调递增如果数据顺序混乱先用sort排序。6.2 把参数结果导出成table直接用于论文论文中通常需要表格列出多组拟合参数。把每次拟合的结果存成行向量最后合并成MATLAB table并使用writetable输出csv或xlsxT table(p(1), p(2), p(3), p(4), SSE, R2, ... VariableNames, {sigma, gamma, eta, amplitude, SSE, R2}); writetable(T, fit_summary.csv);这样生成的csv文件可以直接导入到Excel或LaTeX表格中避免手动复制粘贴出错。对于有多个样本的数据可以在for循环中不断用T_temp拼接。6.3 扩展为多峰Voigt拟合的思路把单峰Voigt扩展成多峰的方法很简单模型函数内部循环累加。下面给出一个支持可变峰数量的模型函数框架function y multi_voigt(x, p, npeaks) y zeros(size(x)); for k 1:npeaks idx (k-1)*4 1; sigma p(idx); gamma p(idx1); eta p(idx2); amp p(idx3); center p(idx4); y y amp * myvoigt(x - center, sigma, gamma, eta); end end调用时令初始参数向量长度为npeaks*5依次排列所有峰的参数。需要注意峰位置的初值必须通过findpeaks预先确定否则优化算法很容易让峰位置乱跑。多峰拟合对初值要求更高建议用MultiStart辅助寻找全局解。如果峰之间重叠严重参数相关性会增强拟合结果的可解释性下降此时可考虑将部分参数固定比如所有峰的eta使用同一个值减少自由参数数量。这样一来这个毕业设计已经从简单的单峰拟合延伸到了真正的分析工具在答辩时可以明显提升项目复杂度评价。本文还有配套的精品资源点击获取
