Auto-contrast

tilemaker/analysis/histogram.py

@@ -41,15 +41,8 @@ def build(
 
         timing_start = perf_counter()
 
-        start = layer.vmin * 4
-        end = layer.vmax * 4
-        bins = 128
-
-        log = log.bind(start=start, end=end, bins=bins)
-
-        edges = np.linspace(start, end, bins + 1)
-        counts = np.zeros(bins)
-
+        # Get the two tiles in the top (coarsest) layer into memory
+        tdata = []
         for tile_x in [0, 1]:
             tile, pushable = tiles.pull(
                 PullableTile(
@@ -63,9 +56,43 @@ def build(
             )
 
             tiles.push(pushable)
+            tdata.append(tile.data)
+
+        tdata = np.asarray(tdata)
+
+        # If vmin or vmax is auto, determine it automatically
+        # from the distribution of pixels in the top-level maps
+        if layer.vmin == "auto" or layer.vmax == "auto":
+            quantile = 0.01  # TODO: Raise this to a configuration/UI option
+            # The following logic is copied from pixell/enplot.py
+            vals = np.sort(tdata[np.isfinite(tdata)])
+            n = len(vals)
+            if n == 0:
+                raise ValueError("No finite values in the map.")
+            i = min(n - 1, int(round(n * quantile)))
+            v1, v2 = vals[i], vals[n - 1 - i]
+            # Avoid division by zero later, in case min and max are the same
+            if v2 == v1:
+                (v1, v2) = (v1 - 1, v2 + 1)
+
+        if layer.vmin == "auto":
+            layer.vmin = v1.item()
+        if layer.vmax == "auto":
+            layer.vmax = v2.item()
+
+        # Construct the histogram
+        start = layer.vmin * 4
+        end = layer.vmax * 4
+        bins = 128
+
+        log = log.bind(start=start, end=end, bins=bins)
 
-            if tile.data is not None:
-                counts += np.histogram(tile.data, bins=edges)[0]
+        edges = np.linspace(start, end, bins + 1)
+        counts = np.zeros(bins)
+
+        for tile_x in [0, 1]:
+            if tdata[tile_x] is not None:
+                counts += np.histogram(tdata[tile_x], bins=edges)[0]
 
         timing_end = perf_counter()
         log = log.bind(dt=timing_end - timing_start)

tilemaker/metadata/definitions.py

@@ -192,8 +192,8 @@ class Layer(AuthenticatedModel):
     number_of_levels: int | None = None
     tile_size: int | None = None
 
-    vmin: float | None = None
-    vmax: float | None = None
+    vmin: float | Literal["auto"] | None = 'auto'
+    vmax: float | Literal["auto"] | None = 'auto'
     cmap: str | None = None
 
     def model_post_init(self, _):

tilemaker/metadata/generation.py

@@ -5,7 +5,7 @@
 import os
 from hashlib import md5
 from pathlib import Path
-from typing import Any
+from typing import Any, Literal
 
 import structlog
 from astropy import units
@@ -88,8 +88,8 @@ def layers_from_fits(
                         description="Unknown layer",
                         quantity=None,
                         units=None,
-                        vmin=-1.0,
-                        vmax=1.0,
+                        vmin='auto',
+                        vmax='auto',
                         cmap="viridis",
                         index=None,
                     )
@@ -129,8 +129,8 @@ class ProtoLayer(BaseModel):
     description: str | None = None
     quantity: str | None = None
     units: str | None = None
-    vmin: float | None = None
-    vmax: float | None = None
+    vmin: float | Literal["auto"] | None = 'auto'
+    vmax: float | Literal["auto"] | None = 'auto'
     cmap: str | None = None
     index: int | None = None