Clarification of rotations

[henry2004y]

Hi,

I have a Python script in which there is a call to scipy.ndimage.rotate for rotating a 2D array, similar to the following

using PyCall
ndimage = pyimport("scipy.ndimage")
# test array
img = Float64[1 2 3 4; 5 6 7 8; 9 10 11 12; 13 14 15 16]
img_90 = ndimage.rotate(img, 90, reshape=false)

I came across this package when searching for the same kind of operation on a plain 2D array in Julia. However, soon I got confused by the rotation methods provided here, especially imrotate.

using ImageTransformations
using Rotations: RotMatrix
using CoordinateTransformations: recenter

# test array
img = Float64[1 2 3 4; 5 6 7 8; 9 10 11 12; 13 14 15 16]

# define transformation
trfm = recenter(RotMatrix(pi/2), center(img))
imgw = warp(img, trfm) # this is what I expected, same as in scipy

imgr = imrotate(img, pi/2) # this gives a 6x6 OffsetArray with values that I don't quite understand

where imrotate returned

6×6 OffsetArray(::Matrix{Float64}, 0:5, 0:5) with eltype Float64 with indices 0:5×0:5:
 NaN  NaN    NaN        NaN        NaN        NaN
 NaN  NaN    NaN          4.99995  NaN        NaN
 NaN   14.0    9.99998    5.99999  NaN        NaN
 NaN  NaN     11.0        7.00002    3.00003  NaN
 NaN  NaN     12.0      NaN        NaN        NaN
 NaN  NaN    NaN        NaN        NaN        NaN

I tried the example given in the doc for imrotate, and it looked just fine.

Can you kindly point out what I missed here? What is the closest possible translation from scipy.ndimage.rotate?

[johnnychen94]

There are two issues here:

• There are some numerical errors on pi/2 because you can't represent pi arbitrarily precisely. (see also #145)
• The default imrotate output is the minimal array that encloses all pixels in the original image

There are two options to workaround this:

• use the Base function rotr90(img) (also rotl90 and rot180)
• clip the output by passing axes: imrotate(img, pi/2, axes(img))

I'm going to close this issue as it's very much duplicated of #145

[henry2004y]

Rotating the 2D array by 90 degree is just for testing purposes: the actual angle is in principle any value in [0, 2pi). What puzzles me is the behavior of imrotate, especially its not NaN structure after rotation. As a first time user to this function, I believe I am not the only one to be confused.

Even with axes:

julia> imrotate(img, pi/2, axes(img))
4×4 Matrix{Float64}:
 NaN    NaN          4.99995  NaN
  14.0    9.99998    5.99999  NaN
 NaN     11.0        7.00002    3.00003
 NaN     12.0      NaN        NaN

I am not sure if using any kind of low order interpolation (e.g. Closest()) will avoid the NaN fillings.

Since my goal is to reproduce the result in scipy, I guess I should go with warp then?

[timholy]

I'm reopening because the differences are due to the strategy for handling roundoff error in https://github.com/JuliaImages/ImageTransformations.jl/pull/58#discussion_r226829801:

julia> using Debugger

julia> @enter imrotate(img, π/2, axes(img))
In #imrotate#9(kwargs, , img, θ, inds) at /home/tim/.julia/packages/ImageTransformations/O9GYH/src/warp.jl:245
 245  function imrotate(img::AbstractArray{T}, θ::Real, inds::Union{Tuple, Nothing} = nothing; kwargs...) where T
 246      # TODO: expose rotation center as a keyword
>247      θ = floor(mod(θ,2pi)*typemax(Int16))/typemax(Int16) # periodic discretezation
 248      tform = recenter(RotMatrix{2}(θ), center(img))
 249      # Use the `nothing` trick here because moving the `autorange` as default value is not type-stable
 250      inds = isnothing(inds) ? autorange(img, inv(tform)) : inds
 251      warp(img, tform, inds; kwargs...)
 252  end

About to run: (*)(2, π)
1|debug> n
In #imrotate#9(kwargs, , img, θ, inds) at /home/tim/.julia/packages/ImageTransformations/O9GYH/src/warp.jl:245
 245  function imrotate(img::AbstractArray{T}, θ::Real, inds::Union{Tuple, Nothing} = nothing; kwargs...) where T
 246      # TODO: expose rotation center as a keyword
 247      θ = floor(mod(θ,2pi)*typemax(Int16))/typemax(Int16) # periodic discretezation
>248      tform = recenter(RotMatrix{2}(θ), center(img))
 249      # Use the `nothing` trick here because moving the `autorange` as default value is not type-stable
 250      inds = isnothing(inds) ? autorange(img, inv(tform)) : inds
 251      warp(img, tform, inds; kwargs...)
 252  end

About to run: (Core.apply_type)(RotMatrix, 2)
1|julia> θ
1.5707876827295755

1|julia> π/2
1.5707963267948966

That's a pretty large difference from what the user entered.

[henry2004y]

Thanks! Another tiny issue I want to mention is the fillvalue keyword:

julia> imrotate(img, pi/2, axes(img), fillvalue=0.0)
4×4 Matrix{Float64}:
 NaN    NaN          4.99995  NaN
  14.0    9.99998    5.99999  NaN
 NaN     11.0        7.00002    3.00003
 NaN     12.0      NaN        NaN

Seemingly contrary to the method doc

• fillvalue: the value that used to fill the new region. The default value is NaN if possible, otherwise is 0.

I don't observe any changes to the output in this case.

[timholy]

I can't replicate that; in a revised version of #148 I added these tests:

@test  any(isnan, imrotate(img, π/3))
@test !any(isnan, imrotate(img, π/3; fillvalue=0.0))
@test !any(isnan, imrotate(img, π/3, fillvalue=0.0))
@test  any(isnan, imrotate(img, π/3, axes(img)))
@test !any(isnan, imrotate(img, π/3, axes(img); fillvalue=0.0))
@test !any(isnan, imrotate(img, π/3, axes(img), fillvalue=0.0))

and they all pass.

[henry2004y]

Sorry I just found that I was accidentally using an earlier version 0.8.13, probably due to some restrictions of other packages.

[timholy]

We really need to get Images 0.25 out!

julia> using RegistryCompatTools

help?> RegistryCompatTools
search: RegistryCompatTools

  RegistryCompatTools makes it easier to discover packages that need [compat] updates.

    •  held_back_by("SomePkg") - which packages are holding back SomePkg?

    •  held_back_packages()["SomePkg"] - which packages does SomePkg hold back?

    •  find_julia_packages_github() - discover which packages I can push to

julia> held_back_by("ImageTransformations")
15-element Vector{String}:
 "ADI"
 "Asimov"
 "Augmentor"
 "BasicTextRender"
 "CameraCalibrations"
 "DataAugmentation"
 "EasyML"
 "HCIToolbox"
 "Images"                    # not good...
 "MIRT"
 "MakieGallery"
 "Photometry"
 "UNet"
 "WordCloud"
 "YOLO"