diff --git a/CHANGELOG.md b/CHANGELOG.md index 0614810..94f0114 100644 --- a/CHANGELOG.md +++ b/CHANGELOG.md @@ -1,5 +1,17 @@ This document contains the Spec2nii release history in reverse chronological order. +0.6.1 (Wednesday 18 January 2023) +--------------------------------- +- Fixed conjugation issue introduced by new nifti-mrs package dependency +- SPAR/SDAT pipeline now handles HYPER special case. +- Data/list pipeline now handles HYPER special case. +- Fixed issue with XA Twix PatientSex and TxOffset attributes. +- Reenable Bruker conversion. + +0.6.0 (Wednesday 11th January 2023) +----------------------------------- +- NIfTI-MRS creation/handling/verification now performed by nifti-mrs package. + 0.5.0 (Friday 2nd December 2022) -------------------------------- - NIfTI-MRS version 0.6 diff --git a/requirements.yml b/requirements.yml index c287c92..7f5bb67 100644 --- a/requirements.yml +++ b/requirements.yml @@ -4,5 +4,6 @@ dependencies: - pydicom - pyMapVBVD>=0.5.2 - scipy - - brukerapi>=0.1.4.0 + - brukerapi>=0.1.8 - pandas + - nifti-mrs>=0.1.3 diff --git a/setup.py b/setup.py index c2636bd..597837b 100644 --- a/setup.py +++ b/setup.py @@ -27,7 +27,6 @@ 'spec2nii.Siemens', 'spec2nii.GSL', 'spec2nii.dcm2niiOrientation', - 'spec2nii.nifti_mrs', 'spec2nii.Philips', 'spec2nii.GE'], install_requires=install_requires, diff --git a/spec2nii/GE/ge_pfile.py b/spec2nii/GE/ge_pfile.py index aa0f6b3..7ae8c95 100644 --- a/spec2nii/GE/ge_pfile.py +++ b/spec2nii/GE/ge_pfile.py @@ -43,9 +43,11 @@ # 3rd party modules import numpy as np +from nifti_mrs.create_nmrs import gen_nifti_mrs_hdr_ext +from nifti_mrs.hdr_ext import Hdr_Ext + # Our modules from spec2nii.GE.ge_read_pfile import Pfile -from spec2nii import nifti_mrs from spec2nii.nifti_orientation import NIFTIOrient from spec2nii import __version__ as spec2nii_ver @@ -104,7 +106,7 @@ def _process_svs_pfile(pfile): out_nmrs = [] for dd, mm in zip(data, meta): - out_nmrs.append(nifti_mrs.NIfTI_MRS(dd, orientation.Q44, dwelltime, mm)) + out_nmrs.append(gen_nifti_mrs_hdr_ext(dd, dwelltime, mm, orientation.Q44, no_conj=True)) return out_nmrs, fname_suffix @@ -228,7 +230,7 @@ def fft_and_shift(x, axis): meta = _populate_metadata(pfile) orientation = NIFTIOrient(_calculate_affine_mrsi(pfile)) - return [nifti_mrs.NIfTI_MRS(data, orientation.Q44, dwelltime, meta), ], ['', ] + return [gen_nifti_mrs_hdr_ext(data, dwelltime, meta, orientation.Q44, no_conj=True), ], ['', ] def _calculate_affine_mrsi(pfile): @@ -297,8 +299,9 @@ def _populate_metadata(pfile, water_suppressed=True): 'Use the --override_nucleus option to specify a different nuclide.') nucleus = "1H" - meta = nifti_mrs.hdr_ext(spec_frequency, - nucleus) + meta = Hdr_Ext( + spec_frequency, + nucleus) # Standard defined metadata # # 5.1 MRS specific Tags diff --git a/spec2nii/Philips/philips.py b/spec2nii/Philips/philips.py index 9c2861d..ee85594 100644 --- a/spec2nii/Philips/philips.py +++ b/spec2nii/Philips/philips.py @@ -6,27 +6,32 @@ from ast import literal_eval import numpy as np +from nifti_mrs.create_nmrs import gen_nifti_mrs_hdr_ext +from nifti_mrs.hdr_ext import Hdr_Ext from spec2nii.nifti_orientation import NIFTIOrient, calc_affine -from spec2nii import nifti_mrs from spec2nii import __version__ as spec2nii_ver -def read_sdat_spar_pair(sdat_file, spar_file, shape=None, tags=None): - """_summary_ - - _extended_summary_ - - :param sdat_file: _description_ - :type sdat_file: _type_ - :param spar_file: _description_ - :type spar_file: _type_ - :param shape: _description_, defaults to None - :type shape: _type_, optional - :param tags: _description_, defaults to None - :type tags: _type_, optional - :return: _description_ - :rtype: _type_ +def read_sdat_spar_pair(sdat_file, spar_file, shape=None, tags=None, fileout=None, special=None): + """Read and convert SPAR/SDAT pairs from Philips scanners + + :param sdat_file: Path to .SDAT file + :type sdat_file: pathlib.Path + :param spar_file: Path to .SPAR file + :type spar_file: pathlib.Path + :param shape: List of dimension shapes to reshape array, defaults to None + :type shape: list of ints, optional + :param tags: List of higher dimension tags, defaults to None + :type tags: List of strings, optional + :param fileout: File name string passed by user, defaults to None/stem of file + :type fileout: str, optional + :param special: Identifier for special-cased sequence, defaults to None + :type special: str, optional + :return: List of NIFTI-MRS objects + :rtype: list of nifti_mrs.nifti_mrs.NIFTI_MRS + :return: List of file name parts + :rtype: list of strings """ spar_params = read_spar(spar_file) data = read_sdat(sdat_file, @@ -62,7 +67,21 @@ def read_sdat_spar_pair(sdat_file, spar_file, shape=None, tags=None): affine = np.diag(np.array([10000, 10000, 10000, 1])) orientation = NIFTIOrient(affine) - return [nifti_mrs.NIfTI_MRS(data, orientation.Q44, dwelltime, meta), ] + # name of output + if fileout is not None: + mainStr = fileout + else: + mainStr = sdat_file.stem + + # Special cases + if spar_params['scan_id'].lower() == 'hyper'\ + or (special is not None and + special.lower() == 'hyper'): + return _special_case_hyper(data, dwelltime, meta, orientation.Q44, mainStr) + else: + # Normal case + return [gen_nifti_mrs_hdr_ext(data, dwelltime, meta, orientation.Q44, no_conj=True), ],\ + [mainStr, ] def read_spar(filename): @@ -143,8 +162,9 @@ def spar_to_nmrs_hdrext(spar_dict): # Extract required metadata and create hdr_ext object cf = float(spar_dict["synthesizer_frequency"]) / 1E6 - obj = nifti_mrs.hdr_ext(cf, - spar_dict["nucleus"]) + obj = Hdr_Ext( + cf, + spar_dict["nucleus"]) def set_standard_def(nifti_mrs_key, location, key, cast=None): try: @@ -266,3 +286,37 @@ def _vax_to_ieee_single_float(data): # may want to raise an exception here ... return f + + +def _special_case_hyper(data, dwelltime, meta, orientation, fout_str): + # Reorganise the data. This unfortunately makes hardcoded assumptions about the size of each part. + data_short_te = data[:, :, :, :, :32] + data_edited = data[:, :, :, :, 32:] + data_edited = data_edited.T.reshape((56, 4, data.shape[3], 1, 1, 1)).T + + meta_short_te = meta.copy() + meta_edited = meta.copy() + + edit_pulse_1 = 1.9 + edit_pulse_2 = 4.58 + edit_pulse_off = 4.18 + dim_info = "HERCULES j-difference editing, four conditions" + dim_header = {"EditCondition": ["A", "B", "C", "D"]} + edit_pulse_val = { + "A": {"PulseOffset": [edit_pulse_1, edit_pulse_2], "PulseDuration": 0.02}, + "B": {"PulseOffset": [edit_pulse_off, edit_pulse_2], "PulseDuration": 0.02}, + "C": {"PulseOffset": edit_pulse_1, "PulseDuration": 0.02}, + "D": {"PulseOffset": edit_pulse_off, "PulseDuration": 0.02}} + + meta_edited.set_dim_info( + 0, + 'DIM_EDIT', + dim_info, + dim_header) + + meta_edited.set_dim_info(1, 'DIM_DYN') + meta_edited.set_standard_def("EditPulse", edit_pulse_val) + + return [gen_nifti_mrs_hdr_ext(data_short_te, dwelltime, meta_short_te, orientation, no_conj=True), + gen_nifti_mrs_hdr_ext(data_edited, dwelltime, meta_edited, orientation, no_conj=True)],\ + [fout_str + '_hyper_short_te', fout_str + '_hyper_edited'] diff --git a/spec2nii/Philips/philips_data_list.py b/spec2nii/Philips/philips_data_list.py index 0683d92..b7346ae 100644 --- a/spec2nii/Philips/philips_data_list.py +++ b/spec2nii/Philips/philips_data_list.py @@ -7,10 +7,19 @@ import pandas as pd import numpy as np +from nibabel.nicom.dicomwrappers import wrapper_from_file +from pydicom.errors import InvalidDicomError + +from nifti_mrs.create_nmrs import gen_nifti_mrs_hdr_ext from spec2nii.Philips.philips import read_spar, spar_to_nmrs_hdrext, _philips_orientation +from spec2nii.Philips.philips_dcm import \ + _enhanced_dcm_svs_to_orientation,\ + _is_new_format,\ + _extractDicomMetadata_new,\ + _extractDicomMetadata_old from spec2nii.nifti_orientation import NIFTIOrient -from spec2nii import nifti_mrs + index_options = ['STD', 'REJ', 'PHC', 'FRC', 'NOI', 'NAV'] index_headers = ['typ', 'mix', 'dyn', 'card', 'echo', 'loca', @@ -29,19 +38,39 @@ class TooManyDimError(Exception): pass -def read_data_list_pair(data_file, list_file, spar_file): +def read_data_list_pair(data_file, list_file, aux_file, special_case=None): df, num_dict, coord_dict, os_dict = _read_list(list_file) sorted_data_dict = _read_data(data_file, df) - spar_params = read_spar(spar_file) + if aux_file.suffix.lower() == '.spar': + spar_params = read_spar(aux_file) - # Dwelltime - dwelltime = 1.0 / float(spar_params["sample_frequency"]) + dwelltime = 1.0 / float(spar_params["sample_frequency"]) + orientation = NIFTIOrient(_philips_orientation(spar_params)) + + # Meta + base_meta = spar_to_nmrs_hdrext(spar_params) + else: + try: + dcm_obj = wrapper_from_file(aux_file) + except (InvalidDicomError, FileNotFoundError) as exc: + raise ValueError('Data/List auxiliary file must be a valid Philips .SPAR or DICOM file.') from exc + + orientation = _enhanced_dcm_svs_to_orientation(dcm_obj) + + if _is_new_format(dcm_obj): + base_meta = _extractDicomMetadata_new(dcm_obj) + dwelltime = 1.0 / dcm_obj.dcm_data.SpectralWidth + else: + base_meta = _extractDicomMetadata_old(dcm_obj) + dwelltime = 1.0 / dcm_obj.dcm_data[('2005', '1357')].value - # Orientation - affine = _philips_orientation(spar_params) - orientation = NIFTIOrient(affine) + base_meta.set_standard_def( + 'OriginalFile', + [data_file.name, + list_file.name, + aux_file.name]) data_out = [] name_out = [] @@ -49,12 +78,7 @@ def read_data_list_pair(data_file, list_file, spar_file): data = sorted_data_dict[data_type] name = data_type - # Meta - meta = spar_to_nmrs_hdrext(spar_params) - meta.set_standard_def('OriginalFile', - [data_file.name, - list_file.name, - spar_file.name]) + meta = base_meta.copy() kept_ind = [] for ii, sha in zip(indices, data.shape[1:]): @@ -67,24 +91,63 @@ def read_data_list_pair(data_file, list_file, spar_file): f' Dimensions are {kept_ind} with shape {out_data.shape[1:]}.' ' NIFTI-MRS can only handle three dynamic dimensions. Unsure how to proceed.') - unknown_counter = 0 - for idx, ki in enumerate(kept_ind): - if ki in defaults: - meta.set_dim_info(idx, - defaults[ki], - info=f'data/list dim {ki}') - else: - meta.set_dim_info(idx, - f'DIM_USER_{unknown_counter}', - info=f'data/list dim {ki}') - unknown_counter += 1 + # Special cases + if data_type == 'STD_0'\ + and (special_case == 'hyper' or 'hyper' in meta['ProtocolName'].lower()): + out_hyper, meta_hyper = _special_case_hyper(out_data, meta) + # Handle the main acqusition of the HYPER (short TE + editing) sequence + + # Insert spatial dimensions + out_shortte = out_hyper[0].reshape((1, 1, 1) + out_hyper[0].shape) + out_edited = out_hyper[1].reshape((1, 1, 1) + out_hyper[1].shape) + + # Apply conjugate + out_shortte = out_shortte.conj() + out_edited = out_edited.conj() + + data_out.append( + gen_nifti_mrs_hdr_ext(out_shortte, dwelltime, meta_hyper[0], orientation.Q44, no_conj=True)) + name_out.append('hyper_short_te') + + data_out.append( + gen_nifti_mrs_hdr_ext(out_edited, dwelltime, meta_hyper[1], orientation.Q44, no_conj=True)) + name_out.append('hyper_edited') + + continue + + elif data_type == 'STD_1'\ + and (special_case == 'hyper' or 'hyper' in meta['ProtocolName'].lower()): + # Handle the water ref acqusition of the HYPER sequence + + meta.set_dim_info( + 1, + 'DIM_USER_0', + info='HYPER water reference') + + name = 'hyper_water_ref' + + else: + # Default case: assign default DIM Tags. + unknown_counter = 0 + for idx, ki in enumerate(kept_ind): + if ki in defaults: + meta.set_dim_info( + idx, + defaults[ki], + info=f'data/list dim {ki}') + else: + meta.set_dim_info( + idx, + f'DIM_USER_{unknown_counter}', + info=f'data/list dim {ki}') + unknown_counter += 1 # Insert spatial dimensions out_data = out_data.reshape((1, 1, 1) + out_data.shape) # Apply conjugate out_data = out_data.conj() - data_out.append(nifti_mrs.NIfTI_MRS(out_data, orientation.Q44, dwelltime, meta)) + data_out.append(gen_nifti_mrs_hdr_ext(out_data, dwelltime, meta, orientation.Q44, no_conj=True)) name_out.append(name) return data_out, name_out @@ -132,6 +195,18 @@ def _read_data(data_file, df): output_dict[f'{tt}_{mix}'][ind] = raw[offset:(offset + dsize)] offset += dsize + # Slight hack - remove all zero data if there is any (skipped indices) + for key in output_dict: + if key == 'NOI': + continue + keyparts = key.split('_') + cdf = df[df.typ == keyparts[0]][df.mix == int(keyparts[1])] + used_df = cdf[indices].loc[:, cdf[indices].max() > 0] + used_col_indices = [indices.index(col) for col in used_df.columns] + for idx, col_idx in enumerate(used_col_indices): + used_indices_in_col = used_df.iloc[:, idx].unique() + output_dict[key] = np.take(output_dict[key], used_indices_in_col, axis=col_idx + 1) + return output_dict @@ -161,7 +236,7 @@ def _read_list(list_file): r'(\d+)\s+\d+\s+\d+\s+number_of_([a-z0-9_]+)\s+:\s+(\d+)', line) if matched[2] not in num_dict: - num_dict[matched[2]] = np.zeros((n_mixes), dtype=np.int) + num_dict[matched[2]] = np.zeros((n_mixes), dtype=int) num_dict[matched[2]][int(matched[1])] = int(matched[3]) elif '_range' in line: matched = re.search( @@ -169,7 +244,7 @@ def _read_list(list_file): line, flags=re.IGNORECASE) if matched[3] not in coord_dict: - coord_dict[matched[3]] = np.zeros((n_mixes, n_echoes, 2), dtype=np.int) + coord_dict[matched[3]] = np.zeros((n_mixes, n_echoes, 2), dtype=int) coord_dict[matched[3]][int(matched[1]), int(matched[2]), 0] = int(matched[4]) coord_dict[matched[3]][int(matched[1]), int(matched[2]), 1] = int(matched[5]) elif '_oversample_factor' in line: @@ -187,3 +262,37 @@ def _read_list(list_file): df[index_headers[1:]] = df[index_headers[1:]].apply(pd.to_numeric) return df, num_dict, coord_dict, os_dict + + +# Special cases +def _special_case_hyper(data, meta): + + data_short_te = data[:, :, :32] + data_edited = data[:, :, 32:] + data_edited = data_edited.T.reshape((56, 4, data.shape[1], data.shape[0])).T + + meta_short_te = meta.copy() + meta_edited = meta.copy() + + edit_pulse_1 = 1.9 + edit_pulse_2 = 4.58 + edit_pulse_off = 4.18 + dim_info = "HERCULES j-difference editing, four conditions" + dim_header = {"EditCondition": ["A", "B", "C", "D"]} + edit_pulse_val = { + "A": {"PulseOffset": [edit_pulse_1, edit_pulse_2], "PulseDuration": 0.02}, + "B": {"PulseOffset": [edit_pulse_off, edit_pulse_2], "PulseDuration": 0.02}, + "C": {"PulseOffset": edit_pulse_1, "PulseDuration": 0.02}, + "D": {"PulseOffset": edit_pulse_off, "PulseDuration": 0.02}} + + meta_edited.set_dim_info( + 1, + 'DIM_EDIT', + dim_info, + dim_header) + + meta_edited.set_dim_info(2, 'DIM_DYN') + meta_edited.set_standard_def("EditPulse", edit_pulse_val) + + return [data_short_te, data_edited],\ + [meta_short_te, meta_edited] diff --git a/spec2nii/Philips/philips_dcm.py b/spec2nii/Philips/philips_dcm.py index 53ae2a0..ffc3773 100644 --- a/spec2nii/Philips/philips_dcm.py +++ b/spec2nii/Philips/philips_dcm.py @@ -2,12 +2,16 @@ Author: William Clarke Copyright (C) 2020 University of Oxford """ +from datetime import datetime + import numpy as np import nibabel.nicom.dicomwrappers + +from nifti_mrs.create_nmrs import gen_nifti_mrs_hdr_ext +from nifti_mrs.hdr_ext import Hdr_Ext + from spec2nii.dcm2niiOrientation.orientationFuncs import dcm_to_nifti_orientation from spec2nii.nifti_orientation import NIFTIOrient -from spec2nii import nifti_mrs -from datetime import datetime from spec2nii import __version__ as spec2nii_ver @@ -155,9 +159,9 @@ def all_equal(lst): meta_used.set_dim_info(0, tag) # Create NIFTI MRS object. - nifti_mrs_out.append(nifti_mrs.NIfTI_MRS(combined_data, or_used.Q44, dt_used, meta_used)) + nifti_mrs_out.append(gen_nifti_mrs_hdr_ext(combined_data, dt_used, meta_used, or_used.Q44, no_conj=True)) if ref_list[0] is not None: - nifti_mrs_out.append(nifti_mrs.NIfTI_MRS(combined_ref, or_used.Q44, dt_used, meta_ref_used)) + nifti_mrs_out.append(gen_nifti_mrs_hdr_ext(combined_ref, dt_used, meta_ref_used, or_used.Q44, no_conj=True)) else: for idx, (dd, rr, oo, dt, mm, mr, ff) in enumerate(zip(data_list, ref_list, @@ -169,12 +173,12 @@ def all_equal(lst): # Add original files to nifti meta information. mm.set_standard_def('OriginalFile', [str(ff.name), ]) fnames_out.append(f'{mainStr}_{idx:03}') - nifti_mrs_out.append(nifti_mrs.NIfTI_MRS(dd, oo.Q44, dt, mm)) + nifti_mrs_out.append(gen_nifti_mrs_hdr_ext(dd, dt, mm, oo.Q44, no_conj=True)) if rr is not None: mr.set_standard_def('OriginalFile', [str(ff.name), ]) fnames_out.append(f'{mainStr}_ref_{idx:03}') - nifti_mrs_out.append(nifti_mrs.NIfTI_MRS(rr, oo.Q44, dt, mr)) + nifti_mrs_out.append(gen_nifti_mrs_hdr_ext(rr, dt, mr, oo.Q44, no_conj=True)) return nifti_mrs_out, fnames_out @@ -265,7 +269,7 @@ def _process_philips_fid(img, verbose): return spec_data_main, spec_data_ref, currNiftiOrientation, dwelltime, meta, meta_r -def _enhanced_dcm_svs_to_orientation(img, verbose): +def _enhanced_dcm_svs_to_orientation(img, verbose=False): '''Convert the Volume Localization Sequence (0018,9126) enhanced DICOM tag to a 4x4 affine format. @@ -325,7 +329,8 @@ def _enhanced_dcm_svs_to_orientation(img, verbose): from scipy.spatial.transform import Rotation rot = Rotation.from_euler('xyz', [-angle_lr, -angle_ap, angle_hf], degrees=True) - imageOrientationPatient = rot.as_matrix() + # THIS NEEDS TESTING!!! CODE WILL PASS BUT I DON"T TRUST IT AT ALL. + imageOrientationPatient = rot.as_matrix()[:, :2].T imagePositionPatient = [-shift_lr, -shift_ap, shift_hf] except KeyError: @@ -357,8 +362,19 @@ def _extractDicomMetadata_old(dcmdata, water_suppressed=True): """ # Extract required metadata and create hdr_ext object - obj = nifti_mrs.hdr_ext(dcmdata.dcm_data.TransmitterFrequency, - dcmdata.dcm_data.ResonantNucleus) + try: + frequency = dcmdata.dcm_data.TransmitterFrequency + except AttributeError: + frequency = float(dcmdata.dcm_data[('2001', '1083')].value) + + try: + nucleus = dcmdata.dcm_data.ResonantNucleus + except AttributeError: + nucleus = dcmdata.dcm_data[('2001', '1087')].value + + obj = Hdr_Ext( + frequency, + nucleus) # private_mrs_tags = dcmdata[('2005', '140f')][0] @@ -369,21 +385,33 @@ def set_if_present(tag, val): # Standard metadata # # 5.1 MRS specific Tags # 'EchoTime' - echo_time = float(dcmdata.dcm_data.PerFrameFunctionalGroupsSequence[0] - .MREchoSequence[0].EffectiveEchoTime) * 1E-3 + try: + echo_time = float(dcmdata.dcm_data.PerFrameFunctionalGroupsSequence[0] + .MREchoSequence[0].EffectiveEchoTime) * 1E-3 + except AttributeError: + echo_time = float(dcmdata.dcm_data[('2005', '1310')].value) obj.set_standard_def('EchoTime', echo_time) + # 'RepetitionTime' - rep_tim = float(dcmdata.dcm_data.PerFrameFunctionalGroupsSequence[0] - .MRTimingAndRelatedParametersSequence[0].RepetitionTime) / 1E3 - obj.set_standard_def('RepetitionTime', rep_tim) + try: + rep_time = float(dcmdata.dcm_data.PerFrameFunctionalGroupsSequence[0] + .MRTimingAndRelatedParametersSequence[0].RepetitionTime) / 1E3 + except AttributeError: + rep_time = dcmdata.dcm_data.RepetitionTime + obj.set_standard_def('RepetitionTime', rep_time) + # 'InversionTime' # Not known # 'MixingTime' # Not known # 'ExcitationFlipAngle' - fa = float(dcmdata.dcm_data.PerFrameFunctionalGroupsSequence[0] - .MRTimingAndRelatedParametersSequence[0].RepetitionTime) + try: + fa = float(dcmdata.dcm_data.PerFrameFunctionalGroupsSequence[0] + .MRTimingAndRelatedParametersSequence[0].RepetitionTime) + except AttributeError: + fa = dcmdata.dcm_data.FlipAngle obj.set_standard_def('ExcitationFlipAngle', fa) + # 'TxOffset' # Not known # 'VOI' @@ -417,7 +445,8 @@ def set_if_present(tag, val): # # 5.3 Sequence information # 'SequenceName' - obj.set_standard_def('SequenceName', dcmdata.dcm_data.PulseSequenceName) + if 'SequenceName' in dcmdata.dcm_data: + obj.set_standard_def('SequenceName', dcmdata.dcm_data.PulseSequenceName) # 'ProtocolName' obj.set_standard_def('ProtocolName', dcmdata.dcm_data.ProtocolName) # # 5.4 Sequence information @@ -461,8 +490,9 @@ def _extractDicomMetadata_new(dcmdata, water_suppressed=True): """ # Extract required metadata and create hdr_ext object - obj = nifti_mrs.hdr_ext(dcmdata.dcm_data.TransmitterFrequency, - dcmdata.dcm_data.ResonantNucleus) + obj = Hdr_Ext( + dcmdata.dcm_data.TransmitterFrequency, + dcmdata.dcm_data.ResonantNucleus) def set_if_present(tag, val): if val: diff --git a/spec2nii/Siemens/dicomfunctions.py b/spec2nii/Siemens/dicomfunctions.py index b472827..3f4129d 100644 --- a/spec2nii/Siemens/dicomfunctions.py +++ b/spec2nii/Siemens/dicomfunctions.py @@ -4,15 +4,18 @@ """ from datetime import datetime import re +import warnings import numpy as np import nibabel.nicom.dicomwrappers -import warnings from nibabel.nicom import csareader as csar -from spec2nii.dcm2niiOrientation.orientationFuncs import dcm_to_nifti_orientation -from spec2nii import nifti_mrs from mapvbvd.read_twix_hdr import parse_buffer + +from nifti_mrs.create_nmrs import gen_nifti_mrs_hdr_ext +from nifti_mrs.hdr_ext import Hdr_Ext + +from spec2nii.dcm2niiOrientation.orientationFuncs import dcm_to_nifti_orientation from spec2nii import __version__ as spec2nii_ver @@ -189,11 +192,11 @@ def not_equal(lst): # Create NIFTI MRS object. try: - nifti_mrs_out.append(nifti_mrs.NIfTI_MRS(combined_data, or_used.Q44, dt_used, meta_used)) + nifti_mrs_out.append(gen_nifti_mrs_hdr_ext(combined_data, dt_used, meta_used, or_used.Q44, no_conj=True)) except np.linalg.LinAlgError: - warnings.warn("""The quaternion passes to NIfTI_MRS was singular. - Most likely your slice position is not well defined. I have set it to None.""") - nifti_mrs_out.append(nifti_mrs.NIfTI_MRS(combined_data, None, dt_used, meta_used)) + warnings.warn("""The quaternion passed to NIfTI_MRS was singular. + Most likely your slice position is not well defined. I have set it to default.""") + nifti_mrs_out.append(gen_nifti_mrs_hdr_ext(combined_data, dt_used, meta_used, no_conj=True)) # If there are any identical names then append an index seen = np.unique(fnames_out) @@ -381,8 +384,9 @@ def extractDicomMetadata_xa(dcmdata): dcm_hdrs1 = dcmdata.dcm_data.SharedFunctionalGroupsSequence[0] # Extract required metadata and create hdr_ext object - obj = nifti_mrs.hdr_ext(dcmdata.dcm_data.TransmitterFrequency, - dcmdata.dcm_data.ResonantNucleus) + obj = Hdr_Ext( + dcmdata.dcm_data.TransmitterFrequency, + dcmdata.dcm_data.ResonantNucleus) # Standard defined metadata def set_standard_def(nifti_mrs_key, location, key, cast=None): @@ -480,8 +484,9 @@ def extractDicomMetadata_vx(dcmdata): """ # Extract required metadata and create hdr_ext object - obj = nifti_mrs.hdr_ext(dcmdata.csa_header['tags']['ImagingFrequency']['items'][0], - dcmdata.csa_header['tags']['ImagedNucleus']['items'][0]) + obj = Hdr_Ext( + dcmdata.csa_header['tags']['ImagingFrequency']['items'][0], + dcmdata.csa_header['tags']['ImagedNucleus']['items'][0]) # Standard defined metadata def set_standard_def(nifti_mrs_key, location, key, cast=None): diff --git a/spec2nii/Siemens/rda.py b/spec2nii/Siemens/rda.py index ec8fece..a401983 100644 --- a/spec2nii/Siemens/rda.py +++ b/spec2nii/Siemens/rda.py @@ -8,7 +8,9 @@ import numpy as np -from spec2nii import nifti_mrs +from nifti_mrs.create_nmrs import gen_nifti_mrs_hdr_ext +from nifti_mrs.hdr_ext import Hdr_Ext + from spec2nii.dcm2niiOrientation.orientationFuncs import dcm_to_nifti_orientation from spec2nii import __version__ as spec2nii_ver @@ -102,7 +104,7 @@ def convert_rda(rda_path, fname_out, verbose): else: name = rda_path.stem - return [nifti_mrs.NIfTI_MRS(data_cmplx, currNiftiOrientation.Q44, dwelltime, meta), ], [name, ] + return [gen_nifti_mrs_hdr_ext(data_cmplx, dwelltime, meta, currNiftiOrientation.Q44, no_conj=True), ], [name, ] def extractRdaMetadata(hdr): @@ -111,12 +113,13 @@ def extractRdaMetadata(hdr): :param hdr: Header read from rda file :type hdr: dict :return: Complete nifti_mrs header object - :rtype: spec2nii.nifti_mrs.hdr_ext + :rtype: nifti_mrs.nifti_mrs.Hdr_Ext """ # Extract required metadata and create hdr_ext object - obj = nifti_mrs.hdr_ext(_locale_float(hdr['MRFrequency']), - hdr['Nucleus']) + obj = Hdr_Ext( + _locale_float(hdr['MRFrequency']), + hdr['Nucleus']) # Standard defined metadata def set_standard_def(nifti_mrs_key, location, key, cast=None): diff --git a/spec2nii/Siemens/twixfunctions.py b/spec2nii/Siemens/twixfunctions.py index c77ac75..b5d3e76 100644 --- a/spec2nii/Siemens/twixfunctions.py +++ b/spec2nii/Siemens/twixfunctions.py @@ -7,10 +7,11 @@ from datetime import datetime import numpy as np +from nifti_mrs.create_nmrs import gen_nifti_mrs_hdr_ext +from nifti_mrs.hdr_ext import Hdr_Ext import spec2nii.GSL.gslfunctions as GSL from spec2nii.dcm2niiOrientation.orientationFuncs import dcm_to_nifti_orientation -from spec2nii import nifti_mrs from spec2nii import __version__ as spec2nii_ver @@ -527,7 +528,7 @@ def assemble_nifti_mrs(data, dwellTime, orientation, meta_obj, dim_tags=None): for idx, dt in zip(range(data.ndim - 4), dim_tags): meta_obj.set_dim_info(idx, dt) - return nifti_mrs.NIfTI_MRS(data, orientation.Q44, dwellTime, meta_obj) + return gen_nifti_mrs_hdr_ext(data, dwellTime, meta_obj, orientation.Q44, no_conj=True) def twix2DCMOrientation(mapVBVDHdr, force_svs=False, verbose=False): @@ -690,8 +691,9 @@ def extractTwixMetadata_xa(mapVBVDHdr, original_file): """ # Extract required metadata and create hdr_ext object - obj = nifti_mrs.hdr_ext(mapVBVDHdr['Dicom'][('lFrequency')] / 1E6, - mapVBVDHdr['MeasYaps'][('sTXSPEC', 'asNucleusInfo', '0', 'tNucleus')].strip('"')) + obj = Hdr_Ext( + mapVBVDHdr['Dicom'][('lFrequency')] / 1E6, + mapVBVDHdr['MeasYaps'][('sTXSPEC', 'asNucleusInfo', '0', 'tNucleus')].strip('"')) # Standard defined metadata # # 5.1 MRS specific Tags @@ -707,7 +709,9 @@ def extractTwixMetadata_xa(mapVBVDHdr, original_file): # 'ExcitationFlipAngle' obj.set_standard_def('ExcitationFlipAngle', float(mapVBVDHdr['MeasYaps'][('adFlipAngleDegree', '0')])) # 'TxOffset' - obj.set_standard_def('TxOffset', empty_str_to_0float(mapVBVDHdr['MeasYaps'][('sSpecPara', 'dDeltaFrequency')])) + obj.set_standard_def( + 'TxOffset', + empty_str_or_val_to_0float(mapVBVDHdr['MeasYaps'], ('sSpecPara', 'dDeltaFrequency'))) # 'VOI' # 'WaterSuppressed' # TO DO @@ -756,14 +760,15 @@ def extractTwixMetadata_xa(mapVBVDHdr, original_file): # 'PatientDoB' obj.set_standard_def('PatientDoB', str(mapVBVDHdr['Meas'][('PatientBirthDay')])) # 'PatientSex' - if mapVBVDHdr['Meas'][('PatientSex')] == 0.0: + patient_sex = mapVBVDHdr['Meas'][('PatientSex')] + if patient_sex == 1.0: sex_str = 'M' - elif mapVBVDHdr['Meas'][('PatientSex')] == 1.0: + elif patient_sex == 2.0: sex_str = 'F' - elif mapVBVDHdr['Meas'][('PatientSex')] == 2.0: + elif patient_sex == 3.0: sex_str = 'O' else: - raise ValueError('Meas, PatientSex, should be 0, 1, or 2.') + raise ValueError(f'Meas, PatientSex, should be 1, 2, or 3, but is {patient_sex}') obj.set_standard_def('PatientSex', sex_str) @@ -802,8 +807,9 @@ def extractTwixMetadata_vx(mapVBVDHdr, original_file): """ # Extract required metadata and create hdr_ext object - obj = nifti_mrs.hdr_ext(mapVBVDHdr['Meas'][('Frequency')] / 1E6, - mapVBVDHdr['Meas'][('ResonantNucleus')]) + obj = Hdr_Ext( + mapVBVDHdr['Meas'][('Frequency')] / 1E6, + mapVBVDHdr['Meas'][('ResonantNucleus')]) # Standard defined metadata # # 5.1 MRS specific Tags @@ -822,7 +828,7 @@ def extractTwixMetadata_vx(mapVBVDHdr, original_file): # 'ExcitationFlipAngle' obj.set_standard_def('ExcitationFlipAngle', float(mapVBVDHdr['Meas'][('FlipAngle')])) # 'TxOffset' - obj.set_standard_def('TxOffset', empty_str_to_0float(mapVBVDHdr['Meas'][('dDeltaFrequency')])) + obj.set_standard_def('TxOffset', empty_str_or_val_to_0float(mapVBVDHdr['Meas'], ('dDeltaFrequency'))) # 'VOI' # 'WaterSuppressed' # TO DO @@ -898,11 +904,14 @@ def extractTwixMetadata_vx(mapVBVDHdr, original_file): return obj -def empty_str_to_0float(value): - if value == '': +def empty_str_or_val_to_0float(param_dict, key): + try: + if param_dict[key] == '': + return 0.0 + else: + return param_dict[key] + except KeyError: return 0.0 - else: - return value def _try_int(value): diff --git a/spec2nii/anonymise.py b/spec2nii/anonymise.py index 65d66b4..194c1fc 100644 --- a/spec2nii/anonymise.py +++ b/spec2nii/anonymise.py @@ -4,14 +4,13 @@ Copyright (C) 2021 University of Oxford """ -import json import re import pprint -import nibabel as nib - -from spec2nii.nifti_mrs.definitions import standard_defined -from spec2nii import nifti_mrs +from nifti_mrs.definitions import standard_defined +from nifti_mrs.create_nmrs import gen_nifti_mrs_hdr_ext +from nifti_mrs.nifti_mrs import NIFTI_MRS +from nifti_mrs.hdr_ext import Hdr_Ext def anon_nifti_mrs(args): @@ -24,11 +23,10 @@ def anon_nifti_mrs(args): :rtype: [str,] """ # Load data - nifti_mrs_img = nib.load(args.file) + nifti_mrs_img = NIFTI_MRS(args.file) # Extract header extension - hdr_ext_codes = nifti_mrs_img.header.extensions.get_codes() - hdr_ext = json.loads(nifti_mrs_img.header.extensions[hdr_ext_codes.index(44)].get_content()) + hdr_ext = nifti_mrs_img.hdr_ext.to_dict() # Loop through fields. Remove those which are marked for anonymisation # Either those fields which are standard-defined and are marked for anonymisation, @@ -75,10 +73,11 @@ def iter_json(in_dict): pp.pprint(anon_hdr) # Make new NIfTI-MRS image - anon_out = nifti_mrs.NIfTI_MRS(nifti_mrs_img.get_fdata(dtype=nifti_mrs_img.get_data_dtype()), - nifti_mrs_img.affine, - nifti_mrs_img.header['pixdim'][4], - anon_hdr) + anon_out = gen_nifti_mrs_hdr_ext( + nifti_mrs_img[:], + nifti_mrs_img.dwelltime, + Hdr_Ext.from_header_ext(anon_hdr), + affine=nifti_mrs_img.getAffine('voxel', 'world')) # Process output name. if args.fileout: diff --git a/spec2nii/bruker.py b/spec2nii/bruker.py index 22d7cfa..9bf2d0e 100644 --- a/spec2nii/bruker.py +++ b/spec2nii/bruker.py @@ -18,8 +18,10 @@ from brukerapi.mergers import FrameGroupMerger from brukerapi.exceptions import FilterEvalFalse +from nifti_mrs.create_nmrs import gen_nifti_mrs_hdr_ext +from nifti_mrs.hdr_ext import Hdr_Ext + from spec2nii.nifti_orientation import NIFTIOrient -from spec2nii import nifti_mrs from spec2nii import __version__ as spec2nii_ver # Default dimension assignments. @@ -41,10 +43,12 @@ def read_bruker(args): # for all Bruker datasets compliant all queries for data, orientation, dwelltime, meta, name in yield_bruker(args): imageOut.append( - nifti_mrs.NIfTI_MRS(data, - orientation.Q44, - dwelltime, - meta) + gen_nifti_mrs_hdr_ext( + data, + dwelltime, + meta, + orientation.Q44, + no_conj=True) ) fileoutNames.append(name) @@ -222,8 +226,9 @@ def _2dseq_meta(d, dump=False): # Extract required metadata and create hdr_ext object cf = d.SpectrometerFrequency - obj = nifti_mrs.hdr_ext(cf, - d.ResonantNucleus) + obj = Hdr_Ext( + cf, + d.ResonantNucleus) # # 5.1 MRS specific Tags # 'EchoTime' @@ -310,8 +315,9 @@ def _fid_meta(d, dump=False): # Extract required metadata and create hdr_ext object cf = d.SpectrometerFrequency - obj = nifti_mrs.hdr_ext(cf, - d.ResonantNucleus) + obj = Hdr_Ext( + cf, + d.ResonantNucleus) # # 5.1 MRS specific Tags # 'EchoTime' diff --git a/spec2nii/jmrui.py b/spec2nii/jmrui.py index 9f651d4..498d6fe 100644 --- a/spec2nii/jmrui.py +++ b/spec2nii/jmrui.py @@ -2,11 +2,15 @@ Author: William Clarke Copyright (C) 2020 University of Oxford """ -import numpy as np import re -from spec2nii.nifti_orientation import NIFTIOrient -from spec2nii import nifti_mrs from datetime import datetime + +import numpy as np + +from nifti_mrs.create_nmrs import gen_nifti_mrs_hdr_ext +from nifti_mrs.hdr_ext import Hdr_Ext + +from spec2nii.nifti_orientation import NIFTIOrient from spec2nii import __version__ as spec2nii_ver @@ -85,10 +89,12 @@ def jmrui_mrui(args): nifti_orientation = NIFTIOrient(affine) - img_out = [nifti_mrs.NIfTI_MRS(data, - nifti_orientation.Q44, - dwelltime, - meta), ] + img_out = [gen_nifti_mrs_hdr_ext( + data, + dwelltime, + meta, + nifti_orientation.Q44, + no_conj=True), ] # File names if args.fileout: @@ -162,8 +168,9 @@ def jmrui_hdr_to_obj_mrui(header, str_info): header['transmitter_frequency'], header['magnetic_field']) - meta = nifti_mrs.hdr_ext(header['transmitter_frequency'], - nucleus) + meta = Hdr_Ext( + header['transmitter_frequency'], + nucleus) # meta.set_standard_def('ManufacturersModelName', header['Spectrometer']) # meta.set_standard_def('PatientName', header['NameOfPatient']) @@ -206,10 +213,12 @@ def jmrui_txt(args): nifti_orientation = NIFTIOrient(affine) - img_out = [nifti_mrs.NIfTI_MRS(data, - nifti_orientation.Q44, - dwelltime, - meta), ] + img_out = [gen_nifti_mrs_hdr_ext( + data, + dwelltime, + meta, + nifti_orientation.Q44, + no_conj=True), ] # File names if args.fileout: @@ -310,8 +319,9 @@ def jmrui_hdr_to_obj(header): header['TransmitterFrequency'], header['MagneticField']) - meta = nifti_mrs.hdr_ext(float(header['TransmitterFrequency']), - nucleus) + meta = Hdr_Ext( + float(header['TransmitterFrequency']), + nucleus) if 'Spectrometer' in header: meta.set_standard_def('ManufacturersModelName', header['Spectrometer']) diff --git a/spec2nii/nifti_mrs/__init__.py b/spec2nii/nifti_mrs/__init__.py deleted file mode 100644 index 4afc882..0000000 --- a/spec2nii/nifti_mrs/__init__.py +++ /dev/null @@ -1,2 +0,0 @@ -from spec2nii.nifti_mrs.hdr_ext import hdr_ext -from spec2nii.nifti_mrs.nifti_mrs import get_mrs_class diff --git a/spec2nii/nifti_mrs/definitions.py b/spec2nii/nifti_mrs/definitions.py deleted file mode 100644 index 148bd7f..0000000 --- a/spec2nii/nifti_mrs/definitions.py +++ /dev/null @@ -1,259 +0,0 @@ -'''Definitions of NIfTI-MRS standard meta data and dimension tags. - -Type fields should either be generic python types: float, int, str -or a tuple indicating an array type and element type : (list, float) or (list, str) - -Copyright Will Clarke, University of Oxford, 2021 -''' - -# Define nifti-mrs version number here. -# First element is major version, second is minor -nifti_mrs_version = [0, 6] - -# Possible dimension tags and descriptions -dimension_tags = {"DIM_COIL": "For storage of data from each individual receiver coil element.", - "DIM_DYN": "For storage of each individual acquisition transient. E.g. for post-acquisition B0 drift correction.", - "DIM_INDIRECT_0": "The indirect detection dimension - necessary for 2D (and greater) MRS acquisitions.", - "DIM_INDIRECT_1": "The indirect detection dimension - necessary for 2D (and greater) MRS acquisitions.", - "DIM_INDIRECT_2": "The indirect detection dimension - necessary for 2D (and greater) MRS acquisitions.", - "DIM_PHASE_CYCLE": "Used for increments of phase-cycling, for example in dephasing unwanted coherence order pathways, or TPPI for 2D spectra.", - "DIM_EDIT": "Used for edited MRS techniques such as MEGA or HERMES.", - "DIM_MEAS": "Used to indicate multiple repeats of the full sequence contained within the same original data file.", - "DIM_USER_0": "User defined dimension.", - "DIM_USER_1": "User defined dimension.", - "DIM_USER_2": "User defined dimension.", - "DIM_ISIS": "Dimension for storing image-selected in vivo spectroscopy (ISIS) acquisitions."} - -# Required metadata fields -required = {'SpectrometerFrequency': - {'doc': 'Precession frequency in MHz of the nucleus being addressed for each spectral axis.', - 'type': (list, float)}, - 'ResonantNucleus': - {'doc': 'Must be one of the DICOM recognised nuclei “1H”, “3HE”, “7LI”, “13C”, “19F”, “23NA”, “31P”, “129XE” or one named in the specified format. I.e. Mass number followed by the chemical symbol in uppercase.', - 'type': (list, str)}} - -# Defined metadata fields -# # 5.1 MRS specific Tags -# 'EchoTime' -# 'RepetitionTime' -# 'InversionTime' -# 'MixingTime' -# 'AcquisitionStartTime' -# 'ExcitationFlipAngle' -# 'TxOffset' -# 'VOI' -# 'WaterSuppressed' -# 'WaterSuppressionType' -# 'SequenceTriggered' -# # 5.2 Scanner information -# 'Manufacturer' -# 'ManufacturersModelName' -# 'DeviceSerialNumber' -# 'SoftwareVersions' -# 'InstitutionName' -# 'InstitutionAddress' -# 'TxCoil' -# 'RxCoil' -# # 5.3 Sequence information -# 'SequenceName' -# 'ProtocolName' -# # 5.4 Sequence information -# 'PatientPosition' -# 'PatientName' -# 'PatientID' -# 'PatientWeight' -# 'PatientDoB' -# 'PatientSex' -# # 5.5 Provenance and conversion metadata -# 'ConversionMethod' -# 'ConversionTime' -# 'OriginalFile' -# # 5.6 Spatial information -# 'kSpace' -# # 5.7 Editing Pulse information structure -# 'EditCondition' -# 'EditPulse' -# # 5.8 Processing Provenance -# 'ProcessingApplied' - -# These fields are optional but must not be redefined. -# Format is a dict of tuples containing (type, unit string, doc string, anonymisation state) -standard_defined = { - # 5.1 MRS specific Tags - 'EchoTime': - (float, - 's', - 'Time from centroid of excitation to start of FID or centre of echo. Units: Seconds', - False), - 'RepetitionTime': - (float, - 's', - 'Sequence repetition time. Units: Seconds', - False), - 'InversionTime': - (float, - 's', - 'Inversion time. Units: Seconds', - False), - 'MixingTime': - (float, - 's', - 'Mixing time in e.g. STEAM sequence. Units: Seconds', - False), - 'AcquisitionStartTime': - (float, - 's', - 'Time, relative to EchoTime, that the acquisition starts. Positive values indicate a time after the EchoTime, negative indicate before the EchoTime, a value of zero indicates no offset. Units: Seconds', - False), - 'ExcitationFlipAngle': - (float, - 'degrees', - 'Nominal excitation pulse flip-angle', - False), - 'TxOffset': - (float, - 'ppm', - 'Transmit chemical shift offset from SpectrometerFrequency', - False), - 'VOI': - ((list, list, float), - None, - 'VoI localisation volume for MRSI sequences. Stored as a 4 x 4 affine using identical conventions to the xform NIfTI affine matrix. Not defined for data stored with a single spatial voxel', - False), - 'WaterSuppressed': - (bool, - None, - 'Boolian value indicating whether data was collected with (True) or without (False) water suppression.', - False), - 'WaterSuppressionType': - (str, - None, - 'Type of water suppression used.', - False), - 'SequenceTriggered': - (bool, - None, - 'Boolian value indicating whether the sequence is triggered. If triggered the repetition time might not be constant.', - False), - # 5.2 Scanner information - 'Manufacturer': - (str, - None, - 'Manufacturer of the device. DICOM tag (0008,0070).', - False), - 'ManufacturersModelName': - (str, - None, - "Manufacturer's model name of the device. DICOM tag (0008,1090).", - False), - 'DeviceSerialNumber': - (str, - None, - "Manufacturer's serial number of the device. DICOM tag (0018,1000).", - False), - 'SoftwareVersions': - (str, - None, - "Manufacturer's designation of the software version. DICOM tag (0018,1020)", - False), - 'InstitutionName': - (str, - None, - "Institution’s Name. DICOM tag (0008,0080).", - False), - 'InstitutionAddress': - (str, - None, - "Institution’s address. DICOM tag (0008,0081).", - False), - 'TxCoil': - (str, - None, - "Name or description of transmit RF coil.", - False), - 'RxCoil': - (str, - None, - "Name or description of receive RF coil.", - False), - # 5.3 Sequence information - 'SequenceName': - (str, - None, - "User defined name. DICOM tag (0018,0024).", - False), - 'ProtocolName': - (str, - None, - "User-defined description of the conditions under which the Series was performed. DICOM tag (0018,1030).", - False), - # 5.4 Sequence information - 'PatientPosition': - (str, - None, - "Patient position descriptor relative to the equipment. DICOM tag (0018,5100). Must be one of the DICOM defined code strings e.g. HFS, HFP.", - False), - 'PatientName': - (str, - None, - "Patient's full name. DICOM tag (0010,0010).", - True), - 'PatientID': - (str, - None, - "Patient identifier. DICOM tag (0010,0020).", - True), - 'PatientWeight': - (float, - 'kg', - "Weight of the Patient in kilograms. DICOM tag (0010,1030).", - False), - 'PatientDoB': - (str, - None, - "Date of birth of the named Patient. YYYYMMDD. DICOM tag (0010,0030).", - True), - 'PatientSex': - (str, - None, - "Sex of the named Patient. ‘M’, ‘F’, ‘O’. DICOM tag (0010,0040)", - False), - # 5.5 Provenance and conversion metadata - 'ConversionMethod': - (str, - None, - "Description of the process or program used for conversion. May include additional information like software version.", - False), - 'ConversionTime': - (str, - None, - "Time and date of conversion. ISO 8601 compliant format", - False), - 'OriginalFile': - ((list, str), - None, - "Name and extension of the original file(s)", - True), - # 5.6 Spatial information - 'kSpace': - ((list, bool), - None, - "Three element list, corresponding to the first three spatial dimensions. If True the data is stored as a dense k-space representation.", - False), - # 5.7 Editing Pulse information structure - 'EditConditon': - ((list, str), - None, - "List of strings that index the entries of the EditPulse structure that are used in this data acquisition. Typically used in dynamic headers (dim_N_header).", - False), - 'EditPulse': - (dict, - None, - "Structure defining editing pulse parameters for each condition. Each condition must be assigned a key.", - False), - # 5.8 Processing Provenance - 'ProcessingApplied': - (dict, - None, - "Describes and records the processing steps applied to the data.", - False)} diff --git a/spec2nii/nifti_mrs/hdr_ext.py b/spec2nii/nifti_mrs/hdr_ext.py deleted file mode 100644 index a8e3ccc..0000000 --- a/spec2nii/nifti_mrs/hdr_ext.py +++ /dev/null @@ -1,127 +0,0 @@ -from .definitions import dimension_tags, standard_defined -import json - - -class hdr_ext: - """Class to hold meta data stored in a NIfTI MRS header extension. - Required fields must be passed to initialise, - Default dimension information automatically generated, but may be modified by set_dim_info method. - Standard defined meta-data and user-defined data can be added using set_standard_def and - set_user_def respectively. - """ - def __init__(self, spec_frequency, resonant_nucleus): - """Initialise class object with the two required bits of meta-data - Set default dimension information. - Inputs: - spec_frequency: Spectrometer frequency in MHz, - resonant_nucleus: Resonant nucleus string e.g. '1H' - """ - if isinstance(spec_frequency, float): - self.SpectrometerFrequency = [spec_frequency, ] - elif isinstance(spec_frequency, (list, tuple))\ - and isinstance(spec_frequency[0], float): - self.SpectrometerFrequency = spec_frequency - else: - raise ValueError('spec_frequency must be a float or array of floats.') - - if isinstance(resonant_nucleus, str): - self.ResonantNucleus = [resonant_nucleus, ] - elif isinstance(resonant_nucleus, (list, tuple))\ - and isinstance(resonant_nucleus[0], str): - self.ResonantNucleus = resonant_nucleus - else: - raise ValueError('resonant_nucleus must be a string or array of strings.') - - self.dim_info = [{"tag": "DIM_COIL", "info": None, "hdr": None}, - {"tag": "DIM_DYN", "info": None, "hdr": None}, - {"tag": "DIM_INDIRECT_0", "info": None, "hdr": None}] - - self.standard_data = {} - self.user_data = {} - - def set_dim_info(self, dim, tag, info=None, hdr=None): - """Set information associated with the optional data dimensions. - Inputs: - dim: May be (0,1,2) or ("5th","6th","7th") - tag: Must be one of the defined dimension tag strings. - info: Optional, free-form use string. - hdr: Dict containing relevant header value names and values. - """ - if tag not in dimension_tags: - raise ValueError("tag must be one of the defined dimension tag.") - - new_info = {"tag": tag, - "info": info, - "hdr": hdr} - - if dim in (0, 1, 2): - self.dim_info[dim] = new_info - elif dim in ("5th", "6th", "7th"): - if dim == "5th": - self.dim_info[0] = new_info - elif dim == "6th": - self.dim_info[1] = new_info - elif dim == "7th": - self.dim_info[2] = new_info - else: - raise ValueError('dim must be 0,1,2 or "5th","6th","7th".') - - def set_standard_def(self, key, value): - """Add a single standard-defined bit of meta-data to the object.""" - if key not in standard_defined: - raise ValueError("key must be one of the standard-defined keys.") - - self.standard_data[key] = value - - def set_user_def(self, all_keys=None, key=None, value=None, doc=None): - """Add user defined meta data keys to the header extension. - Pass dict as kwarg all_keys to set all key/value pairs, or - add keys and values one at a time using key, value and doc. - """ - - if all_keys is not None: - self.user_data = all_keys - else: - if isinstance(value, dict): - self.user_data[key] = value - self.user_data[key].update({'Description': doc}) - else: - self.user_data[key] = {'Value': value, - 'Description': doc} - - def get_json(self, dimensions=7): - """Generate json from properties.""" - - # Required meta-data - json_obj = {'SpectrometerFrequency': self.SpectrometerFrequency, - 'ResonantNucleus': self.ResonantNucleus} - - # Dimension information - if dimensions < 4: - raise ValueError('dimensions must be 4 or greater') - elif dimensions == 4: - pass - else: - update_dict = {} - for idx in range(5, dimensions + 1): - update_dict[f'dim_{idx}'] = self.dim_info[idx - 5]['tag'] - if self.dim_info[idx - 5]['info'] is not None: - update_dict[f'dim_{idx}_info'] = self.dim_info[idx - 5]['info'] - if self.dim_info[idx - 5]['hdr'] is not None: - update_dict[f'dim_{idx}_header'] = self.dim_info[idx - 5]['hdr'] - - json_obj.update(update_dict) - - # Add standard defined - json_obj.update(self.standard_data) - - # Add user defined - json_obj.update(self.user_data) - - return json.dumps(json_obj) - - def __str__(self) -> str: - return self.get_json() - - def __repr__(self) -> str: - return str(self) diff --git a/spec2nii/nifti_mrs/nifti_mrs.py b/spec2nii/nifti_mrs/nifti_mrs.py deleted file mode 100644 index e5fc00b..0000000 --- a/spec2nii/nifti_mrs/nifti_mrs.py +++ /dev/null @@ -1,80 +0,0 @@ -import json - -import nibabel as nib -from .validator import validate_nifti_mrs -from spec2nii.nifti_mrs import hdr_ext -from spec2nii.nifti_mrs.definitions import nifti_mrs_version - - -def get_mrs_class(nifti=2): - '''Generate the NIfTI_MRS class using either NIFTI2 (default) or NIFTI1 formatting.''' - - base_class = nib.nifti1.Nifti1Image if nifti == 1 else nib.nifti2.Nifti2Image - - class NIfTI_MRS(base_class): - """Class to contain a NIfTI MRS dataset. Derived from nibabel's Nifti2Image.""" - - def __init__(self, dataobj, affine, dwelltime, header_ext, header=None, extra=None, file_map=None): - """Initialise a NIfTI MRS object which can be saved to disk. - :param dataobj: Four to seven dimensional complex numpy array. - :param affine: None or (4,4) array-like - homogeneous affine giving relationship between voxel coordinates and - world coordinates. Affine can also be None. In this case, - ``obj.affine`` also returns None, and the affine as written to disk - will depend on the file format. - :param float dwelltime: Spectroscopic time-domain dwelltime in seconds. - :param header_ext : nifti_mrs.hdr_ext object or dict. - :param header: None or mapping or header instance, optional - metadata for this image format - :param extra: None or mapping, optional - metadata to associate with image that cannot be stored in the - metadata of this image type - :param file_map : mapping, optional - mapping giving file information for this image format - """ - super().__init__(dataobj, affine, header=header, extra=extra, file_map=file_map) - - # Add the header extensions - if isinstance(header_ext, hdr_ext): - json_s = header_ext.get_json(dimensions=dataobj.ndim) - elif isinstance(header_ext, dict): - json_s = json.dumps(header_ext) - else: - raise ValueError('header_ext must be a spec2nii.nifti_mrs.hdr_ext object or dict') - - extension = nib.nifti1.Nifti1Extension(44, json_s.encode('UTF-8')) - self.header.extensions.append(extension) - - # Set intent_name - self.set_version_info( - nifti_mrs_version[0], - nifti_mrs_version[1]) - - # Set the dwell time of the time dimension - self.set_dwell_time(dwelltime) - - def set_dwell_time(self, dwelltime): - """Set dwelltime (pixdim[4])""" - pixDim = self.header['pixdim'] - pixDim[4] = dwelltime - self.header['pixdim'] = pixDim - - def set_version_info(self, major, minor): - """puts mrs_v{major}_{minor} into intent_name""" - self.header['intent_name'] = f'mrs_v{major}_{minor}'.encode() - - def validate(self): - """Run NIfTI MRS validation.""" - validate_nifti_mrs(self) - - def save(self, filename): - """ Utility function """ - nib.save(self, filename) - - @property - def hdr_ext(self): - '''Return MRS JSON header extension as python dict''' - hdr_ext_codes = self.header.extensions.get_codes() - return json.loads(self.header.extensions[hdr_ext_codes.index(44)].get_content()) - - return NIfTI_MRS diff --git a/spec2nii/nifti_mrs/validator.py b/spec2nii/nifti_mrs/validator.py deleted file mode 100644 index 26f0730..0000000 --- a/spec2nii/nifti_mrs/validator.py +++ /dev/null @@ -1,146 +0,0 @@ -import json -from .definitions import dimension_tags, standard_defined -from numpy import asanyarray, iscomplexobj - - -class Error(Exception): - """Base class for other exceptions""" - pass - - -class headerExtensionError(Error): - """Raised if problems with header extension are found.""" - pass - - -class niftiHeaderError(Error): - """Raised if problems with nifti header are found.""" - pass - - -class niftiDataError(Error): - """Raised if problems with nifti data are found.""" - pass - - -def validate_nifti_mrs(nifti_mrs): - """Validate a full NIfTI MRS image.""" - - # Validate data - data = asanyarray(nifti_mrs.dataobj) - validate_nifti_data(data) - - # Validate nifit header - validate_nifti_header(nifti_mrs.header) - - # Validate header extension - validate_hdr_ext(nifti_mrs.header.extensions[0].get_content(), - data.ndim) - - -def validate_nifti_data(nifti_img_data): - """Validate the data inside a nibabel nifti image - 1. Check data is complex - 2. Check number of dimensions is at least 4 but less than 8. - """ - - # 1. Check for complexity - if not iscomplexobj(nifti_img_data): - raise niftiDataError('Data must be complex.') - - # 2. Check for between 4 and 7 dimensions - if nifti_img_data.ndim < 4\ - or nifti_img_data.ndim > 7: - raise niftiDataError('Data must have between 4 and 7 dimensions.' - f' It has {nifti_img_data.ndim}.') - - -def validate_nifti_header(nifti_header): - """Validate the header of a nibabel nifti image - Check data type is complex - Check orientation data. - Check dwell time - Check intent name - """ - pass - - -def validate_hdr_ext(header_ex, data_dimensions): - """ Validate the header extension - 1. Check that it is json formatted string. - 2. Check that it contains the required meta-data - 3. Check that it contains any required dimension information. - 4. Check that standard-defined data is of correct type. - Inputs: - - """ - # 1. Check that header_ext is json - try: - json_dict = json.loads(header_ex) - except json.JSONDecodeError as exc: - raise headerExtensionError("Header extension is not json deserialisable.") from exc - - # 2. Check the two required bits of meta-data - if "SpectrometerFrequency" in json_dict: - if not isinstance(json_dict["SpectrometerFrequency"], (list, tuple))\ - and not isinstance(json_dict["SpectrometerFrequency"][0], float): - - raise headerExtensionError("SpectrometerFrequency must be list of floats.") - else: - raise headerExtensionError("Header extension must contain SpectrometerFrequency.") - - if "ResonantNucleus" in json_dict: - if not isinstance(json_dict["ResonantNucleus"], (list, tuple))\ - and not isinstance(json_dict["ResonantNucleus"][0], str): - - raise headerExtensionError("ResonantNucleus must be list of strings.") - else: - raise headerExtensionError("Header extension must contain ResonantNucleus.") - - # 3. Dimension information - for ddx in range(5, 8): - if data_dimensions > (ddx - 1): - if f"dim_{ddx}" in json_dict: - if json_dict[f"dim_{ddx}"] not in dimension_tags: - raise headerExtensionError(f"'dim_{ddx}' must be a defined tag.") - - if f"dim_{ddx}_info" in json_dict\ - and not isinstance(json_dict[f"dim_{ddx}_info"], str): - raise headerExtensionError(f"'dim_{ddx}_info' must be a string.") - - if f"dim_{ddx}_header" in json_dict\ - and not isinstance(json_dict[f"dim_{ddx}_header"], dict): - raise headerExtensionError(f"'dim_{ddx}_header' must be a dict.") - else: - raise headerExtensionError(f" With {data_dimensions} dimensions the header extension" - f" must contain 'dim_{ddx}'.") - - # 4. Check standard-defined data types - for key in json_dict: - if key in standard_defined and not check_type(json_dict[key], standard_defined[key][0]): - raise headerExtensionError(f'{key} must be a {standard_defined[key][0]}. ' - f'{key} is a {type(json_dict[key])}, with value {json_dict[key]}.') - - print('Header extension validated!') - - -def check_type(value, json_type): - '''Checks that values is of type json_type - json_type may be a tuple to handle array types - e.g. (list, float) indicates a list of floats. - ''' - if isinstance(json_type, tuple): - while len(json_type) > 1: - if not check_type(value, json_type[0]): - return False - try: - # TO DO: check more than the first element! - value = value[0] - except TypeError: - return False - json_type = json_type[1:] - return check_type(value, json_type[0]) - else: - if isinstance(value, json_type): - return True - return False diff --git a/spec2nii/other.py b/spec2nii/other.py index afd3cbd..7861943 100644 --- a/spec2nii/other.py +++ b/spec2nii/other.py @@ -6,9 +6,7 @@ import json import pprint -import nibabel as nib - -from spec2nii import nifti_mrs +from nifti_mrs.nifti_mrs import NIFTI_MRS def dump_headers(args): @@ -17,18 +15,15 @@ def dump_headers(args): :param args: Command line arguments """ # Load file - nifti_mrs_img = nib.load(args.file) + nifti_mrs_img = NIFTI_MRS(args.file) # Print NIfTI header print('NIfTI header: ') print(nifti_mrs_img.header) - hdr_ext_codes = nifti_mrs_img.header.extensions.get_codes() - hdr_ext = json.loads(nifti_mrs_img.header.extensions[hdr_ext_codes.index(44)].get_content()) - print('\nNIfTI-MRS header extension: ') pp = pprint.PrettyPrinter() - pp.pprint(hdr_ext) + pp.pprint(nifti_mrs_img.hdr_ext.to_dict()) def extract_hdr_ext(args): @@ -36,9 +31,7 @@ def extract_hdr_ext(args): :param args: Command line arguments """ - nifti_mrs_img = nib.load(args.file) - - hdr_ext_codes = nifti_mrs_img.header.extensions.get_codes() + nifti_mrs_img = NIFTI_MRS(args.file) if args.outdir: args.outdir.mkdir(exist_ok=True, parents=True) @@ -47,29 +40,25 @@ def extract_hdr_ext(args): out_json = args.file.parent / (args.file.with_suffix('').with_suffix('').name + '.json') with open(out_json, 'w') as fp: - json.dump(json.loads(nifti_mrs_img.header.extensions[hdr_ext_codes.index(44)].get_content()), fp, indent=4) + json.dump(nifti_mrs_img.hdr_ext.to_dict(), fp, indent=4) def insert_hdr_ext(args): - """Function for anonymising input NIfTI-MRS files. + """Function for inserting a new header into NIfTI-MRS files. :param args: Command line arguments parsed in spec2nii.py - :return: List of anonymised images - :rtype: [nib.nifti2.Nifti2Image,] + :return: Modified NIfTI-MRS file + :rtype: [nifti_mrs.nifti_mrs.NIFTI_MRS,] :return: List of output names :rtype: [str,] """ # Load data - nifti_mrs_img = nib.load(args.file) + nifti_mrs_img = NIFTI_MRS(args.file) with open(args.json_file) as jf: new_hdr = json.load(jf) - # Make new NIfTI-MRS image - mod_out = nifti_mrs.NIfTI_MRS(nifti_mrs_img.get_fdata(dtype=nifti_mrs_img.get_data_dtype()), - nifti_mrs_img.affine, - nifti_mrs_img.header['pixdim'][4], - new_hdr) + nifti_mrs_img.hdr_ext = new_hdr # Process output name. if args.fileout: @@ -77,4 +66,4 @@ def insert_hdr_ext(args): else: fname_out = [args.file.with_suffix('').with_suffix('').name, ] - return [mod_out, ], fname_out + return [nifti_mrs_img, ], fname_out diff --git a/spec2nii/other_formats.py b/spec2nii/other_formats.py index 4b79649..1c2c8e4 100644 --- a/spec2nii/other_formats.py +++ b/spec2nii/other_formats.py @@ -3,12 +3,15 @@ Author: William Clarke Copyright (C) 2020 University of Oxford """ -import numpy as np -from spec2nii.nifti_orientation import NIFTIOrient -from spec2nii import nifti_mrs from datetime import datetime from os.path import basename, splitext + +import numpy as np +from nifti_mrs.create_nmrs import gen_nifti_mrs_hdr_ext +from nifti_mrs.hdr_ext import Hdr_Ext + from spec2nii import __version__ as spec2nii_ver +from spec2nii.nifti_orientation import NIFTIOrient class InsufficentHeaderInformationError(Exception): @@ -27,8 +30,9 @@ def text(args): # Interpret required arguments (frequency and bandwidth) dwelltime = 1.0 / args.bandwidth - meta = nifti_mrs.hdr_ext(args.imagingfreq, - args.nucleus) + meta = Hdr_Ext( + args.imagingfreq, + args.nucleus) meta.set_standard_def('ConversionMethod', 'spec2nii') conversion_time = datetime.now().isoformat(sep='T', timespec='milliseconds') @@ -44,10 +48,12 @@ def text(args): nifti_orientation = NIFTIOrient(affine) - img_out = [nifti_mrs.NIfTI_MRS(data, - nifti_orientation.Q44, - dwelltime, - meta), ] + img_out = [gen_nifti_mrs_hdr_ext( + data, + dwelltime, + meta, + nifti_orientation.Q44, + no_conj=True), ] # File names if args.fileout: @@ -89,8 +95,9 @@ def lcm_raw(args): "There is no 'centralFrequency' key in the header information, " "please use the optional '-i' argument to provide a frequency.") - meta = nifti_mrs.hdr_ext(spec_frequency, - args.nucleus) + meta = Hdr_Ext( + spec_frequency, + args.nucleus) meta.set_standard_def('ConversionMethod', f'spec2nii v{spec2nii_ver}') conversion_time = datetime.now().isoformat(sep='T', timespec='milliseconds') @@ -106,10 +113,12 @@ def lcm_raw(args): nifti_orientation = NIFTIOrient(affine) - img_out = [nifti_mrs.NIfTI_MRS(data, - nifti_orientation.Q44, - dwelltime, - meta), ] + img_out = [gen_nifti_mrs_hdr_ext( + data, + dwelltime, + meta, + nifti_orientation.Q44, + no_conj=True), ] # File names if args.fileout: diff --git a/spec2nii/spec2nii.py b/spec2nii/spec2nii.py index f777eb9..296bdde 100644 --- a/spec2nii/spec2nii.py +++ b/spec2nii/spec2nii.py @@ -28,7 +28,6 @@ import json from nibabel.nifti2 import Nifti2Image -from spec2nii import nifti_mrs from spec2nii import __version__ as spec2nii_ver # There are case specific imports below @@ -121,6 +120,11 @@ def add_common_parameters(subparser): nargs='+', default=None, help="Specify the shape of higher (5th-7th) dimensions. Applies numpy style reshaping.") + parser_philips.add_argument( + "--special", + type=str, + default=None, + help="Identify special case sequence. Options: 'hyper'.") parser_philips = add_common_parameters(parser_philips) parser_philips.set_defaults(func=self.philips) @@ -128,10 +132,15 @@ def add_common_parameters(subparser): parser_p_dl = subparsers.add_parser('philips_dl', help='Convert from Philips data/list format.') parser_p_dl.add_argument('data', help='.data file', type=Path) parser_p_dl.add_argument('list', help='.list file', type=Path) - parser_p_dl.add_argument('spar', help='.SPAR file', type=Path) + parser_p_dl.add_argument('aux', help='Auxiliary file: .SPAR or DICOM file', type=Path) # for idx in range(5, 8): # parser_p_dl.add_argument(f"-d{idx}", f"--dim{idx}", type=str, help=f"Specify dim {idx} loop counter.") # parser_p_dl.add_argument(f"-t{idx}", f"--tag{idx}", type=str, help=f"Specify dim {idx} NIfTI MRS tag.") + parser_p_dl.add_argument( + "--special", + type=str, + default=None, + help="Identify special case sequence. Options: 'hyper'.") parser_p_dl = add_common_parameters(parser_p_dl) parser_p_dl.set_defaults(func=self.philips_dl) @@ -278,17 +287,12 @@ def add_common_parameters(subparser): self.outputDir = args.outdir - if args.nifti1: - nifti_mrs.NIfTI_MRS = nifti_mrs.get_mrs_class(nifti=1) - else: - nifti_mrs.NIfTI_MRS = nifti_mrs.get_mrs_class(nifti=2) - args.func(args) if self.imageOut: self.implement_overrides(args) self.validate_output() - self.write_output(args.json) + self.write_output(args.json, args.nifti1) self.validate_write(args.verbose) if args.verbose: print(f'Please cite {cite_str}.') @@ -318,10 +322,20 @@ def implement_overrides(self, args): def validate_output(self): """Run NIfTI MRS validation on output.""" + import nifti_mrs.validator as validate + # Currently this repeats the validation before the save. + # But useful here to do exception handling for f_out, nifti_mrs_img in zip(self.fileoutNames, self.imageOut): - nifti_mrs_img.validate() + try: + validate.validate_nifti_mrs(nifti_mrs_img) + except ( + validate.headerExtensionError, + validate.niftiDataError, + validate.niftiHeaderError) as exc: - def write_output(self, write_json=False): + raise Spec2niiError(f'Generated file {f_out} failed validation.') from exc + + def write_output(self, write_json=False, nifti1=False): """Write any NIfTI MRS objects stored. If write_json is true also write meta-data as sidecar. """ @@ -329,7 +343,20 @@ def write_output(self, write_json=False): for f_out, nifti_mrs_img in zip(self.fileoutNames, self.imageOut): out = self.outputDir / (f_out + '.nii.gz') - nifti_mrs_img.save(out) + + if nifti1: + # If nifti1 is requested just remake the file. + # This is more than a bit hacky, but avoids passing the option to every end-point. + from nifti_mrs.create_nmrs import gen_nifti_mrs_hdr_ext + gen_nifti_mrs_hdr_ext( + nifti_mrs_img[:], + nifti_mrs_img.dwelltime, + nifti_mrs_img.hdr_ext, + nifti_mrs_img.getAffine('voxel', 'world'), + nifti_version=1)\ + .save(out) + else: + nifti_mrs_img.save(out) if write_json: out_json = self.outputDir / (f_out + '.json') @@ -448,15 +475,13 @@ def philips(self, args): # philips specific imports from spec2nii.Philips.philips import read_sdat_spar_pair - self.imageOut = read_sdat_spar_pair(args.sdat, args.spar, args.shape, args.tags) - - # name of output - if args.fileout: - mainStr = args.fileout - else: - mainStr = args.sdat.stem - - self.fileoutNames.append(mainStr) + self.imageOut, self.fileoutNames = read_sdat_spar_pair( + args.sdat, + args.spar, + args.shape, + args.tags, + args.fileout, + args.special) # Philips data/list (+SPAR) handler def philips_dl(self, args): @@ -468,7 +493,7 @@ def philips_dl(self, args): # 'tags': (args.tag5, args.tag6, args.tag7)} # self.imageOut, file_names = read_data_list_pair(args.data, args.list, args.spar, overrides) - self.imageOut, file_names = read_data_list_pair(args.data, args.list, args.spar) + self.imageOut, file_names = read_data_list_pair(args.data, args.list, args.aux, args.special) # name of output if args.fileout: @@ -525,16 +550,8 @@ def raw(self, args): # Bruker 2dseq files with FG_COMPLEX def bruker(self, args): from spec2nii.bruker import read_bruker - from warnings import warn - from brukerapi.exceptions import MissingProperty - warn('Bruker conversion must currently be run with numpy<1.20.0') - - try: - self.imageOut, self.fileoutNames = read_bruker(args) - except MissingProperty: - raise Spec2niiError( - 'Bruker conversion must currently be run with numpy<1.20.0. ' - 'The underlying brukerapi package requires updating for ongoing compatibility.') + + self.imageOut, self.fileoutNames = read_bruker(args) # Varian parser def varian(self, args): diff --git a/spec2nii/uih.py b/spec2nii/uih.py index f8a7654..8f1664e 100644 --- a/spec2nii/uih.py +++ b/spec2nii/uih.py @@ -3,13 +3,17 @@ Author: William Clarke Copyright (C) 2020 University of Oxford """ +from datetime import datetime +from warnings import warn + import numpy as np import nibabel.nicom.dicomwrappers + +from nifti_mrs.create_nmrs import gen_nifti_mrs_hdr_ext +from nifti_mrs.hdr_ext import Hdr_Ext + from spec2nii.dcm2niiOrientation.orientationFuncs import dcm_to_nifti_orientation from spec2nii.nifti_orientation import NIFTIOrient -from spec2nii import nifti_mrs -from datetime import datetime -from warnings import warn from spec2nii import __version__ as spec2nii_ver @@ -94,7 +98,13 @@ def all_equal(lst): meta_used.set_dim_info(0, tag) # Create NIFTI MRS object. - nifti_mrs_out.append(nifti_mrs.NIfTI_MRS(combined_data, or_used.Q44, dt_used, meta_used)) + nifti_mrs_out.append( + gen_nifti_mrs_hdr_ext( + combined_data, + dt_used, + meta_used, + or_used.Q44, + no_conj=True)) else: for idx, (dd, oo, dt, mm, ff) in enumerate(zip(data_list, orientation_list, @@ -104,7 +114,13 @@ def all_equal(lst): # Add original files to nifti meta information. mm.set_standard_def('OriginalFile', [str(ff.name), ]) fnames_out.append(f'{mainStr}_{idx:03}') - nifti_mrs_out.append(nifti_mrs.NIfTI_MRS(dd, oo.Q44, dt, mm)) + nifti_mrs_out.append( + gen_nifti_mrs_hdr_ext( + dd, + dt, + mm, + oo.Q44oo.Q44, + no_conj=True)) return nifti_mrs_out, fnames_out @@ -192,8 +208,9 @@ def extractDicomMetadata(dcmdata): """ # Extract required metadata and create hdr_ext object - obj = nifti_mrs.hdr_ext(dcmdata.dcm_data.TransmitterFrequency, - dcmdata.dcm_data.ResonantNucleus) + obj = Hdr_Ext( + dcmdata.dcm_data.TransmitterFrequency, + dcmdata.dcm_data.ResonantNucleus) def set_if_present(tag, val): if val: diff --git a/spec2nii/varian_importer.py b/spec2nii/varian_importer.py index bdfa6ee..4bbe709 100644 --- a/spec2nii/varian_importer.py +++ b/spec2nii/varian_importer.py @@ -4,16 +4,19 @@ Copyright (C) 2021 University of Oxford and University of Aarhus """ -import numpy as np -from spec2nii.nifti_orientation import NIFTIOrient, calc_affine -from spec2nii import nifti_mrs +import re +import warnings from datetime import datetime from os.path import basename, splitext -from spec2nii import __version__ as spec2nii_ver -import re +import numpy as np + +from nifti_mrs.create_nmrs import gen_nifti_mrs_hdr_ext +from nifti_mrs.hdr_ext import Hdr_Ext + +from spec2nii import __version__ as spec2nii_ver +from spec2nii.nifti_orientation import NIFTIOrient, calc_affine import spec2nii.varian as v -import warnings def read_varian(args): @@ -87,7 +90,7 @@ def read_varian(args): orientation = NIFTIOrient(affine) # create object - meta = nifti_mrs.hdr_ext(imagingfreq, nucleus) + meta = Hdr_Ext(imagingfreq, nucleus) meta.set_standard_def('ConversionMethod', f'spec2nii v{spec2nii_ver}') meta.set_standard_def('EchoTime', echotime) meta.set_standard_def('RepetitionTime', repetition_time) @@ -133,16 +136,16 @@ def read_varian(args): else: fname_out = [splitext(basename(args.file))[0], ] - return [nifti_mrs.NIfTI_MRS(data, orientation.Q44, dwelltime, meta), ], fname_out + return [gen_nifti_mrs_hdr_ext(data, dwelltime, meta, orientation.Q44, no_conj=True), ], fname_out def _varian_orientation_1d(params): '''Calculate 1d slice orientation from parameters struct extract voxel positions if available NB: 0d spect -- thk etc not defined -   1d slice -- thk, and then pss along the usual (theta, phi, psi) angles (and orient) -   SVS voxel - vorient, vpsi, vphi, vtheta: define angle of voxel; -   - pos1, pos2, pos3 : define position of voxel + 1d slice -- thk, and then pss along the usual (theta, phi, psi) angles (and orient) + SVS voxel - vorient, vpsi, vphi, vtheta: define angle of voxel; + - pos1, pos2, pos3 : define position of voxel - vox1, vox2, vox3, thkunit : define size of voxel (no thkunit is mm) ''' warnings.warn('Not yet implemented') diff --git a/tests/spec2nii_test_data b/tests/spec2nii_test_data index dc3905d..f3adaec 160000 --- a/tests/spec2nii_test_data +++ b/tests/spec2nii_test_data @@ -1 +1 @@ -Subproject commit dc3905d2efdc778b390a0ad1193a19c4c03623b7 +Subproject commit f3adaec956be2480bcef63bb997fe8364ea2bf3c diff --git a/tests/test_bruker.py b/tests/test_bruker.py index 4760988..c663132 100644 --- a/tests/test_bruker.py +++ b/tests/test_bruker.py @@ -9,7 +9,6 @@ import json import numpy as np -import pytest from .io_for_tests import read_nifti_mrs @@ -18,7 +17,6 @@ data_path = bruker_path / '20201208_105201_lego_rod_1_3' -@pytest.mark.skip(reason="Awaiting fix of brukerapi package") def test_fid(tmp_path): subprocess.check_call(['spec2nii', 'bruker', @@ -64,14 +62,13 @@ def test_fid(tmp_path): # Img 9 - svs img = read_nifti_mrs(tmp_path / 'fid_FID_9.nii.gz') - assert img.shape == (1, 1, 1, 1980, 1) + assert img.shape == (1, 1, 1, 1980) assert np.iscomplexobj(img.dataobj) assert np.isclose(1 / img.header['pixdim'][4], 4401.41) hdr_ext_codes = img.header.extensions.get_codes() hdr_ext = json.loads(img.header.extensions[hdr_ext_codes.index(44)].get_content()) - assert hdr_ext['dim_5'] == 'DIM_DYN' assert np.isclose(hdr_ext['SpectrometerFrequency'][0], 400.32251) assert hdr_ext['ResonantNucleus'][0] == '1H' assert hdr_ext['OriginalFile'][0] == str(data_path.absolute() / '9' / 'fid') @@ -81,14 +78,13 @@ def test_fid(tmp_path): # Img 10 - svs img = read_nifti_mrs(tmp_path / 'fid_FID_10.nii.gz') - assert img.shape == (1, 1, 1, 1980, 1) + assert img.shape == (1, 1, 1, 1980) assert np.iscomplexobj(img.dataobj) assert np.isclose(1 / img.header['pixdim'][4], 4401.41) hdr_ext_codes = img.header.extensions.get_codes() hdr_ext = json.loads(img.header.extensions[hdr_ext_codes.index(44)].get_content()) - assert hdr_ext['dim_5'] == 'DIM_DYN' assert np.isclose(hdr_ext['SpectrometerFrequency'][0], 400.32251) assert hdr_ext['ResonantNucleus'][0] == '1H' assert hdr_ext['OriginalFile'][0] == str(data_path.absolute() / '10' / 'fid') @@ -96,7 +92,6 @@ def test_fid(tmp_path): assert 'acqp' in hdr_ext -@pytest.mark.skip(reason="Awaiting fix of brukerapi package") def test_2dseq(tmp_path): subprocess.check_call(['spec2nii', 'bruker', diff --git a/tests/test_nifti_options.py b/tests/test_nifti_options.py index 4dcec03..20aa50d 100644 --- a/tests/test_nifti_options.py +++ b/tests/test_nifti_options.py @@ -3,7 +3,8 @@ import subprocess import nibabel as nib from pathlib import Path - +import numpy as np +from nifti_mrs.nifti_mrs import NIFTI_MRS # Data paths siemens_path = Path(__file__).parent / 'spec2nii_test_data' / 'Siemens' @@ -30,9 +31,13 @@ def test_process(tmp_path): assert type(img_nifti2) is nib.nifti2.Nifti2Image assert type(img_nifti1) is nib.nifti1.Nifti1Image + img_nifti2 = NIFTI_MRS(tmp_path / 'nifti_2.nii.gz') + img_nifti1 = NIFTI_MRS(tmp_path / 'nifti_1.nii.gz') -def test_generation(tmp_path): - from spec2nii.nifti_mrs import get_mrs_class + assert type(img_nifti2.header) is nib.nifti2.Nifti2Header + assert type(img_nifti1.header) is nib.nifti1.Nifti1Header - assert get_mrs_class(nifti=2).__bases__[0] is nib.nifti2.Nifti2Image - assert get_mrs_class(nifti=1).__bases__[0] is nib.nifti1.Nifti1Image + # Check parameters stored are the same + assert np.allclose(img_nifti2[:], img_nifti1[:]) + assert np.isclose(img_nifti2.dwelltime, img_nifti1.dwelltime) + assert img_nifti2.hdr_ext.current_keys == img_nifti1._hdr_ext.current_keys diff --git a/tests/test_other.py b/tests/test_other.py index 1fc3012..1265a31 100644 --- a/tests/test_other.py +++ b/tests/test_other.py @@ -59,7 +59,7 @@ def test_text(affine_file, tmp_path): assert converted.shape == (1, 1, 1, 4096) assert np.iscomplexobj(converted.dataobj) assert np.allclose(np.loadtxt(affine_file), converted.affine) - assert np.allclose(converted.dataobj, original.dataobj) + assert np.allclose(converted.dataobj[:], original.dataobj[:]) assert (tmp_path / "outdir" / 'text.json').exists() @@ -80,7 +80,7 @@ def test_raw(affine_file, tmp_path): assert converted.shape == (1, 1, 1, 4096) assert np.iscomplexobj(converted.dataobj) assert np.allclose(np.loadtxt(affine_file), converted.affine) - assert np.allclose(converted.dataobj, original.dataobj) + assert np.allclose(converted.dataobj[:], original.dataobj[:]) assert (tmp_path / 'raw.json').exists() @@ -103,5 +103,5 @@ def test_raw_nohdr(affine_file, tmp_path): assert converted.shape == (1, 1, 1, 4096) assert np.iscomplexobj(converted.dataobj) assert np.allclose(np.loadtxt(affine_file), converted.affine) - assert np.allclose(converted.dataobj, original.dataobj) + assert np.allclose(converted.dataobj[:], original.dataobj[:]) assert (tmp_path / 'raw.json').exists() diff --git a/tests/test_philips_data_list.py b/tests/test_philips_data_list.py new file mode 100644 index 0000000..379b513 --- /dev/null +++ b/tests/test_philips_data_list.py @@ -0,0 +1,64 @@ +'''Tests for Philips SDAT/SPAR format conversion. + +Copyright William Clarke, University of Oxford 2021 +Subject to the BSD 3-Clause License. +''' + +import subprocess +from pathlib import Path +import json + +import numpy as np + +from .io_for_tests import read_nifti_mrs + +# Data paths +philips_path = Path(__file__).parent / 'spec2nii_test_data' / 'philips' +list_path = philips_path / 'hyper' / 'raw_226.list' +data_path = philips_path / 'hyper' / 'raw_226.data' +aux_path = philips_path / 'hyper' / 'converted_dcm.dcm' + + +def test_svs_hyper(tmp_path): + + subprocess.check_call(['spec2nii', 'philips_dl', + '-f', 'svs', + '-o', tmp_path, + '-j', + '--special', 'hyper', + str(data_path), + str(list_path), + str(aux_path)]) + + img_1 = read_nifti_mrs(tmp_path / 'svs_hyper_short_te.nii.gz') + img_2 = read_nifti_mrs(tmp_path / 'svs_hyper_edited.nii.gz') + img_3 = read_nifti_mrs(tmp_path / 'svs_NOI.nii.gz') + + assert img_1.shape == (1, 1, 1, 2048, 15, 32) + assert np.iscomplexobj(img_1.dataobj) + assert 1 / img_1.header['pixdim'][4] == 2000.0 + + assert img_2.shape == (1, 1, 1, 2048, 15, 4, 56) + assert np.iscomplexobj(img_2.dataobj) + assert 1 / img_2.header['pixdim'][4] == 2000.0 + + assert img_3.shape == (1, 1, 1, 2048, 15) + assert np.iscomplexobj(img_2.dataobj) + assert 1 / img_2.header['pixdim'][4] == 2000.0 + + hdr_ext_codes = img_1.header.extensions.get_codes() + hdr_ext = json.loads(img_1.header.extensions[hdr_ext_codes.index(44)].get_content()) + + assert hdr_ext['SpectrometerFrequency'][0] == 127.74876 + assert hdr_ext['ResonantNucleus'][0] == '1H' + assert hdr_ext['OriginalFile'][0] == data_path.name + assert hdr_ext['dim_5'] == 'DIM_COIL' + assert hdr_ext['dim_6'] == 'DIM_DYN' + + hdr_ext_codes = img_2.header.extensions.get_codes() + hdr_ext = json.loads(img_2.header.extensions[hdr_ext_codes.index(44)].get_content()) + + assert hdr_ext['dim_5'] == 'DIM_COIL' + assert hdr_ext['dim_6'] == 'DIM_EDIT' + assert hdr_ext['dim_6_header'] == {'EditCondition': ['A', 'B', 'C', 'D']} + assert hdr_ext['dim_7'] == 'DIM_DYN' diff --git a/tests/test_philips_sdat_spar.py b/tests/test_philips_sdat_spar.py index cfe91bf..bad56e2 100644 --- a/tests/test_philips_sdat_spar.py +++ b/tests/test_philips_sdat_spar.py @@ -20,6 +20,9 @@ svs_edit_path_sdat = philips_path / 'HERCULES_spar_sdat' / 'HERCULES_Example_noID.sdat' svs_edit_path_spar = philips_path / 'HERCULES_spar_sdat' / 'HERCULES_Example_noID.spar' +hyper_path_sdat = philips_path / 'hyper' / 'HBCD_HYPER_r5712_Export_WIP_HYPER_5_2_raw_act.SDAT' +hyper_path_spar = philips_path / 'hyper' / 'HBCD_HYPER_r5712_Export_WIP_HYPER_5_2_raw_act.SPAR' + def test_svs(tmp_path): @@ -69,3 +72,39 @@ def test_svs_edit(tmp_path): assert hdr_ext['OriginalFile'][0] == svs_edit_path_sdat.name assert hdr_ext['dim_5'] == 'DIM_EDIT' assert hdr_ext['dim_6'] == 'DIM_DYN' + + +def test_svs_hyper(tmp_path): + + subprocess.check_call(['spec2nii', 'philips', + '-f', 'svs', + '-o', tmp_path, + '-j', + str(hyper_path_sdat), + str(hyper_path_spar)]) + + img_1 = read_nifti_mrs(tmp_path / 'svs_hyper_short_te.nii.gz') + img_2 = read_nifti_mrs(tmp_path / 'svs_hyper_edited.nii.gz') + + assert img_1.shape == (1, 1, 1, 2048, 32) + assert np.iscomplexobj(img_1.dataobj) + assert 1 / img_1.header['pixdim'][4] == 2000.0 + + assert img_2.shape == (1, 1, 1, 2048, 4, 56) + assert np.iscomplexobj(img_2.dataobj) + assert 1 / img_2.header['pixdim'][4] == 2000.0 + + hdr_ext_codes = img_1.header.extensions.get_codes() + hdr_ext = json.loads(img_1.header.extensions[hdr_ext_codes.index(44)].get_content()) + + assert hdr_ext['SpectrometerFrequency'][0] == 127.74876 + assert hdr_ext['ResonantNucleus'][0] == '1H' + assert hdr_ext['OriginalFile'][0] == hyper_path_sdat.name + assert hdr_ext['dim_5'] == 'DIM_DYN' + + hdr_ext_codes = img_2.header.extensions.get_codes() + hdr_ext = json.loads(img_2.header.extensions[hdr_ext_codes.index(44)].get_content()) + + assert hdr_ext['dim_5'] == 'DIM_EDIT' + assert hdr_ext['dim_5_header'] == {'EditCondition': ['A', 'B', 'C', 'D']} + assert hdr_ext['dim_6'] == 'DIM_DYN' diff --git a/tests/test_twix.py b/tests/test_twix.py index d58c9e0..2c58882 100644 --- a/tests/test_twix.py +++ b/tests/test_twix.py @@ -86,7 +86,7 @@ def test_XA20(tmp_path): hdr_ext_codes = img_ref.header.extensions.get_codes() hdr_ext = json.loads(img_ref.header.extensions[hdr_ext_codes.index(44)].get_content()) - assert img_ref.shape == (1, 1, 1, 2080, 42, 1) + assert img_ref.shape == (1, 1, 1, 2080, 42) assert np.iscomplexobj(img_ref.dataobj) assert hdr_ext['dim_5'] == 'DIM_COIL' @@ -113,7 +113,7 @@ def test_XA30(tmp_path): hdr_ext_codes = img_ref.header.extensions.get_codes() hdr_ext = json.loads(img_ref.header.extensions[hdr_ext_codes.index(44)].get_content()) - assert img_ref.shape == (1, 1, 1, 4096, 44, 1) + assert img_ref.shape == (1, 1, 1, 4096, 44) assert np.iscomplexobj(img_ref.dataobj) assert hdr_ext['dim_5'] == 'DIM_COIL'