python - 使用 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 在评论中建议的那样，并行处理是一种选择。但这在python中可能会令人头疼。您可能需要多处理，而不是多线程。
另一个想法是使用 atlite ，一个专为这种空间建模设计的python包。但是它的太阳能建模功能远不如 pvlib 详细，因此它可能对您的情况没有用处。
另一个注意事项:我不知道 WRF 数据是间隔平均值还是瞬时值，但如果您关心准确性，则应该以不同的方式处理它们以进行换位。看到这个 example .
编辑添加:再次查看您的代码后，可能会有另一个显着的加速。您在打电话吗get_POA_irradiance对于位置和时间戳的单一组合？如果是这样，那是不必要的，而且非常慢。传递每个位置的完整时间序列会快得多，即标量纬度/经度但矢量辐照度。