企业官网建设流程全解析

风电项目的第一步，往往不是选风机，而是搞清楚场址的风到底怎么样。我在好几个前期评估项目里拿到过从测风塔导出的历史风力数据，Matlab里一读，缺失值、野值、塔影效应、采集器时钟漂移，各种问题扎堆出现。这篇内容我会从“拿到测风塔原始数据”开始，讲怎么用 Matlab 完成从导入、清洗、计算风功率密度到威布尔拟合、发电量估算的完整流程，把里面值得注意的坑和经验一并写清楚。

如果你正在做风能资源评估、可再生能源前期选址，或者只是拿到了一个测风塔数据集不知道如何下手，这篇内容就是按我实际工作的思路整理的，你可以直接照着跑。

1. 测风塔数据从采集器到 Matlab：先搞明白你手里这堆数是什么

1.1 典型的测风塔配置和输出结构

测风塔（也叫气象塔）是风资源评估最重要的数据来源。一个规范配置的测风塔，通常会在多个高度安装风速计和风向标，比如 10m、30m、50m、70m，甚至更高。除了风速风向，绝大多数塔还会配温度传感器、气压传感器和湿度传感器，这些数据直接关系到空气密度修正，后续计算风功率密度时根本绕不开。

数据采集器一般按 10 分钟的平均间隔记录一组数据。每 10 分钟一组，一天就是 144 组，一年下来一个通道大约 52560 个数据点。原始记录里通常不只有平均风速，还有最大风速、最小风速、标准方差，甚至每个 10 分钟窗口里的逐秒采样瞬时值。不同采集器的存储格式差别很大，有的是直接 CSV，有的是 NRG、Campbell 等采集器特有的格式，还有的会导出为 Excel 表格。

首次拿到数据时，我建议先建一个“数据地图”：

做完这张表，你就知道自己手里有哪些通道，哪些高度可用于风切变分析，哪些通道能用来算湍流强度。没有这张地图，后面每一步都会很被动。

1.2 从原创采集文件到可分析表的常规转换

不同采集器导出的文件，第一行往往是变量名，中间还可能穿插着注释行或设备日志。直接 readtable 经常会把时间列读成文本，风速列读成 cell，气压列带个“NaN”字符串，总之第一步不会太干净。

我的做法是分两步走：

1. 先用人工过一遍文件头，把设备日志、双标题行清理掉，只保留表头和数据行。
2. 然后用 readtable 统一读取，再用 detectImportOptions 把各列类型显式指定好。

% 读取测风塔CSV数据 opts = detectImportOptions('tower_data.csv'); % 显式指定各列类型 opts = setvartype(opts, {'time','WS70_avg','WD70_avg','WS50_avg','WD50_avg','Temp','Press'}, ... {'datetime','double','double','double','double','double','double'}); T = readtable('tower_data.csv', opts); % 检查读取结果 head(T)

有一个非常容易掉的坑：时间列在 Excel 里显示为“2023-05-01 00:10”，但读取进来之后格式可能是“01/05/2023 00:10”，如果直接用小时索引处理，顺序完全错乱。因此读取后立刻统一时间列格式，并给时间列排序：

T = sortrows(T, 'time'); T.time = dateshift(T.time, 'start', 'minute'); % 对齐到整分钟

很多数据采集器还会把没有风速记录的时间留成空格或 -9999，这一步顺手统一换成 NaN，方便后续清洗。做完这些，数据才算真正进到 Matlab 的“操作台”上。

1.3 时间戳连续性检查：采集器时钟漂移不得不防

测风塔在野外常年运行，太阳能供电、冬季低温，采集器的实时时钟经常出现漂移。有些塔几个月后会慢几分钟，甚至慢到半小时。这种漂移对于长期平均风速影响不大，但对逐时变化、逐日变化和极端风事件分析，会直接造成“时间错位”。

我从实际项目里总结了一个快速检查法：用 diff(T.time) 检查相邻两行时间差是否为 10 分钟，把所有非 10 分钟的间隔全列出来。

timeGap = minutes(diff(T.time)); badGap = find(abs(timeGap - 10) > 0.01); if isempty(badGap) disp('时间连续性正常'); else disp(['发现 ', num2str(length(badGap)), ' 处间隔异常']); disp(T.time(badGap(1:min(end,5)))); end

如果异常间隔的区域是小范围且只有一分钟级别的偏移，一般可以在后续清洗时做重采样对齐。如果漂移严重，就必须联系测风塔运维单位，调取现场维护日志，确认哪些时间段不可靠。这部分工作是数据代表性的底线，不查清楚，后面所有分析结果都会存疑。

2. 数据清洗实操：野值、缺测、塔影效应逐个干掉

2.1 风速负值、超量程和逻辑错误怎么判定

测风数据的清洗不能只靠“看着不对劲就删”。我一般先做三道硬过滤，再做一道逻辑过滤。

硬过滤：

• 风速小于 0，直接判定为野值。理论上风速不可能为负，出现负值往往是传感器信号短路或采集器故障。
• 风速大于传感器量程。大多数测风杯的极限量程是 60m/s 或 75m/s，超过量程上限的读数直接剔除。
• 风向小于 0 或大于 360，同样剔除。

逻辑过滤需要结合通道关系。比如 70m 高度风速低于 10m 高度风速且相差特别悬殊，不一定是错，但如果是稳定低于 20% 以上，很可能是低层传感器结冰或高层传感器轴承磨损，这类数据要单独标记，不能盲目删除。

下面这一段清洗代码可以直接套用：

% 硬过滤 T.WS70_avg(T.WS70_avg < 0 | T.WS70_avg > 60) = NaN; T.WS50_avg(T.WS50_avg < 0 | T.WS50_avg > 60) = NaN; T.WD70_avg(T.WD70_avg < 0 | T.WD70_avg > 360) = NaN; % 逻辑检查：70m 风速应为塔内最高，若低于 10m 风速 30% 以上，标记为可疑 idx_logic = (T.WS70_avg < 0.7 .* T.WS10_avg) & ~isnan(T.WS10_avg); T.WS70_avg(idx_logic) = NaN;

不过要注意，逻辑过滤不能做得太激进。强逆温层结、低空急流、局地山谷风环流，都可能让上层风速低于下层，这种情况虽然少见但不是错误。我在项目里一般只标记不删除，后续通过时间切片确认是否为天气系统造成的正常现象。

2.2 卡死值和跳变值的识别技巧

测风杯轴承卡死是野外塔最常见的机械故障之一。卡死时风速通道会持续输出同一个值，持续几小时甚至几天，乍一看曲线很正常，但标准差接近 0。判断卡死值的方法很简单：计算连续相同值的持续时间，超过 30 分钟（3 个连续采样点）就该怀疑了。

% 找出连续卡死 3 个点以上的区间 runLen = 1; badRun = false(height(T), 1); for i = 2:height(T) if T.WS70_avg(i) == T.WS70_avg(i-1) runLen = runLen + 1; else if runLen >= 3 badRun(i-runLen:i-1) = true; end runLen = 1; end end % 卡死区间置为 NaN T.WS70_avg(badRun) = NaN;

跳变值则是相邻两个点风速突然变化超过一定阈值。10 分钟平均风速一般不会瞬间变化超过 5m/s 到 8m/s，除非是阵风锋面或者雷暴出流，这才可能发生剧烈突变。所以跳变值也不能一删了之，我通常是先把突变点标记出来，再结合气压、温度变化判断是否属于真实天气过程。

一个比较实用的替代思路是：先画全年的风速时间序列散点图，把肉眼可见的“尖刺”区间放大看，确认是孤立点还是连续变化。孤立点直接剔除，连续变化则保留并记录说明。

2.3 塔影效应和表头朝向剔除：风向数据里藏着的系统性偏差

测风塔的塔身、拉线、横臂，都会对气流产生扰动。当风向吹向塔体的某一侧时，塔身下游的尾流会让测风杯的读数系统性偏低，这就是“塔影效应”。如果直接把所有风向扇区的风速都纳入统计，平均风速会被轻微拉低，风功率密度会被低估。

一条最常规的做法是：根据测风塔安装日志里的传感器朝向，把风向在塔身方位 ±30° 范围内的数据剔除。比如测风塔的横臂朝向为 90°（正东），当风向在 60°–120° 区间时，风速数据不可用。

% 假设传感器朝向为 90° 和 270° 两侧，塔身两侧 60°–120° 与 240°–300° 为受扰扇区 towerShadowSec = [60 120; 240 300]; isShadow = (T.WD70_avg >= towerShadowSec(1,1) & T.WD70_avg <= towerShadowSec(1,2)) | ... (T.WD70_avg >= towerShadowSec(2,1) & T.WD70_avg <= towerShadowSec(2,2)); T.WS70_avg(isShadow) = NaN;

前面的逻辑过滤会引发一个问题：大量数据被置成 NaN 后，剩余数据量可能不够。行业里对有效数据率的要求一般是 90% 以上，低于这个水平，年度代表性就要打问号。所以塔影剔除法务必要记录剔除比例，如果超过 25%，说明塔的布局本身有问题，不能硬着头皮继续分析。

2.4 缺测数据的插补：短期线性插值、长期需要相关回归

插补策略取决于缺失数据的时长。缺 1–2 个 10 分钟点，直接用前后线性插值完全可以，误差很小。缺 2–3 小时，线性插值还行，但要谨慎。缺一整天甚至更久，唯一的可靠办法是用邻近高度的风速相关性做回归，或者用测风塔附近长期气象站数据做拟合。

% 短时缺测：线性插值 T.WS70_avg = fillmissing(T.WS70_avg, 'linear', 'SamplePoints', T.time); % 长时间缺测：如果 70m 缺得太多，用 50m 回归到 70m 的幂律关系 validIdx = ~isnan(T.WS70_avg) & ~isnan(T.WS50_avg); p = polyfit(log(T.WS50_avg(validIdx)), log(T.WS70_avg(validIdx)), 1); T.WS70_avg(isnan(T.WS70_avg)) = exp(polyval(p, log(T.WS50_avg(isnan(T.WS70_avg)))));

需要特别留意：插补数据不代表真实测量值，在最终报告中必须单独标注插补比例。能作为依据的气象标准里，插补数据比例过高时，风功率密度的置信度会明显下降，这一点在后面的不确定性分析里会继续使用。

3. 平均风速、风功率密度与湍流强度：三个核心指标一次算清

3.1 逐时平均风速与有效样本数量

风速的统计基础是 10 分钟平均序列，但做资源评估时，通常还会按小时、月、年三个尺度统计。小时平均可以滤掉阵风脉动的短时干扰，月平均可以看到季节变化，年平均则是场址长期资源水平的基本盘。

% 按小时平均 T.hourly = dateshift(T.time, 'start', 'hour'); hourlyWS = groupsummary(T, 'hourly', 'mean', 'WS70_avg'); % 年平均 annualWS = groupsummary(T, 'year', 'mean', 'WS70_avg');

这里有个非常关键的原理：平均风速必须在“小时风速”和“10分钟风速算术平均”之间选择正确的口径。规范做法是用 10 分钟平均风速序列直接参与后续风功率密度计算，因为风功率密度与风速立方成正比，而风速立方平均不等于平均风速的立方。要是直接用年平均风速去估算发电量，结果会严重偏低。

3.2 风功率密度公式与空气密度修正

风功率密度表示单位面积上气流的功率，公式是：

[ P = \frac{1}{2} \rho V^3 ]

实际计算时，用 10 分钟平均风速序列逐个点求风速三次方，再取平均，最后乘上 0.5 和空气密度。Matlab 写法：

% 计算空气密度：根据实测气压和温度 Rd = 287.05; % 干空气气体常数 J/(kg·K) T.K = T.Temp + 273.15; T.rho = T.Press ./ (Rd .* T.K); % Press 单位 Pa rho_mean = mean(T.rho, 'omitnan'); % 风功率密度 W/m2 P_wd = 0.5 * rho_mean * mean(T.WS70_avg.^3, 'omitnan');

空气密度往往是新手会忽略的地方。标准空气密度 1.225 kg/m³ 只适用于海平面、15°C 的标准情况。高原风电场空气密度可能只有 1.05 甚至更低，直接套 1.225，风功率密度会被高估 10% 以上。如果你手上的测风塔没有气压传感器，可以用海拔粗略估算：

[ \rho \approx 1.225 \times e^{-0.0001 \times z} ]

其中 z 是海拔，单位米。这个近似公式在 3000m 以下误差可以控制在 2%–3% 左右，前期预评估够用，施工图阶段必须要实测气压。

3.3 湍流强度计算：直接决定风机等级选型

湍流强度反映风速脉动的剧烈程度，定义为标准差除以平均风速：

[ I = \frac{\sigma V}{\overline{V}} ]

计算时需要 10 分钟窗口内的高频采样记录，或者采集器输出的风速标准差通道。如果没有标准差通道，只能用相邻 10 分钟平均值的变差来近似，但这会低估真实湍流强度，不太推荐。

% 用采集器输出的标准差通道计算湍流强度 T.I70 = T.WS70_std ./ T.WS70_avg; T.I70(T.I70 > 0.6) = NaN; % 明显错误值剔除 I70_mean = mean(T.I70, 'omitnan');

为什么要专门算湍流强度？因为风机选型时，IEC 等级直接和湍流强度挂钩，比如 IEC 1A 类风机要求 15m/s 风速下的湍流强度不超过 0.16。如果场址湍流强度太高，风机载荷安全性不达标，要么降低等级选型，要么就得避开高湍流区域排布机位。这个参数对风机排布方案的影响比很多人想象中大得多。

4. 威布尔分布拟合与风玫瑰图：判断风资源质量的两把尺子

4.1 威布尔分布参数估计：最小二乘与极大似然

风速频率分布是描述风资源最重要的统计特征。威布尔分布有两个参数：形状参数 k 和尺度参数 c。k 值决定了风速分布的“尖陡”程度，k 越小，风速分布越矮胖，说明风速变化范围大；k 越大，风速越集中，比如 2–3 左右。c 值近似代表年平均风速换算后的水平。

实际拟合思路：

• 把所有有效风速从小到大排序，用中位秩公式计算经验累积频率。
• 对威布尔累积分布函数做线性化变换，用最小二乘拟合直线，斜率和截距反推 k 和 c。

v = sort(T.WS70_avg(~isnan(T.WS70_avg))); n = length(v); F = (1:n)' / (n + 1); % 中位秩 x = log(v); y = log(-log(1 - F)); p = polyfit(x, y, 1); k_fit = p(1); c_fit = exp(-p(2) / p(1));

如果想省事，Matlab 自带的 wblfit 可以直接做极大似然估计：

params = wblfit(v); c_mle = params(1); k_mle = params(2);

两种方法各有侧重，最小二乘比较直观可控，极大似然统计效率更高。我建议两个都算一遍对比，如果差值超过 5%，说明数据里还有未清理的异常段或者双峰分布（比如季风区冬夏风向差异极大），这时候单纯拟合单一威布尔分布就不够精确，要考虑分季节或分风向扇区拟合。

4.2 风玫瑰图：主导风向比平均风速更能说明问题

风玫瑰图的核心是展示不同方向的频率和风速贡献。选址阶段主要关注两点：主风向是否集中，主风向上的平均风速是否够高。如果主风向是东南风，但强风速全在西北方向，机位排布和道路设计就得有不同的侧重。

Matlab 里用 polarhistogram 和 polaraxes 就能画：

dir_rad = deg2rad(T.WD70_avg(~isnan(T.WD70_avg))); figure; polaraxes; polarhistogram(dir_rad, 16, 'Normalization', 'probability'); title('70m 高度风玫瑰图');

实际项目里我一般画两张图：一张是风向频率玫瑰，一张是风速贡献玫瑰。频率玫瑰告诉你哪边风多，风速贡献玫瑰告诉你哪边风带来的能量多，后者对排布的影响更直接。很多风电场排布方案失败，就是只看了频率玫瑰，没看能量玫瑰。

4.3 风切变指数与轮毂高度风速推算

测风塔各高度层的平均风速可以用来计算风切变指数 α：

[ \alpha = \frac{\ln(V_2 / V_1)}{\ln(z_2 / z_1)} ]

alpha = log(mean(T.WS70_avg,'omitnan') / mean(T.WS30_avg,'omitnan')) / log(70 / 30);

α 的值直接决定从已知高度推算轮毂高度风速的可信度。平坦开阔地形 α 大约在 0.10–0.20 之间，复杂山地可达 0.3 甚至更高。如果测风塔最高层是 70m，而风机轮毂高度是 100m，就需要用 α 外推。

hub_h = 100; ref_h = 70; V_hub = mean(T.WS70_avg, 'omitnan') * (hub_h / ref_h)^alpha;

需要警示一点：α 不是恒定常数，它会随风速区间变化，夜间稳定层结下 α 可能异常偏大。所以在做外推时，最好分风速段统计 α，再取轮毂高度对应风速段的值，而不是用全年平均 α 一把梭。否则 100m 高度风速容易被高估。

5. 发电量估算与机位排布：资源评估最终要回答投资方的问题

5.1 功率曲线法的 AEP 理论计算

资源评估做完了，最终要落到“这个场址一年能发多少电”上。最常用的方法是功率曲线法：把风速频率分布（威布尔分布）与风机功率曲线逐风速点相乘，再乘以全年小时数。

% 假设功率曲线以表格格式存储，列1为风速，列2为功率kW powerCurve = readtable('turbine_power_curve.csv'); v_pc = powerCurve.V; P_pc = powerCurve.kW; % 用威布尔概率密度函数生成风速分布 v_bin = 0:0.5:30; f_v = wblpdf(v_bin, c_fit, k_fit); % 确保功率曲线覆盖 v_bin 范围 P_interp = interp1(v_pc, P_pc, v_bin, 'linear', 0); % 计算年发电量 (单位 kWh) AEP_theoretical = 8760 * sum(f_v .* P_interp);

注意两个细节。第一，功率曲线的风速一般是标准空气密度下的值，如果场址空气密度低于 1.225，需要做密度修正。工程上常见的做法是按实际密度与标准密度的比值等比修正功率曲线，密度越低，同等风速下功率越低。第二，风机的切入切出风速区间会直接影响有效发电小时数，如果风速分布的高风速占比大，切入切出后的截断效应更明显。

5.2 综合折减系数：把理论值拉回现实

理论 AEP 是“风机全年无故障、无停机、无尾流损失、无控制损失、电网永远消纳”的理想值。现实里必须乘一个折减系数。我常用的一组参考值是：

综合下来，年净发电量大约要乘 0.75–0.85。如果某个项目算出来的综合折减超过 90%，要么你的可利用率数据异常优秀，要么你漏了某个大项，要回去重新核对。

netFactor = 0.95 * 0.90 * 0.97 * 0.98 * 0.98 * 0.90; AEP_net = AEP_theoretical * netFactor;

5.3 从资源评估到机位排布：两个容易被忽视的图

机位排布图纸上通常看两样东西，一个是风速分布色块图，一个是湍流强度色块图。风速高的区域优先布置风机，这是直觉；湍流强度高的区域尽量避开，这是安全。如果高湍流区和高风速区重叠，就需要逐机位做载荷安全评估，不能拍脑袋决定。

我处理过的项目中，出现过两个典型情况。第一个，场址西南角年平均风速高出东部 15%，但西南角紧邻断崖，湍流强度超过 0.25，最终那排机位全部后移了 200m，发电量下降但载荷安全性达标。第二个，某个平坦湖岸场址主风向与岸线平行，排布时按正方形网格布置，结果尾流损失达到 14%，改成顺主风向错列排布后降到 8%，年发电量提升明显。这些都说明，资源评估的最终交付物不只是报告，还包括真正能指导排布的分析图和数据表。

行业里的成熟做法是用 CFD 或尾流模型做多轮迭代优化，但在没有商业软件的情况下，用 Matlab 基于威布尔分布、风速分布和功率曲线做一版初步的静态估算，已经能给项目可行性判断提供足够的方向性参考。

我在实际项目里最大的体会是：测风塔数据的处理流程尽量标准化，清洗规则、插补方法、塔影扇区、折减系数全都用脚本固化下来。这样每次拿到新项目的数据，只需改文件路径和塔结构参数，就能快速产出结果。风电项目的周期往往很短，前期评估要快、要准、要可追溯，一套可复用的 Matlab 处理流程比临时写脚本省下不止一半时间。

如果你手头也有一套测风塔数据，建议先按第 1 节的方法建好数据地图，再按清洗、插补、指标计算、分布的流程走一遍。遇到什么奇怪的数据现象，欢迎来交流。