View source: R/radiationtools.R
shortwaveveg | R Documentation |
shortwaveveg
is used to downscale the flux density of shortwave radiation
received at the surface of the Earth, accounting for both topographic and
canopy effects.
shortwaveveg(
dni,
dif,
julian,
localtime,
lat = NA,
long = NA,
dtm = array(0, dim = c(1, 1)),
slope = NA,
aspect = NA,
svv = 1,
alb = 0.23,
fr,
albr = 0.23,
ha = 0,
reso = 1,
merid = NA,
dst = 0,
shadow = TRUE,
x,
l,
difani = TRUE
)
dni |
a single numeric value, SpatRaster object, two-dimensional array or matrix of coarse-resolution direct radiation perpendicular to the solar beam (MJ m-2 hr-1). |
dif |
a single numeric value, SpatRaster object, two-dimensional array or matrix of coarse-resolution diffuse radiation horizontal ot the surface (MJ m-2 hr-1). |
julian |
a single integer representing the Julian as returned by |
localtime |
a single numeric value representing local time (decimal hour, 24 hour clock). |
lat |
an optional single numeric value representing the mean latitude of the location for which downscaled radiation is required (decimal degrees, -ve south of equator). |
long |
an optional single numeric value representing the mean longitude of the location for which downscaled radiation is required (decimal degrees, -ve west of Greenwich meridian). |
dtm |
an optional SpatRaster object, two-dimensional array or matrix of elevations (m), orientated as if derived using |
slope |
an optional single value, SpatRaster object, two-dimensional array or matrix of slopes (º). If an array or matrix, then orientated as if derived using |
aspect |
an optional single value, SpatRaster object, two-dimensional array or matrix of aspects (º). If an array or matrix, then orientated as if derived using |
svv |
an optional SpatRaster object, two-dimensional array or matrix of values representing the proportion of isotropic radiation received by a surface partially obscured by topography relative to the full hemisphere underneath vegetation as returned by |
alb |
an optional single value, SpatRaster object, two-dimensional array or matrix of values representing albedo(s) as returned by |
fr |
a SpatRaster object, two-dimensional array or matrix of fractional canopy cover as returned by |
albr |
an optional single value, SpatRaster object, two-dimensional array or matrix of values representing the albedo(s) of adjacent surfaces as returned by |
ha |
an optional SpatRaster object, two-dimensional array or matrix of values representing the mean slope to the horizon (decimal degrees) of surrounding surfaces from which radiation is reflected for each cell of |
reso |
a single numeric value representing the spatial resolution of |
merid |
an optional numeric value representing the longitude (decimal degrees) of the local time zone meridian (0 for GMT). Default is |
dst |
an optional numeric value representing the time difference from the timezone meridian (hours, e.g. +1 for BST if |
shadow |
an optional logical value indicating whether topographic shading should be considered (False = No, True = Yes). |
x |
a SpatRaster object, two-dimensional array or matrix of numeric values representing the ratio of vertical to horizontal projections of leaf foliage as returned by |
l |
a SpatRaster object, two-dimensional array or matrix of leaf area index values as returned by |
difani |
an optinional logical indicating whether to treat a proportion of the diffuse radiation as anistropic (see details). |
If slope
is unspecified, and dtm
is a SpatRaster slope
and aspect
are calculated from
the SpatRaster. If slope
is unspecified, and dtm
is not a SpatRaster, the slope and aspect
are set to zero. If lat
is unspecified, and dtm
is a SpatRaster with a coordinate reference
system defined, lat
and long
are calculated from the SpatRaster. If lat
is unspecified,
and dtm
is not a SpatRaster, or a SpatRaster without a coordinate reference system defined, an
error is returned. If dtm
is specified, then the projection system used must be such that
units of x, y and z are identical. Use terra::project()
to convert the projection to a
Universal Transverse Mercator type projection system. If dtm
is a SpatRaster object, a SpatRaster
object is returned. If dni
or dif
are SpatRaster
objects, two-dimensional arrays or matrices, then it is assumed that they have been
derived from coarse-resolution data by interpolation, and have the same extent as dtm
.
If no value for ha
is provided, the mean slope to the horizon is assumed
to be 0. If no value for svv
is provided, then the entire hemisphere is
assumed to be in view. If values of albg
and albr
are not specified,
then a default value of 0.23, typical of well-watered grass is assumed. If
single values of albg
and albr
are given, then the entire area is
assumed to have the same albedo. If dtm
is specified, then the projection
system used must be such that the units of x, y and z are identical. Use
terra::project()
to convert the projection to a Universal Transverse
Mercator type projection system. If no value for dtm
is provided,
radiation is downscaled by deriving values on the inclined surfaces
specified in slope
and aspect
and topographic shadowing is ignored. If
single values are provided for slope
and aspect
, the entire extent
covered by fr
is assumed to have the same slope and aspect. Only single
values of lat
and long
are taken as inputs. Under partially
cloudy conditions, a proportion of diffuse radiation is typically anistropic
(directional). If difani
is TRUE (the default), then the assumption is made that
hourly direct radiation transmission can define the portions of the diffuse
radiation to be treated as anistropic and isotropic. If dtm
covers a large
extent, the dtm
is best divided into blocks and seperate calculations
performed on each block. Since horizon angles, topographic shading and
sky view correction factors may be influenced by locations beyond the extent of dtm
, it is best to ensure
dtm
covers a larger extent than that for which radiation values are
needed, and to ensure sub-divided blocks overlap in extent. Calculations are faster
if values for all inputs are provided.
a SpatRaster object, two-dimensional array of numeric values representing net shortwave radiation (MJ per metre squared per hour).
The function terra::terrain()
can be used to derive slopes and aspects from dtm
(see example).
Function shortwavetopo()
returns net shortwave radiation, or components thereof, above the canopy.
library(terra)
# =================================
# Extract data for 2010-05-24 11:00
# =================================
dni <- microvars$dni[564]
dif <- microvars$dif[564]
# ==========================
# Calculate input paramaters
# ==========================
x <- leaf_geometry(rast(veg_hgt))
l <- lai(aerial_image[,,3], aerial_image[,,4])
l <- lai_adjust(l, rast(veg_hgt))
fr <- canopy(l)
alb <- albedo(aerial_image[,,1], aerial_image[,,2], aerial_image[,,3],
aerial_image[,,4])
sv <- skyviewveg(rast(dtm1m), l, x)
jd <- julday(2010, 5, 24)
ha <- mean_slope(rast(dtm1m))
# ===============================================================
# Calculate and plot net shortwave radiation for 2010-05-24 11:00
# ===============================================================
netshort1m <- shortwaveveg(dni, dif, jd, 11, dtm = rast(dtm1m), svv = sv, alb = alb,
fr = fr, ha = ha, x = x, l = l)
plot(mask(netshort1m, rast(dtm1m)), main = "Net shortwave radiation")
Add the following code to your website.
For more information on customizing the embed code, read Embedding Snippets.