How to change only array for Dicom file with Simple ITK in python

Question

I have a bunch of medical images in dicom that I want to correct for bias field inhomogeneity using SimpleITK in Python. The workflow is straightforward: I want to (1) open the dicom image, (2) create a binary mask of the object in the image, (3) apply N4 bias field correction to the masked image, (4) write back the corrected image in dicom format. Note that no spatial transformation is applied to the image, but only intensity transformation, so that I could copy all spatial information and all meta data (except for date/hour of creation and instance number) from the original to the corrected image.

I have written this function to achieve my goal:

def n4_dcm_correction(dcm_in_file):
    metadata_to_set = ["0008|0012", "0008|0013", "0020|0013"]
    filepath = PurePath(dcm_in_file)
    root_dir = str(filepath.parent)
    file_name = filepath.stem
    dcm_reader = sitk.ImageFileReader()
    dcm_reader.SetFileName(dcm_in_file)
    dcm_reader.LoadPrivateTagsOn()
    inputImage = dcm_reader.Execute()
    metadata_to_copy = [k for k in inputImage.GetMetaDataKeys() if k not in metadata_to_set]
    maskImage = sitk.OtsuThreshold(inputImage,0,1,200)
    filledImage = sitk.BinaryFillhole(maskImage)
    floatImage = sitk.Cast(inputImage,sitk.sitkFloat32)
    corrector = sitk.N4BiasFieldCorrectionImageFilter();
    output = corrector.Execute(floatImage, filledImage)
    output.CopyInformation(inputImage)
    for k in metadata_to_copy:
        print("key is: {}; value is {}".format(k, inputImage.GetMetaData(k)))
        output.SetMetaData(k, inputImage.GetMetaData(k))
     output.SetMetaData("0008|0012", time.strftime("%Y%m%d"))
     output.SetMetaData("0008|0013", time.strftime("%H%M%S"))
     output.SetMetaData("0008|0013", str(float(inputImage.GetMetaData("0008|0013")) + randint(1, 999)))
    out_file = "{}/{}_biascorrected.dcm".format(root_dir, file_name)
    writer = sitk.ImageFileWriter()
    writer.KeepOriginalImageUIDOn()
    writer.SetFileName(out_file)
    writer.Execute(sitk.Cast(output, sitk.sitkUInt16))
    return

n4_dcm_correction("/path/to/my/dcm/image.dcm")

As much as the bias correction part works (the bias is removed), the writing part is a mess. I would expect my output dicom to have the exact same metadata of the original one, however they are all missing, notably the patient name, the protocol name and the manufacturer. Similalry, something is very wrong with the spatial information, since if I try to convert the dicom to the nifti format with dcm2niix, the directions are reversed: superior is down and inferior is up, forward is back and backward is front. What step am I missing ?

Answer 1

I suspect you are working with a MRI series, not a single file. Likely this example does what you want, read-modify-write a volume stored in a set of files.

If the example did not resolve your issue, please post to the ITK discourse which is the primary location for ITK/SimpleITK related discussions.