Repository navigation
Convolution fails when kernel sums to zero #1647
Description
The convolve (convolution.convolve) function always normalizes the kernel by its sum prior to performing the convolution:
# Because the Cython routines have to normalize the kernel on the fly, we
# explicitly normalize the kernel here, and then scale the image at the
# end if normalization was not requested.
kernel_sum = kernel_internal.sum()
kernel_internal /= kernel_sum
and scales the result after the convolution:
# If normalization was not requested, we need to scale the array (since
# the kernel was normalized prior to convolution)
if not normalize_kernel:
result *= kernel_sum
Obviously, this procedure fails when kernel_sum = 0. This can occur when using sharpening or edge-detection filters (e.g. kern = [[0, -1, 0], [-1, 4, -1], [0, -1, 0]]). Scipy's ndimage.convolve handles such a zero-sum kernel just fine.
Activity
astrofrog commentedon Oct 22, 2013
cc @cdeil @adonath - in case you have time to fix this (I'm very busy this week, but could look in a few days)
adonath commentedon Oct 23, 2013
I've taken a deeper look into this issue and unfortunately it's not only the problem @larrybradley pointed out.
Due to the way the actual convolution routines are implented now, it is not possible to use kernels which sum to zero.
The
NaNinterpolation requires a renormalization of the convolved value, because not all values of the array are taken into account, whenNaNvalues are present in the array.To fix this I need to understand better how the interpolation is supposed to work. E.g. what should be the expected result for the following example:
astrofrog commentedon Oct 23, 2013
I actually don't know what that should return... The NaN interpolation may not be appropriate with kernels that have negative values actually, I'm not sure... I'm going to have to think about this.
cdeil commentedon Oct 23, 2013
@adonath and I have discussed this a bit today. I think that
NaNregions in a spectrum or image are like boundaries ... there's different things one can do that make sense in different cases (e.g. treat as missing values, interpolate in some way, ...).Currently in the docs the "proper treatment" of
NaNvalues is highlighted as one of the major improvements of the astropy convolution routines over the ones innumpyandscipy, but what that means should be described more clearly in the convolve docstring, which currently states:I'd much prefer a formula and example, because I don't know what "NaN values are replaced with interpolated values using the kernel as an interpolation function" exactly means.
If we don't interpolate, but treat
NaNas missing values, I believe the problem with kernel re-normalisation disappears and all you have to do is apply a correction factor (number of pixels in the kernel) / (number of non-NaN pixels in the kernel)? I can give examples and code in case it's not clear.In any case I think that fixing the original issue (convolve can't handle zero-sum kernels) should take precedence over possible changes in the behaviour in and near NaN values in the input array.
cdeil commentedon Oct 23, 2013
In my last comment I forgot to mention why a correction factor is needed at all near
NaNvalues if they are treated as missing values. I'll do it with an example:[1, 1, 1, NaN, 1, 1, 1]convolved with[1, 1, 1]should give[3, 3, 3, 3, 3, 3, 3]and not[3, 3, 2, 2, 2, 3, 3]like you would get if you simply replaceNaNwith0in the input array.astrofrog commentedon Oct 30, 2013
I think we are going to have to delay this until 0.3.1. It is a bug, but I just don't think we have time to reach a consensus and fix it before Friday.
@cdeil @adonath - if you do think you can fix the original issue here and/or fix the treatment of NaN values, please do open a pull request! I just can't do it this week as I am swamped with other things.
astrofrog commentedon May 21, 2014
@keflavich - this is related to the other issue (can't find the number right now) because if we turn off normalization by default, this will no longer be an issue.
larrybradley commentedon Nov 21, 2014
I think we should treat
NaNsimply as missing values. Relatedly, I'd also like to seeconvolution.convolvetake amaskas input, where masked pixels would also be excluded (in the same way asNaN). If the user wants to interpolate overNaN, then they should do that before runningconvolve.cdeil commentedon Sep 10, 2015
I haven't imported
astropy.convolvein a long time and am busy with other things, so please un-assign me.keflavich commentedon Jun 19, 2017
Closed by #5782
[util] prevent convolve from trying to interpolate nan values