MutualInformationImageRegistration

MutualInformationImageRegistration performs image registration using mutual information.

Here's a complete example on how to use this package:

using MutualInformationImageRegistration, MutualInformationImageRegistration.FastHistograms, Random

# Create the container used to hold intermediate variables for registration
mi = MutualInformationContainer(
    create_fast_histogram(
        FastHistograms.FixedWidth(),
        FastHistograms.Arithmetic(),
        FastHistograms.NoParallelization(),
        [(0x00, 0xff, 4), (0x00, 0xff, 4)],
    ),
)

# Create the full image that the smaller images to register will be pulled from
full_image = rand(UInt8, 500, 300)
view(full_image, 300:330, 200:220) .= 0xff

# The fixed image is the image that the other images are registered against
fixed = full_image[(300-10):(330+10), (200-10):(220+10)]

# The max shift is the maximum search range
MAX_SHIFT = 11
# Padding is used to grow the bbox for higher quality registration
padding = [-10, -10, 10, 10]

# A buffer is needed to hold intermediate data
buffer = Array{UInt8}(undef, (size(fixed) .+ (MAX_SHIFT * 2))...)

# Introduce a random shift to the moving bbox
expected_x = rand(-5:5)
expected_y = rand(-5:5)

# Register the image given by the bbox (called the moving bbox) against the fixed image
shift, mm, mms = register!(
    mi,
    full_image,
    fixed,
    [300, 200, 330, 220] .+ padding .+ [expected_x, expected_y, expected_x, expected_y],
    MAX_SHIFT,
    MAX_SHIFT,
    buffer,
)

# The shift we get out should be equal and opposite of the shift we applied
shift == (-expected_x, -expected_y)

# output

true

Also look at the non-exported traits NoParallelization and SIMD, which can be used when constructing a MutualInformationContainer. NoParallelization is the default; other traits can be used to customize how the mutual information computation is run.

A trait for the ways the mutual information calculation can be parallelized.

No threading nor vectorization.

SIMD vectorization.

Compute the MI. The hist inside mi must already be incremeted.

mutual_information!(
    mi::MutualInformationContainer,
    fixed,
    buffer,
    ::Any,
    moving_bbox,
    range_x,
    range_y,
    prev_mis::AbstractArray{Float32,2};
    get_buffer_crop,
    kwargs...
)

Calculates the mutual information of two images at all shifts within the range_x and range_y. Warm-starts the evaluation using previous results (prev_mis; the return value from a previous call of this function) and using the previously set and filtered contents of the buffer.

mutual_information!(
    mi::MutualInformationContainer,
    fixed,
    buffer,
    full_image,
    moving_bbox,
    range_x,
    range_y,
    ::Missing;
    set_buffer!,
    get_buffer_crop,
    prefilter_frame_crop! = x -> nothing,
)

Calculates the mutual information of two images at all shifts within the range_x and range_y. The fixed image must already be filtered. This will set the buffer and filter its contents using prefilter_frame_crop!.

mutual_information!(mi::MutualInformationContainer, x, y)

Computes the mutual information between the two variables x and y. The histogram within mi must be of the correct type to handle the formats of x and y.

register!(
    mi::MutualInformationContainer,
    full_image::AbstractArray{T,2},
    fixed::AbstractArray{T,2},
    moving_bbox::AbstractVector{Int},
    range_x::AbstractVector{Int},
    range_y::AbstractVector{Int},
    buffer::AbstractArray{T,2},
    prev_mis::Union{Missing,AbstractArray{Float32,2}};
    set_buffer! = (buffer, current_frame, moving_bbox) -> set_buffer!(buffer, current_frame, moving_bbox, maximum(range_x), maximum(range_y)),
    get_buffer_crop = (buffer, moving_bbox, shift_x, shift_y) -> get_buffer_crop(buffer, moving_bbox, shift_x, shift_y, maximum(range_x), maximum(range_y)),
    prefilter_frame_crop! = x -> nothing,
) where {T<:Integer}

Calculates the shift that best aligns the moving_bbox to the fixed image within the full_image. At a high level, this attempts to best match the view of moving_bbox inside full_image to the fixed image. This only considers shifts along the x and y axes; no rotation is considered. All combinations of shifts within range_x and range_y are considered.

Adding some padding to moving_bbox is a good idea to improve registration stability. E.g. my_bbox .+ [-10, -10, 10, 10]. You will need to determine the best padding value for your data.

The parameter buffer is required for temporary storage between calls to this function because it is assumed that you will call this function in a loop to register many similar images. Generally, you can define buffer as Array{T}(undef, (size(fixed) .+ (MAX_SHIFT * 2))...). The buffer must have a size which is at least the size of the fixed image expanded by the maximum and minimum x and y shifts on each respective axis. The parameters set_buffer! and get_buffer_crop are used to write to and read from this buffer, respectively. There is typically no need to change these from their defaults.

The parameter prev_mis is used to memoize the mutual information calculations if this function is being called "incrementally" with an expanding shift horizon. If you are calling this function directly (and not in a loop to gradually expand a maximum shift horizon), you should set this to missing, which would cause all combinations of shifts within range_x and range_y to be considered.

The parameter prefilter_frame_crop! can be specified if you want to apply image filtering before computing the mutual information between the two images. This function must mutate the image it is given with whatever filtering operation you implement. Also, the fixed image must have the filtering pre-applied. See this package's tests for an example.

register!(
    mi::MutualInformationContainer,
    full_image::AbstractArray{T,2},
    fixed::AbstractArray{T,2},
    moving_bbox::AbstractVector{Int},
    max_shift_x::Int,
    max_shift_y::Int,
    buffer::AbstractArray{T,2};
    set_buffer! = (buffer, current_frame, moving_bbox) -> set_buffer!(buffer, current_frame, moving_bbox, max_shift_x, max_shift_y),
    get_buffer_crop = (buffer, moving_bbox, shift_x, shift_y) -> get_buffer_crop(buffer, moving_bbox, shift_x, shift_y, max_shift_x, max_shift_y),
    prefilter_frame_crop! = x -> nothing,
    start_shift_x = 3,
    start_shift_y = 3,
    expand_border = 1,
    expand_increment = 1,
) where {T<:Integer}

Calculates the shift that best aligns the moving_bbox to the fixed image within the full_image. At a high level, this attempts to best match the view of moving_bbox inside full_image to the fixed image. This only considers shifts along the x and y axes; no rotation is considered. This function incrementally expands the maximum shift horizon to evaluate as few shifts as possible. At a minimum, all shifts within ± start_shift_x and ± start_shift_y are considered. If the optimal shift falls within expand_border of the horizon, the horizon will be expanded by expand_increment and all new shifts will be evaluated. This process repeats until the optimal shift does not fall within expand_border of the horizon, or until the horizon has reached max_shift_x and max_shift_y.

Adding some padding to moving_bbox is a good idea to improve registration stability. E.g. my_bbox .+ [-10, -10, 10, 10]. You will need to determine the best padding value for your data.

The parameter buffer is required for temporary storage. Generally, you can define buffer as Array{T}(undef, (size(fixed) .+ (MAX_SHIFT * 2))...). The buffer must have a size which is at least the size of the fixed image expanded by the maximum and minimum x and y shifts on each respective axis. The parameters set_buffer! and get_buffer_crop are used to write to and read from this buffer, respectively. There is typically no need to change these from their defaults.

The parameter prefilter_frame_crop! can be specified if you want to apply image filtering before computing the mutual information between the two images. This function must mutate the image it is given with whatever filtering operation you implement. Also, the fixed image must have the filtering pre-applied. See this package's tests for an example.