Homec4science
Diffusion solar_potential (master) Browse

solarPOT/Potential_hourly/EW_scenarios/compute_tilted_irrad_batch_flat_EW.py
67c48c8e30edmaster

Last Change
Raw File

compute_tilted_irrad_batch_flat_EW.py
View Options

import pandas as pd
import xarray as xr
import numpy as np
from pvlib.tools import cosd
import os
import time
import sys
def read_csv_rows( fp, start, batch_size, usecols = None, OFFSET = 0 ):
field_names = pd.read_csv( fp, nrows = 0).columns
return pd.read_csv(fp, skiprows = start + OFFSET, nrows = batch_size, header = None, names = field_names, usecols = usecols)
## ===== INPUTS ========== ##
POA_IRRAD_FP = sys.argv[1] # name of the tilted radiation file
ROOFS_FP = sys.argv[2] # name of the file with roof characteristics
OUTFILE_FP = sys.argv[3] # name of the output file
## ===== CONSTANTS ===== ##
CHUNKSIZE = 10000
FLAT_TILT = 15
FLAT_DIST = 1
FLAT_DEF = 10 # tilt angle below which roofs are considered flat (default 0)
PARTLY_SHADING_FP = '/work/hyenergy/output/solar_potential/geographic_potential/available_area/shading_stats/CH_partly_shading_replaced_debiased.nc'
SVF_FP = '/work/hyenergy/output/solar_potential/geographic_potential/tilted_irradiance/SVF_CH_replaced.csv'
MUTUAL_SHADE_FP = '/work/hyenergy/output/solar_potential/technical_potential/flat_roofs/mutual_shading.nc'
print('\nShading file path: %s' %PARTLY_SHADING_FP)
print('POA irradiance file path: %s' %POA_IRRAD_FP)
print('Skyview factor file path: %s' %SVF_FP)
print('Mutual shading file path: %s' %MUTUAL_SHADE_FP)
print('Roof surfaces file path: %s' %ROOFS_FP)
print('Output file path: %s' %OUTFILE_FP)
## =========================================================================================== ##
## ===== LOAD DATA ===== ##
tt = time.time()
t0 = time.time()
irrad = xr.open_mfdataset( POA_IRRAD_FP, chunks={'DF_UID': CHUNKSIZE},
combine = 'by_coords', lock=False ).sortby('DF_UID')
shading = xr.open_mfdataset( PARTLY_SHADING_FP, chunks={'DF_UID': CHUNKSIZE},
combine = 'by_coords', lock=False ).sel(DF_UID = irrad.DF_UID.values)
svf = pd.read_csv( SVF_FP )
mutual_shade = xr.open_dataset( MUTUAL_SHADE_FP )
roofs = pd.read_csv( ROOFS_FP, usecols = ['DF_UID', 'NEIGUNG', 'panel_tilt'] )
print('Loaded all data in %.3f seconds' %(time.time()-tt))
flat_roofs = roofs[(roofs.DF_UID.isin(irrad.DF_UID.values)) & (roofs.NEIGUNG <= FLAT_DEF)]
## ===== COMPUTE MEAN IRRADIANCE ON ROOFTOPS ===== ##
print('\nComputing mean irradiance on rooftops ...')
tt = time.time()
tilted_irrad = xr.merge( [ shading, irrad.drop('time') ] ).fillna(0)
mutual_shade_flat = mutual_shade.sel(tilt = FLAT_TILT,
distance = FLAT_DIST,
timestamp = irrad.timestamp[irrad.timestamp.isin(mutual_shade.timestamp)]
).drop(['tilt', 'distance'])
tilted_irrad['svf'] = svf[svf.DF_UID.isin(irrad.DF_UID.values)].set_index('DF_UID')['SVF_unbiased'].to_xarray()
tilted_irrad['is_flat'] = tilted_irrad.DF_UID.isin( flat_roofs.DF_UID.values ).astype(int)
# compute tilted irradiance for tilted and flat surfaces
# by multiplying mutual shading with is_flat (0, 1)
tilted_irrad['tilted_irradiance'] = ( tilted_irrad.partly_shading * tilted_irrad.poa_direct
* (1 - mutual_shade_flat.shading * tilted_irrad.is_flat)
+ tilted_irrad.svf * tilted_irrad.poa_sky_diffuse
* (1 - mutual_shade_flat.SVF_shade * tilted_irrad.is_flat)
+ tilted_irrad.poa_ground_diffuse )
print(tilted_irrad)
print('Finished in %.3f seconds' %(time.time()-tt))
# Compute yearly sum
tilted_irrad['yearly_kWh_m2'] = tilted_irrad.tilted_irradiance.sum( dim = ['timestamp','day'] ) / 1000
## === SAVE OUTPUT FILE ==== ##
tt = time.time()
print('\nSaving output ...')
irrad_save = tilted_irrad[[ 'tilted_irradiance', 'yearly_kWh_m2' ]]
# reduce storage requirement
irrad_save['tilted_irradiance'] = irrad_save['tilted_irradiance'].astype('float32')
irrad_save['yearly_kWh_m2' ] = irrad_save['yearly_kWh_m2' ].astype('float32')
irrad_save['time'] = irrad.time
print(irrad_save)
irrad_save.to_netcdf( OUTFILE_FP )
print('Saved output in %s (%.3fs)' %( OUTFILE_FP, time.time()-tt ))
irrad.close()
shading.close()
print('\nDONE - total execution time: %.3f seconds' %(time.time()-t0))
c4science · Help