Intensity normalization test.ipynb

{
 "cells": [
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "# MRI intensity normalization\n",
    "\n",
    "Intensity normalization of multi-channel MRI images using the method proposed by [Nyul et al. 2000](http://ieeexplore.ieee.org/lpdocs/epic03/wrapper.htm?arnumber=836373).\n",
    "In the original paper, the authors suggest a method where a set of standard histogram landmarks are learned from a set of MRI images. These landmarks are then used to equalize the histograms of the images to normalize. In both learning and transformation, the histograms are used to find the intensity landmarks. \n",
    "\n",
    "## Ackwoledgements: \n",
    "* The Python implementation is based on the awesome implementation available here [Reinhold et al. 2019](https://github.com/jcreinhold/intensity-normalization). \n",
    "* For this particular tutorial, we use a very small subset from the [Calgary-Campinas dataset](https://sites.google.com/view/calgary-campinas-dataset/home/download)."
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "## Train the standard histogram:\n",
    "To train the standard histogram, we just have to create a list of the input images to process. Optionally, we can also provide the brainmasks:"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 1,
   "metadata": {},
   "outputs": [],
   "source": [
    "import os\n",
    "import numpy as np\n",
    "import nibabel as nib\n",
    "from nyul import nyul_train_standard_scale\n",
    "  \n",
    "DATA_DIR = 'data_examples'\n",
    "T1_name = 'T1.nii.gz'\n",
    "MASK_name = 'brainmask.nii.gz'\n",
    "\n",
    "# generate training scans\n",
    "train_scans = [os.path.join(DATA_DIR, folder, T1_name)\n",
    "               for folder in os.listdir(DATA_DIR)]\n",
    "mask_scans = [os.path.join(DATA_DIR, folder, MASK_name)\n",
    "              for folder in os.listdir(DATA_DIR)]"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Then, the train the standard histogram. By default, the parameters are set as follows:\n",
    "* Minimum percentile to consider `i_min=1`\n",
    "* Maximum percentile to consider `i_max=99`\n",
    "* Minimum percentile on the standard histogram `i_s_min=1`\n",
    "* Maximum percentile on the standard histogram `i_s_max=100`\n",
    "* Middle percentile lower bound `l_percentile=10`\n",
    "* Middle percentile upped bound `u_percentile=90`\n",
    "* number of deciles `step=10`\n",
    "\n"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 2,
   "metadata": {},
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "processing scan  data_examples/CC0358/T1.nii.gz\n",
      "processing scan  data_examples/CC0357/T1.nii.gz\n",
      "processing scan  data_examples/CC0006/T1.nii.gz\n",
      "processing scan  data_examples/CC0356/T1.nii.gz\n",
      "processing scan  data_examples/CC0140/T1.nii.gz\n",
      "processing scan  data_examples/CC0139/T1.nii.gz\n",
      "processing scan  data_examples/CC0002/T1.nii.gz\n",
      "processing scan  data_examples/CC0001/T1.nii.gz\n",
      "processing scan  data_examples/CC0226/T1.nii.gz\n",
      "processing scan  data_examples/CC0228/T1.nii.gz\n"
     ]
    }
   ],
   "source": [
    "standard_scale, perc = nyul_train_standard_scale(train_scans, mask_scans)"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "## Save the standard histogram:\n",
    "\n",
    "Save the histogram to apply it to unseen images afterwards:"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 3,
   "metadata": {},
   "outputs": [],
   "source": [
    "standard_path = 'histograms/standard_test.npy'\n",
    "np.save(standard_path, [standard_scale, perc])"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "## Apply intensity normalization to new images:\n",
    "\n",
    "Finally, the learned histogram can be applied to new images. Here, we just use the same images before and after normalization as an example. "
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 5,
   "metadata": {},
   "outputs": [
    {
     "data": {
      "text/plain": [
       "Text(0.5, 1.0, 'Image 1 Normalized')"
      ]
     },
     "execution_count": 5,
     "metadata": {},
     "output_type": "execute_result"
    },
    {
     "data": {
      "image/png": "iVBORw0KGgoAAAANSUhEUgAAAbgAAAEoCAYAAAAqrOTwAAAABHNCSVQICAgIfAhkiAAAAAlwSFlzAAALEgAACxIB0t1+/AAAADh0RVh0U29mdHdhcmUAbWF0cGxvdGxpYiB2ZXJzaW9uMy4xLjEsIGh0dHA6Ly9tYXRwbG90bGliLm9yZy8QZhcZAAAgAElEQVR4nO3df5SUV53n8fdHMJpJTCABGeRHOipHF7MrJr0BR3eGmJE0OA44J+YkOzv0RFacDdmN53iOEucHMRkdcmY1JqPiomGAWQ1homwYJRIW43HjCtIoEyCYpU3IoZEACYSoGeOQfPeP5xZ50lQ3TVdXV9ftz+ucOvXU97n1PPc+Kfqbe+s+txQRmJmZ5eZVja6AmZlZPTjBmZlZlpzgzMwsS05wZmaWJSc4MzPLkhOcmZllyQnOzOpC0p9Kerj0+peS3jjA5/iepP88kMe0fDjB2bAgaZ+k3290PXoj6SxJ96W6hqSZpyn/PUm/ljSpFPt9SfvqXdf+iIhzI+LxRtfDhg8nOLOh5WHgPwFP9bH8r4C/HIgTSxoxEMcxGyqc4GzYSUNnP5B0h6RnJT0u6XdSfL+kw5LaS+XfJ+knkp5L+2/pdrz5kp6U9Iykvyz3FiW9StJiST9L+9dKuqBavSLiNxHx+Yh4GHixj825C7hO0pt6aOu/ST29ZyXtlvSHpX0rJS2TtEHSr4ArUuxLkh5IQ4o/kPTbkj4v6Zikn0p6R+kYlbb9QtKjkj7QU0VTr/TNkt6Qjl15PC8pSuU+JGlPOt9GSReV9r031eG4pC8A6uN1smHICc6Gq+nAI8CFwNeBNcC/B95M0YP6gqRzU9lfAfOBUcD7gP8iaR6ApKnAl4A/BsYD5wMTSuf5r8A84PeANwDHgC8OYDsOAF8BPtV9h6RXA/8EPAi8PtXla5LeUir2H4FPA6+j6D0CXAP8BTAGeAH4IfDj9Po+4HOl9/8M+A8U7f4U8D8lje+twhHx8zRceW5EnAuso7j+SJoLfBL4I2As8H+Ae9K+McA3S3X7GfCuXq+ODWtOcDZcPRERfx8RLwL3ApOAWyPihYh4EPgNRbIjIr4XETsj4qWIeITiD+7vpeNcDfxTRDwcEb8B/gooL/D6Z8CfR0RXRLwA3AJcLWnkALblb4D3S3pbt/gM4Fxgaeodfhf4FnBdqcz9EfGD1LZfp9i6iNieXq8Dfh0Rq0vX6mQPLiL+MSWslyLiXmAvcHlfKy7pE8BbgQ+l0J8BfxMReyLiBPAZYFrqxc0BdkfEfRHxr8Dn6ftQrg1DTnA2XB0qbf8LQER0j50LIGm6pIckHZF0nOKP8JhU7g3A/sqbIuJ54JnScS4C1qUhwmeBPRTDj+MGqiERcQT4AnBrt11vAPZHxEul2JO8soe5n1N1vw5VrwucHJ7dUWrfJbx8bXolaTZwEzAvIv4lhS8C7iwd7yjFMOQETr3W0UP9zQAnOLO++DqwHpgUEecDX+bl734OAhMrBSWdTTHsWbEfmB0Ro0qP10bEgQGu498CVwCXlWI/ByZJKv87n0wxrFnR758TSb2qrwA3AhdGxChgF334XiwNk64CromIcpLaD3yk2/U6OyL+L8W1Ls8YVfm1WXdOcGan9zrgaET8WtLlFN9bVdxHMTz4O5LOohiCLP+B/zLw6cpECUlj0/dMVUl6jaTXppdnSXpt+kPeq4h4Fvgs8PFSeCvwPPBxSa9Otx28n/R91wA4hyJBHkl1v56iB9crSecB91MM3T7cbfeXgZsrw62Szpf0wbTv28DbJP1RGuL9b8BvD0hLLEtOcGandwNwq6RfUHzHtrayIyJ2U0zeWEPRw/glcJhicgbAnRS9vwfT+7dQTHDpyWMUw4ATgI1p+6JeypfdSWn2ZfpO8P3AbOBpiskw8yPip308Xq8i4lGKpPpDimHMfwv8oA9vvRR4C3BHeTZlOuY64HZgjaTnKHqEs9O+p4EPAksphoGn9PF8NkzJP3hqNnDSzMtngSkR8USj62M2nLkHZ1YjSe+X9FuSzgH+O7AT2NfYWpmZE5xZ7eZSTOj4OcWw2bXhoRGzhvMQpZmZZem0PThJK9LSRbtKsVskHUj3v+yQNKe072ZJnZIek3RVKd6WYp2SFpfiF0vamuL3pploldlk96b4VkktpzuHmZlZRV+GKFcCbVXid0TEtPTYACeXLboWeFt6z5ckjVCxiOsXKWZDTaVYO29qOs7t6VhvpljGaEGKLwCOpfgdqVyP5zizZpuZWe5Ou1xQRHy/3Hs6jbnAmrQk0ROSOnl52Z7Oyk9lSFoDzJW0B3gPL99XtIriPqJl6Vi3pPh9FGsDqpdz/LC3io0ZMyZaWvraDDMzawbbt29/OiLGVttXy3p4N0qaD3QAH4uIYxT37mwpleni5WWB9neLT6dY8eHZtOZc9/ITKu+JiBNpiaQLT3OOHrW0tNDR0dH31pmZ2ZAn6cme9vV3FuUy4E3ANIqbWz/bz+PUlaSFkjokdRw5cqTR1TEzs0HUrwQXEYci4sW0iOtXeHkY8gCvXBtuYor1FH8GGFVaWb0Sf8Wx0v7zU/mejlWtnssjojUiWseOrdqDNTOzTPUrwXX7vacPUCynA8WSRNemGZAXU9wT9CNgGzAlzZg8i2KSyPp0r9BDFD85AtBOsUZd5ViVH528GvhuKt/TOczMzE467Xdwku4BZgJjJHUBS4CZkqZRLLS6D/gIFOvySVoLPAqcABal35BC0o0Ua+uNAFakNfwAPkGx7txfAz8B7k7xu4F/SJNIjlIkxV7P0Sgti799Smzf0vc1oCZmZlYxbG70bm1tjXpNMnGCMzNrDEnbI6K12j4v1WVmZllygjMzsyw5wZmZWZac4MzMLEtOcGZmliUnODMzy5ITnJmZZckJzszMsuQEZ2ZmWXKCMzOzLDnBmZlZlpzgzMwsS05wZmaWJSc4MzPLkhOcmZllyQnOzMyy5ARnZmZZcoIzM7MsOcGZmVmWnODMzCxLTnBmZpYlJzgzM8uSE5yZmWXJCc7MzLLkBGdmZllygjMzsyw5wZmZWZac4MzMLEunTXCSVkg6LGlXKXaBpE2S9qbn0SkuSXdJ6pT0iKRLS+9pT+X3SmovxS+TtDO95y5J6u85zMzMKvrSg1sJtHWLLQY2R8QUYHN6DTAbmJIeC4FlUCQrYAkwHbgcWFJJWKnMh0vva+vPOczMzMpOm+Ai4vvA0W7hucCqtL0KmFeKr47CFmCUpPHAVcCmiDgaEceATUBb2ndeRGyJiABWdzvWmZzDzMzspP5+BzcuIg6m7aeAcWl7ArC/VK4rxXqLd1WJ9+ccp5C0UFKHpI4jR470sWlmZpaDmieZpJ5XDEBdBvwcEbE8IlojonXs2LF1qJmZmQ1V/U1whyrDgun5cIofACaVyk1Msd7iE6vE+3MOMzOzk/qb4NYDlZmQ7cD9pfj8NNNxBnA8DTNuBGZJGp0ml8wCNqZ9z0makWZPzu92rDM5h5mZ2UkjT1dA0j3ATGCMpC6K2ZBLgbWSFgBPAtek4huAOUAn8DxwPUBEHJV0G7Atlbs1IioTV26gmKl5NvBAenCm5zAzMys7bYKLiOt62HVllbIBLOrhOCuAFVXiHcAlVeLPnOk5zMzMKrySiZmZZckJzszMsuQEZ2ZmWXKCMzOzLDnBmZlZlpzgzMwsS05wZmaWJSc4MzPLkhOcmZllyQnOzMyy5ARnZmZZcoIzM7MsOcGZmVmWnODMzCxLTnBmZpYlJzgzM8uSE5yZmWXJCc7MzLLkBGdmZllygjMzsyw5wZmZWZac4MzMLEtOcGZmliUnODMzy5ITnJmZZckJzszMsuQEZ2ZmWaopwUnaJ2mnpB2SOlLsAkmbJO1Nz6NTXJLuktQp6RFJl5aO057K75XUXopflo7fmd6r3s5hZmZWMRA9uCsiYlpEtKbXi4HNETEF2JxeA8wGpqTHQmAZFMkKWAJMBy4HlpQS1jLgw6X3tZ3mHGZmZkB9hijnAqvS9ipgXim+OgpbgFGSxgNXAZsi4mhEHAM2AW1p33kRsSUiAljd7VjVzmFmZgbUnuACeFDSdkkLU2xcRBxM208B49L2BGB/6b1dKdZbvKtKvLdzmJmZATCyxve/OyIOSHo9sEnST8s7IyIkRY3n6FVv50hJdyHA5MmT61kNs1doWfztU2L7lr6vATUxG75q6sFFxIH0fBhYR/Ed2qE0vEh6PpyKHwAmld4+McV6i0+sEqeXc3Sv3/KIaI2I1rFjx/a3mWZm1oT63YOTdA7wqoj4RdqeBdwKrAfagaXp+f70lvXAjZLWUEwoOR4RByVtBD5TmlgyC7g5Io5Kek7SDGArMB/4u9Kxqp3DrG7cKzNrLrUMUY4D1qWZ+yOBr0fEdyRtA9ZKWgA8CVyTym8A5gCdwPPA9QApkd0GbEvlbo2Io2n7BmAlcDbwQHpAkdiqncOsqVRLmuDEaTYQ+p3gIuJx4O1V4s8AV1aJB7Coh2OtAFZUiXcAl/T1HGaDracEVWtZM6udVzIxM7MsOcGZmVmWnODMzCxLtd4HZ5Ylf19m1vyc4GxYG6qJrK/18mxLs555iNLMzLLkBGdmZllygjMzsyw5wZmZWZac4MzMLEueRWnWxLyWpVnPnOBs2BiqtwTUg3/5wMxDlGZmliknODMzy5KHKC1Lw2k4sq88bGnDjXtwZmaWJSc4MzPLkhOcmZllyd/BmQ1jvo/OcuYEZ2an8IQUy4ETnDU1z5Y0s574OzgzM8uSe3Bm1icetrRm4wRnTcPDkUOPJ6nYUOYhSjMzy5J7cGY24DycaUOBE5wNSR6OzI+Tng02Jzgza5gz+R8ZJ0M7U02d4CS1AXcCI4CvRsTSBlfJzpB7atZXToZ2ppo2wUkaAXwReC/QBWyTtD4iHm1szawnTmY2WOrxWXPSbD5Nm+CAy4HOiHgcQNIaYC4wZBOc/8CbNa/h9O+3WjJvxu9QmznBTQD2l153AdPLBSQtBBaml7+U9Fid6jIGePoV5769TmdqjFPal5mc25dz28Dtq4u+/v0agL9zA9G+i3ra0cwJ7rQiYjmwvN7nkdQREa31Pk+juH3NK+e2gdvX7Ordvma+0fsAMKn0emKKmZmZNXWC2wZMkXSxpLOAa4H1Da6TmZkNEU07RBkRJyTdCGykuE1gRUTsblB16j4M2mBuX/PKuW3g9jW7urZPEVHP45uZmTVEMw9RmpmZ9cgJzszMsuQEVwNJbZIek9QpaXGj69MfklZIOixpVyl2gaRNkvam59EpLkl3pfY+IunSxtW8byRNkvSQpEcl7ZZ0U4pn0UZJr5X0I0n/nNr3qRS/WNLW1I5700QsJL0mve5M+1saWf++kDRC0k8kfSu9zqZtAJL2SdopaYekjhTL5fM5StJ9kn4qaY+kdw5m25zg+qm0VNhsYCpwnaSpja1Vv6wE2rrFFgObI2IKsDm9hqKtU9JjIbBskOpYixPAxyJiKjADWJT+O+XSxheA90TE24FpQJukGcDtwB0R8WbgGLAglV8AHEvxO1K5oe4mYE/pdU5tq7giIqaV7gnL5fN5J/CdiHgr8HaK/46D17aI8KMfD+CdwMbS65uBmxtdr362pQXYVXr9GDA+bY8HHkvb/wO4rlq5ZnkA91OsX5pdG4HfAn5MsaLP08DIFD/5WaWYdfzOtD0ylVOj695LmyamP4LvAb4FKJe2ldq4DxjTLdb0n0/gfOCJ7v8NBrNt7sH1X7WlwiY0qC4DbVxEHEzbTwHj0nZTtzkNWb0D2EpGbUxDeDuAw8Am4GfAsxFxIhUpt+Fk+9L+48CFg1vjM/J54OPAS+n1heTTtooAHpS0PS0vCHl8Pi8GjgB/n4aYvyrpHAaxbU5w1qso/leq6e8lkXQu8A3goxHxXHlfs7cxIl6MiGkUvZ3Lgbc2uEoDQtIfAIcjYnuj61Jn746ISymG6BZJ+t3yzib+fI4ELgWWRcQ7gF/x8nAkUP+2OcH1X85LhR2SNB4gPR9O8aZss6RXUyS3r0XEN1M4qzYCRMSzwEMUw3ajJFUWcii34WT70v7zgWcGuap99S7gDyXtA9ZQDFPeSR5tOykiDqTnw8A6iv9JyeHz2QV0RcTW9Po+ioQ3aG1zguu/nJcKWw+0p+12iu+tKvH5abbTDOB4aahhSJIk4G5gT0R8rrQrizZKGitpVNo+m+L7xT0Uie7qVKx7+yrtvhr4bvq/6CEnIm6OiIkR0ULx7+u7EfHHZNC2CknnSHpdZRuYBewig89nRDwF7Jf0lhS6kuLnzAavbY3+IrKZH8Ac4P9RfOfx542uTz/bcA9wEPhXiv/jWkDxvcVmYC/wv4ELUllRzBz9GbATaG10/fvQvndTDIE8AuxIjzm5tBH4d8BPUvt2AX+V4m8EfgR0Av8IvCbFX5ted6b9b2x0G/rYzpnAt3JrW2rLP6fH7srfkYw+n9OAjvT5/F/A6MFsm5fqMjOzLHmI0szMsuQEZ2ZmWXKCMzOzLDnBmZlZlpzgzMwsS05wZmaWJSc4MzPLkhOcmZllyQnOzMyyNPL0RfIwZsyYaGlpaXQ1zMxsAG3fvv3piBhbbd+wSXAtLS10dHQ0uhpmZjaAJD3Z0z4PUZqZWZac4MzMLEtOcGZmliUnODMzy9KwmWRiza9l8bdPie1b+r4G1MTMmoF7cGZmliX34GxIqtZbMzM7E+7BmZlZltyDs6bWU0/P382ZmXtwZmaWpX4nOEmTJD0k6VFJuyXdlOIXSNokaW96Hp3iknSXpE5Jj0i6tHSs9lR+r6T2UvwySTvTe+6SpN7OYWZmVlFLD+4E8LGImArMABZJmgosBjZHxBRgc3oNMBuYkh4LgWVQJCtgCTAduBxYUkpYy4APl97XluI9ncPMzAyoIcFFxMGI+HHa/gWwB5gAzAVWpWKrgHlpey6wOgpbgFGSxgNXAZsi4mhEHAM2AW1p33kRsSUiAljd7VjVzmFmZgYM0CQTSS3AO4CtwLiIOJh2PQWMS9sTgP2lt3WlWG/xripxejlH93otpOgtMnny5DNslQ0G3w5gZvVS8yQTSecC3wA+GhHPlfelnlfUeo7e9HaOiFgeEa0R0Tp2bNWfCzIzs0zVlOAkvZoiuX0tIr6ZwofS8CLp+XCKHwAmld4+McV6i0+sEu/tHGZmZkBtsygF3A3siYjPlXatByozIduB+0vx+Wk25QzgeBpm3AjMkjQ6TS6ZBWxM+56TNCOda363Y1U7h5mZGVDbd3DvAv4E2ClpR4p9ElgKrJW0AHgSuCbt2wDMATqB54HrASLiqKTbgG2p3K0RcTRt3wCsBM4GHkgPejmHmZkZUEOCi4iHAfWw+8oq5QNY1MOxVgArqsQ7gEuqxJ+pdg4zM7MKL9VlWfJP65iZl+oyM7MsOcGZmVmWnODMzCxLTnBmZpYlJzgzM8uSE5yZmWXJtwnYoGn0wsq+dcBseHEPzszMsuQEZ2ZmWXKCMzOzLDnBmZlZlpzgzMwsS05wZmaWJd8mYMNaT7cu+PYBs+ZXUw9O0gpJhyXtKsVukXRA0o70mFPad7OkTkmPSbqqFG9LsU5Ji0vxiyVtTfF7JZ2V4q9JrzvT/pZa2mFmZvmpdYhyJdBWJX5HRExLjw0AkqYC1wJvS+/5kqQRkkYAXwRmA1OB61JZgNvTsd4MHAMWpPgC4FiK35HKmZmZnVRTgouI7wNH+1h8LrAmIl6IiCeATuDy9OiMiMcj4jfAGmCuJAHvAe5L718FzCsda1Xavg+4MpU3MzMD6jfJ5EZJj6QhzNEpNgHYXyrTlWI9xS8Eno2IE93irzhW2n88lX8FSQsldUjqOHLkyMC0zMzMmkI9JpksA24DIj1/FvhQHc5zWhGxHFgO0NraGo2ow3DV6HUnzcwGvAcXEYci4sWIeAn4CsUQJMABYFKp6MQU6yn+DDBK0shu8VccK+0/P5U3MzMD6pDgJI0vvfwAUJlhuR64Ns2AvBiYAvwI2AZMSTMmz6KYiLI+IgJ4CLg6vb8duL90rPa0fTXw3VTezMwMqHGIUtI9wExgjKQuYAkwU9I0iiHKfcBHACJit6S1wKPACWBRRLyYjnMjsBEYAayIiN3pFJ8A1kj6a+AnwN0pfjfwD5I6KSa5XFtLO8y680/rmDW/mhJcRFxXJXx3lVil/KeBT1eJbwA2VIk/zstDnOX4r4EPnlFlzcxsWPFSXWZmliUnODMzy5LXorSa+HYAMxuq3IMzM7MsuQdnVgeehWnWeO7BmZlZltyDM+sj98rMmosTnPWZJ5TUxgnSbHB5iNLMzLLkHpxZDdyrNRu63IMzM7MsuQdnVblnMjh6us7+bs6sdu7BmZlZltyDM/fWzCxL7sGZmVmWav3B0xXAHwCHI+KSFLsAuBdoofjB02si4pgkAXcCc4DngT+NiB+n97QDf5EO+9cRsSrFLwNWAmdT/F7cTRERPZ2jlraYDSW+Z86sdrUOUa4EvgCsLsUWA5sjYqmkxen1J4DZwJT0mA4sA6anZLUEaKX4FfDtktanhLUM+DCwlSLBtQEP9HIO64WHIs1sOKlpiDIivg8c7RaeC6xK26uAeaX46ihsAUZJGg9cBWyKiKMpqW0C2tK+8yJiS0QERRKdd5pzmJmZAfX5Dm5cRBxM208B49L2BGB/qVxXivUW76oS7+0cZmZmQJ1nUabvy6JR55C0EFgIMHny5HpWY8jxcKSZDXf16MEdSsOLpOfDKX4AmFQqNzHFeotPrBLv7RyvEBHLI6I1IlrHjh1bU6PMzKy51KMHtx5oB5am5/tL8RslraGYZHI8Ig5K2gh8RtLoVG4WcHNEHJX0nKQZFJNM5gN/d5pzDEvurZmZnarW2wTuAWYCYyR1UcyGXAqslbQAeBK4JhXfQHGLQCfFbQLXA6REdhuwLZW7NSIqE1du4OXbBB5ID3o5h5mZGVBjgouI63rYdWWVsgEs6uE4K4AVVeIdwCVV4s9UO4eZmVmFl+pqIh6KNDPrOye4IcrJzLrz6iZmZ8YJ7gz4p03MzJqHE9wQ4N6amdnAc4IbAB46MjMbepzg6sS9MjOzxvLvwZmZWZac4MzMLEtOcGZmliUnODMzy5InmZg1Md+badYz9+DMzCxLTnBmZpYlJzgzM8uSE5yZmWWpbglO0j5JOyXtkNSRYhdI2iRpb3oeneKSdJekTkmPSLq0dJz2VH6vpPZS/LJ0/M70XtWrLWZm1nzq3YO7IiKmRURrer0Y2BwRU4DN6TXAbGBKeiwElkGRECl+JXw6cDmwpJIUU5kPl97XVue2mJlZExns2wTmAjPT9irge8AnUnx1+tXvLZJGSRqfym6KiKMAkjYBbZK+B5wXEVtSfDUwD3hg0FpiNoR5AXCz+vbgAnhQ0nZJC1NsXEQcTNtPAePS9gRgf+m9XSnWW7yrStzMzAyobw/u3RFxQNLrgU2SflreGREhKep4flJiXQgwefLkep7KzMyGmLr14CLiQHo+DKyj+A7tUBp6JD0fTsUPAJNKb5+YYr3FJ1aJd6/D8ohojYjWsWPHDkSzzMysSdQlwUk6R9LrKtvALGAXsB6ozIRsB+5P2+uB+Wk25QzgeBrK3AjMkjQ6TS6ZBWxM+56TNCPNnpxfOpaZmVndhijHAevSzP2RwNcj4juStgFrJS0AngSuSeU3AHOATuB54HqAiDgq6TZgWyp3a2XCCXADsBI4m2JyiSeYmJnZSSomLuavtbU1Ojo6ajqGf6XbcuOZldbsJG0v3Yr2Cl7JxMzMsuQEZ2ZmWXKCMzOzLDnBmZlZlvyL3mbDmH8R3HLmHpyZmWXJCc7MzLLkIUozO4V/jcBy4B6cmZllyT04M+sT9+qs2TjBmVm/ncnydU6GNtic4MxsULgHaIPNCc7MGsY9QKsnJzgzawpOhnamnODMLDt9TYZOhHlr6gQnqQ24ExgBfDUilja4SmbWRNwrzFvTJjhJI4AvAu8FuoBtktZHxKONrZmZ5ajWHzx2ghx8TZvggMuBzoh4HEDSGmAu4ARnZkOOE+Tga+YENwHYX3rdBUxvUF3MzOqq1gTZaI1I0M2c4E5L0kJgYXr5S0mPNbI+wBjg6QbXIQe+jrXzNaydr+EZ0O1VwwNxDS/qaUczJ7gDwKTS64kpdlJELAeWD2aleiOpIyJaG12PZufrWDtfw9r5Gtau3tewmRdb3gZMkXSxpLOAa4H1Da6TmZkNEU3bg4uIE5JuBDZS3CawIiJ2N7haZmY2RDRtggOIiA3AhkbX4wwMmeHSJufrWDtfw9r5GtaurtdQEVHP45uZmTVEM38HZ2Zm1iMnuEEiqU3SY5I6JS1udH2ahaR9knZK2iGpI8UukLRJ0t70PLrR9RxqJK2QdFjSrlKs6nVT4a702XxE0qWNq/nQ0cM1vEXSgfR53CFpTmnfzekaPibpqsbUeuiQNEnSQ5IelbRb0k0pPmifQye4QVBaVmw2MBW4TtLUxtaqqVwREdNK04kXA5sjYgqwOb22V1oJtHWL9XTdZgNT0mMhsGyQ6jjUreTUawhwR/o8TkvzAEj/nq8F3pbe86X07344OwF8LCKmAjOARek6Ddrn0AlucJxcViwifgNUlhWz/pkLrErbq4B5DazLkBQR3weOdgv3dN3mAqujsAUYJWn84NR06OrhGvZkLrAmIl6IiCeATop/98NWRByMiB+n7V8AeyhWoBq0z6ET3OCotqzYhAbVpdkE8KCk7WllGoBxEXEwbT8FjGtM1ZpOT9fNn88zc2MaQltRGh73NeyFpBbgHcBWBvFz6ARnQ927I+JSiuGLRZJ+t7wzimnAngp8hnzd+m0Z8CZgGnAQ+GxjqzP0SToX+Abw0Yh4rryv3p9DJ7jBcdplxay6iDiQng8D6yiGfQ5Vhi7S8+HG1bCp9HTd/Pnso4g4FBEvRsRLwFd4eRjS17AKSa+mSG5fi4hvpvCgfQ6d4AaHlxXrB0nnSHpdZRuYBeyiuHbtqVg7cH9jath0erpu64H5aRbbDOB4aQjJSrp9J/QBis8jFNfwWkmvkXQxxUSJHw12/YYSSQLuBvZExOdKuwbtc9jUK5k0Cy8r1m/jgHXFvxNGAl+PiO9I2gaslbQAeBK4poF1HJIk3QPMBMZI6gKWAGyCc/IAAAB8SURBVEupft02AHMoJkY8D1w/6BUegnq4hjMlTaMYVtsHfAQgInZLWkvxe5QngEUR8WIj6j2EvAv4E2CnpB0p9kkG8XPolUzMzCxLHqI0M7MsOcGZmVmWnODMzCxLTnBmZpYlJzgzM8uSE5yZmWXJCc7MzLLkBGdmZln6/1tkKd+h3CahAAAAAElFTkSuQmCC\n",
      "text/plain": [
       "<Figure size 432x288 with 2 Axes>"
      ]
     },
     "metadata": {
      "needs_background": "light"
     },
     "output_type": "display_data"
    }
   ],
   "source": [
    "from nyul import nyul_apply_standard_scale\n",
    "import matplotlib.pyplot as plt\n",
    "\n",
    "image_1 = nib.load(train_scans[0]).get_data()\n",
    "mask_1 = nib.load(mask_scans[0]).get_data()\n",
    "image_1_norm = nyul_apply_standard_scale(image_1, standard_path, input_mask=mask_1)\n",
    "\n",
    "fig, axs = plt.subplots(2, 1, constrained_layout=True)\n",
    "f1 = axs[0].hist(image_1.flatten(), bins=64, range=(-10,600))\n",
    "f2 = axs[1].hist(image_1_norm.flatten(), bins=64, range=(-10,200))\n",
    "axs[0].set_title('Image 1 Original')\n",
    "axs[0].set_title('Image 1 Normalized')\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.7.3"
  }
 },
 "nbformat": 4,
 "nbformat_minor": 2
}