Translate images by subpixel amounts

Hello,

I’m trying to translate an image by arbitrary amounts, and am currently doing translations via JuliaImages and CoordinateTransformation packages as such:

t = Translation(y_shift, x_shift)
img = warp(img, t, indices_spatial(img), 0)

However, this requires that y_shift and x_shift are integers. I am wondering if there is support for non-integer translation via interpolation. Is there a way to use warp with any affine transformation? Is there a better way than warp?

Thanks!

1 Like

For 1D signals you can do the translation in the frequency domain, which corresponds to a circular shift. With the proper zero-padding you can do a non-circular shift.

using FFTW: rfft, irfft
using Plots: plot
using DSP: resample
julia> x = randn(50);

julia> function fftshift(x, d)
           L = length(x)
           irfft(rfft(x) .* exp.(-im * d .* range(0, π, length=L÷2+1)), L)
       end
fftshift (generic function with 1 method)

# zero-pad by more than your shift amount
julia> xpad = [x;0];

julia> x2 = fftshift(xpad, 0.4);

# upsample for plotting so you can see between samples
julia> plot([resample(xpad, 8) resample(x2, 8)])

You should be able to do something similar using a 2D FFT.

4 Likes

Not for me (just curious, where did you get that idea?):

julia> t = Translation(0.2, -1.4)

julia> imgs = warp(img, t, indices_spatial(img));

If you want something beyond linear interpolation, try e.g.,

imgs = warp(float.(img), t, indices_spatial(img), Quadratic(Flat(OnCell())), 0)

The float may be necessary if img itself has a restrictive eltype, to ensure that you can represent the output appropriately.

3 Likes

Hi @tim.holy, thanks for the response!

I can’t even reproduce the exact error I was getting from my program anymore :sweat_smile:. However, my first thought was that it was because integer coordinate transformations were required, and adding x_shift, y_shift = round(x_shift), round(y_shift), translations began to work.

Of course, I didn’t just take this as truth and did some testing in the REPL. Turns out that my toy tests were bad because my “images” were integer arrays, and that only seems to work with certain translations. For example

julia> a = [1 2 3; 4 5 6];

julia> t = Translation(.75, .75);

julia> b = warp(a, t, indices_spatial(a), 0); # works

julia> t = Translation(.9, .9);

julia> b = warp(a, t, indices_spatial(a), 0); # doesn't work
ERROR: InexactError: Int64(4.6000000000000005)
Stacktrace:
 [1] Type at ./float.jl:703 [inlined]
 [2] convert at ./number.jl:7 [inlined]
 [3] setindex! at ./array.jl:768 [inlined]
 [4] setindex! at ./multidimensional.jl:488 [inlined]
 [5] warp!(::Array{Int64,2}, ::Interpolations.FilledExtrapolation{Int64,2,Interpolations.BSplineInterpolation{Int64,2,Array{Int64,2},BSpline{Linear},Tuple{Base.OneTo{Int64},Base.OneTo{Int64}}},BSpline{Linear},Int64}, ::Translation{StaticArrays.SArray{Tuple{2},Float64,1,2}}) at /Users/John/.juliapro/JuliaPro_v1.2.0-1/packages/ImageTransformations/7wC0C/src/warp.jl:89
 [6] warp at /Users/John/.juliapro/JuliaPro_v1.2.0-1/packages/ImageTransformations/7wC0C/src/warp.jl:84 [inlined]
 [7] warp(::Array{Int64,2}, ::Translation{StaticArrays.SArray{Tuple{2},Float64,1,2}}, ::Tuple{Base.OneTo{Int64},Base.OneTo{Int64}}, ::Int64) at /Users/John/.juliapro/JuliaPro_v1.2.0-1/packages/ImageTransformations/7wC0C/src/warp.jl:96
 [8] top-level scope at none:0

julia> b = warp(float.(a), t, indices_spatial(a), 0); # now works!

So, my testing in the REPL falsely led me to believe that integer shifts were required. I wish I could reproduce the error I was having in my program, however, because the image was already float values.

I do have another question tho. I am looking at the docs and it says there are descriptions of the boundary condition methods under /doc/latex, but I cannot find that in the repo, and I can’t understand the different interpolation methods such as Quadratic(Flat(OnCell())).

The behavior I am looking for is:

a = [1 2; 3 4]
t = Translation(.5, .5) 

^ just noticed that this also does not work without float casting, and this was probably the example that I first tried and saw InexactError (2.5)

b = warp(a, t, indices_spatial(a), 0)
[2.5 1.5; 1.75 1]

However, the edge pixels become zeros. Is there a method to produce the behavior above?

Sorry, I meant to respond to this @ssfrr. I really appreciate it. Although I was not able to implement a solution using FFT myself.

Somehow I stumbled across this again. Fixed in https://github.com/JuliaImages/ImageTransformations.jl/pull/88 (released as version 0.8.3)

3 Likes