diff --git a/REQUIRE b/REQUIRE index 67db770..e570459 100644 --- a/REQUIRE +++ b/REQUIRE @@ -1,4 +1,4 @@ -julia 0.5 +julia 0.6 Images Interpolations AffineTransforms 0.2.1 diff --git a/src/CachedInterpolations.jl b/src/CachedInterpolations.jl index f5b3ea0..ce44514 100644 --- a/src/CachedInterpolations.jl +++ b/src/CachedInterpolations.jl @@ -48,7 +48,7 @@ A single `CachedInterpolation` represents a "movable" i_2)`. An array of these objects thus implements an array-of-arrays interface. Create them with `cachedinterpolators`. """ -type CachedInterpolation{T,N,M,O,K} <: AbstractInterpolation{T,N,BSpline{Quadratic{InPlace}},OnCell} +mutable struct CachedInterpolation{T,N,M,O,K} <: AbstractInterpolation{T,N,BSpline{Quadratic{InPlace}},OnCell} # Note: M = N+K coefs::Array{T,M} # tiled array of 3x3x... buffers parent::Array{T,M} # the overall array (`P` in the documentation above) @@ -57,12 +57,12 @@ type CachedInterpolation{T,N,M,O,K} <: AbstractInterpolation{T,N,BSpline{Quadrat end const splitN = Base.IteratorsMD.split -Base.size{T,N}(itp::CachedInterpolation{T,N}) = splitN(size(itp.parent), Val{N})[1] -Base.size{T,N}(itp::CachedInterpolation{T,N}, d) = d <= N ? size(itp.parent, d) : 1 +Base.size(itp::CachedInterpolation{T,N}) where {T,N} = splitN(size(itp.parent), Val{N})[1] +Base.size(itp::CachedInterpolation{T,N}, d) where {T,N} = d <= N ? size(itp.parent, d) : 1 -Base.indices{T,N,M,O}(itp::CachedInterpolation{T,N,M,O}) = +Base.indices(itp::CachedInterpolation{T,N,M,O}) where {T,N,M,O} = map((x,o)->x-o, splitN(indices(itp.parent), Val{N})[1], O) -Base.indices{T,N,M,O}(itp::CachedInterpolation{T,N,M,O}, d) = +Base.indices(itp::CachedInterpolation{T,N,M,O}, d) where {T,N,M,O} = d <= N ? indices(itp.parent, d) - O[d] : Base.OneTo(1) """ @@ -87,7 +87,7 @@ coordinates within a tile. For example, one can mimic a `CenterIndexedArray` when `size(parent, d)` is odd for the first `N` dimensions and you supply `origin = (div(size(parent,1)+1, 2), ...)`. """ -function cachedinterpolators{T,M}(parent::Array{T,M}, N, origin=ntuple(d->0,N)) +function cachedinterpolators(parent::Array{T,M}, N, origin=ntuple(d->0,N)) where {T,M} 0 <= N <= M || error("N must be between 0 and $M") length(origin) == N || throw(DimensionMismatch("length(origin) = $(length(origin)) is inconsistent with $N interpolating dimensions")) sz3 = ntuple(d->d<=N ? 3 : size(parent,d), M)::NTuple{M,Int} @@ -100,7 +100,7 @@ function cachedinterpolators{T,M}(parent::Array{T,M}, N, origin=ntuple(d->0,N)) end # function-barriered to circumvent type-instability in sztiles -@noinline function cachedinterpolators{T,N,M,K}(buffer::Array{T,M}, parent::Array{T,M}, origin::NTuple{N,Int}, center::NTuple{N,Int}, sztiles::NTuple{K,Int}) +@noinline function cachedinterpolators(buffer::Array{T,M}, parent::Array{T,M}, origin::NTuple{N,Int}, center::NTuple{N,Int}, sztiles::NTuple{K,Int}) where {T,N,M,K} itps = Array{CachedInterpolation{T,N,M,origin}}(sztiles) for tileindex in CartesianRange(sztiles) itps[tileindex] = CachedInterpolation{T,N,M,origin,K}(buffer, parent, center, tileindex) @@ -108,7 +108,7 @@ end itps end -@generated function Base.getindex{T,N,M,O}(itp::CachedInterpolation{T,N,M,O}, xs::Number...) +@generated function Base.getindex(itp::CachedInterpolation{T,N,M,O}, xs::Number...) where {T,N,M,O} length(xs) == N || error("Must use $N indexes") ibsyms = [Symbol("ib_", d) for d = 1:N] ipsyms = [Symbol("ip_", d) for d = 1:N] @@ -139,18 +139,18 @@ end end end -immutable CoefsWrapper{N,A} +struct CoefsWrapper{N,A} coefs::A end -Base.size{N}(itp::CoefsWrapper{N}) = splitN(size(itp.coefs), Val{N})[1] -Base.size{N}(itp::CoefsWrapper{N}, d) = d <= N ? size(itp.coefs, d) : 1 +Base.size(itp::CoefsWrapper{N}) where {N} = splitN(size(itp.coefs), Val{N})[1] +Base.size(itp::CoefsWrapper{N}, d) where {N} = d <= N ? size(itp.coefs, d) : 1 # FIXME: this function cheats dangerously, because it does _not_ # update the cache. This is equivalent to the assumption you've called # getindex for the current `(x_1, x_2, ...)` location before calling # gradient!. If this is not true, you'll get wrong answers. -@generated function Interpolations.gradient!{T,N,M,O,K}(g::AbstractVector, A::CachedInterpolation{T,N,M,O,K}, ys::Vararg{Number,N}) +@generated function Interpolations.gradient!(g::AbstractVector, A::CachedInterpolation{T,N,M,O,K}, ys::Vararg{Number,N}) where {T,N,M,O,K} IT = Tuple{ntuple(d->BSpline{Quadratic{InPlace}}, N)..., ntuple(d->NoInterp, K)...} BS = Interpolations.BSplineInterpolation{T,M,Array{T,M},IT,OnCell,0} ex = Interpolations.gradient_impl(BS) diff --git a/src/CenterIndexedArrays.jl b/src/CenterIndexedArrays.jl index 10f0d88..739c5bb 100644 --- a/src/CenterIndexedArrays.jl +++ b/src/CenterIndexedArrays.jl @@ -11,7 +11,7 @@ export CenterIndexedArray ## SymRange, an AbstractUnitRange that's symmetric around 0 # These are used as indices for CenterIndexedArrays -immutable SymRange <: AbstractUnitRange{Int} +struct SymRange <: AbstractUnitRange{Int} n::Int # goes from -n:n end @@ -34,9 +34,9 @@ Base.intersect(r::SymRange, s::SymRange) = SymRange(min(last(r), last(s))) s end -Base.promote_rule{UR<:AbstractUnitRange}(::Type{SymRange}, ::Type{UR}) = +Base.promote_rule(::Type{SymRange}, ::Type{UR}) where {UR<:AbstractUnitRange} = UR -Base.promote_rule{T2}(::Type{UnitRange{T2}}, ::Type{SymRange}) = +Base.promote_rule(::Type{UnitRange{T2}}, ::Type{SymRange}) where {T2} = UnitRange{promote_type(T2, Int)} function Base.convert(::Type{SymRange}, r::AbstractUnitRange) first(r) == -last(r) || error("cannot convert $r to a SymRange") @@ -53,26 +53,26 @@ A `CenterIndexedArray` is one for which the array center has indexes CenterIndexedArray(A) "converts" `A` into a CenterIndexedArray. All the sizes of `A` must be odd. """ -immutable CenterIndexedArray{T,N,A<:AbstractArray} <: AbstractArray{T,N} +struct CenterIndexedArray{T,N,A<:AbstractArray} <: AbstractArray{T,N} data::A halfsize::NTuple{N,Int} - function (::Type{CenterIndexedArray{T,N,A}}){T,N,A<:AbstractArray}(data::A) + function CenterIndexedArray{T,N,A}(data::A) where {T,N,A<:AbstractArray} new{T,N,A}(data, _halfsize(data)) end end -CenterIndexedArray{T,N}(A::AbstractArray{T,N}) = CenterIndexedArray{T,N,typeof(A)}(A) -CenterIndexedArray{T}(::Type{T}, dims) = CenterIndexedArray(Array{T}(dims)) -CenterIndexedArray{T}(::Type{T}, dims...) = CenterIndexedArray(Array{T}(dims)) +CenterIndexedArray(A::AbstractArray{T,N}) where {T,N} = CenterIndexedArray{T,N,typeof(A)}(A) +CenterIndexedArray(::Type{T}, dims) where {T} = CenterIndexedArray(Array{T}(dims)) +CenterIndexedArray(::Type{T}, dims...) where {T} = CenterIndexedArray(Array{T}(dims)) # This is the AbstractArray default, but do this just to be sure -@compat Base.IndexStyle{A<:CenterIndexedArray}(::Type{A}) = IndexCartesian() +Base.IndexStyle(::Type{A}) where {A<:CenterIndexedArray} = IndexCartesian() Base.size(A::CenterIndexedArray) = size(A.data) Base.indices(A::CenterIndexedArray) = map(SymRange, A.halfsize) -function Base.similar{T}(A::CenterIndexedArray, ::Type{T}, inds::Tuple{SymRange,Vararg{SymRange}}) +function Base.similar(A::CenterIndexedArray, ::Type{T}, inds::Tuple{SymRange,Vararg{SymRange}}) where T data = Array{T}(map(length, inds)) CenterIndexedArray(data) end @@ -86,7 +86,7 @@ function _halfsize(A::AbstractArray) map(n->n>>UInt(1), size(A)) end -@inline function Base.getindex{T,N}(A::CenterIndexedArray{T,N}, i::Vararg{Number,N}) +@inline function Base.getindex(A::CenterIndexedArray{T,N}, i::Vararg{Number,N}) where {T,N} @boundscheck checkbounds(A, i...) @inbounds val = A.data[map(offset, A.halfsize, i)...] val @@ -96,16 +96,16 @@ offset(off, i) = off+i+1 const Index = Union{Colon,AbstractVector} -Base.getindex{T}(A::CenterIndexedArray{T,1}, I::Index) = CenterIndexedArray([A[i] for i in _cindex(A, 1, I)]) -Base.getindex{T}(A::CenterIndexedArray{T,2}, I::Index, J::Index) = CenterIndexedArray([A[i,j] for i in _cindex(A,1,I), j in _cindex(A,2,J)]) -Base.getindex{T}(A::CenterIndexedArray{T,3}, I::Index, J::Index, K::Index) = CenterIndexedArray([A[i,j,k] for i in _cindex(A,1,I), j in _cindex(A,2,J), k in _cindex(A,3,K)]) +Base.getindex(A::CenterIndexedArray{T,1}, I::Index) where {T} = CenterIndexedArray([A[i] for i in _cindex(A, 1, I)]) +Base.getindex(A::CenterIndexedArray{T,2}, I::Index, J::Index) where {T} = CenterIndexedArray([A[i,j] for i in _cindex(A,1,I), j in _cindex(A,2,J)]) +Base.getindex(A::CenterIndexedArray{T,3}, I::Index, J::Index, K::Index) where {T} = CenterIndexedArray([A[i,j,k] for i in _cindex(A,1,I), j in _cindex(A,2,J), k in _cindex(A,3,K)]) _cindex(A::CenterIndexedArray, d, I::Range) = convert(SymRange, I) _cindex(A::CenterIndexedArray, d, I::AbstractVector) = error("unsupported, use a range") _cindex(A::CenterIndexedArray, d, ::Colon) = SymRange(A.halfsize[d]) -@inline function Base.setindex!{T,N}(A::CenterIndexedArray{T,N}, v, i::Vararg{Number,N}) +@inline function Base.setindex!(A::CenterIndexedArray{T,N}, v, i::Vararg{Number,N}) where {T,N} @boundscheck checkbounds(A, i...) @inbounds A.data[map(offset, A.halfsize, i)...] = v v diff --git a/src/RegisterCore.jl b/src/RegisterCore.jl index 65906b4..8392a74 100644 --- a/src/RegisterCore.jl +++ b/src/RegisterCore.jl @@ -192,7 +192,7 @@ This representation is efficient for `Interpolations.jl`, because it allows interpolation to be performed on "both arrays" at once without recomputing the interpolation coefficients. """ -immutable NumDenom{T<:Number} +struct NumDenom{T<:Number} num::T denom::T end @@ -206,14 +206,14 @@ NumDenom(n, d) = NumDenom(promote(n, d)...) (*)(n::Number, p::NumDenom) = NumDenom(n*p.num, n*p.denom) (*)(p::NumDenom, n::Number) = n*p (/)(p::NumDenom, n::Number) = NumDenom(p.num/n, p.denom/n) -Base.one{T}(::Type{NumDenom{T}}) = NumDenom(one(T),one(T)) +Base.one(::Type{NumDenom{T}}) where {T} = NumDenom(one(T),one(T)) Base.one(p::NumDenom) = one(typeof(p)) -Base.zero{T}(::Type{NumDenom{T}}) = NumDenom(zero(T),zero(T)) +Base.zero(::Type{NumDenom{T}}) where {T} = NumDenom(zero(T),zero(T)) Base.zero(p::NumDenom) = zero(typeof(p)) -Base.promote_rule{T1,T2<:Number}(::Type{NumDenom{T1}}, ::Type{T2}) = NumDenom{promote_type(T1,T2)} -Base.eltype{T}(::Type{NumDenom{T}}) = T -Base.convert{T}(::Type{NumDenom{T}}, p::NumDenom{T}) = p -Base.convert{T}(::Type{NumDenom{T}}, p::NumDenom) = NumDenom{T}(p.num, p.denom) +Base.promote_rule(::Type{NumDenom{T1}}, ::Type{T2}) where {T1,T2<:Number} = NumDenom{promote_type(T1,T2)} +Base.eltype(::Type{NumDenom{T}}) where {T} = T +Base.convert(::Type{NumDenom{T}}, p::NumDenom{T}) where {T} = p +Base.convert(::Type{NumDenom{T}}, p::NumDenom) where {T} = NumDenom{T}(p.num, p.denom) Base.show(io::IO, p::NumDenom) = print(io, "NumDenom(", p.num, ",", p.denom, ")") function Base.showcompact(io::IO, p::NumDenom) print(io, "NumDenom(") @@ -223,7 +223,7 @@ function Base.showcompact(io::IO, p::NumDenom) print(io, ")") end -@compat const MismatchArray{ND<:NumDenom,N,A} = CenterIndexedArray{ND,N,A} +const MismatchArray{ND<:NumDenom,N,A} = CenterIndexedArray{ND,N,A} maxshift(A::MismatchArray) = A.halfsize @@ -232,7 +232,7 @@ maxshift(A::MismatchArray) = A.halfsize `(num,denom)` into a single `MismatchArray`. This is useful preparation for interpolation. """ -function (::Type{M}){M<:MismatchArray}(num::AbstractArray, denom::AbstractArray) +function (::Type{M})(num::AbstractArray, denom::AbstractArray) where M<:MismatchArray size(num) == size(denom) || throw(DimensionMismatch("num and denom must have the same size")) T = promote_type(eltype(num), eltype(denom)) numdenom = CenterIndexedArray(NumDenom{T}, size(num)) @@ -270,7 +270,7 @@ end `mms = mismatcharrays(nums, denom)`, for `denom` a single array, uses the same `denom` array for all `nums`. """ -function mismatcharrays{A<:AbstractArray,T<:Number}(nums::AbstractArray{A}, denom::AbstractArray{T}) +function mismatcharrays(nums::AbstractArray{A}, denom::AbstractArray{T}) where {A<:AbstractArray,T<:Number} first = true local mms for i in eachindex(nums) @@ -285,7 +285,7 @@ function mismatcharrays{A<:AbstractArray,T<:Number}(nums::AbstractArray{A}, deno mms end -function mismatcharrays{A1<:AbstractArray,A2<:AbstractArray}(nums::AbstractArray{A1}, denoms::AbstractArray{A2}) +function mismatcharrays(nums::AbstractArray{A1}, denoms::AbstractArray{A2}) where {A1<:AbstractArray,A2<:AbstractArray} size(nums) == size(denoms) || throw(DimensionMismatch("nums and denoms arrays must have the same number of apertures")) first = true local mms @@ -304,7 +304,7 @@ end `num, denom = separate(mm)` splits an `AbstractArray{NumDenom}` into separate numerator and denominator arrays. """ -function separate{T}(data::AbstractArray{NumDenom{T}}) +function separate(data::AbstractArray{NumDenom{T}}) where T num = Array{T}(size(data)) denom = similar(num) for I in eachindex(data) @@ -320,7 +320,7 @@ function separate(mm::MismatchArray) CenterIndexedArray(num), CenterIndexedArray(denom) end -function separate{M<:MismatchArray}(mma::AbstractArray{M}) +function separate(mma::AbstractArray{M}) where M<:MismatchArray T = eltype(eltype(M)) nums = Array{CenterIndexedArray{T,ndims(M)}}(size(mma)) denoms = similar(nums) @@ -336,7 +336,7 @@ end `denom < thresh`, and `fillval`'s type determines the type of the output. The default is NaN. """ -function ratio{T}(mm::MismatchArray, thresh, fillval::T) +function ratio(mm::MismatchArray, thresh, fillval::T) where T out = CenterIndexedArray(T, size(mm)) for I in eachindex(mm) nd = mm[I] @@ -345,12 +345,12 @@ function ratio{T}(mm::MismatchArray, thresh, fillval::T) out end ratio(mm::MismatchArray, thresh) = ratio(mm, thresh, convert(eltype(eltype(mm)), NaN)) -@inline ratio{T}(nd::NumDenom{T}, thresh, fillval=convert(T,NaN)) = nd.denom < thresh ? fillval : nd.num/nd.denom +@inline ratio(nd::NumDenom{T}, thresh, fillval=convert(T,NaN)) where {T} = nd.denom < thresh ? fillval : nd.num/nd.denom -ratio{T<:Real}(r::CenterIndexedArray{T}, thresh, fillval=convert(T,NaN)) = r +ratio(r::CenterIndexedArray{T}, thresh, fillval=convert(T,NaN)) where {T<:Real} = r -(::Type{M}){M<:MismatchArray,T}(::Type{T}, dims) = CenterIndexedArray(NumDenom{T}, dims) -(::Type{M}){M<:MismatchArray,T}(::Type{T}, dims...) = CenterIndexedArray(NumDenom{T}, dims) +(::Type{M})(::Type{T}, dims) where {M<:MismatchArray,T} = CenterIndexedArray(NumDenom{T}, dims) +(::Type{M})(::Type{T}, dims...) where {M<:MismatchArray,T} = CenterIndexedArray(NumDenom{T}, dims) function Base.copy!(M::MismatchArray, nd::Tuple{AbstractArray, AbstractArray}) num, denom = nd @@ -394,7 +394,7 @@ arrays. end end -function indmin_mismatch{T<:Number}(r::CenterIndexedArray{T}) +function indmin_mismatch(r::CenterIndexedArray{T}) where T<:Number ind = ind2sub(size(r), indmin(r.data)) indctr = map(d->ind[d]-(size(r,d)+1)>>1, (1:ndims(r)...)) CartesianIndex(indctr) @@ -415,7 +415,7 @@ If you do not wish to highpass-filter along a particular axis, put You may optionally specify the element type of the result, which for `Integer` or `FixedPoint` inputs defaults to `Float32`. """ -function highpass{T}(::Type{T}, data::AbstractArray, sigma) +function highpass(::Type{T}, data::AbstractArray, sigma) where T if any(isinf, sigma) datahp = convert(Array{T,ndims(data)}, data) else @@ -424,7 +424,7 @@ function highpass{T}(::Type{T}, data::AbstractArray, sigma) datahp[datahp .< 0] = 0 # truncate anything below 0 datahp end -highpass{T<:AbstractFloat}(data::AbstractArray{T}, sigma) = highpass(T, data, sigma) +highpass(data::AbstractArray{T}, sigma) where {T<:AbstractFloat} = highpass(T, data, sigma) highpass(data::AbstractArray, sigma) = highpass(Float32, data, sigma) @@ -436,13 +436,8 @@ parent. See also `trimmedview`. """ paddedview(A::SubArray) = _paddedview(A, (), (), A.indexes...) -if VERSION < v"0.6.0-dev" - _paddedview{T,N,P,I}(A::SubArray{T,N,P,I}, newindexes, newsize) = - SubArray(A.parent, newindexes, newsize) -else - _paddedview{T,N,P,I}(A::SubArray{T,N,P,I}, newindexes, newsize) = - SubArray(A.parent, newindexes) -end +_paddedview(A::SubArray{T,N,P,I}, newindexes, newsize) where {T,N,P,I} = + SubArray(A.parent, newindexes) @inline function _paddedview(A, newindexes, newsize, index, indexes...) d = length(newindexes)+1 _paddedview(A, (newindexes..., pdindex(A.parent, d, index)), pdsize(A.parent, newsize, d, index), indexes...) @@ -480,7 +475,7 @@ end # For faster and type-stable slicing -immutable ColonFun end -(::Type{ColonFun})(::Int) = Colon() +struct ColonFun end +ColonFun(::Int) = Colon() end diff --git a/src/RegisterDeformation.jl b/src/RegisterDeformation.jl index 6a24a0e..03e6c28 100644 --- a/src/RegisterDeformation.jl +++ b/src/RegisterDeformation.jl @@ -9,6 +9,7 @@ import Interpolations: AbstractInterpolation, AbstractExtrapolation import ImageTransformations: warp, warp! # to avoid `scale` conflict: using AffineTransforms: AffineTransform, TransformedArray +import CoordinateTransformations: compose using Compat export @@ -35,7 +36,7 @@ export const DimsLike = Union{Vector{Int}, Dims} const InterpExtrap = Union{AbstractInterpolation,AbstractExtrapolation} -@compat const Extrapolatable{T,N} = Union{TransformedArray{T,N},AbstractExtrapolation{T,N}} +const Extrapolatable{T,N} = Union{TransformedArray{T,N},AbstractExtrapolation{T,N}} """ # RegisterDeformation @@ -65,11 +66,11 @@ The major functions/types exported by RegisterDeformation are: """ RegisterDeformation -@compat abstract type AbstractDeformation{T,N} end -Base.eltype{T,N}(::Type{AbstractDeformation{T,N}}) = T -Base.ndims{T,N}(::Type{AbstractDeformation{T,N}}) = N -Base.eltype{D<:AbstractDeformation}(::Type{D}) = eltype(supertype(D)) -Base.ndims{D<:AbstractDeformation}(::Type{D}) = ndims(supertype(D)) +abstract type AbstractDeformation{T,N} end +Base.eltype(::Type{AbstractDeformation{T,N}}) where {T,N} = T +Base.ndims(::Type{AbstractDeformation{T,N}}) where {T,N} = N +Base.eltype(::Type{D}) where {D<:AbstractDeformation} = eltype(supertype(D)) +Base.ndims(::Type{D}) where {D<:AbstractDeformation} = ndims(supertype(D)) Base.eltype(d::AbstractDeformation) = eltype(typeof(d)) Base.ndims(d::AbstractDeformation) = ndims(typeof(d)) @@ -94,11 +95,11 @@ Finally, `ϕ = GridDeformation((u1, u2, ...), ...)` allows you to construct the deformation using an `N`-tuple of shift-arrays, each with `N` dimensions. """ -immutable GridDeformation{T,N,A<:AbstractArray,L} <: AbstractDeformation{T,N} +struct GridDeformation{T,N,A<:AbstractArray,L} <: AbstractDeformation{T,N} u::A knots::NTuple{N,L} - function (::Type{GridDeformation{T,N,A,L}}){T,N,A,L,FV<:SVector}(u::AbstractArray{FV,N}, knots::NTuple{N,L}) + function GridDeformation{T,N,A,L}(u::AbstractArray{FV,N}, knots::NTuple{N,L}) where {T,N,A,L,FV<:SVector} typeof(u) == A || error("typeof(u) = $(typeof(u)), which is different from $A") length(FV) == N || throw(DimensionMismatch("Dimensionality $(length(FV)) must match $N knot vectors")) for d = 1:N @@ -106,30 +107,30 @@ immutable GridDeformation{T,N,A<:AbstractArray,L} <: AbstractDeformation{T,N} end new{T,N,A,L}(u, knots) end - function (::Type{GridDeformation{T,N,A,L}}){T,N,A,L,FV<:SVector}(u::ScaledInterpolation{FV,N}) + function GridDeformation{T,N,A,L}(u::ScaledInterpolation{FV,N}) where {T,N,A,L,FV<:SVector} new{T,N,A,L}(u, u.ranges) end end -@compat const InterpolatingDeformation{T,N,A<:AbstractInterpolation} = GridDeformation{T,N,A} +const InterpolatingDeformation{T,N,A<:AbstractInterpolation} = GridDeformation{T,N,A} # Ambiguity avoidance -function GridDeformation{FV<:SVector,N}(u::AbstractArray{FV,N}, - knots::Tuple{}) +function GridDeformation(u::AbstractArray{FV,N}, + knots::Tuple{}) where {FV<:SVector,N} error("Cannot supply an empty knot tuple") end # With knot ranges -function GridDeformation{FV<:SVector,N,L<:AbstractVector}(u::AbstractArray{FV,N}, - knots::NTuple{N,L}) +function GridDeformation(u::AbstractArray{FV,N}, + knots::NTuple{N,L}) where {FV<:SVector,N,L<:AbstractVector} T = eltype(FV) length(FV) == N || throw(DimensionMismatch("$N-dimensional array requires SVector{$N,T}")) GridDeformation{T,N,typeof(u),L}(u, knots) end # With image spatial size -function GridDeformation{FV<:SVector,N,L<:Integer}(u::AbstractArray{FV,N}, - dims::NTuple{N,L}) +function GridDeformation(u::AbstractArray{FV,N}, + dims::NTuple{N,L}) where {FV<:SVector,N,L<:Integer} T = eltype(FV) length(FV) == N || throw(DimensionMismatch("$N-dimensional array requires SVector{$N,T}")) knots = ntuple(d->linspace(1,dims[d],size(u,d)), N) @@ -137,7 +138,7 @@ function GridDeformation{FV<:SVector,N,L<:Integer}(u::AbstractArray{FV,N}, end # Construct from a plain array -function GridDeformation{T<:Number,N}(u::AbstractArray{T}, knots::NTuple{N}) +function GridDeformation(u::AbstractArray{T}, knots::NTuple{N}) where {T<:Number,N} ndims(u) == N+1 || error("Need $(N+1) dimensions for $N-dimensional deformations") size(u,1) == N || error("First dimension of u must be of length $N") uf = convert_to_fixed(SVector{N,T}, u, tail(size(u))) @@ -145,7 +146,7 @@ function GridDeformation{T<:Number,N}(u::AbstractArray{T}, knots::NTuple{N}) end # Construct from a (u1, u2, ...) tuple -function GridDeformation{N}(u::NTuple{N,AbstractArray}, knots::NTuple{N}) +function GridDeformation(u::NTuple{N,AbstractArray}, knots::NTuple{N}) where N ndims(u[1]) == N || error("Need $N dimensions for $N-dimensional deformations") ua = permutedims(cat(N+1, u...), (N+1,(1:N)...)) uf = convert_to_fixed(ua) @@ -153,10 +154,10 @@ function GridDeformation{N}(u::NTuple{N,AbstractArray}, knots::NTuple{N}) end # When knots is a vector -GridDeformation{V<:AbstractVector}(u, knots::AbstractVector{V}) = GridDeformation(u, (knots...,)) -GridDeformation{R<:Real}(u, knots::AbstractVector{R}) = GridDeformation(u, (knots,)) +GridDeformation(u, knots::AbstractVector{V}) where {V<:AbstractVector} = GridDeformation(u, (knots...,)) +GridDeformation(u, knots::AbstractVector{R}) where {R<:Real} = GridDeformation(u, (knots,)) -function GridDeformation{FV<:SVector}(u::ScaledInterpolation{FV}) +function GridDeformation(u::ScaledInterpolation{FV}) where FV<:SVector N = length(FV) ndims(u) == N || throw(DimensionMismatch("Dimension $(ndims(u)) incompatible with vectors of length $N")) GridDeformation{eltype(FV),N,typeof(u),typeof(u.ranges[1])}(u) @@ -167,14 +168,14 @@ end seqeuential deformations. The last dimension of the array `u` should correspond to time; `ϕs[i]` is produced from `u[..., i]`. """ -function griddeformations{N,T<:Number}(u::AbstractArray{T}, knots::NTuple{N}) +function griddeformations(u::AbstractArray{T}, knots::NTuple{N}) where {N,T<:Number} ndims(u) == N+2 || error("Need $(N+2) dimensions for a vector of $N-dimensional deformations") size(u,1) == N || error("First dimension of u must be of length $N") uf = RegisterDeformation.convert_to_fixed(SVector{N,T}, u, Base.tail(size(u))) griddeformations(uf, knots) end -function griddeformations{N,FV<:SVector}(u::AbstractArray{FV}, knots::NTuple{N}) +function griddeformations(u::AbstractArray{FV}, knots::NTuple{N}) where {N,FV<:SVector} ndims(u) == N+1 || error("Need $(N+1) dimensions for a vector of $N-dimensional deformations") length(FV) == N || throw(DimensionMismatch("Dimensionality $(length(FV)) must match $N knot vectors")) colons = ntuple(d->Colon(), N)::NTuple{N,Colon} @@ -183,7 +184,7 @@ end Base.:(==)(ϕ1::GridDeformation, ϕ2::GridDeformation) = ϕ1.u == ϕ2.u && ϕ1.knots == ϕ2.knots -Base.copy{T,N,A,L}(ϕ::GridDeformation{T,N,A,L}) = (u = copy(ϕ.u); GridDeformation{T,N,typeof(u),L}(u, map(copy, ϕ.knots))) +Base.copy(ϕ::GridDeformation{T,N,A,L}) where {T,N,A,L} = (u = copy(ϕ.u); GridDeformation{T,N,typeof(u),L}(u, map(copy, ϕ.knots))) # # TODO: flesh this out # immutable VoroiDeformation{T,N,Vu<:AbstractVector,Vc<:AbstractVector} <: AbstractDeformation{T,N} @@ -205,15 +206,15 @@ function Interpolations.interpolate!(ϕ::GridDeformation, BC) end Interpolations.interpolate!(ϕ::GridDeformation) = interpolate!(ϕ, InPlace()) -Interpolations.interpolate{ T,N,A<:AbstractInterpolation}(ϕ::GridDeformation{T,N,A}) = error("ϕ is already interpolating") +Interpolations.interpolate(ϕ::GridDeformation{T,N,A}) where { T,N,A<:AbstractInterpolation} = error("ϕ is already interpolating") -Interpolations.interpolate!{T,N,A<:AbstractInterpolation}(ϕ::GridDeformation{T,N,A}) = error("ϕ is already interpolating") +Interpolations.interpolate!(ϕ::GridDeformation{T,N,A}) where {T,N,A<:AbstractInterpolation} = error("ϕ is already interpolating") -function vecindex{T,N,A<:AbstractInterpolation}(ϕ::GridDeformation{T,N,A}, x::SVector{N}) +function vecindex(ϕ::GridDeformation{T,N,A}, x::SVector{N}) where {T,N,A<:AbstractInterpolation} x + vecindex(ϕ.u, x) end -@generated function Base.getindex{T,N,A<:AbstractInterpolation}(ϕ::GridDeformation{T,N,A}, xs::Number...) +@generated function Base.getindex(ϕ::GridDeformation{T,N,A}, xs::Number...) where {T,N,A<:AbstractInterpolation} length(xs) == N || throw(DimensionMismatch("$(length(xs)) indexes is not consistent with ϕ dimensionality $N")) xindexes = [:(xs[$d]) for d = 1:N] ϕxindexes = [:(xs[$d]+ux[$d]) for d = 1:N] @@ -226,7 +227,7 @@ end end # Composition ϕ_old(ϕ_new(x)) -function (ϕ_old::GridDeformation{T1,N,A}){T1,T2,N,A<:AbstractInterpolation}(ϕ_new::GridDeformation{T2,N}) +function (ϕ_old::GridDeformation{T1,N,A})(ϕ_new::GridDeformation{T2,N}) where {T1,T2,N,A<:AbstractInterpolation} uold, knots = ϕ_old.u, ϕ_old.knots if !isa(ϕ_new.u, AbstractInterpolation) ϕ_new.knots == knots || error("If knots are incommensurate, ϕ_new must be interpolating") @@ -257,7 +258,7 @@ end lookup(u::AbstractInterpolation, x, i) = vecindex(u, x) lookup(u, x, i) = u[i] -@generated function knot{N}(knots::NTuple{N}, i::Integer) +@generated function knot(knots::NTuple{N}, i::Integer) where N args = [:(knots[$d][s[$d]]) for d = 1:N] quote s = ind2sub(map(length, knots), i) @@ -265,14 +266,14 @@ lookup(u, x, i) = u[i] end end -@generated function knot{N}(knots::NTuple{N}, I) +@generated function knot(knots::NTuple{N}, I) where N args = [:(knots[$d][I[$d]]) for d = 1:N] :(SVector($(args...))) end arraysize(knots::NTuple) = map(k->(x = extrema(k); convert(Int, x[2]-x[1]+1)), knots) -immutable KnotIterator{K,N} +struct KnotIterator{K,N} knots::K iter::CartesianRange{N} end @@ -299,7 +300,7 @@ Reparametrize the deformation `ϕ` so that it has a grid size `gridsize`. ``` for a 3-dimensional deformation `ϕ`. """ -function regrid{T,N}(ϕ::InterpolatingDeformation{T,N}, sz::Dims{N}) +function regrid(ϕ::InterpolatingDeformation{T,N}, sz::Dims{N}) where {T,N} knots_new = map((r,n)->linspace(first(r), last(r), n), ϕ.knots, sz) u = Array{SVector{N,T},N}(sz) for (i, k) in enumerate(eachknot(knots_new)) @@ -307,7 +308,7 @@ function regrid{T,N}(ϕ::InterpolatingDeformation{T,N}, sz::Dims{N}) end GridDeformation(u, knots_new) end -regrid{T,N}(ϕ::GridDeformation{T,N}, sz::Dims{N}) = regrid(interpolate(ϕ), sz) +regrid(ϕ::GridDeformation{T,N}, sz::Dims{N}) where {T,N} = regrid(interpolate(ϕ), sz) """ `ϕ_c = ϕ_old(ϕ_new)` computes the composition of two deformations, @@ -321,7 +322,7 @@ position `(i,j,...)`. You can use `_, g = compose(identity, ϕ_new)` if you need the gradient for when `ϕ_old` is equal to the identity transformation. """ -function compose{T1,T2,N,A<:AbstractInterpolation}(ϕ_old::GridDeformation{T1,N,A}, ϕ_new::GridDeformation{T2,N}) +function compose(ϕ_old::GridDeformation{T1,N,A}, ϕ_new::GridDeformation{T2,N}) where {T1,T2,N,A<:AbstractInterpolation} u, knots = ϕ_old.u, ϕ_old.knots ϕ_new.knots == knots || error("Not yet implemented for incommensurate knots") unew = ϕ_new.u @@ -349,7 +350,7 @@ end `ϕsi_old` is interpolated ϕs_old: e.g) `ϕsi_old = map(Interpolations.interpolate!, copy(ϕs_old))` """ -function compose{G1<:GridDeformation, G2<:GridDeformation}(ϕsi_old::AbstractVector{G1}, ϕs_new::AbstractVector{G2}) +function compose(ϕsi_old::AbstractVector{G1}, ϕs_new::AbstractVector{G2}) where {G1<:GridDeformation, G2<:GridDeformation} n = length(ϕs_new) length(ϕsi_old) == n || throw(DimensionMismatch("vectors-of-deformations must have the same length, got $(length(ϕsi_old)) and $n")) ϕc1, g1 = compose(first(ϕsi_old), first(ϕs_new)) @@ -363,12 +364,12 @@ function compose{G1<:GridDeformation, G2<:GridDeformation}(ϕsi_old::AbstractVec end -function compose{T,N}(f::Function, ϕ_new::GridDeformation{T,N}) +function compose(f::Function, ϕ_new::GridDeformation{T,N}) where {T,N} f == identity || error("Only the identity function is supported") ϕ_new, fill(eye(similar_type(SArray, T, Size(N, N))), size(ϕ_new.u)) end -function medfilt{D<:AbstractDeformation}(ϕs::AbstractVector{D}, n) +function medfilt(ϕs::AbstractVector{D}, n) where D<:AbstractDeformation nhalf = n>>1 2nhalf+1 == n || error("filter size must be odd") T = eltype(eltype(D)) @@ -380,7 +381,7 @@ function medfilt{D<:AbstractDeformation}(ϕs::AbstractVector{D}, n) _medfilt!(ϕout, ϕs, v, vs) # function barrier due to instability of vs end -@noinline function _medfilt!{N,T}(ϕout, ϕs, v, vs::NTuple{N,T}) +@noinline function _medfilt!(ϕout, ϕs, v, vs::NTuple{N,T}) where {N,T} n = size(v,2) nhalf = n>>1 tmp = Vector{eltype(T)}(N) @@ -429,19 +430,19 @@ where evaluated anywhere. See the Interpolations package. - ϕ is an `AbstractDeformation` """ -type WarpedArray{T,N,A<:Extrapolatable,D<:AbstractDeformation} <: AbstractArray{T,N} +mutable struct WarpedArray{T,N,A<:Extrapolatable,D<:AbstractDeformation} <: AbstractArray{T,N} data::A ϕ::D end # User already supplied an interpolatable ϕ -function WarpedArray{T,N,S,A<:AbstractInterpolation}(data::Extrapolatable{T,N}, - ϕ::GridDeformation{S,N,A}) +function WarpedArray(data::Extrapolatable{T,N}, + ϕ::GridDeformation{S,N,A}) where {T,N,S,A<:AbstractInterpolation} WarpedArray{T,N,typeof(data),typeof(ϕ)}(data, ϕ) end # Create an interpolatable ϕ -function WarpedArray{T,N}(data::Extrapolatable{T,N}, ϕ::GridDeformation) +function WarpedArray(data::Extrapolatable{T,N}, ϕ::GridDeformation) where {T,N} itp = scale(interpolate(ϕ.u, BSpline(Quadratic(Flat())), OnCell()), ϕ.knots...) ϕ′ = GridDeformation(itp, ϕ.knots) WarpedArray{T,N,typeof(data),typeof(ϕ′)}(data, ϕ′) @@ -453,7 +454,7 @@ WarpedArray(data, ϕ::GridDeformation) = WarpedArray(to_etp(data), ϕ) Base.size(A::WarpedArray) = size(A.data) Base.size(A::WarpedArray, i::Integer) = size(A.data, i) -@generated function Base.getindex{T,N}(W::WarpedArray{T,N}, x::Number...) +@generated function Base.getindex(W::WarpedArray{T,N}, x::Number...) where {T,N} length(x) == N || error("Must use $N indexes") getindex_impl(N) end @@ -469,11 +470,7 @@ function getindex_impl(N) end end -if VERSION < v"0.5.0-dev" - getindex!(dest, W::WarpedArray, coords...) = Base._unsafe_getindex!(dest, Base.LinearSlow(), W, coords...) -else - getindex!(dest, W::WarpedArray, coords...) = Base._unsafe_getindex!(dest, W, coords...) -end +getindex!(dest, W::WarpedArray, coords...) = Base._unsafe_getindex!(dest, W, coords...) """ `Atrans = translate(A, displacement)` shifts `A` by an amount @@ -501,7 +498,7 @@ array "spins" around its center. The array of grid points defining `ϕ` has size specified by `gridsize`. The dimensionality of `tform` must match that specified by `arraysize` and `gridsize`. """ -function tform2deformation{T,Tv,N}(tform::AffineTransform{T,Tv,N}, arraysize, gridsize) +function tform2deformation(tform::AffineTransform{T,Tv,N}, arraysize, gridsize) where {T,Tv,N} if length(arraysize) != N || length(gridsize) != N error("Dimensionality mismatch") end @@ -531,16 +528,16 @@ function warp(img::AbstractArray, ϕ::AbstractDeformation) warp!(dest, wimg) end -warp_type{T<:AbstractFloat}(img::AbstractArray{T}) = T -warp_type{T<:Number}(img::AbstractArray{T}) = Float32 -warp_type{C<:Colorant}(img::AbstractArray{C}) = warp_type(img, eltype(eltype(C))) -warp_type{C<:Colorant, T<:AbstractFloat}(img::AbstractArray{C}, ::Type{T}) = C -warp_type{C<:Colorant, T}(img::AbstractArray{C}, ::Type{T}) = base_colorant_type(C){Float32} +warp_type(img::AbstractArray{T}) where {T<:AbstractFloat} = T +warp_type(img::AbstractArray{T}) where {T<:Number} = Float32 +warp_type(img::AbstractArray{C}) where {C<:Colorant} = warp_type(img, eltype(eltype(C))) +warp_type(img::AbstractArray{C}, ::Type{T}) where {C<:Colorant, T<:AbstractFloat} = C +warp_type(img::AbstractArray{C}, ::Type{T}) where {C<:Colorant, T} = base_colorant_type(C){Float32} """ `warp!(dest, src::WarpedArray)` instantiates a `WarpedArray` in the output `dest`. """ -@generated function warp!{T,N}(dest::AbstractArray{T,N}, src::WarpedArray) +@generated function warp!(dest::AbstractArray{T,N}, src::WarpedArray) where {T,N} ϕxindexes = [:(I[$d]+ux[$d]) for d = 1:N] quote size(dest) == size(src) || error("dest must have the same size as src") @@ -586,7 +583,7 @@ An alternative syntax is `warp!(T, io, img, uarray; [nworkers=1])`, where `uarray` is an array of `u` values with `size(uarray)[end] == nimages(img)`. """ -function warp!{T}(::Type{T}, dest::Union{IO,HDF5Dataset,JLD.JldDataset}, img, ϕs; nworkers=1) +function warp!(::Type{T}, dest::Union{IO,HDF5Dataset,JLD.JldDataset}, img, ϕs; nworkers=1) where T n = nimages(img) ssz = size(img, coords_spatial(img)...) if n == 1 @@ -609,12 +606,12 @@ function warp!{T}(::Type{T}, dest::Union{IO,HDF5Dataset,JLD.JldDataset}, img, ϕ nothing end -warp!{T,R<:Real}(::Type{T}, dest::Union{IO,HDF5Dataset,JLD.JldDataset}, img, u::AbstractArray{R}; kwargs...) = warp!(T, dest, img, convert_to_fixed(u); kwargs...) +warp!(::Type{T}, dest::Union{IO,HDF5Dataset,JLD.JldDataset}, img, u::AbstractArray{R}; kwargs...) where {T,R<:Real} = warp!(T, dest, img, convert_to_fixed(u); kwargs...) warp!(dest::Union{HDF5Dataset,JLD.JldDataset}, img, u; nworkers=1) = warp!(eltype(dest), dest, img, u; nworkers=nworkers) -function _warp!{T}(::Type{T}, dest, img, ϕs, nworkers) +function _warp!(::Type{T}, dest, img, ϕs, nworkers) where T n = nimages(img) ssz = size(img, coords_spatial(img)...) wpids = addprocs(nworkers) @@ -665,7 +662,7 @@ function warp_write(dest, destarray, i) dest[colons..., i] = destarray end -function extract1{V<:SVector}(u::AbstractArray{V}, N, ssz) +function extract1(u::AbstractArray{V}, N, ssz) where V<:SVector if ndims(u) == N+1 ϕ = GridDeformation(reshape(u, size(u)[1:end-1]), ssz) else @@ -673,22 +670,22 @@ function extract1{V<:SVector}(u::AbstractArray{V}, N, ssz) end ϕ end -extract1{D<:AbstractDeformation}(ϕs::Vector{D}, N, ssz) = ϕs[1] +extract1(ϕs::Vector{D}, N, ssz) where {D<:AbstractDeformation} = ϕs[1] -function extracti{V<:SVector}(u::AbstractArray{V}, i, ssz) +function extracti(u::AbstractArray{V}, i, ssz) where V<:SVector colons = [Colon() for d = 1:ndims(u)-1] GridDeformation(u[colons..., i], ssz) end -extracti{D<:AbstractDeformation}(ϕs::Vector{D}, i, _) = ϕs[i] +extracti(ϕs::Vector{D}, i, _) where {D<:AbstractDeformation} = ϕs[i] -function checkϕdims{V<:SVector}(u::AbstractArray{V}, N, n) +function checkϕdims(u::AbstractArray{V}, N, n) where V<:SVector ndims(u) == N+1 || error("u's dimensionality $(ndims(u)) is inconsistent with the number of spatial dimensions $N of the image") if size(u)[end] != n error("Must have one `u` slice per image") end nothing end -checkϕdims{D<:AbstractDeformation}(ϕs::Vector{D}, N, n) = length(ϕs) == n || error("Must have one `ϕ` per image") +checkϕdims(ϕs::Vector{D}, N, n) where {D<:AbstractDeformation} = length(ϕs) == n || error("Must have one `ϕ` per image") """ @@ -710,7 +707,7 @@ Optionally specify the element type `T` of `img`. See also [`warpgrid`](@ref). """ -function knotgrid{T,N}(::Type{T}, knots::NTuple{N,Range}) +function knotgrid(::Type{T}, knots::NTuple{N,Range}) where {T,N} @assert all(r->first(r)==1, knots) inds = map(r->Base.OneTo(ceil(Int, last(r))), knots) img = Array{T}(map(length, inds)) @@ -722,7 +719,7 @@ function knotgrid{T,N}(::Type{T}, knots::NTuple{N,Range}) end img end -knotgrid{T}(::Type{T}, ϕ::GridDeformation) = knotgrid(T, ϕ.knots) +knotgrid(::Type{T}, ϕ::GridDeformation) where {T} = knotgrid(T, ϕ.knots) knotgrid(arg) = knotgrid(Bool, arg) @@ -796,27 +793,27 @@ to_etp(img, A::AffineTransform) = TransformedArray(to_etp(img), A) # JLD extensions # Serializer for Array{SArray{T,...}} -immutable ArraySArraySerializer{T,DF,DT} +struct ArraySArraySerializer{T,DF,DT} A::Array{T,DT} end -function JLD.readas{T,DT}(serdata::ArraySArraySerializer{T,1,DT}) +function JLD.readas(serdata::ArraySArraySerializer{T,1,DT}) where {T,DT} sz = size(serdata.A) reinterpret(SVector{sz[1],T}, serdata.A, tail(sz)) end -function JLD.readas{T,DT}(serdata::ArraySArraySerializer{T,2,DT}) +function JLD.readas(serdata::ArraySArraySerializer{T,2,DT}) where {T,DT} sz = size(serdata.A) reinterpret(similar_type(SArray,T,Size(sz[1],sz[2])), serdata.A, tail(tail(sz))) end -function JLD.writeas{FSA<:StaticArray}(A::Array{FSA}) +function JLD.writeas(A::Array{FSA}) where FSA<:StaticArray T = eltype(FSA) DF = ndims(FSA) ArraySArraySerializer{T,DF,DF+ndims(A)}(reinterpret(T, A, (size(FSA)..., size(A)...))) end # Extensions to Interpolations and StaticArrays -@generated function vecindex{N}(A::AbstractArray, x::SVector{N}) +@generated function vecindex(A::AbstractArray, x::SVector{N}) where N args = [:(x[$d]) for d = 1:N] meta = Expr(:meta, :inline) quote @@ -825,7 +822,7 @@ end end end -@generated function vecgradient!{N}(g, itp::AbstractArray, x::SVector{N}) +@generated function vecgradient!(g, itp::AbstractArray, x::SVector{N}) where N args = [:(x[$d]) for d = 1:N] meta = Expr(:meta, :inline) quote @@ -834,13 +831,13 @@ end end end -function convert_to_fixed{T}(u::Array{T}, sz=size(u)) +function convert_to_fixed(u::Array{T}, sz=size(u)) where T N = sz[1] convert_to_fixed(SVector{N, T}, u, tail(sz)) end # Unlike the one above, this is type-stable -function convert_to_fixed{T,N}(::Type{SVector{N,T}}, u::AbstractArray{T}, sz=tail(size(u))) +function convert_to_fixed(::Type{SVector{N,T}}, u::AbstractArray{T}, sz=tail(size(u))) where {T,N} if isbits(T) uf = reinterpret(SVector{N,T}, u, sz) else @@ -850,7 +847,7 @@ function convert_to_fixed{T,N}(::Type{SVector{N,T}}, u::AbstractArray{T}, sz=tai uf end -@generated function copy_ctf!{N,T}(dest::Array{SVector{N,T}}, src::Array) +@generated function copy_ctf!(dest::Array{SVector{N,T}}, src::Array) where {N,T} exvec = [:(src[offset+$d]) for d=1:N] quote for i = 1:length(dest) @@ -861,7 +858,7 @@ end end end -function convert_from_fixed{N,T}(uf::AbstractArray{SVector{N,T}}, sz=size(uf)) +function convert_from_fixed(uf::AbstractArray{SVector{N,T}}, sz=size(uf)) where {N,T} if isbits(T) && isa(uf, Array) u = reinterpret(T, uf, (N, sz...)) else @@ -881,12 +878,10 @@ end # :(SMatrix{Tuple{R,C},T}(($(args...),))) # end -if VERSION >= v"0.6.0-rc2" - include_string(""" +include_string(""" const RowVectorArray{T} = RowVector{T,Vector{T}} Base.reinterpret{T,S,N}(::Type{T}, r::RowVectorArray{S}, dims::Dims{N}) = reinterpret(T, r.vec, dims) """) -end end # module diff --git a/src/RegisterFit.jl b/src/RegisterFit.jl index 5a65b65..a7a7f6f 100644 --- a/src/RegisterFit.jl +++ b/src/RegisterFit.jl @@ -187,7 +187,7 @@ end `tf = uisvalid(u, maxshift)` returns `true` if all entries of `u` are within the allowed domain. """ -function uisvalid{T<:Number}(u::AbstractArray{T}, maxshift) +function uisvalid(u::AbstractArray{T}, maxshift) where T<:Number nd = size(u,1) sztail = Base.tail(size(u)) for j in CartesianRange(sztail), idim = 1:nd @@ -201,7 +201,7 @@ end """ `u = uclamp!(u, maxshift)` clamps the values of `u` to the allowed domain. """ -function uclamp!{T<:Number}(u::AbstractArray{T}, maxshift) +function uclamp!(u::AbstractArray{T}, maxshift) where T<:Number nd = size(u,1) sztail = Base.tail(size(u)) for j in CartesianRange(sztail), idim = 1:nd @@ -210,7 +210,7 @@ function uclamp!{T<:Number}(u::AbstractArray{T}, maxshift) u end -function uclamp!{T<:StaticVector}(u::AbstractArray{T}, maxshift) +function uclamp!(u::AbstractArray{T}, maxshift) where T<:StaticVector uclamp!(reinterpret(eltype(T), u, (length(T), size(u)...)), maxshift) u end @@ -220,7 +220,7 @@ end image `img`. `center` is the centroid of intensity, and `cov` the covariance matrix of the intensity. """ -function principalaxes{T,N}(img::AbstractArray{T,N}) +function principalaxes(img::AbstractArray{T,N}) where {T,N} Ts = typeof(zero(T)/1) psums = pa_init(Ts, size(img)) # partial sums along all but one axis # Use a two-pass algorithm @@ -259,8 +259,8 @@ function principalaxes{T,N}(img::AbstractArray{T,N}) m, cov end -@noinline pa_init{S}(::Type{S}, sz) = [zeros(S, s) for s in sz] -@noinline function pa_centroid{S}(psums::Vector{Vector{S}}) +@noinline pa_init(::Type{S}, sz) where {S} = [zeros(S, s) for s in sz] +@noinline function pa_centroid(psums::Vector{Vector{S}}) where S s = sum(psums[1]) s, S[sum(psums[d] .* (1:length(psums[d]))) for d = 1:length(psums)] / s end @@ -331,13 +331,13 @@ end #### Low-level utilities -@generated function qfit_core!{T,N}(dE::Array{T,2}, V4::Array{T,2}, C::Array{T,4}, mm::Array{NumDenom{T},N}, thresh, umin::NTuple{N,Int}, E0, maxsep::NTuple{N,Int}) +@generated function qfit_core!(dE::Array{T,2}, V4::Array{T,2}, C::Array{T,4}, mm::Array{NumDenom{T},N}, thresh, umin::NTuple{N,Int}, E0, maxsep::NTuple{N,Int}) where {T,N} # The cost of generic matrix-multiplies is too high, so we write # these out by hand. quote - @nexprs $N i->(@nexprs $N j->j(@nexprs $N j->j(umin_d = umin[d]) - @nexprs $N iter1->(@nexprs $N iter2->iter2iter3iter4(@nexprs $N iter2->iter2iter3iter4(abs(i_d-umin[d]) > maxsep[d]) d->(continue) d->nothing nd = @nref $N mm i @@ -348,12 +348,12 @@ end @nexprs $N d->(v2 += v_d*v_d) r = num/den dE0 = r-E0 - @nexprs $N j->(@nexprs $N k->k(@nexprs $N iter2->iter2iter3iter4(@nexprs $N k->k(@nexprs $N iter2->iter2iter3iter4(@nexprs $N j->j(@nexprs $N iter2->iter2iter3iter4(@nexprs $N j->j(@nexprs $N iter2->iter2iter3iter4Expr(:(=), Symbol("x_",d), e), 1:N, uindexes)...) ϕxindexes = [:(I[$d] + y[$d]) for d = 1:N] @@ -145,10 +145,10 @@ function penalty_hindsight_data!_gen{IT}(N, ::Type{IT}, Pad) end end -scaledindexes{IT}(::Type{IT}, N) = +scaledindexes(::Type{IT}, N) where {IT} = map(d->Interpolations.iextract(IT, d) != NoInterp ? :(Interpolations.coordlookup(knots[$d], I[$d])) : :(I[$d]), 1:N) -interprange{IT}(::Type{IT}, N::Integer) = _interprange((), IT, N) +interprange(::Type{IT}, N::Integer) where {IT} = _interprange((), IT, N) function _interprange(out, IT, N) if length(out) < N return _interprange((out..., interprange(Interpolations.iextract(IT, length(out)+1))), IT, N) @@ -158,13 +158,13 @@ end interprange(::Type{NoInterp}) = 0:0 interprange(::Type{BSpline{Constant}}) = 0:0 interprange(::Type{BSpline{Linear}}) = 0:1 -interprange{Q<:Quadratic}(::Type{BSpline{Q}}) = -1:1 -interprange{C<:Cubic}(::Type{BSpline{C}}) = -1:2 +interprange(::Type{BSpline{Q}}) where {Q<:Quadratic} = -1:1 +interprange(::Type{BSpline{C}}) where {C<:Cubic} = -1:2 -coef_gen{IT,N}(::Type{IT}, d::Integer, offsets::CartesianIndex{N}) = +coef_gen(::Type{IT}, d::Integer, offsets::CartesianIndex{N}) where {IT,N} = coef_gen(Interpolations.iextract(IT, d+1), IT, d+1, offsets) -function coef_gen{D<:Degree,IT<:DimSpec{BSpline}, N}(::Type{BSpline{D}}, ::Type{IT}, d::Integer, offsets::CartesianIndex{N}) +function coef_gen(::Type{BSpline{D}}, ::Type{IT}, d::Integer, offsets::CartesianIndex{N}) where {D<:Degree,IT<:DimSpec{BSpline}, N} if d <= N sym = offsets[d] == -1 ? Symbol("cm_",d) : offsets[d] == 0 ? Symbol("c_",d) : @@ -176,10 +176,10 @@ function coef_gen{D<:Degree,IT<:DimSpec{BSpline}, N}(::Type{BSpline{D}}, ::Type{ end end -@generated function penalty_hindsight_data{T,N,T1,T2}(ϕ1::InterpolatingDeformation{T,N}, - ϕ2::InterpolatingDeformation{T,N}, - fixed::AbstractArray{T1,N}, - moving::AbstractInterpolation{T2,N}) +@generated function penalty_hindsight_data(ϕ1::InterpolatingDeformation{T,N}, + ϕ2::InterpolatingDeformation{T,N}, + fixed::AbstractArray{T1,N}, + moving::AbstractInterpolation{T2,N}) where {T,N,T1,T2} uindexes = [:((I[$d]-offsets[$d])/steps[$d] + 1) for d = 1:N] ϕ1xindexes = [:(I[$d] + u1[$d]) for d = 1:N] ϕ2xindexes = [:(I[$d] + u2[$d]) for d = 1:N] diff --git a/src/RegisterMismatch.jl b/src/RegisterMismatch.jl index 4ac578a..43375ee 100644 --- a/src/RegisterMismatch.jl +++ b/src/RegisterMismatch.jl @@ -37,7 +37,7 @@ RegisterMismatch FFTW.set_num_threads(min(Sys.CPU_CORES, 8)) const FFTPROD = [2,3] -type NanCorrFFTs{T<:AbstractFloat,N} +mutable struct NanCorrFFTs{T<:AbstractFloat,N} I0::RCpair{T,N} I1::RCpair{T,N} I2::RCpair{T,N} @@ -52,7 +52,7 @@ computations over domains of size `aperture_width`, computing the mismatch up to shifts of size `maxshift`. The keyword arguments allow you to control the planning process for the FFTs. """ -type CMStorage{T<:AbstractFloat,N} +mutable struct CMStorage{T<:AbstractFloat,N} aperture_width::Vector{Float64} maxshift::Vector{Int} getindexes::Vector{UnitRange{Int}} # indexes for pulling padded data, in source-coordinates @@ -66,7 +66,7 @@ type CMStorage{T<:AbstractFloat,N} ifftfunc!::Function shiftindexes::Vector{Vector{Int}} # indexes for performing fftshift & snipping from -maxshift:maxshift - function (::Type{CMStorage{T,N}}){T,N}(::Type{T}, aperture_width::WidthLike, maxshift::DimsLike; flags=FFTW.ESTIMATE, timelimit=Inf, display=true) + function CMStorage{T,N}(::Type{T}, aperture_width::WidthLike, maxshift::DimsLike; flags=FFTW.ESTIMATE, timelimit=Inf, display=true) where {T,N} blocksize = map(x->ceil(Int,x), aperture_width) length(blocksize) == length(maxshift) || error("Dimensionality mismatch") padsz = padsize(blocksize, maxshift) @@ -95,10 +95,10 @@ type CMStorage{T<:AbstractFloat,N} new{T,N}(Float64[aperture_width...], maxshiftv, getindexes, padded, fixed, moving, buf1, buf2, fftfunc, ifftfunc, shiftindexes) end end -CMStorage{T<:Real}(::Type{T}, aperture_width, maxshift; kwargs...) = CMStorage{T,length(aperture_width)}(T, aperture_width, maxshift; kwargs...) +CMStorage(::Type{T}, aperture_width, maxshift; kwargs...) where {T<:Real} = CMStorage{T,length(aperture_width)}(T, aperture_width, maxshift; kwargs...) -eltype{T,N}(cms::CMStorage{T,N}) = T - ndims{T,N}(cms::CMStorage{T,N}) = N +eltype(cms::CMStorage{T,N}) where {T,N} = T + ndims(cms::CMStorage{T,N}) where {T,N} = N """ `mm = mismatch([T], fixed, moving, maxshift; @@ -111,7 +111,7 @@ normalization scheme (`:intensity` or `:pixels`). `fixed` and `moving` must have the same size; you can pad with `NaN`s as needed. See `nanpad`. """ -function mismatch{T<:Real}(::Type{T}, fixed::AbstractArray, moving::AbstractArray, maxshift::DimsLike; normalization = :intensity) +function mismatch(::Type{T}, fixed::AbstractArray, moving::AbstractArray, maxshift::DimsLike; normalization = :intensity) where T<:Real assertsamesize(fixed, moving) maxshiftv = tovec(maxshift) msz = 2maxshiftv.+1 @@ -190,15 +190,15 @@ in a rectangular grid, you can use an `N`-dimensional array-of-tuples (or array-of-vectors) or an `N+1`-dimensional array with the center positions specified along the first dimension. See `aperture_grid`. """ -function mismatch_apertures{T}(::Type{T}, - fixed::AbstractArray, - moving::AbstractArray, - aperture_centers::AbstractArray, - aperture_width::WidthLike, - maxshift::DimsLike; - normalization = :pixels, - flags = FFTW.MEASURE, - kwargs...) +function mismatch_apertures(::Type{T}, + fixed::AbstractArray, + moving::AbstractArray, + aperture_centers::AbstractArray, + aperture_width::WidthLike, + maxshift::DimsLike; + normalization = :pixels, + flags = FFTW.MEASURE, + kwargs...) where T nd = sdims(fixed) assertsamesize(fixed,moving) (length(aperture_width) == nd && length(maxshift) == nd) || error("Dimensionality mismatch") @@ -233,7 +233,7 @@ function mismatch_apertures!(mms, fixed, moving, aperture_centers, cms; normaliz end # Calculate the components needed to "nancorrelate" -function fftnan!{T<:Real}(out::NanCorrFFTs{T}, A::AbstractArray{T}, fftfunc!::Function) +function fftnan!(out::NanCorrFFTs{T}, A::AbstractArray{T}, fftfunc!::Function) where T<:Real I0 = real(out.I0) I1 = real(out.I1) I2 = real(out.I2) @@ -245,7 +245,7 @@ function fftnan!{T<:Real}(out::NanCorrFFTs{T}, A::AbstractArray{T}, fftfunc!::Fu out end -function _fftnan!{T<:Real}(I0, I1, I2, A::AbstractArray{T}) +function _fftnan!(I0, I1, I2, A::AbstractArray{T}) where T<:Real @inbounds for i in CartesianRange(size(A)) a = A[i] f = !isnan(a) @@ -256,14 +256,14 @@ function _fftnan!{T<:Real}(I0, I1, I2, A::AbstractArray{T}) end end -function fillfixed!{T}(cms::CMStorage{T}, fixed::AbstractArray) +function fillfixed!(cms::CMStorage{T}, fixed::AbstractArray) where T fill!(cms.padded, NaN) X = view(cms.padded, ntuple(d->(1:size(fixed,d))+cms.maxshift[d], ndims(fixed))...) copy!(X, fixed) fftnan!(cms.fixed, cms.padded, cms.fftfunc!) end -function fillfixed!{T}(cms::CMStorage{T}, fixed::SubArray) +function fillfixed!(cms::CMStorage{T}, fixed::SubArray) where T fill!(cms.padded, NaN) X = view(cms.padded, ntuple(d->(1:size(fixed,d))+cms.maxshift[d], ndims(fixed))...) get!(X, parent(fixed), parentindexes(fixed), convert(T, NaN)) @@ -272,7 +272,7 @@ end #### Utilities -Base.isnan{T}(A::Array{Complex{T}}) = isnan(real(A)) | isnan(imag(A)) +Base.isnan(A::Array{Complex{T}}) where {T} = isnan(real(A)) | isnan(imag(A)) function sumsq_finite(A) s = 0.0 for a in A diff --git a/src/RegisterMismatchCommon.jl b/src/RegisterMismatchCommon.jl index 346b436..815e12a 100644 --- a/src/RegisterMismatchCommon.jl +++ b/src/RegisterMismatchCommon.jl @@ -7,13 +7,13 @@ export correctbias!, nanpad, mismatch0, aperture_grid, allocate_mmarrays, defaul const DimsLike = Union{AbstractVector{Int}, Dims} const WidthLike = Union{AbstractVector,Tuple} -mismatch{T<:AbstractFloat}(fixed::AbstractArray{T}, moving::AbstractArray{T}, maxshift::DimsLike; normalization = :intensity) = mismatch(T, fixed, moving, maxshift; normalization=normalization) +mismatch(fixed::AbstractArray{T}, moving::AbstractArray{T}, maxshift::DimsLike; normalization = :intensity) where {T<:AbstractFloat} = mismatch(T, fixed, moving, maxshift; normalization=normalization) mismatch(fixed::AbstractArray, moving::AbstractArray, maxshift::DimsLike; normalization = :intensity) = mismatch(Float32, fixed, moving, maxshift; normalization=normalization) -mismatch_apertures{T<:AbstractFloat}(fixed::AbstractArray{T}, moving::AbstractArray{T}, args...; kwargs...) = mismatch_apertures(T, fixed, moving, args...; kwargs...) +mismatch_apertures(fixed::AbstractArray{T}, moving::AbstractArray{T}, args...; kwargs...) where {T<:AbstractFloat} = mismatch_apertures(T, fixed, moving, args...; kwargs...) mismatch_apertures(fixed::AbstractArray, moving::AbstractArray, args...; kwargs...) = mismatch_apertures(Float32, fixed, moving, args...; kwargs...) -function mismatch_apertures{T}(::Type{T}, fixed::AbstractArray, moving::AbstractArray, gridsize::DimsLike, maxshift::DimsLike; kwargs...) +function mismatch_apertures(::Type{T}, fixed::AbstractArray, moving::AbstractArray, gridsize::DimsLike, maxshift::DimsLike; kwargs...) where T cs = coords_spatial(fixed) aperture_centers = aperture_grid(size(fixed, cs...), gridsize) aperture_width = default_aperture_width(fixed, gridsize) @@ -31,12 +31,12 @@ whenever `i` or `j` is zero. Data are imputed by averaging the adjacent non-suspect values. This function works in-place, overwriting the original `mm`. """ -function correctbias!{ND,N}(mm::MismatchArray{ND,N}, w = correctbias_weight(mm)) +function correctbias!(mm::MismatchArray{ND,N}, w = correctbias_weight(mm)) where {ND,N} T = eltype(ND) mxshift = maxshift(mm) Imax = CartesianIndex(mxshift) Imin = CartesianIndex(map(x->-x,mxshift)::NTuple{N,Int}) - I1 = CartesianIndex(ntuple(d->d>2?0:1, N)::NTuple{N,Int}) # only first 2 dims + I1 = CartesianIndex(ntuple(d->d>2 ? 0 : 1, N)::NTuple{N,Int}) # only first 2 dims for I in eachindex(mm) if w[I] == 0 mms = NumDenom{T}(0,0) @@ -55,14 +55,14 @@ function correctbias!{ND,N}(mm::MismatchArray{ND,N}, w = correctbias_weight(mm)) end "`correctbias!(mms)` runs `correctbias!` on each element of an array-of-MismatchArrays." -function correctbias!{M<:MismatchArray}(mms::AbstractArray{M}) +function correctbias!(mms::AbstractArray{M}) where M<:MismatchArray for mm in mms correctbias!(mm) end mms end -function correctbias_weight{ND,N}(mm::MismatchArray{ND,N}) +function correctbias_weight(mm::MismatchArray{ND,N}) where {ND,N} T = eltype(ND) w = CenterIndexedArray(ones(T, size(mm))) for I in eachindex(mm) @@ -92,20 +92,20 @@ function nanpad(fixed, moving) get(fixed, rng, nanval(T)), get(moving, rng, nanval(T)) end -nanval{T<:AbstractFloat}(::Type{T}) = convert(T, NaN) -nanval{T}(::Type{T}) = convert(Float32, NaN) +nanval(::Type{T}) where {T<:AbstractFloat} = convert(T, NaN) +nanval(::Type{T}) where {T} = convert(Float32, NaN) """ `mm0 = mismatch0(fixed, moving, [normalization])` computes the "as-is" mismatch between `fixed` and `moving`, without any shift. `normalization` may be either `:intensity` (the default) or `:pixels`. """ -function mismatch0{Tf,Tm,N}(fixed::AbstractArray{Tf,N}, moving::AbstractArray{Tm,N}; normalization = :intensity) +function mismatch0(fixed::AbstractArray{Tf,N}, moving::AbstractArray{Tm,N}; normalization = :intensity) where {Tf,Tm,N} size(fixed) == size(moving) || throw(DimensionMismatch("Size $(size(fixed)) of fixed is not equal to size $(size(moving)) of moving")) _mismatch0(zero(Float64), zero(Float64), fixed, moving; normalization=normalization) end -function _mismatch0{T,Tf,Tm,N}(num::T, denom::T, fixed::AbstractArray{Tf,N}, moving::AbstractArray{Tm,N}; normalization = :intensity) +function _mismatch0(num::T, denom::T, fixed::AbstractArray{Tf,N}, moving::AbstractArray{Tm,N}; normalization = :intensity) where {T,Tf,Tm,N} if normalization == :intensity for i in eachindex(fixed, moving) vf = T(fixed[i]) @@ -136,7 +136,7 @@ mismatch between `fixed` and `moving`, without any shift. The mismatch is represented in `mms` as an aperture-wise Arrays-of-MismatchArrays. """ -function mismatch0{M<:MismatchArray}(mms::AbstractArray{M}) +function mismatch0(mms::AbstractArray{M}) where M<:MismatchArray mm0 = eltype(M)(0, 0) cr = eachindex(first(mms)) z = cr.start+cr.stop # all-zeros CartesianIndex @@ -152,7 +152,7 @@ grid of aperture centers. The grid has size `gridsize`, and is constructed for an image of spatial size `ssize`. Along each dimension the first and last elements are at the image corners. """ -function aperture_grid{N}(ssize::Dims{N}, gridsize) +function aperture_grid(ssize::Dims{N}, gridsize) where N if length(gridsize) != N if length(gridsize) == N-1 info("ssize and gridsize disagree; possible fix is to use a :time axis (AxisArrays) for the image") @@ -183,7 +183,7 @@ array-of-vectors) or an `N+1`-dimensional array with the center positions specified along the first dimension. (But you may find the `gridsize` syntax to be simpler.) """ -function allocate_mmarrays{T,C<:Union{AbstractVector,Tuple}}(::Type{T}, aperture_centers::AbstractArray{C}, maxshift) +function allocate_mmarrays(::Type{T}, aperture_centers::AbstractArray{C}, maxshift) where {T,C<:Union{AbstractVector,Tuple}} isempty(aperture_centers) && error("aperture_centers is empty") N = length(first(aperture_centers)) sz = map(x->2*x+1, maxshift) @@ -201,7 +201,7 @@ function allocate_mmarrays{T,C<:Union{AbstractVector,Tuple}}(::Type{T}, aperture mms end -function allocate_mmarrays{T,R<:Real}(::Type{T}, aperture_centers::AbstractArray{R}, maxshift) +function allocate_mmarrays(::Type{T}, aperture_centers::AbstractArray{R}, maxshift) where {T,R<:Real} N = ndims(aperture_centers)-1 mms = Array{MismatchArray{T,N}}(size(aperture_centers)[2:end]) sz = map(x->2*x+1, maxshift) @@ -211,7 +211,7 @@ function allocate_mmarrays{T,R<:Real}(::Type{T}, aperture_centers::AbstractArray mms end -function allocate_mmarrays{T<:Real,N}(::Type{T}, gridsize::NTuple{N,Int}, maxshift) +function allocate_mmarrays(::Type{T}, gridsize::NTuple{N,Int}, maxshift) where {T<:Real,N} mms = Array{MismatchArray{NumDenom{T},N}}(gridsize) sz = map(x->2*x+1, maxshift) for i in eachindex(mms) @@ -220,7 +220,7 @@ function allocate_mmarrays{T<:Real,N}(::Type{T}, gridsize::NTuple{N,Int}, maxshi mms end -immutable ContainerIterator{C} +struct ContainerIterator{C} data::C end @@ -228,11 +228,11 @@ Base.start(iter::ContainerIterator) = start(iter.data) Base.done(iter::ContainerIterator, state) = done(iter.data, state) Base.next(iter::ContainerIterator, state) = next(iter.data, state) -immutable FirstDimIterator{A<:AbstractArray,R<:CartesianRange} +struct FirstDimIterator{A<:AbstractArray,R<:CartesianRange} data::A rng::R - (::Type{FirstDimIterator{A,R}}){A,R}(data::A) = new{A,R}(data, CartesianRange(Base.tail(size(data)))) + FirstDimIterator{A,R}(data::A) where {A,R} = new{A,R}(data, CartesianRange(Base.tail(size(data)))) end FirstDimIterator(A::AbstractArray) = FirstDimIterator{typeof(A),typeof(CartesianRange(Base.tail(size(A))))}(A) @@ -250,9 +250,9 @@ AbstractArray-of-tuples or -AbstractVectors, or may be an `AbstractArray` where each point is represented along the first dimension (e.g., columns of a matrix). """ -each_point{C<:Union{AbstractVector,Tuple}}(aperture_centers::AbstractArray{C}) = ContainerIterator(aperture_centers) +each_point(aperture_centers::AbstractArray{C}) where {C<:Union{AbstractVector,Tuple}} = ContainerIterator(aperture_centers) -each_point{R<:Real}(aperture_centers::AbstractArray{R}) = FirstDimIterator(aperture_centers) +each_point(aperture_centers::AbstractArray{R}) where {R<:Real} = FirstDimIterator(aperture_centers) """ `rng = aperture_range(center, width)` returns a tuple of @@ -289,7 +289,7 @@ end `truncatenoise!(mm, thresh)` zeros out any entries of the MismatchArray `mm` whose `denom` values are less than `thresh`. """ -function truncatenoise!{T<:Real}(mm::AbstractArray{NumDenom{T}}, thresh::Real) +function truncatenoise!(mm::AbstractArray{NumDenom{T}}, thresh::Real) where T<:Real for I in eachindex(mm) if mm[I].denom <= thresh mm[I] = NumDenom{T}(0,0) @@ -298,7 +298,7 @@ function truncatenoise!{T<:Real}(mm::AbstractArray{NumDenom{T}}, thresh::Real) mm end -function truncatenoise!{A<:MismatchArray}(mms::AbstractArray{A}, thresh::Real) +function truncatenoise!(mms::AbstractArray{A}, thresh::Real) where A<:MismatchArray for i = 1:length(denoms) truncatenoise!(mms[i], thresh) end @@ -436,16 +436,9 @@ tovec(v::Tuple) = [v...] # TODO: redesign this whole thing to be safer? using Base: ViewIndex, to_indexes, unsafe_length, index_shape, tail -if VERSION < v"0.6.0-dev" - @inline function extraunsafe_view{T,N}(V::SubArray{T,N}, I::Vararg{ViewIndex,N}) - idxs = unsafe_reindex(V, V.indexes, to_indexes(I...)) - SubArray(V.parent, idxs, map(unsafe_length, (index_shape(V.parent, idxs...)))) - end -else - @inline function extraunsafe_view{T,N}(V::SubArray{T,N}, I::Vararg{ViewIndex,N}) - idxs = unsafe_reindex(V, V.indexes, to_indices(V, I)) - SubArray(V.parent, idxs) - end +@inline function extraunsafe_view(V::SubArray{T,N}, I::Vararg{ViewIndex,N}) where {T,N} + idxs = unsafe_reindex(V, V.indexes, to_indices(V, I)) + SubArray(V.parent, idxs) end unsafe_reindex(V, idxs::Tuple{UnitRange, Vararg{Any}}, subidxs::Tuple{UnitRange, Vararg{Any}}) = diff --git a/src/RegisterMismatchCuda.jl b/src/RegisterMismatchCuda.jl index 42440a7..01a28d6 100644 --- a/src/RegisterMismatchCuda.jl +++ b/src/RegisterMismatchCuda.jl @@ -61,7 +61,7 @@ close() = (for md in mdlist; unload(md); end; empty!(mdlist); empty!(ptxdict)) const FFTPROD = [2,3] -type NanCorrFFTs{T<:AbstractFloat,N} +mutable struct NanCorrFFTs{T<:AbstractFloat,N} I0::Tuple{CudaPitchedArray{T,N},CudaPitchedArray{Complex{T},N}} I1::Tuple{CudaPitchedArray{T,N},CudaPitchedArray{Complex{T},N}} I2::Tuple{CudaPitchedArray{T,N},CudaPitchedArray{Complex{T},N}} @@ -80,7 +80,7 @@ computations over domains of size `aperture_width`, computing the mismatch up to shifts of size `maxshift`. The keyword arguments allow you to control the planning process for the FFTs. """ -type CMStorage{T<:AbstractFloat,N} +mutable struct CMStorage{T<:AbstractFloat,N} aperture_width::Vector{Float64} maxshift::Vector{Int} getindexes::Vector{UnitRange{Int}} # indexes for pulling padded data, in source-coordinates @@ -97,7 +97,7 @@ type CMStorage{T<:AbstractFloat,N} fdshift::Vector{Int} # shift needed to unwrap (fftshift) stream - function CMStorage(::Type{T}, aperture_width::WidthLike, maxshift::DimsLike; stream=null_stream) + function CMStorage{T<:AbstractFloat,N}(::Type{T}, aperture_width::WidthLike, maxshift::DimsLike; stream=null_stream) where {T<:AbstractFloat,N} blocksize = map(x->ceil(Int,x), aperture_width) length(blocksize) == length(maxshift) || error("Dimensionality mismatch") padsz = padsize(blocksize, maxshift) @@ -119,7 +119,7 @@ type CMStorage{T<:AbstractFloat,N} end end # Note: display doesn't do anything -CMStorage{T<:Real}(::Type{T}, blocksize, maxshift; stream=null_stream, display=false) = CMStorage{T,length(blocksize)}(T, blocksize, maxshift; stream=stream) +CMStorage(::Type{T}, blocksize, maxshift; stream=null_stream, display=false) where {T<:Real} = CMStorage{T,length(blocksize)}(T, blocksize, maxshift; stream=stream) function free(cms::CMStorage) free(cms.fixed) @@ -130,8 +130,8 @@ end device(cms::CMStorage) = device(cms.num[1]) -eltype{T,N}(cms::CMStorage{T,N}) = T - ndims{T,N}(cms::CMStorage{T,N}) = N +eltype(cms::CMStorage{T,N}) where {T,N} = T + ndims(cms::CMStorage{T,N}) where {T,N} = N # Some tools from Images sdims(A::CudaPitchedArray) = ndims(A) @@ -152,7 +152,7 @@ normalization scheme (`:intensity` or `:pixels`). This operation is synchronous with respect to the host. """ -function mismatch{T<:Real}(::Type{T}, fixed::AbstractArray, moving::AbstractArray, maxshift::DimsLike; normalization = :intensity) +function mismatch(::Type{T}, fixed::AbstractArray, moving::AbstractArray, maxshift::DimsLike; normalization = :intensity) where T<:Real assertsamesize(fixed, moving) d_fixed = CudaPitchedArray(convert(Array{T}, fixed)) d_moving = CudaPitchedArray(convert(Array{T}, moving)) @@ -162,7 +162,7 @@ function mismatch{T<:Real}(::Type{T}, fixed::AbstractArray, moving::AbstractArra mm end -function mismatch{T}(fixed::AbstractCudaArray{T}, moving::AbstractCudaArray{T}, maxshift::DimsLike; normalization = :intensity) +function mismatch(fixed::AbstractCudaArray{T}, moving::AbstractCudaArray{T}, maxshift::DimsLike; normalization = :intensity) where T assertsamesize(fixed, moving) nd = ndims(fixed) maxshiftv = tovec(maxshift) @@ -199,13 +199,13 @@ in a rectangular grid, you can use an `N`-dimensional array-of-tuples (or array-of-vectors) or an `N+1`-dimensional array with the center positions specified along the first dimension. See `aperture_grid`. """ -function mismatch_apertures{T}(::Type{T}, - fixed::AbstractArray, - moving::AbstractArray, - aperture_centers::AbstractArray, - aperture_width::WidthLike, - maxshift::DimsLike; - kwargs...) +function mismatch_apertures(::Type{T}, + fixed::AbstractArray, + moving::AbstractArray, + aperture_centers::AbstractArray, + aperture_width::WidthLike, + maxshift::DimsLike; + kwargs...) where T assertsamesize(fixed, moving) d_fixed = CudaPitchedArray(convert(Array{T}, sdata(fixed))) d_moving = CudaPitchedArray(convert(Array{T}, moving)) @@ -217,13 +217,13 @@ end # only difference here relative to RegisterMismatch is the lack of the # FFTW keywords -function mismatch_apertures{T}(fixed::AbstractCudaArray{T}, - moving::AbstractCudaArray, - aperture_centers::AbstractArray, - aperture_width::WidthLike, - maxshift::DimsLike; - normalization = :pixels, - kwargs...) +function mismatch_apertures(fixed::AbstractCudaArray{T}, + moving::AbstractCudaArray, + aperture_centers::AbstractArray, + aperture_width::WidthLike, + maxshift::DimsLike; + normalization = :pixels, + kwargs...) where T nd = sdims(fixed) assertsamesize(fixed,moving) (length(aperture_width) == nd && length(maxshift) == nd) || error("Dimensionality mismatch") @@ -234,7 +234,7 @@ function mismatch_apertures{T}(fixed::AbstractCudaArray{T}, mms end -function fillfixed!{T}(cms::CMStorage{T}, fixed::CudaPitchedArray; f_indexes = ntuple(i->1:size(fixed,i), ndims(fixed))) +function fillfixed!(cms::CMStorage{T}, fixed::CudaPitchedArray; f_indexes = ntuple(i->1:size(fixed,i), ndims(fixed))) where T dev = device(cms) device(fixed) == dev || error("Fixed and cms must be on the same device") nd = ndims(cms) @@ -273,7 +273,7 @@ end computes the mismatch as a function of shift, storing the result in `mm`. The `fixed` image has been prepared in `cms`, a `CMStorage` object. """ -function mismatch!{T}(mm::MismatchArray, cms::CMStorage{T}, moving::CudaPitchedArray; normalization = :intensity, m_offset = ntuple(i->0, ndims(cms))) +function mismatch!(mm::MismatchArray, cms::CMStorage{T}, moving::CudaPitchedArray; normalization = :intensity, m_offset = ntuple(i->0, ndims(cms))) where T global ptxdict dev = device(cms) device(moving) == dev || error("Moving and cms must be on the same device") diff --git a/src/RegisterOptimize.jl b/src/RegisterOptimize.jl index b8ed1bb..40673d3 100644 --- a/src/RegisterOptimize.jl +++ b/src/RegisterOptimize.jl @@ -43,7 +43,7 @@ RegisterOptimize # Some conveniences for MathProgBase -@compat abstract type GradOnly <: MathProgBase.AbstractNLPEvaluator end +abstract type GradOnly <: MathProgBase.AbstractNLPEvaluator end function MathProgBase.initialize(d::GradOnly, requested_features::Vector{Symbol}) for feat in requested_features @@ -55,14 +55,14 @@ end MathProgBase.features_available(d::GradOnly) = [:Grad, :Jac] -@compat abstract type GradOnlyBoundsOnly <: GradOnly end +abstract type GradOnlyBoundsOnly <: GradOnly end MathProgBase.eval_g(::GradOnlyBoundsOnly, g, x) = nothing MathProgBase.jac_structure(::GradOnlyBoundsOnly) = Int[], Int[] MathProgBase.eval_jac_g(::GradOnlyBoundsOnly, J, x) = nothing -@compat abstract type BoundsOnly <: MathProgBase.AbstractNLPEvaluator end +abstract type BoundsOnly <: MathProgBase.AbstractNLPEvaluator end MathProgBase.eval_g(::BoundsOnly, g, x) = nothing MathProgBase.jac_structure(::BoundsOnly) = Int[], Int[] @@ -215,14 +215,14 @@ function p2rigid(p, SD) end to_float(A, B) = to_float(typeof(oneunit(eltype(A)) - oneunit(eltype(B))), A, B) -to_float{T<:AbstractFloat}(::Type{T}, A, B) = convert(Array{T}, A), convert(Array{T}, B) -to_float{T}(::Type{T}, A, B) = convert(Array{Float32}, A), convert(Array{Float32}, B) +to_float(::Type{T}, A, B) where {T<:AbstractFloat} = convert(Array{T}, A), convert(Array{T}, B) +to_float(::Type{T}, A, B) where {T} = convert(Array{Float32}, A), convert(Array{Float32}, B) ### ### Rigid registration from raw images, MathProg interface ### -type RigidValue{N,A<:AbstractArray,I<:AbstractExtrapolation,SDT} <: MathProgBase.AbstractNLPEvaluator +mutable struct RigidValue{N,A<:AbstractArray,I<:AbstractExtrapolation,SDT} <: MathProgBase.AbstractNLPEvaluator fixed::A wfixed::A moving::I @@ -230,7 +230,7 @@ type RigidValue{N,A<:AbstractArray,I<:AbstractExtrapolation,SDT} <: MathProgBase thresh end -function RigidValue{T<:Real}(fixed::AbstractArray, moving::AbstractArray{T}, SD, thresh) +function RigidValue(fixed::AbstractArray, moving::AbstractArray{T}, SD, thresh) where T<:Real f = copy(fixed) fnan = isnan.(f) f[fnan] = 0 @@ -253,7 +253,7 @@ function (d::RigidValue)(x) sum(abs2, f-m)/den end -type RigidOpt{RV<:RigidValue,G} <: GradOnlyBoundsOnly +mutable struct RigidOpt{RV<:RigidValue,G} <: GradOnlyBoundsOnly rv::RV g::G end @@ -291,7 +291,7 @@ function initial_deformation(ap::AffinePenalty, cs, Qs) _initial_deformation(ap, cs, Qs) end -function _initial_deformation{T,N}(ap::AffinePenalty{T,N}, cs, Qs) +function _initial_deformation(ap::AffinePenalty{T,N}, cs, Qs) where {T,N} if ap.λ <= 0 return cs2u(SVector{N,T}, cs), true end @@ -320,16 +320,16 @@ function _initial_deformation{T,N}(ap::AffinePenalty{T,N}, cs, Qs) convert_to_fixed(SVector{N,T}, x, size(cs)), isconverged end -cs2u{V}(::Type{V}, cs) = [V((c...)) for c in cs] +cs2u(::Type{V}, cs) where {V} = [V((c...)) for c in cs] -function initial_deformation{T,N,V<:SVector,M<:SMatrix}(ap::AffinePenalty{T,N}, cs::AbstractArray{V}, Qs::AbstractArray{M}) +function initial_deformation(ap::AffinePenalty{T,N}, cs::AbstractArray{V}, Qs::AbstractArray{M}) where {T,N,V<:SVector,M<:SMatrix} Tv = eltype(V) eltype(M) == Tv || error("element types of cs ($(eltype(V))) and Qs ($(eltype(M))) must match") size(M,1) == size(M,2) == length(V) || throw(DimensionMismatch("size $(size(M)) of Qs matrices is inconsistent with cs vectors of size $(size(V))")) _initial_deformation(convert(AffinePenalty{Tv,N}, ap), cs, Qs) end -function to_full{T,N}(ap::AffinePenalty{T,N}, Qs) +function to_full(ap::AffinePenalty{T,N}, Qs) where {T,N} FF = ap.F*ap.F' nA = N*size(FF,1) FFN = zeros(nA,nA) @@ -343,7 +343,7 @@ function to_full{T,N}(ap::AffinePenalty{T,N}, Qs) A end -function prep_b{T}(::Type{T}, cs, Qs) +function prep_b(::Type{T}, cs, Qs) where T n = prod(size(Qs)) N = length(first(cs))::Int b = zeros(T, N*n) @@ -369,39 +369,39 @@ function find_opt(P, b) end # A type for computing multiplication by the linear operator -type AffineQHessian{AP<:AffinePenalty,M<:StaticMatrix,N,Φ} +mutable struct AffineQHessian{AP<:AffinePenalty,M<:StaticMatrix,N,Φ} ap::AP Qs::Array{M,N} ϕ_old::Φ end -function AffineQHessian{T,TQ,N}(ap::AffinePenalty{T}, Qs::AbstractArray{TQ,N}, ϕ_old) +function AffineQHessian(ap::AffinePenalty{T}, Qs::AbstractArray{TQ,N}, ϕ_old) where {T,TQ,N} AffineQHessian{typeof(ap),similar_type(SArray,T,Size(N,N)),N,typeof(ϕ_old)}(ap, Qs, ϕ_old) end -Base.eltype{AP,M,N,Φ}(::Type{AffineQHessian{AP,M,N,Φ}}) = eltype(AP) +Base.eltype(::Type{AffineQHessian{AP,M,N,Φ}}) where {AP,M,N,Φ} = eltype(AP) Base.eltype(P::AffineQHessian) = eltype(typeof(P)) Base.size(P::AffineQHessian, d) = length(P.Qs)*size(first(P.Qs),1) # These compute the gradient of (x'*P*x)/2, where P is the Hessian # for the objective in the doc text for initial_deformation. -function (*){T,N}(P::AffineQHessian{AffinePenalty{T,N}}, x::AbstractVector{T}) +function (*)(P::AffineQHessian{AffinePenalty{T,N}}, x::AbstractVector{T}) where {T,N} u = convert_to_fixed(SVector{N,T}, x, size(P.Qs)) g = similar(u) _A_mul_B!(g, P, u) reinterpret(T, g, size(x)) end -function Base.LinAlg.A_mul_B!{T,N}(dest::AbstractVector{T}, - P::AffineQHessian{AffinePenalty{T,N}}, - x::AbstractVector{T}) +function Base.LinAlg.A_mul_B!(dest::AbstractVector{T}, + P::AffineQHessian{AffinePenalty{T,N}}, + x::AbstractVector{T}) where {T,N} u = convert_to_fixed(SVector{N,T}, x, size(P.Qs)) g = convert_to_fixed(SVector{N,T}, dest, size(P.Qs)) _A_mul_B!(g, P, u) dest end -function _A_mul_B!{T,N}(g, P::AffineQHessian{AffinePenalty{T,N}}, u) +function _A_mul_B!(g, P::AffineQHessian{AffinePenalty{T,N}}, u) where {T,N} gridsize = size(P.Qs) λ = P.ap.λ nspatialgrid = size(P.ap.F, 1) @@ -425,7 +425,7 @@ function _A_mul_B!{T,N}(g, P::AffineQHessian{AffinePenalty{T,N}}, u) end affine_part!(g, P, u) = _affine_part!(g, P.ap, u) -function _affine_part!{T,N}(g, ap::AffinePenalty{T,N}, u) +function _affine_part!(g, ap::AffinePenalty{T,N}, u) where {T,N} local s if ndims(u) == N s = penalty!(g, ap, u) @@ -448,7 +448,7 @@ function _affine_part!{T,N}(g, ap::AffinePenalty{T,N}, u) s end -function initial_deformation{T,N}(ap::AffinePenalty{T,N}, cs, Qs, ϕ_old, maxshift) +function initial_deformation(ap::AffinePenalty{T,N}, cs, Qs, ϕ_old, maxshift) where {T,N} error("This is broken, don't use it") b = prep_b(T, cs, Qs) # In case the grid is really big, solve iteratively @@ -463,12 +463,12 @@ end # type for minimization with composition (which turns the problem into # a nonlinear problem) -type InitialDefOpt{AQH,B} <: GradOnlyBoundsOnly +mutable struct InitialDefOpt{AQH,B} <: GradOnlyBoundsOnly P::AQH b::B end -function find_opt{AP,M,N,Φ<:GridDeformation}(P::AffineQHessian{AP,M,N,Φ}, b, maxshift, x0) +function find_opt(P::AffineQHessian{AP,M,N,Φ}, b, maxshift, x0) where {AP,M,N,Φ<:GridDeformation} objective = InitialDefOpt(P, b) solver = IpoptSolver(hessian_approximation="limited-memory", print_level=0, @@ -491,7 +491,7 @@ end MathProgBase.eval_f(d::InitialDefOpt, x::AbstractVector) = _eval_f(d.P, d.b, x) -function _eval_f{T,N}(P::AffineQHessian{AffinePenalty{T,N}}, b, x::AbstractVector) +function _eval_f(P::AffineQHessian{AffinePenalty{T,N}}, b, x::AbstractVector) where {T,N} gridsize = size(P.Qs) n = prod(gridsize) u = convert_to_fixed(x, (N,gridsize...))# reinterpret(SVector{N,T}, x, gridsize) @@ -511,12 +511,12 @@ function MathProgBase.eval_grad_f(d::InitialDefOpt, grad_f, x) copy!(grad_f, P*x-b) end -function affine_part!{AP,M,N,Φ<:GridDeformation}(g, P::AffineQHessian{AP,M,N,Φ}, u) +function affine_part!(g, P::AffineQHessian{AP,M,N,Φ}, u) where {AP,M,N,Φ<:GridDeformation} ϕ_c, g_c = compose(P.ϕ_old, GridDeformation(u, P.ϕ_old.knots)) penalty!(g, P.ap, ϕ_c, g_c) end -function affine_part!{AP,M,N,Φ<:GridDeformation}(::Void, P::AffineQHessian{AP,M,N,Φ}, u) +function affine_part!(::Void, P::AffineQHessian{AP,M,N,Φ}, u) where {AP,M,N,Φ<:GridDeformation} # Sadly, with GradientNumbers this gives an error I haven't traced # down (might be a Julia bug) # ϕ_c = P.ϕ_old(GridDeformation(u, P.ϕ_old.knots)) @@ -530,8 +530,8 @@ end ### itporder(mmis::Array) = itporder(eltype(mmis)) -itporder{T,N,A}(::Type{CenterIndexedArray{T,N,A}}) = itporder(A) -itporder{AI<:AbstractInterpolation}(::Type{AI}) = Interpolations.itptype(AI) +itporder(::Type{CenterIndexedArray{T,N,A}}) where {T,N,A} = itporder(A) +itporder(::Type{AI}) where {AI<:AbstractInterpolation} = Interpolations.itptype(AI) """ `ϕ, fval, fval0 = optimize!(ϕ, ϕ_old, dp, mmis; [tol=1e-6, print_level=0])` @@ -589,7 +589,7 @@ function _optimize!(ϕ, ϕ_old, dp::DeformationPenalty, mmis, T::Type; tol=1e-6, _optimize!(objective, ϕ, dp, mmis, tol, print_level; kwargs...) end -function optimize!{Tf<:Number, T, N}(ϕ, ϕ_old, dp::AffinePenalty{T,N}, mmis::Array{Tf}) +function optimize!(ϕ, ϕ_old, dp::AffinePenalty{T,N}, mmis::Array{Tf}) where {Tf<:Number, T, N} ND = NumDenom{Tf} mmisr = reinterpret(ND, mmis, tail(size(mmis))) mmisc = cachedinterpolators(mmisr, N, ntuple(d->(size(mmisr,d)+1)>>1, N)) @@ -622,12 +622,12 @@ function _optimize!(objective, ϕ, dp, mmis, tol, print_level; kwargs...) ϕ, fval, fval0 end -function u_as_vec{T,N}(ϕ::GridDeformation{T,N}) +function u_as_vec(ϕ::GridDeformation{T,N}) where {T,N} n = length(ϕ.u) uvec = reinterpret(T, ϕ.u, (n*N,)) end -function vec_as_u{T,N}(g::Array{T}, ϕ::GridDeformation{T,N}) +function vec_as_u(g::Array{T}, ϕ::GridDeformation{T,N}) where {T,N} reinterpret(SVector{N,T}, g, size(ϕ.u)) end @@ -637,7 +637,7 @@ function _copy!(ϕ::GridDeformation, x) end # Sequence of deformations -function u_as_vec{D<:GridDeformation}(ϕs::Vector{D}) +function u_as_vec(ϕs::Vector{D}) where D<:GridDeformation T = eltype(D) N = ndims(D) ngrid = length(first(ϕs).u) @@ -650,7 +650,7 @@ function u_as_vec{D<:GridDeformation}(ϕs::Vector{D}) end -type DeformOpt{D,Dold,DP,M} <: GradOnlyBoundsOnly +mutable struct DeformOpt{D,Dold,DP,M} <: GradOnlyBoundsOnly ϕ::D ϕ_old::Dold dp::DP @@ -697,14 +697,14 @@ function optimize!(ϕs, ϕs_old, dp::AffinePenalty, λt, mmis; kwargs...) _copy!(ϕs, result.minimum), result.f_minimum end -function optimize!{Tf<:Number, T, N}(ϕs, ϕs_old, dp::AffinePenalty{T,N}, λt, mmis::Array{Tf}) +function optimize!(ϕs, ϕs_old, dp::AffinePenalty{T,N}, λt, mmis::Array{Tf}) where {Tf<:Number, T, N} ND = NumDenom{Tf} mmisr = reinterpret(ND, mmis, tail(size(mmis))) mmisc = cachedinterpolators(mmisr, N, ntuple(d->(size(mmisr,d)+1)>>1, N)) optimize!(ϕs, ϕs_old, dp, λt, mmisc) end -type DeformTseriesOpt{D,Dsold,DP,T,M} <: GradOnlyBoundsOnly +mutable struct DeformTseriesOpt{D,Dsold,DP,T,M} <: GradOnlyBoundsOnly ϕs::Vector{D} ϕs_old::Dsold dp::DP @@ -724,7 +724,7 @@ function MathProgBase.eval_grad_f(d::DeformTseriesOpt, grad_f, x) penalty!(grad_f, d.ϕs, d.ϕs_old, d.dp, d.λt, d.mmis) end -function _copy!{D<:GridDeformation,T<:Number}(ϕs::Vector{D}, x::Array{T}) +function _copy!(ϕs::Vector{D}, x::Array{T}) where {D<:GridDeformation,T<:Number} N = ndims(first(ϕs)) L = N*length(first(ϕs).u) length(x) == L*length(ϕs) || throw(DimensionMismatch("ϕs is incommensurate with a vector of length $(length(x))")) @@ -745,7 +745,7 @@ object for that grid, and `mmis` is the array-of-mismatch arrays See also: `auto_λ`. """ -function fixed_λ{T,N}(cs, Qs, knots::NTuple{N}, ap::AffinePenalty{T,N}, mmis; ϕs_old = identity, mu_init=0.1, kwargs...) +function fixed_λ(cs, Qs, knots::NTuple{N}, ap::AffinePenalty{T,N}, mmis; ϕs_old = identity, mu_init=0.1, kwargs...) where {T,N} maxshift = map(x->(x-1)>>1, size(first(mmis))) u0, isconverged = initial_deformation(ap, cs, Qs) if !isconverged @@ -771,7 +771,7 @@ end computes an optimal vector-of-deformations `ϕs` for an image sequence, using an temporal penalty coefficient `λt`. """ -function fixed_λ{T,N,TP,L}(cs::AbstractArray{SVector{N,T}}, Qs::AbstractArray{SMatrix{N,N,T,L}}, knots::NTuple{N}, ap::AffinePenalty{TP,N}, λt, mmis; ϕs_old = identity, mu_init=0.1, kwargs...) +function fixed_λ(cs::AbstractArray{SVector{N,T}}, Qs::AbstractArray{SMatrix{N,N,T,L}}, knots::NTuple{N}, ap::AffinePenalty{TP,N}, λt, mmis; ϕs_old = identity, mu_init=0.1, kwargs...) where {T,N,TP,L} λtT = T(λt) apT = convert(AffinePenalty{T,N}, ap) maxshift = map(x->(x-1)>>1, size(first(mmis))) @@ -790,7 +790,7 @@ function fixed_λ{T,N,TP,L}(cs::AbstractArray{SVector{N,T}}, Qs::AbstractArray{S end # This version re-packs variables as read from the .jld file -function fixed_λ{Tf<:Number,T,N}(cs::Array{Tf}, Qs::Array{Tf}, knots::NTuple{N}, ap::AffinePenalty{T,N}, λt, mmis::Array{Tf}; kwargs...) +function fixed_λ(cs::Array{Tf}, Qs::Array{Tf}, knots::NTuple{N}, ap::AffinePenalty{T,N}, λt, mmis::Array{Tf}; kwargs...) where {Tf<:Number,T,N} csr = reinterpret(SVector{N,Tf}, cs, tail(size(cs))) Qsr = reinterpret(similar_type(SArray,Tf,Size(N,N)), Qs, tail(tail(size(Qs)))) if length(mmis) > 10^7 @@ -841,7 +841,7 @@ knots, mmis, (λmin, λmax))` will perform the analysis on the chosen See also: `fixed_λ`. Because `auto_λ` performs the optimization repeatedly for many different `λ`s, it is slower than `fixed_λ`. """ -function auto_λ{R<:Real,S<:Real,N}(fixed::AbstractArray{R}, moving::AbstractArray{S}, gridsize::NTuple{N}, maxshift::NTuple{N}, λrange; thresh=(0.5)^ndims(fixed)*length(fixed)/prod(gridsize), kwargs...) +function auto_λ(fixed::AbstractArray{R}, moving::AbstractArray{S}, gridsize::NTuple{N}, maxshift::NTuple{N}, λrange; thresh=(0.5)^ndims(fixed)*length(fixed)/prod(gridsize), kwargs...) where {R<:Real,S<:Real,N} T = Float64 local mms try @@ -855,24 +855,24 @@ function auto_λ{R<:Real,S<:Real,N}(fixed::AbstractArray{R}, moving::AbstractArr auto_λ(cs, Qs, knots, mmis, λrange; kwargs...) end -function auto_λ{N}(cs, Qs, knots::NTuple{N}, mmis, λrange; kwargs...) +function auto_λ(cs, Qs, knots::NTuple{N}, mmis, λrange; kwargs...) where N ap = AffinePenalty{Float64,N}(knots, λrange[1]) # default to affine-residual penalty, Ipopt needs Float64 auto_λ(cs, Qs, knots, ap, mmis, λrange; kwargs...) end -function auto_λ{N}(stackidx::Integer, cs, Qs, knots::NTuple{N}, mmis, λrange; kwargs...) +function auto_λ(stackidx::Integer, cs, Qs, knots::NTuple{N}, mmis, λrange; kwargs...) where N cs1 = cs[ntuple(d->Colon(),ndims(cs)-1)..., stackidx]; Qs1 = Qs[ntuple(d->Colon(),ndims(Qs)-1)..., stackidx]; mmis1 = mmis[ntuple(d->Colon(),ndims(mmis)-1)..., stackidx]; auto_λ(cs1, Qs1, knots, mmis1, λrange; kwargs...) end -function auto_λ{Tf<:Number,T,N}(cs::Array{Tf}, Qs::Array{Tf}, knots::NTuple{N}, ap::AffinePenalty{T,N}, mmis::Array{Tf}, λrange; kwargs...) +function auto_λ(cs::Array{Tf}, Qs::Array{Tf}, knots::NTuple{N}, ap::AffinePenalty{T,N}, mmis::Array{Tf}, λrange; kwargs...) where {Tf<:Number,T,N} # Ipopt requires Float64 auto_λ(convert(Array{Float64}, cs), convert(Array{Float64}, Qs), knots, ap, convert(Array{Float64}, mmis), λrange; kwargs...) end -function auto_λ{T,N}(cs::Array{Float64}, Qs::Array{Float64}, knots::NTuple{N}, ap::AffinePenalty{T,N}, mmis::Array{Float64}, λrange; kwargs...) +function auto_λ(cs::Array{Float64}, Qs::Array{Float64}, knots::NTuple{N}, ap::AffinePenalty{T,N}, mmis::Array{Float64}, λrange; kwargs...) where {T,N} csr = reinterpret(SVector{N,Float64}, cs, tail(size(cs))) Qsr = reinterpret(similar_type(SArray,Float64,Size(N,N)), Qs, tail(tail(size(Qs)))) mmisr = reinterpret(NumDenom{Float64}, mmis, tail(size(mmis))) @@ -881,7 +881,7 @@ function auto_λ{T,N}(cs::Array{Float64}, Qs::Array{Float64}, knots::NTuple{N}, auto_λ(csr, Qsr, knots, ap64, mmisc, λrange; kwargs...) end -function auto_λ{T,N}(cs, Qs, knots::NTuple{N}, ap::AffinePenalty{T,N}, mmis, λrange; kwargs...) +function auto_λ(cs, Qs, knots::NTuple{N}, ap::AffinePenalty{T,N}, mmis, λrange; kwargs...) where {T,N} λmin, λmax = λrange gridsize = map(length, knots) uc = zeros(T, N, gridsize...) @@ -997,7 +997,7 @@ function auto_λt(Es, cs, Qs, ap, λtrange) λts, datapenalty end -function auto_λt{Tf<:Number,T,N}(Es, cs::Array{Tf}, Qs::Array{Tf}, ap::AffinePenalty{T,N}, λt) +function auto_λt(Es, cs::Array{Tf}, Qs::Array{Tf}, ap::AffinePenalty{T,N}, λt) where {Tf<:Number,T,N} csr = reinterpret(SVector{N,Tf}, cs, tail(size(cs))) Qsr = reinterpret(similar_type(SArray,Tf,Size(N,N)), Qs, tail(tail(size(Qs)))) auto_λt(Es, csr, Qsr, ap, λt) @@ -1006,7 +1006,7 @@ end ### ### Whole-experiment optimization with a temporal roughness penalty ### -function initial_deformation{T,N,V<:SVector,M<:SMatrix}(ap::AffinePenalty{T,N}, λt, cs::AbstractArray{V}, Qs::AbstractArray{M}) +function initial_deformation(ap::AffinePenalty{T,N}, λt, cs::AbstractArray{V}, Qs::AbstractArray{M}) where {T,N,V<:SVector,M<:SMatrix} Tv = eltype(V) eltype(M) == Tv || error("element types of cs ($(eltype(V))) and Qs ($(eltype(M))) must match") length(V) == N || throw(DimensionMismatch("Dimensionality $N of ap does not match $(length(V))")) @@ -1018,37 +1018,37 @@ function initial_deformation{T,N,V<:SVector,M<:SMatrix}(ap::AffinePenalty{T,N}, convert_to_fixed(SVector{N,Tv}, x, size(cs)), isconverged end -immutable TimeHessian{AQH<:AffineQHessian,T} +struct TimeHessian{AQH<:AffineQHessian,T} aqh::AQH λt::T end -Base.eltype{AQH,T}(::Type{TimeHessian{AQH,T}}) = eltype(AQH) +Base.eltype(::Type{TimeHessian{AQH,T}}) where {AQH,T} = eltype(AQH) Base.eltype(P::TimeHessian) = eltype(typeof(P)) Base.size(P::TimeHessian, d) = size(P.aqh, d) -function (*){AQH}(P::TimeHessian{AQH}, x::AbstractVector) +function (*)(P::TimeHessian{AQH}, x::AbstractVector) where AQH y = P.aqh*x ϕs = vec2vecϕ(P.aqh.Qs, x) penalty!(y, P.λt, ϕs) y end -function Base.LinAlg.A_mul_B!{AQH}(y::AbstractVector, - P::TimeHessian{AQH}, - x::AbstractVector) +function Base.LinAlg.A_mul_B!(y::AbstractVector, + P::TimeHessian{AQH}, + x::AbstractVector) where AQH A_mul_B!(y, P.aqh, x) ϕs = vec2vecϕ(P.aqh.Qs, x) penalty!(y, P.λt, ϕs) y end -function vec2vecϕ{T,N,L}(Qs::Array{SMatrix{N,N,T,L}}, x::AbstractVector{T}) +function vec2vecϕ(Qs::Array{SMatrix{N,N,T,L}}, x::AbstractVector{T}) where {T,N,L} xf = convert_to_fixed(SVector{N,T}, x, size(Qs)) _vec2vecϕ(xf, Base.front(size(Qs))) end -@noinline function _vec2vecϕ{N}(x::AbstractArray, sz::NTuple{N,Int}) +@noinline function _vec2vecϕ(x::AbstractArray, sz::NTuple{N,Int}) where N colons = ntuple(ColonFun, Val{N}) [GridDeformation(view(x, colons..., i), sz) for i = 1:size(x)[end]] end @@ -1124,7 +1124,7 @@ function optimize(tform::AffineTransform, mmis, knots) AffineTransform(convert(Matrix{T}, Siopt), convert(Vector{T}, displacementopt)), result.f_minimum end -function affinepenalty!{N}(g, At, mmis, keep, Xt, gridsize::NTuple{N}, gtmp) +function affinepenalty!(g, At, mmis, keep, Xt, gridsize::NTuple{N}, gtmp) where N u = _calculate_u(At, Xt, gridsize) @assert eltype(u) == eltype(At) val = penalty!(gtmp, u, mmis, keep) @@ -1137,7 +1137,7 @@ function affinepenalty!{N}(g, At, mmis, keep, Xt, gridsize::NTuple{N}, gtmp) val end -function _calculate_u{N}(At, Xt, gridsize::NTuple{N}) +function _calculate_u(At, Xt, gridsize::NTuple{N}) where N Ut = Xt*At u = Ut[:,1:size(Ut,2)-1]' # discard the dummy dimension reinterpret(SVector{N, eltype(u)}, u, gridsize) # put u in the shape of the grid @@ -1188,7 +1188,7 @@ function fit_sigmoid(data) end -type SigmoidOpt{G,H} <: BoundsOnly +mutable struct SigmoidOpt{G,H} <: BoundsOnly data::Vector{Float64} g::G h::H @@ -1228,7 +1228,7 @@ function sigpenalty(x, data) sum(abs2, (data-bottom)/(top-bottom) - 1./(1+exp.(-(collect(1:length(data))-center)/width))) end -@generated function RegisterCore.maxshift{T,N}(A::CachedInterpolation{T,N}) +@generated function RegisterCore.maxshift(A::CachedInterpolation{T,N}) where {T,N} args = [:(size(A.parent, $d)>>1) for d = 1:N] Expr(:tuple, args...) end diff --git a/src/RegisterPenalty.jl b/src/RegisterPenalty.jl index 0aaf606..e03380a 100644 --- a/src/RegisterPenalty.jl +++ b/src/RegisterPenalty.jl @@ -31,12 +31,12 @@ The main exported types/functions are: RegisterPenalty -@compat abstract type DeformationPenalty{T,N} end -Base.eltype{T,N}(::Type{DeformationPenalty{T,N}}) = T -Base.eltype{DP<:DeformationPenalty}(::Type{DP}) = eltype(supertype(DP)) +abstract type DeformationPenalty{T,N} end +Base.eltype(::Type{DeformationPenalty{T,N}}) where {T,N} = T +Base.eltype(::Type{DP}) where {DP<:DeformationPenalty} = eltype(supertype(DP)) Base.eltype(dp::DeformationPenalty) = eltype(typeof(dp)) -Base.ndims{T,N}(::Type{DeformationPenalty{T,N}}) = N -Base.ndims{DP<:DeformationPenalty}(::Type{DP}) = ndims(supertype(DP)) +Base.ndims(::Type{DeformationPenalty{T,N}}) where {T,N} = N +Base.ndims(::Type{DP}) where {DP<:DeformationPenalty} = ndims(supertype(DP)) Base.ndims(dp::DeformationPenalty) = ndims(typeof(dp)) """ @@ -96,7 +96,7 @@ function penalty!(g, ϕ, ϕ_old, dp::DeformationPenalty, mmis::AbstractArray, ke end # Allow it to be called without StaticArrays -function penalty!{T<:Number}(g::Array{T}, ϕ, ϕ_old, dp::DeformationPenalty, mmis::AbstractArray, keep = trues(size(mmis))) +function penalty!(g::Array{T}, ϕ, ϕ_old, dp::DeformationPenalty, mmis::AbstractArray, keep = trues(size(mmis))) where T<:Number gf = RegisterDeformation.convert_to_fixed(g, (ndims(dp), size(ϕ.u)...)) penalty!(gf, ϕ, ϕ_old, dp, mmis, keep) end @@ -108,7 +108,7 @@ evaluates the total penalty for temporal sequence of deformations `ϕs`, using the temporal sequence of mismatch data `mmis`. `λt` is the temporal roughness penalty coefficient. """ -function penalty!{D<:AbstractDeformation}(gs, ϕs::AbstractVector{D}, ϕs_old, dp::DeformationPenalty, λt::Number, mmis::AbstractArray, keep = trues(size(mmis))) +function penalty!(gs, ϕs::AbstractVector{D}, ϕs_old, dp::DeformationPenalty, λt::Number, mmis::AbstractArray, keep = trues(size(mmis))) where D<:AbstractDeformation ntimes = length(ϕs) size(mmis)[end] == ntimes || throw(DimensionMismatch("Number of deformations $ntimes does not agree with mismatch data of size $(size(mmis))")) s = _penalty!(gs, ϕs, ϕs_old, dp, mmis, keep, 1) @@ -120,12 +120,12 @@ function penalty!{D<:AbstractDeformation}(gs, ϕs::AbstractVector{D}, ϕs_old, d end -function penalty!{T<:Number,D<:AbstractDeformation}(gs::Array{T}, ϕs::AbstractVector{D}, ϕs_old, dp::DeformationPenalty, λt::Number, mmis::AbstractArray, keep = trues(size(mmis))) +function penalty!(gs::Array{T}, ϕs::AbstractVector{D}, ϕs_old, dp::DeformationPenalty, λt::Number, mmis::AbstractArray, keep = trues(size(mmis))) where {T<:Number,D<:AbstractDeformation} gf = RegisterDeformation.convert_to_fixed(gs, (ndims(dp), size(first(ϕs).u)..., length(ϕs))) penalty!(gf, ϕs, ϕs_old, dp, λt, mmis, keep) end -function _penalty!{T,N}(gs, ϕs, ϕs_old, dp::DeformationPenalty{T,N}, mmis, keeps, i) +function _penalty!(gs, ϕs, ϕs_old, dp::DeformationPenalty{T,N}, mmis, keeps, i) where {T,N} colons = ntuple(ColonFun, Val{N}) indexes = (colons..., i) #mmi = view(mmis, indexes...) # making these views runs afoul of inference limits @@ -146,7 +146,7 @@ end # Data penalty # ################ -@compat const CenteredInterpolant{T,N,A<:AbstractInterpolation} = Union{MismatchArray{T,N,A}, CachedInterpolation{T,N}} +const CenteredInterpolant{T,N,A<:AbstractInterpolation} = Union{MismatchArray{T,N,A}, CachedInterpolation{T,N}} """ `p = penalty!(g, ϕ, mmis, [keep=trues(size(mmis))])` computes the @@ -172,16 +172,16 @@ where each index `_i` refers to a single aperture, and each `p_i` involves just `mmis[i]` and `ϕ.u[:,i]`. `mmis[i]` must be interpolating, so that it can be evaluated for fractional shifts. """ -function penalty!{M<:CenteredInterpolant}(g, ϕ::AbstractDeformation, mmis::AbstractArray{M}, keep=trues(size(mmis))) +function penalty!(g, ϕ::AbstractDeformation, mmis::AbstractArray{M}, keep=trues(size(mmis))) where M<:CenteredInterpolant penalty!(g, ϕ.u, mmis, keep) end -function penalty!{Tg<:Number,M<:CenteredInterpolant}(g::Array{Tg}, ϕ::AbstractDeformation, mmis::AbstractArray{M}, keep=trues(size(mmis))) +function penalty!(g::Array{Tg}, ϕ::AbstractDeformation, mmis::AbstractArray{M}, keep=trues(size(mmis))) where {Tg<:Number,M<:CenteredInterpolant} gf = RegisterDeformation.convert_to_fixed(g, (ndims(ϕ), size(ϕ.u)...)) penalty!(gf, ϕ, mmis, keep) end -function penalty!{Tu,Dim,M<:CenteredInterpolant}(g, u::AbstractArray{SVector{Dim,Tu}}, mmis::AbstractArray{M}, keep=trues(size(mmis))) +function penalty!(g, u::AbstractArray{SVector{Dim,Tu}}, mmis::AbstractArray{M}, keep=trues(size(mmis))) where {Tu,Dim,M<:CenteredInterpolant} # This "outer" function just handles the chain rule for computing the # total penalty and gradient. The "real" work is done by penalty_nd!. nblocks = length(mmis) @@ -247,14 +247,14 @@ end penalty_nd!(gnd, u::AbstractInterpolation, mmis, keep) = error("ϕ must not be interpolating") -@generated function checkbounds_shift{N}(dx::SVector{N}, mxs) +@generated function checkbounds_shift(dx::SVector{N}, mxs) where N quote @nexprs $N d->(if abs(dx[d]) >= mxs[d]-0.5 return false end) true end end -@generated function _wsum{N}(x::SVector{N}, cnum, cdenom) +@generated function _wsum(x::SVector{N}, cnum, cdenom) where N args = [:(cnum*x[$d].num + cdenom*x[$d].denom) for d = 1:N] quote SVector($(args...)) @@ -280,13 +280,13 @@ When the deformation is defined on a regular grid, `knots` can be an NTuple of knot-vectors; otherwise, it should be an `ndims`-by-`npoints` matrix that stores the knot locations in columns. """ -type AffinePenalty{T,N} <: DeformationPenalty{T,N} +mutable struct AffinePenalty{T,N} <: DeformationPenalty{T,N} F::Matrix{T} # geometry data for the affine-residual penalty λ::T # regularization coefficient - (::Type{AffinePenalty{T,N}}){T,N}(F::Matrix{T}, λ::T, _) = new{T,N}(F, λ) + AffinePenalty{T,N}(F::Matrix{T}, λ::T, _) where {T,N} = new{T,N}(F, λ) - function (::Type{AffinePenalty{T,N}}){T,N}(knots::NTuple{N}, λ) + function AffinePenalty{T,N}(knots::NTuple{N}, λ) where {T,N} gridsize = map(length, knots) C = Matrix{Float64}(prod(gridsize), N+1) i = 0 @@ -300,18 +300,18 @@ type AffinePenalty{T,N} <: DeformationPenalty{T,N} new{T,N}(F, λ) end - function (::Type{AffinePenalty{T,N}}){T,N}(knots::AbstractMatrix, λ) + function AffinePenalty{T,N}(knots::AbstractMatrix, λ) where {T,N} C = hcat(knots', ones(eltype(knots), size(knots, 2))) F, _ = qr(C) new{T,N}(F, λ) end end -AffinePenalty{V<:AbstractVector,N}(knots::NTuple{N,V}, λ) = AffinePenalty{eltype(V),N}(knots, λ) -AffinePenalty{V<:AbstractVector}(knots::AbstractVector{V}, λ) = AffinePenalty{eltype(V),length(knots)}((knots...), λ) -AffinePenalty{T}(knots::AbstractMatrix{T}, λ) = AffinePenalty{T,size(knots,1)}(knots, λ) +AffinePenalty(knots::NTuple{N,V}, λ) where {V<:AbstractVector,N} = AffinePenalty{eltype(V),N}(knots, λ) +AffinePenalty(knots::AbstractVector{V}, λ) where {V<:AbstractVector} = AffinePenalty{eltype(V),length(knots)}((knots...), λ) +AffinePenalty(knots::AbstractMatrix{T}, λ) where {T} = AffinePenalty{T,size(knots,1)}(knots, λ) -Base.convert{T,N}(::Type{AffinePenalty{T,N}}, ap::AffinePenalty) = AffinePenalty{T,N}(convert(Matrix{T}, ap.F), convert(T, ap.λ), 0) +Base.convert(::Type{AffinePenalty{T,N}}, ap::AffinePenalty) where {T,N} = AffinePenalty{T,N}(convert(Matrix{T}, ap.F), convert(T, ap.λ), 0) """ `p = penalty!(g, dp::DeformationPenalty, ϕ_c, [g_c])` computes the @@ -328,11 +328,11 @@ If `g` is non-empty, the gradient of the penalty with respect to `u` will be computed. If you write a custom `penalty!` function for a new `DeformationPenalty`, it is your responsibility to set `g` in-place. """ -function penalty!{T,N}(g, dp::AffinePenalty, ϕ_c::AbstractDeformation{T,N}) +function penalty!(g, dp::AffinePenalty, ϕ_c::AbstractDeformation{T,N}) where {T,N} penalty!(g, dp, ϕ_c.u) end -function penalty!{T,N}(g, dp::AffinePenalty, u::AbstractArray{SVector{N,T},N}) +function penalty!(g, dp::AffinePenalty, u::AbstractArray{SVector{N,T},N}) where {T,N} F, λ = dp.F, dp.λ if λ == 0 if g != nothing && !isempty(g) @@ -377,7 +377,7 @@ for a vector `ϕ` of deformations. `g`, if not `nothing`, should be a single real-valued vector with number of entries corresponding to all of the `u` arrays in all of `ϕs`. """ -function penalty!{D<:GridDeformation}(g, λt::Real, ϕs::Vector{D}) +function penalty!(g, λt::Real, ϕs::Vector{D}) where D<:GridDeformation if g == nothing || isempty(g) return penalty(λt, ϕs) end @@ -403,14 +403,14 @@ function penalty!{D<:GridDeformation}(g, λt::Real, ϕs::Vector{D}) s end -function penalty!{T<:Number, D<:GridDeformation}(g::Array{T}, λt::Real, ϕs::Vector{D}) +function penalty!(g::Array{T}, λt::Real, ϕs::Vector{D}) where {T<:Number, D<:GridDeformation} N = ndims(first(ϕs)) sz = size(first(ϕs).u) gf = RegisterDeformation.convert_to_fixed(g, (N, sz..., length(ϕs))) penalty!(gf, λt, ϕs) end -function penalty{D<:GridDeformation}(λt::Real, ϕs::Vector{D}) +function penalty(λt::Real, ϕs::Vector{D}) where D<:GridDeformation s = zero(eltype(D)) ngrid = length(first(ϕs).u) λt2 = λt/2 @@ -425,11 +425,11 @@ function penalty{D<:GridDeformation}(λt::Real, ϕs::Vector{D}) s end -function RegisterDeformation.convert_to_fixed{T,N,A,L}(::Type{GridDeformation{T,N,A,L}}, g::Array{T}) +function RegisterDeformation.convert_to_fixed(::Type{GridDeformation{T,N,A,L}}, g::Array{T}) where {T,N,A,L} reinterpret(SVector{N,T}, g, (div(length(g), N),)) end -function vec2ϕs{T,N}(x::Array{T}, gridsize::NTuple{N,Int}, n, knots) +function vec2ϕs(x::Array{T}, gridsize::NTuple{N,Int}, n, knots) where {T,N} xr = RegisterDeformation.convert_to_fixed(SVector{N,T}, x, (gridsize..., n)) colons = ntuple(d->Colon(), N)::NTuple{N,Colon} [GridDeformation(view(xr, colons..., i), knots) for i = 1:n] @@ -448,7 +448,7 @@ you want to specify the boundary condition (default is `InPlace()`). `mmis = interpolate_mm!(mms, order)` prepares the array-of-MismatchArrays `mms` for interpolation. """ -function interpolate_mm!{T<:MismatchArray}(mms::AbstractArray{T}, spline::Interpolations.Degree) +function interpolate_mm!(mms::AbstractArray{T}, spline::Interpolations.Degree) where T<:MismatchArray f = x->CenterIndexedArray(interpolate!(x.data, BSpline(spline), gridstyle(spline))) map(f, mms) end @@ -457,7 +457,7 @@ function interpolate_mm!(mm::MismatchArray, spline::Interpolations.Degree) CenterIndexedArray(interpolate!(mm.data, BSpline(spline), gridstyle(spline))) end -interpolate_mm!{order<:Interpolations.Degree}(arg, ::Type{order}) = +interpolate_mm!(arg, ::Type{order}) where {order<:Interpolations.Degree} = interpolate_mm!(arg, itporder(order)) interpolate_mm!(arg) = interpolate_mm!(arg, Quadratic) @@ -467,7 +467,7 @@ itporder(::Type{Linear}) = Linear() gridstyle(::Quadratic) = OnCell() gridstyle(::Linear) = OnGrid() -@generated function Interpolations.gradient!{T,N}(g::AbstractVector, A::CenterIndexedArray{T,N}, i::Number...) +@generated function Interpolations.gradient!(g::AbstractVector, A::CenterIndexedArray{T,N}, i::Number...) where {T,N} length(i) == N || error("Must use $N indexes") args = [:(i[$d]+A.halfsize[$d]+1) for d = 1:N] meta = Expr(:meta, :inline) diff --git a/src/RegisterUtilities.jl b/src/RegisterUtilities.jl index 736803b..54e6a2a 100644 --- a/src/RegisterUtilities.jl +++ b/src/RegisterUtilities.jl @@ -9,7 +9,7 @@ export Counter # # Stolen from Grid.jl. Useful when you want to do more math on the iterator. -immutable Counter +struct Counter max::Vector{Int} end Counter(sz::Tuple) = Counter(Int[sz...]) diff --git a/test/register_hindsight.jl b/test/register_hindsight.jl index 6acac9e..8561240 100644 --- a/test/register_hindsight.jl +++ b/test/register_hindsight.jl @@ -6,7 +6,7 @@ using TestImages using Base.Test #add jitter in sampling location, simulating inconsistencies in piezo position when using OCPI under certain conditions -function jitter{T}(img::Array{T,1}, npix::Float64) +function jitter(img::Array{T,1}, npix::Float64) where T etp = extrapolate(interpolate(img, BSpline(Linear()),OnGrid()), Flat()) out = zeros(eltype(img), size(img)) z_def = Float64[] @@ -29,7 +29,7 @@ function dualgrad_data!(g, ϕ, fixed, moving) for i in CartesianRange(indices(ϕ.u.itp.coefs)) for j = 1:nd temp = ur[j, i] - ur[j, i] = dual(value(temp), 1.0) + ur[j, i] = dual(DualNumbers.value(temp), 1.0) gr[j, i] = epsilon(RegisterHindsight.penalty_hindsight_data(ϕ, fixed, moving)) ur[j, i] = temp end @@ -43,7 +43,7 @@ function dualgrad_reg!(g, ap, ϕ) for i in CartesianRange(indices(ϕ.u.itp.coefs)) for j = 1:nd temp = ur[j, i] - ur[j, i] = dual(value(temp), 1.0) + ur[j, i] = dual(DualNumbers.value(temp), 1.0) gr[j, i] = epsilon(RegisterHindsight.penalty_hindsight_reg(ap, ϕ)) ur[j, i] = temp end diff --git a/test/register_optimize.jl b/test/register_optimize.jl index 4c2d102..f1626d8 100644 --- a/test/register_optimize.jl +++ b/test/register_optimize.jl @@ -36,7 +36,7 @@ function initial_guess_direct(A, cs::Matrix, Qs::Matrix) reinterpret(SVector{2,Float64}, x, size(Qs)) end -function build_Ac_b{Tc,TQ}(A, λt, cs::Array{Tc,3}, Qs::Array{TQ,3}) +function build_Ac_b(A, λt, cs::Array{Tc,3}, Qs::Array{TQ,3}) where {Tc,TQ} n = size(Qs,3) l = size(A,1) b = zeros(l*n) @@ -66,7 +66,7 @@ function build_Ac_b{Tc,TQ}(A, λt, cs::Array{Tc,3}, Qs::Array{TQ,3}) Ac, b end -function initial_guess_direct{Tc,TQ}(A, λt, cs::Array{Tc,3}, Qs::Array{TQ,3}) +function initial_guess_direct(A, λt, cs::Array{Tc,3}, Qs::Array{TQ,3}) where {Tc,TQ} Ac, b = build_Ac_b(A, λt, cs, Qs) x = Ac\b reinterpret(SVector{2,Float64}, x, size(Qs)) @@ -361,7 +361,7 @@ mxrot = pi/90 minwidth_rot = [0.0002] SD = eye(ndims(fixed)) -tfm, mm = qd_rigid(fixed, moving, mxshift, mxrot, minwidth_rot, SD; thresh=thresh, rtol = 1e-7, atol = 1e-9 * thresh) +tfm, mm = qd_rigid(fixed, moving, mxshift, mxrot, minwidth_rot, SD; thresh=thresh, maxevals=1000, rtol = 1e-7, atol = 1e-9 * thresh) @test sum(abs.(tfm0.m - tfm.m)) < 1e-3 @@ -378,9 +378,9 @@ mxrot = [pi/90; pi/90; pi/90] minwidth_rot = fill(0.0002, 3) SD = eye(ndims(fixed)) -tfm, mm = qd_rigid(fixed, moving, mxshift, mxrot, minwidth_rot, SD; thresh=thresh, rtol = 1e-7, atol = 1e-9 * thresh) +tfm, mm = qd_rigid(fixed, moving, mxshift, mxrot, minwidth_rot, SD; thresh=thresh, maxevals=1000, rtol = 1e-7, atol = 1e-9 * thresh) -@test mm < 1e-4 +@test mm < 3e-4 @test sum(abs.(vcat(tfm0.m[:], tfm0.v) - vcat(RotXYZ(tfm.m)[:], tfm.v))) < 0.1 @@ -407,7 +407,7 @@ tfm, mm = qd_rigid(fixed, moving, mxshift, mxrot, minwidth_rot, SD; thresh=thres moving = img2 - tform, mm = qd_rigid(fixed, moving, mxshift, mxrot, minwidth_rot, SD, rtol=0, fvalue=0.01) + tform, mm = qd_rigid(fixed, moving, mxshift, mxrot, minwidth_rot, SD, maxevals=1000, rtol=0, fvalue=0.01) # imgw = warp(img2, tform) # inds2 = intersect.(indices(img), indices(imgw)) diff --git a/test/register_penalty.jl b/test/register_penalty.jl index cbac0fb..8759176 100644 --- a/test/register_penalty.jl +++ b/test/register_penalty.jl @@ -251,7 +251,7 @@ function pfun(x, ϕs, ap, λt, mmis) RegisterPenalty.penalty!(nothing, similarϕ(ϕs, x), identity, ap, λt, mmis) end # This is needed for handling GradientNumbers -function similarϕ{Tϕ,N,A,L,Tx}(ϕs::Vector{GridDeformation{Tϕ,N,A,L}}, x::Array{Tx}) +function similarϕ(ϕs::Vector{GridDeformation{Tϕ,N,A,L}}, x::Array{Tx}) where {Tϕ,N,A,L,Tx} len = N*length(first(ϕs).u) length(x) == len*length(ϕs) || throw(DimensionMismatch("ϕs is incommensurate with a vector of length $(length(x))")) xf = RegisterDeformation.convert_to_fixed(SVector{N,Tx}, x, (size(first(ϕs).u)..., length(ϕs)))