Conversation
…gistration methods. In some cases ordering of transform compositions was also incorrect, this fixes that.
|
@ChantalJuntao after this PR gets merged I will make another one that adds QuadDIRECT-based registration for simple (sub-pixel) translation and for full affine transformations. Just wanted to let you know in case you planned to implement them yourself. |
|
The test failure is the same unrelated failure as in #79 (weird, I don't get that error locally). |
timholy
left a comment
timholy
left a comment
There was a problem hiding this comment.
This looks good, just a couple of suggestions.
src/RegisterMismatchCommon.jl
Outdated
| size(fixed) == size(moving) || throw(DimensionMismatch("Size $(size(fixed)) of fixed is not equal to size $(size(moving)) of moving")) | ||
| num = denom = zero(promote_type(typeof((oneunit(Tf) - oneunit(Tm))^2), typeof(oneunit(Tf)^2+oneunit(Tm)^2))) | ||
| if sizeof(typeof(num)) < 4 | ||
| num = denom = Float64(0) #prevents overflow in most cases |
There was a problem hiding this comment.
Very interesting. Looking down at the num += (vf-vm)^2 lines, it occurs to me that the unsigned types may not play well with this:
julia> x, y = 0xff, 0x0f
(0xff, 0x0f)
julia> x-y
0xf0
julia> y-x
0x10
julia> (y-x)^2
0x00
julia> (x-y)^2
0x00and with N0f8:
julia> x, y = reinterpret(N0f8, 0xff), reinterpret(N0f8, 0x0f)
(1.0N0f8, 0.059N0f8)
julia> x-y
0.941N0f8
julia> y-x
0.063N0f8
julia> (y-x)^2
0.004N0f8
julia> (x-y)^2
0.886N0f8I think it would probably be good to declare num::T, denom::T and then use vf, vm = T(fixed[i]), T(moving[i]) to ensure that the difference gets computed with the higher precision.
I wonder if we should just always use Float64? Is there a substantial cost to doing that?
Thanks for catching this!
There was a problem hiding this comment.
I wonder if we should just always use Float64? Is there a substantial cost to doing that?
I thought there might be, but I didn't check the performance. I also avoided T(fixed[i]) for performance reasons, but on second thought we should probably prioritize safety as you suggest.
test/register_mismatch.jl
Outdated
| @test nd0.denom ≈ nd1.denom | ||
|
|
||
| # Test for denom overflow with mismatch0 | ||
| C = rand(N0f16, 7,9)./3 |
There was a problem hiding this comment.
rand might not be the best choice here since whether it overflows may vary from run-to-run. Can you come up with a concrete array that exhibited the overflow before this commit? Could be just a 2x2 matrix, probably.
test/register_mismatch.jl
Outdated
| # Test for denom overflow with mismatch0 | ||
| C = rand(N0f16, 7,9)./3 | ||
| D = rand(N0f16, 7,9)./3 | ||
| mm = RegisterMismatch.mismatch(C, D, (3,3)) |
There was a problem hiding this comment.
mismatch is more complicated than mismatch0, so it's a little weird using it to compute ground truth for mismatch0. If you have concrete values for C and D, then you could compare against the ground truth, right?
|
|
||
| #rotation + shift, slow because it warps for every rotation and shift | ||
| function slow_mm(tfm, fixed, moving, thresh, SD; initial_tfm = -1) | ||
| function slow_mm(tfm, fixed, moving, thresh, SD; initial_tfm = IdentityTransformation()) |
There was a problem hiding this comment.
Even though it's effectively the same thing I like this better.
|
Thanks for the feedback! I just pushed changes that I think address your comments. |
|
Thanks @Cody-G! Great detective work here. |
2 participants
This allows specifying an initial transform for all QuadDIRECT-based registration methods. After #79 I made sure to get the order of transformations correct! I also found and fixed a bug where the
denomaccumulated in themismatch0function can silently overflow when the element types of input images have insufficient dynamic range. Note that my fix doesn't prevent overflow in all cases...we may want to be even more cautious.