From a4263ad2787bcc386703905105448f078dea8ce6 Mon Sep 17 00:00:00 2001 From: spinicist Date: Wed, 24 Jun 2020 11:33:12 +0100 Subject: [PATCH] BLD Minor polish --- Python/T1_mapping.ipynb | 237 +++++++++++++++++ Python/mupa_notebook.ipynb | 395 ----------------------------- Python/steady_state.ipynb | 265 ------------------- Source/PARMESAN/transient_main.cpp | 3 +- Source/Relaxometry/qimcdespot.cpp | 4 +- Source/qi_main.cpp | 11 +- 6 files changed, 246 insertions(+), 669 deletions(-) create mode 100644 Python/T1_mapping.ipynb delete mode 100644 Python/mupa_notebook.ipynb delete mode 100644 Python/steady_state.ipynb diff --git a/Python/T1_mapping.ipynb b/Python/T1_mapping.ipynb new file mode 100644 index 00000000..0f73b72b --- /dev/null +++ b/Python/T1_mapping.ipynb @@ -0,0 +1,237 @@ +{ + "cells": [ + { + "cell_type": "markdown", + "metadata": {}, + "source": [ + "# T1-Mapping with DESPOT1 / VFA\n", + "\n", + "This notebook deomstrates how T1 mapping with a basic `nipype` pipeline with `QUIT` programs.\n", + "\n", + "This pipeline (http://github.com/spinicist/nanslice) downloads the BrainWeb brain & B1 phantoms from http://brainweb.bic.mni.mcgill.ca. It then replaces the tissue classification labels with values of Proton Density and T1, simulates an SPGR/FLASH image with some added noise, and finally uses that simulated data to fit for T1 and PD again.\n", + "\n", + "In order to display the images, and to do some " + ] + }, + { + "cell_type": "markdown", + "metadata": {}, + "source": [ + "## Imports" + ] + }, + { + "cell_type": "code", + "metadata": {}, + "outputs": [], + "source": [ + "%matplotlib inline\n", + "from QUIT.interfaces.relax import DESPOT1, DESPOT1Sim, DESPOT2\n", + "from nanslice import Layer\n", + "import nanslice.jupyter as ns\n", + "import nibabel as nib\n", + "import numpy as np\n", + "import requests\n", + "import gzip\n", + "import os.path" + ] + }, + { + "cell_type": "code", + "execution_count": null, + "metadata": {}, + "outputs": [], + "source": [ + "## Create Phantom Data\n", + "\n", + "This section downloads the Brainweb phantoms. These are stored MINC format so we load the images using `nibabel` to get the raw data later." + ] + }, + { + "cell_type": "markdown", + "metadata": {}, + "source": [ + "# Data Fetching\n", + "\n", + "We now download the Brainweb brain & T1 phantoms, if they are not already present in the current directory. These are stored in MINC format, so we use " + ] + }, + { + "cell_type": "markdown", + "metadata": {}, + "source": [ + "Display the phantom data as a check" + ] + }, + { + "cell_type": "code", + "metadata": {}, + "outputs": [], + "source": [ + "if not isfile('classes.mnc'):\n", + " print('Downloading classes')\n", + " params = {'download_for_real':'[Start Download!]',\n", + " 'do_download_alias':'phantom_1.0mm_normal_crisp',\n", + " 'format_value':'minc',\n", + " 'who_name': 'Tobias Wood',\n", + " 'who_institution': 'KCL',\n", + " 'who_email': 'tobias.wood@kcl.ac.uk'}\n", + " response = requests.get(url='http://brainweb.bic.mni.mcgill.ca/cgi/brainweb1', params=params)\n", + " minc_file = open('classes.mnc', 'wb')\n", + " minc_file.write(response.content)\n", + "if not isfile('rf20_C.mnc'):\n", + " print('Downloading B1')\n", + " params = {'download_for_real':'[Start Download!]',\n", + " 'do_download_alias':'rf20_C.mnc.gz',\n", + " 'format_value':'minc',\n", + " 'zip_value':'none',\n", + " 'who_name': 'Tobias Wood',\n", + " 'who_institution': 'KCL',\n", + " 'who_email': 'tobias.wood@kcl.ac.uk'}\n", + " response = requests.get(url='https://brainweb.bic.mni.mcgill.ca/cgi/brainweb1', params=params)\n", + " b1_file = open('rf20_C.mnc', 'wb')\n", + " b1_file.write(response.content)\n", + "classes = nanslice.Layer('classes.mnc')\n", + "b0_minc = nanslice.Layer('rf20_B.mnc')\n", + "b1_minc = nanslice.Layer('rf20_C.mnc')" + ] + }, + { + "cell_type": "code", + "metadata": {}, + "outputs": [], + "source": [ + "ns.three_plane(classes)" + ] + }, + { + "cell_type": "markdown", + "metadata": {}, + "source": [ + "Now we define a function that allows us to convert the tissue classes to T1 and PD/M0 values, and use it to create our reference T1, PD and B1 images. These will be saved into the current directory. This function also allows us to subsample the images, which is useful when investigating methods that take longer to fit.\n", + "\n", + "The BrainWeb phantom uses the following mapping:\n", + "\n", + "0. Background\n", + "1. CSF\n", + "2. Grey Matter\n", + "3. White Matter\n", + "4. Fat\n", + "5. Muscle/Skin\n", + "6. Skin\n", + "7. Skull\n", + "8. Glial Matter\n", + "9. Connective\n", + "\n", + "To keep things simple, we only specify values for CSF, GM & WM, and set the other tissue types to zero. The B1 data is used as is from the phantom." + ] + }, + { + "cell_type": "code", + "metadata": {}, + "outputs": [], + "source": [ + "def classes_to_params(M0Vals, T1Vals, subsamp=1):\n", + " class_data = classes.image.get_data().astype('int32')\n", + " M0data = np.choose(class_data[::subsamp,::subsamp,::subsamp], M0vals).astype('float32')\n", + " T1data = np.choose(class_data[::subsamp,::subsamp,::subsamp], T1vals).astype('float32')\n", + " T2data = np.choose(class_data[::subsamp,::subsamp,::subsamp], T2vals).astype('float32')\n", + " B1data = b1_minc.image.get_data().astype('float32')[::subsamp,::subsamp,::subsamp]\n", + " # PDdata = np.array(list(map(PDFunc, classes.image.get_data())))\n", + " M0image = nib.nifti1.Nifti1Image(M0data, affine=classes.image.affine)\n", + " T1image = nib.nifti1.Nifti1Image(T1data, affine=classes.image.affine)\n", + " B1image = nib.nifti1.Nifti1Image(B1data, affine=classes.image.affine)\n", + " nib.save(M0image, 'M0.nii.gz')\n", + " nib.save(T1image, 'T1.nii.gz')\n", + " nib.save(B1image, 'B1.nii.gz')\n", + "\n", + "M0vals = np.array([0, 1, 0.8, 0.7, 0, 0, 0, 0, 0, 0])\n", + "T1vals = np.array([0, 3.0, 1.3, 0.9, 0, 0, 0, 0, 0, 0])\n", + "classes_to_params(M0vals, T1vals, 1)" + ] + }, + { + "cell_type": "markdown", + "metadata": {}, + "source": [ + "## Simulate Image\n", + "\n", + "As well as fitting image data to find parameters, QUIT allows you to simulate the images from the paramters. Here we use the phantom images from the step above to create some simple gradient-echo (also called FLASH, SPGR or FFE depending on vendor) images.\n", + "\n", + "First we define a dictionary that sets the sequence parameters for the type of scan we are simulating. In this case we only need the TR and the flip-angle (FA). The TR is specified in seconds, not milliseconds, and there are multiple values for the flip-angle as we need to simulate multiple images.\n", + "\n", + "We then define and run the simulator. This will return a results object that contains the paths to the output files, which we use to display both volumes in the output image. The noise is defined as an SNR level relative assuming a nominal PD of 1 (this is why the PD of CSF was set to 1 above)." + ] + }, + { + "cell_type": "code", + "metadata": {}, + "outputs": [], + "source": [ + "sequence_parameters = {'SPGR': {'TR': 10e-3, 'FA': [3, 18]}}\n", + "simulator_results = DESPOT1Sim(sequence=d1seq, in_file='sim_spgr.nii.gz', noise=0.001, PD='PD.nii.gz', T1='T1.nii.gz').run()\n", + "\n", + "simulated_spgr = Layer(simulator_results.outputs.out_file, volume=0)\n", + "display(ns.three_plane(spgr))\n", + "spgr.volume = 1\n", + "display(ns.three_plane(spgr))" + ] + }, + { + "cell_type": "markdown", + "metadata": {}, + "source": [ + "## Fit Data and Compare to Reference\n", + "\n", + "Now we have simulated data we can run the DESPOT1/VFA fitting code and compare the results to our phantom reference data." + ] + }, + { + "cell_type": "markdown", + "metadata": {}, + "source": [ + "## B1 Mapping\n", + "\n", + "SPGR images are affected by B1 inhomogeneity. In the above simulation, we did not specify a B1 map and assumed that B1 was flat." + ] + }, + { + "cell_type": "code", + "metadata": {}, + "outputs": [], + "source": [ + "fitting_results = DESPOT1(sequence=d1seq, in_file='sim_spgr.nii.gz').run()\n", + "\n", + "display(ns.compare('T1.nii.gz', fitting_results.outputs.t1_map, title='T1 Comparison'))\n", + "display(ns.compare('PD.nii.gz', fitting_results.outputs.pd_map, title='PD Comparison'))" + ] + } + ], + "metadata": { + "kernel_info": { + "name": "python3" + }, + "kernelspec": { + "display_name": "Python 3", + "language": "python", + "name": "python3" + }, + "language_info": { + "codemirror_mode": { + "name": "ipython", + "version": 3 + }, + "file_extension": ".py", + "mimetype": "text/x-python", + "name": "python", + "nbconvert_exporter": "python", + "pygments_lexer": "ipython3", + "version": "3.7.4" + }, + "nteract": { + "version": "0.8.4" + } + }, + "nbformat": 4, + "nbformat_minor": 4 +} \ No newline at end of file diff --git a/Python/mupa_notebook.ipynb b/Python/mupa_notebook.ipynb deleted file mode 100644 index 0d968309..00000000 --- a/Python/mupa_notebook.ipynb +++ /dev/null @@ -1,395 +0,0 @@ -{ - "cells": [ - { - "cell_type": "markdown", - "execution_count": null, - "metadata": {}, - "outputs": [], - "source": [ - "## Imports and Data Fetching" - ] - }, - { - "cell_type": "code", - "execution_count": null, - "metadata": {}, - "outputs": [], - "source": [ - "%matplotlib inline\n", - "%load_ext autoreload\n", - "%autoreload 1\n", - "%aimport nanslice\n", - "%aimport nanslice.jupyter\n", - "%aimport QUIT.interfaces.rufis\n", - "import nibabel as nib\n", - "import numpy as np\n", - "import requests\n", - "import gzip\n", - "from os.path import isfile" - ] - }, - { - "cell_type": "markdown", - "execution_count": null, - "metadata": {}, - "outputs": [], - "source": [ - "## Setup" - ] - }, - { - "cell_type": "code", - "execution_count": null, - "metadata": {}, - "outputs": [], - "source": [ - "if not isfile('classes.mnc'):\n", - " print('Downloading classes')\n", - " params = {'download_for_real':'[Start Download!]',\n", - " 'do_download_alias':'phantom_1.0mm_normal_crisp',\n", - " 'format_value':'minc',\n", - " 'who_name': 'Tobias Wood',\n", - " 'who_institution': 'KCL',\n", - " 'who_email': 'tobias.wood@kcl.ac.uk'}\n", - " response = requests.get(url='http://brainweb.bic.mni.mcgill.ca/cgi/brainweb1', params=params)\n", - " minc_file = open('classes.mnc', 'wb')\n", - " minc_file.write(response.content)\n", - "if not isfile('rf20_C.mnc'):\n", - " print('Downloading B1')\n", - " params = {'download_for_real':'[Start Download!]',\n", - " 'do_download_alias':'rf20_C.mnc.gz',\n", - " 'format_value':'minc',\n", - " 'zip_value':'none',\n", - " 'who_name': 'Tobias Wood',\n", - " 'who_institution': 'KCL',\n", - " 'who_email': 'tobias.wood@kcl.ac.uk'}\n", - " response = requests.get(url='https://brainweb.bic.mni.mcgill.ca/cgi/brainweb1', params=params)\n", - " b1_file = open('rf20_C.mnc', 'wb')\n", - " b1_file.write(response.content)\n", - "classes = nanslice.Layer('classes.mnc')\n", - "b1_minc = nanslice.Layer('rf20_C.mnc')" - ] - }, - { - "cell_type": "code", - "execution_count": null, - "metadata": {}, - "outputs": [], - "source": [ - "# 0=Background\n", - "# 1=CSF\n", - "# 2=Grey Matter\n", - "# 3=White Matter\n", - "# 4=Fat\n", - "# 5=Muscle/Skin\n", - "# 6=Skin\n", - "# 7=Skull\n", - "# 8=Glial Matter\n", - "# 9=Connective\n", - "def classes_to_params(M0Vals, T1Vals, T2vals, subsamp=1):\n", - " class_data = classes.image.get_data().astype('int32')\n", - " M0data = np.choose(class_data[::subsamp,::subsamp,::subsamp], M0vals).astype('float32')\n", - " T1data = np.choose(class_data[::subsamp,::subsamp,::subsamp], T1vals).astype('float32')\n", - " T2data = np.choose(class_data[::subsamp,::subsamp,::subsamp], T2vals).astype('float32')\n", - " B1data = b1_minc.image.get_data().astype('float32')[::subsamp,::subsamp,::subsamp]\n", - " # PDdata = np.array(list(map(PDFunc, classes.image.get_data())))\n", - " M0image = nib.nifti1.Nifti1Image(M0data, affine=classes.image.affine)\n", - " T1image = nib.nifti1.Nifti1Image(T1data, affine=classes.image.affine)\n", - " T2image = nib.nifti1.Nifti1Image(T2data, affine=classes.image.affine)\n", - " B1image = nib.nifti1.Nifti1Image(B1data, affine=classes.image.affine)\n", - " nib.save(M0image, 'M0.nii.gz')\n", - " nib.save(T1image, 'T1.nii.gz')\n", - " nib.save(T2image, 'T2.nii.gz')\n", - " nib.save(B1image, 'B1.nii.gz')" - ] - }, - { - "cell_type": "code", - "execution_count": null, - "metadata": {}, - "outputs": [], - "source": [ - "prep_pulses = {\n", - " 'inv': {'FAeff': 175.3, 'T_long': 0.0365, 'T_trans': 0.004, 'int_b1_sq': 53390.0},\n", - " 't2-20': {'FAeff': 0.07, 'T_long': 0.0055, 'T_trans': 0.020, 'int_b1_sq': 279600.0},\n", - " 't2-40': {'FAeff': 0.07, 'T_long': 0.0055, 'T_trans': 0.040, 'int_b1_sq': 279600.0},\n", - " 't2-60': {'FAeff': 0.07, 'T_long': 0.0055, 'T_trans': 0.060, 'int_b1_sq': 279600.0},\n", - " 't2-80': {'FAeff': 0.07, 'T_long': 0.0055, 'T_trans': 0.080, 'int_b1_sq': 279600.0},\n", - " 't2-120': {'FAeff': 0.07, 'T_long': 0.0055, 'T_trans': 0.120, 'int_b1_sq': 279600.0},\n", - " 't2-160': {'FAeff': 0.07, 'T_long': 0.0055, 'T_trans': 0.160, 'int_b1_sq': 279600.0},\n", - " 't2-inv': {'FAeff': 180., 'T_long': 0.0055, 'T_trans': 0.020, 'int_b1_sq': 279600.0},\n", - " 't1-20': {'FAeff': 20, 'T_long': 0.036487, 'T_trans': 0.004396, 'int_b1_sq': 0},\n", - " 't1-45': {'FAeff': 45, 'T_long': 0.036487, 'T_trans': 0.004396, 'int_b1_sq': 0},\n", - " 't1-60': {'FAeff': 60, 'T_long': 0.036487, 'T_trans': 0.004396, 'int_b1_sq': 0},\n", - " 't1-75': {'FAeff': 75, 'T_long': 0.036487, 'T_trans': 0.004396, 'int_b1_sq': 0},\n", - " 't1-90': {'FAeff': 90, 'T_long': 0.036487, 'T_trans': 0.004396, 'int_b1_sq': 0},\n", - " 't1-120': {'FAeff': 120, 'T_long': 0.036487, 'T_trans': 0.004396, 'int_b1_sq': 0},\n", - " 't1-180': {'FAeff': 180, 'T_long': 0.036487, 'T_trans': 0.004396, 'int_b1_sq': 0},\n", - " 't1-360': {'FAeff': 360, 'T_long': 0.03, 'T_trans': 0.004, 'int_b1_sq': 0},\n", - " 'null': {'FAeff': 0.0, 'T_long': 0.0, 'T_trans': 0.0, 'int_b1_sq': 0.0},\n", - " 'delay-400': {'FAeff': 0.0, 'T_long': 0.4, 'T_trans': 0.0, 'int_b1_sq': 0.0},\n", - " 'delay-800': {'FAeff': 0.0, 'T_long': 0.8, 'T_trans': 0.0, 'int_b1_sq': 0.0}\n", - "}" - ] - }, - { - "cell_type": "markdown", - "execution_count": null, - "metadata": {}, - "outputs": [], - "source": [ - "## Simulate and Display Image" - ] - }, - { - "cell_type": "code", - "execution_count": null, - "metadata": {}, - "outputs": [], - "source": [ - "M0vals = np.array([0, 1, 0.8, 0.7, 0, 0, 0, 0, 0, 0])\n", - "T1vals = np.array([0, 3.0, 1.3, 0.9, 0, 0, 0, 0, 0, 0])\n", - "T2vals = np.array([0, 0.5, 0.08, 0.06, 0, 0, 0, 0, 0, 0])\n", - "classes_to_params(M0vals, T1vals, T2vals, 2)" - ] - }, - { - "cell_type": "code", - "execution_count": null, - "metadata": {}, - "outputs": [], - "source": [ - "seq = {\n", - " \"RUFIS\": {\n", - " \"TR\": 0.002012,\n", - " \"Tramp\": 0.01,\n", - " \"Tspoil\": 0.05,\n", - " \"Trf\": [28e-6, 28e-6, 28e-6],\n", - " \"spokes_per_seg\": 48,\n", - " \"groups_per_seg\": [1, 1, 1],\n", - " \"FA\": [2, 8, 8],\n", - " \"prep\": [\"null\", \"null\", \"t1-360\"],\n", - " \"prep_pulses\": prep_pulses\n", - " }\n", - "}\n", - "region='45,0,0,1,109,91'\n", - "sim_result = QUIT.interfaces.rufis.SteadyStateSim(sequence=seq, in_file='sim.nii.gz', noise=0.0005, M0='M0.nii.gz', T1='T1.nii.gz', B1='B1.nii.gz', mask='M0.nii.gz', subregion=region).run()\n", - "display(nanslice.jupyter.timeseries(sim_result.outputs.out_file, clim=[0,0.025]))\n", - "mupa = QUIT.interfaces.rufis.SteadyState(sequence=seq, in_file='sim.nii.gz', mask_file='M0.nii.gz', subregion=region).run()\n", - "display(nanslice.jupyter.compare('M0.nii.gz', 'RUFIS_SS_M0.nii.gz', axis='z', clim=[0, 1.25], title='M0'))\n", - "display(nanslice.jupyter.compare('T1.nii.gz', 'RUFIS_SS_T1.nii.gz', axis='z', clim=[0.5, 1.5], title='T1'))\n", - "display(nanslice.jupyter.compare('B1.nii.gz', 'RUFIS_SS_B1.nii.gz', axis='z', clim=[0.5, 1.5], title='B1'))" - ] - }, - { - "cell_type": "markdown", - "execution_count": null, - "metadata": {}, - "outputs": [], - "source": [ - "## Emil's November Data" - ] - }, - { - "cell_type": "code", - "execution_count": null, - "metadata": {}, - "outputs": [], - "source": [ - "seq = {\n", - " 'MUPA': {\n", - " 'TR': 2.34e-3,\n", - " 'Tramp': 10e-3,\n", - " 'spokes_per_seg': 256,\n", - " 'groups_per_seg': [1, 1, 1, 1, 1, 1, 1],\n", - " 'FA': [2, 2, 2, 2, 2, 2, 2],\n", - " 'Trf': [24, 24, 24, 24, 24, 24, 24],\n", - " 'prep': ['inv', 'null', 'null', 'null', 't2-40', 't2-80', 't2-160'],\n", - " 'prep_pulses': prep_pulses}\n", - "}\n", - "region='45,0,0,1,109,91'\n", - "sim_result = QUIT.interfaces.rufis.MUPASim(sequence=seq, in_file='sim.nii.gz', noise=0.0005, M0='M0.nii.gz', T1='T1.nii.gz', T2='T2.nii.gz', B1='B1.nii.gz', mask='M0.nii.gz', subregion=region).run()\n", - "display(nanslice.jupyter.timeseries(sim_result.outputs.out_file, clim=[-0.01,0.02]))\n", - "mupa = QUIT.interfaces.rufis.MUPA(sequence=seq, in_file='sim.nii.gz', mask_file='M0.nii.gz', subregion=region).run()\n", - "display(nanslice.jupyter.compare('M0.nii.gz', 'MUPA_M0.nii.gz', axis='z', clim=[0, 1.25], title='M0'))\n", - "display(nanslice.jupyter.compare('T1.nii.gz', 'MUPA_T1.nii.gz', axis='z', clim=[0.5, 1.5], title='T1'))\n", - "display(nanslice.jupyter.compare('T2.nii.gz', 'MUPA_T2.nii.gz', axis='z', clim=[0.0, 0.15], title='T2'))" - ] - }, - { - "cell_type": "markdown", - "execution_count": null, - "metadata": {}, - "outputs": [], - "source": [ - "## 4 Point" - ] - }, - { - "cell_type": "code", - "execution_count": null, - "metadata": {}, - "outputs": [], - "source": [ - "seq = {\n", - " 'MUPA': {\n", - " 'TR': 2.34e-3,\n", - " 'Tramp': 10e-3,\n", - " 'spokes_per_seg': 512,\n", - " 'groups_per_seg': [1, 1, 1, 1],\n", - " 'FA': [2, 2, 2, 2],\n", - " 'Trf': [24, 24, 24, 24],\n", - " 'prep': ['inv', 'delay-800', 't2-80', 'delay-800'],\n", - " 'prep_pulses': prep_pulses}\n", - "}\n", - "region='45,0,0,1,109,91'\n", - "sim_result = QUIT.interfaces.rufis.MUPASim(sequence=seq, in_file='sim.nii.gz', noise=0.0005, M0='M0.nii.gz', T1='T1.nii.gz', T2='T2.nii.gz', B1='B1.nii.gz', mask='M0.nii.gz', subregion=region).run()\n", - "display(nanslice.jupyter.timeseries(sim_result.outputs.out_file, clim=[-0.01,0.025]))\n", - "mupa = QUIT.interfaces.rufis.MUPA(sequence=seq, in_file='sim.nii.gz', mask_file='M0.nii.gz', subregion=region).run()\n", - "display(nanslice.jupyter.compare('M0.nii.gz', 'MUPA_M0.nii.gz', axis='z', clim=[0, 1.25], title='M0'))\n", - "display(nanslice.jupyter.compare('T1.nii.gz', 'MUPA_T1.nii.gz', axis='z', clim=[0.5, 1.5], title='T1'))\n", - "display(nanslice.jupyter.compare('T2.nii.gz', 'MUPA_T2.nii.gz', axis='z', clim=[0.0, 0.15], title='T2'))" - ] - }, - { - "cell_type": "markdown", - "execution_count": null, - "metadata": {}, - "outputs": [], - "source": [ - "## 4 Point with Groups" - ] - }, - { - "cell_type": "code", - "execution_count": null, - "metadata": {}, - "outputs": [], - "source": [ - "seq = {\n", - " 'MUPA': {\n", - " 'TR': 2.34e-3,\n", - " 'Tramp': 10e-3,\n", - " 'spokes_per_seg': 512,\n", - " 'groups_per_seg': [1, 1, 4, 1],\n", - " 'FA': [2, 2, 2, 2],\n", - " 'Trf': [24, 24, 24, 24],\n", - " 'prep': ['inv', 'delay-800','t2-80', 'delay-800'],\n", - " 'prep_pulses': prep_pulses}\n", - "}\n", - "region='45,0,0,1,109,91'\n", - "sim_result = QUIT.interfaces.rufis.MUPASim(sequence=seq, in_file='sim.nii.gz', noise=0.0005, M0='M0.nii.gz', T1='T1.nii.gz', T2='T2.nii.gz', B1='B1.nii.gz', mask='M0.nii.gz', subregion=region).run()\n", - "display(nanslice.jupyter.timeseries(sim_result.outputs.out_file, clim=[-0.01,0.025]))\n", - "mupa = QUIT.interfaces.rufis.MUPA(sequence=seq, in_file='sim.nii.gz', mask_file='M0.nii.gz', subregion=region).run()\n", - "display(nanslice.jupyter.compare('M0.nii.gz', 'MUPA_M0.nii.gz', axis='z', clim=[-0.01, 1.25], title='M0'))\n", - "display(nanslice.jupyter.compare('T1.nii.gz', 'MUPA_T1.nii.gz', axis='z', clim=[0.5, 1.5], title='T1'))\n", - "display(nanslice.jupyter.compare('T2.nii.gz', 'MUPA_T2.nii.gz', axis='z', clim=[0.0, 0.15], title='T2'))" - ] - }, - { - "cell_type": "markdown", - "execution_count": null, - "metadata": {}, - "outputs": [], - "source": [ - "## T2-Prep Only" - ] - }, - { - "cell_type": "code", - "execution_count": null, - "metadata": {}, - "outputs": [], - "source": [ - "seq = {\n", - " 'MUPA': {\n", - " 'TR': 2.34e-3,\n", - " 'Tramp': 10e-3,\n", - " 'spokes_per_seg': 512,\n", - " 'groups_per_seg': [4, 4, 1, 1],\n", - " 'FA': [2, 2, 2, 2],\n", - " 'Trf': [24, 24, 24, 24],\n", - " 'prep': ['t2-40', 't2-80', 'null', 'delay-800'],\n", - " 'prep_pulses': prep_pulses}\n", - "}\n", - "region='45,0,0,1,109,91'\n", - "sim_result = QUIT.interfaces.rufis.MUPASim(sequence=seq, in_file='sim.nii.gz', noise=0.0005, M0='M0.nii.gz', T1='T1.nii.gz', T2='T2.nii.gz', B1='B1.nii.gz', mask='M0.nii.gz', subregion=region).run()\n", - "display(nanslice.jupyter.timeseries(sim_result.outputs.out_file, clim=[0,0.025]))\n", - "mupa = QUIT.interfaces.rufis.MUPA(sequence=seq, in_file='sim.nii.gz', mask_file='M0.nii.gz', subregion=region).run()\n", - "display(nanslice.jupyter.compare('M0.nii.gz', 'MUPA_M0.nii.gz', axis='z', clim=[-0.01, 1.25], title='M0'))\n", - "display(nanslice.jupyter.compare('T1.nii.gz', 'MUPA_T1.nii.gz', axis='z', clim=[0.5, 1.5], title='T1'))\n", - "display(nanslice.jupyter.compare('T2.nii.gz', 'MUPA_T2.nii.gz', axis='z', clim=[0.0, 0.15], title='T2'))" - ] - }, - { - "cell_type": "code", - "execution_count": null, - "metadata": {}, - "outputs": [], - "source": [ - "FAe = np.degrees(np.arccos(np.exp(-0.00234/1)))\n", - "print(FAe)" - ] - }, - { - "cell_type": "code", - "execution_count": null, - "metadata": {}, - "outputs": [], - "source": [ - "seq = {\n", - " 'MUPA': {\n", - " 'TR': 2.34e-3,\n", - " 'Tramp': 10e-3,\n", - " 'spokes_per_seg': 512,\n", - " 'groups_per_seg': [1, 1, 4],\n", - " 'FA': [3, 3, 3],\n", - " 'Trf': [24, 24, 24],\n", - " 'prep': ['inv', 'null','t2-60'],\n", - " 'prep_pulses': prep_pulses}\n", - "}\n", - "region='45,0,0,1,109,91'\n", - "sim_result = QUIT.interfaces.rufis.MUPASim(sequence=seq, in_file='sim.nii.gz', noise=0.0005, M0='M0.nii.gz', T1='T1.nii.gz', T2='T2.nii.gz', B1='B1.nii.gz', mask='M0.nii.gz', subregion=region).run()\n", - "display(nanslice.jupyter.timeseries(sim_result.outputs.out_file, clim=[-0.01,0.025]))\n", - "mupa = QUIT.interfaces.rufis.MUPA(sequence=seq, in_file='sim.nii.gz', mask_file='M0.nii.gz', subregion=region).run()\n", - "display(nanslice.jupyter.compare('M0.nii.gz', 'MUPA_M0.nii.gz', axis='z', clim=[-0.01, 1.25], title='M0'))\n", - "display(nanslice.jupyter.compare('T1.nii.gz', 'MUPA_T1.nii.gz', axis='z', clim=[0.5, 1.5], title='T1'))\n", - "display(nanslice.jupyter.compare('T2.nii.gz', 'MUPA_T2.nii.gz', axis='z', clim=[0.0, 0.15], title='T2'))" - ] - }, - { - "cell_type": "code", - "execution_count": null, - "metadata": {}, - "outputs": [], - "source": [] - } - ], - "metadata": { - "kernel_info": { - "name": "python3" - }, - "kernelspec": { - "display_name": "Python 3", - "language": "python", - "name": "python3" - }, - "language_info": { - "codemirror_mode": { - "name": "ipython", - "version": 3 - }, - "file_extension": ".py", - "mimetype": "text/x-python", - "name": "python", - "nbconvert_exporter": "python", - "pygments_lexer": "ipython3", - "version": "3.7.4-final" - }, - "nteract": { - "version": "0.8.4" - } - }, - "nbformat": 4, - "nbformat_minor": 4 -} \ No newline at end of file diff --git a/Python/steady_state.ipynb b/Python/steady_state.ipynb deleted file mode 100644 index a1835a4b..00000000 --- a/Python/steady_state.ipynb +++ /dev/null @@ -1,265 +0,0 @@ -{ - "cells": [ - { - "cell_type": "markdown", - "execution_count": null, - "metadata": {}, - "outputs": [], - "source": [ - "## Imports and Data Fetching" - ] - }, - { - "cell_type": "code", - "execution_count": null, - "metadata": {}, - "outputs": [], - "source": [ - "%matplotlib inline\n", - "%load_ext autoreload\n", - "%autoreload 1\n", - "%aimport nanslice\n", - "%aimport nanslice.jupyter\n", - "%aimport QUIT.base\n", - "%aimport QUIT.interfaces.rufis\n", - "import nibabel as nib\n", - "import numpy as np\n", - "import requests\n", - "import gzip\n", - "from os.path import isfile" - ] - }, - { - "cell_type": "markdown", - "execution_count": null, - "metadata": {}, - "outputs": [], - "source": [ - "## Setup" - ] - }, - { - "cell_type": "code", - "execution_count": null, - "metadata": {}, - "outputs": [], - "source": [ - "if not isfile('classes.mnc'):\n", - " print('Downloading classes')\n", - " params = {'download_for_real':'[Start Download!]',\n", - " 'do_download_alias':'phantom_1.0mm_normal_crisp',\n", - " 'format_value':'minc',\n", - " 'who_name': 'Tobias Wood',\n", - " 'who_institution': 'KCL',\n", - " 'who_email': 'tobias.wood@kcl.ac.uk'}\n", - " response = requests.get(url='http://brainweb.bic.mni.mcgill.ca/cgi/brainweb1', params=params)\n", - " minc_file = open('classes.mnc', 'wb')\n", - " minc_file.write(response.content)\n", - "if not isfile('rf20_B.mnc'):\n", - " print('Downloading B0')\n", - " params = {'download_for_real':'[Start Download!]',\n", - " 'do_download_alias':'rf20_B.mnc.gz',\n", - " 'format_value':'minc',\n", - " 'zip_value':'none',\n", - " 'who_name': 'Tobias Wood',\n", - " 'who_institution': 'KCL',\n", - " 'who_email': 'tobias.wood@kcl.ac.uk'}\n", - " response = requests.get(url='https://brainweb.bic.mni.mcgill.ca/cgi/brainweb1', params=params)\n", - " b0_file = open('rf20_B.mnc', 'wb')\n", - " b0_file.write(response.content)\n", - "if not isfile('rf20_C.mnc'):\n", - " print('Downloading B1')\n", - " params = {'download_for_real':'[Start Download!]',\n", - " 'do_download_alias':'rf20_C.mnc.gz',\n", - " 'format_value':'minc',\n", - " 'zip_value':'none',\n", - " 'who_name': 'Tobias Wood',\n", - " 'who_institution': 'KCL',\n", - " 'who_email': 'tobias.wood@kcl.ac.uk'}\n", - " response = requests.get(url='https://brainweb.bic.mni.mcgill.ca/cgi/brainweb1', params=params)\n", - " b1_file = open('rf20_C.mnc', 'wb')\n", - " b1_file.write(response.content)\n", - "classes = nanslice.Layer('classes.mnc')\n", - "b0_minc = nanslice.Layer('rf20_B.mnc')\n", - "b1_minc = nanslice.Layer('rf20_C.mnc')" - ] - }, - { - "cell_type": "code", - "execution_count": null, - "metadata": {}, - "outputs": [], - "source": [ - "# 0=Background\n", - "# 1=CSF\n", - "# 2=Grey Matter\n", - "# 3=White Matter\n", - "# 4=Fat\n", - "# 5=Muscle/Skin\n", - "# 6=Skin\n", - "# 7=Skull\n", - "# 8=Glial Matter\n", - "# 9=Connective\n", - "def classes_to_params(M0Vals, T1Vals, T2vals, subsamp=1):\n", - " class_data = classes.image.get_data().astype('int32')\n", - " M0data = np.choose(class_data[::subsamp,::subsamp,::subsamp], M0vals).astype('float32')\n", - " T1data = np.choose(class_data[::subsamp,::subsamp,::subsamp], T1vals).astype('float32')\n", - " T2data = np.choose(class_data[::subsamp,::subsamp,::subsamp], T2vals).astype('float32')\n", - " f0data = (b0_minc.image.get_data().astype('float32')[::subsamp,::subsamp,::subsamp] - 1) * 500\n", - " B1data = b1_minc.image.get_data().astype('float32')[::subsamp,::subsamp,::subsamp]\n", - " # PDdata = np.array(list(map(PDFunc, classes.image.get_data())))\n", - " M0image = nib.nifti1.Nifti1Image(M0data, affine=classes.image.affine)\n", - " T1image = nib.nifti1.Nifti1Image(T1data, affine=classes.image.affine)\n", - " T2image = nib.nifti1.Nifti1Image(T2data, affine=classes.image.affine)\n", - " f0image = nib.nifti1.Nifti1Image(f0data, affine=classes.image.affine)\n", - " B1image = nib.nifti1.Nifti1Image(B1data, affine=classes.image.affine)\n", - " nib.save(M0image, 'M0.nii.gz')\n", - " nib.save(T1image, 'T1.nii.gz')\n", - " nib.save(T2image, 'T2.nii.gz')\n", - " nib.save(f0image, 'f0.nii.gz')\n", - " nib.save(B1image, 'B1.nii.gz')" - ] - }, - { - "cell_type": "markdown", - "execution_count": null, - "metadata": {}, - "outputs": [], - "source": [ - "## Simulate and Display Image" - ] - }, - { - "cell_type": "code", - "execution_count": null, - "metadata": {}, - "outputs": [], - "source": [ - "M0vals = np.array([0, 1, 0.8, 0.7, 0, 0, 0, 0, 0, 0])\n", - "T1vals = np.array([0, 3.0, 1.3, 0.9, 0, 0, 0, 0, 0, 0])\n", - "T2vals = np.array([0, 0.5, 0.08, 0.06, 0, 0, 0, 0, 0, 0])\n", - "classes_to_params(M0vals, T1vals, T2vals, 1)\n", - "region='90,0,0,1,217,181'\n", - "region2='45,0,0,1,109,91'" - ] - }, - { - "cell_type": "code", - "execution_count": null, - "metadata": {}, - "outputs": [], - "source": [ - "seq = {\n", - " \"TR\": 0.002,\n", - " \"Tramp\": 0.01,\n", - " \"Tspoil\": 0.05,\n", - " \"FA\": [2, 8, 8],\n", - " \"Trf\": [28e-6, 28e-6, 28e-6],\n", - " \"spokes_per_seg\": 48,\n", - " \"groups_per_seg\": [1, 1, 1],\n", - " \"prep_FA\": [0, 0, 360],\n", - " \"prep_df\": [0, 0, 0],\n", - " \"prep_Trf\": 0.001,\n", - " \"prep_p1\": 1.,\n", - " \"prep_p2\": 1.\n", - " }\n", - "inputs={'M0':'M0.nii.gz', 'T1':'T1.nii.gz', 'B1':'B1.nii.gz'}\n", - "sim_result = QUIT.interfaces.rufis.SteadyStateSim(sequence=seq, in_file='sim.nii.gz', noise=0.00025, mask='M0.nii.gz', subregion=region, **inputs).run()\n", - "display(nanslice.jupyter.timeseries(sim_result.outputs.out_file, clim=[0,0.025]))\n", - "mupa = QUIT.interfaces.rufis.SteadyState(sequence=seq, in_file='sim.nii.gz', mask_file='M0.nii.gz', subregion=region).run()\n", - "for name, fname in inputs.items():\n", - " display(nanslice.jupyter.compare(fname, 'SS_' + fname, axis='z', title=name))" - ] - }, - { - "cell_type": "code", - "execution_count": null, - "metadata": {}, - "outputs": [], - "source": [ - "seq = {\n", - " \"TR\": 0.002,\n", - " \"Tramp\": 0.01,\n", - " \"Tspoil\": 0.05,\n", - " \"FA\": [2, 8, 2, 8, 2, 8],\n", - " \"Trf\": [28e-6, 28e-6, 28e-6, 28e-6, 28e-6, 28e-6],\n", - " \"spokes_per_seg\": 48,\n", - " \"groups_per_seg\": [1, 1, 1, 1, 1, 1],\n", - " \"prep_FA\": [0, 180, 180, 360, 540, 540],\n", - " \"prep_df\": [0, -0.25, 0.25, 0, -0.1, 0.1],\n", - " \"prep_Trf\": 0.01,\n", - " \"prep_p1\": 1.,\n", - " \"prep_p2\": 1.\n", - " }\n", - "inputs={'M0':'M0.nii.gz', 'T1':'T1.nii.gz', 'T2':'T2.nii.gz', 'f0': 'f0.nii.gz', 'B1':'B1.nii.gz'}\n", - "sim_result = QUIT.interfaces.rufis.SteadyStateSim(sequence=seq, fitT2=True, in_file='sim.nii.gz', noise=0.00025, mask='M0.nii.gz', subregion=region, **inputs).run()\n", - "display(nanslice.jupyter.timeseries(sim_result.outputs.out_file, clim=[0,0.025]))\n", - "mupa = QUIT.interfaces.rufis.SteadyState(sequence=seq, fitT2=True, in_file='sim.nii.gz', mask_file='M0.nii.gz', subregion=region).run()\n", - "for name, fname in inputs.items():\n", - " display(nanslice.jupyter.compare(fname, 'SS_' + fname, axis='z', title=name))" - ] - }, - { - "cell_type": "code", - "execution_count": null, - "metadata": {}, - "outputs": [], - "source": [ - "seq = {\n", - " \"TR\": 0.002,\n", - " \"Tramp\": 0.01,\n", - " \"Tspoil\": 0.05,\n", - " \"FA\": [2, 8, 2, 2, 2, 2, 2],\n", - " \"Trf\": [28e-6, 28e-6, 28e-6, 28e-6, 28e-6, 28e-6, 28e-6],\n", - " \"spokes_per_seg\": 48,\n", - " \"groups_per_seg\": [1, 1, 1, 1, 1, 1, 1],\n", - " \"prep_FA\": [0, 0, 180, 180, 360, 540, 540],\n", - " \"prep_df\": [0, 0, -0.25, 0.25, 0, -0.1, 0.1],\n", - " \"prep_Trf\": 0.01,\n", - " \"prep_p1\": 1.,\n", - " \"prep_p2\": 1.\n", - " }\n", - "inputs={'M0':'M0.nii.gz', 'T1':'T1.nii.gz', 'T2':'T2.nii.gz', 'f0': 'f0.nii.gz', 'B1':'B1.nii.gz'}\n", - "sim_result = QUIT.interfaces.rufis.SteadyStateSim(sequence=seq, fitT2=True, in_file='sim.nii.gz', noise=0.00025, mask='M0.nii.gz', subregion=region, **inputs).run()\n", - "display(nanslice.jupyter.timeseries(sim_result.outputs.out_file, clim=[0,0.025]))\n", - "mupa = QUIT.interfaces.rufis.SteadyState(sequence=seq, fitT2=True, in_file='sim.nii.gz', mask_file='M0.nii.gz', subregion=region).run()\n", - "for name, fname in inputs.items():\n", - " display(nanslice.jupyter.compare(fname, 'SS_' + fname, axis='z', title=name))" - ] - }, - { - "cell_type": "code", - "execution_count": null, - "metadata": {}, - "outputs": [], - "source": [] - } - ], - "metadata": { - "kernel_info": { - "name": "python3" - }, - "kernelspec": { - "display_name": "Python 3", - "language": "python", - "name": "python3" - }, - "language_info": { - "codemirror_mode": { - "name": "ipython", - "version": 3 - }, - "file_extension": ".py", - "mimetype": "text/x-python", - "name": "python", - "nbconvert_exporter": "python", - "pygments_lexer": "ipython3", - "version": "3.7.4-final" - }, - "nteract": { - "version": "0.8.4" - } - }, - "nbformat": 4, - "nbformat_minor": 4 -} \ No newline at end of file diff --git a/Source/PARMESAN/transient_main.cpp b/Source/PARMESAN/transient_main.cpp index 6ab6c10e..4ba63c6a 100644 --- a/Source/PARMESAN/transient_main.cpp +++ b/Source/PARMESAN/transient_main.cpp @@ -28,7 +28,7 @@ /* * Main */ -int rufis_transient_main(args::Subparser &parser) { +int transient_main(args::Subparser &parser) { args::Positional input_path(parser, "INPUT", "Input MUPA file"); QI_COMMON_ARGS; args::Flag mt(parser, "MT", "Use MT model", {"mt"}); @@ -58,7 +58,6 @@ int rufis_transient_main(args::Subparser &parser) { simulate.Get(), subregion.Get()); } else { - QI::Log(verbose, "NS {}", decltype(model)::NS); using FitType = QI::ScaledNumericDiffFit; FitType fit{model}; auto fit_filter = diff --git a/Source/Relaxometry/qimcdespot.cpp b/Source/Relaxometry/qimcdespot.cpp index 27e3d10f..413bb8e9 100644 --- a/Source/Relaxometry/qimcdespot.cpp +++ b/Source/Relaxometry/qimcdespot.cpp @@ -37,9 +37,7 @@ template struct MCDSRCFunctor { const Eigen::ArrayXd &d, const Eigen::ArrayXd &w) : data(d), - weights(w), fixed(f), model(m) { - assert(data.rows() == model.sequence.size()); - } + weights(w), fixed(f), model(m) {} int inputs() const { return Model::NV; } int values() const { return model.spgr.size() + model.ssfp.size(); } diff --git a/Source/qi_main.cpp b/Source/qi_main.cpp index 3ca8f08d..7f3b2a34 100644 --- a/Source/qi_main.cpp +++ b/Source/qi_main.cpp @@ -7,7 +7,8 @@ int main(int argc, char **argv) { args::ArgumentParser parser("QUIT http://github.com/spinicist/QUIT"); args::Group commands(parser, "COMMANDS"); args::GlobalOptions globals(parser, global_group); - args::Flag version(parser, "VERSION", "Print the version of QUIT", {"version"}); + args::Group top(parser, "TOP"); + args::Flag version(top, "VERSION", "Print the version of QUIT", {"version"}); #define ADD(CMD, HELP) args::Command CMD(commands, #CMD, HELP, &CMD##_main); ADD(newimage, "Create a new image"); @@ -84,11 +85,13 @@ int main(int argc, char **argv) { std::cerr << parser << '\n'; exit(EXIT_SUCCESS); } catch (args::Error e) { + if (version) { + std::cout << QI::GetVersion() << '\n'; + exit(EXIT_SUCCESS); + } std::cerr << parser << '\n' << e.what() << '\n'; exit(EXIT_FAILURE); } - if (version) { - std::cout << QI::GetVersion() << '\n'; - } + exit(EXIT_SUCCESS); } \ No newline at end of file