A complete review of MRV fundamentals is outside the scope of this paper, however a short description is provided herein. Hydrogen protons – in this case within aqueous solutions – align with a strong magnetic field and are perturbed out of alignment by radio-frequency pulses. The perturbation causes detectable signal changes measured about the resonant frequency as the spins realign with the field. The application of temporally varying spatial magnetic field gradients causes the resonant frequency to vary in space and time, allowing the accrued phase of the spins to be sensitized to position and velocity. With 3D pulse sequences, the resultant signals can be processed for both position and velocity information in each coordinate direction. The overall duration of the excitation and readout of the signal is on the order of milliseconds. This section provides scan details for each of the five research teams who participated in the challenge. Table 1 summarizes the key equipment and parameters for the experiments.
3.1 Hanyang University (Hanyang) Scan Details
The dataset of the Hanyang University group was obtained using a 3T human MRI system located at Ochang campus of the Korea Basic Science Institute in Chungju, Korea. Commercial 4D flow MRI sequences embedded in a Philips MR scanner were utilized to obtain 4D pulsatile velocity data. The sequence used in the experiments consists of velocity encoding using a 6-point phase contrast method and 3D Cartesian orthogonal space encoding. The velocity encoding was performed as in Table 1with retrospective gating to measure periodic flows. The 5 to 20 mV TTL signal generated by a DAQ device (USB-6361 multifunctional I/O device, National Instruments) was transmitted to the MRI scanner through the wireless ECG device (MR400 ECG/Spo2 Transmitters, Invivo Expression) to control the scan in synchronization with the flow cycle. The 1.5 second period was divided into 20 temporal phases to obtain velocity data with a time resolution of 75 ms, and no interpolation between temporal phases was applied.
Table 1
Summary of Key Scan Details for Participating Teams
|
Hanyang University
|
Stanford/USMA
|
University of Rostock
|
Seoul National University
|
University of Illinois – Urbana Champaign
|
Scanner Type
|
3.0T Philips Achieva (TX)
|
3.0T GE Healthcare Discovery 750
|
3.0T Siemens Magnetrom Trio
|
3.0T Siemens Magnetom Trio
|
3.0T Siemens - Magnetom Prisma
|
Gradients [mT/m]
|
80
|
50
|
38
|
45
|
80
|
Max Slew Rate [mT/(m*ms)]
|
200
|
150
|
170
|
200
|
200
|
Coil
|
8-Channel knee
|
Single Channel Head
|
2x 3-channel body matrix coil
|
6-Channel Body
|
Two coils: 4-Channel Flex and 18 Channel Body Flex arrays
|
Velocity Encodes (VENC) – x, y, z [cm/s]
|
100, 100, 100
|
110, 120, 84
|
80, 80, 80
|
150,150,80
|
125, 125, 80
|
Voxel resolution – x, y, z [mm]
|
0.5 (isotropic)
|
0.8,0.7,0.7
|
0.75 (isotropic)
|
0.8 (isotropic)
|
0.67 (isotropic)
|
Raw Image Matrix Size
|
512 x 60 x 34
|
320 x124 x 60
|
448 x 140 x 72
|
320 x 48 x 48
|
384 x 144 x 48
|
Echo time, Repetition time (TE, TR) [ms]
|
2.5, 37.5
|
1.8, 5.4
|
3.33, 6.2
|
3.68, 8.21
|
7.37, 11.8375
|
Single Scan Times [mm:ss]
|
154:30
|
20:20 (Flow-On)
2:02 (Flow-Off)
|
84:00 (Flow-On)
4:00 (Flow-Off)
|
28:48
|
98:00
|
# Flow-On/off Scans
|
2/1
|
9/10
|
5/6
|
2/1
|
3/8
|
Gating Technique/ # of injection phases
|
Retrospective, 20 measured phases, 20 reconstructed phases
|
Retrospective, 10 measured phases, 20 reconstructed phases
|
Prospective, 20 measured phases, 20 reconstructed phases
|
Prospective, 22 measured phases, 22 reconstructed phases
|
Prospective, 31 measured phases, 31 reconstructed phases
|
Mainstream Pump max power [W]
|
350
|
370
|
1100
|
250
|
370
|
A 20 mM copper sulfate aqueous solution was used as the working fluid, and air bubbles in the solution were removed as much as possible using a vacuum pump before the experiment. The working fluid was maintained at 20.5 ℃ during the experiments by an icepack. The main flow was driven using a multi-stage centrifugal pump (PP2-20Y, Hwarang System) with a maximum output of 350 W and was monitored by a flowmeter (KTM-800, Kometer) located between the test section and the pump to maintain a constant flowrate. The injection flowrate was monitored by a differential pressure flowmeter (L-2LPM-D5V, Alicat) located between the test section and the pump. The voltage waveform was set to have maximum and minimum values of 1.64 V and 0.34 V according to the guidelines given to each group to try to standardize the experimental setup, and the entire pump system was controlled using the same DAQ device and software with the TTL gating signal for the synchronization of flow and MR scans. The full scans consisted of two flow-on and one flow-off conditions and were taken in order of on-off-on. To reduce the background phase error of the spatially corrected data, the process of averaging two flow-on data and subtracting the flow-off data was performed for each voxel.
The MRI raw data was reconstructed using the software included with the Philips MRI console. Although parallel imaging methods were not applied to reduce scan time, the zero-interpolation-filling (ZIP) was applied by a factor of 2 in the phase and slice direction. ZIP may increase errors due to missing information but can improve spatial resolution while preserving scan times (Zhu et al 2013). In addition, a house code, based on the algorithm of barrel distortion correction, was applied to correct the spatial distortion caused by the relatively large FOV setting in the reconstructed image. For an accurate velocity mask, the optimal threshold magnitude was determined per slice based on the given cross-sectional area. Finally, the wall-bounded divergence-free smoothing filter was applied to reduce the noise in the velocity fields. This filter was modified by the Hanyang University group by adding the no-slip wall boundary condition to the DFS filter developed by (Wang et al 2016).
The uncertainty of the velocity measurement (Table 2) was calculated for each temporal phase and velocity component using the uncertainty equation (Eq. 1) proposed by (Bruschewski et al 2016). For the entire temporal phase, the uncertainties did not change significantly and remained within 1.9% of VENC. However, they are slightly higher at the higher injection flow rates (9th to 14th phase) than the other. This phenomenon is probably caused by strong local turbulence and the results are indicative for each scan team and thus not repeated.
$${{\sigma }}_{\text{u}}=\frac{1}{\sqrt{2NSA}}\sqrt{Var\left\{{u}_{1}\left(ROI\right)-{u}_{2}\left(ROI\right)\right\}} (\text{E}\text{q}.1)$$
3.2 Stanford University/United States Military Academy (Stanford) Scan Details
Data were acquired using a 3T GE System as in Table 1 located at the Richard M. Lucas Center for MRI at Stanford University. The triggering system used a DAQ controlled by Labview (National Instruments, Austin, Texas, United States), and the mainstream flow rate was measured by a Transonic Systems TS410 flow module with an ME20PXL probe. The working fluid for the mainstream and the injector which was also triggered by the DAQ was a 0.06 Mol/Liter Copper Sulfate solution in de-gassed water.
Table 2
Temporal phase uncertainties based on HYU data.
The channel was placed in a standard single channel transmit receive head coil and scanned with sagittal slices with a gradient echo phase-contrast sequence triggered using the DAQ generated trigger connected to the MRI ECG system. The frequency direction was the primary streamwise direction, x, with a 25.6 cm field of view (FOV), and the phase FOV was 34% of this. Voxel size was 0.8 x 0.7 x 0.7 mm in the frequency, phase, and slice directions, respectively with matrix details as in Table 1. Phase-averaged data were acquired with a temporal resolution of 150 ms using a retrospectively gated sequence resulting in 10 phases and a scan time of 20 minutes and 20 seconds. These 10 phases were reconstructed in k-space to give 20 phases for the injection cycle. A total of 9 gated scans with flow through the channel and injector (flow-on) and 10 interleaved all flow-off ungated scans were acquired and combined to produce the final data set. The ungated flow-off scans were book-ended before and after a flow-on scan, and the averaged flow-off velocity field was subtracted from each phase of the flow-on scan to remove effects from system drift and eddy currents and produce a single phase-averaged set of velocity fields. The 9 sets of these fields were then averaged to produce the final data. A divergence-free filtering technique was leveraged in post-processing for the final data set (Schiavazzi et al 2014), which preserves local features while reducing random noise for smoother velocity gradients, helpful in calculating flow quantities where derivatives are used.
The uncertainty was estimated using the technique proposed by (Bruschewski et al 2016) and varies based on temporal phase as a result. The highest uncertainty was estimated to be less than 2% of the VENC and was nearly the same for the vertical and streamwise directions in the region around the jet for the maximum injection phases (9–12).
3.3 University of Rostock (Rostock) Scan Details
Data were acquired using a 3T Siemens Magnetom TRIO System (Siemens Healthineers, Erlangen, Germany) using two 3-channel body matrix coils. The working fluid was degassed, deionized water with 0.02 Mol/L Copper Sulfate. The whole flow system included 200 L of this fluid, and a laboratory chiller was used to keep the fluid temperature at 20.5°C. The mainstream flow was provided by a frequency-controlled 1100 W centrifugal pump and measured by a calibrated SM6020 flowmeter (IFM, Bochum, Germany).
The triggering system consisted of a data acquisition unit controlled by Labview (National Instruments, Austin, Texas, United States) using a clock rate of 10kHz with the digital signal sent every 1500 ms to trigger the scanning. With the same timing, a periodic analog voltage signal was sent to the trigger for the pulsatile pump, which was calibrated in steady-state conditions with a similar flow meter as used for the mainstream.
The minimum and maximum time between two received trigger signals was 1492 ms and 1508 ms. The acquisition window size describes the time available for scanning and was set to 1490 ms. This ensured scanning stopped before the subsequent trigger. Each cycle was equally divided into 20 bins with a size of 74.5 ms, which corresponds to the temporal resolution of the phase-locked acquisition. The sequence was designed such that velocity encoding was performed in the inner loop with a loop duration of 4*TR. The TR was set to 6.2 ms so that the time required for 3 k-space lines fits almost exactly into the bin size (4*3*TR = 74.4ms). In every 1500 ms cycle, three k-space lines were measured for all four velocity encodes and all 20 injection phases. Each k-space contains 100,080 k-space lines. Therefore, the whole cycle was repeated 3,360 times which resulted in a total acquisition time of 84 minutes for each flow-on scan. A total of five flow-on scans were performed.
Before and after these scans, a set of three flow-off scans were measured. These scans were not triggered, and the data were averaged and subtracted from each flow-on scan. Due to the lengthy measurement campaign, magnetic field drifts could become significant in the velocity data. It was found that the mean velocities in the flow-off measurements changed by 0.025 m/s between the first and last scans, which were separated by about eight hours. Assuming that this drift is linear, it can be estimated that the velocity error is approximately +/-0.011 m/s over the seven-hour flow-on scans. These deviations correspond to 1% of the VENC value and can therefore be considered negligible.
The stochastic measurement uncertainty in the flow-on scans were calculated with Eq. (1). The velocity component in the streamwise direction u showed the highest uncertainties with values between 0.015 m/s and 0.023 m/s depending on the injection phase. Relatively large uncertainty values occurred during ramping up and down of the jet. The velocity component v points in the direction of the jet and has an uncertainty between 0.014 m/s and 0.020 m/s. Like the stream-wise direction, there are variations in the uncertainty value that seem to correlate with the injection phases. The phase-averaged uncertainty in the third velocity component is 0.016 m/s and almost no phase dependent variations are visible in this component.
3.4 Seoul National University (SNU) Scan Details
The data set was acquired at Seoul National University Hospital using a phased array torso coil with pertinent details in Table 1. The working fluid was an aqueous 0.06 M copper sulfate solution, and the mean temperature was kept constant at 20.5°C using a chiller. The main flow was driven by a 250 W centrifugal pump (PW-350M, Wilo) and maintained at 22.0 L/min, while being monitored by a paddlewheel flowmeter. The pulsatile flow pump was operated by a function generator (AFG3152C, Tektronix) using a filtered square waveform with 1.5 second period as an input to the function generator. The flow rate of the pulsatile flow was monitored by a second paddlewheel flowmeter connected to the inlet tube of the pulsatile pump.
The number of scans used for the final data were 2 flow-on scans and 1 flow-off scan. Each scan measured a single field of view (FOV) of 256 × 38.4 × 38.4 mm3 in the x, y, and z directions at an isotropic voxel resolution of 0.8 mm3. A distortion correction filter was applied to prevent gradient warp due to the large FOV. A commercial 4D flow sequence was utilized. The sequence is comprised of 3D Cartesian spatial encoding and acquisition, and 4-point velocity encoding for three separate VENC values of 150, 150, and 80 cm/s in the x, y, and z directions. For measurement of the periodic flow, prospective gating was applied. The TTL signal for gating was sent from the function generator to the MRI scanner and synchronized with the waveform signal to the pulsatile pump. The acquisition window size describes the time available for scanning and was set to 1450 ms. Each cycle was equally divided into 22 bins with a size of 65.9 ms, which corresponds to the temporal resolution of the phase-locked acquisition. The sequence was designed such that velocity encoding was performed in the inner loop with a loop duration of 4*TR. The TR was set to 8.21 ms so that the time required for 2 k-space lines fits almost exactly into the bin size (4*2*TR = 65.68ms). Parallel imaging to reduce measurement time was not applied.
The k-space data were reconstructed using the open-source MATLAB code “mapVBVD” (developed by Philipp Ehses), which reads the Siemens raw data file. Each dataset was then converted to the respective velocity field. The final velocity data were obtained by averaging 2 flow-on scans and subtracting 1 flow-off scan to eliminate background phase error. A volumetric mask was constructed using the signal intensity data to distinguish between fluid and solid regions. The mask was superimposed on the velocity data, forcing spurious data in solid regions to be zero. A divergence-free smoothing (DFS) filter was applied to the masked velocity data to smooth out noise and reduce error in the velocity field. The DFS filter is based on (Wang et al 2016), and was modified by Hanyang University for wall-bounded flows with no-slip condition at the walls similar to the recommendations in (Benson et al 2020).
Uncertainty in the velocity measurements was evaluated using Eq. (1) by the approach from (Bruschewski et al 2016), where NSA is the number of signal averages and the factor 2 is due to the subtraction of two data sets. The spatial variance is calculated from the difference between two statistically independent velocity data sets comprising of all voxels within the region of interest (ROI). The uncertainty, which was estimated from the raw data, is 0.037, 0.032, and 0.016 m/s for the maximum pulsatile flow rate and 0.032, 0.027, and 0.014 m/s for the minimum pulsatile flow rate in the x, y, and z directions, respectively. It is estimated that the uncertainty is higher during jet injection because of the increased turbulence, which is a major source of noise. Nevertheless, the uncertainty values are within 2.5% of the VENC for all directions and phases.
3.5 University of Illinois – Urbana Champaign (UIUC) Scan Details
The MRV experiments were performed at the Beckman Institute at the University of Illinois, Urbana-Champaign, with the details in Table 1. The main flow was driven by a Little Giant 370W centrifugal pump that supplied 22.2 L/min when the secondary flow was running steadily at 1 L/min. A Blue-White (F1000RB) paddle flow meter calibrated against a Kobold MIM Electromagnetic flow meter was used to control the main flow rate accurately to \(\pm\)3.4% of the nominal test flow rate. An Omega FLR1605A differential flow meter upstream of the pulsatile pump, with a 0 to 2 L/min range and a full-scale accuracy of 2%, was used to pre-adjust the pump voltage inputs to ensure a steady flow rate of 0 and 1 L/min with the main flow at nominal conditions. The pulsatile pump was independently controlled by a programmable Keysight 33600A waveform generator that served as a synchronizer providing the damped sinusoidal waveform to the pump and a 5V TTL trigger signal to the scanner for prospective gating. There was no delay imposed between both signals. Primary and secondary flow loops fed from a 114 L tank in which the recirculating working fluid temperature was kept at 20.5°C with the help of an external PolyScience chiller. A vacuum degassed 0.06M aqueous Copper Sulfate solution was used as the working fluid.
Two flexible coils, a small 4-channel placed under the test-rig and an 18-channel array wrapped around the top and lateral walls, provided the best signal among the different coil configurations evaluated. The dataset was acquired using a Siemens sequence (fl_pc) with Cartesian sampling, and asymmetric phase encodings with flow compensation. With a data matrix of 384 by 144 by 48 in the x-y-z test section axes respectively, parallel imaging K-space domain based iPAT/GRAPPA with an acceleration factor of 2 in the y-direction allowed for reducing the scan time with negligible impact on the measurements. During each 1500 ms cycle, 31 equally spaced temporal phases were acquired for the same slice/phase encodings. This allowed us to maximize the temporal resolution, limited to the 47.35 ms time (total time for four-point flow encoding and three phase components). No filters or corrections were used during the acquisition or reconstruction.
A total of 3 gated scans with the flow-on, each bracketed by 2 non-gated scans with the flow-off before and 2 after for a total of 8 flow-off scans were performed. Non-gated flow-off scans were used to correct for possible eddy currents and system drift over time, none of which were identified over a total testing time of 6 hours. Given the absence of drifts, all 8 flow-off scans were averaged for background subtraction, limiting the impact of noise in the flow-off data.
Data postprocessing using in-house MATLAB codes consisting of averaging each temporal phase of the 3 flow-on scans following reconstruction, subtracting the averaged flow-off maps, and correcting for biases caused by non-linearities of the magnetic field gradients used for spatial encoding. High-order nonlinear components of the gradient coils’ magnetic fields increase in magnitude away from the isocenter, leading to geometric distortions and phase image biases in the reconstructed data when unaccounted for. We used the spherical harmonic coefficients for our scanner provided by Siemens to model the spatial dependence of the field generated by each gradient coil and used them to correct our data in a two-steps process. We first applied a 3D image distortion correction to map the phase-encoded information to the correct spatial location. The correction reproduces the unwrapping method incorporated in the Github HCP Pipelines distribution (GitHub). The second step corrects the impact of the gradient non-linearities on the first moments used to encode the motion, needed for the true magnitude and direction of the velocity data. This phase distortion correction follows the generalized reconstruction technique described by (Markl et al 2003). At last, the divergence-free filter discussed by (Schiavazzi et al 2014) was applied to each temporal phase to reduce random noise.
The uncertainty in the measurement of each velocity component was estimated from the postprocessed unfiltered data using the approach proposed by (Bruschewski et al 2016), under the assumptions that the noise is constant across the selected region for different realizations, and that its spatial distribution is Gaussian. The voxel-to-voxel difference between two independent velocity maps for the same temporal phase and different scans is used to compute the spatial variance. The temporal phase averaged estimated uncertainties are 2.2%, 2.1% and 2.2 % of the VENC values for each direction x, y, z, respectively.
3.6 Specific Scan Timing Details
As each team used different values and sequences for reporting echo (TE) and repetition times (TR), a short description of the general scan sequences is provided to aid investigators interested in duplicating the technique, and to help in overall understanding of the time windows for each since there are fluid dynamics and economic factors associated with overall scan time durations.
Figure 2(a) graphically displays some of the key details associated with the Hanyang University group, who used the 6-point velocity encoding scheme. In their sequence, the velocity is encoded separately for the three directions, by combining a reference excitation with no bipolar gradient applied with a second excitation using a bipolar gradient. This leads to an inner loop time of 2TR. The inner loop time describes the smallest non-divisible time scale in a sequence. Note that all other groups used the 4-point velocity encoding scheme as shown in Fig. 2(b), and all velocity information was encoded simultaneously with an inner loop time of 4TR. The inner loop is composed of four alternating bipolar gradients. Combining adjacent bipolar gradient pairs determines the velocity component in each coordinate direction. Therefore, each velocity component is again associated with a 2TR time window. The acquisition of all three velocity components is associated with a 4TR time window. Note that the acquisitions for different velocity components are not collocated in time. Therefore, the temporal resolution of four-point encoding methods is 4TR when considering all velocity components together.
Except for the velocity encoding method, the acquisition process depicted in Fig. 2 is representative of all groups. Specifically, within a single injection period, one or several k-space lines are acquired in each temporal phase. The injection cycle of 1500 ms must be repeated multiple times until all k-space data for all velocity encodings and all injection phases are sampled. Some teams used a relatively small TR so that multiple k-pace lines could be acquired during each encoding segment, which reduces the number of cycles necessary to sample all k-space data. This is the reason for the large deviations in the single scan times in Table 1.
In the case of the Hanyang University study, one k-space line was measured in each segment. To obtain a temporal resolution of 75 ms (20 phases) and given an inner loop time of 2TR (6-point encoding), the TR was adjusted to 37.5ms. The UIUC team also acquired one line of k-space in each segment but used a slightly different approach. With an inner loop time of 4TR (4-point encoding) and a temporal resolution of 47.35 ms (31 injection phases), the TR was adjusted to 11.8375 ms.
The teams at Stanford University, University of Rostock and Seoul National University leveraged the capability to adjust the number of k-space lines per segment, which means they acquired multiple lines of k-space in each encoding segment. The team at Stanford measured 10 injection phases (reconstructed to 20 injection phases) which resulted in a time resolution of 150ms. With an inner loop time of 4TR (4-point encoding) and a TR of 5.356 ms, they acquired 7 lines of k-space in each segment. As a result, their scan time was the shortest among all teams. The teams at the University of Rostock and the Seoul National University used a similar approach. They acquired 3 or 2 lines of k-space in each segment, respectively. However, when acquiring multiple k-space lines per segment, the temporal resolution is effectively reduced since the k-space data are spread out over time. For instance, the temporal resolution for one velocity component when acquiring one k-space line per segment is 2TR, while acquiring 3 k-space lines with 4-point velocity encoding reduces the resolution to 4TR*3 (or equivalently the length of a segment).
The research grade sequences used by these groups allow for this segmentation, whereas many commercial systems do not have access to this parameter. The cost of multi-segmenting is primarily one of blurring the temporal resolution slightly as the multiple acquisitions don’t happen at the precise same time.
Another aspect worth mentioning is the gating technique. With prospective gating, the segments are acquired at exactly the same points in time in the cycle, which reduces blur compared to retrospective gating. However, each trigger signal has a certain amount of uncertainty that must be considered when interpreting the expected gating. Usually, the acquisition window is made slightly shorter than the time between two trigger signals to account for these small deviations. For example, the teams here using prospective gating used acquisition windows between 1450 ms and 1490 ms (compared to the 1500 ms cycle). As a result, each phase is slightly shorter than the cycle divided by the number of temporal phases and there is a time shift that might accumulate over a cycle. These differences can lead to temporal discrepancies in the data.
With retrospective gating, the k-space lines are acquired with the prescribed timing dependent on the TRs and number of k-space lines in each segment. The gating signal is sampled simultaneously, and the lines of k-space are sorted into their respective phase intervals during reconstruction. If the TRs are not prescribed so that there is an integral number of k-space lines in each phase interval, then the lines start to span across the phase interval boundaries. This is particularly significant for the last phase which depending on the timing may be under sampled. In the case of this experiment, the Stanford TRs were set to match the phase intervals, so this is not a concern. Of course, there is some uncertainty in the timing of the gating signal which may cause some blurring across phases when k-space lines are misplaced during reconstruction. This effect, in general, is negligible relative to the uncertainties described in previous sections.