F-Praktikum-Anleitung-Fortgeschritten2.ipynb

{
 "cells": [
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "# Versuch V-So 15 – Maschinelles Lernen in der wissenschaftlichen Bildanalyse\n",
    "\n",
    "# Einleitung\n",
    "\n",
    "Ziel dieses Versuchs ist es die physikalische Eigenschaften\n",
    "- *Position*\n",
    "- *lokale Dichte*\n",
    "- *Volumen*\n",
    "- *Intensität*\n",
    "- usw.\n",
    "\n",
    "von Zellen mit Hilfe von tiefen neuronalen Netzwerken aus Mikroskopaufnahmen zu extrahieren. Der Versuch ist im Wesentlich in vier Phasen geteilt:\n",
    "1. Vorbereitung der Trainingsdaten und Visualisierung\n",
    "2. Testen von verschieden Netzwerkdesigns\n",
    "3. Anwedung von einem Design auf noch unbekannte Daten\n",
    "3. Berechnung der physikalischen Eigenschaften\n",
    "\n",
    "# Hinweise\n",
    "\n",
    "Dies ist keine Prüfungssituation!\n",
    "- Zu jedem Zeitpunkt kann eine [Suchmaschiene](https://en.wikipedia.org/wiki/Web_search_engine) euer Vorankommen beschleunigen.\n",
    "- Wenn ein Problem auftritt, dann sucht nach Beispiele und Lösungen auf [Github](https://github.com/), [StackOverflow](https://stackoverflow.com/) und in offizielle [Dokumentationen](https://matplotlib.org/gallery/index.html).\n",
    "- Es erwartet keiner, dass ihr alle Python [Befehle](https://docs.python.org/3/library/) auswendig kennt. Es wird aber erwartet, dass ihr in kurzer Zeit euch neue Befehle heraussuchen könnt.\n",
    "- Wenn ihr Fragen habt, lasst es uns wissen.\n",
    "- Speichert das Notebook, wann immer ihr es für notwendig haltet.\n",
    "\n",
    "# Vorbereitung\n",
    "\n",
    "## Macht euch mit der Maschine vertraut:\n",
    "\n",
    "Die für die Numerik sind folgende Komponenten von Bedeutung:\n",
    "- CPU (Model)\n",
    "- Arbeitsspeicher (Speicherplatz)\n",
    "- Speicher (Model (von nvme0n1))\n",
    "- GPU (Model)\n",
    "\n",
    "diese Komponenten entscheiden meist darüber ob die vorgeschlagenen Rechnungen auch ausgeführt werden können und sollten im Protokoll Erwähnungen finden. Die folgenden Zellen machen es euch einfach alles notwendige herauszufinden"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "!cat /proc/cpuinfo | \\\n",
    "awk -v FS=':' '                                       \\\n",
    "  /^physical id/ { if(nb_cpu<$2)  { nb_cpu=$2 } }     \\\n",
    "  /^cpu cores/   { if(nb_cores<$2){ nb_cores=$2 } }   \\\n",
    "  /^processor/   { if(nb_units<$2){ nb_units=$2 } }   \\\n",
    "  /^model name/  { model=$2 }                         \\\n",
    "                                                      \\\n",
    "  END{                                                \\\n",
    "   nb_cpu=(nb_cpu+1);                                 \\\n",
    "   nb_units=(nb_units+1);                             \\\n",
    "                                                      \\\n",
    "   print \"CPU model:\",model;                          \\\n",
    "   print nb_cpu,\"CPU,\",nb_cores,\"physical cores per CPU, total\",nb_units,\"logical CPU units\" \\\n",
    " }' # Quelle: https://superuser.com/questions/388115/interpreting-output-of-cat-proc-cpuinfo"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "!cat /proc/meminfo | grep MemTotal"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "!cat /sys/class/block/nvme0n1/device/model"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "!nvidia-smi --query-gpu=gpu_name,memory.total --format=csv"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "## Macht euch mit der Daten vertraut\n",
    "\n",
    "Im Gegensatz zu vielen anderen Versuche im Rahmen des Fortgeschrittenenpraktikums, ist ein Großteil der für diesen Versuch notwendigen Daten bereits vorhanden. Ihr findet die Daten unter:\n",
    "\n",
    "**Traininglabels (Ground Truth)** - Existierende Binärmasken.\n",
    "\n",
    "`/DATA/GT/<Experiment Name>/`\n",
    "\n",
    "**Traininginput** - Aufnahmen auf denen die Binärmasken basieren.\n",
    "\n",
    "`/DATA/IN/<Experiment Name>/<Parameter>/`\n",
    "\n",
    "**Anwendungsdatensätze** - Datensätze dir ihr analysieren sollt.\n",
    "\n",
    "`/DATA/validation/noise1000/`\n",
    "\n",
    "**Notebook & Benutzerdaten** - Speicherort für *.ipynb* und allen Daten die ihr im Laufe des Versuchs generiert.\n",
    "\n",
    "`/tf/Notebooks`\n",
    "\n",
    "#### Aufgabe:\n",
    "- Lasst euch den Inhalt von den Verzeichnissen anzeigen.\n",
    "- Stellt je 2 verschiedene Ebenen von 2 verschiedenen Trainingslabels und -Inputpaare in einem Diagrammen dar.\n",
    "\n",
    "Hinweise:\n",
    "\n",
    "- Um eine Zelle auszuführen müsst ihr `Shift+Enter` drücken.\n",
    "\n",
    "- Mit vorangestellten `!` kann man in jupyter beliebige Programme auf dem Computer ausführen.\n",
    "\n",
    "- Unter Linux/ MacOSX kann man sich den Inhalt von Verzeichnissen mit dem Programm `ls` und den darauf folgenden Verzeichnissnamen anzeigen lassen.\n",
    "\n",
    "- Nutzt die unten vorgeschriebene Funktion `loadImage`.\n",
    " \n",
    "- Zum Anzeigen von 2D Matrizen kann die Funktion [`imshow`](https://matplotlib.org/gallery/images_contours_and_fields/image_demo.html#sphx-glr-gallery-images-contours-and-fields-image-demo-py) des `matplotlib.pyplot`-[Modules](https://matplotlib.org/gallery/index.html) verwendet werden."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "%matplotlib notebook\n",
    "\n",
    "import numpy as np # importiert das numpy Modul unter den Namen 'np'\n",
    "import imread # importiert imread ohne eignenen Namen\n",
    "\n",
    "def loadImage(filepath):\n",
    "    \"\"\"Liest tiff-Dateien von der Festplatte als numpy array in den Arbeitsspeicher. \"\"\"\n",
    "    return np.asarray(imread.imload_multi(filepath))"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "#### Lösung"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "### Limitiere CUDA auf eine GPU\n",
    "\n",
    "Im Verlauf dieses Versuchs werden wir Grafikkarten benutzen um Tensoroperationen durchzuführen. Dies soll mithilfe des [*tensorflow*](https://www.tensorflow.org/) Python Moduls passieren. Standardmäßig nutzt *tensorflow*  fast den gesamten verfügbare Grafikspeicher.\n",
    "\n",
    "Die momentane Auslastung der Grafikkarten auf den System lassen sich mit dem Programm *nvidia-smi* überprüfen."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "!nvidia-smi"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Über die [Umgebungsvariable](https://de.wikipedia.org/wiki/Umgebungsvariable) *CUDA_VISIBLE_DEVICES* kann man die Geräte einschränken, die für [CUDA](https://de.wikipedia.org/wiki/CUDA) Berechnungen benutzt werden dürfen.\n",
    "\n",
    "Hinweis: \n",
    "- `%env` ist *jupyter* [*Magie*](https://ipython.readthedocs.io/en/stable/interactive/magics.html) um Umgebungsvariablen zu setzen oder anzeigen zu lassen."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "%env CUDA_DEVICE_ORDER=PCI_BUS_ID\n",
    "%env CUDA_VISIBLE_DEVICES=0"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "### Vordefinierte Funktionen\n",
    "\n",
    "Um schneller eine Überblick über die eingelesen 3D-Bilder als auch über die vorhergesagten Segmentierungen zu erhalten, haben wir die Funktion `zSlicer` vorbereitet. Über den Slider unter der Funktion könnt ihr die aktuell dargestellte z-Ebene verändern. Dabei unterstützt der zSlicer eine beliebige Anzahl an 2D oder 3D Arrays als Input."
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "#### Aufgabe\n",
    "1. Lest eine *.tif*-Datei ein und nutzt den Slider von zSlicer um durch das Volumen zu scrollen."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "%matplotlib notebook\n",
    "import ipywidgets as widgets\n",
    "import matplotlib.pyplot as plt\n",
    "\n",
    "def zSlicer(*args):\n",
    "    \"\"\"Zeichnet interactive imshow Diagramme von zwei 3D oder 2D numpy Arrays.\"\"\"\n",
    "    num_img = len(args)\n",
    "            \n",
    "    num_rows = 1\n",
    "    num_cols = len(args)\n",
    "    \n",
    "    # init plot\n",
    "    fig, axes = plt.subplots(num_rows, num_cols, figsize=[9.5,5])\n",
    "    \n",
    "    val = 0\n",
    "    aspect = 'equal'\n",
    "    \n",
    "    plots = []\n",
    "    for i, img in enumerate(args):\n",
    "        if len(args) == 1:\n",
    "            ax = axes\n",
    "        else:\n",
    "            ax = axes[i]\n",
    "        \n",
    "        if img.ndim > 2:\n",
    "\n",
    "            plots.append([ax.imshow(img[val], aspect=aspect, vmin=np.amin(img[1:]), vmax=np.amax(img[1:])), img])\n",
    "\n",
    "        else:\n",
    "            ax.imshow(img, aspect=aspect, vmin=np.amin(img), vmax=np.amax(img))\n",
    "            \n",
    "\n",
    "    def update_plot(change):\n",
    "        val = change['new']\n",
    "        for plot,img in plots:\n",
    "            plot.set_array(img[val])\n",
    "            \n",
    "        return\n",
    "\n",
    "\n",
    "    slider = widgets.IntSlider(value=val, max=args[-1].shape[0]-1)\n",
    "    display(slider)\n",
    "\n",
    "    slider.observe(update_plot, names='value')\n",
    "        "
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "#### Lösung"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "## Keras Einführung\n",
    "\n",
    "[Keras](https://keras.io/) ist eine high-level Deep Learning Bibliothek, die es euch erlaubt einfach eigene Designs für \n",
    "neuronale Netzwerke zu basteln, ohne das ihr *TensorFlow* oder gar *CUDA* Code schreiben müsst. Dabei kann Keras mit einer Vielzahl von Deep Learning Backends benutzt werden. Wir werden es heute ausschließlich mit dem Tensorflowbackend benutzen.\n",
    "\n",
    "Am Ende dieses Versuchs wollen wir Schritt für Schritt einen [Autoencoder](http://science.sciencemag.org/content/313/5786/504/tab-pdf) erstellen, der dem Design von [UNet](http://arxiv.org/abs/1505.04597) entspricht. Dazu müssen wir uns erst einmal mit den Grundlagen vertraut machen.\n",
    "\n",
    "### Die Keras Model API\n",
    "\n",
    "Keras bezeichnet jedes neuronale Netzwerk unabhängig von der Architektur als [*Model*](https://keras.io/models/model/). Ein Model $M$ besteht aus ein oder mehreren Ebenen (Layern), die nichts anders als hintereinander geschaltet Funktionen darstellen. \n",
    "$$M(x) = (O \\circ \\ldots \\circ L_2 \\circ L_1 \\circ I)(x) $$\n",
    "Die Model API von Keras kommt mit einer Vielzahl von nützlichen Befehlen, die z.B. das Übersetzung in CUDA-Code, das Training oder aber die Anwendung des Models auf neue Probleme sehr vereinfacht.\n",
    "\n",
    "#### Beispiel: Minimales Keras Model\n",
    "\n",
    "Ein minimales Model hat lediglich einen Input $I$ und einen Output $O$. In dem unten gegebenen Fall wird über den Input ein *numpy* Array mit den Form ($64 \\times 64 \\times 1$) als Input angenommen und anschließend die Identität $\\mathbb{1}:x \\mapsto x$ angewandt."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "from tensorflow.keras.layers import Input, Lambda\n",
    "from tensorflow.keras.models import Model\n",
    "\n",
    "input_shape = (64, 64, 1)\n",
    "\n",
    "inputs = Input(input_shape)\n",
    "outputs = Lambda(lambda x: x)(inputs)\n",
    "\n",
    "simple_model = Model(inputs=inputs, outputs=outputs)\n",
    "simple_model.summary()"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Um dieses Model nutzen zu können, müssen wir es für die Ausführung *kompilieren* d.h. unser Python Code wird im Hintergrund in Maschienencode umgewandelt. In Keras muss man dafür die Optimierungsfunktion (hier: [Stochastic gradient descent](https://en.wikipedia.org/wiki/Stochastic_gradient_descent)) und eine Verlustfunktion (hier: Mittlere quadratische Abweichung) angeben."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "from tensorflow.keras.optimizers import SGD\n",
    "\n",
    "simple_model.compile(optimizer=SGD(), loss='mean_squared_error')"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Der nächste Schritt wäre schon der Trainingsprozess. Da aber keine freien Parameter existieren (vgl. Ausgabe oben) können wir uns direkt anschauen, wie wir das Model für Vorhersagen benutzen können.\n",
    "\n",
    "Dazu erstellen wir ein Schachbrettmuster im Format $(64 \\times 64)$ und verändern anschließend die Anzahl der Dimensionen damit der Input die richtige Form für das oben definierte Modell hat. Anschließend nutzen wir das einfache Modell um eine Vorhersage zu treffen."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "shape = (64, 64)\n",
    "test_input = np.indices(shape).sum(axis=0) % 2 # https://stackoverflow.com/a/51715491\n",
    "test_input =  test_input.reshape((1,*input_shape))\n",
    "print(test_input.shape)\n",
    "result = simple_model.predict(test_input)\n",
    "print(result.shape)"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Mit der Funktion `zSlicer` können wir den Input und das Resultat der Vorhersage anzeigen lassen."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "zSlicer(test_input[0, ..., 0], result[0, ..., 0])"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Diese Minimalbeispiel ist nicht wirklich nützlich um Binärmasken aus Mikroskopaufnahmen zu generieren. Deshalb wollen wir nun ein Netzwerk mit mehreren hintereinander geschalteten [Konvolutionen](https://towardsdatascience.com/intuitively-understanding-convolutions-for-deep-learning-1f6f42faee1) erstellen."
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "#### Aufgabe:\n",
    "- Erstellt ein neues Netzwerk (Name ist euch überlassen) mit den oben benutzen Input Objekt und 4 drauf folgende 2D [Konvolutionen](https://keras.io/layers/convolutional/) mit dem Eigenschaften:\n",
    "    - Die ersten drei Konvolutionen sollen 64 Filtern besitzen, die letzte Konvolution 32\n",
    "    - Einer Kernelgröße von 3 und der [hier](https://arxiv.org/abs/1502.01852) beschriebenen Initalisierung\n",
    "    - Einer \"rectified linear unit\" als Aktivierungsfunktionen\n",
    "    - Jeder Layer soll den den gesamten Input falten und nicht nur den zulässigen Input\n",
    "    \n",
    "- Bastelt eine Funktion die als Parameter den `input_shape` mit den oben gegebenen Wert als Standardwert annimmt und ein komplettes (noch nicht kompiliertes) Modell als Rückgabewert.\n",
    "- Lasst euch die Zusammenfassung des Netzwerks ausgeben, kompiliert das Netzwerk mit den oben angegebenen Parametern und lasst euch eine Vorhersage von dem so erstellten Netzwerk ausgeben.\n",
    "- Bestimmt die Form vom Resultat.\n",
    "\n",
    "Hinweis:\n",
    "- Ihr müsst (abgesehen von der geforderten) keine einzige Funktion selbst  schreiben.\n",
    "- Links sind dafür da um benutzt zu werden."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "#### Lösung"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Fällt euch bzgl. der Form des Outputs etwas auf?\n",
    "\n",
    "#### Aufgabe\n",
    "- modifiziert die vierte Konvolution und ergänzt eine fünfte und sechste Konvolution in der oben erstellte Funktion:\n",
    "    - Die vierte Konvolution soll jetzt 64 Filter haben\n",
    "    - Die fünfte Konvolution soll bis auf die Filteranzahl (2) identisch mit den vorherhigen Konvolutionen sein. \n",
    "    - Die neue Ausgabekonvolution soll bis auf die Kernelgröße von 1 und der Aktivierungsfunktion eine `2DConv` mit den Standardwerten darstellen. Die Filteranzahl soll von euch so gewählt werden, dass ihr ähnlich in den Minimalbeispiel 2D-Bilder mit der gleichen Form wie die Inputbilder erhaltet. Bitte verwendet als Aktivierungsfunktion:\n",
    "    $$f(x)= \\frac{1}{1 +\\mathrm{e}^{-x}}$$\n",
    "    \n",
    "    - Überprüft das Resultat eures Netzwerks mit `zSlicer` (wie im Minimalbeispiel)\n",
    "    \n",
    "Hinweise:\n",
    "- Ihr müsst NICHT eine Aktivierungsfunktion selbst definieren.\n",
    "- Benutzt Kopieren und Einfügen wann immer möglich."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "#### Lösung"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Ihr solltet in den letzten Beispielen erkannt haben, dass die Anzahl der freien Parameter (engl. *Trainable params*) immer weiter zugenommen haben. Nun stellt sich die Frage: Wie trainieren wir eigentlich ein solches Netzwerk?\n",
    "\n",
    "# Trainingsdaten einlesen\n",
    "\n",
    "Wie oben bereits erwähnt, liegen die Labeldaten unter `\\DATA\\GT` und zugehörige Inputdaten unter `\\DATA\\IN`. Um diese in den Arbeitsspeicher zu laden muss man für jede einzelne Datei den Dateipfad angeben.\n",
    "\n",
    "Weil ein manuelles Abschreiben aller Dateipfade den Zeitrahmen diese Versuchs sprengen würde, lassen wir uns alle Dateipfade ausgeben, die einem bestimmten Muster entsprechen. Dazu benutzen wir die `glob`-Funktion des gleichnamigen Moduls.\n",
    "\n",
    "#### Beispiel"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "import glob\n",
    "\n",
    "filelist_labels = glob.glob('/DATA/GT/*/*')\n",
    "for filepath in filelist_labels:\n",
    "    print(filepath)"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Das Zeichen `*` bedeutet hier eine beliebige Anzahl beliebiger Zeichen. D.h, wir erhalten mit einem Befehl alle Dateien in der zweiten Unterverzeichnissebene unter `DATA/GT/`."
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "#### Aufgabe\n",
    "\n",
    "- Erstelle eine Pfadliste mit allen (.tif) Dateien in `/DATA/IN`"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "#### Lösung"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "### Einfache Pfadoperationen\n",
    "\n",
    "Operationen zur Modifikation von Dateipfaden oder Pfade insgesammt finden sich im Modul [`os.path`](https://docs.python.org/3.7/library/os.path.html) in der Python Standard Library. Wir benötigen in diesem Versuch:\n",
    "- `basename` - extrahiert einen Dateiname (bzw. Namen des letzten Unterverzeichnisses, wenn ein Verzeichniss angegeben ist).\n",
    "- `dirname` - gibt den Namen des Elternverzeichnis des gegebenen Pfades zurück.\n",
    "- `isfile` - überprüft ob eine Datei existiert und gibt *True* oder *False* zurücl.\n",
    "- `join` - fügt mehrere Zeichenketten zu einen Pfad zusammen"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "#### Aufgabe\n",
    "\n",
    "- Erstelle zwei Dateilisten. Eine für alle Inputdateien und eine für die zugehörigen Labeldateien.\n",
    "- Lasst euch die Anzahl aller vollständigen Input-Label-Paare ausgeben.\n",
    "\n",
    "Hinweis: Es gibt weniger Input Dateien als Label Dateien, *Experiment Name* und *Dateiname* sind bei zueinandergehörigen Dateien identisch."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "#### Lösung"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "### Bilderausschnitte exrahieren\n",
    "\n",
    "Ihr solltet bereits festgestellt haben, dass die Form der 3D Volumen nicht die Form von 2-dimensionalen $64 \\times 64$ Bilder haben, die wir für die Vorhersage und dem Training des Netzwerks einsetzen wollen. Wir müssen aus den Bildstacks einzelne Bilder herausschneiden. Da wir die 3D Bilder als Numpy Arrays einlesen können wir Indexing benutzen."
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "#### Aufgabe \n",
    "- Sucht euch eine beliebiges Input/ Label Paar aus der Pfadliste aus.\n",
    "- Wählt eine Positon (hier definiert durch den $z$, $y$, $x$ Index) [zufällig](http://lmgtfy.com/?q=random+integer+numpy+%20in+interval) aus dem Volumen aus.\n",
    "- Schneidet jeweils ein Bild mit der Form $64\\times64$ heraus und zweigt sie in einem Plot an.\n",
    "- Bringe die jeweiligen Bilder auf die Form (1, 64, 64, 1) und vergleicht das Resultat des Letzten von euch erstellten Netzwerkes mit der tatsächlichen Binärmaske.\n",
    "\n",
    "Hinweise: \n",
    "- Sowohl die Inputs als auch die Labels haben eine Übersichtsebene in z == 0.\n",
    "- Falls ihr keinen Bildausschnitt erwischt mit Zellen, nutzt den zSlicer um euch das Labelbild anzeigen zu lassen und wählt manuell eine Position aus, die Zellen enthält"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "#### Lösung"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Wie erwartet ist die Vorhersage des Netzwerks sehr schlecht ...\n",
    "\n",
    "#### Aufgabe\n",
    "- Erzeugt euch zwei numpy Array mit der Form (100, 64, 64, 1) je aus einem Label und Inputbild (Batch)\n",
    "- Benutzt die oben definierte Funktion `plotBatches` um euch die ersten 10 der 100 Minbatches (Untereinheit von einem Batch) anzeigen zu lassen.\n",
    "- Benutzt die unten vorgegebene Funktion `cropToCells` um die Volumen nach Einlesen von der Festplatte auf einen Würfel zu reduzieren in denen sich Zellen finden lassen."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "def cropToCells(im_x, im_y):\n",
    "    \n",
    "    x_vals = np.any(im_y, axis=(0,1))\n",
    "    y_vals = np.any(im_y, axis=(0,2))\n",
    "    z_vals = np.any(im_y, axis=(1,2))\n",
    "    x_occp = np.where(x_vals)\n",
    "    y_occp = np.where(y_vals)\n",
    "    z_occp = np.where(z_vals)\n",
    "\n",
    "    # cut both volumes accordingly\n",
    "    im_y_ = im_y[np.amin(z_occp):np.amax(z_occp), np.amin(y_occp):np.amax(y_occp), np.amin(x_occp):np.amax(x_occp)]\n",
    "    im_x_ = im_x[np.amin(z_occp):np.amax(z_occp), np.amin(y_occp):np.amax(y_occp), np.amin(x_occp):np.amax(x_occp)]\n",
    "\n",
    "    return im_x_, im_y_\n"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "#### Lösung"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "### Training von einem einfachen Netzwerk\n",
    "\n",
    "#### Aufgabe\n",
    "- [Trainiert](https://keras.io/models/sequential/#fit) das oben definierte Netzwerk mit dem soeben generierten Batch. Dafür\n",
    "    - Kompiliert es neu mit *Stochastic gradient descent* als Optimierungsfunktion, dem Metrikparameter `accuracy` und als Verlustfunktion die `binary_crossentropy`.\n",
    "    - Startet das Training auf 5 Epochen, wobei 10 Prozent der Batches für die Validierung benutzt werden sollen."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "#### Lösung"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Glückwunsch Ihr habt ein erstes einfaches Netzwerk erfolgreich trainiert!\n",
    "\n",
    "#### Aufgabe\n",
    "- Erstellt ein neues Batch, welches auf einem anderen Input/ Label-Paar basiert, als das Trainingsbatch und schaut euch die ersten 10 [Vorhersagen](https://keras.io/models/model/#predict) eures Netzwerks an.\n",
    "- Vergleicht das Ergebnis mit der tatäschlichen Binärmaske."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "#### Lösung"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "# Bildgenerator\n",
    "\n",
    "Wir wollen heute die von uns erstellten Netzwerke nicht nur auf einzelnen Bildern oder einzelnen Experimenten trainieren, sondern auf allen verfügbaren. Der komplette Datensatz lässt sich aber nicht in einem einzelnen *numpy*-Array zusammenfassen, der Arbeitsspeicher reicht nicht aus.\n",
    "\n",
    "Um trotzdem die Netzwerke mit allen verfügbaren Daten zu trainieren, benutzt Keras (übrigens genauso wie Python) Generatoren. Generatoren erzeugen Daten erst wenn sie aufgerufen werden.\n",
    "\n",
    "#### Aufgabe:\n",
    " - Ergänzt den unten teilweise vorgeschriebenen Generator und Normierungsfunktionen\n",
    " - Testet den Generator mit den weiter unten stehenden Test-Code.\n",
    " \n",
    "Hinweis:\n",
    "- Startet mit den Normierungen, ergänzt dann die Funktionen `__init__`, `__getimages__` und zum Schluss `__getitem__`."
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "#### Normalizations"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "import numpy as np\n",
    "\n",
    "def zeroMeanUnitVariance(a, axis=None):\n",
    "    \"\"\"\n",
    "    # Aufgabe: Schreibt eine eigene Norm-Funktion die das Inputarray 'a' auf gegebenen Achsen auf einen\n",
    "    Durchschnittwert von 0 und einer Varianz von 1 normiert.\n",
    "    \n",
    "    Hinweis: Benutzt numpy Funktionen die als Schlüsselwörter die Achse akzeptieren und \n",
    "    für leichtere Implementierung die Dimensionen des Inputs erhalten.\n",
    "    \n",
    "    \"\"\"\n",
    "    \n",
    "    return a\n",
    "\n",
    "def zeroToOne(a, axis=None):\n",
    "    \"\"\"\n",
    "    # Aufgabe: Schreibt eine  Norm-Funktion die das Inputarray 'a' auf gegebenen Achsen auf das Intervall\n",
    "    [0, 1] normiert.\n",
    "    \n",
    "    Hinweis: Benutzt numpy Funktionen!\n",
    "    \n",
    "    \"\"\"\n",
    "    a = a - np.amin(a, axis=axis, keepdims=True)\n",
    "    return a / np.amax(a, axis=axis, keepdims=True)\n",
    "    "
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "from tensorflow.keras.utils import Sequence\n",
    "from tensorflow.keras.utils import to_categorical\n",
    "from tensorflow.image import grayscale_to_rgb\n",
    "import glob\n",
    "import os\n",
    "import numpy as np\n",
    "import imread\n",
    "\n",
    "def loadImage(filename):\n",
    "    return np.asarray(imread.imload_multi(filename))\n",
    "\n",
    "# Code (heavily) modified from https://stanford.edu/~shervine/blog/keras-how-to-generate-data-on-the-fly\n",
    "class DataGenerator(Sequence):\n",
    "    \n",
    "    def __init__(self, gt_path = '/DATA/GT', in_path = '/DATA/IN', setname='set*', snr=1000, snrname='dz400_IMax{}',\n",
    "                 patches_per_img=32, img_per_epoch=1000, batch_size=32, dim=(64,64), n_channels=1, n_classes=1,\n",
    "                 use_fraction=None, norm=zeroMeanUnitVariance, shuffle=True):\n",
    "        \n",
    "        \n",
    "        #self.snr = snr\n",
    "        self.setname = setname\n",
    "        \n",
    "        self.patches_per_img = patches_per_img\n",
    "        self.img_per_epoch  = img_per_epoch\n",
    "        self.batch_size = batch_size\n",
    "        self.dim = dim\n",
    "        self.n_channels = n_channels\n",
    "        self.n_classes = n_classes\n",
    "        self.shuffle = shuffle\n",
    "        self.norm = norm\n",
    "        \n",
    "        \n",
    "        \"\"\"\n",
    "        # Aufgabe: Ergänzt den Programmcode zum Erstellen der Dateilisten\n",
    "        \n",
    "        Zu benutzende Variablen:\n",
    "        gt_path - Elternverzeichnis von Experimentordnern (vgl. Standardwert)\n",
    "        in_path - Großelternverzeichnis von Experimentordnern (vgl. Standardwert)\n",
    "        setname - Name des Experimentordners (vgl. Standardwert, soll * benutzen können)\n",
    "        snr - Signal zu Rauschverhältnis (vgl. STandardwert, Teil des Parameterordnernames von Input)\n",
    "        snrname - Formatstring in dem snr eingesetzt werden soll (vgl. \"Using % and .format() for great good!\")\n",
    "        \n",
    "        \n",
    "        Hinweis: Die Standardwerte sind bereits als Parameter von __init__ gegeben.\n",
    "        \n",
    "        \n",
    "        Zu erstellende Variablen:\n",
    "        \n",
    "        x_filelist - Liste von Dateipfaden (Strings) zu existierenden Inputbildern\n",
    "        y_filelist - Liste von Dateipfaden (Strings) zu existierenden Labelbildern\n",
    "        \n",
    "        \"\"\"\n",
    "        \n",
    "        \n",
    "        # Store filelist in member field\n",
    "        \n",
    "        self.x_filelist = x_filelist\n",
    "        self.y_filelist = y_filelist\n",
    "        \n",
    "        print(len(y_filelist))\n",
    "        \n",
    "        self.list_IDs = np.arange(self.y_filelist.shape[0]*self.patches_per_img)\n",
    "        self.indexes = np.arange(len(self.list_IDs))\n",
    "        \n",
    "    def __len__(self):\n",
    "        'Number of batches per epoch'\n",
    "        return int(np.floor(self.y_filelist.shape[0]*self.patches_per_img / self.batch_size)) \n",
    "\n",
    "    def __getitem__(self, index):\n",
    "        'Generate one batch of data'\n",
    "        # Select image from image list (Tweak: multiple images)\n",
    "        indexes = self.indexes[index*self.batch_size:(index+1)*self.batch_size] # shuffeled id list\n",
    "\n",
    "        list_IDs_temp = [self.list_IDs[k] for k in indexes] # ids in planned number of images per epoch\n",
    "\n",
    "        X = np.empty((self.batch_size, *self.dim, self.n_channels))\n",
    "        Y = np.empty((self.batch_size, *self.dim, self.n_classes), dtype=bool)\n",
    "        \n",
    "        file_ids = np.arange(self.y_filelist.shape[0])\n",
    "        file_ids = np.repeat(file_ids, self.patches_per_img)\n",
    "        \n",
    "        previous_file_id = -1\n",
    "        \n",
    "        for i, file_id in enumerate(file_ids[list_IDs_temp]):\n",
    "            \n",
    "            if previous_file_id != file_id:\n",
    "                im_x, im_y, shape =  self.__getimages__(file_id)\n",
    "                    \n",
    "            # else: work with existing im_x, im_y\n",
    "            \n",
    "            \"\"\"\n",
    "            # Aufgabe: Ergänzt eine Funktion, die aus im_x und im_y die Minibatches für Input und Label in \n",
    "            der unten beschriebenen Form an zufällig Positionen ausschneidet und für das Batch vorbereitet.\n",
    "            \n",
    "            \n",
    "            Zu benutzende Variablen\n",
    "            self.norm - Normierungsfunktion von Input-Minibatch (wie oben definiert)\n",
    "            self.dim - y, x Ausdehnung der Minibatches (wie oben definiert)\n",
    "            im_x - Eingesenes Input Volumen mit den Dimensionen (z, y, x)\n",
    "            im_y - Eingelesenes Label Volumen mit den Dimensionen (z, y, x)\n",
    "            shape - Form der beiden eingelesenen Volumen\n",
    "            \n",
    "            \n",
    "            Zu erstellende Variablen:           \n",
    "            im_x_norm - auf die Form (self.dim[0], self.dim[1], 1) zurechtgeschnittenes und normiertes Minibatch\n",
    "                        von im_x.\n",
    "            im_y_norm - auf die Form (self.dim[0], self.dim[1], 1) zurechtgeschnittenes Minibatch von im_y.\n",
    "            \"\"\"\n",
    "            \n",
    "            # stacks input for RGB networks\n",
    "            im_x_norms = ()\n",
    "            for _ in range(self.n_channels):\n",
    "                im_x_norms = im_x_norms + (im_x_norm,)\n",
    "            im_x_norm = np.concatenate(im_x_norms, axis=-1)\n",
    "                \n",
    "            X[i, ...] = im_x_norm\n",
    "            Y[i, ...] = im_y_norm\n",
    "            \n",
    "            previous_file_id = file_id\n",
    "            \n",
    "\n",
    "        return X, Y\n",
    "    \n",
    "    def __getimages__(self, file_id):\n",
    "        \"\"\" Reads .tif Files and cuts them to a cube where Cells can be found.\"\"\"\n",
    "        \n",
    "        \"\"\"\n",
    "        # Aufgabe: Ergänzt das Einlesen von Dateien\n",
    "        \n",
    "        Zu verwendende Funktionen (Bereits in Notebook definiert)\n",
    "        loadImage() \n",
    "        cropToCells()\n",
    "        \n",
    "        \n",
    "        \n",
    "        Benötigte Variablen:\n",
    "        \n",
    "        im_x - eingelesene Input .tif Datei als numpy Array\n",
    "        im_y - eingelesene Label .tif Datei als numpy Array\n",
    "        im_x_ - auf Würfel mit Zellen zurechtgeschnittenes Input Array  \n",
    "        im_y_ - auf Würfel mit Zellen zurechtgeschnittenes Label Array\n",
    "        shape - Zellwürfel Form\n",
    "        \n",
    "        \"\"\"\n",
    "\n",
    "        if not np.any(np.asarray(shape) - np.array([0, self.dim[0], self.dim[1]]) <= 0): # cutted array large enough\n",
    "            im_y = im_y_\n",
    "            im_x = im_x_\n",
    "        else:\n",
    "            shape = im_y.shape\n",
    "\n",
    "        return im_x, im_y, shape\n",
    "\n",
    "    def on_epoch_end(self):\n",
    "        'Updates indexes after each epoch'\n",
    "        if self.shuffle == True:\n",
    "            # random offset\n",
    "            offset = np.random.randint(0, len(self.list_IDs) )\n",
    "            indexes = np.empty_like(self.list_IDs)\n",
    "            indexes[:-offset] = self.list_IDs[offset:]\n",
    "            indexes[-offset:] = self.list_IDs[:offset]\n",
    "            \n",
    "            # random order of chunks with length of batch size\n",
    "            chunk_order = np.arange(self.__len__())\n",
    "            np.random.shuffle(chunk_order)\n",
    "            self.indexes = [i for index in chunk_order for i in indexes[index*self.batch_size:(index+1)*self.batch_size]]\n",
    "        else:\n",
    "            pass\n",
    "    "
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "#### Lösung"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "#### Testcode 1\n",
    "# Erstelle eine neue DataGenerator Instanz\n",
    "beispiel_generator = DataGenerator()"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "#### Testcode 2\n",
    "# Erstellt ein Beispielbild\n",
    "images = beispiel_generator.__getimages__(1)[:2]\n",
    "# Zeigt die Beispielbilder an\n",
    "zSlicer(*images)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "#### Testcode 3\n",
    "# Erstellt ein Beispieltrainingsset\n",
    "batches = beispiel_generator.__getitem__(0)\n",
    "# Zeigt das Beispieltrainingsset an\n",
    "plotBatches(*batches)"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "## Netzwerktraining mit Generator Tensorboard-Callback\n",
    "\n",
    "Wir haben jetzt die Möglichkeit alle Mikroskopaufnahmen zum Training zu benutzen. \n",
    "\n",
    "#### Aufgabe:\n",
    "- Erstelle ein neues Netzwerk mit dem oben benutzten Design.\n",
    "- Kompiliert das Netzwerk mit *stochastic gradient descent* als Optimierungsfunktion,  `binary_crossentropy` als Verlust Funktion und der `accuracy`-Metrik.\n",
    "- [Trainiert](https://keras.io/models/model/#fit_generator) das Netzwerk mit\n",
    "    - Einer neuen `DataGenerator`-Instanz \n",
    "        - mit `setname='set9*'`\n",
    "        - und `use_fraction=0.5`\n",
    "    - Auf 20 Epochen\n",
    "    - Nutzt dabei Multiprocessing und 6 worker um das Erstellen der Minibatches zu beschleunigen.\n",
    "    - Erstellt eine Variable `history` die den Trainingsprozess als Wert zugewiesen bekommt.\n",
    "\n",
    "Hinweis: Der Traininsprozesses dauert auf der verbauten Grafikkarte ca. 10 min."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "#### Lösung"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "In der Variable `history` haben wir Informationen über den Trainingsprozess gespeichert.\n",
    "Mit diesem Informationen wollen wir ein Diagramm erstellen, das den Verlauf des Trainingsprozesses über mehrere Epochen darstellt.\n",
    "\n",
    "#### Beispiel"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "# Data preparation\n",
    "Y = [history.history['acc'], history.history['loss']]\n",
    "YLabels = ['Genauigkeit', 'Verlust']\n",
    "\n",
    "# Plotting\n",
    "f, axes = plt.subplots(1,2, figsize=(9.5, 5))\n",
    "\n",
    "for ax, y, ylabel in zip(axes, Y, YLabels):\n",
    "    ax.plot(np.arange(1, len(y)+1), y, 'ro', label='6x Conv2D')\n",
    "    ax.set_xlabel('Epoche')\n",
    "    ax.set_ylabel(ylabel)\n",
    "    ax.legend()\n",
    "    ax.grid()\n",
    "\n",
    "plt.tight_layout()"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "### **Bitte sammelt alle Trainigsverläufe in entsprechenden Dateien!**\n",
    "\n",
    "Um euch die Arbeit zu erleichtern, haben eine Funktion vorgeschrieben, welches nicht nur den Trainingsverlauf, sondern auch das fertig trainierte Modell abspeichert. Sicher ist sicher ..."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "from datetime import datetime\n",
    "import csv\n",
    "import os\n",
    "\n",
    "def save_model_history(history, net, foldername):\n",
    "    \n",
    "    folder = os.path.join(foldername, datetime.now().strftime('%Y-%m-%d_%H-%M'))\n",
    "    \n",
    "    if not os.path.exists(folder):\n",
    "        os.makedirs(folder)\n",
    "        \n",
    "    hist_pivot = [dict(zip(history.history, col)) for col in zip(*history.history.values())] # https://stackoverflow.com/a/37489316\n",
    "        \n",
    "    with open(os.path.join(folder, 'history.csv'), 'w') as f:\n",
    "        w = csv.DictWriter(f, history.history.keys())\n",
    "        w.writeheader()\n",
    "        for row in hist_pivot:\n",
    "            w.writerow(row)\n",
    "        \n",
    "    \n",
    "    \n",
    "    model_json = net.to_json()\n",
    "    with open(os.path.join(folder, \"model.json\"), \"w\") as json_file:\n",
    "        json_file.write(model_json)\n",
    "        \n",
    "    net.save_weights(os.path.join(folder, \"weights.h5\"))"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "#### Beispiel"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "save_model_history(history, net, 'simple2conv')"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Tatsächlich gibt es eine noch elegantere Variante den Trainingsprozess im Auge zu behalten: *Tensorboard*.\n",
    "Mit *Tensorboard* kann man alle Trainingsprozesse (auch mehrere gleichzeitig) beobachten.\n",
    "\n",
    "Der Keras Callback `Tensorboard` ermöglicht erstellt dabei die notwendigen log-Dateien. Der Callback lässt sich sehr einfach in das Training integrieren. Dazu müsst ihr denn Callback zunächst mit einem (am besten eindeutigen) Verzeichnisnamen unterhalb von `./logs` erstellen:\n",
    "\n",
    "``` python\n",
    "from datetime import datetime\n",
    "from tensorflow.keras.callbacks import TensorBoard\n",
    "\n",
    "datestr = datetime.now().strftime('%Y-%m-%d_%H-%M')\n",
    "tboard = TensorBoard(log_dir='./logs/{}_6x_Conv2D'.format(datestr), update_freq='batch')\n",
    "```\n",
    "\n",
    "In der Funktion `fit_generator` müsst ihr nur noch das Schlüsselwort `callbacks=[tboard]` ergänzen. Der Callback zeichnet automatisch den Trainingsprozess auf. Der im Hintergrund bereits von uns gestartete Tensorboardprozess sammelt die Informationen und zeigt Sie auf einer eigenen [Website](http://localhost:6006) an.\n",
    "\n",
    "#### Aufgabe\n",
    "- Startet das letzte Training mit einem Tensorboard Callback erneut und beobachtet den Graphen auf der Tensorboard Website"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "#### Lösung"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "# Schritt für Schritt zum Autoencoder\n",
    "\n",
    "Die Genauigkeit aller bisher trainierte Netzwerke ist nicht zufriedenstellend.\n",
    "\n",
    "Unser Modell soll bis zu 98% Genauigkeit erreichen. Die bisher verwendete Architektur ist dafür nicht geeignet.\n",
    "\n",
    "Damit dies möglich wird, muss das Netzwerk *tiefer* werden d.h. wir benötigen mehr Ebenen zwischen Input und Output. Gleichzeitig wollen wir Informationen von benachbarten Pixeln zusammenfassen. Dazu werden Pooling Ebenen wie `MaxPooling2D` benutzt. Im Gegensatz zu Konvolutionen modifiziert Pooling nicht die Filteranzahl. Das Gegenstück zu `MaxPooling2D` stellen Upsampling Ebenen wie `Upsampling2D` dar.\n",
    "\n",
    "#### Aufgabe:\n",
    "- Erstellt basierend auf der bereits erstellten Funktion eine neue Funktion die euch fertige Modelle zurückgibt:\n",
    "    1. `Input` (Parameter wie oben)\n",
    "    2. 2x `Conv2D` (Parameter wie oben) \\[1\\]\n",
    "    3. `MaxPooling2D` (mit `pool_size=(2,2)`)\n",
    "    4. 2x `Conv2D` (Parameter wir oben nur doppelter Filter Anzahl)\n",
    "    5. `Upsampling2D`(`size` wie in `MaxPooling2D`)\n",
    "    6. `Conv2D` (Parameter wie oben nur Kernelgröße von 2) \\[2\\]\n",
    "    7. `Concatenate`(entlang der Filter Dimension mit einer Liste mit \\[1\\] und \\[2\\] als Input)\n",
    "    8. 2x `Conv2D` (Parameter wie oben)\n",
    "    9. 2x `Conv2D` (wie die letzten beiden Konvolutions oben)\n",
    "- Testet das Model auf 10 Epochen mit\n",
    "    - Trainingsdatengenerator mit `setname='set1*'`\n",
    "    - Validierungsgenerator mit `setname='set9*'` und `use_fraction=0.20`"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "#### Lösung"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Das neue Design bleibt hinter den Erwartungen zurück. Durch noch tiefere Netzwerke können wir die Genauigkeit auf den angestrebten Wert steigern.\n",
    "#### Aufgabe\n",
    "- Schreibt eine Funktion für den MaxPooling & Convolution Schritt. (3. und 4.)\n",
    "- Schreibt eine Funktion für den UpSampling und Convolution Schritt (5., 6., 7. und 8.)\n",
    "- Erstellt Netzwerke mit 2, 3, 4 MaxPooling Schritte (und der selben Anzahl an UpSampling Layer)\n",
    "- Fügt zusätzlich in die Netzwerke mit 3, 4 Pooling Ebenen eine Dropout Ebene nach dem 3. und 4. Pooling Schritt ein.\n",
    "- Die 'dropout_rate' von den Dropout Ebenen soll 0.5 betragen."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "#### Lösung"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "In manchen Fällen konviert das Stochastische Gradientenabstiegsverfahren in einem lokalen Minimum, dass sehr weit weg von einem nutzbaren Modell liegt. In vielen Fällen reicht es aus eine andere Optimierungsmethode zu verwenden, wie z.B. [Adam](https://arxiv.org/pdf/1412.6980v8.pdf), welche neben dem Gradienten auch die Momente des Gradienten erster und zweiter Ordnung für die Optimierung benutzt.\n",
    "\n",
    "Der [Dropout](https://arxiv.org/pdf/1207.0580.pdf) Ansatz reduziert die für das Training benutzte Neuronen adaptiv um das Modell robuster zu machen.\n",
    "\n",
    "#### Aufgabe:\n",
    "- Trainiert alle oben angegeben Netzwerke (inkl. des Netzwerks welches nur 6 Konvolutionen benuzt) mit der `Adam` Optimierungsfunktion anstatt `SGD`.\n",
    "- Benutzt dabei eine *Learning Rate* von 1e-4.\n",
    "- 20 Epochen\n",
    "- Vergleicht die Performance der Netzwerke.\n",
    "\n",
    "Hinweis:\n",
    "- Schreibt eine Schachtelung von For-loops die die Aufgaben unmittelbar nacheinander ausführt.\n",
    "- Erstellt euch mit Hilfe der Funktion `create_batches_ram` einmalig ein Trainings - und Validationset und verzichtet bei diesem Benchmark auf die direkte Nutzung des Generators während des Trainingsprozesses. "
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "def create_batches_ram(generator):\n",
    "    tmp_batch = generator.__getitem__(0)[0]\n",
    "    batch_size = tmp_batch.shape[0]\n",
    "    num_batches = generator.__len__()\n",
    "    \n",
    "    X = np.zeros((num_batches*batch_size, 64, 64, 1))\n",
    "    Y = np.zeros((num_batches*batch_size, 64, 64, 1))\n",
    "    \n",
    "    \n",
    "    for i in range(num_batches):\n",
    "        out = generator.__getitem__(i)\n",
    "        X[i*batch_size:(i+1)*batch_size] = out[0]\n",
    "        Y[i*batch_size:(i+1)*batch_size] = out[1]\n",
    "    \n",
    "    return X, Y"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "#### Lösung"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "# Treffe Voraussagen auf noch unbekannte Daten"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Im nächsten Teil des Praktikums wollen wir uns mehr auf die Eigenschaften der Zellen innerhalb der Volumen konzentrieren. \n",
    "\n",
    "- Welche physikalischen Zelleigenschaften können wir aus noch unbekannten Mikroskopaufnahmen extrahieren.\n",
    "\n",
    "Dazu benutzen wir die Mikroskopaufnahmen unter `/DATA/validation/noise1000/`. Im Folgenden wollen wir nur das Netzwerk mit der höchsten Genauigkeit auf den Validierungsdaten verwenden."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "!ls /DATA/validation/noise1000/"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Sowohl der Trainingsgenerator, als auch der Validationgenerator erstellt an zufälligen Positionen im Volumen 2D Minibatches mit der Größe von $64 \\times 64$ Pixel.\n",
    "\n",
    "Diese Volumenausschnitte haben auf den ersten Blick wenig mit einer Vorhersage eines kompletten Volumens zu tun. Dazu muss das Volumen in kleine Teilbilder zerschnitten werden. Dabei ist es wichtig, dass eine festgelegte Reihenfolge eingehalten wird, damit hinterher die Vorhersagen wieder zu einem Volumen zusammengesetzt werden können.\n",
    "\n",
    "### Aufgabe:\n",
    "- Schaut euch die unten definierte Funktion `predictImageStack`an.\n",
    "- Benutzt die Funktion um Zell vorhersagen der oben angegeben Mikroskopbilder anzufertigen.\n",
    "- Berechnet von jeder Zellvorhersage den Durchschnitt in $x$, $y$ und $z$-Richtung. Fällt euch etwas auf?"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "def predictImageStack(img, model, num_channels=1, num_classes=1, batch_size=32, norm=zeroMeanUnitVariance):\n",
    "    print('Image shape: {}'.format(img.shape))\n",
    "    input_shape = model.layers[0].output_shape\n",
    "    print('Model input shape: {}'.format(input_shape))\n",
    "    \n",
    "    in_z, in_y, in_x = img.shape\n",
    "    height, width = input_shape[1:3]\n",
    "    nz = in_z - 1\n",
    "    ny = int(np.floor(in_y / height))\n",
    "    nx = int(np.floor(in_x / width))\n",
    "    \n",
    "    print('Number of tiles in x, y, z: ({:d}, {:d}, {:d})'.format(nx, ny, nz))\n",
    "    \n",
    "    margin_x = (in_x - (width * nx)) / 2\n",
    "    margin_y = (in_y - (height * ny)) / 2\n",
    "    \n",
    "    \n",
    "    img = img[..., None]\n",
    "    imgs = ()\n",
    "    for _ in range(num_channels):\n",
    "        imgs = imgs + (img,)\n",
    "    img = np.concatenate(imgs, axis=-1)\n",
    "    \n",
    "    # cut away overview plane in z and reduce size in x, y to fit model input\n",
    "    img_ = img[1:, int(np.floor(margin_y)):-int(np.ceil(margin_y)), int(np.floor(margin_x)):-int(np.ceil(margin_x)), :]\n",
    "    \n",
    "    print('cut image to shape: {}'.format(img_.shape))\n",
    "    \n",
    "    num_batches = int(nx * ny * nz)\n",
    "    \n",
    "    print('Total number of tiles: {}'.format(num_batches))\n",
    "    \n",
    "    prediction_batches = np.zeros((num_batches, height, width, num_channels))\n",
    "    \n",
    "    print('Prediction batches shape: {}'.format(prediction_batches.shape))\n",
    "    \n",
    "    \n",
    "    print('Fill input batches with image data')\n",
    "    k = 0\n",
    "    for i in range(ny):\n",
    "        for j in range(nx):\n",
    "            for l in range(nz):\n",
    "                prediction_batches[k, :, :, :] = norm(img_[l, height*i:height*(i+1), width*j:width*(j+1)], axis=(0, 1, 2))\n",
    "                k += 1\n",
    "                \n",
    "    num_predictions = num_batches // batch_size\n",
    "    \n",
    "    print('Predict batches')\n",
    "    for i in range(num_predictions):\n",
    "        prediction_batches[i*batch_size:(i+1)*batch_size] = model.predict(prediction_batches[i*batch_size:(i+1)*batch_size])\n",
    "    \n",
    "    if num_predictions*batch_size != num_batches:\n",
    "        prediction_batches[num_predictions*batch_size:] = model.predict(prediction_batches[num_predictions*batch_size:])\n",
    "    \n",
    "    img_prediction = np.zeros((*img_.shape[:3], num_classes))\n",
    "    \n",
    "    print('Constuct predicted image with shape: {}'.format(img_prediction.shape))\n",
    "    \n",
    "    k = 0\n",
    "    for i in range(ny):\n",
    "        for j in range(nx):\n",
    "            for l in range(nz):\n",
    "                img_prediction[l, height*i:height*(i+1), width*j:width*(j+1)] = prediction_batches[k, ..., :num_classes]\n",
    "                k += 1\n",
    "\n",
    "    \n",
    "    return np.squeeze(img_prediction), np.squeeze(img_)\n"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "#### Lösung"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "## Isotrope Auflösung\n",
    "\n",
    "Um wirklich quantitative Aussagen über Zellen zu treffen (Volumen, Intensität, Abstände), müssen wir beachten, dass\n",
    "in unterschiedliche Raumrichtungen unterschiedliche Auflösungen existieren. In den Beispielbildern wird von einem Pixelabstand von $63,2\\,\\mathrm{nm}$ in $x$ und $y$-Richtung und von $400\\,\\mathrm{nm}$ in $z$-Richtung verwendet.\n",
    "\n",
    "Für isotrope Auflösungen wollen müssen wir den Pixelabstand in $z$ Richtung auf ebenfalls $63,2\\,\\mathrm{nm}$ linear interpolieren.\n",
    "\n",
    "Dazu bietet sich die Funktion `affine_transform` im Modul `scipy.ndimage` an. Die benötigten Parameter sind:\n",
    " - das Usprungsvolumen (ohne Übersichtsebene in `[0, :, :]`)\n",
    " - die Transformationsmatrix (eine $3 \\times 3$ Diagonalmatrix, die auf der Hauptdiagonale die Einträg $\\left(d, 1, 1\\right)$ besitzt, wobei $d$ der gewünschte Streckungsfaktor in $z$-Richtung darstellt).\n",
    " - die neue Form des Volumens mit bereits angepasster Anzahl an $z$-Ebenen\n",
    " - `order=1` welches die Interpolation auf die lineare Ordnung beschränkt.\n",
    "\n",
    "#### Aufgabe\n",
    "- Berechnet die lineare Interpolation von allen .tif Dateien in `/DATA/validation/`\n",
    "- Speichert die Interpolationen als `.npy`-Datei auf der Festplatte. \n"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "#### Lösung"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "#### Aufgabe: \n",
    "- Schaut euch die Dateigröße der eben gespeicherten Dateien an.\n",
    "    - Dazu können ihr beispielsweise Programm `du` benutzen.\n",
    "    - Input ist der Verzeichnispfad.\n",
    "    - Die Optionen `-hs` erstellt eine gut lesebare Zusammenfassung. "
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "#### Lösung"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Die Dateien sind bis zu mehrere Gigabyte groß.\n",
    "\n",
    "Bei jeder folgenden Rechenoperation muss die gesamte Datei in den Arbeitsspeicher geladen werden. Deshalb wollen wir uns schnell einen Überblick verschaffen, ob wir tatsächlich alle Voxel (=Pixel in 3D) benötigen."
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "#### Aufgabe\n",
    "- Erstellt eine Maximumsprojektion aller 3 Raumachsen\n",
    "- Lasst euch die Projektionen anzeigen und benutzt nur den Teil des Volumens, der Zellen beinhaltet.\n",
    "- Speichert die so reduzierten Volumen als `.npy`-Dateien ebenfalls auf der Festplatte und vergleicht die Dateigrößen."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "#### Lösung"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Die Dateigröße konnte sehr stark reduziert werden. Dies ist gängig Praxis bei Mikroskopdaten, bei denen nur ein bestimmter Bildausschnitt (Region of Interest) für die Analyse interessant ist. Das Verkleinern des Bildausschnitts auf die ROI nennt man *Cropping*. Dies hat vor allem ressourcetechnische Vorteile:\n",
    "1. Die Dateien brauchen weniger Speicherplatz.\n",
    "2. Die Analyse auf den Daten geht schneller.\n",
    "\n",
    "Zu Bedenken ist aber auch, dass dies die Möglichkeit eröffnet nur Bildausschnitte zu veröffentlichen, die die eigenen Thesen untermauern. Deshalb gilt es als gute wissenschaftliche Praxis die *Primärdaten* ebenfalls aufzuheben.\n",
    "\n",
    "**Bitte speichert alle Ergebnisse, für die Ihr die Grafikkarten zur Vorhersage benutzt habt auf der Festplatte. Achtet dabei bitte auf eine eindeutige Benennung bzw. schreibt euch alle Schritte in eurem Laborbuch auf. Denkt daran, dass ihr auf den Daten auch im Nachhinein noch Analysen für das Protokoll anfertigen müsst.**"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "#### Aufgabe\n",
    "- Nutzt ein Netzwerk um die Zellen in den Volumen zu erkennen.\n",
    "- Speichert als `.npy`-Datei:\n",
    "    - die Zellvorhersagen\n",
    "    - die zuerechtgeschnitten Intensitätsbilder\n",
    "- Benutzt die Funktion `zSlicer` um einen Input mit der Vorhersage zu vergleichen\n",
    "- Erstellt ein Histogramm der Pixelwerte in der Vorhersage\n",
    "\n",
    "Hinweis: Benutzt für das Histogram logarithmische Auftragung auf der y-Achse und mindestens 100 Abschnitte."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "#### Lösung"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "# Auswertung\n",
    "\n",
    "**Die Auswertung ist mit den gespeicherten Werten auch Zuhause oder im Lernzentrum möglich.**\n",
    "\n",
    "\n",
    "Mit der Vorhersage wollen wir jetzt weiterarbeiten. Momentan stehen wir noch vor einigen Problemem:\n",
    "1. Die Vorhersage muss diskretisiert werden (von \\[0, 1\\] zu \\{0, 1\\})\n",
    "2. Wir wollen Objekte, die nicht miteinander verbunden sind mit unterschiedlichen [Label](http://scikit-image.org/docs/dev/auto_examples/segmentation/plot_label.html#sphx-glr-auto-examples-segmentation-plot-label-py) versehen \n",
    "\n",
    "#### Aufgabe:\n",
    "- Berechnet ein Volumen in dem jedes Objekt ein eigenes Label bekommt.\n",
    "\n",
    "Hinweise: \n",
    "- In dem gegebenen Beispiel heißt das mit Label versehene Volumen `label_image`.\n",
    "- Nicht alle Operationen werden in unserem Fall benötigt\n",
    "- Benutzt nach Bedarf `binary_erosion`, `ball`, `binary_dilation` um kleine Objekte aus eurem Datensatz zu filtern"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "#### Lösung"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Diesen segmentierten Daten kann man schon eine Reihe von Parameter extrahieren. Im Detail sind dies:\n",
    "- Zellpositionen (Beispiel)\n",
    "- Mittlere/ Maximale/ Minimale  Helligkeit \n",
    "- Zellvolumen\n",
    "- Lokale Dichte\n",
    "- Zellabstände (optional)\n",
    "- Nemantic Orderparamter (optional)\n",
    "\n",
    "## Zellposition\n",
    "\n",
    "#### Beispiel\n",
    "\n",
    "Sei `w` die oben bestimmte Labelmatrix (bzw. das Labelvolumen) und `img_` die der Segmentierung zugehörige Intensitätswerte, dann lässt sich ein numpy Array `centroids_array` mit allen Schwerpunkten aller Objekten in der Labelmatrix bestimmen via:"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "from skimage.measure import regionprops\n",
    "\n",
    "props = regionprops(w, img_)\n",
    "centroids_list = [props[i].centroid for i in range(len(props))]\n",
    "centroids_array = np.asarray(centroids_list)"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Dabei beinhaltet das Ergebnis von `skimage.measure.regionprops` noch eine ganze Reihe von weiteren [Objekteigenschaften](http://scikit-image.org/docs/dev/api/skimage.measure.html#regionprops)."
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "## Helligkeitswerte\n",
    "\n",
    "#### Aufgabe\n",
    "- Bestimmt die mittlere, maximale und minimale Helligkeit aller Objekte in dem gelabelten Volumen.\n",
    "- Erstellt je ein Histogram für jeden der drei Helligkeitsverteilungen."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "#### Lösung"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "## Zellvolumen\n",
    "\n",
    "Auch das Volumen lässt sich auf diese Art und Weise bestimmen, wobei der entsprechende Parametername von `props` nicht direkt ersichtlich ist.\n",
    "\n",
    "#### Aufgabe\n",
    "- Finde den entsprechenden Parametername heraus und stelle die Verteilungen der einzelnen Volumen in metrischen-Einheiten in einem Diagramm dar"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "#### Lösung"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "## Lokale Dichte\n",
    "\n",
    "Obwohl Position, Helligkeit und Volumen interessante Zellparameter sind, können wir aus den Segmentierung noch wesentlich mehr Informationen extrahieren. In den meisten Veröffentlichungen konzentiert man sich auf einige wenige Eigenschaften. Weil keine andere Arbeitsgruppe diese Eigenschaften bereits untersucht hat, muss man sich selbst überlegen wie man die Messwerte aus den Bildern extrahieren kann.\n",
    "\n",
    "Ein solches Beispiel wäre die lokale Dichte, die man nicht in den `regionprops` Eigenschaften findet. Stattdessen sind eigene Programmierkenntnisse notwendig."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "from numpy.linalg import norm\n",
    "\n",
    "def calculateLocalDensity(w, centroid_array, radius = 60):\n",
    "    \"\"\" Calculates occupied volume in a sphere around centroid\n",
    "    \"\"\"\n",
    "\n",
    "    Z, Y, X = np.meshgrid(np.arange(-radius, radius, 1), np.arange(-radius, radius, 1), np.arange(-radius, radius, 1))\n",
    "    Z, Y, X = np.ravel(Z), np.ravel(Y), np.ravel(X)\n",
    "    sphere = norm([Z, Y, X], axis=0) < radius\n",
    "\n",
    "    w_ = w > 0\n",
    "    w_ = np.pad(w_, radius, 'constant', constant_values=-1)\n",
    "\n",
    "    vol_sphere = np.sum(sphere) \n",
    "\n",
    "    localDensity = np.zeros(len(props))\n",
    "    for i in range(len(props)):\n",
    "        z, y, x = centroid_array[i]\n",
    "        z, y, x = np.round([z, y, x])\n",
    "        w_part = w_[int(z):int(z+2*radius), int(y):int(y+2*radius), int(x):int(x+2*radius)]\n",
    "        vol = w_part.flatten()*sphere\n",
    "        out = np.sum(w_part == -1)\n",
    "        localDensity[i] = np.sum(vol==1)/(vol_sphere - out);\n",
    "\n",
    "    return localDensity"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "# Visualisierung\n",
    "\n",
    "\n",
    "Wirklich interessant werden die einzelnen Parameter erst, wenn man die räumliche Verteilung kennt. Dazu eigenen u.A. 3D Scatterplots.\n",
    "\n",
    "In der *matplotlib* Dokumentation findet ihr ein sehr einfaches [Beispiel](https://matplotlib.org/examples/mplot3d/scatter3d_demo.html) dazu."
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "#### Aufgabe\n",
    "- Erstellt ein 3D Scatterplot für die durchschnittliche Helligkeit und die lokale Dichte.\n",
    "    - Die Markerposition soll durch den jeweiligen Schwerpunkt gegeben sein.\n",
    "    - Die Markerfarbe soll den Parameter darstellen.\n",
    "    \n",
    "Hinweis: Berechnet die lokale Dichte mit der gegebenen Funktion `calculateLocalDensity`. Passt dabei den Radius auf ein sinnvolles Maß an. (Dabei kann ein Blick in die Volumenverteilung nicht schaden)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "#### Lösung"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "# 3D Rendering\n",
    "\n",
    "Zum Abschluss stehen euch zwei Möglichkeiten offen 3D Visualisierungen der Zellen zu erstellen.\n",
    "\n",
    "## *matplotlib*\n",
    "\n",
    "Kann in diesem Notebook gemacht werden, dauert aber lange und die Plots lassen sich nicht wirklich interativ drehen.\n",
    "\n",
    "**Bitte entfernt nur die Anführungszeichen, wenn ihr 5min warten könnt!**"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "from skimage.measure import marching_cubes_lewiner\n",
    "\n",
    "verts, faces, normals, values = marching_cubes_lewiner(dilate, 0.5)\n",
    "\n",
    "\"\"\"\n",
    "fig = plt.figure()\n",
    "ax = fig.add_subplot(111, projection='3d')\n",
    "ax.plot_trisurf(verts[:, 0], verts[:,1], faces, verts[:, 2],\n",
    "                cmap='Spectral', lw=1)\n",
    "\n",
    "plt.show()\n",
    "\"\"\""
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "## *ParaView*\n",
    "\n",
    "ParaView  ist ein 3D Rendering Programm für große Datensätze und de facto Standard für die Visualisierung von aufwendigen Simulation. Die folgende Handreichung soll euch dazu ermuntern ParaView zum Rendern einer 3D Darstellung zu benutzen. (beispielsweise kann man mit ParaView Abbildungen wie z.B. in der [Versuchsbeschreibung](https://www.uni-marburg.de/de/fb13/studium/praktika/praktika-fuer-physiker/v-so-15-maschinelles-lernen) erstellen oder gar ganze Simulationsreihen in einen [Film](https://static-content.springer.com/esm/art%3A10.1038%2Fs41567-018-0356-9/MediaObjects/41567_2018_356_MOESM3_ESM.mp4) umwandeln)."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "# Export to VTI\n",
    "from pyevtk.hl import imageToVTK\n",
    "\n",
    "imageToVTK('{}_unet_prediction'.format(datestr),\n",
    "           pointData={\"prediction\":prediction})"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "## Download\n",
    "\n",
    "Ladet bitte den ParaView Installer von der offiziellen [Website](https://www.paraview.org/download/) herunter. Im folgenden wird die Windowsversion benutzt. Für Max OSX und der Linux Distribution eurer Wahl finden sich entsprechende Versionen im Mac App Store und der Paketverwaltung.\n",
    "\n",
    "## Datentransfer\n",
    "\n",
    "Um die Daten auf euren lokalen Computer herunterladen zu können müsst ihr ein Programm benutzen, dass Daten über eine [SSH](https://de.wikipedia.org/wiki/Secure_Shell)-Verbindung übertragen kann. Unter Mac OSX und Linux könnt ihr euren normalen Dateimanager benutzen. Unter Windows benötigt ihr ein zusätzliches Programm (z.B. [WinSCP](https://winscp.net/eng/download.php)).\n",
    "\n",
    "Für eine Verbindung mit WinSCP benötigt Ihr folgende Daten:\n",
    "* File protocol: SFTP\n",
    "* Host name: ***\n",
    "* User name: ***\n",
    "* Passwort: ***\n",
    "\n",
    "Unter nautilus (Linux Dateimanager) findet ihr das entsprechende Menü unter 'Other Locations' und 'Connect to Server'. In die Addresszeile müsst ihr eintragen:\n",
    "`sftp://<username>@<hostname>/`\n",
    "\n",
    "Bitte ladet die `.vti`-Dateien unter `/home/<username>/Notebooks` herunter\n",
    "\n",
    "## Datenvisualisierung in ParaView\n",
    "\n",
    "* Öffnet die Datei `JJJJ-MM-DD_hh-mm_unet_prediction.vti`.\n",
    "* Nutzt den Threshold Filter im Bereich eurer Labels.\n",
    "* Spielt ein bisschen mit den Einstellungen im Bereich *OSPRay Rendering* bis ihr eine schöne 3D Visualizierung erstellt habt\n",
    "\n",
    "Hinweise: \n",
    "* Für Schatten benötigt man eine Fläche die Schatten abbilden.\n",
    "* *path tracer* ist wesentlich schöner anzusehen, kostet aber auch sehr viel Rechenleistung."
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "# Abschluss\n",
    "\n",
    "1. Speichert ALLE Daten/ Notebooks/ DIAGRAMME\n",
    "2. Ladet ALLE .ipynb, .npy, .vti herunter\n",
    "3. Ladet alle .tif Bilder herunter\n",
    "\n",
    "\n",
    "# Sonstiges\n",
    "\n",
    "- Wenn Ihr weiter mit Neuronalen Netzwerken/ Python Notebooks experimentieren möchtet, aber keinen Computer mit GPU zur Verfügung habt:\n",
    "https://colab.research.google.com/notebooks/welcome.ipynb (Account wird benötigt)\n",
    "\n",
    "- Der Jupyter Notebook Server ist Teil der Anaconda Distribution:\n",
    "https://www.anaconda.com/distribution/ (ist im Lernzentrum installiert)\n",
    "\n",
    "- Wenn ihr eine neuere Nvidia Grafikkarte besitzt und diese für das Training von Neuronalen Netzwerken benutzen wollt: \n",
    "https://www.nvidia.com/de-de/gpu-cloud/ (Account wird benötigt)\n",
    "\n",
    "Dieser Notebook Server läuft auf einem Ubuntu System in einem eigenen Docker Container."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": []
  }
 ],
 "metadata": {
  "anaconda-cloud": {},
  "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
}