diff --git a/KomaMRIBase/src/datatypes/Phantom.jl b/KomaMRIBase/src/datatypes/Phantom.jl index 94ee8fa3f..e7c6726b7 100644 --- a/KomaMRIBase/src/datatypes/Phantom.jl +++ b/KomaMRIBase/src/datatypes/Phantom.jl @@ -212,7 +212,7 @@ Default ss=4 sample spacing is 2 mm. Original file (ss=1) sample spacing is .5 m - `axis`: (`::String`, `="axial"`, opts=[`"axial"`, `"coronal"`, `"sagittal"`]) orientation of the phantom - `ss`: (`::Integer or ::Vector{Integer}`, `=4`) subsampling parameter for all axes if scaler, per axis if 2 element vector [ssx, ssy] - `us`: (`::Integer or ::Vector{Integer}`, `=1`) upsampling parameter for all axes if scaler, per axis if 2 element vector [usx, usy], if used ss is set to ss=1 - +- `tissue_properties`: (`::Dict`, `=Dict()`) phantom tissue properties in SI units considering the available tissues # Returns - `obj`: (`::Phantom`) Phantom struct @@ -223,10 +223,26 @@ julia> obj = brain_phantom2D(; axis="sagittal", ss=1) julia> obj = brain_phantom2D(; axis="axial", us=[1, 2]) +julia> phantom_values = + Dict( + # ρ, T1, T2, T2*, Δw + "CSF" => [1, 2.569, 0.329, 0.058, 0], + "GM" => [0.86, 0.833, 0.083, 0.069, 0], + "WM" => [0.77, 0.500, 0.070, 0.061, 0], + "FAT1" => [0, 0, 0, 0, 0], + "MUSCLE" => [0, 0, 0, 0, 0], + "SKIN/MUSCLE" => [0, 0, 0, 0, 0], + "SKULL" => [0, 0, 0, 0, 0], + "VESSELS" => [0, 0, 0, 0, 0], + "FAT2" => [0, 0, 0, 0, 0], + "DURA" => [0, 0, 0, 0, 0], + "MARROW" => [0, 0, 0, 0, 0]) +julia> obj = brain_phantom2D(; tissue_properties=phantom_values) + julia> plot_phantom_map(obj, :ρ) ``` """ -function brain_phantom2D(; axis="axial", ss=4, us=1) +function brain_phantom2D(; axis="axial", ss=4, us=1, tissue_properties = Dict()) # check and filter input ssx, ssy, ssz, usx, usy, usz = check_phantom_arguments(2, ss, us) @@ -235,75 +251,25 @@ function brain_phantom2D(; axis="axial", ss=4, us=1) data = MAT.matread(path * "/phantom/brain2D.mat") # subsample or upsample the phantom data - class = repeat(data[axis][1:ssx:end, 1:ssy:end]; inner=[usx, usy]) + labels = repeat(data[axis][1:ssx:end, 1:ssy:end]; inner=[usx, usy]) + # to make it compatible with default_brain_tissue_properties + labels = reshape(labels, (size(labels, 1), size(labels, 2), 1)) # Define spin position vectors Δx = .5e-3 * ssx / usx Δy = .5e-3 * ssy / usy - M, N = size(class) + M, N = size(labels) FOVx = (M - 1) * Δx #[m] FOVy = (N - 1) * Δy #[m] x = (-FOVx / 2):Δx:(FOVx / 2) #spin coordinates y = (-FOVy / 2):Δy:(FOVy / 2) #spin coordinates x, y = x .+ y' * 0, x * 0 .+ y' #grid points + x = reshape(x, (size(x, 1), size(x, 2), 1)) + y = reshape(y, (size(y, 1), size(y, 2), 1)) - # Define spin property vectors - T2 = - (class .== 23) * 329 .+ #CSF - (class .== 46) * 83 .+ #GM - (class .== 70) * 70 .+ #WM - (class .== 93) * 70 .+ #FAT1 - (class .== 116) * 47 .+ #MUSCLE - (class .== 139) * 329 .+ #SKIN/MUSCLE - (class .== 162) * 0 .+ #SKULL - (class .== 185) * 0 .+ #VESSELS - (class .== 209) * 70 .+ #FAT2 - (class .== 232) * 329 .+ #DURA - (class .== 255) * 70 #MARROW - T2s = - (class .== 23) * 58 .+ #CSF - (class .== 46) * 69 .+ #GM - (class .== 70) * 61 .+ #WM - (class .== 93) * 58 .+ #FAT1 - (class .== 116) * 30 .+ #MUSCLE - (class .== 139) * 58 .+ #SKIN/MUSCLE - (class .== 162) * 0 .+ #SKULL - (class .== 185) * 0 .+ #VESSELS - (class .== 209) * 61 .+ #FAT2 - (class .== 232) * 58 .+ #DURA - (class .== 255) * 61 #MARROW - T1 = - (class .== 23) * 2569 .+ #CSF - (class .== 46) * 833 .+ #GM - (class .== 70) * 500 .+ #WM - (class .== 93) * 350 .+ #FAT1 - (class .== 116) * 900 .+ #MUSCLE - (class .== 139) * 569 .+ #SKIN/MUSCLE - (class .== 162) * 0 .+ #SKULL - (class .== 185) * 0 .+ #VESSELS - (class .== 209) * 500 .+ #FAT2 - (class .== 232) * 2569 .+ #DURA - (class .== 255) * 500 #MARROW - ρ = - (class .== 23) * 1 .+ #CSF - (class .== 46) * 0.86 .+ #GM - (class .== 70) * 0.77 .+ #WM - (class .== 93) * 1 .+ #FAT1 - (class .== 116) * 1 .+ #MUSCLE - (class .== 139) * 1 .+ #SKIN/MUSCLE - (class .== 162) * 0 .+ #SKULL - (class .== 185) * 0 .+ #VESSELS - (class .== 209) * 0.77 .+ #FAT2 - (class .== 232) * 1 .+ #DURA - (class .== 255) * 0.77 #MARROW - Δw_fat = -220 * 2π - Δw = - (class .== 93) * Δw_fat .+ #FAT1 - (class .== 209) * Δw_fat #FAT2 - T1 = T1 * 1e-3 - T2 = T2 * 1e-3 - T2s = T2s * 1e-3 - + # Get tissue properties + ρ, T1, T2, T2s, Δw = default_brain_tissue_properties(labels, tissue_properties) + # Define and return the Phantom struct obj = Phantom{Float64}(; name="brain2D_" * axis, @@ -336,6 +302,7 @@ Default ss=4 sample spacing is 2 mm. Original file (ss=1) sample spacing is .5 m - `ss`: (`::Integer or ::Vector{Integer}`, `=4`) subsampling parameter for all axes if scaler, per axis if 3 element vector [ssx, ssy, ssz] - `us`: (`::Integer or ::Vector{Integer}`, `=1`) upsampling parameter for all axes if scaler, per axis if 3 element vector [usx, usy, usz] - `start_end`: (`::Vector{Integer}`, `=[160,200]`) z index range of presampled phantom, 180 is center +- `tissue_properties`: (`::Dict`, `=Dict()`) phantom tissue properties in SI units considering the available tissues # Returns - `obj`: (`::Phantom`) 3D Phantom struct @@ -346,19 +313,34 @@ julia> obj = brain_phantom3D(; ss=5) julia> obj = brain_phantom3D(; us=[2, 2, 1]) +julia> phantom_values = + Dict( + # ρ, T1, T2, T2*, Δw + "CSF" => [1, 2.569, 0.329, 0.058, 0], + "GM" => [0.86, 0.833, 0.083, 0.069, 0], + "WM" => [0.77, 0.500, 0.070, 0.061, 0], + "FAT1" => [0, 0, 0, 0, 0], + "MUSCLE" => [0, 0, 0, 0, 0], + "SKIN/MUSCLE" => [0, 0, 0, 0, 0], + "SKULL" => [0, 0, 0, 0, 0], + "VESSELS" => [0, 0, 0, 0, 0], + "FAT2" => [0, 0, 0, 0, 0], + "DURA" => [0, 0, 0, 0, 0], + "MARROW" => [0, 0, 0, 0, 0]) +julia> obj = brain_phantom3D(; tissue_properties=phantom_values) + julia> plot_phantom_map(obj, :ρ) ``` """ -function brain_phantom3D(; ss=4, us=1, start_end=[160, 200]) +function brain_phantom3D(; ss=4, us=1, start_end=[160, 200], tissue_properties=Dict()) # check and filter input ssx, ssy, ssz, usx, usy, usz = check_phantom_arguments(3, ss, us) - # Get data from .mat file path = @__DIR__ data = MAT.matread(path * "/phantom/brain3D.mat") # subsample or upsample the phantom data - class = repeat( + labels = repeat( data["data"][1:ssx:end, 1:ssy:end, start_end[1]:ssz:start_end[2]]; inner=[usx, usy, usz], ) @@ -367,7 +349,7 @@ function brain_phantom3D(; ss=4, us=1, start_end=[160, 200]) Δx = .5e-3 * ssx / usx Δy = .5e-3 * ssy / usy Δz = .5e-3 * ssz / usz - M, N, Z = size(class) + M, N, Z = size(labels) FOVx = (M - 1) * Δx #[m] FOVy = (N - 1) * Δy #[m] FOVz = (Z - 1) * Δz #[m] @@ -377,63 +359,9 @@ function brain_phantom3D(; ss=4, us=1, start_end=[160, 200]) x = 1 * xx .+ 0 * yy .+ 0 * zz y = 0 * xx .+ 1 * yy .+ 0 * zz z = 0 * xx .+ 0 * yy .+ 1 * zz - - # Define spin property vectors - T2 = - (class .== 23) * 329 .+ #CSF - (class .== 46) * 83 .+ #GM - (class .== 70) * 70 .+ #WM - (class .== 93) * 70 .+ #FAT1 - (class .== 116) * 47 .+ #MUSCLE - (class .== 139) * 329 .+ #SKIN/MUSCLE - (class .== 162) * 0 .+ #SKULL - (class .== 185) * 0 .+ #VESSELS - (class .== 209) * 70 .+ #FAT2 - (class .== 232) * 329 .+ #DURA - (class .== 255) * 70 #MARROW - T2s = - (class .== 23) * 58 .+ #CSF - (class .== 46) * 69 .+ #GM - (class .== 70) * 61 .+ #WM - (class .== 93) * 58 .+ #FAT1 - (class .== 116) * 30 .+ #MUSCLE - (class .== 139) * 58 .+ #SKIN/MUSCLE - (class .== 162) * 0 .+ #SKULL - (class .== 185) * 0 .+ #VESSELS - (class .== 209) * 61 .+ #FAT2 - (class .== 232) * 58 .+ #DURA - (class .== 255) * 61 #MARROW - T1 = - (class .== 23) * 2569 .+ #CSF - (class .== 46) * 833 .+ #GM - (class .== 70) * 500 .+ #WM - (class .== 93) * 350 .+ #FAT1 - (class .== 116) * 900 .+ #MUSCLE - (class .== 139) * 569 .+ #SKIN/MUSCLE - (class .== 162) * 0 .+ #SKULL - (class .== 185) * 0 .+ #VESSELS - (class .== 209) * 500 .+ #FAT2 - (class .== 232) * 2569 .+ #DURA - (class .== 255) * 500 #MARROW - ρ = - (class .== 23) * 1 .+ #CSF - (class .== 46) * 0.86 .+ #GM - (class .== 70) * 0.77 .+ #WM - (class .== 93) * 1 .+ #FAT1 - (class .== 116) * 1 .+ #MUSCLE - (class .== 139) * 1 .+ #SKIN/MUSCLE - (class .== 162) * 0 .+ #SKULL - (class .== 185) * 0 .+ #VESSELS - (class .== 209) * 0.77 .+ #FAT2 - (class .== 232) * 1 .+ #DURA - (class .== 255) * 0.77 #MARROW - Δw_fat = -220 * 2π - Δw = - (class .== 93) * Δw_fat .+ #FAT1 - (class .== 209) * Δw_fat #FAT2 - T1 = T1 * 1e-3 - T2 = T2 * 1e-3 - T2s = T2s * 1e-3 + + # Get tissue properties + ρ, T1, T2, T2s, Δw = default_brain_tissue_properties(labels, tissue_properties) # Define and return the Phantom struct obj = Phantom{Float64}(; @@ -482,7 +410,6 @@ function pelvis_phantom2D(; ss=4, us=1) # subsample or upsample the phantom data class = repeat(data["pelvis3D_slice"][1:ssx:end, 1:ssy:end]; inner=[usx, usy]) - # Define spin position vectors Δx = .5e-3 * ssx / usx Δy = .5e-3 * ssy / usy @@ -612,5 +539,56 @@ function check_phantom_arguments(nd, ss, us) ssy = ss[2] end end + return ssx, ssy, ssz, usx, usy, usz end + +""" + ρ, T1, T2, T2s, Δw = default_brain_tissue_properties(labels, tissue_properties = nothing) + +This function returns the default brain tissue properties using a labels identifier Matrix +# Arguments +- `labels` : (`::Matrix`) the labels identifier matrix of the phantom +- `tissue_properties` : (`::Dict`, `=Dict()`) phantom tissue properties in ms and Hz considering the available tissues + +# Returns +- `ρ, T1, T2, T2s, Δw`: (`::Matrix`) matrices of the same size of labels with the tissues properties information + +# Examples +```julia-repl +julia> ρ, T1, T2, T2s, Δw = default_brain_tissue_properties(labels, tissue_properties) + +julia> ρ, T1, T2, T2s, Δw = default_brain_tissue_properties(labels) +``` +""" +function default_brain_tissue_properties(labels, tissue_properties = Dict()) + # Load default tissue properties + default_properties = Dict( + # ρ, T1, T2, T2*, Δw + "CSF" => [1, 2.569, 0.329, 0.058, 0], + "GM" => [0.86, 0.833, 0.083, 0.069, 0], + "WM" => [0.77, 0.500, 0.070, 0.061, 0], + "FAT1" => [1, 0.350, 0.070, 0.058, -220*2π], #-220 Hz + "MUSCLE" => [1, 0.900, 0.047, 0.030, 0], + "SKIN/MUSCLE" => [1, 0.569, 0.329, 0.058, 0], + "SKULL" => [0, 0, 0, 0, 0], + "VESSELS" => [0, 0, 0, 0, 0], + "FAT2" => [0.77, 0.500, 0.070, 0.061, -220*2π], #-220 Hz + "DURA" => [1, 2.569, 0.329, 0.058, 0], + "MARROW" => [0.77, 0.500, 0.070, 0.061, 0]) + + tissue_properties = merge(default_properties, tissue_properties) + props = ["ρ", "T1", "T2", "T2s", "Δw"] + Nproperties = length(props) + # Order: CSF, DURA, FAT1, FAT2, GM, MARROW, MUSCLE, SKIN/MUSCLE, SKULL, vESSELS, WM + tissue_labels = [23, 232, 93, 209, 46, 255, 116, 139, 162, 185, 70] + tissue_texts = sort(collect(keys(default_properties)))# + data_properties = zeros(Nproperties, size(labels)...) + for i=1:Nproperties + for (label, tissue) in zip(tissue_labels, tissue_texts) + data_properties[i, :, :, :] += (labels .== label)*tissue_properties[tissue][i] + end + end + + return (data_properties[i,:,:,:] for i in 1:Nproperties) +end