简介面向气象、环境科学与太阳能资源评估方向研究者的Python工具围绕全国840个气象站点日照时数SSD转日总太阳辐射SR这一常见数据处理任务整理。压缩包共3个文件由2个CSV数据表和1个Python脚本组成整体约275KBCSV分别承载站点输入数据与计算结果脚本负责读取、清洗、公式换算和批量输出结构紧凑便于套用到同类观测数据。已有6348人学习。脚本以经纬度、日期和日照时数为基础结合林格曼方程或经验公式估算日总太阳辐射并在计算前处理缺失值、异常值与日期格式不一致等问题减少手工整理成本。对科研人员它提供了可复核的批量计算模板可替换站点数据与参数后用于气候变化、能源利用和气候模型的前期分析对数据分析学习者则可借pandas读写CSV、逐站计算与结果落盘的完整流程理解气象指标换算和脚本组织方式也方便在此基础上扩展统计汇总或可视化。1. 840 个站点的日照时数摆在一起先别急着套公式前几年接手一个跨区域蒸散估算的活手头只有一张 840 个气象站点的日值表站号、日期、日照时数三列。第一反应是 Ångström-Prescotta 取 0.25、b 取 0.50脚本半小时写完跑完一看高纬那几个站夏季算出来的日总太阳辐射比天文辐射还高还有一批站的辐射与天文辐射之比全年是一条水平线。公式本身没错问题出在输入变量 n 上单位是 0.1 小时没有除以 10缺测码 32766 被当成了 3276.6 小时可照时数 N 这一项干脆被漏掉了。这套 840 站点日照时数转日总太阳辐射的 python 脚本真正价值不在公式而在于把单位、缺测、坐标、可照时数这几道坎提前处理掉了。凡是手上有气象站日值、要做蒸散、作物模型或者光伏资源初筛的都可以照着往下走。2. 日照时数怎么变成 MJ/m²Ångström-Prescott 的四个输入量Ångström-Prescott 关系本身只有一行Rs/Ra a b·(n/N)。很多人就卡在这行式子里以为把 n 代进去、乘一个系数就完了。实际上这个式子有四个输入量n 只是其中之一另外三个——天文辐射 Ra、可照时数 N、经验系数 a 和 b——每一个都能独立把结果带偏。下面按“先算天、再定系数、最后对齐数据”的顺序拆开。2.1 天文辐射 Ra 与可照时数 N 的逐日算法Ra 是大气层顶水平面上一天的累积辐射量只跟纬度和年积日有关跟云、气溶胶、海拔全都没关系。它之所以重要是因为它给了 Rs 一个物理上限地面日总辐射不可能显著超过 Ra。任何一站的 Rs/Ra 超过 1一定是数据或代码出了问题这就是后面排查时最有用的一条硬性限制。计算链是这样的先由年积日得到日地距离订正系数 dr 和太阳赤纬 δ再由纬度和赤纬算出日落时角 ωs最后把这三者代进 Ra 的积分式。可照时数 N 则是 ωs 换算成小时的结果ωs 是弧度乘 24/π 才是小时。这一步是第一个翻车点很多人把 ωs 直接当小时用N 就只有 1.5 左右于是 n/N 动不动大于 1。import numpy as np import pandas as pd GSC 0.0820 # 太阳常数, MJ m-2 min-1 DAYS_IN_YEAR 365.0 # 固定 365, 与通用手册口径保持一致 def ra_and_n(doy, lat_deg): doy : 年积日, 标量或 ndarray, 1-366 lat_deg : 纬度, 北纬为正, 十进制度 返回 : Ra (MJ m-2 d-1), N (h) doy np.asarray(doy, dtypefloat) phi np.radians(np.asarray(lat_deg, dtypefloat)) dr 1.0 0.033 * np.cos(2 * np.pi * doy / DAYS_IN_YEAR) # 日地距离订正 delta 0.409 * np.sin(2 * np.pi * doy / DAYS_IN_YEAR - 1.39) # 太阳赤纬, rad x -np.tan(phi) * np.tan(delta) ws np.arccos(np.clip(x, -1.0, 1.0)) # 日落时角, rad Ra (24 * 60 / np.pi) * GSC * dr * ( ws * np.sin(phi) * np.sin(delta) np.cos(phi) * np.cos(delta) * np.sin(ws) ) N 24.0 / np.pi * ws # 可照时数, h return Ra, N逐项说明一下。dr的振幅只有 ±3.3%看着不起眼但它决定了近日点和远日点的差异在年尺度统计上会被放大。delta用的是通用近似式误差在 0.1° 量级对 Ra 的影响小于 0.2%够用。np.clip(x, -1, 1)是防止 arccos 定义域越界的保底措施中国陆域纬度在 3°N 到 54°N 之间实际不会出现极昼极夜但站点表里一旦混进纬度 90 或者符号写错的记录np.tan会让 x 跑到几千从而得到 NaN而 NaN 会一路传到结果表里很难追查。先把它限制在定义域内再单独做一轮纬度范围检查比事后查找 NaN 要省事。参数方面doy直接用dt.dayofyear闰年自动是 366但分母我固定用 365。这不是疏忽通用手册的公式就是以 365 为分母换成 366 会让你的结果和别人的结果对不上差异在 0.1% 量级没有意义。lat_deg必须是十进制度北纬为正南纬为负如果你的站点表用度分整数存储先别急着往里扔。2.2 系数 a、b默认值本文还有配套的精品资源点击获取
