Skip to content
Merged
Show file tree
Hide file tree
Changes from all commits
Commits
File filter

Filter by extension

Filter by extension


Conversations
Failed to load comments.
Loading
Jump to
Jump to file
Failed to load files.
Loading
Diff view
Diff view
4 changes: 3 additions & 1 deletion pyproject.toml
Original file line number Diff line number Diff line change
Expand Up @@ -128,6 +128,8 @@ locale = "en-us"

[tool.typos.default.extend-words]
writeable = "writeable" # this is a valid numpy array attribute
Reson = "Reson" # required for J. Magn. Reson. Ser. A.
Ser = "Ser" # required for J. Magn. Reson. Ser. A.

# PyTest section
[tool.pytest.ini_options]
Expand All @@ -143,4 +145,4 @@ markers = ["sigpy: tests that require sigpy"]
# Ref: https://github.com/codespell-project/codespell#using-a-config-file
skip = '.git*'
check-hidden = true
ignore-words-list = 'te,fo,nd,maked,datas,iten,itel,tread'
ignore-words-list = 'te,fo,nd,maked,datas,iten,itel,tread,Reson,Ser'
4 changes: 2 additions & 2 deletions src/pypulseq/Sequence/block.py
Original file line number Diff line number Diff line change
Expand Up @@ -486,7 +486,7 @@ def get_block(self, block_index: int, add_IDs: bool = False) -> SimpleNamespace:
phase_shape = np.array([], dtype=float)

adc = SimpleNamespace()
adc.num_samples = lib_data[0]
adc.num_samples = int(lib_data[0])
adc.dwell = lib_data[1]
adc.delay = lib_data[2]
adc.freq_ppm = lib_data[3]
Expand All @@ -495,8 +495,8 @@ def get_block(self, block_index: int, add_IDs: bool = False) -> SimpleNamespace:
adc.phase_offset = lib_data[6]
adc.phase_modulation = phase_shape
adc.dead_time = self.system.adc_dead_time
adc.num_samples = int(adc.num_samples)
adc.type = 'adc'
adc.duration = adc.num_samples * adc.dwell

if add_IDs:
adc.id = raw_block.adc
Expand Down
4 changes: 4 additions & 0 deletions src/pypulseq/check_timing.py
Original file line number Diff line number Diff line change
Expand Up @@ -104,6 +104,10 @@ def div_check(a: float, b: float, event: str, field: str, raster: str):
div_check(e.delay, raster, event=event, field='delay', raster=raster_str)

if hasattr(e, 'duration'):
# the ADC delay must lie on the RF raster, but the duration (num_samples * dwell) on the ADC raster
if getattr(e, 'type', '') == 'adc':
raster = seq.system.adc_raster_time
raster_str = 'adc_raster_time'
div_check(e.duration, raster, event=event, field='duration', raster=raster_str)

if hasattr(e, 'dwell'):
Expand Down
47 changes: 32 additions & 15 deletions src/pypulseq/make_adc.py
Original file line number Diff line number Diff line change
@@ -1,4 +1,5 @@
import itertools
import numbers
from math import ceil, floor, gcd, isclose, prod
from types import SimpleNamespace
from typing import List, Tuple, Union
Expand Down Expand Up @@ -27,20 +28,20 @@ def make_adc(

Parameters
----------
num_samples: int
Number of readout samples.
system : Opts, default=Opts()
System limits. Default is a system limits object initialized to default values.
dwell : float, default=0
ADC dead time in seconds (s) after sampling.
duration : float, default=0
Duration in seconds (s) of ADC readout event with `num_samples` number of samples.
num_samples : int
Number of readout samples. Must be a positive integer.
delay : float, default=0
Delay in seconds (s) of ADC readout event.
duration : float, default=0
Duration in seconds (s) of ADC readout event with `num_samples` number of samples.
dwell : float, default=0
Dwell time in seconds (s) of ADC readout event, i.e. the time between samples.
freq_offset : float, default=0
Frequency offset in hertz (Hz) of ADC readout event.
phase_offset : float, default=0
Phase offset in radians (rad) of ADC readout event.
system : Opts, default=Opts()
System limits. Default is a system limits object initialized to default values.
freq_ppm : float, default=0
PPM frequency offset in parts per million (ppm) of ADC readout event.
phase_ppm : float, default=0
Expand All @@ -57,31 +58,47 @@ def make_adc(
Raises
------
ValueError
If `num_samples` is not a positive integer.
If `dwell`, `duration` or `delay` is negative or not finite.
If neither `dwell` nor `duration` are defined.
If `phase_modulation` does not have `num_samples` elements.
"""
if system is None:
system = Opts.default

if not isinstance(num_samples, numbers.Integral) or isinstance(num_samples, bool) or num_samples <= 0:
raise ValueError('num_samples must be a positive integer.')

if not np.isfinite(dwell) or dwell < 0:
raise ValueError('Dwell must be non-negative.')

if not np.isfinite(duration) or duration < 0:
raise ValueError('Duration must be non-negative.')

if not np.isfinite(delay) or delay < 0:
raise ValueError('Delay must be non-negative.')

if (dwell == 0 and duration == 0) or (dwell > 0 and duration > 0):
raise ValueError('Either dwell or duration must be defined')

if phase_modulation is not None and len(phase_modulation) != num_samples:
raise ValueError('ADC Phase modulation vector must have the same length as the number of samples')

adc = SimpleNamespace()
adc.type = 'adc'
adc.num_samples = num_samples
adc.num_samples = int(num_samples)
adc.dwell = dwell
adc.delay = delay
adc.freq_offset = freq_offset
adc.phase_offset = phase_offset
adc.freq_ppm = freq_ppm
adc.phase_ppm = phase_ppm
adc.dead_time = system.adc_dead_time

if (dwell == 0 and duration == 0) or (dwell > 0 and duration > 0):
raise ValueError('Either dwell or duration must be defined')

if phase_modulation is not None and len(phase_modulation) != num_samples:
raise ValueError('ADC Phase modulation vector must have the same length as the number of samples')
adc.phase_modulation = phase_modulation

if duration > 0:
adc.dwell = duration / num_samples
adc.duration = duration

if dwell > 0:
adc.duration = dwell * num_samples
Expand Down
112 changes: 59 additions & 53 deletions src/pypulseq/make_adiabatic_pulse.py
Original file line number Diff line number Diff line change
Expand Up @@ -40,60 +40,20 @@ def make_adiabatic_pulse(
"""
Make an adiabatic inversion pulse.

Note: some parameters only affect certain pulse types and are ignored for other; e.g. bandwidth is ignored if
type='hypsec'.

hypsec(n=512, beta=800, mu=4.9, dur=0.012)
Design a hyperbolic secant adiabatic pulse. `mu` * `beta` becomes the amplitude of the frequency sweep.

Args:
- n (int): number of samples (should be a multiple of 4).
- beta (float): AM waveform parameter.
- mu (float): a constant, determines amplitude of frequency sweep.
- dur (float): pulse time (s).

Returns
-------
2-element tuple containing
- **a** (*array*): AM waveform.
- **om** (*array*): FM waveform (radians/s).

References
----------
Baum, J., Tycko, R. and Pines, A. (1985). 'Broadband and adiabatic
inversion of a two-level system by phase-modulated pulses'.
Phys. Rev. A., 32:3435-3447.

wurst(n=512, n_fac=40, bw=40000.0, dur=0.002)
Design a WURST (wideband, uniform rate, smooth truncation) adiabatic inversion pulse

Args:
- n (int): number of samples (should be a multiple of 4).
- n_fac (int): power to exponentiate to within AM term. ~20 or greater is typical.
- bw (float): pulse bandwidth.
- dur (float): pulse time (s).

Returns
-------
2-element tuple containing
- **a** (*array*): AM waveform.
- **om** (*array*): FM waveform (radians/s).

References
----------
Kupce, E. and Freeman, R. (1995). 'Stretched Adiabatic Pulses for
Broadband Spin Inversion'.
J. Magn. Reason. Set. A., 117:246-256.
Note: some parameters only affect certain pulse types and are ignored for others; e.g. `beta` and `mu` are only
used if `pulse_type='hypsec'`, while `bandwidth` and `n_fac` are only used if `pulse_type='wurst'`.

Parameters
----------
pulse_type : str
One of 'hypsec' or 'wurst' pulse types.
adiabaticity : int, default=4
bandwidth : int, default=40000
Pulse bandwidth in hertz (Hz).
adiabaticity : float, default=4
Dimensionless adiabaticity factor determining the RF amplitude required
to maintain adiabatic behavior during the frequency sweep.
bandwidth : float, default=40000
Pulse bandwidth in hertz (Hz). Only used if `pulse_type='wurst'`.
beta : float, default=800.0
AM waveform parameter in radians per second (rad/s).
AM waveform parameter in radians per second (rad/s). Only used if `pulse_type='hypsec'`.
delay : float, default=0
Delay in seconds (s).
duration : float, default=10e-3
Expand All @@ -106,10 +66,10 @@ def make_adiabatic_pulse(
Maximum gradient strength (Hz/m).
max_slew : float, default=None
Maximum slew rate (Hz/m/s).
mu : float, default=4.9
Constant determining amplitude of frequency sweep.
n_fac : int, default=40
Power to exponentiate to within AM term. ~20 or greater is typical.
Power to exponentiate to within AM term. ~20 or greater is typical. Only used if `pulse_type='wurst'`.
mu : float, default=4.9
Constant determining amplitude of frequency sweep. Only used if `pulse_type='hypsec'`.
phase_offset : float, default=0
Phase offset in radians (rad).
return_gz : bool, default=False
Expand Down Expand Up @@ -142,23 +102,69 @@ def make_adiabatic_pulse(
ValueError
If invalid pulse type is encountered.
If invalid pulse use is encountered.
If `duration` is not positive or not finite.
If `delay` is negative or not finite.
If `dwell` is not positive or not finite.
If `adiabaticity` is not positive.
If `beta` or `mu` is not positive and `pulse_type='hypsec'`.
If `bandwidth` or `n_fac` is not positive and `pulse_type='wurst'`.
If slice thickness is not provided but slice-selective trapezoid event is expected.

Notes
-----
hypsec
Design a hyperbolic secant adiabatic pulse. `mu` * `beta` becomes the
amplitude of the frequency sweep. [1]_

wurst
Design a WURST (wideband, uniform rate, smooth truncation) adiabatic
inversion pulse. [2]_

References
----------
.. [1] Baum, J., Tycko, R. and Pines, A. (1985). 'Broadband and adiabatic
inversion of a two-level system by phase-modulated pulses'.
Phys. Rev. A., 32:3435-3447.
.. [2] Kupce, E. and Freeman, R. (1995). 'Stretched Adiabatic Pulses for
Broadband Spin Inversion'. J. Magn. Reson. Ser. A., 117:246-256.
"""
if system is None:
system = Opts.default

if return_gz and slice_thickness <= 0:
if return_gz and (not np.isfinite(slice_thickness) or slice_thickness <= 0):
raise ValueError('Slice thickness must be provided')

valid_pulse_types = ['hypsec', 'wurst']
if (not pulse_type) or (pulse_type not in valid_pulse_types):
raise ValueError(f'Invalid type parameter. Must be one of {valid_pulse_types}.Passed: {pulse_type}')
raise ValueError(f'Invalid type parameter. Must be one of {valid_pulse_types}. Passed: {pulse_type}')
valid_rf_use_labels = get_supported_rf_uses()
if use != '' and use not in valid_rf_use_labels:
raise ValueError(f'Invalid use parameter. Must be one of {valid_rf_use_labels}. Passed: {use}')

if not np.isfinite(duration) or duration <= 0:
raise ValueError('RF pulse duration must be positive.')

if not np.isfinite(delay) or delay < 0:
raise ValueError('Delay must be non-negative.')

if not np.isfinite(adiabaticity) or adiabaticity <= 0:
raise ValueError('adiabaticity must be positive.')

if pulse_type == 'hypsec':
if not np.isfinite(beta) or beta <= 0:
raise ValueError('beta must be positive.')
if not np.isfinite(mu) or mu <= 0:
raise ValueError('mu must be positive.')
elif pulse_type == 'wurst':
if not np.isfinite(bandwidth) or bandwidth <= 0:
raise ValueError('bandwidth must be positive.')
if not np.isfinite(n_fac) or n_fac <= 0:
raise ValueError('n_fac must be positive.')

if dwell is None:
dwell = system.rf_raster_time
elif not np.isfinite(dwell) or dwell <= 0:
raise ValueError('Dwell must be positive.')

# WTC and PS - we have no idea why eps is added here. Leaving for now.
n_raw = round(duration / dwell + eps)
Expand Down
43 changes: 29 additions & 14 deletions src/pypulseq/make_arbitrary_grad.py
Original file line number Diff line number Diff line change
Expand Up @@ -34,24 +34,24 @@ def make_arbitrary_grad(
channel : str
Orientation of gradient event of arbitrary shape. Must be one of `x`, `y` or `z`.
waveform : numpy.ndarray
Arbitrary waveform.
first : float
Arbitrary waveform. Must be a non-empty 1-D array.
first : float, default=None
Gradient value at the start of the gradient event. (t=0)
Will default to a linear extrapolated value if not provided.
last : float
last : float, default=None
Gradient value at the end of the gradient event. (t=duration)
Will default to a linear extrapolated value if not provided.
delay : float, default=0
Delay in seconds (s).
max_grad : float, default=None
Maximum gradient strength (Hz/m).
Will default to `system.max_grad` if not provided or 0.
max_slew : float, default=None
Maximum slew rate (Hz/m/s).
Will default to `system.max_slew` if not provided or 0.
system : Opts, default=Opts()
System limits.
Will default to `pypulseq.opts.default` if not provided.
max_grad : float
Maximum gradient strength (Hz/m).
Will default to `system.max_grad` if not provided.
max_slew : float
Maximum slew rate (Hz/m/s).
Will default to `system.max_slew` if not provided.
delay : float, default=0
Delay in seconds (s).
oversampling : bool, default=False
Boolean flag to indicate if gradient is oversampled by a factor of 2.

Expand All @@ -64,21 +64,36 @@ def make_arbitrary_grad(
------
ValueError
If invalid `channel` is passed. Must be one of x, y or z.
If `waveform` is not a non-empty 1-D array.
If `waveform` has less than 2 samples and `first` or `last` are not provided.
If `waveform` has an even number of samples and `oversampling` is enabled.
If `delay` is negative or not finite.
If slew rate is violated.
If gradient amplitude is violated.
"""
if system is None:
system = Opts.default

if channel not in ['x', 'y', 'z']:
raise ValueError(f'Invalid channel. Must be one of x, y or z. Passed: {channel}')

waveform = np.atleast_1d(np.squeeze(np.asarray(waveform)))
if waveform.ndim != 1:
raise ValueError(f'waveform should have ndim=1. Passed ndim={waveform.ndim}')
if waveform.size == 0:
raise ValueError('waveform must not be empty.')
if len(waveform) < 2 and (first is None or last is None):
raise ValueError('waveform must contain at least 2 samples unless `first` and `last` are provided.')

if not np.isfinite(delay) or delay < 0:
raise ValueError('Delay must be non-negative.')

if max_grad is None or max_grad == 0:
max_grad = system.max_grad

if max_slew is None or max_slew == 0:
max_slew = system.max_slew

if channel not in ['x', 'y', 'z']:
raise ValueError(f'Invalid channel. Must be one of x, y or z. Passed: {channel}')

if first is None or last is None:

def extrap(a, b):
Expand Down
Loading
Loading