High level functions for operations on MRI scans

3D Volume information

Minc2.voxel_to_worldFunction
voxel_to_world(grid::GridTransform)

Extract voxel to world affine transform from a GridTransform

  • grid - GridTransform
source
voxel_to_world(grid::InverseGridTransform)

Extract voxel to world affine transform from a InverseGridTransform

  • grid - InverseGridTransform
source
voxel_to_world(hdr::MincHeader)::AffineTransform{Float64}

Give AffineTransform for world to voxel transformation based on header

source
voxel_to_world(h::VolumeHandle,ijk::Vector{Float64})::Vector{Float64}

Convert contignuous 0-based voxel indexes (I,J,K) to world coordinates (X,Y,Z) 0-based

source
voxel_to_world(vol::Volume3D)

Extract voxel to world affine transform from a Volume3D

source
Minc2.world_to_voxelFunction
world_to_voxel(grid::GridTransform)

Extract world to voxel affine transform from a GridTransform

  • grid - GridTransform
source
world_to_voxel(grid::InverseGridTransform)

Extract world to voxel affine transform from a InverseGridTransform

  • grid - InverseGridTransform
source
world_to_voxel(h::VolumeHandle, xyz::Vector{Float64})::Vector{Float64}

Convert world coordinates (X,Y,Z) to contignuous voxel indexes (I,J,K) 0-based

source
world_to_voxel(hdr::MincHeader)::AffineTransform{Float64}

Give AffineTransform for voxel to world transformation

source
world_to_voxel(vol::Volume3D)

Extract world to voxel affine transform from a Volume3D

source
Minc2.arrayFunction
array(grid::GridTransform)

Extract underlying plain array

  • grid - GridTransform
source
array(grid::InverseGridTransform)

Extract underlying plain array

  • grid - InverseGridTransform
source
array(vol::Volume3D)

Extract underlying plain array

source

File IO functions

Loading 3D volumes

Minc2.read_volumeFunction
read_volume(fn::String; store::Type{T}=Float64)::Volume3D{T}

Read Volume3D from minc file

  • fn - filename
  • store - underlying array type
source
Minc2.read_nifti_volumeFunction
read_nifti_volume(fn::AbstractString; store::Type{T}=Float64)::Volume3D{T}

Read Volume3D from .nii or .nii.gz file

source

Saving 3D volumes

Minc2.save_volumeFunction
save_volume(fn::AbstractString, 
    vol::Volume3D{T,N}; 
    store::Type{S}=Float32,
    history=nothing)

Save Volume3D to minc file

  • fn - filename
  • vol - Volume3D to save
  • store - underlying MINC data type, to be used for storage
source
Minc2.save_nifti_volumeFunction
save_nifti_volume(fn::AbstractString, vol::Volume3D{T}; 
    store::Type{S}=Float32, history=nothing)

Save Volume3D into .nii or .nii.gz file

source

Loading transformations

Minc2.read_ants_transformFunction
read_ants_transform(fn::AbstractString; store::Type{T}=Float32)::Minc2.AnyTransform

Read .txt and .nii(.nii.gz) transforms produces by ANTs

source

Saving transformations

Minc2.save_transformsFunction
save_transforms(fname::String, 
    xfm::Union{Vector{XFM}, XFM};
    grid_store::Type{T}=Float32 ) where {T, XFM<:AnyTransform}

Save transformations into .xfm file,

Arguments

  • fname: output file name
  • xfm: vector of transformations to save, or a single transformation
  • grid_store: type of storage for grid files (default: Float32)
source
Minc2.save_itk_nifti_transformFunction
save_itk_nifti_transform(fn::AbstractString,
    xfm::Minc2.GridTransform{Float64,T}; store::Type{S}=Float32)

Write ANTs style warp transform

source

3D Volume creation

Minc2.Volume3DType
Volume3D{T,N}

An abstract 3D volume, could be vector field or volume with time dimension

source
Minc2.empty_volume_likeFunction
empty_volume_like(
        fn::String; 
        store::Type{T}=Float64, history=nothing)::Volume3D{T}

Create an empty Volume3D

  • fn - filename of a minc file that is used for sampling information
  • store - underlying array type
source
empty_volume_like(
    vol::Volume3D{T1,N}; 
    store::Type{T}=Float64, 
    history=nothing)

Create an empty Volume3D

  • vol - Volume3D that is used for sampling information
  • store - underlying array type
  • history - minc history
source
Minc2.full_volume_likeFunction
full_volume_like(
        fn::String,
        fill::T=zero(T); 
        store::Type{T}=Float64, history=nothing)::Volume3D{T}

Create an empty Volume3D

  • fn - filename of a minc file that is used for sampling information
  • store - underlying array type
source
full_volume_like(
    vol::Volume3D{T1,N}; 
    store::Type{T}=Float64, 
    history=nothing)

Create an empty Volume3D

  • vol - Volume3D that is used for sampling information
  • store - underlying array type
  • history - minc history
source

3D Volume manipulations

Minc2.resample_volumeFunction
resample_volume(
    in_vol::Volume3D{T,3};
    like::Union{Volume3D{O,3},Nothing}=nothing,
    tfm::Union{Vector{XFM},XFM,Nothing}=nothing, 
    itfm::Union{Vector{XFM},XFM,Nothing}=nothing, 
    interp::I=nothing, 
    fill=0.0,
    order=1,
    ftol=1.0/80,
    max_iter=10)::Volume3D

Resample Volume3D using transformation

  • in_vol - input Volume3D
  • like - Volume3D that is used for sampling information
  • itfm - inverse of the transformation to apply (i.e from output to input)
  • tfm - transformation to apply (i.e from output to input) (instead of itfm)
  • interp - interpolation method
  • fill - fill value
  • ftol - tolerance, for inverse transformations
  • max_iter - maximum number of iterations, for inverse transformations

Equivalent to mincresample minc command

source
Minc2.resample_volume!Function
resample_volume!(in_vol::AbstractArray{T,3}, 
    out_vol::AbstractArray{T,3}, 
    v2w::AffineTransform{C}, 
    w2v::AffineTransform{C}, 
    itfm::Union{Vector{XFM},XFM};
    interp::I=BSpline(Quadratic(Line(OnCell()))),
    fill=0.0,
    ftol=1.0/80,
    max_iter=10)

Resample 3D array using transformation

  • in_vol - input 3D array
  • out_vol - output 3D array
  • v2w - voxel to world affine transform in the output array
  • w2v - world to voxel affine transform in the input array
  • itfm - inverse of the transformation to apply (i.e from output to input)
  • interp - interpolation method
  • fill - fill value
  • ftol - tolerance, for inverse transformations
  • max_iter - maximum number of iterations, for inverse transformations
source
resample_volume!(
    in_vol::Volume3D{T,3}, 
    out_vol::Volume3D{O,3}; 
    tfm::Union{Vector{XFM},XFM,Nothing}=nothing, 
    itfm::Union{Vector{XFM},XFM,Nothing}=nothing, 
    interp::I=nothing, 
    fill=0.0, 
    order=nothing,
    ftol=1.0/80,
    max_iter=10)::Volume3D{O,3}

Resample Volume3D using transformation

  • in_vol - input Volume3D
  • out_vol - output Volume3D
  • itfm - inverse of the transformation to apply (i.e from output to input)
  • tfm - transformation to apply (i.e from output to input) (instead of itfm)
  • interp - interpolation method
  • fill - fill value
  • ftol - tolerance, for inverse transformations
  • max_iter - maximum number of iterations, for inverse transformations

Equivalent to mincresample minc command

source
Minc2.resample_gridFunction
resample_grid(
    in_grid::Volume3D{T,4}, 
    itfm::Union{Vector{XFM}, XFM}; 
    like::Union{Nothing,Volume3D{L,4}}=nothing)::Volume3D{T,4}

Resample Volume3D that contain 4D array, using transformation, assume 1st dimension is non spatial

  • in_grid - input Volume3D with 4D array describing vector field
  • itfm - inverse of the transformation to apply (i.e from output to input)
  • like - Volume3D that is used for sampling information
source
Minc2.resample_grid!Function
resample_grid!(
    in_grid::Volume3D{T,4},
    out_grid::Volume3D{T,4},
    itfm::Union{Vector{XFM}, XFM}=nothing)::Volume3D{T,4}

Resample Volume3D that contain 4D array, using transformation, assume 1st dimension is non spatial

  • in_grid - input Volume3D with 4D array describing vector field
  • out_grid - output Volume3D with 4D array describing vector field
  • itfm - inverse of the transformation to apply (i.e from output to input)
source
Minc2.crop_volumeFunction
crop_volume(in_vol::Volume3D{T,N},crop;
    fill_val::T=zero(T))::Volume3D{T,N}

Crop (or pad) a volume

  • in_vol - input Volume3D
  • crop - crop specification for each direction, e.g. [(10,10),(20,20),(5,10)] , negative values mean padding
source
Minc2.calculate_jacobianFunction
calculate_jacobian(
    tfm::Union{Vector{XFM},XFM};
    like::Volume3D{T,3}; 
    interp::I=BSpline(Quadratic(Line(OnCell()))),
    ftol=1.0/80,
    max_iter=10)::Volume3D{T,3}

Calculate dense jacobian determinant field for an arbitrary transformation

  • tfm - transformation to use
  • out_vol - output Volume3D
  • interp - interpolation method
  • ftol - tolerance, for inverse transformations
  • max_iter - maximum number of iterations, for inverse transformations

Roughly equivalent to mincblob minc command

source
Minc2.calculate_jacobian!Function
calculate_jacobian!(
    tfm::Union{Vector{XFM},XFM},
    out_vol::AbstractArray{T,3},
    out_v2w::AffineTransform{C};
    interp::I=BSpline(Quadratic(Line(OnCell()))),
    ftol=1.0/80,
    max_iter=10)

Calculate dense jacobian determinant field for an arbitrary transformation

  • tfm - transformation to use
  • out_vol - output 3D array
  • out_v2w - voxel to world affine transform in the output array
  • interp - interpolation method
  • ftol - tolerance, for inverse transformations
  • max_iter - maximum number of iterations, for inverse transformations

Roughly equivalent to mincblob minc command

source
calculate_jacobian!(
    tfm::Union{Vector{XFM},XFM}, 
    out_vol::Volume3D{T,3}; 
    interp::I=BSpline(Quadratic(Line(OnCell()))),
    ftol=1.0/80,
    max_iter=10)::Volume3D{T,3}

Calculate dense jacobian determinant field for an arbitrary transformation

  • tfm - transformation to use
  • out_vol - output Volume3D
  • interp - interpolation method
  • ftol - tolerance, for inverse transformations
  • max_iter - maximum number of iterations, for inverse transformations

Roughly equivalent to mincblob minc command

source
Missing docstring.

Missing docstring for Minc2.GridTransform. Check Documenter's build log for details.

Missing docstring.

Missing docstring for Minc2.InverseGridTransform. Check Documenter's build log for details.

Minc2.tfm_to_gridFunction
tfm_to_grid(
    tfm::Union{Vector{XFM}, XFM},
    ref::G;
    store::Type{T}=Float64,ftol=1.0/80,max_iter=10)::Volume3D{T,4}

Convert arbitrary transformation into Volume3D with 4D array containing vector field

  • tfm - transformation to use
  • v2w - voxel to world affine transform in the output array
  • ftol - tolerance, for inverse transformations
  • max_iter - maximum number of iterations, for inverse transformations
source
Minc2.normalize_tfmFunction
normalize_tfm(tfm::Union{Vector{XFM}, XFM},
    ref::G;
    store::Type{T}=Float64,ftol=1.0/80,max_iter=10)::GridTransform{Float64,T}

Convert arbitrary transformation into a single GridTransform

  • tfm - transformation to use
  • ref - GridTransform that is used for sampling information
  • store - underlying array type
  • ftol - tolerance, for inverse transformations
  • max_iter - maximum number of iterations, for inverse transformations
source
normalize_tfm(tfm::Union{Vector{XFM}, XFM},
    ref::G;
    store::Type{T}=Float64,ftol=1.0/80,max_iter=10)::GridTransform{Float64,T}

Convert arbitrary transformation into a single GridTransform

  • tfm - transformation to use
  • ref - Volume3D that is used for sampling information
  • store - underlying array type
  • ftol - tolerance, for inverse transformations
  • max_iter - maximum number of iterations, for inverse transformations
source