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

DM-10096 Add unit test asserts for SpherePoint, SpherePointList and PointList #15

Merged
merged 3 commits into from
Apr 6, 2017
Merged
Show file tree
Hide file tree
Changes from all commits
Commits
File filter

Filter by extension

Filter by extension

Conversations
Failed to load comments.
Jump to
Jump to file
Failed to load files.
Diff view
Diff view
1 change: 1 addition & 0 deletions .gitignore
Original file line number Diff line number Diff line change
@@ -1,3 +1,4 @@
_build.*
.sconsign.dblite
config.log
.sconf_temp
Expand Down
37 changes: 22 additions & 15 deletions python/lsst/shapelet/tests.py
Original file line number Diff line number Diff line change
@@ -1,6 +1,6 @@
#
# LSST Data Management System
# Copyright 2008-2013 LSST Corporation.
# Copyright 2008-2017 LSST Corporation.
#
# This product includes software developed by the
# LSST Project (http://www.lsst.org/).
Expand Down Expand Up @@ -87,9 +87,11 @@ def compareShapeletFunctions(self, a, b, rtolEllipse=1E-13, rtolCoeff=1E-13,
atolEllipse=1E-14, atolCoeff=1E-14):
self.assertEqual(a.getOrder(), b.getOrder())
self.assertEqual(a.getBasisType(), b.getBasisType())
self.assertClose(a.getEllipse().getParameterVector(), b.getEllipse().getParameterVector(),
rtol=rtolEllipse, atol=atolEllipse)
self.assertClose(a.getCoefficients(), b.getCoefficients(), rtol=rtolCoeff, atol=atolCoeff)
self.assertFloatsAlmostEqual(a.getEllipse().getParameterVector(),
b.getEllipse().getParameterVector(),
rtol=rtolEllipse, atol=atolEllipse)
self.assertFloatsAlmostEqual(a.getCoefficients(), b.getCoefficients(),
rtol=rtolCoeff, atol=atolCoeff)

def simplifyMultiShapeletFunction(self, msf):
keep = []
Expand Down Expand Up @@ -126,29 +128,34 @@ def checkMoments(self, function, x, y, z):
)
imageMoments = lsst.afw.geom.ellipses.Ellipse(quadrupole, dipole)
shapeletMoments = function.evaluate().computeMoments()
self.assertClose(imageMoments.getCenter().getX(), shapeletMoments.getCenter().getX(), rtol=1E-3)
self.assertClose(imageMoments.getCenter().getY(), shapeletMoments.getCenter().getY(), rtol=1E-3)
self.assertClose(imageMoments.getCore().getIxx(), shapeletMoments.getCore().getIxx(), rtol=1E-3)
self.assertClose(imageMoments.getCore().getIyy(), shapeletMoments.getCore().getIyy(), rtol=1E-3)
self.assertClose(imageMoments.getCore().getIxy(), shapeletMoments.getCore().getIxy(), rtol=1E-3)
self.assertFloatsAlmostEqual(imageMoments.getCenter().getX(),
shapeletMoments.getCenter().getX(), rtol=1E-3)
self.assertFloatsAlmostEqual(imageMoments.getCenter().getY(),
shapeletMoments.getCenter().getY(), rtol=1E-3)
self.assertFloatsAlmostEqual(imageMoments.getCore().getIxx(),
shapeletMoments.getCore().getIxx(), rtol=1E-3)
self.assertFloatsAlmostEqual(imageMoments.getCore().getIyy(),
shapeletMoments.getCore().getIyy(), rtol=1E-3)
self.assertFloatsAlmostEqual(imageMoments.getCore().getIxy(),
shapeletMoments.getCore().getIxy(), rtol=1E-3)
integral = numpy.trapz(numpy.trapz(z, gx, axis=1), y, axis=0)
self.assertClose(integral, function.evaluate().integrate(), rtol=1E-3)
self.assertFloatsAlmostEqual(integral, function.evaluate().integrate(), rtol=1E-3)

def checkConvolution(self, f1, f2):
bbox = lsst.afw.geom.Box2I(lsst.afw.geom.Point2I(-50, -50), lsst.afw.geom.Point2I(50, 50))
i1 = lsst.afw.image.ImageD(bbox)
f1.evaluate().addToImage(i1)
self.assertClose(i1.getArray().sum(), f1.evaluate().integrate(), rtol=1E-3)
self.assertFloatsAlmostEqual(i1.getArray().sum(), f1.evaluate().integrate(), rtol=1E-3)
i2 = lsst.afw.image.ImageD(bbox)
f2.evaluate().addToImage(i2)
self.assertClose(i2.getArray().sum(), f2.evaluate().integrate(), rtol=1E-3)
self.assertFloatsAlmostEqual(i2.getArray().sum(), f2.evaluate().integrate(), rtol=1E-3)
fc1 = f1.convolve(f2)
fc2 = f2.convolve(f1)
ic1 = lsst.afw.image.ImageD(bbox)
fc1.evaluate().addToImage(ic1)
ic2 = lsst.afw.image.ImageD(bbox)
fc2.evaluate().addToImage(ic2)
self.assertClose(ic1.getArray(), ic2.getArray())
self.assertFloatsAlmostEqual(ic1.getArray(), ic2.getArray())
out = lsst.afw.image.ImageD(bbox)
if scipy is None:
print("Skipping convolution test; scipy could not be imported.")
Expand All @@ -159,6 +166,6 @@ def checkConvolution(self, f1, f2):
# I can't even make the operation commutative, let alone correct.
scipy.ndimage.convolve(i1.getArray(), i2.getArray(), output=out.getArray(),
mode="constant", cval=0.0)
self.assertClose(out.getArray(), ic1.getArray(), rtol=1E-4, atol=1E-5)
self.assertClose(out.getArray(), ic2.getArray(), rtol=1E-4, atol=1E-5)
self.assertFloatsAlmostEqual(out.getArray(), ic1.getArray(), rtol=1E-4, atol=1E-5)
self.assertFloatsAlmostEqual(out.getArray(), ic2.getArray(), rtol=1E-4, atol=1E-5)
return fc1, fc2
47 changes: 25 additions & 22 deletions tests/profiles.py
Original file line number Diff line number Diff line change
@@ -1,8 +1,6 @@
#!/usr/bin/env python
Copy link
Contributor

Choose a reason for hiding this comment

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

some of the following formatting changes should probably be on their own commit, consider changing


#
# LSST Data Management System
# Copyright 2008-2014 LSST Corporation.
# Copyright 2008-2017 LSST Corporation.
#
# This product includes software developed by the
# LSST Project (http://www.lsst.org/).
Expand All @@ -21,11 +19,12 @@
# the GNU General Public License along with this program. If not,
# see <http://www.lsstcorp.org/LegalNotices/>.
#

from __future__ import absolute_import, division, print_function
import unittest
import numpy
import os

import numpy as np

import lsst.utils.tests
import lsst.afw.geom
import lsst.afw.geom.ellipses as el
Expand All @@ -52,8 +51,8 @@ class ProfileTestCase(lsst.shapelet.tests.ShapeletTestCase):
def testRadii(self):
"""Check RadialProfile definitions of moments and half-light radii.
"""
s = numpy.linspace(-20.0, 20.0, 1000)
x, y = numpy.meshgrid(s, s)
s = np.linspace(-20.0, 20.0, 1000)
x, y = np.meshgrid(s, s)
r = (x**2 + y**2)**0.5
dxdy = (s[1] - s[0])**2
for name in ["gaussian", "exp", "ser2", "luv", "lux"]:
Expand All @@ -62,19 +61,20 @@ def testRadii(self):
# lux and luv don't use the true half-light radius; instead they use the half-light radius
# of the exp and dev profiles they approximate
if not name.startswith("lu"):
self.assertClose(z[r < 1].sum(), 0.5*z.sum(), rtol=0.01)
self.assertFloatsAlmostEqual(z[r < 1].sum(), 0.5*z.sum(), rtol=0.01)
# lhs of this comparison is the moments radius (using a sum approximation to the integral)
self.assertClose(((z*x**2).sum() / z.sum())**0.5, profile.getMomentsRadiusFactor(), rtol=0.01)
self.assertFloatsAlmostEqual(((z*x**2).sum() / z.sum())**0.5,
profile.getMomentsRadiusFactor(), rtol=0.01)

def testGaussian(self):
"""Test that the Gaussian profile's shapelet 'approximation' is actually exact.
"""
profile = lsst.shapelet.RadialProfile.get("gaussian")
r = numpy.linspace(0.0, 4.0, 100)
r = np.linspace(0.0, 4.0, 100)
z1 = profile.evaluate(r)
basis = profile.getBasis(1)
z2 = lsst.shapelet.tractor.evaluateRadial(basis, r, sbNormalize=True)[0, :]
self.assertClose(z1, z2, rtol=1E-8)
self.assertFloatsAlmostEqual(z1, z2, rtol=1E-8)

def testShapeletApproximations(self):
psf0 = lsst.shapelet.ShapeletFunction(0, lsst.shapelet.HERMITE, PSF_SIGMA)
Expand All @@ -88,22 +88,22 @@ def testShapeletApproximations(self):
check1 = lsst.afw.image.ImageD(os.path.join("tests", "data", name + "-approx.fits")).getArray()
xc = check1.shape[1] // 2
yc = check1.shape[0] // 2
xb = numpy.arange(check1.shape[1], dtype=float) - xc
yb = numpy.arange(check1.shape[0], dtype=float) - yc
xg, yg = numpy.meshgrid(xb, yb)
xb = np.arange(check1.shape[1], dtype=float) - xc
yb = np.arange(check1.shape[0], dtype=float) - yc
xg, yg = np.meshgrid(xb, yb)

basis = lsst.shapelet.RadialProfile.get(name).getBasis(nComponents, maxRadius)
builder = lsst.shapelet.MatrixBuilderD.Factory(xg.ravel(), yg.ravel(), basis, psf)()
image1 = numpy.zeros(check1.shape, dtype=float)
image1 = np.zeros(check1.shape, dtype=float)
matrix = image1.reshape(check1.size, 1)
builder(matrix, el.Ellipse(ellipse))
self.assertClose(check1, image1, plotOnFailure=False, rtol=5E-5, relTo=check1.max())
self.assertFloatsAlmostEqual(check1, image1, plotOnFailure=False, rtol=5E-5, relTo=check1.max())
msf = basis.makeFunction(el.Ellipse(ellipse, lsst.afw.geom.Point2D(xc, yc)),
numpy.array([1.0], dtype=float))
np.array([1.0], dtype=float))
msf = msf.convolve(psf)
image2 = numpy.zeros(check1.shape, dtype=float)
image2 = np.zeros(check1.shape, dtype=float)
msf.evaluate().addToImage(lsst.afw.image.ImageD(image2, False))
self.assertClose(check1, image2, plotOnFailure=False, rtol=5E-5, relTo=check1.max())
self.assertFloatsAlmostEqual(check1, image2, plotOnFailure=False, rtol=5E-5, relTo=check1.max())

if name == 'exp':
# check2 is the exact profile, again by GalSim.
Expand All @@ -113,7 +113,8 @@ def testShapeletApproximations(self):
check2 = lsst.afw.image.ImageD(
os.path.join("tests", "data", name + "-exact.fits")
).getArray()
self.assertClose(check2, image1, plotOnFailure=False, rtol=1E-3, relTo=check2.max())
self.assertFloatsAlmostEqual(check2, image1, plotOnFailure=False,
rtol=1E-3, relTo=check2.max())

if CHECK_COMPONENT_IMAGES:
# This was once useful for debugging test failures, and may be again, but it's
Expand All @@ -123,9 +124,10 @@ def testShapeletApproximations(self):
check = lsst.afw.image.ImageD(
os.path.join("tests", "data", "%s-approx-%0d.fits" % (name, n))
).getArray()
image = numpy.zeros(check1.shape, dtype=float)
image = np.zeros(check1.shape, dtype=float)
sf.evaluate().addToImage(lsst.afw.image.ImageD(image, False))
self.assertClose(check, image, plotOnFailure=False, rtol=5E-5, relTo=check1.max())
self.assertFloatsAlmostEqual(check, image, plotOnFailure=False,
rtol=5E-5, relTo=check1.max())


def suite():
Expand All @@ -143,5 +145,6 @@ def run(shouldExit=False):
"""Run the tests"""
lsst.utils.tests.run(suite(), shouldExit)


if __name__ == "__main__":
run(True)
14 changes: 7 additions & 7 deletions tests/testFunctorKeys.py
Original file line number Diff line number Diff line change
@@ -1,9 +1,6 @@
from builtins import range
#!/usr/bin/env python

#
# LSST Data Management System
# Copyright 2008-2014 LSST Corporation.
# Copyright 2008-2017 LSST Corporation.
#
# This product includes software developed by the
# LSST Project (http://www.lsst.org/).
Expand All @@ -22,17 +19,19 @@
# the GNU General Public License along with this program. If not,
# see <http://www.lsstcorp.org/LegalNotices/>.
#

from __future__ import absolute_import, division, print_function
from builtins import range
import unittest
import numpy

import numpy as np

import lsst.utils.tests
import lsst.afw.geom.ellipses
import lsst.shapelet.tests
import lsst.afw.image
import lsst.afw.table

numpy.random.seed(500)
np.random.seed(500)


class FunctorKeyTestCase(lsst.shapelet.tests.ShapeletTestCase):
Expand Down Expand Up @@ -104,5 +103,6 @@ def run(shouldExit=False):
"""Run the tests"""
lsst.utils.tests.run(suite(), shouldExit)


if __name__ == "__main__":
run(True)
22 changes: 11 additions & 11 deletions tests/testHermiteTransformMatrix.py
Original file line number Diff line number Diff line change
@@ -1,8 +1,6 @@
#!/usr/bin/env python

#
# LSST Data Management System
# Copyright 2008-2013 LSST Corporation.
# Copyright 2008-2017 LSST Corporation.
#
# This product includes software developed by the
# LSST Project (http://www.lsst.org/).
Expand All @@ -25,7 +23,8 @@
from __future__ import print_function
from builtins import range
import unittest
import numpy

import numpy as np

try:
import scipy.special
Expand All @@ -36,7 +35,7 @@
import lsst.afw.geom
import lsst.shapelet.tests

numpy.random.seed(500)
np.random.seed(500)


class HermiteTransformMatrixTestCase(lsst.shapelet.tests.ShapeletTestCase):
Expand All @@ -50,13 +49,13 @@ def setUp(self):
def ht(n):
"""return a scipy poly1d for the nth 'alternate' Hermite polynomial (i.e. Hermite polynomial
with shapeley normalization)"""
return (scipy.poly1d([(2**n * numpy.pi**0.5 * scipy.special.gamma(n+1))**(-0.5)])
* scipy.special.hermite(n))
return (scipy.poly1d([(2**n * np.pi**0.5 * scipy.special.gamma(n+1))**(-0.5)]) *
scipy.special.hermite(n))

def testCoefficientMatrices(self):
coeff = self.htm.getCoefficientMatrix()
coeffInv = self.htm.getInverseCoefficientMatrix()
self.assertClose(numpy.identity(self.order+1), numpy.dot(coeff, coeffInv))
self.assertFloatsAlmostEqual(np.identity(self.order+1), np.dot(coeff, coeffInv))
# Both matrices should be lower-triangular
for i in range(0, self.order+1):
for j in range(i+1, self.order+1):
Expand All @@ -68,7 +67,7 @@ def testCoefficientMatrices(self):
return
for n in range(0, self.order+1):
poly = self.ht(n)
self.assertClose(coeff[n, :n+1], poly.c[::-1], atol=1E-15)
self.assertFloatsAlmostEqual(coeff[n, :n+1], poly.c[::-1], atol=1E-15)

def testTransformMatrix(self):
if scipy is None:
Expand All @@ -78,7 +77,7 @@ def testTransformMatrix(self):
s = lsst.afw.geom.LinearTransform.makeScaling(2.0, 1.5)
r = lsst.afw.geom.LinearTransform.makeRotation(0.30*lsst.afw.geom.radians)
transforms = [s, r, s*r*s]
testPoints = numpy.random.randn(10, 2)
testPoints = np.random.randn(10, 2)
for transform in transforms:
m = self.htm.compute(transform)
for testPoint in testPoints:
Expand All @@ -90,7 +89,7 @@ def testTransformMatrix(self):
v2 = 0.0
for j, jnx, jny in lsst.shapelet.HermiteIndexGenerator(self.order):
v2 += m[i, j] * self.ht(jnx)(origPoint.getX()) * self.ht(jny)(origPoint.getY())
self.assertClose(v1, v2, rtol=1E-11)
self.assertFloatsAlmostEqual(v1, v2, rtol=1E-11)


def suite():
Expand All @@ -108,5 +107,6 @@ def run(shouldExit=False):
"""Run the tests"""
lsst.utils.tests.run(suite(), shouldExit)


if __name__ == "__main__":
run(True)