Loading [MathJax]/extensions/tex2jax.js
LSST Applications g0d33ba9806+3d1aa5cd78,g0fba68d861+6c07529581,g1e78f5e6d3+aa4d1c21f1,g1ec0fe41b4+f536777771,g1fd858c14a+c6963eae98,g35bb328faa+fcb1d3bbc8,g4af146b050+58cb980876,g4d2262a081+bfae794ebc,g53246c7159+fcb1d3bbc8,g5a012ec0e7+b20b785ecb,g60b5630c4e+3d1aa5cd78,g67b6fd64d1+4086c0989b,g78460c75b0+2f9a1b4bcd,g786e29fd12+cf7ec2a62a,g7b71ed6315+fcb1d3bbc8,g87b7deb4dc+a8b896a16a,g8852436030+75f93ca278,g89139ef638+4086c0989b,g9125e01d80+fcb1d3bbc8,g94187f82dc+3d1aa5cd78,g989de1cb63+4086c0989b,g9f33ca652e+94dd9a85be,g9f7030ddb1+3f50642ad9,ga2b97cdc51+3d1aa5cd78,ga44b1db4f6+7d2d5e68ea,gabe3b4be73+1e0a283bba,gabf8522325+fa80ff7197,gb1101e3267+96cb2ddcf2,gb58c049af0+f03b321e39,gb89ab40317+4086c0989b,gcf25f946ba+75f93ca278,gd6cbbdb0b4+af3c3595f5,gd9a9a58781+fcb1d3bbc8,gde0f65d7ad+ed90e8109d,ge278dab8ac+d65b3c2b70,ge410e46f29+4086c0989b,gf67bdafdda+4086c0989b,gfe06eef73a+e2ab3e8e4f,v29.0.0.rc5
LSST Data Management Base Package
All Classes Namespaces Files Functions Variables Typedefs Enumerations Enumerator Properties Friends Macros Modules Pages
free_form.py
Go to the documentation of this file.
1# This file is part of scarlet_lite.
2#
3# Developed for the LSST Data Management System.
4# This product includes software developed by the LSST Project
5# (https://www.lsst.org).
6# See the COPYRIGHT file at the top-level directory of this distribution
7# for details of code ownership.
8#
9# This program is free software: you can redistribute it and/or modify
10# it under the terms of the GNU General Public License as published by
11# the Free Software Foundation, either version 3 of the License, or
12# (at your option) any later version.
13#
14# This program is distributed in the hope that it will be useful,
15# but WITHOUT ANY WARRANTY; without even the implied warranty of
16# MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the
17# GNU General Public License for more details.
18#
19# You should have received a copy of the GNU General Public License
20# along with this program. If not, see <https://www.gnu.org/licenses/>.
21
22__all__ = ["FactorizedFreeFormComponent"]
23
24from typing import Callable, cast
25
26import numpy as np
27from lsst.scarlet.lite.detect_pybind11 import get_connected_multipeak, get_footprints # type: ignore
28
29from ..bbox import Box
30from ..component import Component, FactorizedComponent
31from ..detect import footprints_to_image
32from ..image import Image
33from ..parameters import Parameter, parameter
34
35
37 """Implements a free-form component
38
39 With no constraints this component is typically either a garbage collector,
40 or part of a set of components to deconvolve an image by separating out
41 the different spectral components.
42
43 See `FactorizedComponent` for a list of parameters not shown here.
44
45 Parameters
46 ----------
47 peaks: `list` of `tuple`
48 A set of ``(cy, cx)`` peaks for detected sources.
49 If peak is not ``None`` then only pixels in the same "footprint"
50 as one of the peaks are included in the morphology.
51 If `peaks` is ``None`` then there is no constraint applied.
52 min_area: float
53 The minimum area for a peak.
54 If `min_area` is not `None` then all regions of the morphology
55 with fewer than `min_area` connected pixels are removed.
56 """
57
59 self,
60 bands: tuple,
61 spectrum: np.ndarray | Parameter,
62 morph: np.ndarray | Parameter,
63 model_bbox: Box,
64 bg_thresh: float | None = None,
65 bg_rms: np.ndarray | None = None,
66 floor: float = 1e-20,
67 peaks: list[tuple[int, int]] | None = None,
68 min_area: float = 0,
69 ):
70 super().__init__(
71 bands=bands,
72 spectrum=spectrum,
73 morph=morph,
74 bbox=model_bbox,
75 peak=None,
76 bg_rms=bg_rms,
77 bg_thresh=bg_thresh,
78 floor=floor,
79 )
80
81 self.peaks = peaks
82 self.min_area = min_area
83
84 def prox_spectrum(self, spectrum: np.ndarray) -> np.ndarray:
85 """Apply a prox-like update to the spectrum
86
87 This differs from `FactorizedComponent` because an
88 `SedComponent` has the spectrum normalized to unity.
89 """
90 # prevent divergent spectrum
91 spectrum[spectrum < self.floor] = self.floor
92 # Normalize the spectrum
93 spectrum = spectrum / np.sum(spectrum)
94 return spectrum
95
96 def prox_morph(self, morph: np.ndarray) -> np.ndarray:
97 """Apply a prox-like update to the morphology
98
99 This is the main difference between an `SedComponent` and a
100 `FactorizedComponent`, since this component has fewer constraints.
101 """
102 from lsst.scarlet.lite.detect_pybind11 import get_connected_multipeak, get_footprints # type: ignore
103
104 if self.bg_thresh is not None and isinstance(self.bg_rms, np.ndarray):
105 bg_thresh = self.bg_rms * self.bg_thresh
106 # Enforce background thresholding
107 model = self.spectrum[:, None, None] * morph[None, :, :]
108 morph[np.all(model < bg_thresh[:, None, None], axis=0)] = 0
109 else:
110 # enforce positivity
111 morph[morph < 0] = 0
112
113 if self.peaks is not None:
114 footprint = get_connected_multipeak(morph, self.peaks, 0)
115 morph = morph * footprint
116
117 if self.min_area > 0:
118 footprints = get_footprints(morph, 4.0, self.min_area, 0, 0, False)
119 bbox = self.bbox.copy()
120 bbox.origin = (0, 0)
121 footprint_image = footprints_to_image(footprints, bbox)
122 morph = morph * (footprint_image > 0).data
123
124 if np.all(morph == 0):
125 morph[0, 0] = self.floor
126
127 return morph
128
129 def resize(self, model_box: Box) -> bool:
130 return False
131
132 def __str__(self):
133 return (
134 f"FactorizedFreeFormComponent(\n bands={self.bands}\n "
135 f"spectrum={self.spectrum})\n center={self.peak}\n "
136 f"morph_shape={self.morph.shape}"
137 )
138
139 def __repr__(self):
140 return self.__str__()
141
142
144 """Implements a component with no spectral or monotonicty constraints
145
146 This is a FreeFormComponent that is not factorized into a
147 spectrum and morphology with no monotonicity constraint.
148 """
149
151 self,
152 bands: tuple,
153 model: np.ndarray | Parameter,
154 model_bbox: Box,
155 bg_thresh: float | None = None,
156 bg_rms: np.ndarray | None = None,
157 floor: float = 1e-20,
158 peaks: list[tuple[int, int]] | None = None,
159 min_area: float = 0,
160 ):
161 if len(bands) != 1:
162 raise ValueError(f"{type(self)} only supports one band")
163 super().__init__(bands=bands, bbox=model_bbox)
164 self._model = parameter(model)
165 self.bg_rms = bg_rms
166 self.bg_thresh = bg_thresh
167 self.floor = floor
168 self.peaks = peaks
169 self.min_area = min_area
170
171 @property
172 def model(self) -> np.ndarray:
173 return self._model.x
174
175 def get_model(self) -> Image:
176 return Image(self.model, bands=self.bands, yx0=cast(tuple[int, int], self.bbox.origin))
177
178 @property
179 def shape(self) -> tuple:
180 return self.model.shape
181
182 def grad_model(self, input_grad: np.ndarray, model: np.ndarray) -> np.ndarray:
183 return input_grad
184
185 def prox_model(self, model: np.ndarray) -> np.ndarray:
186 if self.bg_thresh is not None and isinstance(self.bg_rms, np.ndarray):
187 bg_thresh = self.bg_rms * self.bg_thresh
188 # Enforce background thresholding
189 model[model < bg_thresh[:, None, None]] = 0
190 else:
191 # enforce positivity
192 model[model < 0] = 0
193
194 if self.peaks is not None:
195 # Remove pixels not connected to one of the peaks
196 model2d = np.sum(model, axis=0)
197 footprint = get_connected_multipeak(model2d, self.peaks, 0)
198 model = model * footprint[None, :, :]
199
200 if self.min_area > 0:
201 # Remove regions with fewer than min_area connected pixels
202 model2d = np.sum(model, axis=0)
203 footprints = get_footprints(model2d, 4.0, self.min_area, 0, 0, False)
204 bbox = self.bbox.copy()
205 bbox.origin = (0, 0)
206 footprint_image = footprints_to_image(footprints, bbox)
207 model = model * (footprint_image > 0).data[None, :, :]
208
209 if np.all(model == 0):
210 # If the model is all zeros, set a single pixel to the floor
211 model[0, 0] = self.floor
212
213 return model
214
215 def resize(self, model_box: Box) -> bool:
216 return False
217
218 def update(self, it: int, grad_log_likelihood: np.ndarray):
219 self._model.update(it, grad_log_likelihood)
220
221 def parameterize(self, parameterization: Callable) -> None:
222 """Convert the component parameter arrays into Parameter instances
223
224 Parameters
225 ----------
226 parameterization: Callable
227 A function to use to convert parameters of a given type into
228 a `Parameter` in place. It should take a single argument that
229 is the `Component` or `Source` that is to be parameterized.
230 """
231 # Update the spectrum and morph in place
232 parameterization(self)
233 # update the parameters
234 self._model.grad = self.grad_model
235 self._model.prox = self.prox_model
236
237 def __str__(self):
238 result = f"FreeFormComponent<bands={self.bands}, shape={self.shape}>"
239 return result
240
241 def __repr__(self):
242 return self.__str__()
A class to represent a 2-dimensional array of pixels.
Definition Image.h:51
__init__(self, tuple bands, np.ndarray|Parameter spectrum, np.ndarray|Parameter morph, Box model_bbox, float|None bg_thresh=None, np.ndarray|None bg_rms=None, float floor=1e-20, list[tuple[int, int]]|None peaks=None, float min_area=0)
Definition free_form.py:69
np.ndarray prox_spectrum(self, np.ndarray spectrum)
Definition free_form.py:84
np.ndarray prox_model(self, np.ndarray model)
Definition free_form.py:185
__init__(self, tuple bands, np.ndarray|Parameter model, Box model_bbox, float|None bg_thresh=None, np.ndarray|None bg_rms=None, float floor=1e-20, list[tuple[int, int]]|None peaks=None, float min_area=0)
Definition free_form.py:160
np.ndarray grad_model(self, np.ndarray input_grad, np.ndarray model)
Definition free_form.py:182
update(self, int it, np.ndarray grad_log_likelihood)
Definition free_form.py:218
None parameterize(self, Callable parameterization)
Definition free_form.py:221
MatrixB get_connected_multipeak(py::EigenDRef< const M > image, const std::vector< std::vector< int > > &centers, const double thresh=0)
Proximal operator to trim pixels not connected to one of the source centers.
std::vector< Footprint > get_footprints(py::EigenDRef< const M > image, const double min_separation, const int min_area, const double peak_thresh, const double footprint_thresh, const bool find_peaks=true, const int y0=0, const int x0=0)
Get all footprints in an image.