tickets/DM-27208

.travis.yml

@@ -26,9 +26,10 @@ install:
   # Use pip to install our rtd required packages
   - pip install -r docs/rtd-pip-requirements
   - python setup.py install
+
 # Run test
 script:
-  - pytest tests
+  - pytest tests --ignore=tests/test_regressions.py --ignore=tests/test_docs.py --ignore=tests/test_multiresolution.py
 
 # Calculate coverage
 #after_success:

MANIFEST.in

@@ -1 +1,2 @@
 include _version.txt
+include scarlet/testing/lookup.db

docs/0-quickstart.ipynb

@@ -49,7 +49,7 @@
     "filters = data[\"filters\"]\n",
     "catalog = data[\"catalog\"]\n",
     "weights = 1/data[\"variance\"]\n",
-    "psfs = scarlet.PSF(data[\"psfs\"])"
+    "psfs = scarlet.ImagePSF(data[\"psfs\"])"
    ]
   },
   {
@@ -113,8 +113,7 @@
    "metadata": {},
    "outputs": [],
    "source": [
-    "from functools import partial\n",
-    "model_psf = scarlet.PSF(partial(scarlet.psf.gaussian, sigma=.8), shape=(None, 8, 8))"
+    "model_psf = scarlet.GaussianPSF(sigma=(0.8,)*len(filters))"
    ]
   },
   {
@@ -200,9 +199,9 @@
     "\n",
     "You now need to define sources that are going to be fit. The full model, which we will call `Blend`, is a collection of those sources. We provide several pre-built source types, e.g.:\n",
     "\n",
-    "* `RandomSource` fit per-band amplitude and non-parametric morphology starting from uniform random draws for both.\n",
+    "* `RandomSource` fits per-band amplitude and non-parametric morphology starting from uniform random draws for both.\n",
     "* `PointSource` fits centers and per-band amplitude using the observed PSF model.\n",
-    "* `ExtendedSource` fits per-band amplitude and a non-parametric morphology (which can be constrained to be symmetric and/or monotonic with respect to the center).\n",
+    "* `ExtendedSource` fits per-band amplitude and a non-parametric morphology (which is constrained to be monotonically decreasing from the center).\n",
     "* `ExtendedSource(K=2)` splits an `ExtendedSource` into `K` components that are initially radially separated and stacked on top of each other like a pyramid.\n",
     "\n",
     "In our example, we assume *prior* knowledge that object 0 is a star, and object 1 should be modeled with two components. Everything else is assumed a single-component galaxy. **We generally recommend `ExtendedSource` as default** if additional information about the source is not available."
@@ -221,7 +220,7 @@
     "    elif k == 1:\n",
     "        new_source = scarlet.ExtendedSource(model_frame, (src['y'], src['x']), observation, K=2)\n",
     "    else:\n",
-    "        new_source = scarlet.ExtendedSource(model_frame, (src['y'], src['x']), observation)\n",
+    "        new_source = scarlet.ExtendedSource(model_frame, (src['y'], src['x']), observation, compact=True)\n",
     "    sources.append(new_source)"
    ]
   },
@@ -334,7 +333,7 @@
     "\n",
     "### Measure Fluxes\n",
     "\n",
-    "The color information in these plots stems from the per-band amplitude, which could be obtained as `source.spectrum`. However, it is more useful to compute per-band fluxes, which integrate over the morphology. The source spectra plots above have done exactly that. The convention of these fluxes is given by the units and ordering of the original data cube. In the case of multi-component sources, the fluxes of all components are combined."
+    "The color information in these plots stems from the per-band amplitude, which are computed from the hyperspectral model by integrating over the morphology. The source spectra plots in the right panels above have done exactly that. The convention of these fluxes is given by the units and ordering of the original data cube. In the case of multi-component sources, the fluxes of all components are combined."
    ]
   },
   {

docs/rtd-pip-requirements

@@ -12,3 +12,4 @@ prompt_toolkit>=2.0.10
 sep
 gitpython
 boto3
+gitpython

docs/tutorials.rst

@@ -8,4 +8,3 @@ Tutorials
    tutorials/point_source
    tutorials/wavelet_model
    tutorials/multiresolution
-

docs/tutorials/display.ipynb

@@ -36,7 +36,7 @@
     "catalog = data[\"catalog\"]\n",
     "filters = data[\"filters\"]\n",
     "weights = 1/data[\"variance\"]\n",
-    "psfs = scarlet.PSF(data[\"psfs\"])\n",
+    "psfs = scarlet.ImagePSF(data[\"psfs\"])\n",
     "\n",
     "observation = scarlet.Observation(\n",
     "    images, \n",

docs/tutorials/multiresolution.ipynb

@@ -55,7 +55,7 @@
     "# Load the HSC PSF data\n",
     "psf_hsc = fits.open('../../data/test_resampling/PSF_HSC.fits')[0].data\n",
     "Np1, Np2 = psf_hsc[0].shape\n",
-    "psf_hsc = scarlet.PSF(psf_hsc)\n",
+    "psf_hsc = scarlet.ImagePSF(psf_hsc)\n",
     "\n",
     "# Load the HST image data\n",
     "hst_hdu = fits.open('../../data/test_resampling/Cut_HST1.fits')\n",
@@ -66,7 +66,7 @@
     "# Load the HST PSF data\n",
     "psf_hst = fits.open('../../data/test_resampling/PSF_HST.fits')[0].data\n",
     "psf_hst = psf_hst[None,:,:]\n",
-    "psf_hst = scarlet.PSF(psf_hst)\n",
+    "psf_hst = scarlet.ImagePSF(psf_hst)\n",
     "\n",
     "# Scale the HST data\n",
     "n1,n2 = np.shape(data_hst)\n",

docs/tutorials/point_source.ipynb

@@ -21,7 +21,7 @@
     "import matplotlib\n",
     "import matplotlib.pyplot as plt\n",
     "# use a better colormap and don't interpolate the pixels\n",
-    "matplotlib.rc('image', cmap='inferno', interpolation='none')\n",
+    "matplotlib.rc('image', cmap='inferno', interpolation='none', origin='lower')\n",
     "\n",
     "import numpy as np\n",
     "import scarlet\n",
@@ -77,12 +77,11 @@
    "metadata": {},
    "outputs": [],
    "source": [
-    "from functools import partial\n",
-    "model_psf = scarlet.PSF(partial(scarlet.psf.gaussian, sigma=0.9), shape=(None,31,31))\n",
+    "model_psf = scarlet.GaussianPSF(sigma=0.9)\n",
     "model_frame = scarlet.Frame(images.shape, psfs=model_psf, channels=filters)\n",
     "\n",
     "observation = scarlet.Observation(images, \n",
-    "                                  psfs=scarlet.PSF(psfs), \n",
+    "                                  psfs=scarlet.ImagePSF(psfs), \n",
     "                                  weights=weights,\n",
     "                                  channels=filters)\n",
     "observation = observation.match(model_frame)"
@@ -169,7 +168,7 @@
    "metadata": {},
    "outputs": [],
    "source": [
-    "%time it, logL = blend.fit(100, e_rel=1e-4)\n",
+    "%time it, logL = blend.fit(200, e_rel=1e-4)\n",
     "print(f\"scarlet ran for {it} iterations to logL = {logL}\")\n",
     "scarlet.display.show_likelihood(blend)\n",
     "plt.show()"
@@ -200,15 +199,17 @@
     "                             norm=norm, \n",
     "                             observation=observation,\n",
     "                             show_rendered=True, \n",
-    "                             show_observed=True)\n",
+    "                             show_observed=True,\n",
+    "                             add_boxes=True\n",
+    "                            )\n",
     "plt.show()"
    ]
   },
   {
    "cell_type": "markdown",
    "metadata": {},
    "source": [
-    "We can see that the model overall performs well, especially in fitting the amplitude and locations of the point sources. It is worth noting that the stellar residuals have red rings and blue cores, which suggest that the PSF model of the observation is not quite correct across all bands. The extended sources are less convincing. The effects of the symmetry constraint are clearly visible, and the color between the two galaxies is not very different, so these two sources remain poorly separated."
+    "We can see that the model overall performs well, especially in fitting the amplitude and locations of the point sources. It is worth noting that the stellar residuals have red rings and blue cores, which suggest that the PSF model of the observation is not quite correct across all bands. The extended sources are less convincing. The colors of the two galaxies are not very different, so these two sources remain poorly separated."
    ]
   }
  ],

docs/tutorials/wavelet_model.ipynb

@@ -48,7 +48,7 @@
     "filters = data[\"filters\"]\n",
     "catalog = data[\"catalog\"]\n",
     "weights = 1/data[\"variance\"]\n",
-    "psfs = scarlet.PSF(data[\"psfs\"])"
+    "psfs = scarlet.ImagePSF(data[\"psfs\"])"
    ]
   },
   {
@@ -161,8 +161,7 @@
    "metadata": {},
    "outputs": [],
    "source": [
-    "from functools import partial\n",
-    "model_psf = scarlet.PSF(partial(scarlet.psf.gaussian, sigma=.8), shape=(None, 8, 8))"
+    "model_psf = scarlet.GaussianPSF(sigma=0.8)"
    ]
   },
   {

scarlet/__init__.py

@@ -13,6 +13,7 @@
 from . import display
 from . import initialization
 from . import measure
+from . import testing
 
 try:
     from ._version import version

scarlet/component.py

@@ -6,6 +6,8 @@
 from .constraint import PositivityConstraint
 from .bbox import Box, overlapped_slices
 from .fft import fast_zero_pad
+from .morphology import Morphology
+from .spectrum import Spectrum
 
 
 class Component(Model):
@@ -136,13 +138,10 @@ class FactorizedComponent(Component):
     """
 
     def __init__(self, frame, spectrum, morphology):
-        from .spectrum import Spectrum
-
         assert isinstance(spectrum, Spectrum)
-        from .morphology import Morphology
-
         assert isinstance(morphology, Morphology)
-        bbox = spectrum.bbox @ morphology.bbox
+
+        bbox = spectrum.bbox @ morphology.bbox[-2:]
 
         super().__init__(frame, children=[spectrum, morphology], bbox=bbox)
 
@@ -163,7 +162,12 @@ def get_model(self, *parameters, frame=None):
             (Channels, Height, Width) image of the model
         """
         spectrum, morphology = self.get_models_of_children(*parameters)
-        model = spectrum[:, None, None] * morphology[None, :, :]
+        if len(morphology.shape) == 2:
+            model = spectrum[:, None, None] * morphology[None, :, :]
+        elif len(morphology.shape) == 3:
+            model = spectrum[:, None, None] * morphology
+        else:
+            raise AttributeError("morphology must be 2D or 3D")
 
         # project the model into frame (if necessary)
         if frame is not None:

scarlet/display.py

@@ -2,6 +2,7 @@
 from astropy.visualization.lupton_rgb import LinearMapping, AsinhMapping
 import matplotlib.pyplot as plt
 from matplotlib.patches import Rectangle, Polygon
+from matplotlib.ticker import MaxNLocator
 from .bbox import Box
 from .component import Component
 
@@ -181,6 +182,7 @@ def show_likelihood(blend, figsize=None, **kwargs):
     fig, ax = plt.subplots(1, 1, figsize=figsize)
     ax.plot(blend.log_likelihood, **kwargs)
     ax.set_xlabel("Iteration")
+    ax.xaxis.set_major_locator(MaxNLocator(integer=True))
     ax.set_ylabel("log-Likelihood")
     return fig
 
@@ -230,7 +232,7 @@ def show_observation(
         psf_image = np.zeros(observation.images.shape)
 
         if observation.frame.psf is not None:
-            psf_model = observation.frame.psf.image.copy()  # get_model()
+            psf_model = observation.frame.psf.get_model()
             # make PSF as bright as the brightest pixel of the observation
             psf_model *= (
                 observation.images.mean(axis=0).max() / psf_model.mean(axis=0).max()