Description
Hello, this seems very similar to the issue posted here: #1014
I notice an issue where I simply read a dicom series (CT image data, Float64) and then re-write it without changes to test the functionality.
If I explicitly set the rescale and intercept values to the original CT dicom values (-1000.035 and 0.09397), then the output series is correct - data ranges from ~-1000 to +2055.
However, if I set the rescale value to 0.01 and the intercept to 0 (per the documentation: https://simpleitk.readthedocs.io/en/master/link_DicomSeriesFromArray_docs.html) then the output series data ranges from zero to some irrelevant value depending on how many significant digits the rescale value is.
Am I missing something? I was under the impression that ITK would internally do the scaling but this seems not to be the case.
Thanks,
Marc
(Code below for reference)
def writeSeries(IMAGEobj):
out_dir = "output_DICOM_"+time.strftime("%Y%m%d%H%M%S")
os.mkdir(out_dir)
#pixel_dtypes = {"int16": np.int16, "float64": np.float64}
pixel_dtype = sitk.GetArrayFromImage(IMAGEobj).dtype # np.int16 #np.float64
# Copy relevant tags from the original meta-data dictionary (private tags are
# also accessible).
tags_to_copy = [
"0010|0010", # Patient Name
"0010|0020", # Patient ID
"0010|0030", # Patient Birth Date
"0020|000d", # Study Instance UID, for machine consumption
"0020|0010", # Study ID, for human consumption
"0008|0020", # Study Date
"0008|0030", # Study Time
"0008|0050", # Accession Number
"0008|0060", # Modality
"0018|5100", # Patient Position (Marc)
]
def writeSlices(series_tag_values, new_img, out_dir, i):
image_slice = new_img[:, :, i]
# Tags shared by the series.
list(
map(
lambda tag_value: image_slice.SetMetaData(tag_value[0], tag_value[1]),
series_tag_values,
)
)
# Slice specific tags.
# Instance Creation Date
image_slice.SetMetaData("0008|0012", time.strftime("%Y%m%d"))
# Instance Creation Time
image_slice.SetMetaData("0008|0013", time.strftime("%H%M%S"))
# Setting the type to CT so that the slice location is preserved and
# the thickness is carried over.
image_slice.SetMetaData("0008|0060", "CT")
# (0020, 0032) image position patient determines the 3D spacing between
# slices.
# Image Position (Patient)
image_slice.SetMetaData(
"0020|0032",
"\\".join(map(str, new_img.TransformIndexToPhysicalPoint((0, 0, i)))),
)
# Instance Number
image_slice.SetMetaData("0020|0013", str(i))
# Write to the output directory and add the extension dcm, to force
# writing in DICOM format.
writer.SetFileName(os.path.join(out_dir, str(i) + ".dcm"))
writer.Execute(image_slice)
####
writer = sitk.ImageFileWriter()
# Use the study/series/frame of reference information given in the meta-data
# dictionary and not the automatically generated information from the file IO
writer.KeepOriginalImageUIDOn()
modification_time = time.strftime("%H%M%S")
modification_date = time.strftime("%Y%m%d")
# Copy some of the tags and add the relevant tags indicating the change.
# For the series instance UID (0020|000e), each of the components is a number,
# cannot start with zero, and separated by a '.' We create a unique series ID
# using the date and time. Tags of interest:
direction = IMAGEobj.GetDirection()
series_tag_values = [
(k, moving_reader.GetMetaData(0, k))
for k in tags_to_copy
if moving_reader.HasMetaDataKey(0, k)
] + [
("0008|0031", modification_time), # Series Time
("0008|0021", modification_date), # Series Date
("0008|0008", "DERIVED\\SECONDARY"), # Image Type
(
"0020|000e",
"1.2.826.0.1.3680043.2.1125." + modification_date + ".1" + modification_time,
), # Series Instance UID
(
"0020|0037",
"\\".join(
map(
str,
(
direction[0],
direction[3],
direction[6],
direction[1],
direction[4],
direction[7],
),
)
),
), # Image Orientation
# (Patient)
("0008|103e", "Created-SimpleITK"), # Series Description
]
if pixel_dtype == np.float64:
print("Processing as float64")
# If we want to write floating point values, we need to use the rescale
# slope, "0028|1053", to select the number of digits we want to keep. We
# also need to specify additional pixel storage and representation
# information.
rescale_slope = 0.1# 0.93977#0.001 # keep three digits after the decimal point
series_tag_values = series_tag_values + [
("0028|1053", str(rescale_slope)), # rescale slope
("0028|1052", "0"),#-1000.0358"), # rescale intercept
("0028|0100", "16"), # bits allocated
("0028|0101", "16"), # bits stored
("0028|0102", "15"), # high bit
("0028|0103", "0"), # pixel representation
]
# Write slices to output directory
list(
map(
lambda i: writeSlices(series_tag_values, IMAGEobj, out_dir, i),
range(IMAGEobj.GetDepth()),
)
);
print("...Done!")
def sitkDCM(dir, meta=False):
print("Reading Dicom directory:", dir)
reader = sitk.ImageSeriesReader()
dicom_names = reader.GetGDCMSeriesFileNames(dir)
reader.SetFileNames(dicom_names)
reader.MetaDataDictionaryArrayUpdateOn() #MM
reader.LoadPrivateTagsOn() #MM
image = reader.Execute()
size = image.GetSize()
print("Image size:", size[0], size[1], size[2])
return image, reader
moving_image, moving_reader = sitkDCM(dir= "zzCase10_LFrontal")
print( sitk.GetArrayFromImage(moving_image).dtype ) ## PRINTS "INT16"
writeSeries(moving_image)
Activity