.. DO NOT EDIT. .. THIS FILE WAS AUTOMATICALLY GENERATED BY SPHINX-GALLERY. .. TO MAKE CHANGES, EDIT THE SOURCE PYTHON FILE: .. "generated/gallery/02-tours/01_fid_and_spin_echo.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_02-tours_01_fid_and_spin_echo.py: ================= FID and spin echo ================= The simplest signal you can acquire is the free induction decay (FID): a pulse tips the magnetisation into the transverse plane and the receiver records what the precessing magnetisation emits. This Tour builds that sequence, checks the flip angle of its pulses with a Bloch simulation, and then adds a refocusing pulse and a pair of crusher gradients so that the acquisition is centred on a spin echo at the echo time you prescribe. You will see how the flip angle sets the FID amplitude, how the echo time is timed from pulse centres rather than block edges, and how the crusher dephasing sets the shortest echo time. **Prerequisites:** lesson 1 of the :doc:`course `. The representation these objects belong to is described in :doc:`/explanations/pulseq-file`. .. GENERATED FROM PYTHON SOURCE LINES 22-32 .. GENERATED FROM PYTHON SOURCE LINES 33-40 Free induction decay -------------------- A free induction decay needs two things: an RF pulse and an acquisition window. The scanner limits, the dead times and the rasters the events are designed against are in the ``system`` object, as in :doc:`lesson 1 `. .. GENERATED FROM PYTHON SOURCE LINES 40-54 .. code-block:: Python import numpy as np import pypulseqpp as pp 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, ) .. GENERATED FROM PYTHON SOURCE LINES 55-61 The pulse is a *block pulse*, a rectangle of constant amplitude, 200 µs long. You ask for a flip angle and pypulseqpp works out the amplitude. You make four of them, at 10, 30, 60 and 90 degrees, to see how the signal follows the flip angle. The ``delay`` starts each pulse after the RF dead time, and ``use="excitation"`` tells an interpreter and the k-space analysis what the pulse is for. .. GENERATED FROM PYTHON SOURCE LINES 61-74 .. code-block:: Python FLIP_ANGLES_DEG = (10.0, 30.0, 60.0, 90.0) pulses = [ pp.make_block_pulse( flip_angle=np.deg2rad(flip_angle), duration=200e-6, delay=system.rf_dead_time, system=system, use="excitation", ) for flip_angle in FLIP_ANGLES_DEG ] .. GENERATED FROM PYTHON SOURCE LINES 75-79 The ADC window records 8192 samples. The time between two samples is the *dwell time*, and its inverse is the receiver bandwidth. You ask for a dwell of 32 µs; :func:`~pypulseqpp.calc_adc_timing` moves it onto the ADC raster and makes the window end on the gradient raster. .. GENERATED FROM PYTHON SOURCE LINES 79-94 .. code-block:: Python dwell, acquisition = pp.calc_adc_timing( 8192, 32e-6, grad_raster_time=system.grad_raster_time, adc_raster_time=system.adc_raster_time, ) adc = pp.make_adc( num_samples=8192, dwell=dwell, delay=system.adc_dead_time, system=system ) print( f"dwell {dwell * 1e6:.0f} us, receiver bandwidth {1e-3 / dwell:.1f} kHz, " f"acquisition {acquisition * 1e3:.1f} ms" ) .. rst-class:: sphx-glr-script-out .. code-block:: none dwell 40 us, receiver bandwidth 25.0 kHz, acquisition 327.7 ms .. GENERATED FROM PYTHON SOURCE LINES 95-98 One repetition is two blocks: the pulse, then the acquisition. The 500 ms delay in the acquisition block is longer than the ADC window, so it sets the block duration and with it the repetition time. .. GENERATED FROM PYTHON SOURCE LINES 98-106 .. code-block:: Python seq = pp.Sequence(system=system) for pulse in pulses: seq.add_block(pulse) seq.add_block(adc, pp.make_delay(500e-3)) ok, errors = seq.check_timing() print(f"timing {ok}, {seq.num_blocks} blocks, {seq.duration()[0]:.2f} s") .. rst-class:: sphx-glr-script-out .. code-block:: none timing True, 8 blocks, 2.00 s .. GENERATED FROM PYTHON SOURCE LINES 107-109 The solid trace of the diagram is the first repetition; the shaded traces are the other three, which differ only in the pulse amplitude. .. GENERATED FROM PYTHON SOURCE LINES 109-111 .. code-block:: Python seq.paper_plot() .. image-sg:: /generated/gallery/02-tours/images/sphx_glr_01_fid_and_spin_echo_001.png :alt: 01 fid and spin echo :srcset: /generated/gallery/02-tours/images/sphx_glr_01_fid_and_spin_echo_001.png :class: sphx-glr-single-img .. GENERATED FROM PYTHON SOURCE LINES 112-119 Flip angle ---------- :func:`~pypulseqpp.sim_rf` runs a Bloch simulation of each pulse as it is stored in the sequence. On resonance, a block pulse rotates the magnetisation by its flip angle :math:`\alpha`, so the transverse magnetisation is :math:`|M_{xy}| = \sin\alpha`: .. GENERATED FROM PYTHON SOURCE LINES 119-135 .. code-block:: Python on_resonance = [] for pulse in pulses: mxy, frequency = pp.sim_rf(pulse)[1:3] on_resonance.append(abs(mxy[np.argmin(abs(frequency))])) .. image-sg:: /generated/gallery/02-tours/images/sphx_glr_01_fid_and_spin_echo_002.png :alt: 01 fid and spin echo :srcset: /generated/gallery/02-tours/images/sphx_glr_01_fid_and_spin_echo_002.png :class: sphx-glr-single-img .. GENERATED FROM PYTHON SOURCE LINES 136-139 The four simulated pulses sit on the curve: the 90 degree pulse puts all the magnetisation in the transverse plane, and the others a fraction :math:`\sin\alpha` of it. .. GENERATED FROM PYTHON SOURCE LINES 141-154 Spin echo --------- A FID starts decaying as soon as the pulse ends, and it decays faster than :math:`T_2` because spins at different positions precess at slightly different frequencies. A 180 degree *refocusing pulse* at :math:`\mathrm{TE}/2` reverses the phase each spin has accumulated, so the spins are in phase again at :math:`\mathrm{TE}`, and the ADC is centred there. The refocusing pulse is rectangular like the excitation and has the same duration, so twice the flip angle is twice the amplitude. The ``use`` tag lets the k-space analysis invert the accumulated gradient integral at the pulse centre, and lets an interpreter identify the pulse. .. GENERATED FROM PYTHON SOURCE LINES 154-178 .. code-block:: Python ECHO_TIME = 24e-3 VOXEL = 5e-3 CRUSHER_CYCLES = 4.0 excitation = pp.make_block_pulse( flip_angle=np.pi / 2, duration=200e-6, delay=system.rf_dead_time, system=system, use="excitation", ) refocusing = pp.make_block_pulse( flip_angle=np.pi, duration=200e-6, delay=system.rf_dead_time, system=system, use="refocusing", ) print( f"peak B1: excitation {abs(excitation.signal).max():.0f} Hz, " f"refocusing {abs(refocusing.signal).max():.0f} Hz" ) .. rst-class:: sphx-glr-script-out .. code-block:: none peak B1: excitation 1250 Hz, refocusing 2500 Hz .. GENERATED FROM PYTHON SOURCE LINES 179-196 Crusher gradients ----------------- A refocusing pulse with a flip angle that is not exactly 180 degrees also produces a FID, which would overlap the echo. Two *crusher* gradients of equal area and equal polarity, one before and one after the refocusing pulse, remove it. The pulse inverts the phase accumulated in the refocused pathway, so the second crusher unwinds what the first wound. A pathway that crosses the pulse without being inverted, such as that FID, sees the two areas add, and is left wound through ``CRUSHER_CYCLES`` cycles across a voxel, which integrates to nothing over it. So you prescribe a crusher by the phase it winds across a voxel: an area :math:`A` turns the magnetisation through :math:`A\,\Delta x` cycles across an extent :math:`\Delta x`. :func:`~pypulseqpp.make_crusher` takes the dephasing and the voxel and solves the shortest gradient that delivers it within the amplitude and slew limits. .. GENERATED FROM PYTHON SOURCE LINES 196-203 .. code-block:: Python crusher = pp.make_crusher(CRUSHER_CYCLES, VOXEL, channel="z", system=system)[0] print( f"crusher: {1e3 * pp.calc_duration(crusher):.2f} ms for " f"{CRUSHER_CYCLES:.0f} cycles across {1e3 * VOXEL:.0f} mm" ) .. rst-class:: sphx-glr-script-out .. code-block:: none crusher: 0.86 ms for 4 cycles across 5 mm .. GENERATED FROM PYTHON SOURCE LINES 204-216 Echo timing ----------- The echo time is measured between pulse centres, not block edges. A block longer than the events in it is a delay, and the two delays below absorb whatever the pulses, their dead and ringdown times and the crushers do not occupy: the refocusing pulse centre goes at :math:`\mathrm{TE}/2` and the acquisition centre at :math:`\mathrm{TE}`. You put each delay on the block-duration raster with :func:`~pypulseqpp.round_to_raster`; a duration off the raster is rounded up when the block is added, which moves the echo by up to one raster period. The rasters are described in :doc:`/explanations/pulseq-file`. .. GENERATED FROM PYTHON SOURCE LINES 216-283 .. code-block:: Python SAMPLES = 512 dwell, acquisition = pp.calc_adc_timing( SAMPLES, 20e-6, grad_raster_time=system.grad_raster_time, adc_raster_time=system.adc_raster_time, ) adc = pp.make_adc( num_samples=SAMPLES, dwell=dwell, delay=system.adc_dead_time, system=system ) print( f"dwell {dwell * 1e6:.0f} us, acquisition {acquisition * 1e3:.2f} ms, " f"half of it {acquisition * 5e2:.2f} ms" ) def centre(pulse): """Time from the start of a block to the centre of the pulse in it.""" return pulse.delay + pulse.shape_dur / 2 def delays(echo_time, crusher_event): """The two delay blocks, on either side of the refocusing pulse.""" crusher_duration = pp.calc_duration(crusher_event) before = ( pp.calc_duration(excitation) - centre(excitation) + crusher_duration + centre(refocusing) ) # The acquisition is centred on the echo, so what has to fit after the # refocusing pulse is the dead time and half the acquisition window. after = ( pp.calc_duration(refocusing) - centre(refocusing) + crusher_duration + system.adc_dead_time + acquisition / 2 ) raster = system.block_duration_raster return ( pp.round_to_raster(echo_time / 2 - before, raster), pp.round_to_raster(echo_time / 2 - after, raster), ) def spin_echo(echo_time, crusher_event): """One repetition of the spin echo, at the given echo time.""" first, second = delays(echo_time, crusher_event) seq = pp.Sequence(system=system) seq.add_block(excitation) seq.add_block(pp.make_delay(first)) seq.add_block(crusher_event) seq.add_block(refocusing) seq.add_block(crusher_event) seq.add_block(pp.make_delay(second)) seq.add_block(adc) return seq seq = spin_echo(ECHO_TIME, crusher) ok, errors = seq.check_timing() print(f"timing {ok}, {seq.num_blocks} blocks, {1e3 * seq.duration()[0]:.2f} ms") .. rst-class:: sphx-glr-script-out .. code-block:: none dwell 20 us, acquisition 10.24 ms, half of it 5.12 ms timing True, 7 blocks, 29.34 ms .. GENERATED FROM PYTHON SOURCE LINES 284-289 To check the timing, read the pulse centres back from the k-space analysis: the refocusing pulse should be at half the echo time after the excitation, and the midpoint of the acquisition window at the echo time. The analysis returns the time of each ADC sample, so the midpoint of the window is the mean of the first and the last of them. .. GENERATED FROM PYTHON SOURCE LINES 289-299 .. code-block:: Python _, _, t_excitation, t_refocusing, t_adc = seq.calculate_kspace() window_centre = (t_adc[0] + t_adc[-1]) / 2 print( f"prescribed TE/2 {1e3 * ECHO_TIME / 2:.3f} ms, " f"refocusing centre {1e3 * (t_refocusing[0] - t_excitation[0]):.3f} ms\n" f"prescribed TE {1e3 * ECHO_TIME:.3f} ms, " f"acquisition midpoint {1e3 * (window_centre - t_excitation[0]):.3f} ms" ) .. rst-class:: sphx-glr-script-out .. code-block:: none prescribed TE/2 12.000 ms, refocusing centre 12.000 ms prescribed TE 24.000 ms, acquisition midpoint 24.010 ms .. GENERATED FROM PYTHON SOURCE LINES 300-307 The refocusing pulse lands on the prescription and the acquisition midpoint lands half a block-raster period from it. Rounding a delay to the raster moves what follows it by at most half a period, so the echo time is realizable only to that resolution. The sequence diagram shows the crushers on the slice axis, which carries no other gradient in this non-selective experiment: .. GENERATED FROM PYTHON SOURCE LINES 307-309 .. code-block:: Python seq.paper_plot() .. image-sg:: /generated/gallery/02-tours/images/sphx_glr_01_fid_and_spin_echo_003.png :alt: 01 fid and spin echo :srcset: /generated/gallery/02-tours/images/sphx_glr_01_fid_and_spin_echo_003.png :class: sphx-glr-single-img .. GENERATED FROM PYTHON SOURCE LINES 310-319 Minimum echo time ----------------- Both delays shrink as the echo time is shortened, and the first of them to reach zero sets the shortest echo time. The crusher pair is the term in that budget you control: more dephasing is a longer gradient at the same amplitude limit, on both sides of the refocusing pulse. You find the shortest echo time for each dephasing, and check that the sequence built at that echo time passes the timing check. .. GENERATED FROM PYTHON SOURCE LINES 319-377 .. code-block:: Python CYCLES = (1.0, 2.0, 4.0, 8.0, 16.0, 32.0) def shortest_echo_time(crusher_event): """The smallest echo time at which neither delay block is negative.""" raster = system.block_duration_raster before, after = (-delay for delay in delays(0.0, crusher_event)) return raster * np.ceil(2 * max(before, after) / raster) shortest = [] for cycles in CYCLES: event = pp.make_crusher(cycles, VOXEL, channel="z", system=system)[0] echo_time = shortest_echo_time(event) ok, _ = spin_echo(echo_time, event).check_timing() shortest.append( { "cycles": cycles, "crusher": pp.calc_duration(event), "echo_time": echo_time, "timing": ok, } ) .. image-sg:: /generated/gallery/02-tours/images/sphx_glr_01_fid_and_spin_echo_004.png :alt: 01 fid and spin echo :srcset: /generated/gallery/02-tours/images/sphx_glr_01_fid_and_spin_echo_004.png :class: sphx-glr-single-img .. rst-class:: sphx-glr-script-out .. code-block:: none cycles crusher shortest TE timing 1 0.40 ms 11.28 ms True 2 0.56 ms 11.60 ms True 4 0.86 ms 12.20 ms True 8 1.46 ms 13.40 ms True 16 2.64 ms 15.76 ms True 32 5.02 ms 20.52 ms True .. GENERATED FROM PYTHON SOURCE LINES 378-404 Half of the acquisition window sits between the refocusing pulse and the echo, so the shortest echo time is bounded below by the acquisition duration whatever the crushers do, and the crusher pair is added to that floor. The duration of the pair grows as the square root of the dephasing while it is slew-limited and in proportion to it once the amplitude limit is reached, so the last doublings add the most to the echo time. Specification ------------- What this Tour built, stated the way you would ask an agent for it: .. code-block:: text Write a pypulseqpp pulse-acquire sequence for a scanner with 32 mT/m and 130 T/m/s, 100 us RF dead time, 20 us ringdown and 10 us ADC dead time. Play four 200 us block pulses at 10, 30, 60 and 90 degrees, each followed by an 8192-sample ADC at 32 us dwell, with a TR of 500 ms, and verify the flip angles with sim_rf against sin(alpha). Then write a spin echo: a 90 degree excitation and a 180 degree refocusing block pulse of 200 us, a z crusher of 4 cycles across 5 mm on either side of the refocusing pulse, and a 512-sample ADC at about 20 us dwell centred on TE = 24 ms. Put the two delays on the block raster, check the timing, read the pulse centres back with calculate_kspace, and find the shortest TE for 1 to 32 crusher cycles. .. rst-class:: sphx-glr-timing **Total running time of the script:** (0 minutes 0.360 seconds) .. _sphx_glr_download_generated_gallery_02-tours_01_fid_and_spin_echo.py: .. only:: html .. container:: sphx-glr-footer sphx-glr-footer-example .. container:: sphx-glr-download sphx-glr-download-jupyter :download:`Download Jupyter notebook: 01_fid_and_spin_echo.ipynb <01_fid_and_spin_echo.ipynb>` .. container:: sphx-glr-download sphx-glr-download-python :download:`Download Python source code: 01_fid_and_spin_echo.py <01_fid_and_spin_echo.py>` .. container:: sphx-glr-download sphx-glr-download-zip :download:`Download zipped: 01_fid_and_spin_echo.zip <01_fid_and_spin_echo.zip>` .. only:: html .. rst-class:: sphx-glr-signature `Gallery generated by Sphinx-Gallery `_