viewshed

Visible area from point on or above terrain

Syntax

``[vis,visR] = viewshed(Z,R,latObs,lonObs)``
``[vis,visR] = viewshed(Z,R,latObs,lonObs,hObs)``
``[vis,visR] = viewshed(Z,R,latObs,lonObs,hObs,hTarget)``
``[vis,visR] = viewshed(Z,R,latObs,lonObs,hObs,hTarget,hObsRef)``
``[vis,visR] = viewshed(Z,R,latObs,lonObs,hObs,hTarget,hObsRef,hTargetRef)``
``[vis,visR] = viewshed(Z,R,latObs,lonObs,hObs,hTarget,hObsRef,hTargetRef,rad)``
``[vis,visR] = viewshed(Z,R,latObs,lonObs,hObs,hTarget,hObsRef,hTargetRef,rad,effectiveRad)``

Description

Specify Coordinates and Heights

example

````[vis,visR] = viewshed(Z,R,latObs,lonObs)` calculates the visible area (the viewshed) for the point specified by `latObs` and `lonObs`. Specify spatially referenced terrain data, in meters, using `Z` and `R`. This function assumes that the Earth is a sphere.```
````[vis,visR] = viewshed(Z,R,latObs,lonObs,hObs)` specifies the height, in meters, of the observer.```
````[vis,visR] = viewshed(Z,R,latObs,lonObs,hObs,hTarget)` specifies the height, in meters, of the target points. The target points are all the points specified by the spatial reference.```
````[vis,visR] = viewshed(Z,R,latObs,lonObs,hObs,hTarget,hObsRef)` references the height of the observer to either the terrain (ground level) or the sphere (mean sea level).```
````[vis,visR] = viewshed(Z,R,latObs,lonObs,hObs,hTarget,hObsRef,hTargetRef)` references the height of the target points to either the terrain (ground level) or the sphere (mean sea level).```

Specify Reference Sphere

````[vis,visR] = viewshed(Z,R,latObs,lonObs,hObs,hTarget,hObsRef,hTargetRef,rad)` specifies the radius in meters of the reference sphere. This syntax is useful for finding line-of-sight visibility for planetary bodies other than Earth.```

example

````[vis,visR] = viewshed(Z,R,latObs,lonObs,hObs,hTarget,hObsRef,hTargetRef,rad,effectiveRad)` specifies a larger radius for the propagation of line-of-sight paths. You can use this syntax to account for the curvature of signal paths due to refraction in the atmosphere.```

Examples

collapse all

Create sample terrain data by using the `peaks` and `georefcells` functions.

```Z = 500*peaks(100); R = georefcells([-0.1 0],[0 0.1],size(Z));```

Calculate the viewshed for a sample point that is `200` meters above the terrain.

```latObs = -0.027; lonObs = 0.05; [vis,visR] = viewshed(Z,R,latObs,lonObs);```

Display the viewshed as a surface on an `axesm`-based map. Plot the sample point so that it appears above the terrain.

```figure axesm("globe","geoid",earthRadius) meshm(vis,visR,size(Z),Z) axis tight camposm(-10,-10,1e6) camupm(0,0) [row,col] = geographicToDiscrete(R,latObs,lonObs); h = Z(row,col); plot3m(latObs,lonObs,h+200,"ko","MarkerSize",7,"MarkerFaceColor","r")```

```colormap(flipud(summer(2))) brighten(0.75) shading interp camlight cb = lcolorbar(["Obscured" "Visible"]); cb.Position = [0.85 0.5 0.02 0.1];```

Display the area that is visible by radar for an aircraft flying above the terrain. You can use the `viewshed` function to model radio wave propagation in the atmosphere as straight lines on a sphere.

Load terrain elevation data for an area around the Korean peninsula. To make the water areas appear consistent on the plot, set all elevations below sea level (`korea5c<0`) to the same value.

```load korea5c korea5c(korea5c<0) = -1;```

Specify the location of the aircraft. Then, calculate visible areas for the aircraft, at `3000` meters above mean sea level, of the surface. Model radio propagation in the atmosphere by specifying an effective radius for the reference sphere.

```latObs = 34.0931; lonObs = 125.6578; hObs = 3000; hTarget = 0; rad = earthRadius; [vis1,vis1R] = viewshed(korea5c,korea5cR,latObs,lonObs,hObs,hTarget,"MSL","AGL",rad,4/3*rad);```

Set up a relief map for the area by adjusting the data aspect ratio, the plot box aspect ratio, and the lighting.

```figure worldmap(korea5c,korea5cR) setm(gca,"geoid",[1 0]) da = daspect; pba = pbaspect; da(3) = 7.5*pba(3)/da(3); daspect(da); camlight(90,5) camlight(0,5) material([0.25 0.8 0])```

Display the terrain data, the viewshed, and the aircraft on the map. Adjust the colormap so that the visible areas appear in blue and the obscured areas appear in red. Note that some ocean areas are obscured by an island.

```meshm(vis1,vis1R,size(korea5c),korea5c) plotm(latObs,lonObs,"ko","MarkerFaceColor","k") colormap([1 0 0; 0 0 1])```

Calculate the viewshed again, this time calculating visible areas for the aircraft at `3000` meters of other aircraft flying above it at `5000` meters. The visible area is larger, with some of the area obscured by the island and mountains.

```hTarget = 5000; [vis2,vis2R] = viewshed(korea5c,korea5cR,latObs,lonObs,hObs,hTarget,"MSL","MSL",rad,4/3*rad); clmo surface meshm(vis2,vis2R,size(korea5c),korea5c)```

Since R2024a

Calculate a viewshed using terrain from multiple adjacent terrain tiles by merging the tiles into one raster.

Read elevation data from two DTED files into the workspace. Merge the data by using the `mergetiles` function.

```[Z1,R1] = readgeoraster("n39_w106_3arc_v2.dt1",OutputType="double"); [Z2,R2] = readgeoraster("n40_w106_3arc_v2.dt1",OutputType="double"); [Z,R] = mergetiles(Z1,R1,Z2,R2);```

Specify the location of an observation point that is `100` meters above the terrain.

```lat = 39.9597; lon = -105.6883; h = 100;```

Calculate the viewshed for the observation point. The output indicates visible areas using `1` values and obscured areas using `0` values. To avoid plotting the obscured areas, replace the `0` values with `NaN` values.

```v = viewshed(Z,R,lat,lon,h); idx = (v == 0); v(idx) = NaN;```

Create a map of the region.

```figure usamap(Z,R)```

Display the elevation data as a surface. To make the elevation data appear beneath the viewshed, specify the surface heights using `0`s and specify the surface colors using the elevations. Change the colormap.

```zrs = zeros(R.RasterSize); geoshow(zrs,R,DisplayType="surface",CData=Z) colormap gray```

Display the viewshed as a blue surface.

`geoshow(v,R,DisplayType="surface",FaceColor="#0072BD")`

Display the observation point using a yellow circle marker.

```geoshow(lat,lon,ZData=h,DisplayType="point", ... Marker="o",MarkerFaceColor="#EDB120",MarkerEdgeColor="k")```

```c = colorbar; c.Label.String = "Elevation (m)";```

The elevation data used in this example is from the US Geological Survey.

Input Arguments

collapse all

Elevation data grid, in meters, specified as an m-by-n array.

Data Types: `single` | `double`

Spatial reference for `Z`, specified as a `GeographicCellsReference` or `GeographicPostingsReference` object. The `RasterSize` property of `R` must be consistent with `size(Z)`.

Latitude of the observer, in degrees, specified as a scalar.

If the latitude is outside the latitude limits of `R`, then the `viewshed` function issues a warning and sets all elements of `vis` to `0`. You can find the latitude limits of `R` by querying its `LatitudeLimits` property.

Data Types: `single` | `double`

Longitude of the observer, in degrees, specified as a scalar.

If the longitude is outside the latitude limits of `R`, then the `viewshed` function issues a warning and sets all elements of `vis` to `0`. You can find the longitude limits of `R` by querying its `LongitudeLimits` property.

Data Types: `single` | `double`

Height of the observer, in meters, specified as a scalar.

Data Types: `single` | `double`

Height of the target points, in meters, specified as a scalar.

Data Types: `single` | `double`

Height reference for the observer, specified as one of these options:

• `"AGL"` — Reference `hObs` to the terrain (ground level).

• `"MSL"` — Reference `hObs` to the sphere (mean sea level).

Data Types: `char` | `string`

Height reference for the target points, specified as one of these options:

• `"AGL"` — Reference the target points to the terrain (ground level).

• `"MSL"` — Reference the target points to the sphere (mean sea level).

Data Types: `char` | `string`

Radius of the reference sphere in meters, specified as a positive scalar.

Data Types: `single` | `double`

Effective radius of the reference sphere in meters, specified as a positive scalar.

This argument enables you to account for the curvature of signal paths due to refraction in the atmosphere. For example, you can treat radio propagation in the atmosphere as straight-line propagation on a sphere with 4/3 the radius of the Earth by specifying `rad` as `6371000` and `effectiveRad` as `4/3*6371000`.

To calculate line-of-sight visibility for a flat Earth, specify this argument as `Inf`.

If you do not specify this argument, then the function uses the value of `rad`.

Data Types: `single` | `double`

Output Arguments

collapse all

Visibility indicator, returned as an m-by-n array containing `0` and `1` values. The size of `vis` matches the size of `Z`.

• A value of `1` indicates that the observer has line-of-sight visibility with the target point.

• A value of `0` indicates that the line of sight between the observer and the target point is obscured by terrain.

Spatial reference for `vis`, returned as a `GeographicCellsReference` or `GeographicPostingsReference` object.

The output `visR` is equivalent to the input `R`.

Version History

Introduced before R2006a

expand all