script_DRR_similarity.ipynb

{
 "cells": [
  {
   "cell_type": "code",
   "execution_count": 1,
   "id": "saved-marking",
   "metadata": {},
   "outputs": [],
   "source": [
    "# Assess the similarity between generated DRRs and dataLoaded DRRs from PolyU\n",
    "import torch\n",
    "import numpy as np\n",
    "import matplotlib.pyplot as plt\n",
    "import os, sys, time, datetime\n",
    "import pathlib\n",
    "import pydicom.fileset\n",
    "import pandas as pd\n",
    "\n",
    "from scipy import ndimage, misc\n",
    "import Metrics\n",
    "import ImageProcessing\n",
    "\n",
    "import skimage.exposure\n"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 12,
   "id": "dramatic-collins",
   "metadata": {},
   "outputs": [],
   "source": [
    "# Dataloader for DRR images\n",
    "def load_DRR_image(csv_file_name = pathlib.Path(\"D:/data/CT-Covid-19-August2020/DRR_output/DRR_directory.csv\"),\n",
    "                   images_root_directory=pathlib.Path(\"D:\\data\\CT-Covid-19-August2020\\DRR_output\"), \n",
    "                   idx=0, image_type='source', flip=False):\n",
    "    \n",
    "    df = pd.read_csv(csv_file_name)\n",
    "    image = np.load(os.path.join( images_root_directory, df[image_type][idx]))\n",
    "    \n",
    "    PixelSpacing_i = df[\"PixelSize_i\"][idx]\n",
    "    PixelSpacing_j = df[\"PixelSize_j\"][idx]\n",
    "    PixelSpacing = [PixelSpacing_i, PixelSpacing_j]\n",
    "    \n",
    "    if flip:\n",
    "        image = np.flip(image, axis=0)\n",
    "    return image, PixelSpacing\n",
    "\n",
    "def zScoreNormalisation(image):\n",
    "    return (image - np.median(image.flatten()))/np.sqrt(np.var(image.flatten()))\n",
    "    \n",
    "def rescalingNormalisation(image):\n",
    "    return (image - np.amin(image))/(np.amax(image) - np.amin(image))\n",
    "\n",
    "\n",
    "\n",
    "\n",
    "\n",
    "# Create a dataloader for CXR images\n",
    "directoryPath = \"D:\\data\\POLYU_COVID19_CXR_CT_Cohort1\"\n",
    "path_to_file = []\n",
    "for root, dirs, files in os.walk(directoryPath):\n",
    "    for name in files:\n",
    "        if name==\"DICOMDIR\":\n",
    "            path_to_file.append(os.path.join(root, name))\n",
    "\n",
    "# Select a CXR directory\n",
    "path_idx = 0\n",
    "a_path = path_to_file[path_idx]\n",
    "\n",
    "fs = pydicom.fileset.FileSet(a_path)\n",
    "\n",
    "for instance in fs:\n",
    "    if \"1.2.840.10008.5.1.4.1.1.1\" in instance.SOPClassUID:  # check if CR or DX image\n",
    "        # CXR\n",
    "        ds = instance.load()\n",
    "        PixelSpacing_CXR = ds.ImagerPixelSpacing\n",
    "        image = ds.pixel_array\n",
    "        image_CXR = np.abs(np.amax(image)-image) # Dark is low intensity, bright is high intensity\n",
    "                \n",
    "        # DRR\n",
    "        image_DRR, PixelSpacing_DRR = load_DRR_image(idx = 632, flip=True)\n",
    "        \n",
    "        \n",
    "        #Z-normalisation\n",
    "        image_DRR = rescalingNormalisation(image_DRR)\n",
    "        image_CXR = rescalingNormalisation(image_CXR)\n",
    "        \n",
    "        # Rescale image in order to best match histograms\n",
    "        \n",
    "        \n",
    "        # Upsample\n",
    "        image_DRR = ndimage.zoom(image_DRR, (PixelSpacing_DRR[0]/PixelSpacing_CXR[0], PixelSpacing_DRR[1]/PixelSpacing_CXR[1]), order=0)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 33,
   "id": "damaged-funeral",
   "metadata": {},
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "1.3.51.0.7.2909126416.24284.5964.39496.11845.57303.9026\n"
     ]
    },
    {
     "data": {
      "application/vnd.jupyter.widget-view+json": {
       "model_id": "7ae72b2e7cd64f42bad082139e4a5785",
       "version_major": 2,
       "version_minor": 0
      },
      "text/plain": [
       "Canvas(toolbar=Toolbar(toolitems=[('Home', 'Reset original view', 'home', 'home'), ('Back', 'Back to previous …"
      ]
     },
     "metadata": {},
     "output_type": "display_data"
    },
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "0.7615501060925516\n"
     ]
    }
   ],
   "source": [
    "print(instance.SOPInstanceUID)\n",
    "\n",
    "%matplotlib widget\n",
    "\n",
    "fig, axs = plt.subplots(1, 2)\n",
    "axs[0].imshow(image_DRR , cmap='bone')\n",
    "axs[1].imshow(image_CXR , cmap='bone')\n",
    "\n",
    "# Quantify similarity\n",
    "image = image_DRR\n",
    "means, variances = Metrics.patchMetrics(image)\n",
    "#print(means)\n",
    "SSIM = Metrics.SSIM(image_CXR, image_DRR, filter_kernel=(101,101), dynamic_range=4)\n",
    "print(SSIM)\n",
    "\n"
   ]
  }
 ],
 "metadata": {
  "kernelspec": {
   "display_name": "Python 3",
   "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.6.13"
  }
 },
 "nbformat": 4,
 "nbformat_minor": 5
}