shadow (version 0.7.1)

radiation: Estimation of Direct and Diffuse Radiation Load on Extruded Polygon Surfaces

Description

This is a wrapper function for calculating total diffuse, direct and total radiation load per unit area on extruded polygon surfaces. The function operates on obstacle geometry and a set of sun positions with associated meteorological estimates for direct and diffuse radiation (see Details below).

Usage

radiation(
  grid,
  obstacles,
  obstacles_height_field,
  solar_pos = solarpos2(obstacles, time),
  time = NULL,
  solar_normal,
  solar_diffuse,
  radius = Inf,
  returnList = FALSE,
  parallel = getOption("mc.cores")
)

Arguments

grid

A 3D SpatialPointsDataFrame layer, such as returned by function surfaceGrid, specifying the locations where radiation is to be estimated. The layer must include an attribute named type, with possible values being "roof" or "facade", expressing surface orientation per 3D point. The layer must also include an attribute named facade_az, specifying facade azimuth (only for "facade" points, for "roof" points the value should be NA). The type and facade_az attributes are automatically created when creating the grid with the surfaceGrid function

obstacles

A SpatialPolygonsDataFrame object specifying the obstacles outline, inducing self- and mutual-shading on the grid points

obstacles_height_field

Name of attribute in obstacles with extrusion height for each feature

solar_pos

A matrix with two columns representing sun position(s); first column is the solar azimuth (in decimal degrees from North), second column is sun elevation (in decimal degrees); rows represent different sun positions corresponding to the solar_normal and the solar_diffuse estimates. For example, if solar_normal and solar_diffuse refer to hourly measurements in a Typical Meteorological Year (TMY) dataset, then solar_pos needs to contain the corresponding hourly sun positions

time

When solar_pos is unspecified, time can be passed to automatically calculate solar_pos based on the time and the centroid of obstacles, using function maptools::solarpos. In such case obstacles must have a defined CRS (not NA). The time value must be a POSIXct or POSIXlt object

solar_normal

Direct Normal Irradiance (e.g. in Wh/m^2), at sun positions corresponding to solar_pos. Must be a vector with the same number of elements as the number of rows in solar_pos

solar_diffuse

Diffuse Horizontal Irradiance (e.g. in Wh/m^2), at sun positions corresponding to solar_pos. Must be a vector with the same number of elements as the number of rows in solar_pos

radius

Effective search radius (in CRS units) for considering obstacles when calculating shadow and SVF. The default is to use a global search, i.e. radius=Inf. Using a smaller radius can be used to speed up the computation when working on large areas. Note that the search radius is not specific per grid point; instead, a buffer is applied on all grid points combined, then "dissolving" the individual buffers, so that exactly the same obstacles apply to all grid points

returnList

Logical, determines whether to return summed radiation over the entire period per 3D point (default, FALSE), or to return a list with all radiation values per time step (TRUE)

parallel

Number of parallel processes or a predefined socket cluster. With parallel=1 uses ordinary, non-parallel processing. Parallel processing is done with the parallel package

Value

If returnList=FALSE (the default), then returned object is a data.frame, with rows corresponding to grid points and four columns corresponding to the following estimates:

• svf Computed Sky View Factor (see function SVF)

• direct Total direct radiation for each grid point

• diffuse Total diffuse radiation for each grid point

• total Total radiation (direct + diffuse) for each grid point

Each row of the data.frame gives summed radiation values for the entire time period in solar_pos, solar_normal and solar_diffuse

If returnList=TRUE then returned object is a list with two elements:

• direct Total direct radiation for each grid point

• diffuse Total diffuse radiation for each grid point

Each of the elements is a matrix with rows corresponding to grid points and columns corresponding to time steps in solar_pos, solar_normal and solar_diffuse

Details

Input arguments for this function comprise the following:

• An extruded polygon obstacles layer (obstacles and obstacles_height_field) inducing shading on the queried grid

• A grid of 3D points (grid) where radiation is to be estimated. May be created from the 'obstacles' layer, or a subset of it, using function surfaceGrid. For instance, in the code example (see below) radiation is estimated on a grid covering just one of four buildings in the build layer (the first building), but all four buildings are taken into account for evaluating self- and mutual-shading by the buildings.

• Solar positions matrix (solar_pos)

• Direct and diffuse radiation meteorological estimate vectors (solar_normal and solar_diffuse)

Given these inputs, the function goes through the following steps:

• Determining whether each grid point is shaded, at each solar position, using inShadow

• Calculating the coefficient of Direct Normal Irradiance reduction, using coefDirect

• Summing direct radiation considering (1) mutual shading, (2) direst radiation coefficient and (3) direct radiation estimates

• Calculating the Sky View Factor (SVF) for each point, using SVF

• Summing diffuse radiation load considering (1) SVF and (2) diffuse radiation estimates

• Summing total (direct + diffuse) radiation load

Examples

# NOT RUN {
# Create surface grid
grid = surfaceGrid(
  obstacles = build[1, ],
  obstacles_height_field = "BLDG_HT",
  res = 2
)

solar_pos = tmy[, c("sun_az", "sun_elev")]
solar_pos = as.matrix(solar_pos)

# Summed 10-hour radiation estimates for two 3D points
rad1 = radiation(
  grid = grid[1:2, ],
  obstacles = build,
  obstacles_height_field = "BLDG_HT",
  solar_pos = solar_pos[8:17, , drop = FALSE],
  solar_normal = tmy$solar_normal[8:17],
  solar_diffuse = tmy$solar_diffuse[8:17],
  returnList = TRUE
)
rad1

# }
# NOT RUN {
# Same, using 'time' instead of 'solar_pos'

rad2 = radiation(
  grid = grid[1:2, ],
  obstacles = build,
  obstacles_height_field = "BLDG_HT",
  time = as.POSIXct(tmy$time[8:17], tz = "Asia/Jerusalem"),
  solar_normal = tmy$solar_normal[8:17],
  solar_diffuse = tmy$solar_diffuse[8:17],
  returnList = TRUE
)
rad2

# Differences due to the fact that 'tmy' data come with their own
# solar positions, not exactly matching those calulated using 'maptools::solarpos'
rad1$direct - rad2$direct
rad1$diffuse - rad2$diffuse

# }
# NOT RUN {
# }
# NOT RUN {
### Warning! The calculation below takes some time.

# Annual radiation estimates for entire surface of one building
rad = radiation(
  grid = grid,
  obstacles = build,
  obstacles_height_field = "BLDG_HT",
  solar_pos = solar_pos,
  solar_normal = tmy$solar_normal,
  solar_diffuse = tmy$solar_diffuse,
  parallel = 3
)

# 3D plot of the results
library(plot3D)
opar = par(mfrow=c(1, 3))

scatter3D(
  x = coordinates(grid)[, 1],
  y = coordinates(grid)[, 2],
  z = coordinates(grid)[, 3],
  colvar = rad$direct / 1000,
  scale = FALSE,
  theta = 55,
  pch = 20,
  cex = 1.35,
  clab = expression(paste("kWh / ", m^2)),
  main = "Direct radiation"
)
scatter3D(
  x = coordinates(grid)[, 1],
  y = coordinates(grid)[, 2],
  z = coordinates(grid)[, 3],
  colvar = rad$diffuse / 1000,
  scale = FALSE,
  theta = 55,
  pch = 20,
  cex = 1.35,
  clab = expression(paste("kWh / ", m^2)),
  main = "Diffuse radiation"
)
scatter3D(
  x = coordinates(grid)[, 1],
  y = coordinates(grid)[, 2],
  z = coordinates(grid)[, 3],
  colvar = rad$total / 1000,
  scale = FALSE,
  theta = 55,
  pch = 20,
  cex = 1.35,
  clab = expression(paste("kWh / ", m^2)),
  main = "Total radiation"
)

par(opar)
# }