Ticket #873 (closed defect: fixed)
Convolution of a MaskedImage with a spatially varying LinearCombinationKernel can give incorrect variance
Reported by: | rowen | Owned by: | rowen |
---|---|---|---|
Priority: | normal | Milestone: | |
Component: | afw | Keywords: | |
Cc: | rhl, becker | Blocked By: | |
Blocking: | Project: | LSST | |
Version Number: | |||
How to repeat: |
An example is supplied on in afw/tickets/786/tests/ConvolveMaskedImage_1.py: look for DISABLEDtestTicket873 and delete the DISABLED to enable the test (you may wish to disable the other tests to speed it up). This reports: AssertionError: variance planes differ: maxDiff=22278.9649909 at position (68, 64); value=45422.5078125 vs. 23143.5428216 (for doNormalize=False, copyEdge=False) However, any spatially varying LCK that has basis kernels with non-overlapping 0s should show the problem (as long as you avoid DeltaFunctionKernels, since the convolution code contains a special case to handle those). |
Description
The code that optimizes the speed of convolution with a spatially varying LinearCombinationKernel? (LCK) works by treating each kernel separately. This is fine for computing the image and mask planes. However, it is not necessarily adequate for computing the variance plane because the kernel is squared for that computation. Thus there are cross terms and I don't know how to handle them correctly.
The code for DC3a used a clever fix from RHL that approximates the cross terms. In limited testing this fix appears to work quite well. However, this fix must be disabled for an LCK consisting of a delta function basis kernel set, because those basis kernels have no cross terms (a change first implemented on afw/tickets/786, and thus not in time for DC3a).
Clearly there must be intermediate cases where the variance will be mis-computed, e.g. basis kernels that have many zero-valued pixels. At an extereme, one can get incorrect variance by supplying delta function kernels implemented as FixedKernel? or AnalyticKernel? (since the code detects a delta function kernel by asking if it is an instance of DeltaFunctionKernel?).
The question is whether we can solve the problem in a general way while still preserving the speed advantage of handling each basis kernel separately.
Change History
comment:3 Changed 10 years ago by rowen
- Component changed from unknown to afw
RHL reports: I think that the correct value for "alpha" when the per-pixel variance is constant over the Kernel footprint is Spearman's correlation coefficient between the two Kernels K and L:
alpha = sum_i (K_i L_i ) ---------------------------------- sqrt(sum_i(K^2_i) sum_i(L^2_i))
For normalised kernels, this is just the dot product between the two kernels.
For delta-function kernels this is 0.0; for more general kernels it isn't 1.0 so I don't understand how we ever got away with the 1.0 "fix" that I put in.
comment:5 Changed 9 years ago by rowen
I believe the math is as follows:
Var_outpix = sum_inpix Var_inpix K_inpix^2 K_inpix = sum_n Cn_output Kn_inpix where Kn_inpix are the spatially invariant basis kernels and Cn_outpix are the spatially varying scalar coefficients of those basis kernels (hence _outpix) = sum_inpix Var_inpix (sum_n Cn_outpix Kn_inpix)^2 = sum_inpix Var_inpix (sum_n_m Cn_outpix Cm_outpix Kn_inpix Km_inpix)
I don't yet see an efficient way to compute this. The problem is that the actual amount of cross correlation depends on the input variance so it has to be computed at each point.
If we assume the variance is constant enough to ignore this as a second-order effect then the problem becomes more tractable; we can compute the cross-correlation terms in advance, e.g. as a matrix:
matrix_m_n = sum_inpix Kn_inpix Km_inpix
However, I doubt we should make this assumption.
My suspicion is that we would be better off exploring other ways of computing spatially varying linear combination kernels -- at least those that don't use delta function kernels. Perhaps ticket #1198.
comment:6 Changed 9 years ago by rowen
- Status changed from assigned to closed
- Resolution set to fixed
Fixed on the trunk in revision 15035. The fix is to convolve using linear interpolation instead of attempting to compute the cross terms for the basis kernels. Hence fixed as a side effect of ticket #1198. However, convolution with a delta function basis kernel still uses the old technique because it is faster for that case (for reasonable size kernels). Run examples/timeSpatiallyVaryingConvolution.py to see this.