From efecfe6be27c12ac67ead02ec4710f5e08c6b26d Mon Sep 17 00:00:00 2001
From: Tobias Knopp <tobias@knoppweb.de>
Date: Fri, 11 Dec 2020 17:28:50 +0100
Subject: [PATCH] add mridataorg example

---
 examples/mridataorg/README.md       | 39 +++++++++++++
 examples/mridataorg/downloadData.jl | 17 ++++++
 examples/mridataorg/example.jl      | 91 +++++++++++++++++++++++++++++
 3 files changed, 147 insertions(+)
 create mode 100644 examples/mridataorg/README.md
 create mode 100644 examples/mridataorg/downloadData.jl
 create mode 100644 examples/mridataorg/example.jl

diff --git a/examples/mridataorg/README.md b/examples/mridataorg/README.md
new file mode 100644
index 00000000..552b65ee
--- /dev/null
+++ b/examples/mridataorg/README.md
@@ -0,0 +1,39 @@
+# Example of mridata.org reconstruction
+
+This folder contains example code for downloading data from the website mridata.org 
+and performing image reconstruction of the undersampled FSE data. The data is stored
+in the ISMRMRD file format
+
+## Installation
+
+In order to use this code one first has to download Julia (version 1.5) and install
+`MRIReco` by executing
+
+```julia
+using Pkg
+Pkg.add("MPIReco")
+```
+
+Load the package by entering
+```julia
+using MRIReco
+```
+
+You can then switch to the directory of this example by entering
+```julia
+dir = joinpath(dirname(pathof(MRIReco)), "..","examples","mridataorg")
+cd(dir)
+```
+
+## Execution
+After installation and switching to the example directory, the example code can be
+executed by entering
+
+```julia
+include("example.jl")
+```
+
+This will first download the ISMRMRD data (about 167 MB) and then perform a reconstruction.
+Parameters of the reconstruction are documented in the Julia script and can be
+changed. After the reconstruction is done, the script opens a PyPlot window
+and show the reconstruction result. 
diff --git a/examples/mridataorg/downloadData.jl b/examples/mridataorg/downloadData.jl
new file mode 100644
index 00000000..ac0ec5d6
--- /dev/null
+++ b/examples/mridataorg/downloadData.jl
@@ -0,0 +1,17 @@
+using HTTP
+
+function download_(filenameServer, filenameLocal)
+  if !isfile(filenameLocal)
+    @info "download $(filenameLocal)..."
+    HTTP.open("GET", filenameServer) do http
+      open(filenameLocal, "w") do file
+        write(file, http)
+      end
+    end
+  end
+end
+
+mkpath("./data/")
+
+# Download data
+download_("http://mridata.org/download/fd28ac5d-ed5a-462f-9064-0a5a13d76291", "./data/ksp_knee_3dfse.h5")
\ No newline at end of file
diff --git a/examples/mridataorg/example.jl b/examples/mridataorg/example.jl
new file mode 100644
index 00000000..9a5da86d
--- /dev/null
+++ b/examples/mridataorg/example.jl
@@ -0,0 +1,91 @@
+using Pkg
+
+isinstalled(pkg::String) = any(x -> x.name == pkg && x.is_direct_dep, values(Pkg.dependencies()))
+
+# Install required packages
+for P in ["HTTP", "PyPlot"]
+  !isinstalled(P) && Pkg.add(P)
+end
+
+# Download data
+include("downloadData.jl")
+
+using PyPlot, MRIReco
+
+##################
+## Data Loading ##
+##################
+f = ISMRMRDFile("data/ksp_knee_3dfse.h5")
+acqRaw = RawAcquisitionData(f)
+
+# kspace-center is not specified in the raw data file.
+# So we let MRIReco estimate it
+acqData = AcquisitionData(acqRaw; estimateProfileCenter=true)
+
+# convert to 2d
+acqData2d = convert3dTo2d(acqData)
+
+# extract slices
+sl = [50,100,150,200]
+acqData2d.kdata = acqData2d.kdata[:,sl,:]
+acqData2d.encodingSize = [274,208,1]
+
+####################################
+## generate coil sensitivity maps ##
+####################################
+@info "Espirit"
+smaps = espirit(acqData2d, (6,6), 30, eigThresh_1=0.04, eigThresh_2=0.98)
+
+
+###############################
+## Reconstruction Parameters ##
+###############################
+params = Dict{Symbol, Any}()
+params[:reco] = "multiCoil"
+params[:reconSize] = Tuple(acqData2d.encodingSize[1:2])
+params[:solver] = "admm"
+params[:regularization] = "L1"
+params[:sparseTrafo] = "Wavelet"
+params[:λ] = 2.e-1
+params[:iterations] = 30
+params[:ρ] = 0.1
+params[:absTol] = 1.e-4
+params[:relTol] = 1.e-3
+params[:tolInner] = 1.e-2
+params[:senseMaps] = smaps
+params[:normalizeReg] = true
+
+############################
+## Perform reconstruction ##
+############################
+
+@info "Perform Reconstruction "
+@time img = fftshift( reconstruction(acqData2d, params).data[:,:,:,1,1], 2)
+
+###########################
+## Visualize the results ##
+###########################
+
+figure(1, figsize=(6,1.9))
+
+for s in [1,2,3,4]
+  subplot(1,4,s)
+  imshow(abs.(img[:,:,s,1,1]))
+  title("slice $(sl[s])")
+  axis("off")
+end
+
+subplots_adjust(left=0.01,bottom=0.01,wspace=0.3,hspace=0.2,right=0.99,top=0.92)
+
+
+
+
+
+
+
+
+
+
+
+
+