Laura

.gitignore

@@ -1,5 +1,5 @@
 *.jl.*.cov
 *.jl.cov
 *.jl.mem
-/Manifest.toml
-/.vscode
+Manifest.toml
+.vscode

Project.toml

@@ -4,16 +4,20 @@ authors = ["Laura Beghini"]
 version = "0.1.0"
 
 [deps]
+DSP = "717857b8-e6f2-59f4-9121-6e50c889abd2"
+DataInterpolations = "82cc6244-b520-54b8-b5a6-8a565e85f1d0"
 FileIO = "5789e2e9-d7fb-5bc7-8068-2c6fae9b9549"
 Images = "916415d5-f1e6-5110-898d-aaa5f9f070e0"
 JLD2 = "033835bb-8acc-5ee8-8aae-3f567f8a3819"
+MRIBase = "f7771a9a-6e57-4e71-863b-6e4b6a2f17df"
 MRICoilSensitivities = "c57eb701-aafc-44a2-a53c-128049758959"
 MRIFiles = "5a6f062f-bf45-497d-b654-ad17aae2a530"
 MRIReco = "bdf86e05-2d2b-5731-a332-f3fe1f9e047f"
 NIfTI = "a3a9e032-41b5-5fc4-967a-a6b7a19844d3"
 PolygonOps = "647866c9-e3ac-4575-94e7-e3d426903924"
 Scratch = "6c6a2e73-6563-6170-7368-637461726353"
 Setfield = "efcf1570-3423-57d1-acb7-fd33fddbac46"
+Statistics = "10745b16-79ce-11e8-11f9-7d13ad32a3b2"
 
 [compat]
 julia = "1.3"

README.md

@@ -2,3 +2,19 @@
 
 [![Build Status](https://github.com/Laura2305/MRINavigator.jl/actions/workflows/CI.yml/badge.svg?branch=main)](https://github.com/Laura2305/MRINavigator.jl/actions/workflows/CI.yml?query=branch%3Amain)
 [![Coverage](https://codecov.io/gh/Laura2305/MRINavigator.jl/branch/main/graph/badge.svg)](https://codecov.io/gh/Laura2305/MRINavigator.jl)
+
+MRINavigator.jl provides multiple navigator-based correction pipelines for magnetic resonance data. These aim at demodulating time dependent field variations. The package was developed with a focus on spinal cord imaging, however it can be used for multiple imaging applications. 
+The corrections are to be applied on the raw data before the image reconstruction. MRIReco.jl can be used to reconstruct the images.
+
+More details can be found in the Online Documentation (soon to be published).
+
+MRIReco reference: https://onlinelibrary.wiley.com/doi/epdf/10.1002/mrm.28792
+
+## How to give credit
+
+soon to be poublished
+
+
+## Community Standards
+
+This project is part of the Julia community and follows the [Julia community standards](https://julialang.org/community/standards/). 

src/AdjustData.jl

@@ -1,9 +1,13 @@
-export OrderSlices!, ExtractNoiseData!, ReverseBipolar!, RemoveRef!, CopyTE!, AdjustSubsampleIndices!
+export OrderSlices!, ExtractNoiseData!, ReverseBipolar!, RemoveRef!, CopyTE!, AdjustSubsampleIndices!, ExtractNavigator, ExtractFlags, selectEcho!, selectSlice!
 
 """
     OrderSlices!(rawData::RawAcquisitionData)
 
-Spatially order the slices in the raw data structure.
+Spatially order the slices in the MRIReco.jl raw data structure.
+The slices are ordered basing on the position coordinates saved in each profile.
+If these are not present the slices can not be ordered.
+
+MRIReco reference: https://onlinelibrary.wiley.com/doi/epdf/10.1002/mrm.28792
 
 # Arguments
 * `rawData::RawAcquisitionData` - raw data structure obtained loading raw data with MRIReco.jl
@@ -31,9 +35,11 @@ end
 """
     flags = ExtractFlags(rawData::RawAcquisitionData) 
 
-Extract the acquisition flags from raw data profiles.
+Extract the acquisition flags from the MRIReco.jl raw data profiles.
 Return a 31 elements vector for each profile.
 
+MRIReco reference: https://onlinelibrary.wiley.com/doi/epdf/10.1002/mrm.28792
+
 # Arguments
 * `rawData::RawAcquisitionData` - raw data structure obtained loading raw data with MRIReco.jl
 """
@@ -54,8 +60,11 @@ end
 """
     noisemat = ExtractNoiseData!(rawData::RawAcquisitionData, flags::Array{Int64})
 
-Extract and return the noise acquisition from the raw data.
-The noise acquisition is one of the profiles with slice = 0, contrast = 0, repetition = 0.
+Extract and return the noise acquisition from the MRIReco.jl raw data.
+The noise acquisition is usually the first profile with slice = 0, contrast = 0, repetition = 0.
+The noise profile should have the 19th flag element qual to 1. Check with ExtractFlags if errors occur.
+
+MRIReco reference: https://onlinelibrary.wiley.com/doi/epdf/10.1002/mrm.28792
 
 # Arguments
 * `rawData::RawAcquisitionData` - raw data structure obtained loading raw data with MRIReco.jl
@@ -87,7 +96,9 @@ end
 """
     ReverseBipolar!(rawData::RawAcquisitionData)
 
-Reflect the raw data profiles for bipolar acquisition.
+Reflect the MRIReco.jl raw data profiles for bipolar acquisition.
+
+MRIReco reference: https://onlinelibrary.wiley.com/doi/epdf/10.1002/mrm.28792
 
 # Arguments
 * `rawData::RawAcquisitionData` - raw data structure obtained loading raw data with MRIReco.jl
@@ -115,32 +126,18 @@ end
 """
     RemoveRef!(rawData::RawAcquisitionData, slices::Union{Vector{Int64}, Nothing}, echoes::Union{Vector{Int64}, Nothing})
 
-Remove reference data that are acquired with the phase stabilization on Siemens scanners.
+Remove reference data that are not useful for the navigator-based crrection from acquisitions with phase stabilization on Siemens scanners.
 Not solid to recalls.
 
+MRIReco reference: https://onlinelibrary.wiley.com/doi/epdf/10.1002/mrm.28792
+
 # Arguments
 * `rawData::RawAcquisitionData` - raw data structure obtained loading raw data with MRIReco.jl
-* `slices::Union{Vector{Int64}, Nothing}` - vector containing the numbers of slices to be loaded with MRIReco.jl. Nothing loads all.
-* `echoes::Union{Vector{Int64}, Nothing}` - vector containing the numbers of echoes to be loaded with MRIReco.jl. Nothing loads all.
 """
-function RemoveRef!(rawData::RawAcquisitionData, slices::Union{Vector{Int64}, Nothing}, echoes::Union{Vector{Int64}, Nothing})
+function RemoveRef!(rawData::RawAcquisitionData)
 
-    numSlices = 0
-    numEchoes = 0
-    if slices === nothing
-        numSlices = rawData.params["enc_lim_slice"].maximum+1
-    else
-        numSlices = size(slices, 1)
-    end
-    if echoes !== nothing
-        if 0 in echoes
-            numEchoes = size(echoes, 1) +1 # the navigator is saved as echo zero
-        else
-            numEchoes = size(echoes, 1)
-        end
-    else
-        numEchoes = size(rawData.params["TE"],1)+1
-    end
+    numSlices = rawData.params["enc_lim_slice"].maximum + 1
+    numEchoes = size(rawData.params["TE"],1) + 1
 
     #Apply this only if using phase stabilizaion
     removeIndx = numSlices*(numEchoes)
@@ -152,11 +149,13 @@ end
 """
     CopyTE!(rawData::RawAcquisitionData, acqData::AcquisitionData)
 
-Copy the TE values from the raw data structor to the acquisition data structor.
+Copy the TE values from the MRIReco.jl raw data structure to the acquisition data structure.
+
+MRIReco reference: https://onlinelibrary.wiley.com/doi/epdf/10.1002/mrm.28792
 
 # Arguments
 * `rawData::RawAcquisitionData` - raw data structure obtained loading raw data with MRIReco.jl
-* `acqData::RawAcquisitionData` - acquisition data structure obtained converting raw data with MRIReco.jl
+* `acqData::AcquisitionData` - acquisition data structure obtained converting raw data with MRIReco.jl
 """
 function CopyTE!(rawData::RawAcquisitionData, acqData::AcquisitionData)
 
@@ -170,11 +169,13 @@ end
 """
     AdjustSubsampleIndices!(acqData::AcquisitionData)
 
-Add subsamples indices in the acquisition data structure.
-Needed when conveting data not acquired at the first repetition.
+Add subsamples indices in the MRIReco.jl acquisition data structure.
+Needed when conveting data not acquired in the first repetition.
+
+MRIReco reference: https://onlinelibrary.wiley.com/doi/epdf/10.1002/mrm.28792
 
 # Arguments
-* `acqData::RawAcquisitionData` - acquisition data structure obtained converting raw data with MRIReco.jl
+* `acqData::AcquisitionData` - acquisition data structure obtained converting raw data with MRIReco.jl
 """
 function AdjustSubsampleIndices!(acqData::AcquisitionData)
 
@@ -184,4 +185,106 @@ function AdjustSubsampleIndices!(acqData::AcquisitionData)
         end
     end
 
+end
+
+
+"""
+    (nav, nav_time) = ExtractNavigator(rawData::RawAcquisitionData, slices::Union{Vector{Int64}, Nothing})
+
+Extract the navigator profiles from the MRIReco.jl raw data structure.
+These are registered with the same indices (contract, slice, encoding step) as the image data for the first echo time.
+Return a navigator array and a navigator time array. The navigator array has four dimensions in order: k-space samples, coils, k-space lines, slices.
+
+MRIReco reference: https://onlinelibrary.wiley.com/doi/epdf/10.1002/mrm.28792
+
+# Arguments
+* `rawData::RawAcquisitionData` - raw data structure obtained loading raw data with MRIReco.jl
+"""
+function ExtractNavigator(rawData::RawAcquisitionData)
+
+    total_num = length(rawData.profiles)
+    numberslices = rawData.params["enc_lim_slice"].maximum +1
+    contrasts = zeros(Int64, total_num)
+    slices = zeros(Int64, total_num)
+    lines = zeros(Int64, total_num)
+    for ii = 1:total_num
+        contrasts[ii] = rawData.profiles[ii].head.idx.contrast
+        slices[ii] = rawData.profiles[ii].head.idx.slice
+        lines[ii] = rawData.profiles[ii].head.idx.kspace_encode_step_1
+    end
+    # keep only the indexes of data saved in the first echo (this includes navigator)
+    contrastsIndx = findall(x->x==0, contrasts)
+    slices = slices[contrastsIndx]
+    lines = lines[contrastsIndx]
+
+    nav = zeros(ComplexF32, size(rawData.profiles[1].data)[1], size(rawData.profiles[1].data)[2],
+        rawData.params["reconSize"][2], numberslices)
+
+    nav_time = zeros(Float64,
+        rawData.params["reconSize"][2], numberslices)
+    #Odd indexes are data first echo, Even indexes are navigator data
+    for ii = 2:2:length(slices)
+        nav[:,:,lines[ii]+1,slices[ii]+1] = rawData.profiles[contrastsIndx[ii]].data
+        nav_time[lines[ii]+1,slices[ii]+1] = rawData.profiles[contrastsIndx[ii]].head.acquisition_time_stamp
+    end
+    #Remove the rows filled with zeroes
+    lines = unique(lines) .+1
+    nav = nav[:,:,lines,:]
+    nav_time = nav_time[lines,:]
+
+    return nav, nav_time
+    #navigator[k-space samples, coils, k-space lines, slices]
+
+end
+
+"""
+    SelectEcho!(acqd, idx_echo)
+
+Extract one or more echoes from the acquisition data structure
+
+# Arguments
+* `acqd::AcquisitionData` - acquisition data structure obtained converting raw data with MRIReco.jl
+* `idx_echo::Vector{Int64}` - vector containing the indexes of the echoes to be selected (starting from 0)
+
+MRIReco reference: https://onlinelibrary.wiley.com/doi/epdf/10.1002/mrm.28792
+"""
+function selectEcho!(acqd::AcquisitionData, idx_echo::Vector{Int64})
+
+    if !isempty(idx_echo)
+        contrasts = size(acqd.kdata)[1]
+        indices = Vector{Int64}(undef, contrasts)
+        for ii=1:contrasts
+            indices[ii] = ii
+        end
+        deleteat!(indices, idx_echo)
+        deleteat!(acqd.subsampleIndices, indices)
+        deleteat!(acqd.traj, indices)
+        acqd.kdata = acqd.kdata[idx_echo,:,:];
+    end
+
+end
+
+"""
+    SelectSlice!(acqd, nav, nav_time, idx_slice)
+
+Extract one or more echoes from the acquisition data structure
+
+# Arguments
+* `acqd::AcquisitionData` - acquisition data structure obtained converting raw data with MRIReco.jl
+* `idx_slice::Vector{Int64}` - vector containing the indexes of the slices to be selected (starting from 0, downer slice)
+
+# Optional arguments with default value = nothing
+* `nav::Union{Array{Complex{T}, 4}, Nothin} = nothing` - navigator profiles obtained with the ExtractNavigator function
+* `nav_time::Union{Array{Complex{Float32}, 2}, Nothing}` - time stamps for the navigator data obtained with ExtractNavigator (in ms from the beginning of the day)
+
+MRIReco reference: https://onlinelibrary.wiley.com/doi/epdf/10.1002/mrm.28792
+"""
+function selectSlice!(acqd::AcquisitionData, idx_slice::Vector{Int64}, nav::Union{Array{Complex{T}, 4}, Nothing} = nothing, nav_time::Union{Array{Float64, 2}, Nothing} = nothing) where {T}
+
+    acqd.kdata = acqd.kdata[:,idx_slice,:]
+    if !isnothing(nav) && !isnothing(nav_time)
+        nav = nav[:,:,:,idx_slice]
+        nav_time = nav_time[:,idx_slice]
+    end
+
 end