Merge pull request #4 from NordicMRspine/aj_code_review

Proposed formatting changes and Function redefinition

.github/workflows/CI.yml

@@ -18,7 +18,7 @@ jobs:
       fail-fast: false
       matrix:
         version:
-          - '1.1'
+          - '1.6'
           - '1.8'
           - 'nightly'
         os:

Project.toml

@@ -15,6 +15,7 @@ MRIFiles = "5a6f062f-bf45-497d-b654-ad17aae2a530"
 MRIReco = "bdf86e05-2d2b-5731-a332-f3fe1f9e047f"
 NIfTI = "a3a9e032-41b5-5fc4-967a-a6b7a19844d3"
 PolygonOps = "647866c9-e3ac-4575-94e7-e3d426903924"
+REPL = "3fa0cd96-eef1-5676-8a61-b3b8758bbffb"
 Scratch = "6c6a2e73-6563-6170-7368-637461726353"
 Setfield = "efcf1570-3423-57d1-acb7-fd33fddbac46"
 Statistics = "10745b16-79ce-11e8-11f9-7d13ad32a3b2"

src/AdjustData.jl

@@ -73,9 +73,11 @@ function ExtractNoiseData!(rawData::RawAcquisitionData)
 
     flags = ExtractFlags(rawData)
     total_num = length(rawData.profiles)
+
     if total_num != size(flags, 1)
         @error "size of flags and number of profiles in rawData do not match"
     end
+
     noisemat = Matrix{typeof(rawData.profiles[1].data)}
 
     for ii=1:total_num
@@ -107,6 +109,7 @@ function ReverseBipolar!(rawData::RawAcquisitionData)
 
     flags = ExtractFlags(rawData)
     total_num = length(rawData.profiles)
+
     if total_num != size(flags, 1)
         @error "size of flags and number of profiles in rawData do not match"
     end
@@ -180,9 +183,11 @@ MRIReco reference: https://onlinelibrary.wiley.com/doi/epdf/10.1002/mrm.28792
 function AdjustSubsampleIndices!(acqData::AcquisitionData)
 
     if isempty(acqData.subsampleIndices[1])
+
         for ii = 1:size(acqData.subsampleIndices)[1]
             acqData.subsampleIndices[ii]=1:size(acqData.kdata[1,1,1])[1]
         end
+
     end
 
 end
@@ -207,26 +212,32 @@ function ExtractNavigator(rawData::RawAcquisitionData)
     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 = 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)
 
-    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,:]
@@ -251,17 +262,20 @@ 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
+        acqd.kdata = acqd.kdata[idx_echo,:,:]
 
+    end
 end
 
 """
@@ -281,10 +295,14 @@ 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}
 
+    # get kdata from slice
     acqd.kdata = acqd.kdata[:,idx_slice,:]
+
     if !isnothing(nav) && !isnothing(nav_time)
+
         nav = nav[:,:,:,idx_slice]
         nav_time = nav_time[:,idx_slice]
+
     end
 
 end

src/CoilSensMap.jl

@@ -1,4 +1,4 @@
-export CompSensit, ResizeSensit, CompRoughMask
+export CompSensit, ResizeSensit!, CompRoughMask
 
 """
     sensit = CompSensit(acq::AcquisitionData, thresh = 0.135)
@@ -17,9 +17,9 @@ function CompSensit(acq::AcquisitionData, thresh = 0.135)
     sensit = espirit(acq,(6,6),30,eigThresh_1=0.02, eigThresh_2=0)
     slices = numSlices(acq)
     coils = numChannels(acq)
+
     # compute mask
     mask = CompRoughMask(acq, slices, thresh)
-
     for ii = 1:slices
 
         mask_slice = mask[:,:,ii]
@@ -55,6 +55,7 @@ function CompRoughMask(acq::AcquisitionData, slices::Int64, thresh)
     I_sum = sqrt.(sum(abs.(img) .^ 2, dims = 5)) .+ eps()
     I_sum = dropdims(I_sum, dims = tuple(findall(size(I_sum) .== 1)...))
     I_max = ones(Float64, slices)
+
     mask = zeros(size(I_sum))
     for ii = 1:slices
         I_max[ii] = maximum(abs.(I_sum[:,:,ii]))
@@ -97,16 +98,19 @@ Add a 3 voxels safety margin.
 function removeBehindBack!(mask_slice::Array{T,2}) where{T}
 
     # remove noisy voxels on the left of the image
-    # corresponding to the back of the subject
+    # corresponding to the back of the subject 
     density = sum(mask_slice, dims = 1)[1,:]
+
     # compute dervative of points density
     lines = div(size(mask_slice, 2),2)
     dder = zeros(Int64, lines)
     for ii = 1:lines
         dder[ii] = density[ii+1] - density[ii]
     end
+
     # put to zero everything behind the subject back
     position = findmax(dder)[2] -3
+
     for jj=1:position
         mask_slice[:,jj].=0
     end
@@ -126,7 +130,9 @@ function homogeneousMask!(mask_slice::Array{T,2}) where{T}
     cartes_index_slice = CartesianIndices(mask_slice)
     Bimask_slice = convert(BitMatrix, mask_slice)
     hull = convexhull(Bimask_slice)
+
     push!(hull, hull[1])
+
     inside = [inpolygon(p, hull; in=true, on=true, out=false) for p in cartes_index_slice]
     mask_slice[inside] .= 1
 
@@ -146,16 +152,18 @@ Image data and reference data must have the same slice center.
 
 MRIReco reference: https://onlinelibrary.wiley.com/doi/epdf/10.1002/mrm.28792
 """
-function ResizeSensit(sensit::Array{Complex{T},4}, acqMap::AcquisitionData, acqData::AcquisitionData) where {T}
+function ResizeSensit!(sensit::Array{Complex{T},4}, acqMap::AcquisitionData, acqData::AcquisitionData) where {T}
 
     # Define the relevant sensit region assuming the same slices center between ref and image data
     (freq_enc_FoV, freq_enc_samples, phase_enc_FoV, phase_enc_samples) = Find_scaling_sensit(acqMap, acqData)
     sizeSensit = size(sensit)
 
     if freq_enc_samples[1] != sizeSensit[1] && freq_enc_samples[2] != sizeSensit[2]
         @warn "The coils sensitivity maps have already been resized, the function cannot be executed."
+
     elseif freq_enc_FoV[1] < freq_enc_FoV[2] || phase_enc_FoV[1] < phase_enc_FoV[2]
         @error "The reference data field of view is smaller than the image data field of view."
+
     else
         freq_enc_FoV_disc = round(Int64, (freq_enc_FoV[1] - freq_enc_FoV[2]) / (freq_enc_FoV[1]/freq_enc_samples[1]) / 2)
         phase_enc_FoV_disc = round(Int64, (phase_enc_FoV[1] - phase_enc_FoV[2]) / (phase_enc_FoV[1]/phase_enc_samples[1]) / 2)
@@ -169,15 +177,19 @@ function ResizeSensit(sensit::Array{Complex{T},4}, acqMap::AcquisitionData, acqD
         for ii in cartes_index
             mask[ii] = 1
         end
+
         # Linear interpolation
         sensit = mapslices(x ->imresize(x, (freq_enc_samples[2], phase_enc_samples[2])), sensit, dims=[1,2])
         mask = mapslices(x ->imresize(x, (freq_enc_samples[2], phase_enc_samples[2])), mask, dims=[1,2])
+
         # Remove interpolated outline
         cartes_index = findall(x -> x!=1, mask)
         for ii in cartes_index
             mask[ii] = 0
         end
+
         sensit = mask .* sensit
+
     end
 end
 
@@ -196,7 +208,7 @@ Image data and reference data must have the same slice center.
 
 MRIReco reference: https://onlinelibrary.wiley.com/doi/epdf/10.1002/mrm.28792
 """
-function Find_scaling_sensit(acqMap::AcquisitionData, acqData::AcquisitionData) where {T}
+function Find_scaling_sensit(acqMap::AcquisitionData{T}, acqData::AcquisitionData{T}) where {T}
 
     freq_enc_FoV = [acqMap.fov[1], acqData.fov[1]]
     freq_enc_samples = [acqMap.encodingSize[1], acqData.encodingSize[1]]

src/NavData.jl

@@ -1,6 +1,7 @@
 export additionalNavInput, navOutput
 
 mutable struct additionalNavInput
+
     numslices::Int64
     numechoes::Int64
     numsamples::Int64
@@ -15,6 +16,7 @@ mutable struct additionalNavInput
     noisemat::Array{Complex{Float32}, 2}
     trace::Union{Matrix{Float64}, Nothing}
     centerline::Union{Vector{Float64}, Nothing}
+
 end
 
 """
@@ -57,14 +59,18 @@ function additionalNavInput(
     numsamples = acqData.encodingSize[1]
     numlines = convert(Int64, size(acqData.kdata[1],1)/numsamples)
     TR = rawData.params["TR"]
+
     ii=1
     while rawData.profiles[ii].head.user_int[8] < rawData.profiles[ii+1].head.user_int[8]
         ii=ii+1
     end
+
+    # set up nav timing
     dt_nav = convert(Float64, rawData.profiles[ii-1].head.sample_time_us) .* 1e-6
     TE_nav = rawData.profiles[ii].head.user_int[8] .* 1e-6 # get TE nav
+
     if !isnothing(acqMap) && !isnothing(acqData)
-    (freq_enc_FoV, freq_enc_samples, phase_enc_FoV, phase_enc_samples) = Find_scaling_sensit(acqMap, acqData)
+        (freq_enc_FoV, freq_enc_samples, phase_enc_FoV, phase_enc_samples) = Find_scaling_sensit(acqMap, acqData)
     end
 
     return additionalNavInput(numslices, numechoes, numsamples, numlines, TR, TE_nav, dt_nav,
@@ -73,8 +79,10 @@ function additionalNavInput(
 end
 
 mutable struct navOutput
+
     navigator::Array{Float64, 4}
     centerline::Union{Array{Float64, 1}, Nothing}
     correlation::Union{Array{Float64, 1}, Matrix{Float64}, Nothing}
     wrapped_points::Union{Array{Int8, 2}, Nothing}
+
 end

src/NavParameters.jl

@@ -35,6 +35,7 @@ SCT reference: https://spinalcordtoolbox.com
 
 """
 function defaultNavParams()
+
     params = Dict{Symbol,Any}()
     params[:slices] = nothing
     params[:echoes] = nothing
@@ -45,8 +46,7 @@ function defaultNavParams()
     params[:trust_SCT] = false
     params[:use_SCT] = false
     params[:corr_type] = "FFT"
-    params[:FFT_interval] = 35 # millimiters
-
+    params[:FFT_interval] = 35 # [millimiters]
 
     return params
   end