ImageGeoms overview¶

This page explains the Julia package ImageGeoms.

This page was generated from a single Julia file: 1-example.jl.

Setup¶

Packages needed here.

In [ ]:

using ImageGeoms
using MIRTjim: jim, prompt
using Plots: scatter, plot!, default; default(markerstrokecolor=:auto)
using InteractiveUtils: versioninfo

The following line is helpful when running this file as a script; this way it will prompt user to hit a key after each figure is displayed.

In [ ]:

isinteractive() ? jim(:prompt, true) : prompt(:draw);

Overview¶

When performing tomographic image reconstruction, one must specify the geometry of the grid of image pixels. (In contrast, for image denoising and image deblurring problems, one works with the given discrete image and no physical coordinates are needed.)

The key parameters of a grid of image pixels are

the size (dimensions) of the grid, e.g., 128 × 128,
the spacing of the pixels, e.g., 1mm × 1mm,
the offset of the pixels relative to the origin, e.g., (0,0)

The data type ImageGeom describes such a geometry, for arbitrary dimensions (2D, 3D, etc.).

There are several ways to construct this structure. The default is a 128 × 128 grid with pixel size \Delta_X = \Delta_Y = 1 (unitless) and zero offset:

In [ ]:

ig = ImageGeom()

Here is a 3D example with non-cubic voxel size:

In [ ]:

ig = ImageGeom( (512,512,128), (1,1,2), (0,0,0) )

To avoid remembering the order of the arguments, named keyword pairs are also supported:

In [ ]:

ig = ImageGeom( dims=(512,512,128), deltas=(1,1,2), offsets=(0,0,0) )

Units¶

The pixel dimensions deltas can (and should!) be values with units.

Here is an example for a video (2D+time) with 12 frames per second:

In [ ]:

using UnitfulRecipes
using Unitful: mm, s

ig = ImageGeom( dims=(640,480,1000), deltas=(1mm,1mm,(1//12)s) )

Methods¶

An ImageGeom object has quite a few methods; axes and axis are especially useful:

In [ ]:

ig = ImageGeom( dims=(7,8), deltas=(3,2), offsets=(0,0.5) )
axis(ig, 2)

or an axis of length n with spacing Δ (possibly with units) and (always unitless but possibly non-integer) offset the axis is a subtype of AbstractRange of the form ( (0:n-1) .- ((n - 1)/2 + offset) ) * Δ

These axes are useful for plotting:

In [ ]:

ig = ImageGeom( dims=(12,10), deltas=(1mm,1mm), offsets=(0.5,0.5) )

In [ ]:

showgrid = (ig) -> begin # x,y grid locations of pixel centers
    x = axis(ig, 1)
    y = axis(ig, 2)
    (xg, yg) = grids(ig)
    scatter(xg, yg, label="", xlabel="x", ylabel="y",
        xlims = maximum(abs, x) * 1.2 .* (-1,1),
        xticks = [x[1], zero(eltype(x)), x[end]],
        ylims = maximum(abs, y) * 1.2 .* (-1,1),
        yticks = [y[1], zero(eltype(y)), y[end]],
        aspect_ratio=1, title="offsets $(ig.offsets)")
end
showgrid(ig)

Unit labels on the axes due to UnitfulRecipes.jl

In [ ]:

prompt();

Offsets (unitless translation of grid)¶

The default offsets are zeros, corresponding to symmetric sampling around origin:

In [ ]:

ig = ImageGeom( dims=(12,10), deltas=(1mm,1mm) )
p = showgrid(ig)

In [ ]:

prompt();

That default for offsets is natural for tomography when considering finite pixel size:

In [ ]:

square = (x,y,Δ) -> plot!(p, label="", color=:black,
    x .+ Δ[1] * ([0,1,1,0,0] .- 0.5),
    y .+ Δ[2] * ([0,0,1,1,0] .- 0.5),
)
square2 = (x,y) -> square(x, y, ig.deltas)
square2.(grids(ig)...)
plot!(p)

In [ ]:

prompt();

In that default geometry, the center (0,0) of the image is at a corner of the middle 4 pixels (for even image sizes). That default is typical for tomographic imaging (e.g., CT, PET, SPECT). One must be careful when using operations like imrotate or fft.

AxisArrays¶

There is a natural connection between ImageGeom and AxisArrays. Note the automatic labeling of units (when relevant) on all axes by MIRTjim.jim.

In [ ]:

using AxisArrays
using Unitful: mm
ig = ImageGeom( dims=(60,48), deltas=(1.5mm,1mm) )
za = AxisArray( ellipse(ig) * 10/mm ; x=axis(ig,1), y=axis(ig,2) )
jim(za, "AxisArray example")

In [ ]:

prompt();

Resizing¶

Often we have a target grid in mind but want coarser sampling for debugging. The downsample method is useful for this.

In [ ]:

ig = ImageGeom( dims = (512,512), deltas = (500mm,500mm) ./ 512 )
ig_down = downsample(ig, 4)

Other times we want to avoid an "inverse crime" by using finer sampling to simulate data; use oversample for this.

In [ ]:

ig_over = oversample(ig, (2,2))

Reproducibility¶

This page was generated with the following version of Julia:

In [ ]:

io = IOBuffer()
versioninfo(io)
split(String(take!(io)), '\n')

And with the following package versions

In [ ]:

import Pkg; Pkg.status()

This notebook was generated using Literate.jl.