Skip to content

Commit

Permalink
Merge 5c93b99 into e9afa8d
Browse files Browse the repository at this point in the history
  • Loading branch information
larrybradley committed May 1, 2015
2 parents e9afa8d + 5c93b99 commit 6ee4303
Show file tree
Hide file tree
Showing 6 changed files with 164 additions and 37 deletions.
2 changes: 2 additions & 0 deletions CHANGES.rst
Expand Up @@ -14,6 +14,8 @@ New Features

- ``astropy.convolution``

- Fix issue with repeated normalizations of ``Kernels``. [#3747]

- ``astropy.coordinates``

- ``astropy.cosmology``
Expand Down
74 changes: 48 additions & 26 deletions astropy/convolution/core.py
Expand Up @@ -45,10 +45,7 @@ class Kernel(object):

def __init__(self, array):
self._array = array
if self._array.sum() == 0:
self._normalization = np.inf
else:
self._normalization = 1. / self._array.sum()
self.calculate_normalization()

@property
def truncation(self):
Expand Down Expand Up @@ -88,42 +85,66 @@ def center(self):
"""
return [axes_size // 2 for axes_size in self._array.shape]

def calculate_normalization(self, mode='integral'):
"""
Calculate the kernel normalization factor.
Parameters
----------
mode : {'integral', 'peak'}
One of the following modes:
* 'integral' (default)
Kernel is normalized such that its integral = 1.
* 'peak'
Kernel is normalized such that its peak = 1.
"""

if mode == 'integral':
if self._array.sum() == 0:
self._normalization = np.inf
else:
self._normalization = 1. / self._array.sum()
elif mode == 'peak':
self._normalization = 1. / self._array.max()
else:
raise ValueError("invalid mode, must be 'integral' or 'peak'")
return self._normalization

@property
def normalization(self):
"""
Kernel normalization factor
The kernel normalization factor.
"""
return self._normalization

def normalize(self, mode='integral'):
"""
Force normalization of filter kernel.
Normalize the filter kernel.
Parameters
----------
mode : {'integral', 'peak'}
One of the following modes:
* 'integral' (default)
Kernel normalized such that its integral = 1.
Kernel is normalized such that its integral = 1.
* 'peak'
Kernel normalized such that its peak = 1.
Kernel is normalized such that its peak = 1.
"""
# There are kernel that sum to zero and
# the user should be warned in this case

np.multiply(self._array, self.calculate_normalization(mode=mode),
self.array)

# Warn the user for kernels that sum to zero
if np.isinf(self._normalization):
warnings.warn('Kernel cannot be normalized because the '
'normalization factor is infinite.',
AstropyUserWarning)
return

if np.abs(self._normalization) > MAX_NORMALIZATION:
warnings.warn("Normalization factor of kernel is "
"exceptionally large > {0}.".format(MAX_NORMALIZATION),
warnings.warn('The kernel normalization factor is exceptionally '
'large, > {0}.'.format(MAX_NORMALIZATION),
AstropyUserWarning)
if mode == 'integral':
self._array *= self._normalization
if mode == 'peak':
np.divide(self._array, self._array.max(), self.array)
self._normalization = 1. / self._array.sum()

@property
def shape(self):
Expand Down Expand Up @@ -228,6 +249,8 @@ class Kernel1D(Kernel):
def __init__(self, model=None, x_size=None, array=None, **kwargs):
# Initialize from model
if array is None:
if self._model is None:
raise TypeError("Must specify either array or model.")

if x_size is None:
x_size = self._default_size
Expand All @@ -246,8 +269,7 @@ def __init__(self, model=None, x_size=None, array=None, **kwargs):
# Initialize from array
elif array is not None:
self._model = None
else:
raise TypeError("Must specify either array or model.")

super(Kernel1D, self).__init__(array)


Expand Down Expand Up @@ -288,6 +310,8 @@ def __init__(self, model=None, x_size=None, y_size=None, array=None, **kwargs):

# Initialize from model
if array is None:
if self._model is None:
raise TypeError("Must specify either array or model.")

if x_size is None:
x_size = self._default_size
Expand All @@ -296,7 +320,7 @@ def __init__(self, model=None, x_size=None, y_size=None, array=None, **kwargs):

if y_size is None:
y_size = x_size
elif x_size != int(x_size):
elif y_size != int(y_size):
raise TypeError("y_size should be an integer")

# Set ranges where to evaluate the model
Expand All @@ -316,8 +340,6 @@ def __init__(self, model=None, x_size=None, y_size=None, array=None, **kwargs):
# Initialize from array
elif array is not None:
self._model = None
else:
raise TypeError("Must specify either array or model.")

super(Kernel2D, self).__init__(array)

Expand Down Expand Up @@ -348,8 +370,8 @@ def kernel_arithmetics(kernel, value, operation):
if operation == "sub":
new_array = add_kernel_arrays_1D(kernel.array, -value.array)
if operation == "mul":
raise Exception("Kernel operation not supported. Maybe you want to"
"use convolve(kernel1, kernel2) instead.")
raise Exception("Kernel operation not supported. Maybe you want "
"to use convolve(kernel1, kernel2) instead.")
new_kernel = Kernel1D(array=new_array)
new_kernel._separable = kernel._separable and value._separable
new_kernel._is_bool = kernel._is_bool or value._is_bool
Expand All @@ -361,8 +383,8 @@ def kernel_arithmetics(kernel, value, operation):
if operation == "sub":
new_array = add_kernel_arrays_2D(kernel.array, -value.array)
if operation == "mul":
raise Exception("Kernel operation not supported. Maybe you want to"
"use convolve(kernel1, kernel2) instead.")
raise Exception("Kernel operation not supported. Maybe you want "
"to use convolve(kernel1, kernel2) instead.")
new_kernel = Kernel2D(array=new_array)
new_kernel._separable = kernel._separable and value._separable
new_kernel._is_bool = kernel._is_bool or value._is_bool
Expand Down
8 changes: 3 additions & 5 deletions astropy/convolution/tests/test_convolve.py
Expand Up @@ -25,12 +25,10 @@ def test_list(self):
"""

x = [1, 4, 5, 6, 5, 7, 8]

y = [0.2, 0.6, 0.2]

z = convolve([1, 4, 5, 6, 5, 7, 8], [0.2, 0.6, 0.2], boundary=None)

assert_array_almost_equal_nulp(z, np.array([ 0. , 3.6, 5. , 5.6, 5.6, 6.8, 0. ]), 10)
z = convolve(x, y, boundary=None)
assert_array_almost_equal_nulp(z,
np.array([0., 3.6, 5., 5.6, 5.6, 6.8, 0.]), 10)

@pytest.mark.parametrize(('dtype_array', 'dtype_kernel'), VALID_DTYPES)
def test_dtype(self, dtype_array, dtype_kernel):
Expand Down
2 changes: 1 addition & 1 deletion astropy/convolution/tests/test_convolve_fft.py
Expand Up @@ -462,7 +462,7 @@ def test_uniform_3x3_withnan(self, boundary, interpolate_nan,
def test_big_fail(self):
""" Test that convolve_fft raises an exception if a too-large array is passed in """

with pytest.raises((ValueError, MemoryError)) as ex:
with pytest.raises((ValueError, MemoryError)):
# while a good idea, this approach did not work; it actually writes to disk
#arr = np.memmap('file.np', mode='w+', shape=(512, 512, 512), dtype=np.complex)
# this just allocates the memory but never touches it; it's better:
Expand Down
3 changes: 1 addition & 2 deletions astropy/convolution/tests/test_discretize.py
Expand Up @@ -11,8 +11,7 @@

from ..utils import discretize_model
from ...modeling.functional_models import (
Gaussian1D, Box1D, MexicanHat1D, Trapezoid1D,
Gaussian2D, Box2D, MexicanHat2D)
Gaussian1D, Box1D, MexicanHat1D, Gaussian2D, Box2D, MexicanHat2D)
from ...modeling.tests.example_models import models_1D, models_2D
from ...modeling.tests.test_models import create_model

Expand Down
112 changes: 109 additions & 3 deletions astropy/convolution/tests/test_kernel_class.py
Expand Up @@ -5,18 +5,19 @@
import itertools

import numpy as np
from numpy.testing import assert_almost_equal
from numpy.testing import assert_almost_equal, assert_allclose

from ...tests.helper import pytest
from ..convolve import convolve, convolve_fft
from ..kernels import (
Gaussian1DKernel, Gaussian2DKernel, Box1DKernel, Box2DKernel,
Trapezoid1DKernel, TrapezoidDisk2DKernel, MexicanHat1DKernel,
Tophat2DKernel, MexicanHat2DKernel, AiryDisk2DKernel, Ring2DKernel,
CustomKernel, Model1DKernel, Model2DKernel)
CustomKernel, Model1DKernel, Model2DKernel, Kernel1D, Kernel2D)

from ..utils import KernelSizeError
from ...modeling.models import Box2D, Gaussian1D, Gaussian2D
from ...utils.exceptions import AstropyWarning, AstropyUserWarning

try:
from scipy.ndimage import filters
Expand Down Expand Up @@ -211,6 +212,31 @@ def test_rmultiply_scalar_type(self, number):
gauss_new = gauss * number
assert type(gauss_new) == Gaussian1DKernel

def test_multiply_kernel1d(self):
"""Test that multipling two 1D kernels raises an exception."""
gauss = Gaussian1DKernel(3)
with pytest.raises(Exception):
gauss * gauss

def test_multiply_kernel2d(self):
"""Test that multipling two 2D kernels raises an exception."""
gauss = Gaussian2DKernel(3)
with pytest.raises(Exception):
gauss * gauss

def test_multiply_kernel1d_kernel2d(self):
"""
Test that multipling a 1D kernel with a 2D kernel raises an
exception.
"""
with pytest.raises(Exception):
Gaussian1DKernel(3) * Gaussian2DKernel(3)

def test_add_kernel_scalar(self):
"""Test that adding a scalar to a kernel raises an exception."""
with pytest.raises(Exception):
Gaussian1DKernel(3) + 1

def test_model_1D_kernel(self):
"""
Check Model1DKernel against Gaussian1Dkernel
Expand Down Expand Up @@ -302,7 +328,7 @@ def test_custom_kernel_odd_error(self):
Check if CustomKernel raises if the array size is odd.
"""
with pytest.raises(KernelSizeError):
custom = CustomKernel([1, 1, 1, 1])
CustomKernel([1, 1, 1, 1])

def test_add_1D_kernels(self):
"""
Expand Down Expand Up @@ -421,4 +447,84 @@ def test_box_kernels_even_size(self, width):
assert np.all([_ % 2 != 0 for _ in kernel_2D.shape])
assert kernel_2D.array.sum() == 1.

def test_kernel_normalization(self):
"""
Test that repeated normalizations do not change the kernel [#3747].
"""

kernel = CustomKernel(np.ones(5))
kernel.normalize()
data = np.copy(kernel.array)

kernel.normalize()
assert_allclose(data, kernel.array)

kernel.normalize()
assert_allclose(data, kernel.array)

def test_kernel_normalization_mode(self):
"""
Test that an error is raised if mode is invalid.
"""
with pytest.raises(ValueError):
kernel = CustomKernel(np.ones(3))
kernel.normalize(mode='invalid')

def test_kernel_normalization_inf(self, recwarn):
"""
Test that a warning is issued when normalizing a zero-sum kernel.
"""
kernel = CustomKernel(np.array([-1, 0, 1]))
kernel.normalize()
w = recwarn.pop(AstropyUserWarning)
assert issubclass(w.category, AstropyWarning)

def test_kernel_normalization_large(self, recwarn):
"""
Test that a warning is issued when the inverse normalization factor
is large.
"""
kernel = CustomKernel(np.ones(3) * 1.e-3)
kernel.normalize()
w = recwarn.pop(AstropyUserWarning)
assert issubclass(w.category, AstropyWarning)

def test_kernel1d_int_size(self):
"""
Test that an error is raised if ``Kernel1D`` ``x_size`` is not
an integer.
"""
with pytest.raises(TypeError):
Gaussian1DKernel(3, x_size=1.2)

def test_kernel2d_int_xsize(self):
"""
Test that an error is raised if ``Kernel2D`` ``x_size`` is not
an integer.
"""
with pytest.raises(TypeError):
Gaussian2DKernel(3, x_size=1.2)

def test_kernel2d_int_ysize(self):
"""
Test that an error is raised if ``Kernel2D`` ``y_size`` is not
an integer.
"""
with pytest.raises(TypeError):
Gaussian2DKernel(3, x_size=5, y_size=1.2)

def test_kernel1d_initialization(self):
"""
Test that an error is raised if an array or model is not
specified for ``Kernel1D``.
"""
with pytest.raises(TypeError):
Kernel1D()

def test_kernel2d_initialization(self):
"""
Test that an error is raised if an array or model is not
specified for ``Kernel2D``.
"""
with pytest.raises(TypeError):
Kernel2D()

0 comments on commit 6ee4303

Please sign in to comment.