{
 "metadata": {
  "name": ""
 },
 "nbformat": 3,
 "nbformat_minor": 0,
 "worksheets": [
  {
   "cells": [
    {
     "cell_type": "heading",
     "level": 1,
     "metadata": {},
     "source": [
      "DKsort ROI Figure Notebook"
     ]
    },
    {
     "cell_type": "markdown",
     "metadata": {},
     "source": [
      "This notebook is used to create the figure displaying the parcellation used in these analyses."
     ]
    },
    {
     "cell_type": "code",
     "collapsed": false,
     "input": [
      "%pylab wx\n",
      "from IPython.display import display, Image\n",
      "import seaborn\n",
      "seaborn.set(context=\"paper\", style=\"nogrid\")"
     ],
     "language": "python",
     "metadata": {},
     "outputs": []
    },
    {
     "cell_type": "code",
     "collapsed": false,
     "input": [
      "import os\n",
      "import os.path as op\n",
      "from surfer import Brain\n",
      "from mayavi import mlab\n",
      "from nipy.labs import viz\n",
      "import nibabel as nib\n",
      "from scipy.ndimage import binary_erosion\n",
      "import lyman\n",
      "import moss"
     ],
     "language": "python",
     "metadata": {},
     "outputs": []
    },
    {
     "cell_type": "code",
     "collapsed": false,
     "input": [
      "views = dict(pfc=dict(lh=[135, 80, 290, [20, 15, -5]],\n",
      "                      rh=[45, 80, 290, [-20, 15, -5]]),\n",
      "             net=dict(lh=[230, 65, 335, [40, 15, -20]],\n",
      "                      rh=[-105, 105, 290, [0, -15, 0]]))"
     ],
     "language": "python",
     "metadata": {},
     "outputs": []
    },
    {
     "cell_type": "code",
     "collapsed": false,
     "input": [
      "opts = dict(width=600, height=600, background=\"white\")\n",
      "mpl.rcParams[\"image.interpolation\"] = \"spline36\""
     ],
     "language": "python",
     "metadata": {},
     "outputs": []
    },
    {
     "cell_type": "markdown",
     "metadata": {},
     "source": [
      "Copy the partially-inflated fsaverage surfaces into the SUBJECTS_DIR."
     ]
    },
    {
     "cell_type": "code",
     "collapsed": false,
     "input": [
      "!cp rois/?h.semi7 $SUBJECTS_DIR/fsaverage/surf/"
     ],
     "language": "python",
     "metadata": {},
     "outputs": []
    },
    {
     "cell_type": "markdown",
     "metadata": {},
     "source": [
      "Use PySurfer to draw the surfaces with the full annotation"
     ]
    },
    {
     "cell_type": "code",
     "collapsed": false,
     "input": [
      "annot_shots = dict(pfc=dict(), net=dict())\n",
      "for hemi in [\"lh\", \"rh\"]:\n",
      "    b = Brain(\"fsaverage\", hemi, \"semi7\", config_opts=opts)\n",
      "    b.add_annotation(\"dksort\", alpha=.8, borders=False)\n",
      "    for view in [\"pfc\", \"net\"]:\n",
      "        mlab.view(*views[view][hemi])\n",
      "        annot_shots[view][hemi] = mlab.screenshot()\n",
      "    b.close()"
     ],
     "language": "python",
     "metadata": {},
     "outputs": []
    },
    {
     "cell_type": "markdown",
     "metadata": {},
     "source": [
      "----"
     ]
    },
    {
     "cell_type": "markdown",
     "metadata": {},
     "source": [
      "Now draw the IFS label on fsaverage"
     ]
    },
    {
     "cell_type": "code",
     "collapsed": false,
     "input": [
      "label_shots = dict()\n",
      "b = Brain(\"fsaverage\", \"rh\", \"semi7\", config_opts=opts)\n",
      "b.add_label(\"rois/rh.yeo17_ifs.label\", color=\"#CC3333\", borders=False, alpha=.8)\n",
      "mlab.view(*views[\"pfc\"][\"rh\"])\n",
      "label_shots[\"fsaverage\"] = mlab.screenshot()\n",
      "b.close()"
     ],
     "language": "python",
     "metadata": {},
     "outputs": []
    },
    {
     "cell_type": "markdown",
     "metadata": {},
     "source": [
      "----"
     ]
    },
    {
     "cell_type": "markdown",
     "metadata": {},
     "source": [
      "Make a temporary directory for files we'll create as we go and move into it"
     ]
    },
    {
     "cell_type": "code",
     "collapsed": false,
     "input": [
      "import os\n",
      "import tempfile\n",
      "tempdir = tempfile.mkdtemp()\n",
      "orig_dir = os.getcwd()\n",
      "os.chdir(tempdir)"
     ],
     "language": "python",
     "metadata": {},
     "outputs": []
    },
    {
     "cell_type": "markdown",
     "metadata": {},
     "source": [
      "Transform the annotation to the first subject"
     ]
    },
    {
     "cell_type": "code",
     "collapsed": false,
     "input": [
      "%%bash\n",
      "mri_surf2surf --srcsubject fsaverage --trgsubject dk02 --hemi rh \\\n",
      "              --sval-annot $SUBJECTS_DIR/fsaverage/label/rh.dksort.annot --tval ./rh.dksort.annot &> /dev/null\n",
      "mri_annotation2label --subject dk02 --hemi rh --annotation ./rh.dksort.annot --outdir . &> /dev/null"
     ],
     "language": "python",
     "metadata": {},
     "outputs": []
    },
    {
     "cell_type": "markdown",
     "metadata": {},
     "source": [
      "Now plot the label on the native subject surface"
     ]
    },
    {
     "cell_type": "code",
     "collapsed": false,
     "input": [
      "b = Brain(\"dk02\", \"rh\", \"inflated\", config_opts=opts)\n",
      "b.add_label(\"rh.IFS.label\", color=\"#CC3333\", borders=False, alpha=.8)\n",
      "view = [50, 80, 380, [-30, 15, -5]]\n",
      "mlab.view(*view)\n",
      "label_shots[\"native\"] = mlab.screenshot()\n",
      "b.close()"
     ],
     "language": "python",
     "metadata": {},
     "outputs": []
    },
    {
     "cell_type": "markdown",
     "metadata": {},
     "source": [
      "----"
     ]
    },
    {
     "cell_type": "markdown",
     "metadata": {},
     "source": [
      "Next transform the fsaverage curvature onto the native space"
     ]
    },
    {
     "cell_type": "code",
     "collapsed": false,
     "input": [
      "%%bash\n",
      "mri_surf2surf --srcsubject fsaverage --trgsubject dk02 --hemi rh \\\n",
      "              --sval $SUBJECTS_DIR/fsaverage/surf/rh.curv \\\n",
      "              --tval ./rh.curv.dk02.mgh &> /dev/null"
     ],
     "language": "python",
     "metadata": {},
     "outputs": []
    },
    {
     "cell_type": "markdown",
     "metadata": {},
     "source": [
      "Plot the reverse-normalized curvature on the native sphere"
     ]
    },
    {
     "cell_type": "code",
     "collapsed": false,
     "input": [
      "b = Brain(\"dk02\", \"rh\", \"sphere\", config_opts=opts)\n",
      "b.add_label(\"rh.IFS.label\", color=\"#CC3333\", borders=False, alpha=.6)\n",
      "curv = nib.load(\"rh.curv.dk02.mgh\").get_data().squeeze()\n",
      "curv = (curv > 0).astype(int)\n",
      "b.add_contour_overlay(curv, 0, .75, 2, 5)\n",
      "b.contour[\"colorbar\"].lut_mode = \"copper\"\n",
      "b.contour[\"colorbar\"].visible = False\n",
      "view = [70, 70, 400, [0, 0, 0], 270]\n",
      "mlab.view(*view)\n",
      "label_shots[\"sphere\"] = mlab.screenshot()\n",
      "b.close()"
     ],
     "language": "python",
     "metadata": {},
     "outputs": []
    },
    {
     "cell_type": "markdown",
     "metadata": {},
     "source": [
      "-----"
     ]
    },
    {
     "cell_type": "markdown",
     "metadata": {},
     "source": [
      "Prepare the data to plot the IFS mask on the mean functional image"
     ]
    },
    {
     "cell_type": "code",
     "collapsed": false,
     "input": [
      "os.chdir(orig_dir)"
     ],
     "language": "python",
     "metadata": {},
     "outputs": []
    },
    {
     "cell_type": "code",
     "collapsed": false,
     "input": [
      "project = lyman.gather_project_info()\n",
      "bet_file = op.join(project[\"analysis_dir\"], \"dksort/dk02/preproc/run_1/functional_mask.nii.gz\")\n",
      "mean_file = op.join(project[\"analysis_dir\"], \"dksort/dk02/preproc/run_1/mean_func.nii.gz\")\n",
      "mask_file = op.join(project[\"data_dir\"], \"dk02/masks/yeo17_ifs.nii.gz\")"
     ],
     "language": "python",
     "metadata": {},
     "outputs": []
    },
    {
     "cell_type": "code",
     "collapsed": false,
     "input": [
      "z = 17\n",
      "mean_data = nib.load(mean_file).get_data()[..., z].T\n",
      "bet_data= nib.load(bet_file).get_data()[..., z].T\n",
      "vmin = moss.percentiles(mean_data[mean_data > 0], 6)\n",
      "vmax = moss.percentiles(mean_data[mean_data > 0], 98)\n",
      "mean_data[~binary_erosion(bet_data, iterations=4)] = np.nan\n",
      "mask_data = nib.load(mask_file).get_data()[..., z].T\n",
      "mask_data = np.where(mask_data, 1, np.nan)"
     ],
     "language": "python",
     "metadata": {},
     "outputs": []
    },
    {
     "cell_type": "markdown",
     "metadata": {},
     "source": [
      "----"
     ]
    },
    {
     "cell_type": "markdown",
     "metadata": {},
     "source": [
      "Now stitch everything together into one image"
     ]
    },
    {
     "cell_type": "markdown",
     "metadata": {},
     "source": [
      "First stick each annotation view into a single array"
     ]
    },
    {
     "cell_type": "code",
     "collapsed": false,
     "input": [
      "annot_panels = [annot_shots[\"pfc\"][\"rh\"],\n",
      "                annot_shots[\"pfc\"][\"lh\"],\n",
      "                annot_shots[\"net\"][\"lh\"][:, 25:-25],\n",
      "                annot_shots[\"net\"][\"rh\"][:, 100:-100]]\n",
      "full_annot = np.concatenate(annot_panels, axis=1)\n",
      "shape = [105]\n",
      "shape.extend(full_annot.shape[1:])\n",
      "full_annot = np.concatenate([full_annot, np.ones(shape) * 255], axis=0)"
     ],
     "language": "python",
     "metadata": {},
     "outputs": []
    },
    {
     "cell_type": "markdown",
     "metadata": {},
     "source": [
      "Find the ROI colors from the lookup table"
     ]
    },
    {
     "cell_type": "code",
     "collapsed": false,
     "input": [
      "lut = np.genfromtxt(\"rois/dksort.lut\", object)\n",
      "rois = lut[1:, 1]\n",
      "colors = np.array(lut[1:, 2:5], np.float) / 255."
     ],
     "language": "python",
     "metadata": {},
     "outputs": []
    },
    {
     "cell_type": "markdown",
     "metadata": {},
     "source": [
      "Plot everything together"
     ]
    },
    {
     "cell_type": "code",
     "collapsed": false,
     "input": [
      "f = figure(figsize=(6.85, 4), facecolor=\"white\")\n",
      "\n",
      "annot_ax = f.add_axes([.01, .45, .98, .55])\n",
      "annot_ax.axis(\"off\")\n",
      "annot_ax.imshow(full_annot / 255)\n",
      "xlim, ylim = annot_ax.get_xlim(), annot_ax.get_ylim()\n",
      "for roi, color in zip(rois, colors):\n",
      "    annot_ax.bar(max(xlim) + 1, 0, color=color, label=roi, edgecolor=color)\n",
      "annot_ax.set_xlim(*xlim)\n",
      "annot_ax.set_ylim(*ylim)\n",
      "annot_ax.legend(loc=\"lower left\", bbox_to_anchor=(0.001, 0.001),\n",
      "                ncol=len(rois), frameon=False, fontsize=8)\n",
      "\n",
      "group_ax = f.add_axes([.02, .1, .26, .35])\n",
      "group_ax.axis(\"off\")\n",
      "group_ax.imshow(label_shots[\"fsaverage\"])\n",
      "\n",
      "sphere_ax = f.add_axes([.28, .1, .26, .35])\n",
      "sphere_ax.axis(\"off\")\n",
      "sphere_ax.imshow(label_shots[\"sphere\"])\n",
      "\n",
      "native_ax = f.add_axes([.54, .1, .26, .35])\n",
      "native_ax.axis(\"off\")\n",
      "native_ax.imshow(label_shots[\"native\"])\n",
      "\n",
      "epi_ax = f.add_axes([.75, .05, .3, .4])\n",
      "epi_ax.axis(\"off\")\n",
      "epi_ax.imshow(mean_data, cmap=\"gray\", vmin=vmin, vmax=vmax)\n",
      "cmap = mpl.colors.ListedColormap([\"#CC3333\"])\n",
      "epi_ax.imshow(mask_data, cmap=cmap, vmin=0, vmax=1.5, interpolation=\"nearest\")\n",
      "\n",
      "arrow_ax = f.add_axes([0, 0, 1, .1])\n",
      "arrow_ax.axis(\"off\")\n",
      "arrow_ax.arrow(.15, .6, .64, 0, fc=\"k\", ec=\"k\", head_width=.17, head_length=.04)\n",
      "text_kws = dict(backgroundcolor=\"w\", color=\"k\", size=8, ha=\"center\", va=\"center\")\n",
      "arrow_ax.text(.15, .6, \" Common Surface \", **text_kws)\n",
      "arrow_ax.text(.41, .6, \" Spherical Registration \", **text_kws)\n",
      "arrow_ax.text(.67, .6, \" Native Surface \", **text_kws)\n",
      "arrow_ax.text(.9, .6, \" Native Volume \", **text_kws)\n",
      "\n",
      "f.text(.02, .94, \"A\", size=12)\n",
      "f.text(.02, .42, \"B\", size=12)\n",
      "\n",
      "f.savefig(\"figures/figure_3.pdf\", dpi=300)\n",
      "f.savefig(\"figures/figure_3.png\", dpi=300)\n",
      "\n",
      "plt.close(f)"
     ],
     "language": "python",
     "metadata": {},
     "outputs": []
    },
    {
     "cell_type": "code",
     "collapsed": false,
     "input": [
      "display(Image(\"figures/figure_3.png\"))"
     ],
     "language": "python",
     "metadata": {},
     "outputs": []
    },
    {
     "cell_type": "code",
     "collapsed": false,
     "input": [],
     "language": "python",
     "metadata": {},
     "outputs": []
    }
   ],
   "metadata": {}
  }
 ]
}