Skip to content
New issue

Have a question about this project? Sign up for a free GitHub account to open an issue and contact its maintainers and the community.

By clicking “Sign up for GitHub”, you agree to our terms of service and privacy statement. We’ll occasionally send you account related emails.

Already on GitHub? Sign in to your account

tickets/DM-15638: Fix DcrModel subfilter order. #101

Merged
merged 2 commits into from Oct 11, 2018
Merged

Conversation

isullivan
Copy link
Contributor

A few lines of code were changed, and several unit tests were added.

The order of the DrModel subfilters was reversed, which was due to the calculated shift due to DCR having the wrong sign. There were several bugs that worked together in this: the rotation calculation was written for a flipped WCS, but if a non-flipped WCS was supplied it would be off by -90 degrees. Then, afwGeom.Extent2D() takes arguments in the order X,Y but had been supplied in the order Y,X, so that had been off by another 90 degrees. Finally, DCR is mostly symmetric around the center of the band, so everything still mostly worked even though the shift in each band had the wrong sign.

Copy link
Contributor

@parejkoj parejkoj left a comment

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

See comments: the tests should be reworked so they check more than one value (shades of xkcd 221).

On the other hand, yay more tests!

if wcs.isFlipped:
cdAngle = (np.arctan2(-cd[0, 1], cd[0, 0]) + np.arctan2(cd[1, 0], cd[1, 1]))/2.
else:
cdAngle = (np.arctan2(cd[0, 1], -cd[0, 0]) + np.arctan2(cd[1, 0], cd[1, 1]))/2.
Copy link
Contributor

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

It took me a bit to distinguish the difference between these two lines, but I don't know of a cleaner way to do it.

Separately, do you have a ticket to lift this code out of here and into somewhere more generic (afw, maybe)? Other places would like to use it.

Copy link
Contributor Author

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

I made a ticket, it is DM-16130.

"""Test that the bluest subfilter always has the largest amplitude.
"""
seed = 6
rng = np.random.RandomState(seed)
Copy link
Contributor

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

Why these different random seeds? Why not just put one np.random.seed(5) in setUp and use np.random.X in the methods?

@@ -134,14 +134,16 @@ def makeDummyWcs(self, rotAngle, pixelScale, crval):
Pixel scale of the projection.
crval : `lsst.afw.geom.SpherePoint`
Coordinates of the reference pixel of the wcs.
flipX : `bool`, optional
Optionally flip the direction of increasing Right Ascension.
Copy link
Contributor

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

Don't say "optionally flip": default True implies it defaults to being flipped, so it's not optional.

afwGeom.Extent2D(0.3886592703, 0.2248919247)]
refShift = [afwGeom.Extent2D(-0.3103517169, -0.5363512808),
afwGeom.Extent2D(0.001092054612, 0.001887293861),
afwGeom.Extent2D(0.2248919247, 0.3886592703)]
Copy link
Contributor

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

Where did these numbers come from? A comment about what they are would be good.

Copy link
Contributor Author

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

There's a comment at the start of the test, but I'll add another here.

pixelScale = 0.2*arcseconds
rotAngle = 360.*rng.rand()*degrees
azimuth = 360.*rng.rand()*degrees
elevation = (45. + rng.rand()*40.)*degrees # Restrict to 45 < elevation < 85 degrees
Copy link
Contributor

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

Given that you set a random seed, this will always have the same value. Maybe lift the internals into a checkDcrSubFilterOrder method, and have a loop over some carefully chosen values (including zenith, north, south, east, west) in the testDcrSub... method that call check...?

if dx < maskOffsetEnd:
refImage.mask.array[:, dx - maskOffsetEnd:] = maskValue
if dy < maskOffsetEnd:
refImage.mask.array[dy - maskOffsetEnd:, :] = maskValue
Copy link
Contributor

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

The above logic, with the dx/dy, and maskOffsetStart/End values defined just above it, suggests to me that you want to test a few different cases. Again, a checkDcrShift that this calls with a few different dx/dy and offsets seems like a good idea.

rng = np.random.RandomState(seed)
refAngle = 0.*degrees
# Any additional arbitrary rotation should fall out of the calculation
cdRotAngle = 360*rng.rand()*degrees
Copy link
Contributor

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

This and the comment above definitely tell me there should be a loop over a range of cdRotAngles.

def testRotationSouthZero(self):
"""Test that an observation pointed due South has zero rotation angle.

An observation pointed due South should have zenith directly to the
Copy link
Contributor

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

Your test isn't of "due south", but rather "-30 deg, along the meridian".

seed = 8
rng = np.random.RandomState(seed)
# Any additional arbitrary rotation should fall out of the calculation
cdRotAngle = 360*rng.rand()*degrees
Copy link
Contributor

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

As above, test a range of rotAngles.

self.assertAnglesAlmostEqual(-cdRotAngle, rotAngle, maxDiff=1e-6*radians)

def testRotationNotFlipped(self):
"""Check the interpretation of rotations in the WCS.
Copy link
Contributor

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

This and testRotationFlipped could be factored out into checkRotation that takes a flipX parameter and tests against the negative cdRotAngle if True. Then you can check a range of angles both flipped and not using the same code block.

@parejkoj
Copy link
Contributor

Thanks for updating the tests. I still think testApplyDcr should be split into a few cases, since it has a couple of if statements that are always the same because dx and maxOffset are always the same.

Sign up for free to join this conversation on GitHub. Already have an account? Sign in to comment
Labels
None yet
Projects
None yet
Development

Successfully merging this pull request may close these issues.

None yet

2 participants