From 29f1952feda40122eff2125e121f96e59e4ac12b Mon Sep 17 00:00:00 2001 From: jmcox9 <61171698+jmcox9@users.noreply.github.com> Date: Wed, 12 Apr 2023 15:40:10 -0400 Subject: [PATCH 1/3] Add files via upload This notebook allows a single .AIM or .ISQ file, specified through the variable filename, to be imported and undergo operations for conversion from Hounsfield Units to Density in g/cc and to a binary image mask given user-specified image thresholds and numeric constants extracted from the image header metadata. These operations are NOT optimized and 2 main issues remain: 1. The image header metadata is incomplete, as noted by other users in 2 separate open issues. This workbook uses a template header to show proof of concept, but a finalized version will require header information extracted from each file separately. 2. Imported images have a high memory requirement for array operations needed to convert them to alternate units or a binary image mask. It may be necessary to optimize memory use to reduce required memory and increase throughput --- ITKIOScanco_UnitConversion.ipynb | 212 +++++++++++++++++++++++++++++++ 1 file changed, 212 insertions(+) create mode 100644 ITKIOScanco_UnitConversion.ipynb diff --git a/ITKIOScanco_UnitConversion.ipynb b/ITKIOScanco_UnitConversion.ipynb new file mode 100644 index 0000000..ba536b3 --- /dev/null +++ b/ITKIOScanco_UnitConversion.ipynb @@ -0,0 +1,212 @@ +{ + "cells": [ + { + "attachments": {}, + "cell_type": "markdown", + "metadata": {}, + "source": [ + "This workbook loads and converts a SCANCO .ISQ or .AIM file in the directory to several \n", + ".nrrd formatted images with the following units: density (g/cc), Hounsfield Units, and a\n", + "binary thresholded image mask \n", + "\n", + "Created by: Jason Cox, 04/07/2023\n", + "\n", + "Citations: McCormick, Matthew, et al. \"ITK: Enabling Reproducible \n", + "Research and Open Science.” Frontiers in Neuroinformatics, vol. 8, \n", + "Frontiers Media SA, 2014, doi:10.3389/fninf.2014.00013." + ] + }, + { + "cell_type": "code", + "execution_count": 20, + "metadata": {}, + "outputs": [], + "source": [ + "import numpy as np\n", + "import cv2\n", + "import re\n", + "import itk\n", + "import glob" + ] + }, + { + "cell_type": "code", + "execution_count": 2, + "metadata": {}, + "outputs": [], + "source": [ + "# ISQ and AIM metadata returned by dict (cell 3) contain the same fields. Therefore, this program could be modified to loop over all files\n", + "# in the directory, if that is the desired implementation. The variable all_filenames is unused in the current version, but can be used in a\n", + "# modified version that implements a loop.\n", + "\n", + "all_AIM_filenames = glob.glob(\"*.AIM\")\n", + "all_ISQ_filenames = glob.glob(\"*.ISQ\")\n", + "all_filenames = all_AIM_filenames + all_ISQ_filenames\n", + "\n", + "del all_AIM_filenames, all_ISQ_filenames\n" + ] + }, + { + "cell_type": "code", + "execution_count": 22, + "metadata": {}, + "outputs": [], + "source": [ + "filename = 'C0005757_SCAP.AIM'\n", + "image = itk.imread(filename)\n", + "metadata = dict(image)\n", + "HUpixels = itk.array_view_from_image(image)\n", + "del image\n", + "# HUpixels = np.array(HUpixels) # Array arithmetic below seems to work without this line, but certain operations were not tested as they require too much memory" + ] + }, + { + "cell_type": "code", + "execution_count": null, + "metadata": {}, + "outputs": [], + "source": [ + "# As discussed below, some required constants are missing from metadata. Compare the output of this cell with the file header at:\n", + "# https://docs.google.com/document/d/1L2GCnHa4jmK87p0pGI9MNN10MWRMwrk2WqPbPeklP_I/edit\n", + "\n", + "print(metadata)" + ] + }, + { + "attachments": {}, + "cell_type": "markdown", + "metadata": {}, + "source": [ + "The next cell provides a template file header from C0005757_SCAP.AIM, but this needs to be read from the .AIM or .ISQ.\n", + "\n", + "Then, the following cell identifies numeric constants in the header necessary for converting the image in HU (units given by itk.imread)\n", + "to density, or for thresholding/masking the image.\n", + "\n", + "The only constants not readable from the header which must be provided by the user are the per-milli values for low and high thresholds. In basic terms, the low/high threshold pixel values are computed as the respective low/high per-milli value multiplied by the highest pixel value in the image (a quanitity that is given by the header)." + ] + }, + { + "cell_type": "code", + "execution_count": 11, + "metadata": {}, + "outputs": [], + "source": [ + "header = '{Global Density: intercept=-1.96600006e+02, Global Calib. default unit type=2 (Density), Global Position Slice 1 [um]=47050, Global HU: mu water=0.51400, Global Reconstruction-Alg.=3, Global Parameter units=[1/cm], Global Average data value=322.51215, Global Minimum data value=-2891.00000, Global Created by=ISQ_TO_AIM (IPL), Global Patient Name=JC-NBPI-Rat-May17, Global Energy [V]=70000, Global Scanner ID=2135, Global Standard deviation=0.32443, Global Time=29-JUN-2022 12:39:11.75, Global Orig-ISQ-Dim-p=1792 1792 1781, Global Original Creation-Date=20-AUG-2018 23:18:59.99, Global Parameter name=Linear Attenuation, Global Orig-ISQ-Dim-um=17918 17918 17808, Global Scanner type=10, Global Scan Distance [um]=20480, Global Index Patient=136, Global Index Measurement=7866, Global Minimum value=-0.70581, Global Maximum data value=32767.00000, Global Intensity [uA]=114, Global Maximum value=7.99976, Global Density: slope=3.55519989e+02, Global Mu_Scaling=4096, Global Reference line [um]=47050, Global No. samples=2048, Global Default-Eval=15, Global Site=4, Global Average value=0.07874, Global Calibration Data=Unknown. Backconverted from AIM., Global Standard data deviation=1328.87610, Global No. projections per 180=1000, Global Original file=DK0:[MICROCT.DATA.00000136.00007866]B0005757.ISQ, Global Scaled by factor=4096, Global Angle-Offset [mdeg]=0, Global Integration time [us]=800000, Global Density: unit=mg HA/ccm}'" + ] + }, + { + "cell_type": "code", + "execution_count": 16, + "metadata": {}, + "outputs": [], + "source": [ + "# Gather numeric constants needed for unit conversion. Currently, the file metadata returned with dict does not contain all of the\n", + "# necessary constants. For an example of a full file header with highlighted variables of interest, see the google doc \n", + "# \"SCANCO micro-CT Units Conversion\" at: https://docs.google.com/document/d/1L2GCnHa4jmK87p0pGI9MNN10MWRMwrk2WqPbPeklP_I/edit\n", + "\n", + "mu_threshLow = 0.250\n", + "mu_threshHigh = 1.0\n", + "\n", + "scancoMU = float(re.search('Global Scaled by factor=(.*?),', header).group(1))\n", + "# scancoMU = metadata['MuScaling']\n", + "maxVal = float(re.search('Global Maximum data value=(.*?),', header).group(1))\n", + "# maxVal = metadata['**********']\n", + "densitySlope = float(re.search('slope=(.*?),', header).group(1))\n", + "# densitySlope = metadata['**********']\n", + "densityIntercept = float(re.search('intercept=(.*?),', header).group(1))\n", + "# densityIntercept = metadata['**********']\n", + "water_mu = float(re.search('water=(.*?),', header).group(1))\n", + "# water_mu = metadata['MuWater']\n", + "\n", + "del header # when the code base is modified such that metadata gives all relevant constants, this line and others that reference \"header\" can be deleted\n", + "del metadata\n", + "\n", + "muLow = mu_threshLow * maxVal\n", + "mu_pixelLow = muLow/scancoMU\n", + "muHigh = mu_threshHigh * maxVal\n", + "mu_pixelHigh = muHigh/scancoMU\n", + "HU_pixelLow = (-1000) + mu_pixelLow * (1000/water_mu)\n", + "HU_pixelHigh = (-1000) + mu_pixelHigh * (1000/water_mu)\n", + "HUvalues = [HU_pixelLow, HU_pixelHigh]" + ] + }, + { + "cell_type": "code", + "execution_count": 7, + "metadata": {}, + "outputs": [], + "source": [ + "# Create unit-converted outputs. MUpixels is not desired as an output due to representing non-standard units. However, creation of MUpixels\n", + "# can be skipped, as RHOpixels can be created by combining the functions needed to convert HU to MU and MU to RHO.\n", + "\n", + "# MUpixels = (HUpixels+1000)/(1000/water_mu)\n", + "# RHOpixels = (MUpixels*densitySlope) + densityIntercept\n", + "RHOpixels = ((((HUpixels+1000)/(1000/water_mu))*densitySlope) + densityIntercept)/1000 # Dividing by 1000 at the end changes units from mg/cc to g/cc, which is more commonly published\n", + "temp = HUpixels.copy()\n", + "temp[temp>HUvalues[1]] = 0\n", + "BWmask = cv2.threshold(temp, HUvalues[0], 1, cv2.THRESH_BINARY)[1]" + ] + }, + { + "cell_type": "code", + "execution_count": 8, + "metadata": {}, + "outputs": [], + "source": [ + "# Create .nrrd file output with units of density (mg/cc or g/cc)\n", + "\n", + "RHOout = itk.image_view_from_array(RHOpixels)\n", + "RHOoutname = filename[0:-4] + '_RHO.nrrd'\n", + "itk.imwrite(RHOout, RHOoutname)" + ] + }, + { + "cell_type": "code", + "execution_count": 9, + "metadata": {}, + "outputs": [], + "source": [ + "# Create a binary image mask .nrrd file output based on user-provided low and high threshold values.\n", + "\n", + "BWout = itk.image_view_from_array(BWmask)\n", + "BWoutname = filename[0:-4] + '_BWmask.nrrd'\n", + "itk.imwrite(RHOout, BWoutname)\n" + ] + }, + { + "cell_type": "code", + "execution_count": 6, + "metadata": {}, + "outputs": [], + "source": [ + "# Create .nrrd file output with Hounsfield Units (equivalent to loading the file in 3D Slicer)\n", + "\n", + "HUout = itk.image_view_from_array(HUpixels)\n", + "HUoutname = filename[0:-4] + '_HU.nrrd'\n", + "itk.imwrite(HUout, HUoutname)" + ] + } + ], + "metadata": { + "kernelspec": { + "display_name": "kitfisto", + "language": "python", + "name": "python3" + }, + "language_info": { + "codemirror_mode": { + "name": "ipython", + "version": 3 + }, + "file_extension": ".py", + "mimetype": "text/x-python", + "name": "python", + "nbconvert_exporter": "python", + "pygments_lexer": "ipython3", + "version": "3.11.2" + }, + "orig_nbformat": 4 + }, + "nbformat": 4, + "nbformat_minor": 2 +} From cb646017bda180142f0bd36c36126b01b83111f7 Mon Sep 17 00:00:00 2001 From: jmcox9 <61171698+jmcox9@users.noreply.github.com> Date: Wed, 5 Jul 2023 14:47:36 -0400 Subject: [PATCH 2/3] Revised units conversion Previous version implemented a binary thresholding operation through cv2.threshold. The new version uses itk.BinaryThresholdImageFilter to create a binarized volume version of the input image. --- ITKIOScanco_UnitConversion.ipynb | 67 +++++++++++++++++++------------- 1 file changed, 41 insertions(+), 26 deletions(-) diff --git a/ITKIOScanco_UnitConversion.ipynb b/ITKIOScanco_UnitConversion.ipynb index ba536b3..c6a2437 100644 --- a/ITKIOScanco_UnitConversion.ipynb +++ b/ITKIOScanco_UnitConversion.ipynb @@ -18,15 +18,15 @@ }, { "cell_type": "code", - "execution_count": 20, + "execution_count": 1, "metadata": {}, "outputs": [], "source": [ "import numpy as np\n", - "import cv2\n", "import re\n", "import itk\n", - "import glob" + "import glob\n", + "import argparse" ] }, { @@ -48,7 +48,7 @@ }, { "cell_type": "code", - "execution_count": 22, + "execution_count": 3, "metadata": {}, "outputs": [], "source": [ @@ -87,7 +87,7 @@ }, { "cell_type": "code", - "execution_count": 11, + "execution_count": 7, "metadata": {}, "outputs": [], "source": [ @@ -96,7 +96,7 @@ }, { "cell_type": "code", - "execution_count": 16, + "execution_count": 11, "metadata": {}, "outputs": [], "source": [ @@ -118,8 +118,8 @@ "water_mu = float(re.search('water=(.*?),', header).group(1))\n", "# water_mu = metadata['MuWater']\n", "\n", - "del header # when the code base is modified such that metadata gives all relevant constants, this line and others that reference \"header\" can be deleted\n", - "del metadata\n", + "# del header # when the code base is modified such that metadata gives all relevant constants, this line and others that reference \"header\" can be deleted\n", + "# del metadata\n", "\n", "muLow = mu_threshLow * maxVal\n", "mu_pixelLow = muLow/scancoMU\n", @@ -132,7 +132,38 @@ }, { "cell_type": "code", - "execution_count": 7, + "execution_count": 12, + "metadata": {}, + "outputs": [], + "source": [ + "PixelType = itk.SS\n", + "Dimension = 3\n", + "\n", + "ImageType = itk.Image[PixelType, Dimension]\n", + "\n", + "reader = itk.ImageFileReader[ImageType].New()\n", + "reader.SetFileName(filename)\n", + "\n", + "thresholdFilter = itk.BinaryThresholdImageFilter[ImageType, ImageType].New()\n", + "thresholdFilter.SetInput(reader.GetOutput())\n", + "\n", + "thresholdFilter.SetLowerThreshold(int(HUvalues[0]))\n", + "thresholdFilter.SetUpperThreshold(int(HUvalues[1]))\n", + "thresholdFilter.SetOutsideValue(int(0))\n", + "thresholdFilter.SetInsideValue(int(1))\n", + "\n", + "BWoutname = filename[0:-4] + '_BWmask.nrrd'\n", + "\n", + "writer = itk.ImageFileWriter[ImageType].New()\n", + "writer.SetFileName(BWoutname)\n", + "writer.SetInput(thresholdFilter.GetOutput())\n", + "\n", + "writer.Update()" + ] + }, + { + "cell_type": "code", + "execution_count": null, "metadata": {}, "outputs": [], "source": [ @@ -141,10 +172,7 @@ "\n", "# MUpixels = (HUpixels+1000)/(1000/water_mu)\n", "# RHOpixels = (MUpixels*densitySlope) + densityIntercept\n", - "RHOpixels = ((((HUpixels+1000)/(1000/water_mu))*densitySlope) + densityIntercept)/1000 # Dividing by 1000 at the end changes units from mg/cc to g/cc, which is more commonly published\n", - "temp = HUpixels.copy()\n", - "temp[temp>HUvalues[1]] = 0\n", - "BWmask = cv2.threshold(temp, HUvalues[0], 1, cv2.THRESH_BINARY)[1]" + "RHOpixels = ((((HUpixels+1000)/(1000/water_mu))*densitySlope) + densityIntercept)/1000 # Dividing by 1000 at the end changes units from mg/cc to g/cc, which is more commonly published\n" ] }, { @@ -160,19 +188,6 @@ "itk.imwrite(RHOout, RHOoutname)" ] }, - { - "cell_type": "code", - "execution_count": 9, - "metadata": {}, - "outputs": [], - "source": [ - "# Create a binary image mask .nrrd file output based on user-provided low and high threshold values.\n", - "\n", - "BWout = itk.image_view_from_array(BWmask)\n", - "BWoutname = filename[0:-4] + '_BWmask.nrrd'\n", - "itk.imwrite(RHOout, BWoutname)\n" - ] - }, { "cell_type": "code", "execution_count": 6, From fb2fe4c25b08b7696e818e6c0797ab45576b3af2 Mon Sep 17 00:00:00 2001 From: jmcox9 <61171698+jmcox9@users.noreply.github.com> Date: Tue, 8 Aug 2023 17:22:26 -0400 Subject: [PATCH 3/3] Unit Conversion Second Revision --- ITKIOScanco_UnitConversion.ipynb | 150 +++++++++++++++++-------------- 1 file changed, 85 insertions(+), 65 deletions(-) diff --git a/ITKIOScanco_UnitConversion.ipynb b/ITKIOScanco_UnitConversion.ipynb index c6a2437..00e4880 100644 --- a/ITKIOScanco_UnitConversion.ipynb +++ b/ITKIOScanco_UnitConversion.ipynb @@ -18,15 +18,27 @@ }, { "cell_type": "code", - "execution_count": 1, + "execution_count": null, + "metadata": {}, + "outputs": [], + "source": [ + "# Install non-default packages\n", + "import sys\n", + "!{sys.executable} -m pip install itk pooch" + ] + }, + { + "cell_type": "code", + "execution_count": 4, "metadata": {}, "outputs": [], "source": [ - "import numpy as np\n", + "# Import required packages\n", + "from pathlib import Path\n", "import re\n", "import itk\n", - "import glob\n", - "import argparse" + "import pooch\n", + "import glob" ] }, { @@ -41,23 +53,28 @@ "\n", "all_AIM_filenames = glob.glob(\"*.AIM\")\n", "all_ISQ_filenames = glob.glob(\"*.ISQ\")\n", - "all_filenames = all_AIM_filenames + all_ISQ_filenames\n", - "\n", - "del all_AIM_filenames, all_ISQ_filenames\n" + "all_filenames = all_AIM_filenames + all_ISQ_filenames\n" ] }, { "cell_type": "code", - "execution_count": 3, + "execution_count": 15, "metadata": {}, "outputs": [], "source": [ - "filename = 'C0005757_SCAP.AIM'\n", - "image = itk.imread(filename)\n", + "# This cell originally operated on the image C0005757_SCAP.AIM, though for the purposes of employing pooch, a sample file stored online at the provided URL is now used.\n", + "# This section will need to be modified to handle general implementation.\n", + "\n", + "# filename = 'C0005757_SCAP.AIM'\n", + "filename = 'C0004255.ISQ'\n", + "file_url = 'https://data.kitware.com/api/v1/file/591e56178d777f16d01e0d20/download'\n", + "file_sha256 = 'c2a3750c75826cb077d92093d43976cc0350198b55edecd681265eebabfb438b'\n", + "file_path = pooch.retrieve(file_url, file_sha256, fname=filename, progressbar=False)\n", + "image = itk.imread(file_path)\n", "metadata = dict(image)\n", - "HUpixels = itk.array_view_from_image(image)\n", - "del image\n", - "# HUpixels = np.array(HUpixels) # Array arithmetic below seems to work without this line, but certain operations were not tested as they require too much memory" + "HU_pixels = itk.array_view_from_image(image)\n", + "# del image\n", + "# HU_pixels = np.array(HU_pixels) # Array arithmetic below seems to work without this line, but certain operations were not tested as they require too much memory" ] }, { @@ -68,8 +85,14 @@ "source": [ "# As discussed below, some required constants are missing from metadata. Compare the output of this cell with the file header at:\n", "# https://docs.google.com/document/d/1L2GCnHa4jmK87p0pGI9MNN10MWRMwrk2WqPbPeklP_I/edit\n", + "# This section does, however, collect some metadata and write it to a .txt file.\n", "\n", - "print(metadata)" + "print(metadata)\n", + "metadata_str = str(metadata)\n", + "\n", + "metadata_outname = Path(filename).stem + \"_metadata.txt\"\n", + "with open(metadata_outname, \"w\") as file:\n", + " file.write(metadata_str)" ] }, { @@ -87,7 +110,7 @@ }, { "cell_type": "code", - "execution_count": 7, + "execution_count": 17, "metadata": {}, "outputs": [], "source": [ @@ -96,7 +119,7 @@ }, { "cell_type": "code", - "execution_count": 11, + "execution_count": 18, "metadata": {}, "outputs": [], "source": [ @@ -104,101 +127,98 @@ "# necessary constants. For an example of a full file header with highlighted variables of interest, see the google doc \n", "# \"SCANCO micro-CT Units Conversion\" at: https://docs.google.com/document/d/1L2GCnHa4jmK87p0pGI9MNN10MWRMwrk2WqPbPeklP_I/edit\n", "\n", - "mu_threshLow = 0.250\n", - "mu_threshHigh = 1.0\n", - "\n", - "scancoMU = float(re.search('Global Scaled by factor=(.*?),', header).group(1))\n", - "# scancoMU = metadata['MuScaling']\n", - "maxVal = float(re.search('Global Maximum data value=(.*?),', header).group(1))\n", - "# maxVal = metadata['**********']\n", - "densitySlope = float(re.search('slope=(.*?),', header).group(1))\n", - "# densitySlope = metadata['**********']\n", - "densityIntercept = float(re.search('intercept=(.*?),', header).group(1))\n", - "# densityIntercept = metadata['**********']\n", + "# Thresholds must be entered by the user, they are not contained in the file header\n", + "mu_threshold_low = 0.250\n", + "mu_threshold_high = 1.0\n", + "\n", + "scanco_mu = float(re.search('Global Scaled by factor=(.*?),', header).group(1))\n", + "# scanco_mu = metadata['MuScaling']\n", + "max_val = float(re.search('Global Maximum data value=(.*?),', header).group(1))\n", + "# max_val = metadata['**********']\n", + "density_slope = float(re.search('slope=(.*?),', header).group(1))\n", + "# density_slope = metadata['**********']\n", + "density_intercept = float(re.search('intercept=(.*?),', header).group(1))\n", + "# density_intercept = metadata['**********']\n", "water_mu = float(re.search('water=(.*?),', header).group(1))\n", "# water_mu = metadata['MuWater']\n", "\n", "# del header # when the code base is modified such that metadata gives all relevant constants, this line and others that reference \"header\" can be deleted\n", "# del metadata\n", "\n", - "muLow = mu_threshLow * maxVal\n", - "mu_pixelLow = muLow/scancoMU\n", - "muHigh = mu_threshHigh * maxVal\n", - "mu_pixelHigh = muHigh/scancoMU\n", - "HU_pixelLow = (-1000) + mu_pixelLow * (1000/water_mu)\n", - "HU_pixelHigh = (-1000) + mu_pixelHigh * (1000/water_mu)\n", - "HUvalues = [HU_pixelLow, HU_pixelHigh]" + "mu_low = mu_threshold_low * max_val\n", + "mu_pixel_low = mu_low/scanco_mu\n", + "mu_high = mu_threshold_high * max_val\n", + "mu_pixel_high = mu_high/scanco_mu\n", + "HU_pixel_low = (-1000) + mu_pixel_low * (1000/water_mu)\n", + "HU_pixel_high = (-1000) + mu_pixel_high * (1000/water_mu)\n", + "HU_values = [HU_pixel_low, HU_pixel_high]" ] }, { "cell_type": "code", - "execution_count": 12, + "execution_count": 25, "metadata": {}, "outputs": [], "source": [ - "PixelType = itk.SS\n", - "Dimension = 3\n", - "\n", - "ImageType = itk.Image[PixelType, Dimension]\n", + "# Create a binarized version of the input image by applying low and high thresholds. Threshold values are accepted as int type variables, so some rounding is necessary.\n", + "lower_threshold = int(round(HU_values[0]))\n", + "upper_threshold = int(round(HU_values[1]))\n", + "inside_value = 1\n", + "outside_value = 0\n", "\n", - "reader = itk.ImageFileReader[ImageType].New()\n", - "reader.SetFileName(filename)\n", + "BW_image = itk.binary_threshold_image_filter(image,lower_threshold=lower_threshold,upper_threshold=upper_threshold,inside_value=inside_value,outside_value=outside_value)\n", "\n", - "thresholdFilter = itk.BinaryThresholdImageFilter[ImageType, ImageType].New()\n", - "thresholdFilter.SetInput(reader.GetOutput())\n", + "BW_outname = Path(filename).stem + \"_BWmask.nrrd\"\n", "\n", - "thresholdFilter.SetLowerThreshold(int(HUvalues[0]))\n", - "thresholdFilter.SetUpperThreshold(int(HUvalues[1]))\n", - "thresholdFilter.SetOutsideValue(int(0))\n", - "thresholdFilter.SetInsideValue(int(1))\n", + "itk.imwrite(BW_image, BW_outname)\n", "\n", - "BWoutname = filename[0:-4] + '_BWmask.nrrd'\n", - "\n", - "writer = itk.ImageFileWriter[ImageType].New()\n", - "writer.SetFileName(BWoutname)\n", - "writer.SetInput(thresholdFilter.GetOutput())\n", - "\n", - "writer.Update()" + "# The following can be enabled to free up memory\n", + "# del image\n", + "# del BW_image" ] }, { "cell_type": "code", - "execution_count": null, + "execution_count": 26, "metadata": {}, "outputs": [], "source": [ - "# Create unit-converted outputs. MUpixels is not desired as an output due to representing non-standard units. However, creation of MUpixels\n", + "# Create unit-converted outputs. MUpixels is not generally desired as an output due to representing non-standard units. However, creation of MUpixels\n", "# can be skipped, as RHOpixels can be created by combining the functions needed to convert HU to MU and MU to RHO.\n", "\n", "# MUpixels = (HUpixels+1000)/(1000/water_mu)\n", "# RHOpixels = (MUpixels*densitySlope) + densityIntercept\n", - "RHOpixels = ((((HUpixels+1000)/(1000/water_mu))*densitySlope) + densityIntercept)/1000 # Dividing by 1000 at the end changes units from mg/cc to g/cc, which is more commonly published\n" + "RHO_pixels = ((((HU_pixels+1000)/(1000/water_mu))*density_slope) + density_intercept)/1000 # Dividing by 1000 at the end changes units from mg/cc to g/cc, which is more commonly published\n" ] }, { "cell_type": "code", - "execution_count": 8, + "execution_count": 29, "metadata": {}, "outputs": [], "source": [ "# Create .nrrd file output with units of density (mg/cc or g/cc)\n", "\n", - "RHOout = itk.image_view_from_array(RHOpixels)\n", - "RHOoutname = filename[0:-4] + '_RHO.nrrd'\n", - "itk.imwrite(RHOout, RHOoutname)" + "RHO_out = itk.image_view_from_array(RHO_pixels)\n", + "for k, v in metadata.items():\n", + " RHO_out[k] = v\n", + "RHO_outname = Path(filename).stem + \"_RHO.nrrd\"\n", + "itk.imwrite(RHO_out, RHO_outname)" ] }, { "cell_type": "code", - "execution_count": 6, + "execution_count": 30, "metadata": {}, "outputs": [], "source": [ "# Create .nrrd file output with Hounsfield Units (equivalent to loading the file in 3D Slicer)\n", "\n", - "HUout = itk.image_view_from_array(HUpixels)\n", - "HUoutname = filename[0:-4] + '_HU.nrrd'\n", - "itk.imwrite(HUout, HUoutname)" + "HU_out = itk.image_view_from_array(HU_pixels)\n", + "for k, v in metadata.items():\n", + " HU_out[k] = v\n", + "HU_outname = Path(filename).stem + \"_HU.nrrd\"\n", + "itk.imwrite(HU_out, HU_outname)" ] } ],