Skip to content
New issue

Have a question about this project? Sign up for a free GitHub account to open an issue and contact its maintainers and the community.

By clicking “Sign up for GitHub”, you agree to our terms of service and privacy statement. We’ll occasionally send you account related emails.

Already on GitHub? Sign in to your account

Real & Imaginary field map images to BIDS #455

Closed
dlevitas opened this issue Dec 1, 2020 · 15 comments
Closed

Real & Imaginary field map images to BIDS #455

dlevitas opened this issue Dec 1, 2020 · 15 comments

Comments

@dlevitas
Copy link

dlevitas commented Dec 1, 2020

I've come across two Philips datasets where the field maps were stored as real/imaginary as opposed to magnitude/phase. Following dcm2niix, they appear as as follows (sorted), with each echo containing a real/imaginary pair.

*_e1_imaginary.nii.gz',
*_e1_real.nii.gz',
*_e2_imaginary.nii.gz',
*_e2_real.nii.gz'

Once fixed, would this be the BIDS specification field map case of 2 magnitude 2 phase images?

*_magnitude1.nii.gz
*_phase1.nii.gz
*_magnitude2.nii.gz
*_phase2.nii.gz

One last question. In these cases I've only ever seen to echoes with a real/imaginary images pair for each echo. This is typically how such field map acquisitions are stored or are there other ways that are common?

@neurolabusc
Copy link
Collaborator

Precisely, you will want to use fslcomplex to convert to the format supported by FSL FUGUE. The types of images stored reflect the settings on the console. I really do not have an overview of what is common. On our equipment, I tend to set up sequences to save data in such a way to require the least post processing. I would strongly recommend centers carefully evaluate their settings prior to data collection. For our own studies, I typically bring on two MRI physicists as consultants to evaluate our proposed sequences prior to acquiring data from participants.

@dlevitas
Copy link
Author

dlevitas commented Dec 1, 2020

Perfect, thanks! I actually stumbled across your nii_complex2MagPhase.m code (which I converted to Python), and a visual inspection of the resulting images looks like what I'd expect for mag/phase images. Is this way equivalent to the FSL FUGUE steps, or would there be some sort of incompatibility?

@neurolabusc
Copy link
Collaborator

Nice, I would convert the phase image to a FLOAT32 datatype and scale the voxel values to radians and then you could use Prelude to unwrap the data.

Alternatively, you could combine the phase and magnitude image as a NIfTI complex datatype, and then use prelude with the --complex option.

Once you have completed this, I would encourage you to share your Python code as a GitHub repository to help others.

@dlevitas
Copy link
Author

dlevitas commented Dec 1, 2020

Based on your code example, the phase image(s) should indeed be FLOAT32 and scaled to radians already; however, I didn't unwrap the phases images with FSL's prelude, as I was under the impression that the resulting images were "ready to go". To clarify, in order to use the mag/phase images in a preprocessing pipeline (e.g. FSL, fMRIPrep) the phase image(s) would also need to be unwrapped with prelude?

@drmclem
Copy link

drmclem commented Dec 1, 2020

HI - just to chip in - there are two very different workflows you can use with the Philips system depending on how you set the scans up. The default B0 mapping protocols for fMRI produce a greyscale field map calibrated in Hz. This map is derived from the phase maps from a dual echo sequence with unwrapping and processing done for you. This would not be able to be stored as a real/imaginary pair.

If the data is from a much older software release where B0 mapping was not avalable, or you have decided to "roll you own" protocol from a dual echo gradient sequence then you can either generate the phase maps directly (as a single image) or less ideally save the images as real and imaginary and do the phase mapping from there yourself. In both of these cases you would need to do the unrapping and scaling yourself based on the echo time and intervals - in this case the data would be the source data prior to generating the B0 map - its not clear from your first post which of these cases it is.

However, I would caution whether this data can be used in this way as the set up for extracting valid phase difference information from a two echo sequence in earlier software was not straight forward. If you have access to a Philips clinical scientist at your site I would suggest you reach out them to check this data is valid.

M

@neurolabusc
Copy link
Collaborator

Visual inspection could resolve this, but I would expect raw phase maps will show wrapping. You would want to use prelude to unwrap.

@dlevitas
Copy link
Author

dlevitas commented Dec 1, 2020

@drmclem The files in question come from part of the ADNI dataset collected back in early 2012 on a Philips Achieva. Our institution doesn't have a Philips scanner either, so I'm unable to reach out to a Philips clinical scientist to verify its validity.

@neurolabusc From looking at the generated phase image in FSLeyes, it appears that it's already unwrapped, based on the picture on the FSL FUGUE page. I've shared the generated NIFTI file to double check.

@neurolabusc
Copy link
Collaborator

@drmclem when Philips generates a field map calibrated in Hz, is there a way to explicitly detect this image type?

  1. With PAR/REC images, what value is used for the image_type_mr - is it 18?
  2. With DICOM images, how does one distinguish these field maps from magnitude/phase/real/imaginary data? For example, is a special value stored in the image type (0008,0008) tag?

If there is a way to distinguish these images, we could have dcm2niix automatically detect and name them, avoiding confusion. It would also be great to have an example of a Philips DICOM dataset that demonstrated this field map (an image of a phantom would be fine).

@neurolabusc
Copy link
Collaborator

Latest development commit appends _fieldmaphz to PAR/REC images where image_type_mr is 18. A similar approach could be applied to DICOM images, if we can work out how to distinguish this image type for that format.

@neurolabusc
Copy link
Collaborator

Thanks to @drmclem for providing an exemplar DICOM dataset where the fieldmap is saved in Hz. The latest commit will create the post-fix _fieldmaphz for the filename.

fieldmap

yarikoptic added a commit to neurodebian/dcm2niix that referenced this issue Apr 6, 2021
* commit 'v1.0.20200427-97-g0587941': (22 commits)
  Overflow for Siemens data [missing protocol name] (rordenlab#466)
  MacOS notarization, minor tweaks
  Increase details for series stacking, enhance merge override mode (rordenlab#463)
  Bump version date
  Retain Philips Scaling Values for images where scaling does not vary but volumes must be separated (rordenlab#461)
  Support huge PAR/REC files (rordenlab#460)
  Prevent dti4D overflow
  Fix PAR/REC ADC detection (rordenlab#462)
  DICOM field map calibrated in Hz given name _fieldmaphz (rordenlab#455)
  PAR/REC field map calibrated in Hz given name _fieldmaphz (rordenlab#455)
  Ignore LocationsInAcquisition (0020,1002) if number of slices is not evenly divisible with this value (rordenlab#459)
  Convert DICOM series where bit allocated (0028,0100) varies across slices (rordenlab#458)
  Clear RepetitionTimeExcitation (rordenlab#457)
  leak (rordenlab#454)
  Single file mode memory allocation (rordenlab#454)
  When -n is specified, only save BIDS for requested series (rordenlab#453)
  Use 1st study time if multiple times provided.
  Apple. M1. use C++ not. C
  Fix CMake for Apple Silicon (note similar change required for CloudFlare zlib)
  Only report b-value for GE data with 0043,1039 if EPI (0018,9018) (rordenlab#449)
  ...
@MauricePasternak
Copy link

@drmclem when Philips generates a field map calibrated in Hz, is there a way to explicitly detect this image type?

Would testing a look at (2005,1346) MRSpectroFrequencyUnit help? I've attached gdcmdump output from, apparently, one of those pre-processed fieldmap images output in Hz units found in older Philips v3 Achieva

image

Philips_Achieva_v3_gdcmdump_output.md

@neurolabusc neurolabusc reopened this Sep 12, 2024
@neurolabusc
Copy link
Collaborator

@MauricePasternak note that dcm2niix already identifies the Philips fieldmap as Hz without needing to inspect this private tag:

"Units": "Hz",`

Therefore, dcm2niix is creating valid and rich BIDS sidecars. Wrappers for dcm2niix that assign filenames and folders to match the BIDS structure can use this information to correctly sort files. The general division of labor is the that dcmniix creates valid BIDS/NIfTI images and sidecars, and wrappers organize these into BIDS directory hierarchy and filenames (a step that usually requires a configuration file that has details regarding the study intention that can not be determined by the DICOM images in isolation).

@MauricePasternak
Copy link

@neurolabusc The provided example does not result in the "Units": "Hz" addition to the json sidecar.

{
	"Modality": "MR",
	"MagneticFieldStrength": 3,
	"ImagingFrequency": 127.765447,
	"Manufacturer": "Philips",
	"ManufacturersModelName": "Achieva",
	"BodyPartExamined": "BRAIN",
	"PatientPosition": "HFS",
	"SoftwareVersions": "3.2.1\\3.2.1.1",
	"MRAcquisitionType": "2D",
	"ScanningSequence": "GR",
	"SequenceVariant": "SS",
	"ScanOptions": "FC",
	"PulseSequenceName": "FFE",
	"ImageType": ["ORIGINAL", "PRIMARY", "VELOCITY MAP", "P", "PCA", "PHASE"],
	"SeriesNumber": 601,
	"AcquisitionTime": "16:11:51.070000",
	"AcquisitionNumber": 6,
	"PhilipsRWVSlope": 0.0991453,
	"PhilipsRWVIntercept": -203,
	"PhilipsRescaleSlope": 0.0991453,
	"PhilipsRescaleIntercept": -203,
	"PhilipsScaleSlope": 10.0737,
	"UsePhilipsFloatNotDisplayScaling": 1,
	"SliceThickness": 3,
	"SpacingBetweenSlices": 4,
	"SAR": 0.139658,
	"EchoTime": 0.004604,
	"RepetitionTime": 0.688001,
	"MTState": false,
	"FlipAngle": 60,
	"CoilString": "Dual coil",
	"PercentPhaseFOV": 100,
	"PercentSampling": 98.6111,
	"PartialFourierDirection": "FREQUENCY",
	"PhaseEncodingSteps": 71,
	"FrequencyEncodingSteps": 72,
	"PhaseEncodingStepsOutOfPlane": 1,
	"AcquisitionMatrixPE": 71,
	"ReconMatrixPE": 80,
	"WaterFatShift": 1.5585,
	"EstimatedEffectiveEchoSpacing": 2.27069e-05,
	"EstimatedTotalReadoutTime": 0.00179385,
	"AcquisitionDuration": 100.448,
	"PixelBandwidth": 279,
	"PhaseEncodingAxis": "j",
	"ImageOrientationPatientDICOM": [
		0.998648,
		0.0222905,
		0.0469672,
		-0.0223638,
		0.999749,
		0.00103639	],
	"InPlanePhaseEncodingDirectionDICOM": "COL",
	"BidsGuess": ["anat","_acq-SSGRFFE_run-601_part-phase_T2starw"],
	"ConversionSoftware": "dcm2niix",
	"ConversionSoftwareVersion": "v1.0.20240202"
}

I assume that BidsGuess is wrong here (fair enough, it is a guess) since the output of running dcm2niix on the DICOM directory series creates output that matches Case 3 in the BIDS specification on fieldmaps. My only concern is that the Units key is missing in this case.

@MauricePasternak
Copy link

Yes, it does look like Unit isn't added when it should be. The same site updated from v3 to v5 of an Achieva scanner and define the ImageType more explicitly:

{
	"Modality": "MR",
	"MagneticFieldStrength": 3,
	"ImagingFrequency": 127.754609,
	"Manufacturer": "Philips",
	"ManufacturersModelName": "Achieva",
	"BodyPartExamined": "BRAIN",
	"PatientPosition": "HFS",
	"SoftwareVersions": "5.3.1\\5.3.1.1",
	"MRAcquisitionType": "2D",
	"ScanningSequence": "RM",
	"SequenceVariant": "SS",
	"ScanOptions": "FC",
	"PulseSequenceName": "FFE",
	"ImageType": ["ORIGINAL", "PRIMARY", "B0 MAP", "B0", "UNSPECIFIED", "FIELDMAPHZ"],
	"SeriesNumber": 501,
	"AcquisitionTime": "14:27:51.020000",
	"AcquisitionNumber": 5,
	"PhilipsRWVSlope": 0.0991453,
	"PhilipsRWVIntercept": -203,
	"PhilipsRescaleSlope": 0.0991453,
	"PhilipsRescaleIntercept": -203,
	"PhilipsScaleSlope": 10.0737,
	"UsePhilipsFloatNotDisplayScaling": 1,
	"SliceThickness": 3,
	"SpacingBetweenSlices": 4,
	"SAR": 0.139658,
	"EchoNumber": 1,
	"MTState": false,
	"FlipAngle": 60,
	"CoilString": "MULTI COIL",
	"PercentPhaseFOV": 100,
	"PercentSampling": 100,
	"PhaseEncodingSteps": 64,
	"FrequencyEncodingSteps": 64,
	"PhaseEncodingStepsOutOfPlane": 1,
	"AcquisitionMatrixPE": 64,
	"ReconMatrixPE": 64,
	"WaterFatShift": 1.50277,
	"EstimatedEffectiveEchoSpacing": 2.74579e-05,
	"EstimatedTotalReadoutTime": 0.00172985,
	"AcquisitionDuration": 90.8161,
	"PixelBandwidth": 289,
	"PhaseEncodingAxis": "j",
	"ImageOrientationPatientDICOM": [
		0.997658,
		0.057506,
		-0.0370271,
		-0.0641677,
		0.974351,
		-0.215693	],
	"InPlanePhaseEncodingDirectionDICOM": "COL",
	"BidsGuess": ["fmap","_acq-SSRMFFE_run-501_magnitude"],
	"ConversionSoftware": "dcm2niix",
	"ConversionSoftwareVersion": "v1.0.20240202"
}

image

With Unit still missing in the sidecar.

@neurolabusc
Copy link
Collaborator

I am unable to replicate as both the current stable release (v1.0.20240202) and latest development build do populate the field with the validation dataset. I do not have access to Philips hardware, so if your want further feedback you will need to share a sample dataset with my institutional email. It might also be worth checking the provenance of your images to see if some tags were modified. I do notice in the two examples you are providing the BidsGuess is phase for one of your examples and magnitude for the other, suggesting it has not detected these volumes as being derived fieldmaps, rather they are detected as different image types.

As an aside, I would suggest you contact the Philips Clinical Scientist associated with your center, as the core issue here is that these are not valid DICOMs. The specification demands that the derived image (fieldmap) should be given a different series number than the raw images that it is derived from. Siemens shows the correct behavior here. If the specification was followed, it would be much easier for tools like dcm2niix (and humans) to parse. While I do try to make dcm2niix tolerant to DICOM as seen in the real world, it does assume that the inputs are valid and truthful.

The Unix command line code below shows that dcm2niix does populate this tag for the validation dataset:

git clone [email protected]:neurolabusc/dcm_qa_fmap.git
cd dcm_qa_fmap 
./batch.sh
grep '"Units": "Hz"' ./Out/Ph_801_WIP_B0_NS_e2_fieldmaphz.json

will generate the following JSON for the current stable release:

{
	"Modality": "MR",
	"MagneticFieldStrength": 3,
	"ImagingFrequency": 127.763573,
	"Manufacturer": "Philips",
	"ManufacturersModelName": "Achieva dStream",
	"InstitutionName": "Vanderbilt Health and williamson",
	"InstitutionalDepartmentName": "Research MRI",
	"InstitutionAddress": "1161 21st Ave South",
	"DeviceSerialNumber": "17240",
	"StationName": "mr3tb",
	"BodyPartExamined": "BRAIN",
	"PatientPosition": "HFS",
	"SoftwareVersions": "5.3.0\\5.3.0.3",
	"MRAcquisitionType": "3D",
	"SeriesDescription": "B0 NS",
	"ProtocolName": "WIP B0 NS",
	"ScanningSequence": "RM",
	"SequenceVariant": "SS",
	"ScanOptions": "OTHER",
	"PulseSequenceName": "FFE",
	"ImageType": ["ORIGINAL", "PRIMARY", "T1", "MIXED", "REAL", "FIELDMAPHZ"],
	"SeriesNumber": 801,
	"AcquisitionTime": "15:25:9.020000",
	"AcquisitionNumber": 8,
	"PhilipsRWVSlope": 0.2442,
	"PhilipsRWVIntercept": -500,
	"PhilipsRescaleSlope": 0.2442,
	"PhilipsRescaleIntercept": -500,
	"PhilipsScaleSlope": 4.095,
	"UsePhilipsFloatNotDisplayScaling": 1,
	"Units": "Hz",
	"SliceThickness": 4,
	"SpacingBetweenSlices": 4,
	"EchoTime": 0.00152,
	"RepetitionTime": 0.0038045,
	"MTState": false,
	"FlipAngle": 12,
	"CoilString": "MULTI COIL",
	"PercentPhaseFOV": 100,
	"PercentSampling": 78.5398,
	"PhaseEncodingSteps": 64,
	"FrequencyEncodingSteps": 64,
	"PhaseEncodingStepsOutOfPlane": 32,
	"AcquisitionMatrixPE": 64,
	"ReconMatrixPE": 64,
	"WaterFatShift": 0.140136,
	"EstimatedEffectiveEchoSpacing": 2.56032e-06,
	"EstimatedTotalReadoutTime": 0.0001613,
	"AcquisitionDuration": 15.18,
	"PixelBandwidth": 3099.81,
	"ImageOrientationPatientDICOM": [
		1,
		0,
		2.74252e-05,
		0,
		1,
		0	],
	"InPlanePhaseEncodingDirectionDICOM": "ROW",
	"BidsGuess": ["fmap","_acq-SSRMFFE_run-801_fieldmap"],
	"ConversionSoftware": "dcm2niix",
	"ConversionSoftwareVersion": "v1.0.20240202"
}

Sign up for free to join this conversation on GitHub. Already have an account? Sign in to comment
Labels
None yet
Projects
None yet
Development

No branches or pull requests

4 participants