Skip to content

Commit

Permalink
export siemens csa header
Browse files Browse the repository at this point in the history
constracti committed May 8, 2019
1 parent 5f91cff commit 53eac54
Showing 3 changed files with 114 additions and 91 deletions.
85 changes: 0 additions & 85 deletions csa2.py

This file was deleted.

98 changes: 95 additions & 3 deletions dicomtools.py
Original file line number Diff line number Diff line change
@@ -2,6 +2,7 @@
import re
import math
import datetime
from collections import OrderedDict

import numpy
import pydicom
@@ -10,8 +11,6 @@
except ImportError:
pass

import csa2

def file_is_valid(filename):
try:
pydicom.dcmread(filename, specific_tags=[])
@@ -48,7 +47,7 @@ def get_affine(dataset1, dataset2):
# number of slices
nslices = 1
if (0x0029, 0x0010) in dataset1 and dataset1[0x0029, 0x0010].value == "SIEMENS CSA HEADER":
csa_image_header_info = csa2.decode(dataset1[0x0029, 0x1010].value)
csa_image_header_info = csa2_decode(dataset1[0x0029, 0x1010].value)
if csa_image_header_info["NumberOfImagesInMosaic"]["Data"]:
nslices = int(csa_image_header_info["NumberOfImagesInMosaic"]["Data"][0])
# shape, zooms & affine
@@ -89,6 +88,99 @@ def get_affine(dataset1, dataset2):
# return
return shape, zooms, affine


########
# csa2 #
########

# http://nipy.org/nibabel/dicom/siemens_csa.html

def csa2_decode(arr):
assert arr[0:4] == b"SV10"
# arr[4:8] # unused
ntags = int.from_bytes(arr[8:12], "little")
# int.from_bytes(arr[12:16], "little") # unused 77
i = 16
hdr = OrderedDict()
for ctags in range(0, ntags):
key = arr[i:i+64].split(bytes(1), 1)[0].decode()
hdr[key] = {}
hdr[key]["VM"] = int.from_bytes(arr[i+64:i+68], "little")
hdr[key]["VR"] = arr[i+68:i+72].split(bytes(1), 1)[0].decode()
hdr[key]["SyngoDT"] = int.from_bytes(arr[i+72:i+76], "little")
nitems = int.from_bytes(arr[i+76:i+80], "little")
# int.from_bytes(arr[80:84], "little") # unused 77 or 205
i += 84
hdr[key]["Data"] = []
for citems in range(nitems):
item_len = int.from_bytes(arr[i:i+4], "little")
# int.from_bytes(arr[i+4:i+8], "little") # unused
# int.from_bytes(arr[i+8:i+12], "little") # unused
# int.from_bytes(arr[i+12:i+16], "little") # unused
if item_len:
hdr[key]["Data"].append(arr[i+16:i+16+item_len].split(bytes(1), 1)[0].decode())
else:
hdr[key]["Data"].append(None)
i += 16 + math.ceil(item_len / 4) * 4
return hdr

def csa2_encode(hdr):
arr = b""
arr += b"SV10"
arr += bytes(4) # unused
arr += len(hdr).to_bytes(4, "little")
arr += bytes(4) # unused 77
for key, val in hdr.items():
arr += key.encode().ljust(64, bytes(1))
arr += val["VM"].to_bytes(4, "little")
arr += val["VR"].encode().ljust(4, bytes(1))
arr += val["SyngoDT"].to_bytes(4, "little")
arr += len(val["Data"]).to_bytes(4, "little")
arr += bytes(4) # unused 77 or 205
for dat in val["Data"]:
item_len = len(dat) + 1 if dat is not None else 0
arr += item_len.to_bytes(4, "little")
arr += bytes(4) # unused
arr += bytes(4) # unused
arr += bytes(4) # unused
if item_len:
arr += dat.encode().ljust(math.ceil(item_len / 4) * 4, bytes(1))
arr += bytes(4) # append one more zero
return arr

def csa2_diff(hdr1, hdr2):
keys1 = iter(sorted(hdr1.keys()))
keys2 = iter(sorted(hdr2.keys()))
key1 = next(keys1, None)
key2 = next(keys2, None)
while key1 and key2:
val1 = hdr1[key1]
val2 = hdr2[key2]
if key1 < key2:
print("< {}: {}".format(key1, val1))
key1 = next(keys1, None)
elif key2 < key1:
print("> {}: {}".format(key2, val2))
key2 = next(keys2, None)
else:
if val1 != val2:
print("< {}: {}".format(key1, val1))
print("> {}: {}".format(key2, val2))
key1 = next(keys1, None)
key2 = next(keys2, None)
while key1:
val1 = hdr1[key1]
print("< {}: {}".format(key1, val1))
key1 = next(keys1, None)
while key2:
val2 = hdr2[key2]
print("> {}: {}".format(key2, val2))
key2 = next(keys2, None)

##########
# linear #
##########

def _prop(dataset, tag):
if type(tag) is str:
return dataset.data_element(tag)
22 changes: 19 additions & 3 deletions nifti2dicom.py
Original file line number Diff line number Diff line change
@@ -43,8 +43,6 @@ def nifti2dicom(path):
if (0x0019, 0x0010) in dataset and dataset[0x0019, 0x0010].value == "SIEMENS MR HEADER":
# (0x0019, 0x100b) SliceMeasurementDuration
dataset.pop((0x0019, 0x100b), None)
# (0x0019, 0x1016) TimeAfterStart
dataset.pop((0x0019, 0x1016), None)
# (0x0019, 0x1029) MosaicRefAcqTimes
dataset.pop((0x0019, 0x1029), None)
# Window Explanation Algorithm not specified
@@ -56,7 +54,7 @@ def nifti2dicom(path):
dataset.pop((0x0028, 0x1055), None)
if (0x0029, 0x0010) in dataset and dataset[0x0029, 0x0010].value == "SIEMENS CSA HEADER":
# (0x0029, 0x1010) CSA Image Header Info
dataset.pop((0x0029, 0x1010), None)
csa_image_header_info = dicomtools.csa2_decode(dataset[0x0029, 0x1010].value)
if (0x0043, 0x0010) in dataset and dataset[0x0043, 0x0010].value == "GEMS_PARM_01":
# (0x0043, 0x1028) Unique image iden
dataset.pop((0x0043, 0x1028), None)
@@ -117,6 +115,9 @@ def nifti2dicom(path):
# (0x0019, 0x1015) SlicePosition_PCS
if (0x0019, 0x1015) in dataset:
dicomtools.linear_float_array((0x0019, 0x1015), dataset, dataset1, dataset2)
# (0x0019, 0x1016) TimeAfterStart
if (0x0019, 0x1016) in dataset:
dicomtools.linear_float((0x0019, 0x1016), dataset, dataset1, dataset2)
elif (0x0019, 0x0010) in dataset and dataset[0x0019, 0x0010].value in ["GEMS_ACQU_01", "GEMS_IDEN_01"]:
# (0x0019, 0x10a2) Raw data run number
if (0x0019, 0x10a2) in dataset:
@@ -139,6 +140,21 @@ def nifti2dicom(path):
# (0x0028, 0x0107) Largest Image Pixel Value
if (0x0028, 0x0107) in dataset:
dataset[0x0028, 0x0107].value = data_slice.max()
if (0x0029, 0x0010) in dataset and dataset[0x0029, 0x0010].value == "SIEMENS CSA HEADER":
# (0x0029, 0x1010) CSA Image Header Info
if csa_image_header_info["Actual3DImaPartNumber"]["Data"]:
csa_image_header_info["Actual3DImaPartNumber"]["Data"][0] = str(f).ljust(8)
elif not csa_image_header_info["MosaicRefAcqTimes"]["Data"]: # TODO linear int on csa_image_header_info
csa_image_header_info["ProtocolSliceNumber"]["Data"][0] = str(f).ljust(8)
# csa_image_header_info["GSWDDataType"] CORONAL
# csa_image_header_info["RFSWDDataType"] CORONAL
# csa_image_header_info["ICE_Dims"]["Data"][0] *
# csa_image_header_info["MosaicRefAcqTimes"]["Data"] FMRI
# csa_image_header_info["SliceMeasurementDuration"]["Data"][0] CORONAL
csa_image_header_info["SlicePosition_PCS"]["Data"][0:3] = ["{:.8f}".format(x) for x in dataset.ImagePositionPatient]
if csa_image_header_info["TimeAfterStart"]["Data"]:
csa_image_header_info["TimeAfterStart"]["Data"][0] = "{:.8f}".format(dataset[0x0019, 0x1016].value)
dataset[0x0029, 0x1010].value = dicomtools.csa2_encode(csa_image_header_info)
if (0x2001, 0x0010) in dataset and dataset[0x2001, 0x0010].value == "Philips Imaging DD 001":
# (0x2001, 0x100a) Slice Number MR
if (0x2001, 0x100a) in dataset:

0 comments on commit 53eac54

Please sign in to comment.