Skip to content
Open
Show file tree
Hide file tree
Changes from all commits
Commits
File filter

Filter by extension

Filter by extension

Conversations
Failed to load comments.
Loading
Jump to
Jump to file
Failed to load files.
Loading
Diff view
Diff view
22 changes: 22 additions & 0 deletions README.md
Original file line number Diff line number Diff line change
Expand Up @@ -6,6 +6,28 @@ Ray tracing using CUDA, accessible from Python.

*Real-time viewshed analysis with GPU-accelerated ray tracing. Green areas are visible from the observer position (blue dot). Run `python examples/playground.py` to try it interactively.*

## Quick Start

```python
import rtxpy # registers the .rtx xarray accessor
import rioxarray

# Load a GeoTIFF DEM as an xarray DataArray
dem = rioxarray.open_rasterio('elevation.tif').squeeze()

# Move data to GPU
dem = dem.rtx.to_cupy()

# Compute hillshade with ray-traced shadows
hillshade = dem.rtx.hillshade(shadows=True)

# Compute viewshed from an observer location (pixel coordinates)
viewshed = dem.rtx.viewshed(x=500, y=300, observer_elev=2)

# Launch interactive 3D terrain explorer
dem.rtx.explore()
```

## Prerequisites

- NVIDIA GPU with RTX support (Maxwell architecture or newer)
Expand Down
518 changes: 162 additions & 356 deletions examples/cell_tower_instancing_demo.ipynb

Large diffs are not rendered by default.

278 changes: 262 additions & 16 deletions examples/oblique_3d_terrain_visualization.ipynb

Large diffs are not rendered by default.

128 changes: 119 additions & 9 deletions examples/playground.py
Original file line number Diff line number Diff line change
Expand Up @@ -3,14 +3,23 @@
This example demonstrates real-time terrain exploration using
GPU-accelerated ray tracing via the xarray accessor's explore mode.

Supports two modes:
python playground.py # single-layer elevation explorer
python playground.py --dataset # multi-layer Dataset with landcover overlay

In Dataset mode, press G to cycle between elevation and landcover coloring.

Requirements:
pip install rtxpy[analysis] matplotlib xarray rioxarray requests
"""

import numpy as np
import requests
import xarray as xr
from pathlib import Path

from xrspatial import slope, aspect, quantile

# Import rtxpy to register the .rtx accessor
import rtxpy

Expand Down Expand Up @@ -116,7 +125,7 @@ def load_terrain():
terrain = terrain[::2, ::2]

# Scale down elevation for visualization (optional)
terrain.data = terrain.data * 0.2
terrain.data = terrain.data * 0.1

# Ensure contiguous array before GPU transfer
terrain.data = np.ascontiguousarray(terrain.data)
Expand All @@ -134,16 +143,73 @@ def load_terrain():
return terrain


def make_landcover(terrain):
"""Derive a synthetic landcover classification from terrain elevation.

Classes (by Z-code):
100 Water / lake bottom — lowest 5 % of elevation
500 Barren rock — steepest 30 % of slopes
1000 Shrubland — default (moderate elevation, moderate slope)
1500 Forest — upper-mid elevation, gentle slope
2000 Alpine — highest 5 % of elevation

Parameters
----------
terrain : xarray.DataArray
Elevation data (numpy or cupy).

Returns
-------
numpy.ndarray
2-D float32 array of landcover class codes, same shape as terrain.
"""
data = terrain.data
if hasattr(data, 'get'):
data = data.get()
else:
data = np.asarray(data)

dy, dx = np.gradient(data)
slope = np.sqrt(dx**2 + dy**2)

p5 = np.nanpercentile(data, 5)
p70 = np.nanpercentile(data, 70)
p95 = np.nanpercentile(data, 95)
sp70 = np.nanpercentile(slope, 70)

lc = np.full(data.shape, 1000.0, dtype=np.float32) # shrubland
lc[data <= p5] = 100.0 # water
lc[slope > sp70] = 500.0 # barren rock
lc[(data > p70) & (slope <= sp70)] = 1500.0 # forest
lc[data >= p95] = 2000.0 # alpine

unique, counts = np.unique(lc, return_counts=True)
names = {100: 'Water', 500: 'Barren', 1000: 'Shrubland',
1500: 'Forest', 2000: 'Alpine'}
total = lc.size
print("Landcover classes:")
for val, cnt in zip(unique, counts):
print(f" {names.get(val, val):12s} ({val:5.0f}): {100*cnt/total:5.1f}%")

return lc


if __name__ == "__main__":
import argparse
parser = argparse.ArgumentParser(description="rtxpy terrain playground")
parser.add_argument("--dataset", action="store_true",
help="Build an xr.Dataset with landcover overlay (G to cycle)")
args = parser.parse_args()

# Load terrain data (downloads if needed)
terrain = load_terrain()

print("\nLaunching explore mode...")
print("Controls:")
print("\nControls:")
print(" W/S/A/D or Arrow keys: Move camera")
print(" Q/E or Page Up/Down: Move up/down")
print(" I/J/K/L: Look around")
print(" +/-: Adjust movement speed")
print(" G: Cycle overlay layers" if args.dataset else " G: Cycle geometry layers")
print(" O: Place observer (for viewshed)")
print(" V: Toggle viewshed (teal glow)")
print(" [/]: Adjust observer height")
Expand All @@ -153,10 +219,54 @@ def load_terrain():
print(" H: Toggle help overlay")
print(" X: Exit\n")

# Launch interactive explore mode
terrain.rtx.explore(
width=1024,
height=768,
render_scale=0.5
)
# Camera: south edge looking north
H, W = terrain.shape
elev_data = terrain.data
if hasattr(elev_data, 'get'):
elev_np = elev_data.get()
else:
elev_np = np.asarray(elev_data)
elev_max = float(np.nanmax(elev_np))
elev_mean = float(np.nanmean(elev_np))
diag = np.sqrt(H**2 + W**2)
start_pos = (W / 2, H * 1.05, elev_max + diag * 0.08)
look_target = (W / 2, H / 2, elev_mean)

if args.dataset:
import cupy as cp

print("Building Dataset with landcover overlay...")
lc = make_landcover(terrain)

coords = {d: terrain.coords[d] for d in terrain.dims}

ds = xr.Dataset({
'elevation': terrain.rename(None),
'landcover': xr.DataArray(cp.asarray(lc), dims=terrain.dims, coords=coords),
'slope': slope(terrain),
'aspect': aspect(terrain),
'quantile': quantile(terrain),
})
print(ds)
print("\nLaunching Dataset explore (press G to cycle elevation ↔ landcover)...\n")
ds.rtx.explore(
z='elevation',
mesh_type='voxel',
width=1024,
height=768,
render_scale=0.5,
start_position=start_pos,
look_at=look_target,
)
else:
print("Launching explore mode...\n")
terrain.rtx.explore(
mesh_type='voxel',
width=1024,
height=768,
render_scale=0.5,
start_position=start_pos,
look_at=look_target,
)

print("Done")
Loading
Loading