FEATURED · 精选文章

基于Matlab的视日轨迹跟踪算法实现与工程化应用

发布时间 / 2026/9/13 15:58:33
来源 / 创域科博编辑部
栏目 / 资讯中心
基于Matlab的视日轨迹跟踪算法实现与工程化应用 简介这是一份基于MATLAB实现的视日轨迹跟踪算法仿真资源主要面向新能源、自动化等专业的本科与硕士生以及从事光伏支架控制、太阳能利用等相关研究的工程技术人员。资源围绕太阳位置计算与视日运动轨迹模拟展开能够为课程设计、毕业设计或科研课题提供可直接复现的算法参考与实验起点。压缩包内共4个文件核心是一个MATLAB程序文件可直接运行仿真流程另有两个png格式的效果展示图便于直观核对输出结果一个txt格式说明文档对版本与运行环境等作了简要交代。整包大小约459KB轻量、紧凑下载后可快速部署到本机环境。目前已有158人浏览学习此资源使用时可结合源码、说明和截图快速理清视日轨迹跟踪的建模思路与实现步骤并在此基础上进行参数修改、功能扩展或算法改进对于深入理解光伏视日运动规律和开发相关控制策略都很有帮助。1. 视日轨迹跟踪算法在Matlab里的工程化落地视日轨迹跟踪算法解决的问题很直接给定所在地经纬度、当前日期与时刻算出太阳此刻在哪。与常见的光敏传感器闭环追踪方案相比这种开环前馈计算不依赖天气和硬件反馈在光伏双轴支架、塔式光热定日镜场、建筑遮阳与采光分析中都有直接用途。这个Matlab工程包内置了完整的main.m运行脚本基于2014a及以上版本均可直接运行输出太阳高度角与方位角的全天曲线。我拿到源码后的第一反应是核对它的时间基准与坐标系定义因为这两处是后续所有计算正确与否的分水岭也是论文里最容易被评委追问的细节。2. 赤纬角、时角与真太阳时时间基准先立住视日轨迹算法的本质是把太阳在黄道上的运动投影到观测者所在的水平坐标系。这一步要先后经过赤道坐标系和地平坐标系中间涉及三个关键输入赤纬角太阳直射点纬度、时角太阳相对当地子午线的角位移、以及一个干净的时间基准。任何一项有偏差最终的高度角和方位角都会整体偏移所以这一章先把天文基础打牢。2.1 赤纬角的两个模型Cooper与Spencer的取舍赤纬角的定义是地心天球坐标系中太阳中心相对天赤道的角距离取值范围在-23.45°到23.45°之间。工程中最常用的简化模型是Cooper方程δ 23.45 × sin(2π/365 × (284 N))其中N为年积日1月1日为1。这个模型的推导基础是把地球公转近似为椭圆与自转轴的固定夹角误差来源是忽略黄经与真近点角的差异全年最大偏差约0.35°出现在春秋分附近。如果用来驱动双轴跟踪支架Cooper模型的偏差换算成跟踪角度误差可能让光斑偏移一个塔式定日镜的镜面宽度此时需要用Spencer模型做更高精度的逼近% 年积日计算 N datenum(year, month, day) - datenum(year, 1, 1) 1; % Spencer赤纬角(弧度)B为回归年比例系数 B 2*pi*(N-1)/365; delta 0.006918 - 0.399912*cos(B) 0.070257*sin(B) ... - 0.006758*cos(2*B) 0.000907*sin(2*B) ... - 0.002697*cos(3*B) 0.001480*sin(3*B);这段代码把年积日换算成弧度相位B再代入7项傅里叶级数。Spencer模型的精度在0.01°量级相比Cooper模型提升了一个数量级代价只是多5次三角函数运算对现代处理器来说可以忽略。工程包里选择Spencer还是Cooper可以从仿真时间跨度来推断单日曲线用哪个都行全年逐时扫描就必须用Spencer否则累加误差会在年发电量估算中放大。两个模型的适用边界可以这样对照| 模型 | 最大赤纬误差 | 适用场景 | 计算成本 | | Cooper | ±0.35° | 快速估算、教学演示 | 1次sin | | Spencer | ±0.02° | 跟踪控制、全年仿真 | 5次sin/cos |提示赤纬角公式里的365与366闰年差异很小实际影响低于0.01°多数工程实现不区分平闰年直接用365即可。2.2 时差方程与经度修正从钟表时间到真太阳时视日轨迹算法要求的时间是真太阳时日常使用的钟表时间是平太阳时。两者差值由两个因素叠加一是地球公转速率不均匀导致的时差方程Equation of Time二是观测点经度与所在时区中心经度的偏移。工程里常见的错误是把钟表时间直接当太阳时使用这样在北京东经116.4°时正午会出现约14分钟的偏移换算成时角就是3.6°足以让高度角误差达到0.5°以上。% 时差方程(分钟)公式与Spencer同源 EoT 229.18 * (0.000075 0.001868*cos(B) - 0.032077*sin(B) ... - 0.014615*cos(2*B) - 0.040849*sin(2*B)); % 真太阳时(小数小时) solar_time hour minute/60 ((lon - 15*tz)*4 EoT) / 60; % 时角(度)每1小时对应15度 H_deg 15 * (solar_time - 12);代码里lon代入的是东经正值、西经负值tz代入的是UTC偏移小时数。经度修正项(lon - 15×tz)×4的含义是经度每偏离时区中央经线1度地方时间差4分钟东侧时间领先所以为正值。EoT的单位是分钟除以60换成小时后再叠加。这个修正完成后H_deg在日出时约为负值正午时为0日落时转为正值方向约定与后续方位角公式严格配套。2.3 高度角与方位角球面三角形的两种解法有了赤纬角δ和时角H高度角α由球面余弦定理直接给出。这里需要把经纬度和角度全部统一为弧度再计算否则Matlab的sin/cos函数会因为单位混用给出完全错误的结果。lat_r lat * pi/180; % 纬度转弧度 dec_r delta; % delta已是弧度 H_r H_deg * pi/180; % 时角转弧度 % 高度角(度) sin_alt sin(lat_r)*sin(dec_r) cos(lat_r)*cos(dec_r)*cos(H_r); alt asin(sin_alt) * 180/pi; % 方位角(度从正北顺时针) az_r atan2(sin(H_r), cos(H_r)*sin(lat_r) - tan(dec_r)*cos(lat_r)); az mod(az_r*180/pi 180, 360);方位角的求解有两个常见路线。第一条路是通过余弦定理求出相对南方的角度再用上午/下午判断东西侧需要额外的逻辑分支上面代码用的是atan2单一个表达式分母的综合相位天然对应了太阳在天空的象限省去判断分支。atan2的优势是值域覆盖-pi到pi不会出现acos那种在0°、180°附近对微小扰动过于敏感的情况。mod(az180,360)把参考零点从南方移到正北得到的是北偏东为正的标准气象方位角。这套函数封装完成之后剩下的就是main.m里的逐时刻扫描与可视化。我在第3章把它展开到完整可运行的脚本。3. main.m 实现拆解从函数封装到全天曲线工程包的main.m脚本承担三个任务设定站点与日期参数、逐时刻调用太阳位置函数、绘制高度角与方位角随时间的曲线。下面给出的实现保留了原工程的核心流程并补齐了日出日落界面的处理是实际工作中最常用的写法。3.1 输入参数区与函数封装站点参数放在脚本顶部集中管理改一个城市只需要改四个数字。这里用datenum统一处理日期跨月、跨年循环时不用手工计算每月天数也不容易漏闰年。clc; clear; close all; % 站点北京 lat 39.9042; % 纬度(°)北纬为正 lon 116.4074; % 经度(°)东经为正 tz 8; % 时区(小时)UTC8 % 日期与扫描步长 year 2024; month 6; day 21; % 夏至日 step_min 10; % 扫描步长(分钟) % 生成时间序列 t_start datenum(year, month, day, 0, 0, 0); t_stop datenum(year, month, day, 23, 59, 59); t_vec t_start : step_min/(24*60) : t_stop;这里的datenum把时间编码为连续浮点数1代表1天。step_min/(24×60)把分钟换算成天数增量t_vec就得到从当天0点到23:59的全部时间戳。用浮点时间的好处是调用函数时直接解析时分秒不需要单独维护每个小时刻的索引。3.2 逐时刻扫描与数组预分配太阳位置函数需要独立的year、month、day、hour、minute参数所以循环里要从t_vec反解这些量。Matlab的datevec可以一次完成n length(t_vec); alt zeros(n, 1); az zeros(n, 1); for i 1:n [y, mo, d, h, mi, ~] datevec(t_vec(i)); [alt(i), az(i)] solar_position(lat, lon, y, mo, d, h, mi, tz); enddatevec返回的第6个元素是秒用波浪线丢弃。数组alt、az预先用zeros分配避免循环里动态扩容导致的性能下降——这在扫描全年、步长1分钟的场景下约52万次循环能省下10倍以上的时间。solar_position函数内部按第2章的公式逐步求值其输入输出约定与绘图脚本完全解耦单独维护函数文件的好处是后续切换Cooper/Spencer模型时只动一处。function [alt_deg, az_deg] solar_position(lat_deg, lon_deg, ... year, month, day, hour, minute, tz) % 输入均为原始单位输出高度角、方位角(°从北顺时针) N datenum(year, month, day) - datenum(year, 1, 1) 1; B 2*pi*(N-1)/365; delta 0.006918 - 0.399912*cos(B) 0.070257*sin(B) ... - 0.006758*cos(2*B) 0.000907*sin(2*B) ... - 0.002697*cos(3*B) 0.001480*sin(3*B); EoT 229.18*(0.000075 0.001868*cos(B) - 0.032077*sin(B) ... - 0.014615*cos(2*B) - 0.040849*sin(2*B)); solar_time hour minute/60 ((lon_deg - 15*tz)*4 EoT)/60; H_deg 15 * (solar_time - 12); lat_r lat_deg * pi/180; dec_r delta; H_r H_deg * pi/180; sin_alt sin(lat_r)*sin(dec_r) cos(lat_r)*cos(dec_r)*cos(H_r); alt_deg asin(sin_alt) * 180/pi; az_r atan2(sin(H_r), cos(H_r)*sin(lat_r) - tan(dec_r)*cos(lat_r)); az_deg mod(az_r*180/pi 180, 360); end这个函数里纬度和经度、时区作为普通数值参数传入便于后续用数组或表格驱动批量仿真。一个需要注意的边界是当高度角为负太阳在地平线下时方位角数值仍有意义但水平坐标系下不含大气层的几何计算会给出一个虚拟太阳方向画图时通常用逻辑索引过滤掉。3.3 绘图与时间轴格式化主程序的绘图部分需要同时呈现高度角与方位角两条曲线并把横轴格式化为时刻而非浮点日期数figure(Color,w,Position,[100 100 900 400]); subplot(1,2,1); plot(t_vec, alt, b-, LineWidth, 1.2); hold on; plot(t_vec, zeros(size(t_vec)), k--); xlabel(时刻); ylabel(太阳高度角(°)); title(北京 2024-06-21 高度角); datetick(x, HH:MM); xlim([t_start t_stop]); ylim([-10 90]); grid on; subplot(1,2,2); plot(t_vec, az, r-, LineWidth, 1.2); xlabel(时刻); ylabel(太阳方位角(°)); title(北京 2024-06-21 方位角); datetick(x, HH:MM); xlim([t_start t_stop]); ylim([0 360]); grid on;datetick(x,HH:MM)是Matlab里处理时间轴的关键函数它把datenum浮点数横轴重新标记为人类可读的时分格式。xlim与t_start、t_stop绑定可以防止datetick在数据范围外多画出多余刻度。黑色虚线标记地平线高度角曲线与虚线的两个交点就是日出与日落时刻从图上可以直接读出来。如果实测中发现交点位置与天文年历不同问题基本都出在第2章的时间修正项上。4. 仿真结果判读与精度校验误差藏在这些地方模型跑出曲线只是第一步判断曲线对不对、误差来自哪里是工程落地前必须做的功课。这一章给出三条校验路径正午峰值对照、全天形态比对、边界时刻验证。4.1 正午高度角的特征检验太阳高度角在真太阳时正午达到峰值。对北半球中纬度的观测者正午高度角的理论值由90°-纬度赤纬角直接给出。以北京39.9042°N2024年夏至日为例| 节气 | 日期(近似) | 赤纬角(°) | 正午高度角理论值(°) | | 春分 | 3月20日 | 0 | 50.1 | | 夏至 | 6月21日 | 23.45 | 73.6 | | 秋分 | 9月23日 | 0 | 50.1 | | 冬至 | 12月22日 | -23.45 | 26.6 |跑完main.m后把夏至日曲线峰值与73.6°对比如果偏差超过0.5°优先检查时区tz与经度lon是否输错。这里有一类常见误用有人把经度代入正值但忘了东八区的tz8或者把时区当成0UTC而经度还是东经116°真太阳时直接偏晚8小时高度角曲线整体从正午向右平移峰值不再出现在12:00附近而出现在19:00附近这是最容易肉眼识别的故障。4.2 误差来源的量化拆解即使参数全部正确计算值与理想几何模型之间仍然存在系统性偏差。这套算法产出的是几何太阳位置不是视太阳位置两者的主要差异来自大气折射。高度角大于10°时折射造成的抬高量小于0.1°可以忽略但在日出日落附近折射会造成约0.5°以上的表观抬升工程上判断昼夜分界时一般把高度角修正为-0.833°再判零也就是把折射补偿近似固定为0.833°。时差方程与赤纬角的模型误差也会叠加入最终结果。Spencer模型里EoT的误差约在±0.5分钟以内对应时角误差约0.125°赤纬角0.02°的误差乘以cos(φ)的影响系数后最终高度角误差在0.02°量级。三者叠加后的综合误差不超过0.2°满足绝大多数光伏跟踪支架0.5°以内的转角控制需求。如果发现峰值偏差超过这个范围我一般先用行星历表抽查某个特定时刻再回过来核对公式里的经纬度符号。4.3 日出日落时刻与负高度角过滤日出日落时间的解析解可以直接从高度角方程反推更简单的做法是在已有曲线上做阈值检测。对于10分钟步长线性插值就能把边界时刻误差控制在几分钟内如果步长是60分钟日出日落附近曲线曲率高线性插值会产生较大偏差这时就要用fzero做精确求根。% 精确计算日出时刻在高度角定义中嵌入折射补偿-0.833° refraction_offset -0.833 * pi/180; % 定义高度角函数输入为儒略日输出为高度角-补偿值 f (t) sin_alt_function(t, lat, lon, tz) - sin(refraction_offset); % 日出在0~12点之间搜索 sunrise fzero(f, [t_start, t_start 0.5]); % 日落在12~24点之间搜索 sunset fzero(f, [t_start 0.5, t_stop]);fzero需要函数在区间两端异号所以日出区间必须选在高度角从负到正的范围内日落区间选在正到负。这里的sin_alt_function复用solar_position函数只需把datenum分解出时分秒再调用。对比几何法按0°判零与折射修正法按-0.833°判零后者得到的白昼时长在春秋季平均长4-6分钟这在实际光热电站的运行策略里会直接影响早晨吸热器启动时刻的选择。5. 进阶从视日轨迹到双轴跟踪与全年辐照评估视日轨迹计算的最终价值在于驱动机构和评估产能。这一章给出两个最直接的进阶用法。5.1 双轴支架转角解算对双轴跟踪系统最常见的结构是高度角-方位角型AZ-EL和俯仰-滚转型Tilt-Roll。前者直接使用本算法输出的alt和az作为两个转轴的指令角后者需要先转换到天顶角θz 90° - alt再把天顶角分解到支架的倾斜轴和旋转轴上。需要注意视日轨迹跟踪是纯几何位置跟踪没有考虑阴天散射辐照的优势方向因此云量较大地区应考虑低成本时控策略代替全时跟踪以节省驱动能耗。% AZ-EL双轴转角指令 theta_z 90 - alt; % 天顶角 azimuth_cmd az; % 旋转轴指令 % 检查角度变化率防止跳变 delta_az [0; diff(az)]; delta_az mod(delta_az 180, 360) - 180; % 折回处理delta_az的处理是为了应对方位角从355°跨到5°时差分出现的大幅跳变。直接用diff会得到一个约-350°的错误变化率折回处理后变成10°这个值才是伺服电机的真实角速度需求。省掉这一步在仿真报告里会出现瞬时上千度每秒的转速尖刺实际电机选型时会误判。5.2 全年逐时仿真与年发电量粗估把main.m的单日计算扩展到全年可以得到365×24的位置矩阵。结合一个简单的晴空辐照模型例如Hottel模型就能在没装传感器的情况下估算双轴跟踪相对固定角度安装的发电增益。% 全年逐时太阳位置矩阵 days 1:365; hours_of_day 0:23; alt_matrix zeros(length(days), 24); az_matrix zeros(length(days), 24); for d 1:length(days) [y, mo, da] datevec(datenum(year, 1, 1) d - 1); for h 1:24 [alt_matrix(d,h), az_matrix(d,h)] solar_position(... lat, lon, y, mo, da, h, 0, tz); end end % 过滤掉高度角小于5°的低光照时段 valid alt_matrix 5;这个矩阵的生产成本在Matlab里约在0.5秒以下。后续处理时把晴空直射辐照IDNI乘以cos(入射角)再按valid掩码累加就能得到月均发电量相对值。相比用PVsyst做全年仿真这套自研流程胜在可控和透明——每个系数、每个滤波条件和辐照模型都能在论文或者技术报告里写清楚对于硕士课题和预研项目已经足够。矩阵里的高度角变化率也可以顺带统计出来用来校核支架电机的最大角速度需求这一步在招标技术参数表里经常作为硬性指标出现。本文还有配套的精品资源点击获取
RELATED — 相关阅读

相关资讯

LATEST — 最新资讯

最新发布

TODAY — 本日精选

新闻

WEEKLY — 本周精选

新闻

MONTHLY — 本月精选

新闻