search

rordenlab / dcm2niix

dcm2nii DICOM to NIfTI converter: compiled versions available from NITRC
https://www.nitrc.org/plugins/mwiki/index.php/dcm2nii:MainPage
Other
900 stars 229 forks source link

Suggestions for GE dicom header info to use in BIDS sidecar #163

Closed mick-d closed 6 years ago

mick-d commented 6 years ago

Hi,

I needed to use the BIDS sidecar created by dcm2niix for one project, and unfortunately some information about DICOM data from a GE Signa (HDxt) scanner was missing. I am opening an issue here as recommended by Chris Gorgolewski on this NeuroStars post to provide information which might be useful for the missing information I could find and the missing information I could not solve. Please find the original post below for convenience.

I could not get any of the following fields (among others) with dcm2niix (including latest version I compiled myself at the time of writing):

I could manually get EffectiveEchoSpacing by looking at DICOM field 0043 102c (cf GE doc) and SliceTiming by looking at DICOM field “Trigger Time”.

I could compute TotalReadoutTime with (<DICOM field “IMG Columns”> - 1) * (1 / EffectiveEchoSpacing).

However I could not find a way to get BandwidthPerPixelPhaseEncode (the DICOM field “ACQ Pixel Bandwidth” seems to be something different) and could only get part of the PhaseEncodingDirection: the direction but not the sign (direction was “j” as DICOM field “ACQ Phase Encoding Direction” was ROW).

neurolabusc commented 6 years ago

I do not have access to a GE hardware. We would need a set of validation images (ideally ones we could share publicly, so it would help other conversion tools). As an example @mharms provided publicly available Siemens data with careful attention to deriving effective echo spacing. Note his sample images systematically vary PhaseResolution, PhaseOversampling, PercentPhaseFOV, EchoTrainLength, PhaseEncodingSteps, AcquisitionMatrixPE, ReconMatrixPE, and BandwidthPerPixelPhaseEncoding.

When I wrote the BIDS export for GE data, I tried to include all the values I could reliably determine. However, in my sample datasets I was not able to determine the BIDS tags you are requesting. At the time, Paul Morgan provided my with several GE datasets designed to elucidate these values, but I was unable to confidently work these out.

One other suggestion - you may want to try Xiangrui Li's dicm2nii. I think he devoted more time than me in attempting to decode GE datasets. You may find that it provides a richer set of BIDS data than dcm2niix. If this is the case, you may want to use dicm2nii instead of dcm2niix, or alternatively since dicm2nii and dcm2niix are both open source, you can decipher the method used by dicm2nii and provide a pull request for dcm2niix. This would help all GE users who rely on dcm2niix. You should find dicm2nii and dcm2niix have relatively similar functions, as we worked together to developer our tools.

mick-d commented 6 years ago

Many thanks for the feedback @neurolabusc. I did not know about dicm2nii. I will look into it and get back to you.

neurolabusc commented 6 years ago

Please test the latest pre-release. The notes describe the features extracted from the GE protocol data block as well as the limitations. @mick-d please test with your data and report any issues.

mick-d commented 6 years ago

@neurolabusc Will do and report here. Thank you!

leej3 commented 6 years ago

I did a little more testing with this. And sat down with one of our physicists to try figure out some of the inconsistencies I was observing. Here's a summary of what we figured out:

Many EPI sequences collected at out site use an interface at the scanner console that writes to the GE proprietary "Protocol Data Block". The code code works for these nicely. This includes a stock GE epi sequence. I've linked to scans collected on a phantom in both directions below.

It seems that the phase encoding polarity cannot be extracted for multiphase EPI sequence variant used for the ADNI protocol (http://adni.loni.usc.edu) though. The sequence does not use the same console interface to specify the image polarity and so VIEWORDER is not written to the data block. Currently the code seems to "find" a VIEWORDER. I'm not sure why. If the information is buried in there somewhere else that would be great.

It appears product DTI (including ADNI scans) leaves out SLICEORDER and incorrectly specifies VIEWORDER. Again dcm2niix specifies j- for these scans. Perhaps NA or something would be more appropriate? I've again linked to scans collected in both directions below.

@roopchansinghv can add some adni epis in both directions if that would help

Thanks

John

ADNI DTIs , both directions: https://drive.google.com/file/d/1wLCgISdE4guL_ROb4U_Lb77vVvSo7GMb/view?usp=sharing

Stock GE EPI, both directions: https://drive.google.com/file/d/1ULITxqJ5SPn3OMKO1oV0daHaY4DEah4r/view?usp=sharing

neurolabusc commented 6 years ago

I agree, the Protocol Data Block in the sample ADNI images reports VIEWORDER "1" regardless of phase encoding direction. Is your proposed solution that we ignore VIEWORDER (and leave phase encoding direction unknown in the BIDS header) if the PHASEACCEL tag is present and greater than 1? In the example ADNI dataset PHASEACCEL "2.00". Since I do not have GE hardware available, I appreciate any advice or examples that test the robustness of dcm2niix.

By the way, since the Protocol Data Block (0025,101B) is gzip compressed, you can not view it with most DICOM header viewing utilities (e.g. dcmdump, Horos, etc). If you run dcm2niix in extreme verbosity mode -v 2, it will display the text fields of the Protocol Data Block.

leej3 commented 6 years ago

Thanks. We had indeed been using -v 2. It has been very helpful. We hadn’t realized PHASEACCEL was a part of this too.

Maybe rather than stating unknown it would be better to explicitly state there is no standardized way of fully specifying phase encoding direction from the Dicom header. The more people that know about it the better. People can changes series description-names/use different sequences etc.

I’ll be away from the computer for a few days but after that I can help map out the various situations in which this occurs. Thanks for the help.

On Mar 31, 2018, at 6:56 AM, Chris Rorden notifications@github.com wrote:

I agree, the Protocol Data Block in the sample ADNI images reports VIEWORDER "1" regardless of phase encoding direction. Is your proposed solution that we ignore VIEWORDER (and leave phase encoding direction unknown in the BIDS header) if the PHASEACCEL tag is present and greater than 1? In the example ADNI dataset PHASEACCEL "2.00". Since I do not have GE hardware available, I appreciate any advice or examples that test the robustness of dcm2niix.

By the way, since the Protocol Data Block (0025,101B) is gzip compressed, you can not view it with most DICOM header viewing utilities (e.g. dcmdump, Horos, etc). If you run dcm2niix in extreme verbosity mode -v 2, it will display the text fields of the Protocol Data Block.

— You are receiving this because you commented. Reply to this email directly, view it on GitHub, or mute the thread.

neurolabusc commented 6 years ago

Since I suspect this will take a few weeks to sort out, I think we should push ahead with a new release. I would suggest that for the imminent release for GE images we insert the BIDS tags "ProbableSliceTiming" and "ProbablePhaseEncodingDirection", since we can often determine these but in a few corner cases our current methods fail.

@chrisfilo - do you think this is a good solution for the next stable release, or is it best to omit these tags entirely if we are not yet able to predict these values with 100% accuracy? I will use whichever approach you use. The question is whether to include tags that are usually correct and useful, but where we know there are some instances where they fail (e.g. phase direction for ADNI sample) or suspect they will fail (slice timing for multi-band).

By the way, I know that dicm2nii uses the ProtocolDataBlock SLICEORDER only when the timing can not be determined using the RTIA_Timer (0021,105E). However, none of the sample I have seen provide information in 0021,105E, so I am vary to add support for a method where I have no sample datasets.

Finally, I am set the label for this as an enhancement. These are wish list features rather than known bugs. I am eager to have a new stable release as the handling of enhanced Philips DICOM has improved a lot since December.

chrisgorgo commented 6 years ago

Thanks for looking into this!

Ideally I would keep the real field names you are fairly confident values are correct and print a warning (on stderr) when a known edge case is discovered.

neurolabusc commented 6 years ago

@chrisfilo - thanks for your rapid comments. The issue is that I know we get the wrong answer with the ADNI example, and I do not know how to detect it (adjusting PHASEACCEL is a testable hypothesis). I also know that computing slicetiming based on SLICEORDER will give the wrong value if the acquisition is sparse or hyperband (e.g. SMS/multiband). We get SLICEORDER from decoding a proprietary DICOM element, and I have very few examples of GE datasets. Therefore, the issue is how we identify BIDS tags that we suspect might be unreliable, or if we only include tags that we are confident of.

mharms commented 6 years ago

Personally, I think it is dangerous to start including tags as part of the DICOM conversion that have suspect reliability. At a minimum, I think the tag should include “Possible” or “Probable”, as Chris R. suggested, although I think even that is risky, because users may just ignore the caution implied by the “Possible/Probable” tag.

-- Michael Harms, Ph.D.

Associate Professor of Psychiatry Washington University School of Medicine Department of Psychiatry, Box 8134 660 South Euclid Ave. Tel: 314-747-6173 St. Louis, MO 63110 Email: mharms@wustl.edu

From: Chris Rorden notifications@github.com Reply-To: rordenlab/dcm2niix reply@reply.github.com Date: Monday, April 2, 2018 at 6:03 PM To: rordenlab/dcm2niix dcm2niix@noreply.github.com Cc: "Harms, Michael" mharms@wustl.edu, Mention mention@noreply.github.com Subject: Re: [rordenlab/dcm2niix] Suggestions for GE dicom header info to use in BIDS sidecar (#163)

@chrisfilohttps://github.com/chrisfilo - thanks for your rapid comments. The issue is that I know we get the wrong answer with the ADNI example, and I do not know how to detect it (adjusting PHASEACCEL is a testable hypothesis). I also know that computing slicetiming based on SLICEORDER will give the wrong value if the acquisition is sparse or hyperband (e.g. SMS/multiband). We get SLICEORDER from decoding a proprietary DICOM element, and I have very few examples of GE datasets. Therefore, the issue is how we identify BIDS tags that we suspect might be unreliable, or if we only include tags that we are confident of.

— You are receiving this because you were mentioned. Reply to this email directly, view it on GitHubhttps://github.com/rordenlab/dcm2niix/issues/163#issuecomment-378073316, or mute the threadhttps://github.com/notifications/unsubscribe-auth/ADswBgGXmgC4mOnw6PtHxyl3XjFXKsx0ks5tkq4sgaJpZM4SKb17.

The materials in this message are private and may contain Protected Healthcare Information or other information of a sensitive nature. If you are not the intended recipient, be advised that any unauthorized use, disclosure, copying or the taking of any action in reliance on the contents of this information is strictly prohibited. If you have received this email in error, please immediately notify the sender via telephone or return mail.

neurolabusc commented 6 years ago

Thanks @mharms for your input. If @chrisfilo agrees the next major release will comment out the line #define myReadGeProtocolBlock, and therefore will not attempt to decode GE's proprietary Protocol Block. Advanced users and developers will still be able to compile their own version. From what I have seen, GE changes the included fields and usage between versions. Indeed, the ADNI sample provided by @leej3 suggests that for some variations of GE the phase direction is simply not recorded in the header. I spent a lot of time testing this out, so psychologically I am in an escalation of commitment spiral. However, my sense is we tell GE users there is currently no reliable way to extract slice timing and phase encoding and perhaps they can lobby the vendor for proper support.

chrisgorgo commented 6 years ago

Sounds very reasonable to me.

On Tue, Apr 3, 2018, 6:29 AM Chris Rorden notifications@github.com wrote:

Thanks @mharms https://github.com/mharms for your input. If @chrisfilo https://github.com/chrisfilo agrees the next major release will comment out the line #define myReadGeProtocolBlock, and therefore will not attempt to decode GE's proprietary Protocol Block. Advanced users and developers will still be able to compile their own version. From what I have seen, GE changes the included fields and usage between versions. Indeed, the ADNI sample provided by @leej3 https://github.com/leej3 suggests that for some variations of GE the phase direction is simply not recorded in the header. I spent a lot of time testing this out, so psychologically I am in an escalation of commitment spiral. However, my sense is we tell GE users there is currently no reliable way to extract slice timing and phase encoding and perhaps they can lobby the vendor for proper support.

— You are receiving this because you were mentioned. Reply to this email directly, view it on GitHub https://github.com/rordenlab/dcm2niix/issues/163#issuecomment-378250120, or mute the thread https://github.com/notifications/unsubscribe-auth/AAOkp6iDxVlAzlA0W4WuLFM6t08mK8Ubks5tk3krgaJpZM4SKb17 .

mick-d commented 6 years ago

@neurolabusc Apologies for the delay, the data was not mine and it took time to get it back for testing.

The results are that “EffectiveEchoSpacing”, “TotalReadoutTime”, “SliceTiming” and “PhaseEncodingDirection” are now included, although "PhaseEncodingDirection" includes direction ("i") but no sign ("+/-"). There is not "BandwidthPerPixelPhaseEncode" but there is "PixelBandwidth". I can also confirm that the first 4 fields are the same as the Matlab tool dicm2nii, while neither "BandwidthPerPixelPhaseEncode" or "PixelBandwidth" is output by that tool.

Thanks a lot.

neurolabusc commented 6 years ago

Great. Note you can recompile with the compiler flag -DmyReadGeProtocolBlock and it will report PhaseEncodingDirection with the polarity (e.g. j+ or j-). Just be aware that we know of explicit cases where this value is wrong (the ADNI file linked above). @xiangruili (of dicm2nii fame), Paul Morgan and I have looked exhaustively at the GE headers and we have concluded there is no simple and reliable method to determine these features. In theory, a system that new the GE software version and a few other features might be able to report some of these elements some of the time. However, this would require some knowledge of when these features are used.

Therefore, I am closing this issue. It seems like the current release of dcm2niix reports everything we are able to reliably extract from GE data. I strongly encourage you to contact GE and ask them to explicitly export and document the parameters that are important to you in private fields of the DICOM image (rather than in their proprietary protocol block).

mharms commented 6 years ago

If I'm reading this thread correctly, the current release (without any recompiling) will report a PhaseEncodingDirection for GE data, but the polarity will always be positive? If so, I think that is dangerous, because users will have no indication that polarity is NOT included. If "direction" is to be reported, but without polarity, I think that needs its own specific tag, so as to avoid any mixing up with the existing PhaseEncodingDirection.

xiangruili commented 6 years ago

I agree it is bad to include possibly wrong parameter. But the problem is that GE includes wrong parameter in dicom, so they should fix it eventually.

If we know certain flags that indicate the possibly wrong parameter, the converters can omit the PE polarity and leave it to users. But this kind of flags may not exist.

Xiangrui Li, Ph.D. Center for Cognitive and Behavioral Brain Imaging The Ohio State University Phone: 614-292-1847

From: Michael Harms Date: 2018-04-10 16:54 To: rordenlab/dcm2niix CC: xiangruili; Mention Subject: Re: [rordenlab/dcm2niix] Suggestions for GE dicom header info to use in BIDS sidecar (#163) If I'm reading this thread correctly, the current release (without any recompiling) will report a PhaseEncodingDirection for GE data, but the polarity will always be positive? If so, I think that is dangerous, because users will have no indication that polarity is NOT included. If "direction" is to be reported, but without polarity, I think that needs its own specific tag, so as to avoid any mixing up with the existing PhaseEncodingDirection. — You are receiving this because you were mentioned. Reply to this email directly, view it on GitHub, or mute the thread.

mharms commented 6 years ago

If we can't confidently know both the direction and polarity for GE data, then I don't think there should be a PhaseEncodingDirection tag. To have such a tag in that context is going to lead to all sorts of errors. There could be a differently named tag (e.g., PhaseEncodingAxis), but I think it is really problematic to have a PhaseEncodingDirection tag that doesn't convey polarity (because there is no indication from the name that polarity is NOT being conveyed).

neurolabusc commented 6 years ago

@mharms - I always appreciate your perspective. @chrisfilo would you be willing to make a decision on this. The background is that for Siemens scanners we can determine both phaseencoding dimension as well as phaseencoding polarity. Therefore, we can complete populate the JSON header, e.g. "PhaseEncodingDirection": "j-". In contrast, we are only reliably able to determine phase encoding dimension for GE systems, so for these systems we only specify "PhaseEncodingDirection": "j" or "PhaseEncodingDirection": "i". My previous thought was that we give the user all the information we feel confident about. However, @mharms would prefer that we use a separate tag if we are unable to determine both the dimension and the polarity. I can see logic in both approaches, so I am not lobbying for one approach. Can you provide any guidance on this.

mick-d commented 6 years ago

Thanks for the feedback @neurolabusc. For the flag "-DmyReadGeProtocolBlock", I got a warning that flag was not used when running cmake -DmyReadGeProtocolBlock=ON .. :


(`CMake Warning:
  Manually-specified variables were not used by the project:

    myReadGeProtocolBlock
`)

How should one use this flag?
neurolabusc commented 6 years ago

To build from the command line (while in the console folder) you would run:

g++ -dead_strip -O3 -I. main_console.cpp nii_foreign.cpp nii_dicom.cpp jpg_0XC3.cpp ujpeg.cpp nifti1_io_core.cpp nii_ortho.cpp nii_dicom_batch.cpp  -o dcm2niix -DmyDisableOpenJPEG -DmyReadGeProtocolBlock
mick-d commented 6 years ago

Thanks @neurolabusc. I got the same json file than before when using this custom built (i.e. still without polarity).

Also not sure if this is a problem, but it seems the build command misinterpreted the first option -dead-strip as a debug flag followed by specific flags:

cc1plus: warning: unrecognized gcc debugging option: e
cc1plus: warning: unrecognized gcc debugging option: d
cc1plus: warning: unrecognized gcc debugging option: _
cc1plus: warning: unrecognized gcc debugging option: s
cc1plus: warning: unrecognized gcc debugging option: t
cc1plus: warning: unrecognized gcc debugging option: r
cc1plus: warning: unrecognized gcc debugging option: i
cc1plus: warning: unrecognized gcc debugging option: e
cc1plus: warning: unrecognized gcc debugging option: d
cc1plus: warning: unrecognized gcc debugging option: _
cc1plus: warning: unrecognized gcc debugging option: s
cc1plus: warning: unrecognized gcc debugging option: t
cc1plus: warning: unrecognized gcc debugging option: r
cc1plus: warning: unrecognized gcc debugging option: i
cc1plus: warning: unrecognized gcc debugging option: e
cc1plus: warning: unrecognized gcc debugging option: d
cc1plus: warning: unrecognized gcc debugging option: _
cc1plus: warning: unrecognized gcc debugging option: s
cc1plus: warning: unrecognized gcc debugging option: t
cc1plus: warning: unrecognized gcc debugging option: r
cc1plus: warning: unrecognized gcc debugging option: i
cc1plus: warning: unrecognized gcc debugging option: e
cc1plus: warning: unrecognized gcc debugging option: d
cc1plus: warning: unrecognized gcc debugging option: _
cc1plus: warning: unrecognized gcc debugging option: s
cc1plus: warning: unrecognized gcc debugging option: t
cc1plus: warning: unrecognized gcc debugging option: r
cc1plus: warning: unrecognized gcc debugging option: i
cc1plus: warning: unrecognized gcc debugging option: e
cc1plus: warning: unrecognized gcc debugging option: d
cc1plus: warning: unrecognized gcc debugging option: _
cc1plus: warning: unrecognized gcc debugging option: s
cc1plus: warning: unrecognized gcc debugging option: t
cc1plus: warning: unrecognized gcc debugging option: r
cc1plus: warning: unrecognized gcc debugging option: i
cc1plus: warning: unrecognized gcc debugging option: e
cc1plus: warning: unrecognized gcc debugging option: d
cc1plus: warning: unrecognized gcc debugging option: _
cc1plus: warning: unrecognized gcc debugging option: s
cc1plus: warning: unrecognized gcc debugging option: t
cc1plus: warning: unrecognized gcc debugging option: r
cc1plus: warning: unrecognized gcc debugging option: i
cc1plus: warning: unrecognized gcc debugging option: e
cc1plus: warning: unrecognized gcc debugging option: d
cc1plus: warning: unrecognized gcc debugging option: _
cc1plus: warning: unrecognized gcc debugging option: s
cc1plus: warning: unrecognized gcc debugging option: t
cc1plus: warning: unrecognized gcc debugging option: r
cc1plus: warning: unrecognized gcc debugging option: i
cc1plus: warning: unrecognized gcc debugging option: e
cc1plus: warning: unrecognized gcc debugging option: d
cc1plus: warning: unrecognized gcc debugging option: _
cc1plus: warning: unrecognized gcc debugging option: s
cc1plus: warning: unrecognized gcc debugging option: t
cc1plus: warning: unrecognized gcc debugging option: r
cc1plus: warning: unrecognized gcc debugging option: i
nii_dicom_batch.cpp:635:0: warning: "myReadGeProtocolBlock" redefined
  #define myReadGeProtocolBlock
 ^
<command-line>:0:0: note: this is the location of the previous definition
nii_dicom_batch.cpp: In function ‘int geProtocolBlock(const char*, int, int, int, int*, int*)’:
nii_dicom_batch.cpp:701:121: warning: format ‘%d’ expects argument of type ‘int’, but argument 4 has type ‘size_t {aka long unsigned int}’ [-Wformat=]
   printMessage("GE Protocol Block %s bytes %d compressed, %d uncompressed @ %d\n", filename, geLength, unCmpSz, geOffset);
                                                                                                                         ^
neurolabusc commented 6 years ago

You can ignore the warnings about dead strip (or just remove -dead_strip), that is a compiler script that has a few specific calls for old versions of gcc and for clang/llvm. It sounds like you are attempting to compile an older version of dcm2niix, as #define myReadGeProtocolBlock is commented out in the latest code. You can either re-download the repository or if you cloned the repository you can run git pull from the command line while in your dcm2niix folder to update your source code.

The old code automatically parsed the GE Protocol Block by defining #define myReadGeProtocolBlock. However, since the ADNI dataset provided above shows that the Protocol Block is not a reliable way to infer Phase Encoding Polarity, the new code requires you to explicitly compile for this experimental feature.

leej3 commented 6 years ago

Thanks again for all the help with this. It is disappointing that so many sequences from GE have this issue.

I'm sharing a full session collected by @roopchansinghv with more informative series descriptions. Hopefully this helps in the future if any further testing is required, or perhaps other GE issues. It contains GE product DTI (forward and reverse blip), multiphase EPI (forward and reversed blip), and epiRT (forward and reversed blip). As determined above, phase encoding direction should only be inferred in the last case (with epiRT).

https://drive.google.com/file/d/12KRejaoQ1FyidUQAhPze9DGawJYZqFbs/view?usp=sharing

mick-d commented 6 years ago

Thank you for the help @neurolabusc. However, even with the latest version I did not get the polarity. Below is an extract from the resulting json file (created by the latest dcm2niix code)

"PhaseEncodingDirection": "i",
"ProbablePhaseEncodingDirection": "i",
"ConversionSoftwareVersion": "v1.0.20180403 GCC5.4.0"
xiangruili commented 6 years ago

Thanks for the testing images.

I compared all the dicom tags for following two series: Axial_rsfMRI_multiphase_EPI Axial_rsfMRI_multiphase_EPI_reverse_blip

I believe I have narrowed down the possible difference (which could related to "reverse_blip") to the following tag: (0043,102A) 'OB' 'UserDefineData'. Unfortunately, it is not human-readable.

I am not a GE guy. I know, if we have a dicom file for Siemens, the scanner can import the sequence from the dicom file. Is this true for GE? In other words, can GE system import correct "blip" info from the above two dicom files? If so, we can try to decode "UserDefineData". Otherwise I think it is end of the game for this purpose.

Xiangrui Li, Ph.D. Center for Cognitive and Behavioral Brain Imaging The Ohio State University Phone: 614-292-1847

From: john lee Date: 2018-04-12 09:18 To: rordenlab/dcm2niix CC: xiangruili; Mention Subject: Re: [rordenlab/dcm2niix] Suggestions for GE dicom header info to use in BIDS sidecar (#163) Thanks again for all the help with this. It is disappointing that so many sequences from GE have this issue. I'm sharing a full session collected by @roopchansinghv with more informative series descriptions. Hopefully this helps in the future if any further testing is required, or perhaps other GE issues. It contains GE product DTI (forward and reverse blip), multiphase EPI (forward and reversed blip), and epiRT (forward and reversed blip). As determined above, phase encoding direction should only be inferred in the last case (with epiRT). https://drive.google.com/file/d/12KRejaoQ1FyidUQAhPze9DGawJYZqFbs/view?usp=sharing — You are receiving this because you were mentioned. Reply to this email directly, view it on GitHub, or mute the thread.

neurolabusc commented 6 years ago

@mick-d according to the BIDS spec the possiblevalues:“i”,“j”,“k”,“i-”,“j-”,“k-”. http://bids.neuroimaging.io/bids_spec.pdf so in your case i suggests positive blip direction. As @mharms pointed out, this can cause confusion for GE data where we are not able to determine phase encoding polarity. @chrisfilo may want to weight in here. Options would be to have the option to report i+ (known positive), i- (known negative) and i (known dimension but polarity unknown), or to simply exclude the phaseEncodingDirection tag when either dimension or direction are unknown.

mharms commented 6 years ago

I see that @chrisfilo hasn't had a chance to weigh in yet, but following up on the previous comment from @neurolabusc, the BIDS Spec entry for PhaseEncodingDirection further contains:

The polarity of the phase encoding is assumed to go from zero index to maximum index unless ‘-’ sign is present (then the order is reversed - starting from the highest index instead of zero).

So, the BIDS Spec explicitly indicates that PhaseEncodingDirection by definition must also encode for polarity. And we can't use i by itself as an indication of "known axis but polarity unknown" because the BIDS Spec has already defined that the lack of a sign can be interpreted as "positive" polarity. Thus, if we don't know the polarity for GE, I just don't see how we can use PhaseEncodingDirection as a tag under the current BIDS Spec.

chrisgorgo commented 6 years ago

Hi, Sorry for the late reply - I was travelling.

I agree with @mharms we should not use the PhaseEncodingDirection field if we don't know the polarity. I also agree that using a different field (such as PhaseEncodingAxis) is a better solution for the troublesome GE data.

In retrospect the standard should've probably encode axis and polarity in separate fields. Live and learn - we need to maintain backwards compatibility for now.

neurolabusc commented 6 years ago

As @mharms suggested, PhaseEncodingDirection is provided when both dimension and polarity of image is known, PhaseEncodingAxis is used when only the dimension is known. These parameters are not stored for 3D files (as those do not require slice time correction). For both 2D and 3D files, InPlanePhaseEncodingDirectionDICOM is reported as COL or ROW. Note that when converting 3D images to NIfTI they are rotated to standard space, so for these COL and ROW do not necessarily map to i and j.

leej3 commented 6 years ago

Some updates on this. One of our physicists was looking into some of the proprietary GE blocks and realized that the most reliable way of inferring phase encoding polarity was contained in a byte stream used for reconstruction (0043 102a). These files are organized according to a header file (see here for an example: https://github.com/ScottHaileRobertson/GE-MRI-Tools/blob/master/GePackage/%2BGE/%2BPfile/%2BHeader/%2BRDB15/rdbm.h ). The attached code (that requires exiftool to run) parses two relevant flags set within this bytestream to infer the polarity for all EPI images. I wonder could you incorporate something like this? It should be more reliable than the previous user block (which worked for some functional scans and no diffusion weighted scans).

In using this code it turns out that the above test set that I linked to was incorrect, specifically for the axial_epi_epirt_bottom_up. It appears that certain settings on the scanner break the previous way of inferring the polarity in a way that we could not detect through VIEWORDER etc. My metric for it being incorrect is that it does not match the output of this code and a visual inspection of the distortion suggests this is the case.

C file that requires exiftool:

#include <stdio.h>
#include <stdlib.h>
#include <string.h>

int main(int argc, char *argv[])
{
  if (argc < 2) {
    fprintf(stderr, "Syntax: %s [-v] [-f] image1.dcm ...\n", argv[0]);
    exit(1);
  }

  int verbose = 0;
  int force   = 0;

  int argi;
  for (argi=1; argi<argc; ++argi) {

    if (!strcmp(argv[argi], "-v")) {
      verbose = !verbose;
      continue;
    }
    if (!strcmp(argv[argi], "-f")) {
      force = !force;
      continue;
    }

    /*
    ** -- Read the Dicom Tag stream --
    ** 
    ** Dicom tag: (0043,102a) "UserDefinedData"
    */
    char *fname = argv[argi];
    char cmd[256];
    snprintf(cmd, 256, "exiftool -b -s -UserDefinedData %s", fname);

    if (verbose) fprintf(stderr, "Command: %s\n", cmd);

    char data[16384];
    FILE *fp = NULL;

    if ((fp = popen(cmd, "r")) == NULL) {
      fprintf(stderr, "%s", fname);
      perror(cmd);
      exit(1);
    }
    int size = fread(data, 1, 16384, fp);
    if (ferror(fp)) {
      fprintf(stderr, "%s", fname);
      perror("Error reading tag stream data");
      exit(1);
    }
    pclose(fp);

    if (size <= 0) {
      fprintf(stderr, "%s: Incomplete Tag stream data\n", fname);
      continue;
    }

    /*
    ** -- Interpret the Dicom Tag stream --
    */
    short unsigned hdr_offset = *(short unsigned *)(data + 0x0018);

    if (verbose) fprintf(stderr, "hdr_offset: %d 0x%x\n", hdr_offset, hdr_offset);

    char const *hdr = data + hdr_offset;

    /* Determine Tag stream format
    */
    float version = *(float const *)hdr;
    if (!force) {
      if (version < 5.0 || version > 40.0) {
    fprintf(stderr, "Warning: %s: File format incorrect\n", fname);
    exit(1);
      }
    }
    if (verbose) fprintf(stderr, "Version: %f\n", version);

    /* Skip ahead 19 newly added 4-byte fields in Version 26.0
    */
    if (version >= 26.0) hdr += (19*4);

    /* Checking if EPI
    */
    if (!force) {
      int check1 = *(short const *)(hdr + 0x36a);
      int check2 = *(short const *)(hdr + 0x372);
      if (check1 == 0 || check2 == 0) {
    fprintf(stderr, "%s: Warning: Data is not EPI\n", fname);
    continue;
      }
    }

    /*
    ** -- Get the Phase encoding polarity/direction --
    */

    /* Check for PE polarity
    */
    int flag1 = *(short const *)(hdr + 0x0030) & 0x0004;

    if (verbose) fprintf(stderr, "flag1: %d\n", flag1);

    /* Check for ky direction (view order)
    */
    int flag2 = *(int const *)(hdr + 0x0394);

    if (verbose) fprintf(stderr, "flag2: %d\n", flag2);

    /*
    ** -- Output results --
    */
    switch (flag2) {
    case 0:
    case 2:
      if ((flag1 && !flag2) || (!flag1 && flag2)) {
    printf("%s: Bottom up\n", fname);
      }
      else {
    printf("%s: Top down\n", fname);
      }
      break;

    case 1:
      if (flag1) {
    printf("%s: Center out, polarity reversed\n", fname);
      }
      else {
    printf("%s: Center out, polarity normal\n", fname);
      }
      break;

    default:
      printf("%s: Unknown Ky encoding direction", fname);
    }

  }  
  return 0;
}
neurolabusc commented 6 years ago

@leej3 this sounds promising and very interesting. If we can get this reliable it will be easy to implement. However, when I compile your sample program, all the examples provided above by others report "Warning: Data is not EPI", even though they are clearly EPI images. Do you have examples where this works, or does your code need to have the checks modified.

leej3 commented 6 years ago

You can use -f to force it to read out something.

It’s ok if those checks don’t work. If you can implement a more general check for whether the scan is an EPI go ahead. The checks were just trying to pre-empt some obvious errors.

John

On May 31, 2018, at 7:16 PM, Chris Rorden notifications@github.com wrote:

@leej3 this sounds promising and very interesting. If we can get this reliable it will be easy to implement. However, when I compile your sample program, all the examples provided above by others report "Warning: Data is not EPI", even though they are clearly EPI images. Do you have examples where this works, or does your code need to have the checks modified.

— You are receiving this because you were mentioned. Reply to this email directly, view it on GitHub, or mute the thread.

neurolabusc commented 6 years ago

@leej3 and @mick-d - perhaps you can check out the latest commit and provide any insight or suggestions.

neurolabusc commented 6 years ago

@agt24 @leej3 and @mick-d would any of you be able to acquire a minimal phantom dataset for validating GE parameters. Smaller files (low resolution scans with few slices) are ideal. I would like to have something similar to the Siemens dcm_qa dataset that we validate dcm2niix against with every new commit, so smaller is faster and easier to host on Github. A minimal dataset might be EPI scans that demonstrate ROW versus COLumn phase encoding, each with positive and negative phase encoding. It would also be nice to have EPI volumes with ascending, descending and interleaved slice order (for slice time correction). The dataset would need to be public (BSD, MIT, CreativeCommons, etc). I do not have access to GE hardware, but I think this would really help the community. I do host a number of GE datasets on the dcm2niix wiki but none of these fit the bill. Also, are there other special situations where it would be nice to have a GE reference dataset (e.g. fieldmaps, multi-echo images in general, phase/real/imaginary data)?

neurolabusc commented 6 years ago

@leej3 your current program does not discriminate between the Stock GE EPI, both directions sample you uploaded. Can this be detected?

$ ./fx ge_epi_2.5_mm_forward-00001.dcm Command: exiftool -b -s -UserDefinedData ge_epi_2.5_mm_forward-00001.dcm hdr_offset: 3232 0xca0 Version: 26.002001 @ 3232 4182 checks 0 0 ge_epi_2.5_mm_forward-00001.dcm: Warning: Data is not EPI flag1: 0 flag2: 0 ge_epi_2.5_mm_forward-00001.dcm: Top down

$ ./fx ge_epi_2.5_mm_reverse-00001.dcm Command: exiftool -b -s -UserDefinedData ge_epi_2.5_mm_reverse-00001.dcm hdr_offset: 3232 0xca0 Version: 26.002001 @ 3232 4182 checks 0 0 ge_epi_2.5_mm_reverse-00001.dcm: Warning: Data is not EPI flag1: 0 flag2: 0 ge_epi_2.5_mm_reverse-00001.dcm: Top down

$ ./fx axial_dti_24dirs_flipped_rkv-00040.dcm Command: exiftool -b -s -UserDefinedData axial_dti_24dirs_flipped_rkv-00040.dcm hdr_offset: 3232 0xca0 Version: 26.002001 @ 3232 4182 checks 16448 0 axial_dti_24dirs_flipped_rkv-00040.dcm: Warning: Data is not EPI flag1: 4 flag2: 0 axial_dti_24dirs_flipped_rkv-00040.dcm: Bottom up

$ ./fx axial_dti_b_1000_rkv-00019.dcm Command: exiftool -b -s -UserDefinedData axial_dti_b_1000_rkv-00019.dcm hdr_offset: 3232 0xca0 Version: 26.002001 @ 3232 4182 checks 16576 0 axial_dti_b_1000_rkv-00019.dcm: Warning: Data is not EPI flag1: 0 flag2: 0 axial_dti_b_1000_rkv-00019.dcm: Top down

leej3 commented 6 years ago

This looks great. Thank you. Just to start addressing some of your points/questions:

None of the options are related to spiral. Only the first two are likely to be seen with bold and diffusion scans at this time. The latter two, which we included for completeness, are for scans with multiple shots when the readout starts in the center of k-space . The slice order is probably encoded in this reconstruction block. Just to clarify, are you defining slice order as described here?

http://dbic.dartmouth.edu/wiki/index.php/Slice_Acquisition_Order

@jad11nih has done this work and will look into it and report back and also change the case to if/else statements. Also we’re seeing some issues with the file and offset length. We’ll get back to you on that too.

neurolabusc commented 6 years ago

Yes, it would be great to detect the order of 2D EPI acquisitions in a 3D volume - this would help for slice time correction. The ultimate output would be the BIDS SliceTiming tag, which would not only reveal the order of slices, but also discriminate whether there is a pause between the acquisition of the last 2D slice of a 3D volume and beginning of the next volume (as seen with Sparse imaging, or on Philips systems with prospect. motion corr. enabled):

"SliceTiming": [
        3.04,
        0,
        3.14,
        0.1,
        3.2425,
..
        5.9775,
        2.9375  ],

If it helps, you can insert the commented code below into the latest commit of dcm2niix and it will save the entire contents of 0043,102A as a binary file - so you do not need exiftool and the complete length of the header may be a clue regarding which data blocks GE has inserted.

                if (lLength < 4224) {
                    printMessage(" GE header too small to be valid\n");
                    break;
                }
                //debug code to export binary data
                // FILE *pFile = fopen("ge.bin", "wb");
        // fwrite(&buffer[lPos], 1, lLength, pFile);
                // fclose (pFile);

A couple comments would really help understand your solution, for example I think int flag2 = *(int const *)(hdr + 0x0394); refers to the following line from rdbm.h short rdb_hdr_kydir; /* Ky traversal dircetion 0: top-down, 1:center out */ and if so, I would use the same naming and comments.

Again, thanks for all your work. Feels like we are very close.

mharms commented 6 years ago

I’ve been following this tangentially. Is the goal to ultimately use the same json field names that we are using for Siemens data? And if so, have you figured out the relationship between the GE and Siemen’s polarity conventions?

Cheers, -MH

-- Michael Harms, Ph.D.

Associate Professor of Psychiatry Washington University School of Medicine Department of Psychiatry, Box 8134 660 South Euclid Ave. Tel: 314-747-6173 St. Louis, MO 63110 Email: mharms@wustl.edu From: Chris Rorden notifications@github.com Reply-To: rordenlab/dcm2niix reply@reply.github.com Date: Friday, June 1, 2018 at 10:40 AM To: rordenlab/dcm2niix dcm2niix@noreply.github.com Cc: "Harms, Michael" mharms@wustl.edu, Mention mention@noreply.github.com Subject: Re: [rordenlab/dcm2niix] Suggestions for GE dicom header info to use in BIDS sidecar (#163)

Yes, it would be great to detect the order of 2D EPI acquisitions in a 3D volume - this would help for slice time correctionhttp://www.mccauslandcenter.sc.edu/crnl/tools/stc. The ultimate output would be the BIDS SliceTiming tag, which would not only reveal the order of slices, but also discriminate whether there is a pause between the acquisition of the first slice and beginning of the next volume (as seen with Sparse imaging, or on Philips systems with prospect. motion corr. enabled):

"SliceTiming": [

           3.04,

           0,

           3.14,

           0.1,

           3.2425,

..

           5.9775,

           2.9375  ],

If it helps, you can insert the commented code below into the latest commit of dcm2niix and it will save the entire contents of 0043,102A as a binary file - so you do not need exiftool and the complete length of the header may be a clue regarding which data blocks GE has inserted.

           if (lLength < 4224) {

                   printMessage(" GE header too small to be valid\n");

                   break;

           }

           //debug code to export binary data

           // FILE *pFile = fopen("ge.bin", "wb");

           // fwrite(&buffer[lPos], 1, lLength, pFile);

           // fclose (pFile);

A couple comments would really help understand your solution, for example I think int flag2 = (int const )(hdr + 0x0394); refers to the following line from rdbm.hhttps://github.com/ScottHaileRobertson/GE-MRI-Tools/blob/master/GePackage/%2BGE/%2BPfile/%2BHeader/%2BRDB15/rdbm.h short rdb_hdr_kydir; / Ky traversal dircetion 0: top-down, 1:center out / and if so, I would use the same naming and comments.

Again, thanks for all your work. Feels like we are very close.

— You are receiving this because you were mentioned. Reply to this email directly, view it on GitHubhttps://github.com/rordenlab/dcm2niix/issues/163#issuecomment-393919835, or mute the threadhttps://github.com/notifications/unsubscribe-auth/ADswBtYD3f1tyuOZn9-z-zrwESEVNXZXks5t4WBdgaJpZM4SKb17.

The materials in this message are private and may contain Protected Healthcare Information or other information of a sensitive nature. If you are not the intended recipient, be advised that any unauthorized use, disclosure, copying or the taking of any action in reliance on the contents of this information is strictly prohibited. If you have received this email in error, please immediately notify the sender via telephone or return mail.

leej3 commented 6 years ago

How did you arrive at 4224 for the minimum block length? We have files that are producing 3600 for which the code still works.

neurolabusc commented 6 years ago

Latest commit now includes three checks for minimum block lengths - first to make sure there is enough space to read the offset, then enough space to read the version, and finally once the version is known we can have a good idea of final size.

leej3 commented 6 years ago

That's it for now. Slice timing/order will take a little more digging.

Also, the previous failed tests may work now. There was a misplaced & 0x0004 on flag2.

neurolabusc commented 6 years ago

Fixed the flag, but still does not distinguish ge_epi_2.5_mm_reverse from ge_epi_2.5_mm_forward samples.

leej3 commented 6 years ago

agreed. will keep you updated

mick-d commented 6 years ago

@neurolabusc I am not sure the collaborator can but i will try, thanks!

neurolabusc commented 6 years ago

@leej3, @mick-d, @nikadon, @jdkent @brussj @mick-d @mick-d @vmagnotta Since each of you works at a center with GE hardware, I would appreciate it if you could test the latest code for dcm2niix. In theory, this seems to detect phase encoding polarity as well as slice timing.

The slice timing should match how slices are saved to disk in the NIfTI file. For example, as a convenience dcm2niix stores axial slices on disk in the order foot -> head. If the data is flipped, so are the slice times, so the SliceTiming BIDS vector should be descending if you acquired in the images in the (descending) order head -> foot and ascending if you acquire slices in the (ascending) order foot -> head.

Ideally (and for Siemens) I would like to be able to populate the BIDS phase encoding direction tag:

If anyone wants to acquire a small validation dataset that could be shared publicly I would be grateful. I ideally would like a very compact dataset that can be used to validate dcm2niix everytime we submit a patch. For this, the ideal functional data has very few slices, very low resolution in plane, and just a few time-points. It helps if the head has a clear left-right fiducial (e.g. saline bag). I would suggest a sparse design where the TA is rapid (~2 sec) but the TR is long (~10s). This gives good SNR (T1 effects) and also allows the participant to intentionally rotate their head during the last volume. This head movement ensures that the slice timing is correctly encoded.

ihnorton commented 6 years ago

I just debugged conversion of some multi-shell GE data (from a 3.0T DISCOVERY MR750 scanner), which encodes the b-value by varying the length of gradient stored in the header. As far as I can tell, dcm2niix does not currently understand such encoding, as it reports the same b-value for all volumes in this dataset (in other words, 0043,1039 is the same in all volumes, but needs to be scaled by norm(g)^2).

I will be suggesting to fix this in DWIConvert, and wanted to give a heads-up here as something to consider double-checking if test data acquisitions are being done. If I get any sharable data I will pass it on -- I'm asking for permission to share the data mentioned above, but usual issues apply.

Also, in the current BIDS JSON output, dcm2nii does not appear to be saving the raw contents of the DICOM tags from which it chose to extract b-value and gradient information. Doing so in the future would be immensely helpful for working around/debugging issues like the above, for example when re-conversion is difficult.

brussj commented 6 years ago

Chris-

We'll (V. Magnotta, J. Kent, myself) work to gather some data this Friday and will get back to you when we have something to share with the group.

-Joel

From: Chris Rorden notifications@github.com Sent: Saturday, August 25, 2018 10:55:14 AM To: rordenlab/dcm2niix Cc: Bruss, Joel E; Mention Subject: [External] Re: [rordenlab/dcm2niix] Suggestions for GE dicom header info to use in BIDS sidecar (#163)

@leej3https://github.com/leej3, @mick-dhttps://github.com/mick-d, @nikadonhttps://github.com/nikadon, @jdkenthttps://github.com/jdkent @brussjhttps://github.com/brussj @mick-dhttps://github.com/mick-d @mick-dhttps://github.com/mick-d @vmagnottahttps://github.com/vmagnotta Since each of you works at a center with GE hardware, I would appreciate it if you could test the latest code for dcm2niix. In theory, this seems to detect phase encoding polarity as well as slice timing.

The slice timing should match how slices are saved to disk in the NIfTI file. For example, as a convenience dcm2niix stores axial slices on disk in the order foot -> head. If the data is flipped, so are the slice times, so the SliceTiming BIDS vector should be descending if you acquired in the images in the (descending) order head -> foot and ascending if you acquire slices in the (ascending) order foot -> head.

Ideally (and for Siemens) I would like to be able to populate the BIDS phase encoding direction tag:

If anyone wants to acquire a small validation dataset that could be shared publicly I would be grateful. I ideally would like a very compact dataset that can be used to validate dcm2niix everytime we submit a patch. For this, the ideal functional data has very few slices, very low resolution in plane, and just a few time-points. It helps if the head has a clear left-right fiducial (e.g. saline bag). I would suggest a sparse design where the TA is rapid (~2 sec) but the TR is long (~10s). This gives good SNR (T1 effects) and also allows the participant to intentionally rotate their head during the last volume. This head movement ensures that the slice timing is correctly encodedhttps://www.mccauslandcenter.sc.edu/crnl/tools/stc.

— You are receiving this because you were mentioned. Reply to this email directly, view it on GitHubhttps://github.com/rordenlab/dcm2niix/issues/163#issuecomment-415978630, or mute the threadhttps://github.com/notifications/unsubscribe-auth/AXjibItS66dCJMAoUugRdxWXg8zSd77bks5uUXNigaJpZM4SKb17.

Notice: This UI Health Care e-mail (including attachments) is covered by the Electronic Communications Privacy Act, 18 U.S.C. 2510-2521 and is intended only for the use of the individual or entity to which it is addressed, and may contain information that is privileged, confidential, and exempt from disclosure under applicable law. If you are not the intended recipient, any dissemination, distribution or copying of this communication is strictly prohibited. If you have received this communication in error, please notify the sender immediately and delete or destroy all copies of the original message and attachments thereto. Email sent to or from UI Health Care may be retained as required by law or regulation. Thank you.

leej3 commented 6 years ago

Hi @brussj. At my site, we were hoping to make an attempt at this too next week. It'll be great to get some data for this out there. I've put down some notes of some of things we wish to test. Feel free to take any of the ideas, add your own, especially after the scan session if you get an "if only we had..." feeling. Thanks.

https://docs.google.com/document/d/1ii--eAuvvP-RIuJNvJjkz-1h0LgyfsTATkPwkuKH7uc/edit?usp=sharing

neurolabusc commented 6 years ago

@ihnorton can you please test the most recent code on your sample dataset. @nikadon acquired a sample GE multi-shell dataset, but it does not exhibit the behavior you describe. His system is also a 3.0T DISCOVERY MR750 (running 24_LX_MR_Software_release:DV24.0_R02_1607.b). If you can share a dataset or help nikadon generate a sample that shows this behavior it would be nice to have a validation dataset.

Hopefully my solution works even though I do not have a sample image. If this does not work, could you tune the code commented bVal scaled by norm(g)^2? The GE gradient vectors are stored as text with limited precision, so I added a bit of tolerance for some slop in the gradient length.