Package 'viewscape'

Title: Viewscape Analysis
Description: Provides tools for viewscape analysis from one or multiple viewpoints using a digital surface or elevation model. Core functionality includes computing viewsheds, quantifying visual magnitude, calculating a suite of viewscape configuration metrics (extent, depth, relief, Sky View Factor, skyline variation, and patch-based landscape structure), Shannon Diversity Index and land cover feature proportions within the visible area, pairwise intervisibility networks, panoramic view generation, and visualizing results as rasters or polygons. Viewscape configuration metrics follow the methods of Tabrizian et al. (2020) <doi:10.1016/j.landurbplan.2019.103704>. The viewshed algorithm is based on Franklin & Ray (1994) <https://api.semanticscholar.org/CorpusID:10680920> and Wang et al. (2000) <https://api.semanticscholar.org/CorpusID:131687018>. Visual magnitude is derived from Chamberlain & Meitner (2013) <doi:10.1016/j.landurbplan.2013.01.003>. Sky View Factor is computed following Oke (1981) <doi:10.1002/joc.3370010304> as implemented in the 'shadow' package (Dorman et al. 2019) <doi:10.32614/RJ-2019-024>.
Authors: Xiaohao Yang [aut, cre, cph], Nathan Fox [aut], Derek Van Berkel [aut], Mark Lindquist [aut]
Maintainer: Xiaohao Yang <[email protected]>
License: GPL-3
Version: 3.0.1
Built: 2026-10-18 08:59:54 UTC
Source: https://github.com/billbillbilly/viewscape

Help Index


calculate_diversity

Description

The calculate_diversity function is designed to calculate landscape diversity metrics within a viewshed. It takes as input a land cover raster, a viewshed object representing the observer's line of sight, and an optional parameter to compute class proportions.

Usage

calculate_diversity(viewshed, land, proportion = FALSE)

Arguments

viewshed

Viewshed object.

land

Raster. The raster of land use/land cover representing different land use/cover classes.

proportion

logical (Optional), indicating whether to return class proportions along with the Shannon Diversity Index (SDI). (default is FALSE).

Value

List. a list containing the Shannon Diversity Index (SDI) and, if the proportion parameter is set to TRUE, a table of class proportions within the viewshed.

Examples

library(viewscape)
# Load a viewpoint
test_viewpoint <- sf::read_sf(system.file("test_viewpoint.shp", package = "viewscape"))
# load dsm raster
dsm <- terra::rast(system.file("test_dsm.tif", package ="viewscape"))
#Compute viewshed
output <- viewscape::compute_viewshed(dsm = dsm,
                                      viewpoints = test_viewpoint,
                                      offset_viewpoint = 6, r = 1600)
# load landuse raster
test_landuse <- terra::rast(system.file("test_landuse.tif",
                                        package ="viewscape"))
diversity <- viewscape::calculate_diversity(output,
                                            test_landuse)

calculate_feature

Description

The calculate_feature function is designed to extract specific feature-related information within a viewshed. It allows you to compute the proportion of the feature that is present in the viewshed.

Usage

calculate_feature(viewshed, feature, type, exclude_value)

Arguments

viewshed

Viewshed object.

feature

Raster. Land cover or land use

type

Numeric. The input type of land cover raster. type=1: percentage raster (that represents the percentage of area in each cell). type=2: binary raster (that only uses two values to represent whether the feature exists in each cell).

exclude_value

Numeric. the value of those cells need to be excluded in the analysis. If type = 2, exclude_value is reqired.

Value

Numeric. The canopy area in the viewshed.

Examples

library(viewscape)
# Load a viewpoint
test_viewpoint <- sf::read_sf(system.file("test_viewpoint.shp", package = "viewscape"))
# load dsm raster
dsm <- terra::rast(system.file("test_dsm.tif", package ="viewscape"))
#Compute viewshed
viewshed <- viewscape::compute_viewshed(dsm = dsm,
                                        viewpoints = test_viewpoint,
                                        offset_viewpoint = 6)
# load canopy raster
test_canopy <- terra::rast(system.file("test_canopy.tif",
                                       package ="viewscape"))
# calculate the percentage of canopy coverage
test_canopy_proportion <- viewscape::calculate_feature(viewshed = viewshed,
                                                       feature = test_canopy,
                                                       type = 2,
                                                       exclude_value = 0)

calculate_viewmetrics

Description

The calculate_viewmetrics function is designed to compute a set of configuration metrics based on a given viewshed object and optionally, digital surface models (DSM) and digital terrain models (DTM) for terrain analysis. The function calculates various metrics that describe the visibility characteristics of a landscape from a specific viewpoint.The metrics include:

  1. Extent: The total area of the viewshed, calculated as the number of visible grid cells multiplied by the grid resolution

  2. Depth: The furthest visible distance within the viewshed from the viewpoint

  3. Vdepth: The standard deviation of distances to visible points, providing a measure of the variation in visible distances

  4. Horizontal: The total visible horizontal or terrestrial area within the viewshed

  5. Relief: The standard deviation of elevations of the visible ground surface

  6. Sky View Factor (SVF): The proportion of visible sky hemisphere at the viewpoint, computed by casting rays in 36 azimuth directions and calculating SVF = mean(cos^2(max obstruction angle)). Ranges from 0 (fully enclosed) to 1 (open sky).

  7. Skyline: The standard deviation of the vertical viewscape, including visible canopy and buildings, when specified

  8. Number of patches: Visible fragmentation measured by total visible patches with the viewscape

  9. Mean shape index: Visible patchiness based on average perimeter-to-area ratio for all viewscape patches (vegetation and building)

  10. Edge density: A measure of visible complexity based on the length of patch edges per unit area

  11. Patch size: Total average size of a patches over the entire viewscape area

  12. Patch density: Visible landscape granularity based on measuring patch density

  13. Shannon diversity index: The abundance and evenness of land cover/use in a viewshed

  14. Proportion of object: Proportion of a single type of land use or cover in a viewshed

Usage

calculate_viewmetrics(viewshed, dsm, dtm, masks = list())

Arguments

viewshed

Viewshed object.

dsm

Raster, Digital Surface Model for the calculation of

dtm

Raster, Digital Terrain Model

masks

List, a list including rasters of canopy and building footprints. For example of canopy raster, the value for cells without canopy should be 0 and the value for cells with canopy can be any number.

Value

List

References

Tabrizian, P., Baran, P.K., Berkel, D.B., Mitásová, H., & Meentemeyer, R.K. (2020). Modeling restorative potential of urban environments by coupling viewscape analysis of lidar data with experiments in immersive virtual environments. Landscape and Urban Planning, 195, 103704.

Dorman, M., Vulkan, A., Erell, E., & Kloog, I. (2019). shadow: Geometric Shadow Calculations. The R Journal, 11(1), 287–309. https://github.com/michaeldorman/shadow

Oke, T. R. (1981). Canyon geometry and the nocturnal urban heat island: comparison of scale model and field observations. Journal of Climatology, 1(3), 237–254.

Examples

# Load in DSM
test_dsm <- terra::rast(system.file("test_dsm.tif",
                                    package ="viewscape"))
# Load DTM
test_dtm <- terra::rast(system.file("test_dtm.tif",
                                    package ="viewscape"))

# Load canopy raster
test_canopy <- terra::rast(system.file("test_canopy.tif",
                                       package ="viewscape"))

# Load building footprints raster
test_building <- terra::rast(system.file("test_building.tif",
                                         package ="viewscape"))

# Load in the viewpoint
test_viewpoint <- sf::read_sf(system.file("test_viewpoint.shp",
                                          package = "viewscape"))

# Compute viewshed
output <- viewscape::compute_viewshed(dsm = test_dsm,
                                      viewpoints = test_viewpoint,
                                      offset_viewpoint = 6, r = 1600)

# calculate metrics given the viewshed, canopy, and building footprints
test_metrics <- viewscape::calculate_viewmetrics(output,
                                                 test_dsm,
                                                 test_dtm,
                                                 list(test_canopy, test_building))

compute_viewshed

Description

The compute_viewshed function is designed for computing viewsheds, which are areas visible from specific viewpoints, based on a Digital Surface Model (DSM). It provides flexibility for single or multi-viewpoint analyses and allows options for parallel processing, raster output, and plotting. Optionally, a binary canopy layer can be supplied so that tree crowns stop acting as solid ground-to-crown walls, which is what makes a viewpoint beneath a canopy usable (see details).

Usage

compute_viewshed(
  dsm,
  viewpoints,
  offset_viewpoint = 1.7,
  offset_height = 0,
  r = NULL,
  refraction_factor = 0.13,
  method = "view_tree",
  parallel = FALSE,
  workers = 1,
  raster = FALSE,
  plot = FALSE,
  canopy = NULL,
  dem = NULL,
  cbh = NULL
)

Arguments

dsm

Raster, the digital surface model/digital elevation model

viewpoints

sf point(s) or vector including x,y coordinates of a viewpoint or a matrix including several viewpoints with x,y coordinates

offset_viewpoint

numeric, observer height above the surface in metres (default is 1.7 m). Automatically converted to the CRS unit of the DSM.

offset_height

numeric, target height above the surface in metres (default is 0 m). Automatically converted to the CRS unit of the DSM.

r

Numeric (optional), analysis radius in metres (default is 1000 m). Automatically converted to the CRS unit of the DSM.

refraction_factor

Number, indicating the refraction factor. The refraction factor adjusts the effect of atmospheric refraction on the apparent curvature of the Earth. In most standard applications, a refraction factor of 0.13 is used, and so does this function. However, the appropriate refraction factor may vary depending on environmental conditions.

method

Character, The algorithm for computing a viewshed: "view_tree", "plane", and "los" (see details). "view_tree" is used as default.

parallel

Logical, (default is FALSE) indicating if parallel computing should be used to compute viewsheds of multiview points. When it is TRUE, arguements 'raster' and 'plot' are ignored

workers

Numeric, indicating the number of CPU cores that will be used for parallel computing. It is required if 'parallel' is 'TRUE'.

raster

Logical, (default is FALSE) if it is TRUE, the raster of viewshed will be returned. The default is FALSE

plot

Logical, (default is FALSE) if it is TRUE, the raster of viewshed will be displayed

canopy

Raster (optional), a binary tree-canopy layer aligned with 'dsm', in which cells greater than 0 are canopy and cells that are 0 or NA are not. When it is supplied, canopy cells stop acting as ground-to-crown obstructions and 'dem' becomes required (see details).

dem

Raster (optional), the bare-earth digital elevation model, aligned with 'dsm'. Required when 'canopy' is supplied. It is used both to replace the canopy-top elevation in the obstruction surface and to place an observer who stands underneath a crown.

cbh

Numeric or Raster (optional), crown base height: the height of the lowest foliage above bare earth. Canopy cells whose crown base sits at or below the observer's eye level are kept opaque, because foliage that reaches down to eye level does obstruct the view; the rest are cleared. A single number is interpreted in metres and converted to the CRS unit of the DSM; a raster is taken to be in the vertical unit of 'dsm' already, and its NA cells are treated as 0. When it is NULL (the default) every canopy cell is cleared. Only used when 'canopy' is supplied.

Details

For method, "view_tree" is the view-tree algorithm (Wang et al. 2022). It converts the DEM into a tree data structure based on spatial occlusion relationships among cells, then traverses the tree depth-first, propagating a minimum visible slope (MVS) from parent to child nodes to determine visibility. The method achieves precision comparable to "los" while being approximately 40% faster than R2-class methods.

"plane" is the reference plane algorithm introduced by Wang et al. (2000) and "los" is the line of sight algorithm (Franklin & Ray, 1994). The reference plane algorithm can be more time-efficient than the line of sight algorithm, whereas the accuracy of the line of sight is better.

A DSM stores one elevation per cell, so a tree crown is represented by its top surface only and every viewshed algorithm necessarily treats that crown as opaque from the ground up. This misrepresents an observer standing beneath a canopy, in a city street as much as in a forest, where the line of sight passes under the crown. Supplying 'canopy' together with 'dem' replaces the elevation of canopy cells with the bare-earth elevation, so crowns stop blocking sight lines while non-canopy cells, including buildings and other built structures, keep their DSM elevation and stay opaque.

The choice is deliberately all-or-nothing per cell. In a 2.5D surface a cell obstructs everything below its elevation, so the open space beneath a crown cannot be encoded as an elevation at all: writing dem + cbh into a canopy cell would wall off exactly the gap the observer looks through, and would block more than the canopy top does for the descending sight lines that reach distant ground. Representing a crown as open below the foliage and solid above requires testing the height of each sight line against a height band as it passes through the cell, which a single elevation matrix cannot express.

'cbh' is therefore used as a per-cell opacity classifier rather than as an elevation. Canopy cells whose crown base sits at or below the observer's eye level keep their DSM elevation, since foliage that hangs down to eye level genuinely obstructs the view: hedges, shrubs and low crowns stay solid. Canopy cells whose lowest foliage is above eye level are cleared. This uses the observer's eye height as a stand-in for the height of the sight line, which is exact on flat ground for targets at eye level and approximate elsewhere.

The observer is placed from a separate surface. 'offset_viewpoint' is measured above bare earth when the viewpoint falls on a canopy cell and above the DSM everywhere else, so an observer under a tree stands on the ground rather than on top of the crown, while an observer on a roof or a bridge deck is unaffected. This is what makes a viewpoint beneath a canopy usable at all: the elevation the observer is offset from no longer depends on whether a crown happens to overhang.

Two limitations are worth stating. Trunks are not represented, so sight lines pass through stems; where stem positions are known they can be burned into 'canopy' as non-canopy cells to restore them. And a cleared crown is fully transparent rather than attenuating, so visibility through a deep stand of high-crowned trees is overestimated. Where that matters, computing the viewshed twice, once with 'canopy' and once without, brackets the true value.

Value

Raster or list. For single-viewpoint analysis, the function returns either a raster (raster is TRUE) or a viewshed object. Value 1 means visible while value 0 means invisible. For multi-viewpoint analysis, a list of viewsheds is returned.

References

Franklin, W. R., & Ray, C. (1994, May). Higher isn't necessarily better: Visibility algorithms and experiments. In Advances in GIS research: sixth international symposium on spatial data handling (Vol. 2, pp. 751-770). Edinburgh: Taylor & Francis.

Wang, J., Robinson, G. J., & White, K. (2000). Generating viewsheds without using sightlines. Photogrammetric engineering and remote sensing, 66(1), 87-90.

Wang, Z., Xiong, L., Guo, Z., Zhang, W., & Tang, G. (2022). A view-tree method to compute viewsheds from digital elevation models. International Journal of Geographical Information Science. https://doi.org/10.1080/13658816.2022.2094385

See Also

fov_mask() visual_magnitude()

Examples

# Load a viewpoint
test_viewpoint <- sf::read_sf(system.file("test_viewpoint.shp", package = "viewscape"))
# load dsm raster
dsm <- terra::rast(system.file("test_dsm.tif", package ="viewscape"))
#Compute viewshed
output <- viewscape::compute_viewshed(dsm = dsm,
                                      viewpoints = test_viewpoint,
                                      offset_viewpoint = 6, r = 1600)

compute_viewshed_understory

Description

Experimental. Computes a viewshed for an observer standing beneath a tree canopy, where the sight line passes under the crowns. Unlike a conventional 2.5D viewshed, in which a canopy cell obstructs everything below the canopy top, each crown is modelled as a prism between its base and its top and blocks only sight lines that actually pass through that height band. Lines that stay below the crown base pass underneath and lines that have risen above the crown top pass over.

Usage

compute_viewshed_understory(
  dsm,
  dem,
  canopy,
  viewpoints,
  cbh = 3,
  offset_viewpoint = 1.7,
  offset_height = 0,
  r = NULL,
  canopy_removal_r = 0,
  clear_mask = NULL,
  prefilter = FALSE,
  method = "view_tree"
)

Arguments

dsm

Raster, the digital surface model. Supplies the crown tops and the elevation of every non-canopy obstruction.

dem

Raster, the bare-earth digital elevation model, aligned with 'dsm'.

canopy

Raster, a binary tree-canopy layer aligned with 'dsm', in which cells greater than 0 are canopy and cells that are 0 or NA are not.

viewpoints

sf point, or a vector of x,y coordinates. One viewpoint per call; loop over points for several.

cbh

Numeric or Raster (default is 3), crown base height above bare earth: the height of the lowest foliage. A single number is interpreted in metres and converted to the CRS unit of the DSM; a raster is taken to be in the vertical unit of 'dsm' already. Where the crown base would reach the canopy top, the whole column is foliage and the cell becomes opaque from the ground up, so hedges and shrubs need no special handling. That rule reads the canopy height model, so it also fires wherever 'cbh' simply exceeds the crowns present: past that point raising 'cbh' makes cells opaque rather than opening them, and visibility falls. A warning is issued when it affects more than a tenth of the canopy cells. Check 'cbh' against dsm - dem before trusting a value.

offset_viewpoint

Numeric, observer height in metres above the canopy-cleared surface (default is 1.7 m): bare earth when the viewpoint falls under a crown, the DSM elsewhere. An observer beneath a tree therefore stands on the ground rather than on the canopy top, while one on a roof or a bridge deck keeps that elevation.

offset_height

Numeric, target height above the surface in metres (default is 0 m).

r

Numeric, analysis radius in metres (default is 1000 m).

canopy_removal_r

Numeric (default is 0), radius in metres around the viewpoint within which crowns are cleared to bare earth and contribute no obstruction at all. The default leaves only the observer's own cell exempt, which is the current behaviour; a positive value is a near-field correction (see details). Cleared cells are still reported as canopy in the class raster, since the crowns are physically there.

clear_mask

Raster, sf or SpatVector (optional), open space where crowns are known not to obstruct: a road corridor, a path, a clearing. Crown prisms are suppressed on these cells exactly as within 'canopy_removal_r', and the two combine as a union. A raster is treated as open where cells are greater than 0 and has to line up with 'dsm'; a geometry is transformed and rasterized onto 'dsm' internally. See road_corridor_mask() for building one from road centrelines, and the details below for why this is preferable to a radius.

prefilter

Logical (default is FALSE), whether to run compute_viewshed() on the canopy-cleared surface first and skip the surface targets it rejects. This is faster on large windows, but it substitutes that algorithm's occlusion semantics for this one's, so results can differ from the unfiltered run; the filtered result is always a subset of it (see details).

method

Character, the algorithm used for the pass-1 filter: "view_tree", "plane" or "los". Ignored when 'prefilter' is FALSE.

Details

The test is carried out in slope space. For a target at slope s=(zT−zO)/dTs = (z_T - z_O)/d_T the sight line reaches height zO+s dz_O + s\,d at distance dd, and a canopy cell at distance dcd_c with crown base bb and top tt blocks the target when

b≤zO+s dc≤t.b \le z_O + s\,d_c \le t.

Written per canopy cell as an interval of blocked slopes, [(b−zO)/dc, (t−zO)/dc][(b - z_O)/d_c,\ (t - z_O)/d_c], this is what an azimuth sweep accumulates, so the same test carries over to a sweep implementation.

Two consequences are worth noting. The critical pitch is not a constant: the slope at which a given crown starts to block scales as 1/dc1/d_c, so the same tree blocks steeply near the observer and shallowly far away. And visibility is not monotone in pitch, because a sufficiently steep line clears the crown top; a single pitch threshold cannot express either behaviour.

Non-canopy cells keep the usual rule and stay opaque below their DSM elevation, so buildings, walls and terrain are unaffected. Bare earth beneath a crown still obstructs, which is what makes the method behave correctly on slopes.

At the extremes the crown model switches itself off. A crown base of 0 makes the prism span the whole column, and a crown base above the canopy top collapses to the same thing by the rule described under 'cbh'. Either way the result is bit-for-bit identical to running this function on a surface with no canopy at all, which is the check that isolates the band test from everything else. It is close to but not identical with compute_viewshed() on the raw DSM, for the reason given below.

Because canopy cells are targets in their own right, an observer under a crown sees the crowns around and above them. That is reported separately in the class raster rather than mixed into surface visibility, since the ground under a tree and the tree itself are different observations about the same cell.

'canopy_removal_r' addresses where the prism model is weakest. The blocked band a crown subtends grows as the crown gets closer: with an eye at 1.7 m, a crown base at 3 m and a top at 14 m, a canopy cell 1 m away blocks pitches from 52 to 85 degrees, one 2 m away blocks 33 to 80, and one 5 m away blocks 15 to 68. A handful of neighbouring cells therefore wall off almost everything above the horizontal, even though those cells are modelled as a continuous solid slab where a real crown is a third to a half gaps, and even though an observer beneath a canopy is normally standing in an opening: a road, a path, the space between trunks. Clearing crowns within a small radius is a blunt but explicit correction for that, and for a viewpoint snapped to a road it has a direct reading, namely that the road corridor is open even where overhanging crowns are mapped across it.

In practice the radius is almost entirely a canopy-visibility parameter rather than a ground-visibility one. On a 250 m window under 32% canopy cover, raising it from 0 to 40 m moved visible surface area by less than a tenth of a percent while multiplying visible canopy roughly fiftyfold. The reason is that surface targets are reached by descending or near-horizontal sight lines that were already passing below the crown base, so the near-field crowns were never blocking them; those crowns only intercept the rising lines, which are the ones aimed at other crowns and at the sky. Structural visibility is therefore fairly insensitive to the choice, and canopy exposure is very sensitive to it, which is worth stating explicitly whenever the two are reported together.

It remains a free parameter, so it wants reporting and a sensitivity analysis rather than tuning. Where the open space is actually known, 'clear_mask' is the better instrument, because a road corridor is not a disc. A radius clears equally in every direction, including straight into the crowns flanking the road where nothing should be cleared, whereas a corridor clears along the street and leaves the flanking canopy intact. It should also be applied along the whole length of the corridor within 'r' rather than only near the viewpoint: a tree-lined street is a tunnel, the crowns arch over it from both sides, and the space below crown base is open for the entire run of the street, so the crowns that matter are the ones overhanging the road a hundred metres out rather than the ones directly overhead.

Note what the canopy classes do and do not count. A canopy cell is aimed at its crown top, so under a closed canopy the sight line has to climb through the neighbouring foliage to get there and almost nothing beyond the adjacent ring of crowns qualifies. That is the enclosure effect behaving correctly, but it means classes 2 and 3 answer "is that crown's top visible", not "is there foliage in that direction". For the latter, what matters is the first cell at which a sight line enters a crown, which is the quantity a hemispherical (azimuth by pitch) formulation reports and which this function does not yet compute.

The sight line is modelled as a corridor half a cell wide rather than an infinitely thin ray: at each step along its major axis the march covers the whole range the corridor spans on the minor axis. Sampling a single cell per step would let a one-cell obstruction slip through on a shallow diagonal, and the resulting shadow would open and close along the same bearing depending on how the target's distance happened to round.

This is point-to-point visibility, evaluated independently for every target, which is not what compute_viewshed() computes. The view-tree and reference plane algorithms propagate a horizon between cells, so an isolated obstruction casts a shadow that persists and widens behind it, whereas exact geometry lets a sight line reopen past a narrow obstacle it genuinely passes beside. On the same surface with no canopy this function is the more conservative of the two at grazing incidence: it will not report a cell visible that compute_viewshed() rejects, but it does reject some that compute_viewshed() accepts, mostly along the far side of flat roofs and other runs of equal elevation. That is a difference in occlusion semantics rather than in the canopy model, and it is why 'prefilter' defaults to FALSE.

Known limitations: a crown is a vertical prism over its cell footprint, so rounded crown margins are wrong at roughly cell scale, and an isolated single-cell crown casts a less persistent shadow than a contiguous patch; the prism is fully opaque rather than attenuating, so clipping 1 m and 30 m of foliage are treated alike; trunks are absent, since the column below the crown base is empty; and curvature and refraction corrections are not applied.

Value

A list with three elements:

viewshed

a Viewshed object whose 'visible' matrix is 1 where the ground or built surface is visible, for compatibility with calculate_viewmetrics() and visualize_viewshed().

class

a SpatRaster of visibility classes: 0 invisible, 1 visible surface, 2 visible canopy only, 3 both the surface and the crown of a canopy cell.

candidate

the pass-1 viewshed used as the filter, or NULL when 'prefilter' is FALSE. Useful as the upper bound to compare against.

See Also

compute_viewshed() calculate_viewmetrics()

Examples

# a flat scene with a band of high-crowned trees to the east
dem <- terra::rast(nrows = 61, ncols = 61, xmin = 0, xmax = 61,
                   ymin = 0, ymax = 61, crs = "EPSG:32617", vals = 0)
cnp <- terra::setValues(dem, 0)
m <- terra::as.matrix(cnp, wide = TRUE)
m[, 35:39] <- 1
cnp <- terra::setValues(cnp, m)
dsm <- dem + cnp * 10

out <- compute_viewshed_understory(dsm = dsm, dem = dem, canopy = cnp,
                                   viewpoints = c(30.5, 30.5),
                                   cbh = 3, r = 60)
terra::plot(out$class)

fov_mask

Description

The fov_mask function is designed to subset a viewshed based on its viewpoint and the field of view

Usage

fov_mask(viewshed, fov)

Arguments

viewshed

viewshed object, generated by compute_viewshed()

fov

Vector, specifying the field of view with two angles in degrees (e.g. c(-45, 45)) for masking a viewshed based on its viewpoint. See details.

Details

Angles are measured clockwise from north, ranging from -180 (exclusive) to 180 (inclusive). The two angles may be given in either order. When fov[1] < fov[2] the sector spans clockwise from fov[1] to fov[2]. When fov[1] > fov[2] the sector wraps around ±180 (south) — for example c(135, -135) captures a 90° south-facing arc.

-45 0 45
-90 v 90
-135 180 135

Here, 'v' represents the viewpoint. North is 0°, east is 90°, west is -90°, south is ±180°.

Value

viewshed object

See Also

compute_viewshed()

Examples

# Load a viewpoint
test_viewpoint <- sf::read_sf(system.file("test_viewpoint.shp", package = "viewscape"))
# load dsm raster
dsm <- terra::rast(system.file("test_dsm.tif", package ="viewscape"))
# Compute viewshed
viewshed <- viewscape::compute_viewshed(dsm,
                                        viewpoints = test_viewpoint,
                                        offset_viewpoint = 6)
# subset viewshed using the field of view (north-facing 120° arc)
out <- viewscape::fov_mask(viewshed, c(-60, 60))

intervis_network

Description

Computes an intervisibility network among a set of viewpoints, determining which pairs of viewpoints have a direct line of sight to each other. Results can be returned as a binary adjacency matrix or as an sf collection of lines suitable for mapping the network.

Usage

intervis_network(
  viewpoints,
  dsm,
  offset_viewpoint = 1.7,
  r = NULL,
  refraction_factor = 0.13,
  method = "view_tree",
  output = "matrix"
)

Arguments

viewpoints

sf point(s), or a matrix/data.frame of x,y coordinates. At least 2 viewpoints are required.

dsm

SpatRaster. A digital surface / elevation model.

offset_viewpoint

Numeric. Observer height above the surface in metres. Either a single value applied to all viewpoints, or a vector of per-viewpoint heights. Default is 1.7 m. Automatically converted to the CRS unit of the DSM.

r

Numeric (optional). Analysis radius in metres. Pairs separated by more than r are treated as non-visible. Default is 1000 m. Automatically converted to the CRS unit of the DSM.

refraction_factor

Numeric. Atmospheric refraction coefficient. Default is 0.13.

method

Character. Viewshed algorithm: "view_tree" (default), "plane", or "los". See compute_viewshed for details.

output

Character. Return type: "matrix" (default) returns an N x N integer adjacency matrix where row i, column j is 1 if viewpoint i can see viewpoint j, and 0 otherwise (NA on the diagonal); "lines" returns an sf data frame of LINESTRINGs connecting visible pairs, with columns from, to, and mutual (TRUE when both directions are visible).

Value

An N x N integer matrix (output = "matrix") or an sf data frame of lines (output = "lines").

See Also

compute_viewshed() fov_mask()

Examples

library(viewscape)
# Load viewpoints (requires at least 2)
test_viewpoints <- sf::read_sf(system.file("test_viewpoints.shp", package = "viewscape"))
# Load DSM
dsm <- terra::rast(system.file("test_dsm.tif", package = "viewscape"))
# Compute intervisibility matrix
mat <- intervis_network(test_viewpoints, dsm, offset_viewpoint = 1.7)
# Compute intervisibility as lines
net <- intervis_network(test_viewpoints, dsm, offset_viewpoint = 1.7,
                        output = "lines")

pano_view

Description

Generate an equirectangular or cylindrical panoramic view from a DSM.

For method = "equirectangular":

  • No semantic: automatically downloads ESRI WorldImagery tiles via greenSD and projects the satellite RGB colours onto the panoramic rays, producing a 3-band colour panorama (R, G, B, 0–255).

  • With semantic: projects the supplied land-use / land-cover raster onto the panoramic rays, returning a 2-layer SpatRaster (depth + land-cover codes).

For method = "cylindrical" the function always returns a single-layer distance panorama (with an optional second semantic layer).

Usage

pano_view(
  dsm = NULL,
  vpt = NULL,
  h = 6,
  semantic = NULL,
  method = "equirectangular",
  sky_threshold = 3,
  step_size = 0.5,
  max_dist = 500,
  pano_dim = c(256, 512),
  ares = 0.0174533,
  heading = 0,
  color_source = NULL,
  sky_color = c(135, 206, 235),
  satellite_zoom = 17L,
  plot = FALSE,
  legend = TRUE,
  axes = TRUE
)

Arguments

dsm

A SpatRaster object (single layer) of surface elevation.

vpt

A POINT sf object or numeric coordinate pair representing the viewpoint (in projected coordinates).

h

Numeric. Height of the viewpoint above the ground surface (in metres). Default: 6.

semantic

A SpatRaster object (single layer) of land-use / land-cover class values. When supplied with method = "equirectangular", the panorama encodes the land-cover code at each ray-hit cell (depth layer plus semantic layer). Binary 0/1 rasters are treated as occupied vertical structures. When supplied with method = "cylindrical", the return value has two layers: depth and semantic class.

method

Character. "equirectangular" (default) or "cylindrical".

sky_threshold

Numeric. Elevation buffer (metres) below the ray position that is still treated as terrain (default: 3.0). Must be less than h; if not, it is automatically clamped to h - 0.01.

step_size

Numeric. Ray marching step size in DSM coordinate units (default: 0.5).

max_dist

Numeric. Maximum ray distance (default: 500).

pano_dim

Integer vector of length 2: panorama height and width in pixels (default: c(256, 512)).

ares

Numeric. Cylindrical vertical angular resolution in radians per pixel. Ignored when method = "equirectangular" (default: 0.0174533).

heading

Numeric. Camera heading in degrees (default: 0 = north). See Details.

color_source

Optional SpatRaster with exactly 3 bands (R, G, B, values 0–255). Overrides the automatic greenSD satellite download when method = "equirectangular" and semantic = NULL. Useful when you already have a local satellite image.

sky_color

Numeric vector of length 3 (R, G, B, 0–255). Colour painted on sky pixels in the colour panorama. Default: sky blue c(135, 206, 235).

satellite_zoom

Integer. Tile zoom level passed to greenSD::get_tile_green() when auto-fetching satellite imagery (default: 17, approximately 1.2 m/px at 42 degrees N).

plot

Logical. Whether to plot the result (default: FALSE).

legend

Logical. Whether to display a legend when plotting (default: TRUE).

axes

Logical. Whether to display axes when plotting (default: TRUE).

Details

For heading:

  • heading = 0 → facing north (default)

  • heading = 90 → facing east

  • heading = 180 → facing south

  • heading = 270 → facing west

The equirectangular colour panorama requires the greenSD package (install with devtools::install_github("billbillbilly/greenSD")). It calls greenSD::get_tile_green(bbox, zoom, provider = "esri") and uses the raw RGB tile (⁠[[2]]⁠) as the colour source.

Value

  • method = "equirectangular", no semantic: 3-band SpatRaster (R, G, B, 0–255). Plot with terra::plotRGB().

  • method = "equirectangular", with semantic: 2-layer SpatRaster (depth in metres, land-cover code).

  • method = "cylindrical", no semantic: single-layer distance SpatRaster (NA = sky).

  • method = "cylindrical", with semantic: 2-layer SpatRaster (depth, semantic class).

References

Lu, X., Li, Z., Cui, Z., Oswald, M. R., Pollefeys, M., & Qin, R. (2020). Geometry-aware satellite-to-ground image synthesis for urban areas. In Proceedings of the IEEE/CVF Conference on Computer Vision and Pattern Recognition (pp. 859-867).

Examples

dsm <- terra::rast(system.file("test_dsm.tif", package = "viewscape"))
vpt <- sf::read_sf(system.file("test_viewpoint.shp", package = "viewscape"))
# cylindrical distance panorama (no internet needed)
result <- viewscape::pano_view(dsm, vpt, h = 6, method = "cylindrical")

pano_view_understory

Description

Experimental. Renders a hemispherical panorama from a viewpoint beneath a tree canopy, and summarises it as solid-angle fractions. Each crown is modelled as a prism between its base and its top, so a ray that stays below the foliage passes underneath and one that enters the foliage stops there. This is the first-hit counterpart of compute_viewshed_understory(): rather than asking whether a cell is visible, it asks what the observer sees in each direction, which is what canopy exposure, sky view and any comparison against a segmented street-level image need.

Usage

pano_view_understory(
  dsm,
  dem,
  canopy,
  vpt,
  h = 1.7,
  cbh = 3,
  canopy_removal_r = 0,
  clear_mask = NULL,
  built = NULL,
  built_threshold = 2,
  foliage_k = NULL,
  t_min = 0.001,
  method = "equirectangular",
  sky_threshold = 0,
  step_size = 0.5,
  max_dist = 500,
  pano_dim = c(256, 512),
  ares = 0.0174533,
  heading = 0,
  plot = FALSE
)

Arguments

dsm

Raster, the digital surface model. Supplies the crown tops and the elevation of every non-canopy obstruction.

dem

Raster, the bare-earth digital elevation model, aligned with 'dsm'.

canopy

Raster, a binary tree-canopy layer aligned with 'dsm', in which cells greater than 0 are canopy and cells that are 0 or NA are not.

vpt

A POINT sf object, or a numeric x,y coordinate pair.

h

Numeric, observer height in metres above the canopy-cleared surface (default is 1.7 m): bare earth when the viewpoint falls under a crown, the DSM elsewhere. Note this differs from pano_view(), whose default of 6 m suits a camera mast rather than an eye.

cbh

Numeric or Raster (default is 3), crown base height above bare earth. Where the crown base would reach the canopy top the whole column is foliage and the cell becomes opaque from the ground up, so hedges and shrubs need no special handling. That rule reads the canopy height model, so it also fires wherever 'cbh' simply exceeds the crowns present: past that point raising 'cbh' makes cells opaque rather than opening them, and sky view falls. A warning is issued when it affects more than a tenth of the canopy cells. Check 'cbh' against dsm - dem before trusting a value.

canopy_removal_r

Numeric (default is 0), radius in metres within which crowns are treated as transparent. See compute_viewshed_understory().

clear_mask

Raster, sf or SpatVector (optional), known-open ground where crowns do not obstruct, such as a road corridor from road_corridor_mask().

built

Raster (optional), a binary layer marking built structures, used only to separate the built class from the ground class. When it is NULL, non-canopy cells standing more than 'built_threshold' above bare earth are taken to be built. Canopy cells are never counted as built either way, since a canopy cell's surface is bare earth by construction.

built_threshold

Numeric (default is 2), height in metres above bare earth at which a non-canopy cell is classified as built. Ignored when 'built' is supplied.

foliage_k

Numeric (optional), foliage extinction coefficient per metre. When it is NULL (the default) crowns are opaque and a ray stops at the first foliage it meets. When it is positive, crowns attenuate instead: optical depth accumulates along the path travelled inside the crown and each direction is partitioned between sky, foliage and surface rather than labelled with one of them (see details). Reasonable values are roughly 0.25 to 1.0 per metre.

t_min

Numeric (default is 0.001), transmittance below which a direction is treated as shut and the march stops. Only used when 'foliage_k' is set; it bounds the cost of rays that would otherwise run the full 'max_dist' through dense canopy.

method

Character, "equirectangular" (default) or "cylindrical". Solid-angle metrics are only returned for the equirectangular projection, which covers the whole sphere; see details.

sky_threshold

Numeric (default is 0), elevation tolerance in metres applied when testing the ray against the surface, as in pano_view(), where a positive value stops a coarse step from skipping a thin ridge or parapet. It has to be smaller than 'h', and the default is 0 for a reason: the test counts a hit as soon as the surface is within 'sky_threshold' of the ray, so a tolerance of 3 m with an eye at 1.7 m reports the ground at the first step in every direction and returns a panorama with no sky in it. The crown band is tested without tolerance, since a prism has thickness of its own.

step_size

Numeric (default is 0.5), ray marching step in DSM units.

max_dist

Numeric (default is 500), maximum ray distance in metres.

pano_dim

Integer vector of length 2, panorama height and width in pixels (default is c(256, 512)).

ares

Numeric, cylindrical vertical angular resolution in radians per pixel (default is 0.0174533, about one degree).

heading

Numeric, camera heading in degrees (default is 0 = north).

plot

Logical (default is FALSE), whether to plot the depth and class panoramas side by side.

Details

Pixels are not equal in solid angle. In an equirectangular panorama a row at latitude ϕ\phi subtends cos⁡ϕ\cos\phi times the solid angle of a row at the horizon, so counting pixels overweights the zenith and the nadir badly. Every fraction reported here is weighted by cos⁡ϕ\cos\phi, which is what makes 'svf' comparable with a sky view factor computed by any other means. The cylindrical projection samples latitude proportionally to tan⁡ϕ\tan\phi and does not cover the poles at all, so no fractions are returned for it; use it for a picture, not for a metric.

The classes are deliberately the ones a segmented street-level image gives you: sky, tree, building, road or ground. That makes the panorama directly comparable to a semantic segmentation of a matching photograph, which is the most defensible way to validate an under-canopy visibility model, since no amount of internal consistency proves the crown geometry is right.

'canopy_enclosure' is restricted to the horizon band because that is where foliage acts as a wall rather than as a ceiling, and the two are different experiences of the same trees. 'canopy_overhead' captures the ceiling.

With 'foliage_k' set, a crown attenuates rather than blocks. Optical depth accumulates along the path the ray travels inside the prism,

τ=k L,T=e−τ,\tau = k\,L, \qquad T = e^{-\tau},

the Beer-Lambert form standard in canopy radiative transfer, and each direction is partitioned instead of labelled: a ray reaching a surface at distance SS contributes T(S)T(S) to that surface and 1−T(S)1 - T(S) to foliage, while one that escapes contributes TT to sky and the remainder to foliage. The weights sum to 1 in every direction, so 'sky' becomes a sky view factor with vegetation transmissivity rather than a count of unobstructed pixels. Grazing a crown edge lets most of a direction through; looking through a whole crown lets almost none.

For a turbid-medium crown the coefficient would follow k=G⋅LADk = G \cdot LAD, where LADLAD is leaf area density in square metres of leaf per cubic metre and GG is the mean projection of unit leaf area, about 0.5 for a spherical leaf angle distribution, which puts broadleaf crowns at roughly 0.25 to 1.0 per metre. Do not take the value from that relation. The prism fills the entire column from crown base to canopy top with foliage, whereas a real crown is largely a shell and a DSM-derived canopy top over a patch is the upper envelope of several crowns rather than the lid of a solid block. The model therefore overstates the path length through leaves, sometimes by a large factor, and the coefficient that reproduces observed sky has to absorb that: on a 32% cover test scene, k = 0.5 per metre left a sky view factor of 0.006 and only k near 0.1 produced the 0.2 to 0.3 range a tree-lined street plausibly has.

Treat 'foliage_k' as an effective parameter to be calibrated, not a leaf-area measurement. It is one number, it acts monotonically on sky fraction, and the 'sky' layer is directly comparable to the soft sky mask of a segmented street-level photograph, so it can be fitted per species or per season given enough matched points. Reporting the fitted value alongside the metrics is more honest than adopting a physically derived one that the geometry does not support.

The same limitations apply as to compute_viewshed_understory(): a crown is a vertical prism over its cell footprint, the prism is fully opaque rather than attenuating, trunks are absent, and no curvature or refraction correction is applied. Left opaque, the crown model bites harder here than in the viewshed, because a panorama is mostly made of rising rays and those are the ones foliage intercepts, so 'canopy' is an upper bound and 'sky' a lower one; 'foliage_k' is what turns those bounds into an estimate.

Value

A list with two elements:

pano

a six-layer SpatRaster over the extent 0-360 by -90 to 90. 'depth' is the distance in map units to the first opaque obstruction, NA where none was met; note that with 'foliage_k' set, foliage is not opaque, so depth refers to the surface behind it. 'class' is the dominant contributor per direction: NA sky, 1 ground, 2 built, 3 canopy. 'sky', 'canopy', 'ground' and 'built' are the per-direction weights, which sum to 1 everywhere and are 0/1 indicators when crowns are opaque. The 'sky' layer is the one to compare against a soft sky mask from a segmented photograph.

metrics

a named numeric vector of solid-angle fractions, or NULL for the cylindrical projection. 'sky', 'ground', 'built' and 'canopy' are fractions of the full sphere and sum to 1; 'svf' is the sky fraction of the upper hemisphere; 'canopy_overhead' is the canopy fraction of the upper hemisphere; 'canopy_enclosure' is the canopy fraction of the band within 30 degrees of the horizon; and 'mean_depth' is the mean distance to the first hit, ignoring sky.

See Also

compute_viewshed_understory() road_corridor_mask() pano_view()

Examples

dem <- terra::rast(nrows = 201, ncols = 201, xmin = 0, xmax = 201,
                   ymin = 0, ymax = 201, crs = "EPSG:32617", vals = 0)
cnp <- terra::setValues(dem, 0)
m <- terra::as.matrix(cnp, wide = TRUE)
m[95:105, ] <- 1
cnp <- terra::setValues(cnp, m)
dsm <- dem + cnp * 14

out <- pano_view_understory(dsm = dsm, dem = dem, canopy = cnp,
                            vpt = c(100.5, 100.5), cbh = 3, max_dist = 100)
round(out$metrics, 3)

road_corridor_mask

Description

Experimental. Builds a binary open-space raster from road centrelines by buffering them to a corridor. Road networks are usually distributed as centrelines rather than as paved polygons, and the corridor is what compute_viewshed_understory() needs for its 'clear_mask': the strip of ground over which mapped canopy is really an overhanging crown rather than foliage the observer has to look through.

Usage

road_corridor_mask(roads, dsm, half_width = 6, cap = "ROUND")

Arguments

roads

sf or SpatVector, road centrelines. Polygons are accepted too and are used as they are, without buffering.

dsm

Raster whose geometry, resolution and CRS the output matches.

half_width

Numeric or Character (default is 6), half the corridor width in metres, i.e. the distance from the centreline to the edge. Supply a single number for a constant width, a numeric vector with one value per feature, or the name of a column in 'roads' holding one value per feature. Converted to the CRS unit of the DSM.

cap

Character (default is "ROUND"), how corridor ends are closed: "ROUND", "FLAT" or "SQUARE", passed to sf::st_buffer().

Details

Widths are best taken from the network itself. Many municipal centreline layers carry a right-of-way or pavement width field, which beats any constant; failing that, a per-class lookup is a reasonable second choice and the vector form of 'half_width' is meant for it:

hw <- c(motorway = 12, primary = 9, secondary = 8,
        residential = 6, service = 4)
m  <- road_corridor_mask(roads, dsm, half_width = hw[roads$highway])

Overlapping buffers are dissolved, so intersections come out as the blobs the union produces. That is roughly right, since intersections genuinely are more open than the streets feeding them.

Build the corridor across the whole analysis window rather than only near the viewpoint. A tree-lined street is a tunnel whose interior is open along its entire length, so the crowns that matter for a long view down the street are the ones overhanging the road far away, not the ones directly overhead.

Before reaching for this, check whether the city publishes a paved-surface or impervious-surface polygon layer. That is the road polygon this function approximates, and it is often filed under land cover rather than transportation.

Value

A single-layer SpatRaster aligned with 'dsm', 1 inside the corridor and 0 outside.

See Also

compute_viewshed_understory()

Examples

dsm <- terra::rast(nrows = 61, ncols = 61, xmin = 0, xmax = 61,
                   ymin = 0, ymax = 61, crs = "EPSG:32617", vals = 0)
line <- sf::st_sf(geometry = sf::st_sfc(
  sf::st_linestring(cbind(c(0, 61), c(30.5, 30.5))), crs = "EPSG:32617"))
corridor <- road_corridor_mask(line, dsm, half_width = 5)
terra::plot(corridor)

An S4 class to represent the viewshed

Description

A viewshed object contains a 'matrix' of visible and invisible area, resolution, extent, and crs

Slots

visible

matrix

resolution

vector

extent

numeric

crs

character


visual_magnitude

Description

Visual Magnitude quantifies the extent of a visible region as perceived by an observer. It is derived from the surface's slope and angle features, alongside the observer's relative distance from the area (Chamberlain & Meitner).

Usage

visual_magnitude(viewshed, dsm)

Arguments

viewshed

Viewshed object.

dsm

SpatRaster. A digital surface / elevation model

Value

SpatRaster

References

Chamberlain, B. C., & Meitner, M. J. (2013). A route-based visibility analysis for landscape management. Landscape and Urban Planning, 111, 13-24.

See Also

compute_viewshed()

Examples

# Load a viewpoint
test_viewpoint <- sf::read_sf(system.file("test_viewpoint.shp", package = "viewscape"))
# load dsm raster
dsm <- terra::rast(system.file("test_dsm.tif", package ="viewscape"))
# Compute viewshed
viewshed <- viewscape::compute_viewshed(dsm = dsm,
                                        viewpoints = test_viewpoint,
                                        offset_viewpoint = 6)
# Compute visual magnitude
vm <- viewscape::visual_magnitude(viewshed, dsm)

visualize_viewshed

Description

The visualize_viewshed function is designed for the visualization of a viewshed analysis, providing users with various options for visualizing the results. The function works with a viewshed object and offers multiple plotting and output types.

Usage

visualize_viewshed(viewshed, plottype = "", outputtype = "")

Arguments

viewshed

Viewshed object

plottype

Character, specifying the type of visualization ("polygon" or "raster").

outputtype

Character, specifying the type of output object ("raster" or "polygon").

Value

Visualized viewshed as either a raster or polygon object, depending on the outputtype specified.

Examples

# Load a viewpoint
test_viewpoint <- sf::read_sf(system.file("test_viewpoint.shp", package = "viewscape"))
# load dsm raster
dsm <- terra::rast(system.file("test_dsm.tif", package ="viewscape"))
#Compute viewshed
viewshed <- compute_viewshed(dsm = dsm,
                             viewpoints = test_viewpoint,
                             offset_viewpoint = 6)
# Visualize the viewshed as polygons
visualize_viewshed(viewshed, plottype = "polygon")
# Visualize the viewshed as a raster
visualize_viewshed(viewshed, plottype = "raster")
# Get the visualized viewshed as a polygon object
polygon_viewshed <- visualize_viewshed(viewshed,
                                       plottype = "polygon",
                                       outputtype = "polygon")