Skip to content

Commit

Permalink
Merge pull request #1584 from pierotofy/hillcons
Browse files Browse the repository at this point in the history
More consistent hillshading across zoom levels
  • Loading branch information
pierotofy authored Jan 11, 2025
2 parents 773a4e4 + c6225ff commit 743f3f4
Show file tree
Hide file tree
Showing 3 changed files with 8 additions and 14 deletions.
14 changes: 4 additions & 10 deletions app/api/hillshade.py
Original file line number Diff line number Diff line change
Expand Up @@ -71,12 +71,6 @@ def hillshade(self, elevation, vert_exag=1, dx=1, dy=1, fraction=1.):
A 2d array of illumination values between 0-1, where 0 is
completely in shadow and 1 is completely illuminated.
"""

# Because most image and raster GIS data has the first row in the array
# as the "top" of the image, dy is implicitly negative. This is
# consistent to what `imshow` assumes, as well.
dy = -dy

# compute the normal vectors from the partial derivatives
e_dy, e_dx = np.gradient(vert_exag * elevation, dy, dx)

Expand Down Expand Up @@ -115,19 +109,19 @@ def shade_normals(self, normals, fraction=1.):
intensity = normals.dot(self.direction.astype(np.float32))

# Apply contrast stretch
imin, imax = intensity.min(), intensity.max()
# imin, imax = np.nanmin(intensity), np.nanmax(intensity)
intensity *= fraction

# Rescale to 0-1, keeping range before contrast stretch
# If constant slope, keep relative scaling (i.e. flat should be 0.5,
# fully occluded 0, etc.)
if (imax - imin) > 1e-6:
# if (imax - imin) > 1e-6:
# Strictly speaking, this is incorrect. Negative values should be
# clipped to 0 because they're fully occluded. However, rescaling
# in this manner is consistent with the previous implementation and
# visually appears better than a "hard" clip.
intensity -= imin
intensity /= (imax - imin)
# intensity -= imin
# intensity /= (imax - imin)
intensity = np.clip(intensity, 0, 1)

return intensity
4 changes: 2 additions & 2 deletions app/api/tiler.py
Original file line number Diff line number Diff line change
Expand Up @@ -463,9 +463,9 @@ def get(self, request, pk=None, project_pk=None, tile_type="", z="", x="", y="",
if tile.data.shape[0] != 1:
raise exceptions.ValidationError(
_("Cannot compute hillshade of non-elevation raster (multiple bands found)"))
delta_scale = (maxzoom + ZOOM_EXTRA_LEVELS + 1 - z) * 4
delta_scale = (maxzoom + ZOOM_EXTRA_LEVELS + 1 - z) ** 2
dx = src.dataset.meta["transform"][0] * delta_scale
dy = -src.dataset.meta["transform"][4] * delta_scale
dy = src.dataset.meta["transform"][4] * delta_scale
ls = LightSource(azdeg=315, altdeg=45)

# Remove elevation data from edge buffer tiles
Expand Down
4 changes: 2 additions & 2 deletions app/raster_utils.py
Original file line number Diff line number Diff line change
Expand Up @@ -233,9 +233,9 @@ def write_band(arr, dst, band_num):

intensity = None
if hillshade is not None and hillshade > 0:
delta_scale = (ZOOM_EXTRA_LEVELS + 1) * 4
delta_scale = ZOOM_EXTRA_LEVELS ** 2
dx = src.meta["transform"][0] * delta_scale
dy = -src.meta["transform"][4] * delta_scale
dy = src.meta["transform"][4] * delta_scale
ls = LightSource(azdeg=315, altdeg=45)
intensity = ls.hillshade(arr[0], dx=dx, dy=dy, vert_exag=hillshade)
intensity = intensity * 255.0
Expand Down

0 comments on commit 743f3f4

Please sign in to comment.