pmelchior / scarlet

hyperspectral galaxy modeling and deblending
MIT License
49 stars 22 forks source link

ResolutionRenderer fails for rotated images #276

Closed charlotteaward closed 1 year ago

charlotteaward commented 1 year ago

When creating a Frame from observations which are rotated relative to one another the ResolutionRenderer produces the following error when generating the unrotated coordinates:

Frame.from_observations(observations, coverage='intersection', model_psf=model_psf) File "/mnt/c/Users/davro/scarlet/scarlet/scarlet/frame.py", line 285, in from_observations obs.match(model_frame) File "/mnt/c/Users/davro/scarlet/scarlet/scarlet/observation.py", line 109, in match self.renderer = ResolutionRenderer(self, model_frame) File "/mnt/c/Users/davro/scarlet/scarlet/scarlet/renderer.py", line 322, in init Y_unrot = ( ValueError: cannot reshape array of size 100 into shape (100,100)

Changing lines 322-333 in renderer.py to this instead:

       Y_unrot = (
            (coord_hr[:, 0] - center_y) * self.angle[0]
            - (coord_hr[:, 1] - center_x) * self.angle[1]
        ).reshape(lr_shape[0])
        X_unrot = (
            (coord_hr[:, 1] - center_x) * self.angle[0]
            + (coord_hr[:, 0] - center_y) * self.angle[1]
        ).reshape(lr_shape[1])

        self.Y_unrot = Y_unrot
        self.X_unrot = X_unrot

allows the code to run without errors, but models do not have the correct angles when rendered (among other rendering issues), e.g:

scene_ZTF18abxxohm_2_i0

scene_ZTF18abxxohm_1_r0

Unless the original error is arising from something to do with the input images themselves (perhaps related to their wcs format or large 1" pixel scale), I think a fix may be needed. Since coord_hr usually has size (N, 2) and lr_shape is usually (NxN), I'm not sure how the reshape lines are supposed to work without producing the ValueError - any ideas greatly appreciated!

pmelchior commented 1 year ago

The first problem is clearly a bug, and it comes from the fact that Y_rot is a list of coordinates, and not an image. I think that the reshape operation is completely redundant.

For the second one, it looks like the render method rotates the image in the opposite direction. I suspect that something is wrong about the relative angle between the two frames.

charlotteaward commented 1 year ago

I have been able to fix the problem by making additional changes in renderer.py. Lines 322-348 should be the following for interpolation.get_angles to work as intended:

            Y_unrot = (
                (coord_hr[:, 0] - center_y) * self.angle[0]
                - (coord_hr[:, 1] - center_x) * self.angle[1]
            ).reshape(lr_shape[0])
            X_unrot = (
                (coord_hr[:, 1] - center_x) * self.angle[0]
                + (coord_hr[:, 0] - center_y) * self.angle[1]
            ).reshape(lr_shape[1])

            # Removing redundancy
            self.Y_unrot = Y_unrot
            self.X_unrot = X_unrot

            if self.small_axis:
                self.shifts = np.array(
                    [self.Y_unrot * self.angle[0], -self.Y_unrot * self.angle[1]]
                )
                self.other_shifts = np.array(
                    [self.angle[1] * self.X_unrot, self.angle[0] * self.X_unrot,]
                )
            else:
                self.shifts = np.array(
                    [self.angle[1] * self.X_unrot, self.angle[0] * self.X_unrot,]
                )
                self.other_shifts = np.array(
                    [self.Y_unrot * self.angle[0], -self.Y_unrot * self.angle[1],]
                )

The scenes become: scene_ZTF18abxxohm_2_i0

scene_ZTF18abxxohm_0_g0

There are some additional artifacts in the residual related to the PSF sampling and a saturated star but I think the angle issue is fixed.

pmelchior commented 1 year ago

OK, great. Can you put your changes into a PR, pls?

charlotteaward commented 1 year ago

Yep will do! I noticed one other thing: renderer.py also fails here: https://github.com/pmelchior/scarlet/blob/11d1363a30de893e3a635f2aae48161602afa23a/scarlet/renderer.py#L274

File "/home/cwar4677/anaconda3/lib/python3.9/site-packages/autograd/numpy/numpy_wrapper.py", line 94, in stack raise ValueError('all input arrays must have the same shape') ValueError: all input arrays must have the same shape

Do we expect the resolution renderer to work for rectangular images? If so, there may be an easy fix by padding to make images square.

pmelchior commented 1 year ago

Good catch! Don't pad, this needs to work with non-square images. I'll check.

pmelchior commented 1 year ago

On second look (and reviewing eq. 20 in Remy's paper), the code seem to require that the observations have square frames. The model frame doesn't need to be, I think.

This is unfortunately an undocumented requirement for ResolutionRenderer, which we should make explicit.