1 of 66

Skope User Meeting 2024

Workshops Presentations

The 7th Skope User Meeting was held at the Max Planck Institute for Biological Cybernetics in Tübingen (Germany) on September 3, 2024. Measurements were performed on the institute’s 9.4 T MRI scanner.

Title slide

2 of 66

Overview

  • Gradient system characterization by impulse response measurements with a Skope Field Camera
    • This workshop showcased a rapid and highly sensitive approach for characterizing the dynamic performance of an MR gradient system. Leveraging the open-source framework Pulseq (available at pulseq.github.io), we generated gradient waveforms and captured the resulting fields using a Field Camera. This method allows valuable insights into system behavior, aiding applications such as pre-emphasis optimization, quality assurance, and image reconstruction.
  • Gradient-transfer function estimation
    • In this workshop, we delved into the step-by-step derivation of a gradient-transfer function (GTF) based on data collected during Workshop 1. We investigated gradient dynamics and had a look at subtle features such as eddy currents and mechanically induced field oscillations.
  • Image reconstruction of measured and predicted trajectories
    • The final workshop in this series highlighted the reconstruction of both field-monitored and GTF-predicted trajectories. The derived scanner model was applied to predict the trajectories of an imaging sequence based on nominal gradient waveforms. We demonstrated how to implement gridding reconstruction in MATLAB using a sparse matrix representation. Additionally, we investigated the strengths and limitations of GTF-based trajectory prediction.

Demo sessions

The User Meeting comprised three workshop sessions. The first session was performed on the 9.4 T MRI scanner. The two subsequent sessions evaluated and used the data acquired in the first session.

Text 100%

3 of 66

Theory

Gradient system characterization

This section gives background information on the theory of gradient characterization based on impulse response measurements.

Title slide

4 of 66

System Characterization

Gradient transfer function

Output

Input

GTF

To obtain a Gradient Transfer Function (GTF), we need to probe the gradient system. As input, we select gradient waveforms that sufficiently sample the relevant frequency range of the scanner's response. By comparing input and output, the GTF can be calculated.

Text 100%

5 of 66

Probing a Linear Time-Invariant System

  •  

Theory and methods

measured output

nominal input

GIRF

GTF

ij(t)

i1(t)

The GTF calculation is usually performed using a least-squares fit.

Text 100%

6 of 66

Application Examples

Gradient system characterization

Title slide

7 of 66

Why use a GIRF?

Characterizing real-world system performance

Characterize gradient development

Improve gradient safety

Improve Sequence performance

Monitor hardware interactions

ISMRM 2024 #4073

Le Ster et al

Comparison at 7T of the Impulse

head gradient and whole-body SC72

gradient transfer functions

ISMRM 2024 #0946

McCreedy et al

Safe Spirals for Your Scanner

ISMRM 2024 #4775

Loecher et al

Removing Background Velocity Errors in

PC-MRI with Optimized Spoiler

Gradient Waveforms

Boulant et al

Commissioning of the Iseult CEA 11.7 T

whole-body MRI: current status,

gradient–magnet interaction tests

and first imaging experience

DOI: 10.1007/s10334-023-01063-5

The characterization of real-world systems has many applications. Some recent examples are shown on the current slide.

Text 100%

8 of 66

Pulseq

Open-source sequence development tool

Pulseq was used to create the MR sequences that were played out for the impulse response measurements.

Title slide

9 of 66

Pulseq

  • Framework to easily implement pulse sequences

  • Available on GitHub
    • git clone git@github.com:pulseq/pulseq.git
  • Interpreter sequences available for Siemens and GE
  • Sequences used in this demonstration are also available on GitHub:

Open-source framework for pulse sequence programming

Pulseq is an open-source framework for pulse sequence programming. MATLAB or Python can be used to create text-based sequence files that are interpreted by the scanner software. Pulseq is available for Siemens and GE. The sequences used in this workshop can be downloaded from Skope’s GitHub page.

Text 100%

10 of 66

Measurement Setup

Field Monitoring

This section describes the setup used for performing the impulse response measurements.

Title slide

11 of 66

Hardware Setup

Field monitoring system

11

The hardware setup consisted of a Clip-on Camera mounted on a Plexiglas frame. The frame was positioned in front of a head coil that was connected to the scanner. The trigger output of the 9.4T scanner was connected to the Skope Camera Acquisition System. skope-fx was used to set the acquisition parameters for field monitoring and acquiring the data.

Frame used for field measurements.

Text 100%

12 of 66

Impulse Response Measurements

MR Sequences

This section describes the MR sequences that were used during Workshop 1 to perform the impulse response measurements. Note that the sequences are prototypes and might not be optimally suited for the intended task.

Title slide

13 of 66

Blip sequence

  • Properties
    • 20 blips per axis
    • TR 1s
    • 5 averages (innermost loop)

  • Acquisition settings
    • CameraAcqDuration 70 ms
    • CameraAqDelay 0 ms
    • CameraInterleaveTR 400 ms
    • CameraNrDynamics 300
    • CameraNrSyncDynamics 0

Impulse-response measurements (skope_gtf.seq)

Since real-word systems cannot play out delta functions, the most commonly used shapes for impulse response measurements are gradient blips. The used Pulseq sequence played out out 20 different blips with varying amplitude on each physical gradient axis. To improve SNR, each waveform was played out five times.

Text 100%

14 of 66

Chirp sequence

  • Properties
    • 1 sweep per axis
    • TR 1s
    • 100 averages (innermost loop)

  • Acquisition settings
    • CameraAcqDuration 70 ms
    • CameraAqDelay 0 ms
    • CameraInterleaveTR 400 ms
    • CameraNrDynamics 300
    • CameraNrSyncDynamics 0

Impulse-response measurements (skope_sweep.seq)

An alternative to the gradient blips are chirp pulses. As for the gradient blips, multiple averages were acquired to improve SNR.

Note that the acquisition parameters for the field monitoring experiments can also be found in the headers of the Pulseq sequences.

Text 100%

15 of 66

(Variable) Frequency Sweeps

  • Chirps
    • Broad-band input pulses
    • Essentially flat energy distribution in the frequency
  • Variable sweep speed
    • Account for noise impulse response measurement device (i.e. field camera)
  • Design and evaluation code available on GitHub:

Chirps

Chirp pulses allow to achieve a flat energy distribution in the frequency domain. The sweeping speed can be modulated to account for the varying noise level of the measuring device. The used shape was designed with the code available on GitHub.

Both types of measurement can be used in combination for the calculation of the GTF.

Text 100%

16 of 66

Linearity check

  • Properties
    • 20 blips per axis
    • Full and half amplitude
    • TR 1s
    • 1 average

  • Acquisition settings
    • CameraAcqDuration 70 ms
    • CameraAqDelay 0 ms
    • CameraInterleaveTR 400 ms
    • CameraNrDynamics 120
    • CameraNrSyncDynamics 0

Linearity check (skope_gtf_linearityCheck.seq)

Since the core assumption for the creation of a GIRF is the linearity and time-invariance of the system, it is worthwhile to validate this expectation via a measurement. The linearity check entails playing out two gradient blips with the same length but half the amplitude/slew rate.

Text 100%

17 of 66

Blip sequence

  • Play out gradients with the same width but half the amplitude/slew rate
  • System is linear if outputs are scaled correctly

Linearity check (skope_gtf_linearityCheck.seq)

1

1/2

1

x2

1/2

The system is linear if the measured amplitudes scale correctly.

Text 100%

18 of 66

Gradient Transfer �Function Estimation

Processing steps

This section describes the in Workshop 2 presented steps to generate a Gradient Transfer Function.

Title slide

19 of 66

Processing steps

Workshop 2

Data�pre-processing

Field-drift correction

Expansion to spherical harmonics

Calculation of actual gradient-waveform

Simulation of nominal gradient waveform

Workshop 1�Data acquisition

Workshop 3

Image reconstruction

GTF calibration

Workshop 2�Gradient-transfer function estimation

Calculation of �k-space

The estimation involves 4 steps

  1. k-space calculation from the acquired raw data
  2. Calculation of the actual gradient waveform from the k-space
  3. Simulation of the nominal gradient waveform
  4. Least-squares fit of the nominal and measured data to generate the GTF

Text 100%

20 of 66

Processing steps

Workshop 2

  •  

Data�pre-processing

Field-drift correction

Expansion to spherical harmonics

Calculation of actual gradient-waveform

Simulation of nominal gradient waveform

GTF calibration

Calculation of �k-space

To calculate the k-space from the raw-data, several processing steps are required.

Text 100%

21 of 66

Processing steps

Workshop 2

Data�pre-processing

Field-drift correction

Expansion to spherical harmonics

Calculation of �k-space

  • Perform linear fit on first data samples
    • Remove drift & offset via linear regression on early samples

Calculation of actual gradient-waveform

Simulation of nominal gradient waveform

GTF calibration

Gradient

Accrued Phase

fitting�region

Measured phase (with drift)

Corrected phase

During the GTF measurements, field-drift is likely to occur (e.g. because of gradient heating). The resulting drift must be accounted for. For this purpose, a 2 ms gradient-free interval was added after the trigger event and before the gradient blips.

Text 100%

22 of 66

Processing steps

Workshop 2

Data�pre-processing

Field-drift correction

Expansion to spherical harmonics

Calculation of �k-space

Higher Order k-space� k0 kx ky kz kxy kyz … kx³-3xy²

Calculation of actual gradient-waveform

Simulation of nominal gradient waveform

GTF calibration

Given the corrected and averaged phase evolution of each probe, a spherical harmonics expansion can be performed to obtain the (higher-order) k-space.

Text 100%

23 of 66

Processing steps

Workshop 2

  •  

Data�pre-processing

Field-drift correction

Expansion to spherical harmonics

Calculation of actual gradient-waveform

Simulation of nominal gradient waveform

GTF calibration

Calculation of �k-space

The actual gradient can be obtained by temporal derivative of the measured k-space.

The nominal gradient can be obtained from the PulSeq file.

Text 100%

24 of 66

Gradient Input and Output

Nominal gradient input

This slide presents both the nominal and measured gradient waveforms. The scanner’s gradient raster time is 10 µs. To align with the field camera’s temporal resolution of 1 µs, the resolution of the input waveforms has been enhanced. This adjustment results in the stair-step appearance of the input waveforms.

Text 100%

25 of 66

Nominal Gradient Waveforms

https://github.com/pulseq/pulseq/blob/master/src/ExternalSequence.cpp

  • mexSequenceSimulator based on ExternalSequence.cpp
  • Mex provides
    • Gradient waveform [Nx3] [mT/m]
    • Transmit pulse center [Mx1] [s]
    • Receiver events [3xL]
      • First data sample [s]
      • Number of samples
      • Dwell time [s]
    • Labels: [SLC, SEG, REP, AVG, ECO, PHS, SET, LIN, PAR]
    • Flags: [NAV, REV, SMS, PMC]
    • Pulseq definitions

Pulseq

A mex-tool allows for generating the nominal gradient waveform from a .seq file.

Text 100%

26 of 66

Processing steps

Workshop 2

  •  

Data�pre-processing

Field-drift correction

Expansion to spherical harmonics

Calculation of actual gradient-waveform

Simulation of nominal gradient waveform

GTF calibration

Calculation of �k-space

A least squares find is used to generate the GTF H from the actual and nominal gradient data.

Text 100%

27 of 66

Linear & Time-invariant System assumption

Workshop 2

  • Is my system time-invariant?
    • No (most likely)
    • But depending on your application, this might not be a problem for you

  • Is my system linear?
    • Let’s check!
    • Play gradients with same width, but half the amplitude/slew rate
    • Expectation for linear system: ½ -Blip output scaled by 2 should be similar to 1-Blip output

1

1/2

1

1/2

Inputs

x2

Output

reality

expectation for linear system

x2

The underlying assumption for the use of a GTF is that the system is linear and time-invariant.

Text 100%

28 of 66

Trajectory Prediction

GTF application

After having created the GTF, it can be used to predict trajectories.

Title slide

29 of 66

Probing a Linear Time-Invariant System

  •  

Theory and methods

measured output

nominal input

GIRF

GTF

B0

X

Y

Z

B0

X

Y

Z

B0

X

Y

Z

XY

ZY

Output channels ->

Input channels ->

Frequency samples ->

H

As described above, the GTF is represented by a three-dimensional matrix. One dimension represents the physical inputs (X,Y,Z), another the frequency or time domain and yet another the output fields. The number of output fields is determined by the impulse response measurements. The field cameras can easily measure up to 16 field terms.

Text 100%

30 of 66

Prediction of Gradient Waveforms

Nominal input

Prediction

range

Filter run-in samples

GIRF length

Once the nominal waveforms are known, they can be convolved with the GTF. In case of long waveforms, it is more efficient to perform the prediction for a reduced range only. Note that sufficient samples before the prediction range need to be included in the prediction. These run-in samples must be discarded and cannot be used for image reconstruction.

Text 100%

31 of 66

Prediction of Gradient Waveforms

  • Convolve nominal gradient waveforms with GIRF
  • Integrate waveforms to obtain k-space trajectory

Nominal input

 

The k-space trajectories are calculated based on the predicted gradient waveforms via integration.

Text 100%

32 of 66

Predict of k-Space Trajectories

GIRF-based trajectory prediction

k-space needs to be set to zero

  • at center of RF pulse
  • or at trigger location to imitate field camera measurement

Trigger

RF center

The integrated waveforms need to be set to zero at the center of the transmit pulses. If the field camera acquisition should be mimicked, the trajectories must be set to zero at the trigger locations.

Text 100%

33 of 66

Predict of k-Space Trajectories

GIRF-based trajectory prediction

Trigger

RF center

Filter run-in samples

Once the run-in samples have been discarded, the trajectory can be interpolated to the ADC raster time. The ADC properties (dwell time, location, number of samples) can be extracted from the Pulseq definitions.

Text 100%

34 of 66

Data Cropping

GIRF-based trajectory prediction

ADC event

Echo 1

ADC event

Echo 2

In this example, the readout consisted of two contrasts (i.e. two echo times), which will provide the necessary information to create static off-resonance maps.

Text 100%

35 of 66

Nominal Gradient Simulation

MATLAB code

The evaluation of the data was done in MATLAB. Once input and output folders are defined, the nominal gradient waveforms can be calculated. Note that Pulseq flips the polarity of the x-axis, which is undone in line 64. The unit of the gradient amplitude is also changed from mT/m to T/m.

Text 100%

36 of 66

Gradient and k-space Prediction

MATLAB code

The predicted and measured k-space trajectories match well. However, small differences can be observed. There are many possible reasons for the seen discrepancies:

  • The GFT and the GRE data were acquired on different days
  • The scanner was not in the same state for both measurements (cold vs warm)

Text 100%

37 of 66

Conversion to MRD format

MATLAB code

For image reconstruction, it is advantageous to store the data in a vendor-independent data format. The provided MATLAB function calls the siemens_to_ismrmrd application with a more elegant interface.

runSiemensConverter.m

Text 100%

38 of 66

MRD Format

  • Vendor-independent
  • Standardized MR raw data format
    • Extendible XML header
    • Defined by XSL style sheets
    • Acquisition header
    • with fixed memory layout
  • Available for Windows and Linux
  • Description and source code:
    • https://ismrmrd.github.io/
  • Converters for major MRI vendors:
    • https://github.com/ismrmrd
  • Viewer

Open-source data format

MRD is a data format proposed by the ISMRM. It consists of a flexible XML header and a fixed layout for the acquisitions. The data conversion can be adjusted via parameter maps and style sheets. The parameter maps translate vendor into MRD parameters, while the stylesheets define the structure of the MRD XML header.

Text 100%

39 of 66

Addition of Trajectory to MRD File

  • On the previous slides, only a single acquisition was predicted
  • For image reconstruction, the entire trajectory needs to be available

Open-source data format

Note the timeShift parameter that is used to shift the predicted data with respect to the nominal gradient input. During GTF creation, data delays (Camera hardware delay, scanner trigger delay, etc) were not considered. This delay is important to align the ADC events correctly with the predicted data.

Text 100%

40 of 66

Addition of Trajectory to MRD File

Open-source data format

While adding the k-space data to the file, the MRD header is populated with additional data. Parameters such as the trajectory description are optional and only needed if the data is reconstructed with skope-i.

addTrajectoryToMrdFile.m

Text 100%

41 of 66

Addition of Trajectory to MRD File

Open-source data format

Pulseq sequences do not always add the correct labels and flags to the data. Therefore, the acquisition headers are updated as well.

addTrajectoryToMrdFile.m

Text 100%

42 of 66

Addition of Trajectory to MRD File

Open-source data format

Once the MRD file is updated, the data corresponding to a specific slice can be read from the file.

addTrajectoryToMrdFile.m

Text 100%

43 of 66

Image Reconstruction

Predicted and measured trajectories

This section describes the gridding image reconstruction process.

Title slide

44 of 66

Image reconstruction

Cartesian Imaging

Data sample

Fourier Transform

If the acquired data is located on a Cartesian grid, image reconstruction is straightforward and only requires an FFT.

Text 100%

45 of 66

Image reconstruction

Non-Cartesian Imaging

Data sample

?

If this is not the case, the image reconstruction becomes a bit more complicated.

Text 100%

46 of 66

Gridding Reconstruction

Non-Cartesian Imaging

Kaiser-Bessel function

Data sample to grid

Data sample

The most commonly used method for non-Cartesian image reconstruction is Kaiser-Bessel gridding. Each data point is interpolated to its surrounding grid locations. The closer the data point is to a grid location, the higher its contribution is to that k-space position.

Text 100%

47 of 66

Bookkeeping

Gridding Reconstruction

1

Index of non-Cartesian data point (i)

Index of position in Cartesian k-space (J)

Gridding kernel values (S)

1

21

0.5

1

27

0.1

1

22

0.4

1

….

….

….

….

2

21

0.3

2

30

0.05

1

2

3

4

5

6

7

8

9

10

11

12

13

14

15

16

17

18

19

20

21

22

23

24

25

26

27

28

29

30

31

32

33

34

35

36

The entire gridding process can be represented in a large table, which summarizes how each data point affects the underlaying grid.

Text 100%

48 of 66

Gridding Operation

Sparse matrix implementation

Non-Cartesian data

Cartesian data

Gridding weights

The aforementioned table can be represented as a matrix holding the gridding weights. The matrix is very sparse since each data point can only affect a very limited number of grid points. Image reconstruction simplifies in this way to a matrix vector multiplication, where the vector represents the non-Cartesian data.

Text 100%

49 of 66

Gridding Operation

Sparse matrix implementation

Fourier Transform

The resulting vector can be rearranged and Fourier-transformed to obtain the image. For iterative image reconstruction algorithms, the inverse operation also needs to be performed. This can easily be done by using the transposed sparse matrix.

Text 100%

50 of 66

How Sparse is the Matrix?

Very!

Number of data samples ->

Number of image pixels ->

This figure demonstrates how sparse the matrix is. Due to its sparsity, it is possible to hold the matrix in memory.

Text 100%

51 of 66

Image Reconstruction

MATLAB

The first step for the image reconstruction in MATLAB consists in setting up the gridding operator, i.e. the sparse matrix holding the gridding weights.

Text 100%

52 of 66

Image Reconstruction

MATLAB

As shown above, the image reconstruction only consists of a matrix multiplication for each coil channel. Note that a deapodization needs to be applied to the final image to eliminate the effect of the Kaiser-Bessel interpolation.

reconstructImage.m

Text 100%

53 of 66

Image Reconstruction

MATLAB

The rooted sum-of-squares algorithm can be used to combine the single-channel images.

Text 100%

54 of 66

Image Reconstruction

MATLAB

The same exact procedure can be used to process the EPI data.

Text 100%

55 of 66

Image Reconstruction

MATLAB

  • The GTF does not include any delays between the gradients and the receiver chain
  • Data delays that are caused by the receiver chain must either be calculated based on the sequence protocol (dwell time, bandwidth, etc) if possible or calibrated in some other way
  • In case of field monitoring, data delays between the camera and the scanner data are estimated via a synchronization pre-scan. In case of GTF-predicted trajectories, there is no such pre-scan that could be used for this purpose.

In contrast to GRE, EPI is very sensitive to timing delays as the latter create ghosting artifacts. In this example, the delay was adjusted manually to provide images with the least amount of ghosting artefacts.

Text 100%

56 of 66

Comparison with skope-i

Gradient echo images

  • Iterative gridding reconstruction

The MRD files can also be used by skope-i to perform image reconstruction. This image shows the five GRE slices as reconstructed by Skope’s image reconstruction framework which uses an iterative conjugate-gradient algorithm.

Text 100%

57 of 66

Comparison with skope-i

B0 and sensitivity maps

  • B0 maps for all slices

The two echoes of the GRE can be used to create static off-resonance maps. Sensitivity maps are created based on the first echo.

  • Sensitivity maps (magnitude and phase) for one slice and all 32 channels

Text 100%

58 of 66

Eddy current simulation

B0 component

One detail that was not considered during the previous reconstructions is the scanner’s B0 eddy current compensation. When switching gradients, a B0 eddy current is created. This eddy current is compensated for by the scanner by adjusting the receiver frequency.

Text 100%

59 of 66

Eddy current simulation

B0 component

 

 

Gradient

Slew rate

Exponential ECC model

ECC B0 Term

The correction is based on the nominal gradient waveforms and a sum of decaying exponentials.

Text 100%

60 of 66

Eddy current simulation

B0 component

 

Phase

Phase during ADC

 

Once the B0 modulation is known, the corresponding phase is added by the scanner to the raw data. Since this correction cannot be switched off, it needs to be undone before the k0 measured by the camera can be used for image reconstruction. If the correction is not undone, k0 should be excluded from image reconstruction.

Text 100%

61 of 66

Comparison with skope-i

EPI images

  • Iterative gridding reconstruction without B0 correction and SENSE

Performing an image reconstruction without SENSE and B0 correction produces the EPI images shown on this slide.

Text 100%

62 of 66

Comparison with skope-i

EPI images

  • Iterative gridding reconstruction with B0 correction and SENSE

Adding B0 correction (and SENSE) improves image quality greatly by reducing distortions.

Text 100%

63 of 66

Scripts and Demo Data

Additional resources

63

Should you have any further questions feel free to reach out to your Skope representative.

Text 100%

64 of 66

Acknowledgments

  • We thank Praveen Valsala, Jonas Bause, Dario Bosch, Tina Schröder and Klaus Scheffler for hosting the Skope User Meeting in Tübingen.
  • We thank Johanna Vannesjo and Felix Breuer for the GTF and image reconstruction code.

Text 100%

65 of 66

Disclaimer

This presentation is intended exclusively for the recipient.

The Skope products (NeuroCamTM, Dynamic Field Camera, Clip-on Camera, skope-i, skope-dm) are not medical devices. The statements made regarding the products have not been evaluated by the Food and Drug Administration. The safety and efficacy of the products have not been confirmed by FDA- approved research. The products are not intended to diagnose, treat, cure, or prevent any disease.​

Under no circumstances is the content or parts of the content authorized to be forwarded or to be made available in any form, be it orally or written to third parties. Moreover, the content or parts of the content may not be used in any form without the permission by Skope Magnetic Resonance Technologies AG.

Copyright © 2011-2024 Skope Magnetic Resonance Technologies AG. All data and information contained in this presentation are legally not binding and shall not create any warranties or liabilities whatsoever of Skope.

Text 100%

66 of 66

Text 100%