{ "cells": [ { "attachments": {}, "cell_type": "markdown", "metadata": {}, "source": [ "## GRE walkthrough\n", "\n", "This jupyter notebook is a walkthrough to construct a Gradient Recalled Echo (GRE) pulse sequence.\n", "\n", "\n", "\n", "1. First, import the necessary packages." ] }, { "cell_type": "code", "execution_count": 2, "metadata": {}, "outputs": [], "source": [ "import math\n", "\n", "import numpy as np\n", "\n", "from pypulseq.Sequence.sequence import Sequence\n", "from pypulseq.calc_duration import calc_duration\n", "from pypulseq.make_adc import make_adc\n", "from pypulseq.make_delay import make_delay\n", "from pypulseq.make_sinc_pulse import make_sinc_pulse\n", "from pypulseq.make_trapezoid import make_trapezoid\n", "from pypulseq.opts import Opts" ] }, { "attachments": {}, "cell_type": "markdown", "metadata": {}, "source": [ "2. Construct a `Sequence` object. A `Sequence` object is analogous to a pulse sequence. The following steps will demonstrate creating pulse sequence events and adding them to this `Sequence` object. The MR scanner's hardware limits are defined in an `Opts` object." ] }, { "cell_type": "code", "execution_count": 3, "metadata": {}, "outputs": [], "source": [ "seq = Sequence()\n", "fov = 256e-3 # field of view\n", "Nx = 256 # number of frequency encodes\n", "Ny = 256 # number of phase encodes\n", "alpha = 10 # RF flip\n", "slice_thickness = 3e-3\n", "TE = 7.38e-3 # echo time\n", "TR = 100e-3 # repetition time\n", "\n", "rf_spoiling_inc = 117\n", "\n", "sys = Opts(max_grad=28, grad_unit='mT/m', max_slew=150, slew_unit='T/m/s', rf_ringdown_time=20e-6, rf_dead_time=100e-6,\n", " adc_dead_time=10e-6)" ] }, { "attachments": {}, "cell_type": "markdown", "metadata": {}, "source": [ "3. Now, construct the radio-frequency (RF), gradient, delay and ADC readout events." ] }, { "cell_type": "code", "execution_count": 4, "metadata": {}, "outputs": [ { "ename": "TypeError", "evalue": "cannot unpack non-iterable types.SimpleNamespace object", "output_type": "error", "traceback": [ "\u001b[0;31m---------------------------------------------------------------------------\u001b[0m", "\u001b[0;31mTypeError\u001b[0m Traceback (most recent call last)", "\u001b[1;32m/Users/bilal/workspace/pypulseq_bilal/doc/walkthrough/gre_walkthrough.ipynb Cell 6\u001b[0m in \u001b[0;36m\u001b[0;34m()\u001b[0m\n\u001b[0;32m----> 1\u001b[0m rf, gz, gzr \u001b[39m=\u001b[39m make_sinc_pulse(flip_angle\u001b[39m=\u001b[39malpha \u001b[39m*\u001b[39m math\u001b[39m.\u001b[39mpi \u001b[39m/\u001b[39m \u001b[39m180\u001b[39m, duration\u001b[39m=\u001b[39m\u001b[39m4e-3\u001b[39m, slice_thickness\u001b[39m=\u001b[39mslice_thickness,\n\u001b[1;32m 2\u001b[0m apodization\u001b[39m=\u001b[39m\u001b[39m0.5\u001b[39m, time_bw_product\u001b[39m=\u001b[39m\u001b[39m4\u001b[39m, system\u001b[39m=\u001b[39msys)\n\u001b[1;32m 4\u001b[0m delta_k \u001b[39m=\u001b[39m \u001b[39m1\u001b[39m \u001b[39m/\u001b[39m fov\n\u001b[1;32m 5\u001b[0m gx \u001b[39m=\u001b[39m make_trapezoid(channel\u001b[39m=\u001b[39m\u001b[39m'\u001b[39m\u001b[39mx\u001b[39m\u001b[39m'\u001b[39m, flat_area\u001b[39m=\u001b[39mNx \u001b[39m*\u001b[39m delta_k, flat_time\u001b[39m=\u001b[39m\u001b[39m6.4e-3\u001b[39m, system\u001b[39m=\u001b[39msys)\n", "\u001b[0;31mTypeError\u001b[0m: cannot unpack non-iterable types.SimpleNamespace object" ] } ], "source": [ "rf, gz, gzr = make_sinc_pulse(flip_angle=alpha * math.pi / 180, duration=4e-3, slice_thickness=slice_thickness,\n", " apodization=0.5, time_bw_product=4, system=sys, return_gz=True)\n", "\n", "delta_k = 1 / fov\n", "gx = make_trapezoid(channel='x', flat_area=Nx * delta_k, flat_time=6.4e-3, system=sys)\n", "adc = make_adc(num_samples=Nx, duration=gx.flat_time, delay=gx.rise_time, system=sys)\n", "gx_pre = make_trapezoid(channel='x', area=-gx.area / 2, duration=2e-3, system=sys)\n", "gz_reph = make_trapezoid(channel='z', area=-gz.area / 2, duration=2e-3, system=sys)\n", "phase_areas = (np.arange(Ny) - Ny / 2) * delta_k\n", "\n", "gx_spoil = make_trapezoid(channel='x', area=2 * Nx * delta_k, system=sys)\n", "gz_spoil = make_trapezoid(channel='z', area=4 / slice_thickness, system=sys)\n", "\n", "delay_TE = math.ceil((TE - calc_duration(gx_pre) - gz.fall_time - gz.flat_time / 2 - calc_duration(\n", " gx) / 2) / seq.grad_raster_time) * seq.grad_raster_time\n", "delay_TR = math.ceil((TR - calc_duration(gx_pre) - calc_duration(gz) - calc_duration(\n", " gx) - delay_TE) / seq.grad_raster_time) * seq.grad_raster_time\n", "\n", "assert np.all(delay_TR >= calc_duration(gx_spoil, gz_spoil))" ] }, { "attachments": {}, "cell_type": "markdown", "metadata": {}, "source": [ "4. Add the constructed pulse sequence events to the `Sequence` object created earlier.\n", "\n", "The functionality of the `add_block()` method can be tested by adding only one event to the `Sequence` object and obtaining the duration of the added events. For example, add only the `rf` object to `seq`. Then, call `seq.duration()` to obtain the duration of all the events added to `seq` (in this case, only the `rf` event so far)." ] }, { "cell_type": "code", "execution_count": null, "metadata": {}, "outputs": [], "source": [ "rf_phase = 0\n", "rf_inc = 0\n", "\n", "for i in range(Ny): # We have Ny phase encodes\n", " rf.phase_offset = rf_phase / 180 * np.pi\n", " adc.phase_offset = rf_phase / 180 * np.pi\n", " rf_inc = divmod(rf_inc + rf_spoiling_inc, 360.0)[1]\n", " rf_phase = divmod(rf_phase + rf_inc, 360.0)[1]\n", "\n", " seq.add_block(rf, gz)\n", " gy_pre = make_trapezoid(channel='y', area=phase_areas[i], duration=2e-3, system=sys)\n", " seq.add_block(gx_pre, gy_pre, gz_reph)\n", " seq.add_block(make_delay(delay_TE))\n", " seq.add_block(gx, adc)\n", " gy_pre.amplitude = -gy_pre.amplitude\n", " seq.add_block(make_delay(delay_TR), gx_spoil, gy_pre, gz_spoil)" ] }, { "attachments": {}, "cell_type": "markdown", "metadata": {}, "source": [ "5. Visualise the constructed pulse sequence by calling the `plot()` command. The `plot()` command visualises the ADC and RF events in one window, and the gradient events in another window. For GRE, the `plot()` command should display two plots that look like this:\n", "\n", "![GRE Plot 1](./gre_1.png)\n", "![GRE Plot 2](./gre_2.png)" ] }, { "cell_type": "code", "execution_count": null, "metadata": {}, "outputs": [], "source": [ "seq.plot()" ] }, { "attachments": {}, "cell_type": "markdown", "metadata": {}, "source": [ "6. Export the pulse sequence as a `.seq` file. This file can be executed on MR scanners (Siemens/[GE](https://toppemri.github.io)/\n", "[Bruker](https://github.com/pulseq/bruker_interpreter)) to acquire data." ] }, { "cell_type": "code", "execution_count": null, "metadata": {}, "outputs": [], "source": [ "seq.write('gre_pypulseq.seq')" ] } ], "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.10.4" } }, "nbformat": 4, "nbformat_minor": 2 }