{
"cells": [
{
"cell_type": "markdown",
"metadata": {},
"source": [
"# k-Space lab in Python\n",
"\n",
"> **_NOTE:_** If you use ChatGPT or Github Copilot, please make sure that you can answer any questions about the code.\n",
"\n",
"This is an interactive Jupyter Notebook. The notebook consists of markdown cells (like this) and code cells. To edit a code cell, you just click in it and start typing. To execute code you can either click the Run button, or use the keyboard shortcut Shift + Enter.\n",
"\n",
"If you're running the lab using Google Colab, there is a button above to save a copy of the notebook to your Google Drive.\n",
"\n",
"This lab is making use of the NumPy (tutorial) and PyPlot (tutorial) packages.\n",
"\n",
"A few tips for numpy:\n",
"* Numpy uses zero-based indexing, e.g. ``my_nparray[0]`` is the first element.\n",
"* Arrays can be sliced using the ``[start:stop:step]`` notation, e.g. ``my_nparray[::2]`` for every second element or ``my_nparray[0:20]`` for the 20 first elements.\n",
"* Before modifying a variable, it is a good idea to make a copy, e.g. ``my_kspace = kspace_res512_FOV24.copy()``, otherwise you might end up with unwanted results.\n",
"* If you want to initialize a zero matrix with ``np.zeros``, you need to specify that it's a complex datatype using ``dtype=np.complex64``\n",
"* The command ``np.log`` may give undesired results for arrays containing zeros. You can use ``np.log(x+1)`` for those instances.\n",
"* ``ij`` or ``0+ij`` can be used to represent complex numbers. ``np.pi`` can be used for pi.\n",
"* In numpy convention, the first dimension (axis=0) corresponds to rows and the second dimension (axis=1) corresponds to columns.\n",
"\n",
"A few tips for matplotlib:\n",
"* Use ``plt.imshow`` to display images, and use the option ``cmap='gray'`` to display them in black and white.\n",
"* When the question calls for multiple plots, you can either use ``plt.figure()`` to create multiple plots or use ``fig, axes = plt.subplots(2,2)`` to create a 2 by 2 grid.\n",
"* If your figures appear too small, you can use ``figsize=(20,20)`` to make the plotting area larger, e.g. ``plt.figure(figsize=(20,20))`` or ``fig, axes = plt.subplots(2,2, figsize=(20,20))``\n",
"\n",
"Finally, if you need to reset your workspace, run the first cell again!"
]
},
{
"cell_type": "code",
"execution_count": null,
"metadata": {},
"outputs": [],
"source": [
"!wget -nc https://raw.githubusercontent.com/fyrdahl/kspace-lab/master/data.h5\n",
"\n",
"import numpy as np\n",
"import matplotlib.pyplot as plt\n",
"import h5py\n",
"\n",
"from scipy.signal import convolve2d\n",
"\n",
"with h5py.File('data.h5', 'r') as F:\n",
" kspace_res512_FOV24 = np.array(F['kspace_res512_FOV24'])\n",
" kspace_res512_FOV48 = np.array(F['kspace_res512_FOV48'])\n",
"\n",
"fftc = lambda f: np.fft.fftshift(np.fft.fft(np.fft.ifftshift(f)))\n",
"ifftc = lambda F: np.fft.fftshift(np.fft.ifft(np.fft.ifftshift(F)))\n",
"fft2c = lambda f: np.fft.fftshift(np.fft.fft2(np.fft.ifftshift(f)))\n",
"ifft2c = lambda F: np.fft.fftshift(np.fft.ifft2(np.fft.ifftshift(F)))"
]
},
{
"cell_type": "markdown",
"metadata": {},
"source": [
"## Task 1\n",
"\n",
"Use ``kspace_res512_FOV24``\n",
"\n",
"1. Show the magnitude, phase, real, and imaginary part of the k-space. To better visualize the k-space, you could log the matrix using the ``np.log`` command."
]
},
{
"cell_type": "code",
"execution_count": null,
"metadata": {},
"outputs": [],
"source": [
"fig, axes = plt.subplots(2,2)\n",
"axes[0][0].imshow(np.abs(np.log(kspace_res512_FOV24)), cmap='gray')\n",
"axes[0][0].set_title('Magnitude')\n",
"axes[0][0].set_axis_off()\n",
"\n",
"# Plot phase here\n",
"axes[0][1].imshow(..., cmap='gray')\n",
"axes[0][1].set_title('Phase')\n",
"axes[0][1].set_axis_off()\n",
"\n",
"# Plot real part here\n",
"axes[1][0].imshow(..., cmap='gray')\n",
"axes[1][0].set_title('Real')\n",
"axes[1][0].set_axis_off()\n",
"\n",
"# Plot imaginary part here\n",
"axes[1][1].imshow(..., cmap='gray')\n",
"axes[1][1].set_title('Imaginary')\n",
"axes[1][1].set_axis_off()\n",
"\n",
"plt.show()"
]
},
{
"cell_type": "markdown",
"metadata": {},
"source": [
"2. Do an inverse FT of the k-space data along one dimension (ifftc).\n",
" * Display the magnitude and phase of the result.\n",
" * Describe what happened, which domain is the data in – spatial or frequency?\n",
" * What are the unit on the x and y axis?"
]
},
{
"cell_type": "code",
"execution_count": null,
"metadata": {},
"outputs": [],
"source": []
},
{
"cell_type": "markdown",
"metadata": {},
"source": [
"3. Do an inverse FT of the k-space data in both the frequency and phase encoding direction.\n",
" * Display the magnitude and phase of the result.\n",
" * Which domain is the data in – spatial or frequency?\n",
" * Which dimension is most likely to be the phase encoding direction?"
]
},
{
"cell_type": "code",
"execution_count": null,
"metadata": {},
"outputs": [],
"source": []
},
{
"cell_type": "markdown",
"metadata": {},
"source": [
"## Task 2\n",
"\n",
"Use ``kspace_res512_FOV48``\n",
"\n",
"1. Decrease the FOV in the spatial domain to 24x24 cm, and 16x16 cm by removing parts of k-space and show the results"
]
},
{
"cell_type": "code",
"execution_count": null,
"metadata": {},
"outputs": [],
"source": []
},
{
"cell_type": "markdown",
"metadata": {},
"source": [
"## Task 3\n",
"\n",
"Use ``kspace_res512_FOV24``\n",
"\n",
"1. Decrease the spatial resolution by replacing k-space data with zeros before performing the ifft2c. The new images\n",
"should have a 512x256, 256x256, 64x64, 512x64 and 64x512 resolution. Display the magnitude of k-space and the magnitude image in the spatial domain. Remember that you should use ``np.log(x+1)`` when taking the \"log\" a matrix containing zeros."
]
},
{
"cell_type": "code",
"execution_count": null,
"metadata": {},
"outputs": [],
"source": []
},
{
"cell_type": "markdown",
"metadata": {},
"source": [
"2. What is the difference between the 512x64 and the 64x512 image? Describe the artifact that is appearing."
]
},
{
"cell_type": "markdown",
"metadata": {},
"source": [
"## Task 4\n",
"\n",
"Use ``kspace_res512_FOV24``\n",
"\n",
"1. Set the k-space element ``[kx,ky]`` to 200000 for these elements: ``[100,100]``, ``[200,200]``, ``[250,250]``, ``[255,255]``, ``[257,257]``. Display the image for _each_ k-space adjustment. Describe the artifact.\n",
"\n",
"Remember that numpy arrays are mutable. If you want to create an image for _each_ k-space location, you need to use ``.copy()`` to create copies of the original array."
]
},
{
"cell_type": "code",
"execution_count": null,
"metadata": {},
"outputs": [],
"source": []
},
{
"cell_type": "markdown",
"metadata": {},
"source": [
"## Task 5\n",
"\n",
"Use ``kspace_res512_FOV24``\n",
"\n",
"Change the k-space phase without altering the magnitude. Remember that complex\n",
"numbers can be written as ``z = mag * np.exp(1j*phase)``\n",
"\n",
"1. Increase the phase by 2 radians in the following k-space rows: 1–20, 201–220, 237–256 and 247–266. Display the magnitude image in the spatial domain for each k-space adjustment"
]
},
{
"cell_type": "code",
"execution_count": null,
"metadata": {},
"outputs": [],
"source": []
},
{
"cell_type": "markdown",
"metadata": {},
"source": [
"## Task 6\n",
"\n",
"Use ``kspace_res512_FOV24``\n",
"\n",
"1. Filter the spatial image with a boxcar filter. \n",
"2. Use two kernel sizes: 7x7 and 4x4.\n",
"3. Do an inverse FT and display the magnitude of k-space for each kernel size. Describe the results.\n",
"\n",
"To define a boxcar filter, you can use ``filter = np.ones((n,n))``, then use ``convolve2d(image,filter)`` to perform a 2D convolution. "
]
},
{
"cell_type": "code",
"execution_count": null,
"metadata": {},
"outputs": [],
"source": []
}
],
"metadata": {
"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.8.2-final"
}
},
"nbformat": 4,
"nbformat_minor": 2
}