wishlist: TotalReadoutTime for sidecar BIDS .json

[yarikoptic]

BIDS 1.0.1 added
TotalReadoutTime: defined as the time (in seconds) from the center of the first echo to the center of the last echo (aka “FSL definition” - see here and here how to calculate it). This parameter is required if a corresponding multiple phase encoding directions fieldmap (see 8.9.4) data is present.

• https://fsl.fmrib.ox.ac.uk/fsl/fslwiki/topup/Faq#How_do_I_know_what_phase-encode_vectors_to_put_into_my_--datain_text_file.3F
• http://lcni.uoregon.edu/kb-articles/kb-0003

I wondered how life would be better if sidecar file had it automagically computed for me and BIDS validator shut itself up ;)

[neurolabusc]

Can you provide a sample modern DICOM and a recipe for how you compute this using only the DICOM and CSA headers? I think the links you provide are obsolete (as well as this) as modern Siemens scans no longer seem to create the private tag 0019,1028. I note that dicm2nii attempts to compute this value by reading the dreaded ASCII portion of the CSA header, which I would really rather not do (for both performance and reliability reasons). In addition, with my own recent images, the values generated by dicm2nii do not look plausible. I now use TOPUP, so have not done the FUGUE undistortion for years. However, my recollection was that this value was used simply to save the unwrapped fieldmap estimate in calibrated units, and had no influence on the actual deformations estimated or applied.

I do not know of a reliable way to extract this information from modern Siemens DICOM images. Looking at the available tools I think the method for determining this has changed with different software versions, so even if we have a solution that works for one software version (D13), I am not convinced it would work for prior (B17) or more recent (E11) tools. Therefore, I think one would need to go through an archive of old datasets to work out how the recipe has changed. Indeed, since we have to rely on private tags, any software we create may prove brittle and generate bogus values for future releases.

The situation appears simpler for Philips, though I do not have example datasets.

In brief, this is not a standard DICOM property, and the private vendor specific methods for inferring this value appear to have changed over time. Therefore, computation of this value does not seem to be a task for the DICOM to NIfTI stage (as it is not part of the DICOM standard). @chrisfilo (BIDS maintainer for dcm2niix) may want to comment, but I assume the intention is that this value is provided by the physicist rather than extracted from the DICOM data. Further, if my hunch is correct that this only changes the absolute values (which are scaled based on the data) such that it does not impact the unwarping, I would presume this is an optional field.

[chrisgorgo]

If this value can be automatically extracted by the dcm2niix converter this would be much more convenient and (potentially) less error prone that getting it from the local MR physicist. This field is not mandatory in BIDS, but since it's necessary if you are planning to use FSL TOPUP to do susceptibility distortion correction we highly recommend including it.

…

On Wed, May 3, 2017 at 1:01 AM, Chris Rorden ***@***.***> wrote: Can you provide a sample modern DICOM and a recipe for how you compute this using only the DICOM and CSA headers? I think the links you provide are obsolete (as well as this <https://groups.google.com/forum/#!topic/7t_trt/fDR35kgU4mU>) as modern Siemens scans no longer seem to create the private tag 0019,1028. I note that dicm2nii attempts to compute this value by reading the dreaded ASCII portion of the CSA header, which I would really rather not do (for both performance and reliability reasons). In addition, with my own recent images, the values generated by dicm2nii do not look plausible. I now use TOPUP, so have not done the FUGUE undistortion for years. However, my recollection was that this value was used simply to save the unwrapped fieldmap estimate in calibrated units, and had no influence on the actual deformations estimated or applied. I do not know of a reliable way to extract this information from modern Siemens DICOM images. Looking at the available tools I think the method for determining this has changed with different software versions, so even if we have a solution that works for one software version (D13), I am not convinced it would work for prior (B17) or more recent (E11) tools. Therefore, I think one would need to go through an archive of old datasets to work out how the recipe has changed. Indeed, since we have to rely on private tags, any software we create may prove brittle and generate bogus values for future releases. The situation appears simpler for Philips <https://www.jiscmail.ac.uk/cgi-bin/webadmin?A2=fsl;162ab1a3.1308>, though I do not have example datasets. In brief, this is not a standard DICOM property, and the private vendor specific methods for inferring this value appear to have changed over time. Therefore, computation of this value does not seem to be a task for the DICOM to NIfTI stage (as it is not part of the DICOM standard). @chrisfilo <https://github.com/chrisfilo> (BIDS maintainer for dcm2niix) may want to comment, but I assume the intention is that this value is provided by the physicist rather than extracted from the DICOM data. Further, if my hunch is correct that this only changes the absolute values (which are scaled based on the data) such that it does not impact the unwarping, I would presume this is an optional field. — You are receiving this because you were mentioned. Reply to this email directly, view it on GitHub <#98 (comment)>, or mute the thread <https://github.com/notifications/unsubscribe-auth/AAOkp4dKmUvp18YjZec5aFv6OxIWfXgbks5r2DRqgaJpZM4NOkcc> .

[neurolabusc]

@chrisfilo
1.) looking at recent Siemens Prisma data (D13), I do think the previous methods for extracting these parameters are now broken.
2.) I note that the modern CSA header lists "RealDwellTime" and "SliceMeasurementDuration", which might allow us to estimate this, though we need to see how these are influenced by iPAT and partial Fourier. I think someone needs to sit down at a scanner console and vary a few parameters (bandwidth, lines, iPAT, Fourier) to validate these. The bigger concern I outlined is how to get a general solution that works across several generations of Siemens scanners, as I assume dicm2nii's solution used to work once upon a time. I really do not want to include features into dcm2niix which are not robust.
3.) In typical situations where one applies TOPUP, I think the readout time is identical for all acquisitions and therefore you don't neccessarily have to specify a valid value in this column.

[chrisgorgo]

2. This might help: https://lcni.uoregon.edu/kb-articles/kb-0003
3. Correct. It only really matters if you using files with different total readout time values

[neurolabusc]

As I noted earlier, the information from lcni is obsolete, as recent Siemens DICOMs do not generate the tag 0019, 1028.

[dangom]

Isn't the TotalReadoutTime in EPI equal to the EchoTrainLength (0018,0091)? This parameter should already take PF and Parallel Imaging into account.

[neurolabusc]

No, ETL is "Number of lines in k-space acquired per excitation per image", and does not tell you how long it took to acquire those lines. While it is part of the equation, it is not the whole one, and we do need to validate that partial fourier and in plane parallel imaging (SENSE/GRAPPA) are taken into account. I do think someone needs to rigorously investigate this with sample datasets, and hopefully find a robust solution that can be used across vendors and across software updates.

n.b. For future reference, in practice Siemens EPI scans always record the ETL as 1, so you should look at EPI factor. http://mriquestions.com/echo-planar-imaging.html

[dangom]

You are right. The TotalReadoutTime would be the EchoSpacing * ETL. Interestingly enough, the EchoSpacing is given on the parameter card when setting up the sequence, but isn't written to the DICOM.

[neurolabusc]

@dangom - for completeness the FSL (and therefore BIDS) definition is TotalReadoutTIme = EchoSpacing * (ETL-1) (fencepost problem).

[jonclayden]

I'm also trying to tackle this issue at the moment, since I'm building dcm2niix into our local pipeline and we use TOPUP routinely in research. The ASCII part of the DICOM files from our Prisma (D13) contains the lines

sFastImaging.lEPIFactor	 = 	128
sFastImaging.lEchoSpacing	 = 	570

(with the latter in microseconds — thanks @jmtyszka for the correction!), and this is what I have relied on in the past. I don't know how stable this has been/will be across software versions, though.

[jmtyszka]

@jonclayden - these fields are from the Series Shadow Data within the Siemens CSA header of the DICOM file. There's a documentation attempt here: http://gdcm.sourceforge.net/wiki/index.php/Gdcmdump#SIEMENS_CSA_Header, but the lEchoSpacing field doesn't appear in the GDCM dictionary and could be a D13 addition. If you're using product EPI sequences, the EPIFactor and EchoSpacing fields might be stable during that software version (D13), but there's no guarantee they won't change in the future.

[Minor point: lEchoSpacing is in microseconds, not milliseconds]

[neurolabusc]

@jmtyszka, Yyou can get some archival Siemens EPI data here.
http://www.cabiatl.com/mricro/mricron/dcm2nii.html
This does confirm what I mentioned earlier - there is probably a way to get this fairly reliably by looking at the ASCII portion of the CSA header. I guess we can play with adding support for this, but I really hate delving into this part of the header.

In any case, surveying those archival files it does look like we can find something there... the B13 image was acquired in 2007.

B12 (ep2d_bold) field does not exist
B13 (ep2d_diff) sFastImaging.lEchoSpacing = 800
B15 (ep2d_diff) sFastImaging.lEchoTrainDuration = 700
B17 (ep2d_bold) sFastImaging.lEchoTrainDuration = 700
B17 (ep2d_diff) sFastImaging.lEchoTrainDuration = 700

[jmtyszka]

Thanks Chris! Very cool that you still have a stash of older VB DICOMs and I agree, the CSA header can be a bit like the shifting sands of the Sahara. I have to look a bit deeper at how these fields are populated by a Siemens sequence, particularly whether it's consistent between product and 3rd party sequences like the CMRR multiband package and some of the Siemens WIP releases. Will update when I have a chance.

[neurolabusc]

I have just uploaded some new code that attempts to extract data from the ASCII portion of the CSA header. Unlike @jonclayden, neither of the two D13 systems I have data from provide EchoSpacing in this field (perhaps sequence differences - we use the CMRR multiband sequences). If successful, the software populates the BIDS with the following values
"EchoSpacing": 800,
"EchoTrainDuration": 700,
"EPIFactor": 128,
However, for my sample data I only get all 3 for B13, whereas B15, B17 and D13 omit the "EchoSpacing".

Note that since scanning ASCII data is slow the new code is only invoked for Siemens data, with a CSA Series Header and where the ScanningSequence is "EP" (EPI), and only once for each series.

@jmtyszka - I do have a deeper concern with peaking into the CSA header. At one stage (hopefully now fixed) Siemens did not zero-pad the binary portions of the CSA header, and this meant that random data from previous images could appear in the CSA portions. This makes scanning the CSA for ASCII treacherous, as your key words can get randomly intermixed with real data. I will try to refine my code later to skip the binary portions, but it does lead to ugly and hard to maintain code.

[jonclayden]

Thanks, Chris! This does indeed pull out all three values from my test data. I'll poke around with some other data when I get a chance.

[neurolabusc]

I have just added a small update to restore Windows compatibility and to restrict the key searching to the ASCII portion of the series header (to avoid the known issue that the uninitialized bytes in the BINARY portion often include fragments from prior ASCII data).

Unfortunately, this work confirms my earlier suspicion that the needed information is often not stored in the DICOM images, and one needs to sit down on the console to determine TotalReadoutTime. Unless someone else has insight, I will close this issue.

[neurolabusc]

I think the latest release (v1.0.20170609) solves this. I would appreciate if all of you (@yarikoptic, @chrisfilo, @jonclayden, @jmtyszka and @dangom) could test out the latest release.

I got some help from Edward J. Auerbach (Center for Magnetic Resonance Research), though I take credit for any errors. In particular, for the "FSL" definition of TotalReadoutDuration we need to deal with the fence post issue. So if you have 90 reconstructed lines with AccelFactPE = 2, we need to know the time when we started to acquire the 88th line (as the 89th was interpolated).

EffectiveEchoSpacing = 1.0 / (BandwidthPerPixelPhaseEncode * PhaseEncodingLines)

TotalReadoutDuration = EffectiveEchoSpacing * (NumberOfPhaseEncodingSteps - AccelFactPE)

TrueEchoSpacing = EffectiveEchoSpacing * AccelFactPE

CSAImageHeader has the tag named "BandwidthPerPixelPhaseEncode"
(0051,1011) AccelFactPE (1 if not present: no in slice acceleration)
(0018,0089) NumberOfPhaseEncodingSteps accounts for iPAT and partial Fourier.
(0018,1310) PhaseEncodingLines (large of 3rd and 4th item, depends on PE direction

I no longer need to read data from the dreaded ASCII portion of the CSASeriesHeader, so that troublesome code is no longer used. This means the superfluous JSON tags "EchoSpacing", "EchoTrainDuration", and "EPIFactor" are not created. If anyone does want these tags, you can compile the software with the -DmyReadAsciiCsa compiler directive

[chrisgorgo]

Thanks for putting so much effort into getting this right. It's very much appreciated.

One things that I am concern with is the change of the name from TotalReadoutTime to TotalReadoutDuration which makes it incompatible with current version of BIDS. What was the motivation behind this change?

[neurolabusc]

@chrisfilo - I will change this to the BIDS name once I get feedback that this is the correct solution across sequences, versions, etc. At the moment, it LOOKS like the right solution to me, but I have only tested it on a tiny subset of images. Once the code gets validated, we can change the name and distribute via NITRC (as part of the MRIcroGL package).

[chrisgorgo]

Excellent - thanks!

…

On Sun, Jun 11, 2017 at 10:31 AM, Chris Rorden ***@***.***> wrote: @chrisfilo <https://github.com/chrisfilo> - I will change this to the BIDS name once I get feedback that this is the correct solution across sequences, versions, etc. At the moment, it LOOKS like the right solution to me, but I have only tested it on a tiny subset of images. Once the code gets validated, we can change the name and distribute via NITRC (as part of the MRIcroGL package). — You are receiving this because you were mentioned. Reply to this email directly, view it on GitHub <#98 (comment)>, or mute the thread <https://github.com/notifications/unsubscribe-auth/AAOkpy5gHnZWivbVzf1h5w7-sd_mpx2sks5sDCRkgaJpZM4NOkcc> .

[yarikoptic]

fwiw, I intend to test on some test sequences we collected for this one, but didn't have a chance yet and the local physicist, who could verify, just became a father again... So might take a few more days

[dangom]

I tested with three in-plane accelerated datasets (from a Siemens VD13D Prisma) with the following specs:

Echo Spacing 0.58ms rBW 2204 Hz/px
Echo Spacing 0.58ms rBW 2380 Hz/px
Echo Spacing 0.59ms (free Echo Spacing) rBW 2380 Hz/px

All tests passed.

[neurolabusc]

The latest release now creates TotalReadoutTime . I acquired 6 series of images on our Prisma D13: (a) no slice acceleration, (b) 6/8 partial fourier and (c) iPAT = 2 for both the Siemens Product as well as the CMRR multi-band sequences. The values derived from the bandwidth values stored in the DICOM match the traditional method based on using the echo spacing value stored in the protocol PDF. I am including both the images and the PDF in a new github project dcm_qa which is a very simple dcm2niix validator. The validator also includes B13 data that demonstrate slice timing and the DICOM-specific lossless JPEG transfer syntax. This repository should make it easier to spot regressions in future releases.

[chrisgorgo]

The validator is an excellent idea!

…

On Wed, Jun 21, 2017 at 9:45 AM, Chris Rorden ***@***.***> wrote: Closed #98 <#98>. — You are receiving this because you were mentioned. Reply to this email directly, view it on GitHub <#98 (comment)>, or mute the thread <https://github.com/notifications/unsubscribe-auth/AAOkpyM0o9XHVtdSu1cV7xypIUVMcxpsks5sGR5pgaJpZM4NOkcc> .

[mharms]

See issue #130 for relevant updates to this issue, including nuances regarding the calculation of EffectiveEchoSpacing and TotalReadoutTime when using Phase Resolution < 100% and/or Interpolation.