diff --git a/src/RegisterQD.jl b/src/RegisterQD.jl index 6eef47c..4e7c320 100644 --- a/src/RegisterQD.jl +++ b/src/RegisterQD.jl @@ -35,7 +35,7 @@ function qd_rigid(fixed, moving, mxshift::VecLike, mxrot::Union{Number,VecLike}, end function qd_affine(fixed, moving, mxshift, linmins, linmaxs, SD; - thresh=0.5*sum(abs2.(fixed[.!(isnan.(fixed))])), + thresh=0.5*sum(_abs2.(fixed[.!(isnan.(fixed))])), initial_tfm=IdentityTransformation(), print_interval=100, kwargs...) diff --git a/src/affine.jl b/src/affine.jl index d6d5224..f907ba9 100644 --- a/src/affine.jl +++ b/src/affine.jl @@ -84,7 +84,7 @@ which is why this is "coarse" optimization. function qd_affine_coarse(fixed, moving, mxshift, linmins, linmaxs; SD=I, initial_tfm=IdentityTransformation(), - thresh=0.1*sum(abs2.(fixed[.!(isnan.(fixed))])), + thresh=0.1*sum(_abs2.(fixed[.!(isnan.(fixed))])), minwidth=default_lin_minwidths(moving), maxevals=5e4, kwargs...) @@ -117,7 +117,7 @@ As a consequence, any large translations must be supplied with reasonable accura function qd_affine_fine(fixed, moving, linmins, linmaxs; SD=I, initial_tfm=IdentityTransformation(), - thresh=0.1*sum(abs2.(fixed[.!(isnan.(fixed))])), + thresh=0.1*sum(_abs2.(fixed[.!(isnan.(fixed))])), minwidth_mat=default_lin_minwidths(fixed)./10, maxevals=5e4, kwargs...) @@ -184,7 +184,7 @@ overlap between the two images; with non-zero `thresh`, it is not permissible to function qd_affine(fixed, moving, mxshift, linmins, linmaxs; presmoothed=false, SD=I, - thresh=0.5*sum(abs2.(fixed[.!(isnan.(fixed))])), + thresh=0.5*sum(_abs2.(fixed[.!(isnan.(fixed))])), initial_tfm=IdentityTransformation(), print_interval=100, kwargs...) diff --git a/src/rigid.jl b/src/rigid.jl index 116bc30..83b95ef 100644 --- a/src/rigid.jl +++ b/src/rigid.jl @@ -86,7 +86,7 @@ end function qd_rigid_coarse(fixed, moving, mxshift, mxrot, minwidth_rot; SD=I, initial_tfm=IdentityTransformation(), - thresh=0.1*sum(abs2.(fixed[.!(isnan.(fixed))])), + thresh=0.1*sum(_abs2.(fixed[.!(isnan.(fixed))])), kwargs...) #note: if a trial rotation results in image overlap < thresh for all possible shifts then QuadDIRECT throws an error f(x) = rigid_mm_fast(x, mxshift, fixed, moving, thresh, SD; initial_tfm=initial_tfm) @@ -104,7 +104,7 @@ end function qd_rigid_fine(fixed, moving, mxrot, minwidth_rot; SD=I, initial_tfm=IdentityTransformation(), - thresh=0.1*sum(abs2.(fixed[.!(isnan.(fixed))])), + thresh=0.1*sum(_abs2.(fixed[.!(isnan.(fixed))])), kwargs...) f(x) = rigid_mm_slow(x, fixed, moving, thresh, SD; initial_tfm=initial_tfm) upper_shft = fill(2.0, ndims(fixed)) @@ -179,7 +179,7 @@ function qd_rigid(fixed, moving, mxshift::VecLike, mxrot::Union{Number,VecLike}; presmoothed=false, SD=I, minwidth_rot=default_minwidth_rot(CartesianIndices(fixed), SD), - thresh=0.1*sum(abs2.(fixed[.!(isnan.(fixed))])), + thresh=0.1*sum(_abs2.(fixed[.!(isnan.(fixed))])), initial_tfm=IdentityTransformation(), print_interval=100, kwargs...) diff --git a/src/translations.jl b/src/translations.jl index 49b5393..c615d9d 100644 --- a/src/translations.jl +++ b/src/translations.jl @@ -16,7 +16,7 @@ end function qd_translate_fine(fixed, moving; initial_tfm=IdentityTransformation(), minwidth=fill(0.01, ndims(fixed)), - thresh=0.1*sum(abs2.(fixed[.!(isnan.(fixed))])), + thresh=0.1*sum(_abs2.(fixed[.!(isnan.(fixed))])), kwargs...) f(x) = translate_mm_slow(x, fixed, moving, thresh; initial_tfm=initial_tfm) upper = fill(1.0, ndims(fixed)) @@ -50,13 +50,13 @@ If you have a good initial guess at the solution, pass it with the `initial_tfm` `thresh` enforces a certain amount of sum-of-squared-intensity overlap between the two images; with non-zero `thresh`, it is not permissible to "align" the images by shifting one entirely out of the way of the other. -If the `crop` keyword arg is `true` then `fixed` is cropped by `mxshift` (after the optional `initial_tfm`) on all sides +If the `crop` keyword arg is `true` then `fixed` is cropped by `mxshift` (after the optional `initial_tfm`) on all sides so that there will be complete overlap between `fixed` and `moving` for any evaluated shift. This avoids edge effects that can occur due to normalization when the transformed `moving` doesn't fully overlap with `fixed`. """ function qd_translate(fixed, moving, mxshift; presmoothed=false, - thresh=0.1*sum(abs2.(fixed[.!(isnan.(fixed))])), + thresh=0.1*sum(_abs2.(fixed[.!(isnan.(fixed))])), initial_tfm=IdentityTransformation(), minwidth=fill(0.01, ndims(fixed)), print_interval=100, crop=false, kwargs...) fixed, moving = float(fixed), float(moving) diff --git a/src/util.jl b/src/util.jl index 7fe010a..7959cc4 100644 --- a/src/util.jl +++ b/src/util.jl @@ -240,3 +240,8 @@ function _analyze(f, lower, upper; kwargs...) splits = ([[lower[i]; lower[i]+(upper[i]-lower[i])/2; upper[i]] for i=1:length(lower)]...,) QuadDIRECT.analyze(f, splits, lower, upper; kwargs...) end + +# abs2 was removed in https://github.com/JuliaGraphics/ColorVectorSpace.jl/pull/131 +# this is the current recommended replacement +# TODO: follow https://github.com/JuliaGraphics/ColorVectorSpace.jl/issues/157 and adjust for changes +_abs2(c) = mapreducec(v->float(v)^2, +, float(zero(eltype(c))), c) diff --git a/test/qd_random.jl b/test/qd_random.jl index 82d7c44..e49b700 100644 --- a/test/qd_random.jl +++ b/test/qd_random.jl @@ -5,6 +5,7 @@ using RegisterQD.Images using RegisterQD.CoordinateTransformations using RegisterQD.Rotations using RegisterQD.RegisterMismatch +using RegisterQD: _abs2 using Random #import BlockRegistration, RegisterOptimize @@ -23,7 +24,7 @@ using Test, TestImages itp = interpolate(newfixed, BSpline(Linear())) etp = extrapolate(itp, NaN) fixed = etp(Base.axes(moving)...) #often the warped array has one-too-many pixels in one or more dimensions due to extrapolation - thresh = 0.1 * sum(abs2.(fixed[.!(isnan.(fixed))])) + thresh = 0.1 * sum(_abs2.(fixed[.!(isnan.(fixed))])) mxshift = (10,10) tfm, mm = qd_translate(fixed, moving, mxshift; maxevals=1000, thresh=thresh, rtol=0) @@ -41,7 +42,7 @@ using Test, TestImages itp = interpolate(newfixed, BSpline(Linear())) etp = extrapolate(itp, NaN) fixed = etp(axes(moving)...) #often the warped array has one-too-many pixels in one or more dimensions due to extrapolation - thresh = 0.1 * sum(abs2.(fixed[.!(isnan.(fixed))])) + thresh = 0.1 * sum(_abs2.(fixed[.!(isnan.(fixed))])) mxshift = (5,5,5) tfm, mm = qd_translate(fixed, moving, mxshift; maxevals=1000, thresh=thresh, rtol=0) @@ -56,7 +57,7 @@ using Test, TestImages itp = interpolate(newfixed, BSpline(Linear())) etp = extrapolate(itp, NaN) fixed = etp(axes(moving)...) #often the warped array has one-too-many pixels in one or more dimensions due to extrapolation - thresh = 0.1 * sum(abs2.(fixed[.!(isnan.(fixed))])) + thresh = 0.1 * sum(_abs2.(fixed[.!(isnan.(fixed))])) mxshift = (10,10) mxrot = pi/90 minwidth_rot = [0.0002] @@ -73,7 +74,7 @@ using Test, TestImages itp = interpolate(newfixed, BSpline(Linear())) etp = extrapolate(itp, NaN) fixed = etp(axes(moving)...) #often the warped array has one-too-many pixels in one or more dimensions due to extrapolation - thresh = 0.1 * sum(abs2.(fixed[.!(isnan.(fixed))])) + thresh = 0.1 * sum(_abs2.(fixed[.!(isnan.(fixed))])) mxshift = (5,5,5) mxrot = [pi/90; pi/90; pi/90] minwidth_rot = fill(0.0002, 3) @@ -95,7 +96,7 @@ using Test, TestImages itp = interpolate(newfixed, BSpline(Linear())) etp = extrapolate(itp, NaN) fixed = etp(axes(moving)...) #often the warped array has one-too-many pixels in one or more dimensions due to extrapolation - thresh = 0.5 * sum(abs2.(fixed[.!(isnan.(fixed))])) + thresh = 0.5 * sum(_abs2.(fixed[.!(isnan.(fixed))])) mxshift = (5,5) SD = SDiagonal(@SVector(ones(ndims(fixed)))) @@ -118,7 +119,7 @@ using Test, TestImages #inds = intersect.(axes(moving), axes(newfixed)) #fixed = newfixed[inds...] #moving = moving[inds...] - #thresh = 0.1 * (sum(abs2.(fixed[.!(isnan.(fixed))]))+sum(abs2.(moving[.!(isnan.(moving))]))); + #thresh = 0.1 * (sum(_abs2.(fixed[.!(isnan.(fixed))]))+sum(_abs2.(moving[.!(isnan.(moving))]))); #mxshift = (10,10,10) #SD = eye(ndims(fixed)); @@ -141,7 +142,7 @@ using Test, TestImages #inds = intersect.(axes(moving), axes(newfixed)) #fixed = newfixed[inds...] #moving = moving[inds...] - #thresh = 0.5 * sum(abs2.(fixed[.!(isnan.(fixed))])); + #thresh = 0.5 * sum(_abs2.(fixed[.!(isnan.(fixed))])); #mxshift = (5,5,5) #SD = eye(ndims(fixed)); #@test RegisterOptimize.aff(vcat(tfm00.translation[:], tfm00.linear[:]), fixed, SD) == tfm0