使用 WRF (netCDF) 数据优化阵列平面 (POA) 辐照度计算

Question

我需要使用 python 的 pvlib 包 (https://pvlib-python.readthedocs.io/en/stable/) 计算阵列平面 (POA) 辐照度。为此，我想使用 WRF 模型（GHI、DNI、DHI）的输出数据。输出数据为 netCDF 格式，我使用 netCDF4 包打开，然后使用 wrf-python 包提取必要的变量。

这样我得到了一个 xarray.Dataset ，其中包含我将使用的变量。然后我使用 xarray.Dataset.to_dataframe() 方法将其转换为 pandas 数据帧，然后使用 dataframe.values 将数据帧转换为 numpy 数组。然后我做一个循环，在每次迭代中，我使用函数 irradiance.get_total_irradiance (https://pvlib-python.readthedocs.io/en/stable/auto_examples/plot_ghi_transposition.html) 作为网格点计算 POA。

到目前为止，我就是这样做的，但是我在 WRF 域中有超过 160000 个网格点，数据是每小时一次，跨越 365 天。这提供了非常大量的数据。我相信如果 pvlib 可以直接与 xarray.dataset 一起使用它可能会更快。但是，我只能这样做，将数据转换为 numpy.array 并循环遍历行。谁能告诉我如何优化这个计算？因为我开发的代码非常耗时。

如果有人可以帮助我，我将不胜感激。也许是对代码的改进，或者是从 WRF 数据计算 POA 的另一种方法......

我正在提供我目前构建的代码：

from pvlib import location
from pvlib import irradiance

import os

import matplotlib.pyplot as plt
import numpy as np
import pandas as pd
import xarray as xr
import netCDF4
import wrf

获取 WRF 数据

variaveis = ['T2',
             'U10',
             'V10',
             'SWDDNI',
             'SWDDIF',
             'SWDOWN']

netcdf_data = netCDF4.Dataset('wrfout_d02_2003-11-01_00_00_00')

first = True

for v in variaveis:

    var = wrf.getvar(netcdf_data, v, timeidx=wrf.ALL_TIMES)
    
    if first:
        met_data = var
        first = False
    else:
        met_data = xr.merge([met_data, var])

met_data = xr.Dataset.reset_coords(met_data, ['XTIME'], drop=True)
met_data['T2'] = met_data['T2'] - 273.15

WS10 = (met_data['U10']**2 + met_data['V10']**2)**0.5
met_data['WS10'] = WS10

df = met_data[['SWDDIF', 
               'SWDDNI', 
               'SWDOWN', 
               'T2', 
               'WS10']].to_dataframe().reset_index().drop(columns=['south_north', 
                                                                   'west_east'])

df.rename(columns={'SWDOWN': 'ghi',
                   'SWDDNI':'dni', 
                   'SWDDIF':'dhi', 
                   'T2':'temp_air', 
                   'WS10':'wind_speed',
                   'XLAT': 'lat',
                   'XLONG': 'lon',
                   'Time': 'time'}, inplace=True)
df.set_index(['time'], inplace=True)

df = df[df.ghi>0]
df.index = df.index.tz_localize('America/Recife')

获取 POA 辐照度的函数

def get_POA_irradiance(lon, lat, date, dni, dhi, ghi, tilt=10, surface_azimuth=0):

    site_location = location.Location(lat, lon, tz='America/Recife')

    # Get solar azimuth and zenith to pass to the transposition function
    solar_position = site_location.get_solarposition(times=date)
    
    # Use the get_total_irradiance function to transpose the GHI to POA
    POA_irradiance = irradiance.get_total_irradiance(
        surface_tilt = tilt,
        surface_azimuth = surface_azimuth,
        dni = dni,
        ghi = ghi,
        dhi = dhi,
        solar_zenith = solar_position['apparent_zenith'],
        solar_azimuth = solar_position['azimuth'])
    
    # Return DataFrame with only GHI and POA
    
    return pd.DataFrame({'lon': lon,
                         'lat': lat,
                         'GHI': ghi,
                         'POA': POA_irradiance['poa_global']}, index=[date])

在数组的每一行（时间）循环

array = df.reset_index().values
    
list_poa = []
    
def loop_POA():   
    for i in tqdm(range(len(array) - 1)):
        POA = get_POA_irradiance(lon=array[i,6], 
                                 lat=array[i,7], 
                                 dni=array[i,2], 
                                 dhi=array[i,1], 
                                 ghi=array[i,3], 
                                 date=str(array[i,0]))
        list_poa.append(POA)
    
    return list_poa

poa_final = pd.concat(lista)

Answer 1

感谢您提出一个好问题和使用 pvlib！没错，pvlib 旨在对单个位置进行建模，而不是为与 xarray 数据集一起使用而设计的，尽管某些函数可能会巧合地与它们一起使用。

我强烈怀疑您看到的大部分运行时间都是用于太阳位置计算。您可以切换到更快的方法（请参阅method 选项here），因为默认的太阳位置方法非常准确，但在计算批量位置时也很慢。安装 numba 会有所帮助，但它对你来说仍然可能太慢，所以你可能会检查其他模型（星历，pyephem）。还有一些快速但低精度的方法，但您需要稍微更改代码才能使用它们。请参阅“低精度太阳位置计算的相关性和解析表达式”下的列表here。

就像 Michael Delgado 在 cmets 中建议的那样，并行处理是一种选择。但这在 python 中可能是一个令人头疼的问题。您可能需要多处理，而不是多线程。

另一个想法是使用atlite，这是一个专为这种空间建模而设计的python包。但它的太阳能建模功能不如 pvlib 详细，因此它可能对您的情况没有用处。

另一个注意事项：我不知道 WRF 数据是区间平均值还是瞬时值，但如果您关心准确性，则应该以不同的方式处理它们以进行转置。看到这个example。

编辑添加：再次查看您的代码后，可能会有另一个显着的加速。您是否为位置和时间戳的单一组合致电get_POA_irradiance？如果是这样，那是不必要的并且非常缓慢。传递每个位置的完整时间序列会快得多，即标量纬度/经度但矢量辐照度。