Dwipreproc -rpe_pair: Why equal number of up/down b0 encodings?

Hi all,

using dwipreproc RC1-44, also for topup.

Have dMRI full data in AP (e.g. 10x b0 & multishell bXXX data) and reverse phase b0 data (e.g. 5x b0 & no bXXX data).

Before I used something like:
dwipreproc -cuda -rpe_pair DTI_b0_AP.mif DTI_b0_PA.mif -force -nthreads $ncpu AP DTI_dn.mif DTI_dn_pp_with_APPA.mif

thus supplying AP and PA data for topup.

I know now we use a concatenated AP+PA.

BUT WHY, does the number of AP and PA b0 images need to be identical?

Now I get the error below.

Sure, I can make identical number of AP and PA b0 images, but than I trow away data…

Thx,

Stefan.

PS.

The error:

Treanus-MBP:DTI_mif stefan$ dwipreproc S_004_055.mif S_004_055_PP.mif -pe_dir AP -rpe_pair -se_epi S_004_055_all_b0.mif 
dwipreproc: 
dwipreproc: Note that this script makes use of commands / algorithms that have relevant articles for citation; INCLUDING FROM EXTERNAL SOFTWARE PACKAGES. Please consult the help page (-help option) for more information.
dwipreproc: 
dwipreproc: Generated temporary directory: /Users/stefan/Google_Cumulus/DICOM_DTI/DTI_mif/dwipreproc-tmp-0GPWW2/
Command:  mrconvert /Users/stefan/Google_Cumulus/DICOM_DTI/DTI_mif/S_004_055.mif /Users/stefan/Google_Cumulus/DICOM_DTI/DTI_mif/dwipreproc-tmp-0GPWW2/dwi.mif
Command:  mrconvert /Users/stefan/Google_Cumulus/DICOM_DTI/DTI_mif/S_004_055_all_b0.mif /Users/stefan/Google_Cumulus/DICOM_DTI/DTI_mif/dwipreproc-tmp-0GPWW2/topup_in.mif
dwipreproc: Changing to temporary directory (/Users/stefan/Google_Cumulus/DICOM_DTI/DTI_mif/dwipreproc-tmp-0GPWW2/)
dwipreproc: [ERROR] If using -rpe_pair option, image provided using -se_epi must contain an even number of volumes

Hi Stefan,

It’s a good question. There is indeed nothing preventing you from using all AP-PA images for field map estimation, and topup allows you to do this easily. The dwipreproc script is designed as a simple wrapper around topup and eddy to deal with the most common use cases. We never intended to replicate all its functionality, so if you require something more involved we recommend using the FSL tools directly.

That said, what I used to do with these scans was to average all AP b0’s and all PA b0’s and then use these averages in the dwipreproc script. I don’t expect the topup field estimate to differ much either way. However, in cases with large head motion it might be better to select (individual) AP-PA images acquired closely in time. Preferably, these are also acquired at the start of the scan (otherwise I think there’s an eddy parameter to specify the reference).

Cheers,
Daan

Thx @dchristiaens,

doing an average of all AP and PA b0 now.
But keeping in mind your consideration about large bulk head motion…

Groetjes,

S.

Hi Stefan,

I might initially go on a bit of a tangent on this one; I promise it’ll come around to answer the question, but I’m hoping it’ll give anyone curious a better understanding of how I made this all work.

With 3.0_RC1, dwipreproc has been basically re-written from scratch, largely due to the requirements of being able to handle any arbitrary input data without user interaction for the BIDS apps (now published). In such a use case, we can no longer rely on a user informing the script of the phase encoding design of the acquisition, therefore such information must be determined from the data.

In 3,0_RC1, in addition to a diffusion weighting table, it is also possible for image volumes to have stored in the header a phase encoding table. This captures the direction of phase encoding, and the total readout time, for each volume, just like the topup “configuration file”. When importing from DICOM, MRtrix3 will now read the appropriate entries, generate this table, and store it in the image header.

(Note that while BIDS defines fields for storing phase encoding direction and readout time, it doesn’t have the capacity for storing a 4D image where the phase encoding direction / readout time varies between volumes)

Now, when you run dwipreproc with the new -rpe_header option, it examines the phase encoding table in the image header, and sets up e.g. the topup / eddy calls appropriately based on this information. This occurs quite differently to how the dwipreproc script used to work, with the -rpe_none / -rpe_pair / -rpe_all options passing through substantially different code branches. But I couldn’t remove these old command-line options from the new version of the script, since it can’t be guaranteed that the phase-encoding table will be available in all cases (and also people have become accustomed to using them).

The way I resolved this in the newest version of the script is: When the user uses one of the -rpe_none / -rpe_pair / -rpe_all options in conjunction with the -pe_dir option (and possibly also the -readout_time option), the script uses this information to construct the phase-encoding table based on the user’s inputs and the number of image volumes. This then means that once the phase-encoding table is constructed, these three older -rpe_* “operating modes” of the script in fact pass through the very same code path as the -rpe_header option.

However, this is only possible if the script has enough information to construct that phase-encoding table…


OK, promised I’d get there in the end. The reason -rpe_pair requires an even number of volumes in the -se_epi image is: If the number of volumes were odd, it would be impossible for the script to unambiguously determine the phase-encoding table for those volumes. E.g. for an image with 5 volumes, is it two A>>P and three P>>A images, or three A>>P and two P>>A images?

The alternative - assuming you have a Siemens scanner and don’t use any non-MRtrix3 commands between DICOM import and dwipreproc invocation - is to use the -rpe_header option (this will happily accept any -se_epi image, as long as there is phase-encoding contrast with which to estimate the inhomogeneity field. Note that the output of dwipreproc with the -rpe_header option should always be visually checked, since I can’t guarantee that the phase-encoding import is bulletproof; but I’m hoping that over time more people will make use of this feature.

Rob

Hi @rsmith, Rob,

I’ll look into the -rpe_header option.

At present working with GE and Philips data.

Concerning: “assuming you have a Siemens scanner and don’t use any non-MRtrix3 commands between DICOM import and dwipreproc invocation”, and being lazy, and not (yet) reading the source code, how do you determine the phase-encoding contrast from the (dicom) header? Should be possible on GE/PMS too, I guess?

Thx,

S

Currently it reads 0018 1312 and a couple of Siemens CSA entries. Should absolutely be possible to import from other vendors, but I only have access to Siemens data and so have not had the information necessary to configure the import for data from other vendors. If anybody can get me data on the appropriate fields, I’d be happy to update the import code.

Dear Rob
I have a data set from Philips which I think (and would appreciate if you can also confirm) has different readout times for rev (3 vols) and dwis (102 vols).
So I was wondering whether it is possible to pass different readout times to dwipreproc using -read_time option?
Please let me know if you would like to check this data set.
Thanks and cheers,
Hamed

Hamed,

No, you can’t manually specify different readout times for different volumes. dwipreproc's only purpose for existence is to simplify the interaction with FSL’s topup and eddy tools in the most common use cases; if the command-line interface were sufficiently complex to handle any arbitrary acquisition scheme, it would be no more simple than using the FSL tools themselves directly, and it would therefore cease to serve a purpose.

One benefit of the new -rpe_header option is that dwipreproc can now handle such more complex acquisitions without necessitating a more complex user interface. But there’s a good chance that it won’t work with Philips data without an update to the DICOM import code, which requires either that a developer with Philips data make the necessary changes, or somebody provide me with example data & protocol so that I can locate whatever mechanisms Philips use to encode this information.

Rob

A post was split to a new topic: Dwipreproc: Averaging of b=0 volumes