healpix_analyse.decomp ====================== .. py:module:: healpix_analyse.decomp .. autoapi-nested-parse:: Exactly reconstructing multiscale decomposition for HEALPix maps. ``HealPixDecomp`` builds a local Laplacian pyramid from matched :class:`HealPixDown` and :class:`HealPixUp` operators. Every detail band is a fine map minus the prediction obtained by downsampling and upsampling it. The same stored Up operator is reused during synthesis, making reconstruction algebraically exact up to floating-point round-off, including on masked maps. Classes ------- .. autoapisummary:: healpix_analyse.decomp.HealPixPyramid healpix_analyse.decomp.HealPixWeightedPyramid healpix_analyse.decomp.HealPixDecomp Module Contents --------------- .. py:class:: HealPixPyramid Bases: :py:obj:`Sequence`\ [\ :py:obj:`ArrayLike`\ ] List-like multiscale coefficients returned by :meth:`HealPixDecomp.compute`. ``bands[:-1]`` are detail maps from fine to coarse and ``bands[-1]`` is the final low-pass map. ``cell_ids[i]`` and ``levels[i]`` describe the pixel domain of ``bands[i]``. The object behaves like a read-only sequence, so existing list-oriented code can iterate over it directly. .. py:attribute:: bands :type: tuple[ArrayLike, ...] .. py:attribute:: cell_ids :type: tuple[numpy.ndarray, ...] .. py:attribute:: levels :type: tuple[int, ...] .. py:property:: details :type: tuple[ArrayLike, ...] All band-pass/detail maps, ordered from fine to coarse. .. py:property:: coarse :type: ArrayLike Final low-pass scaling map. .. py:property:: images :type: tuple[ArrayLike, ...] Alias for :attr:`bands`. .. py:class:: HealPixWeightedPyramid Paired analysis of confidence-weighted data and confidence weights. Returned by :meth:`HealPixDecomp.compute_weighted`. ``q`` is the pyramid obtained by running the *same* linear analysis (:class:`HealPixDown`/ :class:`HealPixUp` chain) on ``weights * x`` (with missing values replaced by zero), and ``m`` is the pyramid obtained by running that identical analysis on ``weights`` alone. Because the analysis operator ``W`` is data-independent, ``q`` and ``m`` share ``cell_ids``/``levels`` band for band. The detail bands of ``m`` are *signed correction terms*, not per-band confidences: they are not bounded in ``[0, 1]`` and must not be clipped, thresholded, or divided into individually. Only a full synthesis of both channels followed by a single division -- performed by :meth:`HealPixDecomp.invert` when given a ``HealPixWeightedPyramid``, or by a kernel-pyramid convolution applied identically to both channels -- recovers a valid normalized result. .. py:attribute:: q :type: HealPixPyramid .. py:attribute:: m :type: HealPixPyramid .. py:property:: cell_ids :type: tuple[numpy.ndarray, ...] .. py:property:: levels :type: tuple[int, ...] .. py:class:: HealPixDecomp(level: int, cell_ids=None, Jmax: int = -1, *, ellipsoid: str = 'WGS84', radius_deg: Optional[float] = None, sigma_deg: Optional[float] = None, weight_norm: str = 'l1', up_norm: str = 'col_l1', dtype: torch.dtype = torch.float32, device: Optional[Union[str, torch.device]] = None) Bases: :py:obj:`torch.nn.Module` Local, exactly reconstructing HEALPix multiscale decomposition. :Parameters: * **level** (:py:class:`int`) -- Fine input Grid4Earth/HEALPix level, with ``nside = 2**level``. * **cell_ids** (:py:class:`array-like`, *optional*) -- NESTED identifiers of a masked/partial map, in the same order as the input data. ``None`` selects the complete sphere. * **Jmax** (:py:class:`int`, *default* ``-1``) -- Number of Down operations. ``-1`` uses every available scale down to level zero. ``0`` returns a pyramid containing only the input map. * **ellipsoid** (:py:class:`str`, *default* ``"WGS84"``) -- Geometry used by every analysis filter. * **radius_deg, sigma_deg** (:py:class:`float`, *optional*) -- Gaussian smoothing parameters passed to every ``HealPixDown``. When omitted, each level uses its own resolution-dependent defaults. * **weight_norm** (``{"l1", "l2", "none"}``, *default* ``"l1"``) -- Down-filter normalisation. ``l1`` preserves constant maps locally. * **up_norm** (``{"adjoint", "col_l1", "diag_l2"}``, *default* ``"col_l1"``) -- Matched synthesis normalisation. Exact pyramid reconstruction does not depend on this choice because the prediction is stored as a residual, but it changes the interpretation of individual bands. * **dtype, device** (**optional**) -- PyTorch execution options. CUDA is selected when available. .. admonition:: Notes At scale ``j`` the analysis is ``coarse[j+1] = Down[j](coarse[j])`` ``detail[j] = coarse[j] - Up[j](coarse[j+1])``. Synthesis applies the exact reverse recursion ``coarse[j] = detail[j] + Up[j](coarse[j+1])``. Every ``HealPixUp`` is constructed with ``paired_down=...`` and therefore reuses the exact sparse matrix and fine-cell domain of its analysis step. .. py:attribute:: level .. py:attribute:: Jmax :value: -1 .. py:attribute:: n_scales .. py:attribute:: n_bands .. py:attribute:: ellipsoid :value: '' .. py:attribute:: weight_norm :value: '' .. py:attribute:: up_norm :value: '' .. py:attribute:: dtype :value: Ellipsis .. py:attribute:: partial .. py:attribute:: n_pixels .. py:attribute:: down_layers .. py:attribute:: up_layers .. py:attribute:: levels .. py:attribute:: cell_ids_per_scale .. py:attribute:: sizes .. py:property:: device :type: torch.device Current device of the decomposition operators. .. py:property:: cell_ids :type: tuple[numpy.ndarray, ...] Canonical NESTED identifiers for every returned band. .. py:method:: compute(data: ArrayLike) -> HealPixPyramid Decompose a map into fine-to-coarse details and one coarse map. The returned object is list-like and contains ``n_scales + 1`` images. Arbitrary leading dimensions are supported, for example ``[B,C,N]``. .. py:method:: compute_weighted(data: ArrayLike, weights: Optional[ArrayLike] = None) -> HealPixWeightedPyramid Decompose a map that may contain NaN/missing values. Non-finite entries of ``data`` (and, when given, non-finite or implicitly-absent entries of ``weights``) are treated as missing. ``weights`` defaults to a 0/1 confidence mask built from ``isfinite(data)`` when omitted, i.e. plain unweighted missing-data handling. Internally this runs the *same* analysis operator ``W`` (the stored Down/Up chain) on ``q = weights * data`` (missing values zeroed) and on ``m = weights`` (also zeroed at missing positions), and returns both resulting pyramids paired in a :class:`HealPixWeightedPyramid`. Reconstruct with :meth:`invert`, which performs the single division ``S(q) / S(m)`` after synthesis, per the normalized-convolution convention -- never divide the two pyramids band by band. The boolean missing/valid decision itself carries no gradient (as with any masking based on ``isfinite``); the surviving finite values of ``data``/``weights`` remain differentiable. .. py:method:: forward(data: ArrayLike) -> HealPixPyramid Alias for :meth:`compute`, enabling normal ``nn.Module`` usage. .. py:method:: invert(decomposition: Union[HealPixPyramid, HealPixWeightedPyramid, Sequence[ArrayLike]], *, restore_mask: bool = True, eps: float = 1e-08) -> ArrayLike Reconstruct the original map from multiscale coefficients. Given a plain :class:`HealPixPyramid` (or a bare sequence of bands), this is the exact algebraic reverse of :meth:`compute`, unaffected by ``restore_mask``/``eps``. Given a :class:`HealPixWeightedPyramid` (from :meth:`compute_weighted`), both channels are synthesized with the *same* Up chain and combined with a single division ``y = S(q) / S(m)`` -- never per-band. Output pixels whose synthesized weight ``S(m)`` falls at or below ``eps`` have no reliable support; when ``restore_mask`` is True (default) they are set to NaN, otherwise to 0. .. py:method:: expand(decomposition: Union[HealPixPyramid, Sequence[ArrayLike]]) -> tuple[ArrayLike, ...] Lift every band to the fine grid so their direct sum reconstructs. ``compute`` deliberately returns compact native-resolution bands. This method applies the stored synthesis operators to each band independently, restores the original fine-cell order, and returns ``n_bands`` fine-resolution components. Their elementwise sum equals :meth:`invert` up to floating-point round-off. .. py:method:: extra_repr() -> str Return the extra representation of the module. To print customized extra information, you should re-implement this method in your own modules. Both single-line and multi-line strings are acceptable.