.. DO NOT EDIT. .. THIS FILE WAS AUTOMATICALLY GENERATED BY SPHINX-GALLERY. .. TO MAKE CHANGES, EDIT THE SOURCE PYTHON FILE: .. "generated/gallery/01-course/06_radial_mprage.py" .. LINE NUMBERS ARE GIVEN BELOW. .. only:: html .. note:: :class: sphx-glr-download-link-note :ref:`Go to the end ` to download the full example code. .. rst-class:: sphx-glr-example-title .. _sphx_glr_generated_gallery_01-course_06_radial_mprage.py: ================ 6. Radial MPRAGE ================ MPRAGE is the T1-weighted 3D scan of almost every brain protocol. An inversion pulse flips the longitudinal magnetisation; while it recovers, at a speed set by each tissue's T1, a fast train of gradient echoes reads it out; then the magnetisation recovers before the next inversion. White and grey matter recover at different speeds, and the train catches them while they differ. In this lesson your radial gradient echo of lesson 5 becomes the readout train of a 3D MPRAGE, played as a *stack of stars*: the spokes are in-plane, and a phase encoding along z picks the partition. Radial sampling adds something a Cartesian MPRAGE does not have: every spoke passes through the centre of k-space, so every part of the train carries the contrast of its own moment, and the train can be split into several images of different contrast. **Learning objectives** - Add an adiabatic inversion pulse and its crusher, and set the inversion time. - Turn the 2D radial readout into a stack of stars, with a partition encoding merged into the slab rephaser. - Play one train of golden-angle spokes per inversion, and label every readout with its spoke and partition. - See how the signal evolves along the train, and why golden-angle spokes let you reconstruct it as several contrasts. Previous: :doc:`05_radial_sampling`. Next: :doc:`07_hardware_and_safety_checks`, where you check your MPRAGE against the scanner's limits. .. GENERATED FROM PYTHON SOURCE LINES 33-45 .. GENERATED FROM PYTHON SOURCE LINES 46-61 Inversion-recovery train ------------------------ Each inversion is one period of the scan, and it plays, in order: 1. the inversion pulse, then a crusher that dephases whatever it left in the transverse plane; 2. a wait, so the first spoke comes at the inversion time (TI) after the inversion; 3. the train: one gradient echo per spoke, all in the same partition; 4. a recovery delay, until the next inversion. Your protocol: a 220 mm field of view in-plane and a 128 matrix, 16 partitions of 5 mm, 64 spokes per inversion, TI 900 ms and an inversion every 2.5 s. .. GENERATED FROM PYTHON SOURCE LINES 61-89 .. code-block:: Python import numpy as np import pypulseqpp as pp from pypulseqpp import sequences system = pp.Opts( max_grad=32.0, grad_unit="mT/m", max_slew=130.0, slew_unit="T/m/s", rf_dead_time=100e-6, rf_ringdown_time=20e-6, adc_dead_time=10e-6, ) FOV = 220e-3 MATRIX = 128 PARTITIONS = 16 PARTITION_THICKNESS = 5e-3 SPOKES_PER_TRAIN = 64 FLIP_ANGLE_DEG = 8.0 TI = 900e-3 T_INVERSION = 2.5 SLAB = PARTITIONS * PARTITION_THICKNESS delta_k = 1 / FOV delta_kz = 1 / SLAB .. GENERATED FROM PYTHON SOURCE LINES 90-98 Inversion pulse --------------- An adiabatic pulse inverts the magnetisation uniformly even where the transmit field is weak or strong, which a 180 degree sinc would not. :func:`~pypulseqpp.make_adiabatic_pulse` makes a hyperbolic secant, 10 ms long. ``use="inversion"`` tells the scanner and the checks what it is. The crusher after it is the spoiler of lesson 3, along z. .. GENERATED FROM PYTHON SOURCE LINES 98-109 .. code-block:: Python rf_inv = pp.make_adiabatic_pulse( pulse_type="hypsec", duration=10e-3, bandwidth=40e3, adiabaticity=4, delay=system.rf_dead_time, use="inversion", system=system, ) gz_crush = pp.make_crusher(4.0, PARTITION_THICKNESS, channel="z", system=system)[0] .. GENERATED FROM PYTHON SOURCE LINES 110-115 Readout train ------------- The excitation now selects the whole slab, and the readout is the merged radial readout of lesson 5. .. GENERATED FROM PYTHON SOURCE LINES 115-156 .. code-block:: Python rf, gz, gz_reph = pp.make_sinc_pulse( flip_angle=np.deg2rad(FLIP_ANGLE_DEG), duration=1e-3, slice_thickness=SLAB, apodization=0.5, time_bw_product=4.0, delay=system.rf_dead_time, system=system, use="excitation", return_gz=True, ) dwell, readout_time = pp.calc_adc_timing( MATRIX, 10e-6, grad_raster_time=system.grad_raster_time, adc_raster_time=system.adc_raster_time, ) raster = system.grad_raster_time plateau = readout_time + raster * np.ceil(2 * system.adc_dead_time / raster) G = MATRIX * delta_k / readout_time _, times, amplitudes = pp.make_extended_trapezoid_area( channel="x", grad_start=0.0, grad_end=G, area=-G * plateau / 2, system=system ) plateau_start = times[-1] g_read = pp.make_extended_trapezoid( "x", times=np.append(times, plateau_start + plateau), amplitudes=np.append(amplitudes, G), system=system, ) adc = pp.make_adc( MATRIX, dwell=dwell, delay=plateau_start + (plateau - readout_time) / 2, system=system, ) g_spoil, _, _ = pp.make_extended_trapezoid_area( channel="x", grad_start=G, grad_end=0.0, area=2 * MATRIX * delta_k, system=system ) .. GENERATED FROM PYTHON SOURCE LINES 157-170 Stack of stars -------------- Each partition :math:`p` sits at :math:`k_z = (p - P/2)\,\Delta k_z`. The partition encoding plays on z while the prewinder plays on x and y, in the same block as the slab rephaser, and a block holds one gradient per channel. So you design one z gradient per partition that does both: the rephaser's area plus the partition's :math:`k_z`, as long as the prewinder. After the readout, a z gradient in the spoiler block brings :math:`k_z` back to zero. These are 16 different gradients, and you still design them before the loop, once each. .. GENERATED FROM PYTHON SOURCE LINES 170-180 .. code-block:: Python kz = (np.arange(PARTITIONS) - PARTITIONS // 2) * delta_kz gz_encode = [ pp.make_trapezoid("z", area=gz_reph.area + k, duration=plateau_start, system=system) for k in kz ] gz_rewind = [ pp.make_trapezoid("z", area=-k, duration=pp.calc_duration(g_spoil), system=system) for k in kz ] .. GENERATED FROM PYTHON SOURCE LINES 181-188 Inversion time -------------- TI runs from the centre of the inversion pulse to the centre of the first excitation. What lies between them is the rest of the inversion block, the crusher and the start of the excitation block; the wait block makes up the difference. .. GENERATED FROM PYTHON SOURCE LINES 188-210 .. code-block:: Python inversion_to_end = pp.calc_duration(rf_inv) - ( rf_inv.delay + pp.calc_rf_center(rf_inv)[0] ) start_to_excitation = rf.delay + pp.calc_rf_center(rf)[0] ti_wait = pp.round_to_raster( TI - inversion_to_end - pp.calc_duration(gz_crush) - start_to_excitation, system.block_duration_raster, ) spoke_tr = ( pp.calc_duration(rf, gz) + pp.calc_duration(g_read, adc, gz_encode[0]) + pp.calc_duration(g_spoil, gz_rewind[0]) ) train = SPOKES_PER_TRAIN * spoke_tr played = pp.calc_duration(rf_inv) + pp.calc_duration(gz_crush) + ti_wait + train recovery = pp.round_to_raster(T_INVERSION - played, system.block_duration_raster) print( f"spoke TR {1e3 * spoke_tr:.2f} ms, train {1e3 * train:.0f} ms, " f"recovery {1e3 * recovery:.0f} ms" ) .. rst-class:: sphx-glr-script-out .. code-block:: none spoke TR 4.06 ms, train 260 ms, recovery 1336 ms .. GENERATED FROM PYTHON SOURCE LINES 211-219 Scan loop --------- One inversion per partition. The spokes keep advancing by the golden angle across the whole scan, not only within a train, so every partition gets its own set of angles and together they fill the plane evenly. The labels say which spoke (``LIN``) and which partition (``PAR``) each readout is, and RF spoiling runs as in lesson 3. .. GENERATED FROM PYTHON SOURCE LINES 219-247 .. code-block:: Python GOLDEN_ANGLE = np.pi * (3 - np.sqrt(5)) / 2 seq = pp.Sequence(system) labels = sequences.Labels() n = 0 for partition in range(PARTITIONS): seq.add_block(rf_inv) seq.add_block(gz_crush) seq.add_block(pp.make_delay(ti_wait)) for spoke in range(SPOKES_PER_TRAIN): angle = (n * GOLDEN_ANGLE) % np.pi phase = (np.deg2rad(117.0) * spoke * (spoke + 1) / 2) % (2 * np.pi) rf.phase_offset = adc.phase_offset = phase seq.add_block(rf, gz, *labels(LIN=spoke, PAR=partition)) seq.add_block( *pp.rotate(g_read, angle=angle, axis="z"), adc, gz_encode[partition] ) seq.add_block(*pp.rotate(g_spoil, angle=angle, axis="z"), gz_rewind[partition]) n += 1 seq.add_block(pp.make_delay(recovery)) seq.set_definition("Name", "radial_mprage") seq.set_definition("FOV", [FOV, FOV, SLAB]) seq.set_definition("TI", TI) ok, errors = seq.check_timing() print(f"timing ok: {ok}, scan time {seq.duration()[0]:.1f} s") .. rst-class:: sphx-glr-script-out .. code-block:: none timing ok: True, scan time 40.0 s .. GENERATED FROM PYTHON SOURCE LINES 248-250 One inversion period, drawn over 2.5 s: the inversion, the wait, the train, the recovery. The train is a solid band at this scale. .. GENERATED FROM PYTHON SOURCE LINES 250-252 .. code-block:: Python seq.paper_plot(time_range=(0, T_INVERSION)) .. image-sg:: /generated/gallery/01-course/images/sphx_glr_06_radial_mprage_001.png :alt: 06 radial mprage :srcset: /generated/gallery/01-course/images/sphx_glr_06_radial_mprage_001.png :class: sphx-glr-single-img .. GENERATED FROM PYTHON SOURCE LINES 253-255 And the first spokes of the train, where you can see the partition encoding on z next to the prewinder: .. GENERATED FROM PYTHON SOURCE LINES 255-258 .. code-block:: Python start = pp.calc_duration(rf_inv) + pp.calc_duration(gz_crush) + ti_wait seq.paper_plot(time_range=(start, start + 3 * spoke_tr)) .. image-sg:: /generated/gallery/01-course/images/sphx_glr_06_radial_mprage_002.png :alt: 06 radial mprage :srcset: /generated/gallery/01-course/images/sphx_glr_06_radial_mprage_002.png :class: sphx-glr-single-img .. GENERATED FROM PYTHON SOURCE LINES 259-268 Signal evolution ---------------- The simulation below follows the longitudinal magnetisation of three tissues through the scan: inverted, recovering for TI, tipped by 8 degrees and recovering for one spoke TR at a time along the train, then recovering until the next inversion. After a few inversions it repeats exactly; the figure shows that steady period, with the signal of every spoke (:math:`M_z \sin\alpha` just before each pulse). .. GENERATED FROM PYTHON SOURCE LINES 268-312 .. image-sg:: /generated/gallery/01-course/images/sphx_glr_06_radial_mprage_003.png :alt: 06 radial mprage :srcset: /generated/gallery/01-course/images/sphx_glr_06_radial_mprage_003.png :class: sphx-glr-single-img .. GENERATED FROM PYTHON SOURCE LINES 313-330 A Cartesian MPRAGE puts the centre of k-space in the middle of the train, and that single moment fixes the contrast of the image. Here every spoke passes through the centre, so each spoke has the contrast of its own moment. Cut the train into bins, as in the figure, and each bin is an image with its own inversion time. Because the angles advance by the golden angle, the spokes of every bin cover the plane evenly, just more sparsely: each bin is an undersampled radial image, which a constrained or subspace reconstruction recovers. This is how the same scan gives you a series of contrasts, for example to fit T1, without being accelerated by design. If you need a shorter scan, the stack of stars has one more axis to undersample: you can skip partitions along z and recover them with parallel imaging, exactly as on a Cartesian phase-encoding axis. The labels tell the reconstruction which partition and which spoke each readout is; the time of each spoke in the train, and so its bin, follows from ``LIN``: .. GENERATED FROM PYTHON SOURCE LINES 330-335 .. code-block:: Python counters = seq.evaluate_labels(evolution="adc") print(f"{len(counters['LIN'])} readouts") print("PAR of the first train:", np.unique(counters["PAR"][:SPOKES_PER_TRAIN])) print("LIN of the first train:", counters["LIN"][:6], "...") .. rst-class:: sphx-glr-script-out .. code-block:: none 1024 readouts PAR of the first train: [0] LIN of the first train: [0 1 2 3 4 5] ... .. GENERATED FROM PYTHON SOURCE LINES 336-356 The scan is the shipped :doc:`mprage_stack_of_stars3D ` in miniature. In lesson 9 you rebuild it from modules. Specification ------------- What this lesson built, stated the way you would ask an agent for it: .. code-block:: text Write a pypulseqpp 3D radial MPRAGE as a stack of stars: 220 mm FOV, 128 matrix, 16 partitions of 5 mm, 64 golden-angle spokes per inversion (one partition per inversion, the angle advancing across the whole scan), 8 degree slab-selective sinc, 10 us dwell, merged prewinder and readout, TI 900 ms from the inversion centre to the first excitation, one inversion every 2.5 s with a 10 ms hyperbolic-secant pulse and a z crusher. Combine the slab rephaser and the partition encoding in one z gradient per partition, rewind kz in the spoiler block, label LIN and PAR, RF spoil at 117 degrees, and check the timing. .. rst-class:: sphx-glr-timing **Total running time of the script:** (0 minutes 1.312 seconds) .. _sphx_glr_download_generated_gallery_01-course_06_radial_mprage.py: .. only:: html .. container:: sphx-glr-footer sphx-glr-footer-example .. container:: sphx-glr-download sphx-glr-download-jupyter :download:`Download Jupyter notebook: 06_radial_mprage.ipynb <06_radial_mprage.ipynb>` .. container:: sphx-glr-download sphx-glr-download-python :download:`Download Python source code: 06_radial_mprage.py <06_radial_mprage.py>` .. container:: sphx-glr-download sphx-glr-download-zip :download:`Download zipped: 06_radial_mprage.zip <06_radial_mprage.zip>` .. only:: html .. rst-class:: sphx-glr-signature `Gallery generated by Sphinx-Gallery `_