Select bounding box for each peak about scarlet HOT 4 CLOSED

I have an idea about how to implement number 1, outlined in this ipython notebook, which uses a real HSC image to illustrate how this might work. It doesn't actually solve for the solution, but shows how one might go about choosing the bounding boxes for each object.

The basic idea is that the SED in neighboring pixels changes very little, and in fact differ very little from a peak except in blended regions. So we can mask out the image except for the region that encloses a peak, and then allow that region to grow by some reasonable factor to ensure that all of the flux is enclosed. This works well for the complicated blend I tested it on but more testing needs to be done to test the validity of the method. For example, some of the objects are given significantly larger regions because nearby sources have similar colors, but the important thing is that (at least in this example) none of the objects appear to be missing any flux. It would be interesting to see how this method works on spiral galaxies that have multiple colors.

from scarlet.

That looks quite good already. I don't fully understand the logic that decides the size/shape of the bounding box. However, this is not really what this ticket is about. Even if we had bounding boxes, we wouldn't have a mechanism to use them.

from scarlet.

Here's a different thought that should solve all of the problems we've had about separable objects
run bSDMM with a extended set of variables [S1,S2,...,SK,A], where Si is the spatial part of i-th object. This way they can be treated independently (in whatever boxes we want to have) and the convergence can be checked for any of them separately. This last part should be very good for diagnostic purposes.

This requires a few changes.

1. deblender calls proxmin.bsdmm directly.
2. a new class for building the model and calculating likelihood gradients.
3. initialization of constraints for S needs to be split.

Second item will be most important. There are several redundant computations for the likelihood gradients. By pulling those together in a class, we can store the results of previous calculations and control when updates to internal variables (such as the fully assembled S matrix, or the A matrix) are done. It would also be the interface for the user to inspect and render the deblender outputs.

So, the action of the deblend method would be to take a list/tree of peaks (with constraints, i.e. a PeakCatalog) and to turn it into a list/tree of children with flags plus a mixing matrix. I think we can alter the currently half-useful Deblend class to incorporate the necessary functions (). I would however rename the class to Blend because it is the model of the blended scene. And that should be returned from a call to deblend.

Lastly, we already have a class Step_AS that acts in a similar way by storing A and S to have quicker estimates of the Lipschitz constant. We should fold that class in with Blend.

from scarlet.

Another piece of logic for item 2), and how we deal with the boxes:

• We treat each component separately: [S1,...,Sk]
• Put it into the full matrix S
• Compute the residual R=(Y-APS)
• Compute updated S' = S - 1/L_S A^T(Y-APS)
• Slice out the Sk from S' and report those individually as result of prox_f(Xk), which would then go into the prox_g.

The question is: which of the elements do we need to store/update for the A update: A - 1/L_a (Y-APS) (PS)^T? Note that the convolved sources appear in the last term. Formally, all occurrence of S needs to be replace with S', so we'd have to compute the model twice per iteration. I think that the differences in gradient direction (with and without update of S) will be rather modest, at least after a few iterations. Given that the convolution operation to build the entire model is expensive, this suggest that we keep either PS or R fixed for the updates of S and A.

The simplest thing would this be to use the old S matrix for the A update. This is a deviation from the coordinate descent method, but we'd only have to build the model APS once.

from scarlet.