{"cells":[{"metadata":{"id":"_HAFr8fgGY4w"},"cell_type":"markdown","source":"![](https://aapmchallenges.blob.core.windows.net/public/logos/logo_2.jpg)\n# Table of contents\n\n* [Introduction](#Introduction)\n    - [Motivation](#Motivation)\n    - [Reconstruction versus Windows](#Reconstruction-versus-Windows)\n    - [Reusable code](#Reusable-code)\n\n\n* [DICOM-CT data 3D visualizations](#DICOM-CT-data-3D-visualizations)\n    - [Reading data](#Reading-data)\n    - [Windows and LUTs](#Windows-and-LUTs)\n    - [3D reconstruction with manual windowing and MarchingCubes algorithm](#3D-reconstruction-with-manual-windowing-and-MarchingCubes-algorithm)\n    - [Lungs segmentation](#Lungs-segmentation)\n    - [Other stuff: Vtk to numpy support](#Other-stuff:-Vtk-to-numpy-support)\n\n\n* [Gentle introduction to VTK](#Gentle-introduction-to-VTK)\n    - [VTK in standalone Kaggle notebooks](#VTK-in-standalone-Kaggle-notebooks)\n    - [Cylinder = Hello World example](#Cylinder-=-Hello-World-example)\n    - [Cylinder rotation with animation](#Cylinder-rotation-with-animation)\n    - [Kitchen example](#Kitchen-example)\n    - [Frog example](#Frog-example) \n \n \n* [Appendix: Slicer3D in Kaggle notebooks](#Appendix:-Slicer3D-in-Kaggle-notebooks)\n    - [Installation](#Installation)\n    - [MRHead-Example](#MRHead-Example)\n    - [CTChest-Example](#CTChest-Example)\n    - [OSIC training data Example](#OSIC-training-data-Example)"},{"metadata":{},"cell_type":"markdown","source":"# Introduction"},{"metadata":{"id":"7pcIltEAGY4y"},"cell_type":"markdown","source":"## Motivation\nThe main idea behind this notebook is 3D visualization of DICOM-CT data. I had an impression, that visualization of medical data of this type has great possibilities. On the other hand many popular libraries require x-server, which is not available on Kaggle/Colab. I had earlier, short but very good experience with Slicer3D software. In this program, visualization is based on VTK, which I would like to present. In this notebook you will find: a short introduction to VTK, DICOM-CT visualization and instructions to run Slicer from the command line.\n\n**Warning: This notebooks presents images containing blood visualizations. If you are afraid of the sight of blood, you read on your own responsibility.**"},{"metadata":{},"cell_type":"markdown","source":"## Reconstruction vs Windows\n* source: https://www.youtube.com/watch?v=KZld-5W99cI\n\nThere is a thing that we want to clarify here: the difference between a CT reconstruction algorithms versus the different windows that we display information on.\n\n### Reconstruction \n***Reconstruction*** require the image processing of the raw data commonly by the technologist, but sometimes by the radiologist using specialized software programs the data can be processed in different ways to spending depending on what data we wish to obtain from it or what tissues we wish to display. So typical reconstructions might include a soft tissue, a bone and a lung reconstruction. We can also reconstruct in multiple planes as well as the coronal and sagittal planes. We can do some specialized planes along the organ for any it along the plane of any organ of interest. We can do 3d reconstructions such as surface  reconstruction; maximum intensity projection reconstructions. \nReconstructions are done usually by the technologists or by some separate type of software on a different workstation by the radiologist to allow the raw CT data to be manipulated into different projections or different slices depending on what our clinical question is. Example of this procedure is shown in the diagram below:\n\n![](https://d3i71xaburhd42.cloudfront.net/d4fb82cbcfda18ebe327e7003ec81ac777666299/3-Figure1-1.png)\n\n### Windowing\nNow, the different windows are the different ways that we can visually display the digital information from any set of CT slices and our wide variety of windows common windows that are used at:\n* lung windows\n* soft tissue windows \n* bone windows \n* brain windows \n* blood windows \nand many others.\n\nExample windows:\n![](https://www.stepwards.com/wp-content/uploads/2019/12/Screen-Shot-2019-12-28-at-4.20.41-PM.png)\n\nIn my humble opinion, this topic is quite widely developed (also among the Kaggle community). You can find many interesting materials on this topic, my favorite notebook is the following: https://www.kaggle.com/redwankarimsony/rsna-str-pe-gradient-sigmoid-windowing\nIt is worth noting that this notebook is mainly focused on reconstruction."},{"metadata":{},"cell_type":"markdown","source":"## Reusable code\nProgramming in VTK is a separate section in this notebook and if you would like to use my work just copy the code below. I have presented more details about it in the second section (Gentle introduction to VTK):"},{"metadata":{"trusted":true,"_kg_hide-input":true,"id":"jyPjZRYUGY40"},"cell_type":"code","source":"! pip install pyvirtualdisplay -q\n! apt-get install -y xvfb >> /dev/null\n\nfrom IPython.display import Image\nimport imageio\nimport os\nimport shutil\nimport matplotlib.pyplot as plt \nimport matplotlib.image as mpimg \nimport gc\nfrom vtk.util import numpy_support\nimport numpy\nfrom pyvirtualdisplay import Display\n\ndisp = Display().start()\nimport vtk\ndisp.stop()\n\nN =  18\ndefault_width = 512\ndefault_height = 512\n\ndef vtk_show(renderer, width = default_width, height = default_height, filename = \"\"):\n\n    renderWindow = vtk.vtkRenderWindow()\n    \n    renderWindow.SetOffScreenRendering(1)\n    renderWindow.AddRenderer(renderer)\n    renderWindow.SetSize(width, height)\n    renderWindow.Render()\n     \n    windowToImageFilter = vtk.vtkWindowToImageFilter()\n    windowToImageFilter.SetInput(renderWindow)\n    windowToImageFilter.Update()\n     \n    writer = vtk. vtkPNGWriter()\n    \n    if filename == \"\":\n        writer.SetWriteToMemory(1)\n        writer.SetInputConnection(windowToImageFilter.GetOutputPort())\n        writer.Write()    \n        return bytes(memoryview(writer.GetResult()))\n    else:\n        writer.SetFileName(filename+\".png\")\n        writer.SetInputConnection(windowToImageFilter.GetOutputPort())\n        writer.Write()    \n        return None\n    \ndef vtk_render_gif(renderer, N, name, Roll = False, Azimuth = False, Elevation = False, Actor = None, RotateX = False, RotateY = False, RotateZ = False, Zoom = 0, Dolly = 0, standard = True, width = default_width, height = default_height):    \n    if standard:\n        renderer.ResetCamera()\n        camera = renderer.MakeCamera()\n        renderer.ResetCameraClippingRange()\n        camera.SetPosition(0,0,0)\n    os.makedirs(name,exist_ok=True)\n    \n    if Zoom != 0:\n        renderer.GetActiveCamera().Zoom(Zoom)\n        \n    if Dolly != 0:\n        renderer.GetActiveCamera().Dolly(Dolly)\n        \n    #tmpN = 1\n    if N >0: # render gif\n        for fi in range(N):\n            if Roll:\n                renderer.GetActiveCamera().Roll(360//N) \n            if Azimuth:\n                renderer.GetActiveCamera().Azimuth(360//N) \n            if Elevation:\n                renderer.GetActiveCamera().Elevation(360//N)\n            if Actor is not None:\n                if RotateX:\n                    Actor.RotateX(360//N)\n                if RotateY:\n                    Actor.RotateY(360//N)\n                if RotateZ:\n                    Actor.RotateZ(360//N)                    \n            vtk_show(renderer,filename = name + \"/shot\"+str(fi), width = width, height = height)\n        # render gif and cleanup\n        img_list = []\n        for fi in range(N):\n            img_list.append(mpimg.imread(name + '/shot' + str(fi) + '.png'))\n        shutil.rmtree(name)\n        imageio.mimsave(name + \".gif\", img_list, duration=0.5)\n\n    #if N == 1: # render png\n       #vtk_show(renderer,filename = name + \".gif\")\n\ndef CreateLut():\n    colors = vtk.vtkNamedColors()\n\n    colorLut = vtk.vtkLookupTable()\n    colorLut.SetNumberOfColors(17)\n    colorLut.SetTableRange(0, 16)\n    colorLut.Build()\n\n    colorLut.SetTableValue(0, 0, 0, 0, 0)\n    colorLut.SetTableValue(1, colors.GetColor4d(\"salmon\"))  # blood\n    colorLut.SetTableValue(2, colors.GetColor4d(\"beige\"))  # brain\n    colorLut.SetTableValue(3, colors.GetColor4d(\"orange\"))  # duodenum\n    colorLut.SetTableValue(4, colors.GetColor4d(\"misty_rose\"))  # eye_retina\n    colorLut.SetTableValue(5, colors.GetColor4d(\"white\"))  # eye_white\n    colorLut.SetTableValue(6, colors.GetColor4d(\"tomato\"))  # heart\n    colorLut.SetTableValue(7, colors.GetColor4d(\"raspberry\"))  # ileum\n    colorLut.SetTableValue(8, colors.GetColor4d(\"banana\"))  # kidney\n    colorLut.SetTableValue(9, colors.GetColor4d(\"peru\"))  # l_intestine\n    colorLut.SetTableValue(10, colors.GetColor4d(\"pink\"))  # liver\n    colorLut.SetTableValue(11, colors.GetColor4d(\"powder_blue\"))  # lung\n    colorLut.SetTableValue(12, colors.GetColor4d(\"carrot\"))  # nerve\n    colorLut.SetTableValue(13, colors.GetColor4d(\"wheat\"))  # skeleton\n    colorLut.SetTableValue(14, colors.GetColor4d(\"violet\"))  # spleen\n    colorLut.SetTableValue(15, colors.GetColor4d(\"plum\"))  # stomach\n\n    return colorLut\n\ndef CreateTissueMap():\n    tissueMap = dict()\n    tissueMap[\"blood\"] = 1\n    tissueMap[\"brain\"] = 2\n    tissueMap[\"duodenum\"] = 3\n    tissueMap[\"eyeRetina\"] = 4\n    tissueMap[\"eyeWhite\"] = 5\n    tissueMap[\"heart\"] = 6\n    tissueMap[\"ileum\"] = 7\n    tissueMap[\"kidney\"] = 8\n    tissueMap[\"intestine\"] = 9\n    tissueMap[\"liver\"] = 10\n    tissueMap[\"lung\"] = 11\n    tissueMap[\"nerve\"] = 12\n    tissueMap[\"skeleton\"] = 13\n    tissueMap[\"spleen\"] = 14\n    tissueMap[\"stomach\"] = 15\n\n    return tissueMap\n\ntissueMap = CreateTissueMap()\n\ncolorLut = CreateLut()\n\ndef CreateTissue(reader, ThrIn, ThrOut, color = \"skeleton\", isoValue = 127.5):\n    selectTissue = vtk.vtkImageThreshold()\n    selectTissue.ThresholdBetween(ThrIn,ThrOut)\n    selectTissue.ReplaceInOn()\n    selectTissue.SetInValue(255)\n    selectTissue.ReplaceOutOn()\n    selectTissue.SetOutValue(0)\n    selectTissue.Update()\n    selectTissue.SetInputConnection(reader.GetOutputPort())\n\n    gaussianRadius = 5\n    gaussianStandardDeviation = 2.0\n    gaussian = vtk.vtkImageGaussianSmooth()\n    gaussian.SetStandardDeviations(gaussianStandardDeviation, gaussianStandardDeviation, gaussianStandardDeviation)\n    gaussian.SetRadiusFactors(gaussianRadius, gaussianRadius, gaussianRadius)\n    gaussian.SetInputConnection(selectTissue.GetOutputPort())\n\n    #isoValue = 127.5\n    mcubes = vtk.vtkMarchingCubes()\n    mcubes.SetInputConnection(gaussian.GetOutputPort())\n    mcubes.ComputeScalarsOff()\n    mcubes.ComputeGradientsOff()\n    mcubes.ComputeNormalsOff()\n    mcubes.SetValue(0, isoValue)\n\n    smoothingIterations = 5\n    passBand = 0.001\n    featureAngle = 60.0\n    smoother = vtk.vtkWindowedSincPolyDataFilter()\n    smoother.SetInputConnection(mcubes.GetOutputPort())\n    smoother.SetNumberOfIterations(smoothingIterations)\n    smoother.BoundarySmoothingOff()\n    smoother.FeatureEdgeSmoothingOff()\n    smoother.SetFeatureAngle(featureAngle)\n    smoother.SetPassBand(passBand)\n    smoother.NonManifoldSmoothingOn()\n    smoother.NormalizeCoordinatesOn()\n    smoother.Update()\n\n    normals = vtk.vtkPolyDataNormals()\n    normals.SetInputConnection(smoother.GetOutputPort())\n    normals.SetFeatureAngle(featureAngle)\n\n    stripper = vtk.vtkStripper()\n    stripper.SetInputConnection(normals.GetOutputPort())\n\n    mapper = vtk.vtkPolyDataMapper()\n    mapper.SetInputConnection(stripper.GetOutputPort())\n\n    actor = vtk.vtkActor()\n    actor.SetMapper(mapper)\n    actor.GetProperty().SetColor( colorLut.GetTableValue(tissueMap[color])[:3])\n    actor.GetProperty().SetSpecular(.5)\n    actor.GetProperty().SetSpecularPower(10)\n    \n    return actor\n\ndef render_lungs(workdir, datadir, patient):\n    PathDicom = datadir + patient\n    reader = vtk.vtkDICOMImageReader()\n    reader.SetDirectoryName(PathDicom)\n    reader.Update()    \n    disp = Display().start()\n    renderer = vtk.vtkRenderer()\n    actor = CreateTissue(reader,-2000,-300,\"lung\", isoValue = 170)\n    renderer.AddActor(actor)\n    renderer.SetBackground(1.0, 1.0, 1.0)\n\n    renderer.ResetCamera()\n    renderer.ResetCameraClippingRange()\n    camera = renderer.GetActiveCamera()\n    camera.Elevation(120)\n    camera.Elevation(120)\n    renderer.SetActiveCamera(camera)\n\n    name = workdir + patient + '_lungs'\n\n    vtk_render_gif(renderer, 1, name, Dolly = 1.5,width = 400, height = 400)\n    disp.stop()\n    gc.collect()","execution_count":null,"outputs":[]},{"metadata":{"id":"TK8idxDjGY6z"},"cell_type":"markdown","source":"# DICOM-CT data 3D visualizations\nLet's move on to the main topic of this notebook. 3D visualization of medical data is a very broad topic. We will start with reading the data:\n\nSources:\n* https://gist.github.com/somada141/38d313a65581341f23fd\n* https://pyscience.wordpress.com/2014/09/11/surface-extraction-creating-a-mesh-from-pixel-data-using-python-and-vtk/\n\n## Reading data\nReading data in VTK is really easy and pleasant, but unfortunately also problematic. This is because this library is focused on visualization, not DICOM processing. This problem can be solved by i.e. importing data into VTK directly from the numpy array, but here we will not focus on this and we will move on using the standard input for DICOM files.\n\nFor the following example, we can load data with only three (!) lines:"},{"metadata":{"trusted":true,"id":"NZcioyS4GY60"},"cell_type":"code","source":"## supporting lines =  tidying up\nworkdir = '/kaggle/working/patients/'\nos.makedirs(workdir, exist_ok = True)\ndatadir = \"/kaggle/input/osic-pulmonary-fibrosis-progression/train/\"\npatients = os.listdir(datadir)\npatients.sort()\npatient = patients[17]\n\n## vtk reading dicom\nreader = vtk.vtkDICOMImageReader()\nreader.SetDirectoryName(datadir + patient)\nreader.Update()","execution_count":null,"outputs":[]},{"metadata":{"id":"eVKy9v2HGY7x"},"cell_type":"markdown","source":"## Windows and LUTs\nChoosing the right parameters for Windowing is a job for radiologists. Therefore - regardless of the 3D rendering algorithm, we will display our data as points in 3D space with standard windowing. The data displayed below use standard LUTs color palettes for Slicer:\n\nSources: \n* https://github.com/pletzer/dicom-data-reader/blob/master/python/readVTK.py\n* https://radiopaedia.org/articles/windowing-ct\n* https://www.slicer.org/wiki/Documentation/4.1/SlicerApplication/LookupTables/GenericAnatomyColors\n* https://github.com/Slicer/Slicer/blob/47b845f91ef332e6008af4a166b294912de6d052/Modules/Loadable/VolumeRendering/Logic/vtkSlicerVolumeRenderingLogic.cxx#L458-L540\n* https://dgobbi.github.io/vtk-dicom/doc/api/image_display.html"},{"metadata":{},"cell_type":"markdown","source":"### Axial view:"},{"metadata":{"trusted":true,"_kg_hide-input":true,"id":"zzFt2IbBGY7y"},"cell_type":"code","source":"windowing = {}\nwindowing['lungs'] = [1500,-600,64,123,147]\nwindowing['mediastinum'] = [350,50,255,244,209]\nwindowing['bones'] = [300,400,177,122,101]\nwindowing['blood'] = [5,80,216,101,79]\n\npatient = \"ID00012637202177665765362\"\nreader.SetDirectoryName(datadir+patient)\nreader.Update()\n\nimageData = reader.GetOutput()\nvolumeMapper = vtk.vtkSmartVolumeMapper()\nvolumeMapper.SetInputData(imageData)\nvolumeProperty = vtk.vtkVolumeProperty()\nvolumeProperty.SetInterpolationType(vtk.VTK_LINEAR_INTERPOLATION)\n\nfor cur_windowing in windowing:\n    cur_w = windowing[cur_windowing]\n    opacity_function = vtk.vtkPiecewiseFunction()\n    opacity_function.AddPoint(cur_w[1]-cur_w[0]/2,   0.0)\n    opacity_function.AddPoint(cur_w[1],   1.0)\n    opacity_function.AddPoint(cur_w[1]+cur_w[0]/2,   0.0)\n    volumeProperty.SetScalarOpacity(opacity_function)\n\n    color_function = vtk.vtkColorTransferFunction()\n    color_function.SetColorSpaceToDiverging()\n    color_function.AddRGBPoint(cur_w[1]-cur_w[0]/2,0,0,0)\n    color_function.AddRGBPoint(cur_w[1],cur_w[2],cur_w[3],cur_w[4])\n    color_function.AddRGBPoint(cur_w[1]+cur_w[0]/2, 0,0,0)\n\n    volumeProperty.SetColor(color_function)\n\n    volume = vtk.vtkVolume()\n    volume.SetMapper(volumeMapper)\n    volume.SetProperty(volumeProperty)\n\n    disp = Display().start()\n    renderer = vtk.vtkRenderer();\n    volumeMapper.SetRequestedRenderModeToRayCast()\n    renderer.AddViewProp(volume)\n    \n    renderer.ResetCamera()\n    renderer.SetBackground(1,1,1);\n    renderer.ResetCamera()\n    renderer.ResetCameraClippingRange()\n    camera = renderer.MakeCamera()\n\n    camera.SetPosition(0,0,0)\n    camera = renderer.GetActiveCamera()\n    camera.Dolly(1.5)\n\n\n    camera.Roll(360)\n    name = workdir + patient + cur_windowing + '_top'\n    vtk_render_gif(renderer, N = 1 ,name =  name, standard = False)\n    \n    name = workdir + patient + cur_windowing + '_front'\n    camera.Elevation(240)\n    camera.Elevation(20)\n    vtk_render_gif(renderer, N = 1 ,name =  name, standard = False)\n    \n    disp.stop()\n    \nplt.rcParams[\"figure.figsize\"] = (40,40)\nidp = 0 \nfor cur_windowing in windowing:\n    idp += 1\n    plt.subplot(len(windowing),4, idp)\n    try:\n        im = mpimg.imread( workdir + patient + cur_windowing+'_top.gif') \n        plt.imshow(im) \n        plt.title('Windowing: ' + cur_windowing, fontsize =20)\n    except:\n        pass","execution_count":null,"outputs":[]},{"metadata":{},"cell_type":"markdown","source":"### Coronal view:"},{"metadata":{"_kg_hide-input":true,"trusted":true,"id":"50Bz5vHoGY70"},"cell_type":"code","source":"plt.rcParams[\"figure.figsize\"] = (40,40)\nidp = 0 \nfor cur_windowing in windowing:\n    idp += 1\n    plt.subplot(len(windowing),4, idp)\n    try:\n        im = mpimg.imread( workdir + patient + cur_windowing+'_front.gif') \n        plt.imshow(im) \n        plt.title('Windowing: ' + cur_windowing, fontsize =20)\n    except:\n        pass","execution_count":null,"outputs":[]},{"metadata":{"id":"fDEqAOzsGY7Q"},"cell_type":"markdown","source":"## 3D reconstruction with manual windowing and MarchingCubes algorithm\nFinally, 3D visualization. We assume that the reader has knowledge about windowing and 3D rendering algorithms (at this point we do not want to discuss the details related to it). We will start with the MarchingCubes algorithm with standard windowing for lungs, skeleton and blood:"},{"metadata":{"id":"QOhoVB2zGY7U"},"cell_type":"markdown","source":"### Front view (skeleton + lungs + blood)"},{"metadata":{"_kg_hide-input":true,"trusted":true,"id":"XFerCPmhGY7V"},"cell_type":"code","source":"%%time\n\npatient = patients[17]\nreader = vtk.vtkDICOMImageReader()\nreader.SetDirectoryName(datadir+patient)\nreader.Update()\n\ndisp = Display().start()\nrenderer = vtk.vtkRenderer()\nrenderer.AddActor(CreateTissue(reader,-900,-400,\"lung\"))\nrenderer.AddActor(CreateTissue(reader,0,120,\"blood\"))\nrenderer.AddActor(CreateTissue(reader,100,2000,\"skeleton\"))\nrenderer.SetBackground(1.0, 1.0, 1.0)\n\nrenderer.ResetCamera()\nrenderer.ResetCameraClippingRange()\ncamera = renderer.GetActiveCamera()\ncamera.Elevation(120)\ncamera.Roll(180)\nrenderer.SetActiveCamera(camera)\n\nname = workdir + patient + \"_front\"\nvtk_render_gif(renderer, 1, name, Dolly = 1.5)\ndisp.stop()\n\nImage(filename=name + \".gif\", format='png')    ","execution_count":null,"outputs":[]},{"metadata":{"id":"V7yRjsCGGY7j"},"cell_type":"markdown","source":"### Back view (skeleton + lungs + blood)"},{"metadata":{"_kg_hide-input":true,"trusted":true,"id":"_RoRfY0sGY7k"},"cell_type":"code","source":"%%time\ndisp = Display().start()\nrenderer = vtk.vtkRenderer()\nrenderer.AddActor(CreateTissue(reader,-900,-400,\"lung\"))\nrenderer.AddActor(CreateTissue(reader,0,120,\"blood\"))\nrenderer.AddActor(CreateTissue(reader,100,2000,\"skeleton\"))\n\nrenderer.SetBackground(1.0, 1.0, 1.0)\n\nrenderer.ResetCamera()\nrenderer.ResetCameraClippingRange()\ncamera = renderer.GetActiveCamera()\ncamera.Elevation(120)\ncamera.Elevation(120)\ncamera.Roll(180)\nrenderer.SetActiveCamera(camera)\n\nname = workdir + patient + \"_back\"\nvtk_render_gif(renderer, 1, name, Dolly = 1.5)\ndisp.stop()\n\nImage(filename=name + \".gif\", format='png')    ","execution_count":null,"outputs":[]},{"metadata":{"id":"anEql2NpGY7n"},"cell_type":"markdown","source":"## Lungs segmentation\nOne of the reasons why we do such visualizations is to check if the segmentation went well. After adjustment: WindowWidth, WindowCenter and IsoLevel (for the MarchingCubes algorithm), the selected lungs may look like this:"},{"metadata":{"trusted":true,"_kg_hide-input":true,"id":"tdNyTggXGY7o"},"cell_type":"code","source":"%%time\ndisp = Display().start()\nrenderer = vtk.vtkRenderer()\nactor = CreateTissue(reader,-2000,-300,\"lung\", isoValue = 170)\nrenderer.AddActor(actor)\nrenderer.SetBackground(1.0, 1.0, 1.0)\n\nrenderer.ResetCamera()\nrenderer.ResetCameraClippingRange()\ncamera = renderer.GetActiveCamera()\ncamera.Elevation(120)\ncamera.Elevation(120)\nrenderer.SetActiveCamera(camera)\n\nname = workdir + patient + '_lungs'\n\nvtk_render_gif(renderer, 1, name, Dolly = 1.5)\ndisp.stop()\n\nImage(filename=name + \".gif\", format='png')","execution_count":null,"outputs":[]},{"metadata":{"id":"cFWlD8EvGY7s"},"cell_type":"markdown","source":"Awesome! But that does not mean that such a combination of parameters will work for all patients. We will perform lung rendering for the pre-selected patients from the OSIC training set:"},{"metadata":{"trusted":true,"_kg_hide-input":true,"id":"j_Ci_AjoGY7t"},"cell_type":"code","source":"%%time\ns_patients = [patients[i] for i in [1,2,4,6,8,12,13,14,15,16]]\nfor patient in s_patients:\n    try:\n        render_lungs(workdir, datadir, patient)\n        print(patient + ' completed render lungs')\n    except:\n        print(patient + ' failed render lungs')\n\nplt.rcParams[\"figure.figsize\"] = (40,120)\nidp = 0 \nfor patient in s_patients:\n    idp += 1\n    plt.subplot(10,2, idp)\n    try:\n        im = mpimg.imread( workdir + patient + '_lungs.gif') \n        plt.imshow(im) \n        plt.title('OSIC PatientID: '+patient, fontsize=20)\n    except:\n        pass","execution_count":null,"outputs":[]},{"metadata":{"id":"m8cLwAnGGY7w"},"cell_type":"markdown","source":"We pre-selected this patients, because it turns out that for some patients WindowWidth and WindowCenter are inappropriate. And for some of them (probably) the import went wrong. It is worth adding that it is possible to import data from the numpy array (the method presented above was the export to numpy)."},{"metadata":{"id":"8oGZ5H4nGY64"},"cell_type":"markdown","source":"## Other stuff: Vtk to numpy support\nThis data from reader instance can be exported to a numpy array:"},{"metadata":{"trusted":true,"id":"jXINkagoGY67"},"cell_type":"code","source":"from vtk.util import numpy_support\nimport numpy\n\nreader = vtk.vtkDICOMImageReader()\nreader.SetDirectoryName(datadir + patient)\nreader.Update()\n# Load dimensions using `GetDataExtent`\n_extent = reader.GetDataExtent()\nConstPixelDims = [_extent[1]-_extent[0]+1, _extent[3]-_extent[2]+1, _extent[5]-_extent[4]+1]\n\n# Load spacing values\nConstPixelSpacing = reader.GetPixelSpacing()\n\n# Get the 'vtkImageData' object from the reader\nimageData = reader.GetOutput()\n# Get the 'vtkPointData' object from the 'vtkImageData' object\npointData = imageData.GetPointData()\n# Ensure that only one array exists within the 'vtkPointData' object\nassert (pointData.GetNumberOfArrays()==1)\n# Get the `vtkArray` (or whatever derived type) which is needed for the `numpy_support.vtk_to_numpy` function\narrayData = pointData.GetArray(0)\n\n# Convert the `vtkArray` to a NumPy array\nArrayDicom = numpy_support.vtk_to_numpy(arrayData)\n# Reshape the NumPy array to 3D using 'ConstPixelDims' as a 'shape'\nArrayDicom = ArrayDicom.reshape(ConstPixelDims, order='F')\nArrayDicom.shape","execution_count":null,"outputs":[]},{"metadata":{"id":"sdDzsHGeGY6_"},"cell_type":"markdown","source":"### Checking if everything is OK ;)\nNow, the DICOM import by VTK includes the conversion to HU units. To check if everything is okay, we will import data manually and compare the averages (for performance we do not compare every element of the matrix).\n\nSource: \n* https://www.kaggle.com/allunia/pulmonary-dicom-preprocessing"},{"metadata":{"_kg_hide-input":true,"trusted":true,"id":"Hg7aDNfIGY7A"},"cell_type":"code","source":"## source https://www.kaggle.com/allunia/pulmonary-dicom-preprocessing\n\nbasepath = \"\"\nfrom os import listdir\nimport pydicom\nimport numpy as np\n\ndef load_scans(dcm_path):\n    if basepath == \"../input/osic-pulmonary-fibrosis-progression/\":\n        # in this competition we have missing values in ImagePosition, this is why we are sorting by filename number\n        files = listdir(dcm_path)\n        file_nums = [np.int(file.split(\".\")[0]) for file in files]\n        sorted_file_nums = np.sort(file_nums)[::-1]\n        slices = [pydicom.dcmread(dcm_path + \"/\" + str(file_num) + \".dcm\" ) for file_num in sorted_file_nums]\n    else:\n        # otherwise we sort by ImagePositionPatient (z-coordinate) or by SliceLocation\n        slices = [pydicom.dcmread(dcm_path + \"/\" + file) for file in listdir(dcm_path)]\n        slices.sort(key = lambda x: float(x.ImagePositionPatient[2]))\n    return slices\ndef transform_to_hu(slices):\n    images = np.stack([file.pixel_array for file in slices])\n    images = images.astype(np.int16)\n\n    #images = set_outside_scanner_to_air(images)\n    \n    # convert to HU\n    for n in range(len(slices)):\n        \n        intercept = slices[n].RescaleIntercept\n        slope = slices[n].RescaleSlope\n        \n        if slope != 1:\n            images[n] = slope * images[n].astype(np.float64)\n            images[n] = images[n].astype(np.int16)\n            \n        images[n] += np.int16(intercept)\n    \n    return np.array(images, dtype=np.int16)\n","execution_count":null,"outputs":[]},{"metadata":{"trusted":true,"id":"wLzQUMrpGY7G"},"cell_type":"code","source":"scans = load_scans(datadir + patient)\nhu_scans = transform_to_hu(scans)\nnp.mean(hu_scans) == np.mean(ArrayDicom)","execution_count":null,"outputs":[]},{"metadata":{"id":"ryq7jLBZGY7K"},"cell_type":"markdown","source":"In any case, conversion can be done anyway:"},{"metadata":{"trusted":true,"_kg_hide-input":true,"id":"RkrMXf3BGY7L"},"cell_type":"code","source":"### if needed:\n#shiftScale = vtk.vtkImageShiftScale()\n#shiftScale.SetScale(reader.GetRescaleSlope())\n#shiftScale.SetShift(reader.GetRescaleOffset())\n#shiftScale.SetInputConnection(reader.GetOutputPort())\n#shiftScale.Update()","execution_count":null,"outputs":[]},{"metadata":{"id":"wVMl5vVCGY5C"},"cell_type":"markdown","source":"# Gentle introduction to VTK\nThe Visualization Toolkit (VTK) is an open-source, freely available software system for 3D computer graphics, modeling, image processing, volume rendering, scientific visualization, and 2D plotting.\n\nSource:\n* https://vtk.org/"},{"metadata":{"id":"U5AlIh9IGY5E"},"cell_type":"markdown","source":"## VTK in standalone Kaggle notebooks\nSince there is no x-server on Kaggle, it is not possible to use the VTK library directly. With **pyvirtualdisplay** which imitates x-server, VTK usage is available.\n\nSources: \n* https://pyscience.wordpress.com/2014/09/03/ipython-notebook-vtk/\n* https://pyscience.wordpress.com/2014/11/16/volume-rendering-with-python-and-vtk/\n* https://nbviewer.jupyter.org/urls/bitbucket.org/somada141/pyscience/raw/master/20141029_VolumeRendering/Material/VolumeRendering.ipynb"},{"metadata":{"id":"BEWVU0EaGY5G"},"cell_type":"markdown","source":"In order to use pyvirt you just need to install it via pip and then turn it on (disp.start()) and off (disp.stop()) x-server:"},{"metadata":{"trusted":true,"id":"_btbnoZ-GY5H"},"cell_type":"code","source":"! pip install pyvirtualdisplay -q\nfrom pyvirtualdisplay import Display\n\ndisp = Display().start()\nimport vtk\ndisp.stop()","execution_count":null,"outputs":[]},{"metadata":{"id":"9ZTyfY7iGY5Q"},"cell_type":"markdown","source":"After a few days of using VTK, the main disadvantage seems to be that it is overdesigned (in my humble opinion). For example, this is how we can read the major number of the installed VTK library (and this is just the beginning) :)"},{"metadata":{"trusted":true,"id":"I2tXIvaaGY5R"},"cell_type":"code","source":"vtk.vtkVersion().GetVTKMajorVersion()","execution_count":null,"outputs":[]},{"metadata":{"id":"oM7TXcT4GY5Z"},"cell_type":"markdown","source":"Using a virtual x-server requires that we save the images immediately to the memory/disk. We will use the following function for this:"},{"metadata":{"trusted":true,"id":"Sdp5n5tPGY5a"},"cell_type":"code","source":"def vtk_show(renderer, width = default_width, height = default_height):\n\n    renderWindow = vtk.vtkRenderWindow()\n    \n    renderWindow.SetOffScreenRendering(1)\n    renderWindow.AddRenderer(renderer)\n    renderWindow.SetSize(width, height)\n    renderWindow.Render()\n     \n    windowToImageFilter = vtk.vtkWindowToImageFilter()\n    windowToImageFilter.SetInput(renderWindow)\n    windowToImageFilter.Update()\n     \n    writer = vtk.vtkPNGWriter()\n    writer.SetWriteToMemory(1)\n    writer.SetInputConnection(windowToImageFilter.GetOutputPort())\n    writer.Write()    \n    \n    return bytes(memoryview(writer.GetResult()))","execution_count":null,"outputs":[]},{"metadata":{"id":"urwmYMMJGY5i"},"cell_type":"markdown","source":"## Cylinder = Hello World example\nThis example creates a minimal visualization program, demonstrating VTK's basic rendering and pipeline creation. Notice, that standard pipeline looks like this: Source -> Mapper -> Actor -> Renderer\n\nSource: \n* https://lorensen.github.io/VTKExamples/site/Python/GeometricObjects/CylinderExample/"},{"metadata":{"id":"CG7nzwRTGY5k"},"cell_type":"markdown","source":"We will start with virtual x-server initialization:"},{"metadata":{"trusted":true,"id":"M1iLOdaPGY5l"},"cell_type":"code","source":"disp = Display().start()","execution_count":null,"outputs":[]},{"metadata":{"id":"AmOSC_Q0GY5s"},"cell_type":"markdown","source":"and background color definition ;)"},{"metadata":{"trusted":true,"id":"S3ukJtZtGY5t"},"cell_type":"code","source":"colors = vtk.vtkNamedColors()\n# Set the background color.\nbkg = map(lambda x: x / 255.0, [26, 51, 102, 255])\ncolors.SetColor(\"BkgColor\", *bkg)","execution_count":null,"outputs":[]},{"metadata":{"id":"P6Jd-bi3GY5x"},"cell_type":"markdown","source":"### Source\nThis creates a polygonal cylinder model with eight circumferential facets:"},{"metadata":{"trusted":true,"id":"M03jOwLtGY5y"},"cell_type":"code","source":"cylinder = vtk.vtkCylinderSource()\ncylinder.SetResolution(8)","execution_count":null,"outputs":[]},{"metadata":{"id":"XbOF_82UGY54"},"cell_type":"markdown","source":"### Mapper\nThe mapper is responsible for pushing the geometry into the graphics library. It may also do color mapping, if scalars or other attributes are defined:"},{"metadata":{"trusted":true,"id":"5QIK3_0jGY55"},"cell_type":"code","source":"cylinderMapper = vtk.vtkPolyDataMapper()\ncylinderMapper.SetInputConnection(cylinder.GetOutputPort())","execution_count":null,"outputs":[]},{"metadata":{"id":"VvE-EQdtGY5_"},"cell_type":"markdown","source":"### Actor\nThe actor is a grouping mechanism: besides the geometry (mapper), it also has a property, transformation matrix, and/or texture map. Here we set its color and rotate it:"},{"metadata":{"trusted":true,"id":"fFxAzUpqGY5_"},"cell_type":"code","source":"cylinderActor = vtk.vtkActor()\ncylinderActor.SetMapper(cylinderMapper)\ncylinderActor.GetProperty().SetColor(colors.GetColor3d(\"Tomato\"))\ncylinderActor.RotateX(30.0)\ncylinderActor.RotateY(-45.0)","execution_count":null,"outputs":[]},{"metadata":{"id":"yhkoMse2GY6G"},"cell_type":"markdown","source":"### Renderer\nCreate the graphics structure. The renderer renders into the render window. The render window interactor captures mouse events and will perform appropriate camera or actor manipulation depending on the nature of the events."},{"metadata":{"trusted":true,"id":"ZwRXcr55GY6H"},"cell_type":"code","source":"ren = vtk.vtkRenderer()\n# useless lines: in kaggle notebook\n#renWin = vtk.vtkRenderWindow() \n#renWin.AddRenderer(ren)\n#iren = vtk.vtkRenderWindowInteractor()\n#iren.SetRenderWindow(renWin)","execution_count":null,"outputs":[]},{"metadata":{"id":"bTEbY_lVGY6O"},"cell_type":"markdown","source":"Add the actors to the renderer, set the background and size:"},{"metadata":{"trusted":true,"id":"ZL2lOyh1GY6P"},"cell_type":"code","source":"ren.AddActor(cylinderActor)\nren.SetBackground(colors.GetColor3d(\"BkgColor\"))\n\nren.ResetCamera()\nren.GetActiveCamera().Zoom(1.5)\n#renWin.Render() # useless line in kaggle nb","execution_count":null,"outputs":[]},{"metadata":{"id":"TsOxyd_6GY6T"},"cell_type":"markdown","source":"Finally, gather results, stop dummy x-server and display image:"},{"metadata":{"trusted":true,"id":"EVn4p2CdGY6T"},"cell_type":"code","source":"img = vtk_show(ren)\ndisp.stop()\nImage(img)","execution_count":null,"outputs":[]},{"metadata":{"id":"AK6QZu53GY6Z"},"cell_type":"markdown","source":"## Cylinder rotation with animation\nAnimation of the above example is simply a rotation of the actor (or camera if you want) and saving to a list of result images, i.e. to a .gif file:"},{"metadata":{"_kg_hide-input":true,"trusted":true,"id":"OKTIsmxVGY6b"},"cell_type":"code","source":"def vtk_show(renderer, width = default_width, height = default_height, filename = \"\"):\n\n    renderWindow = vtk.vtkRenderWindow()\n    \n    renderWindow.SetOffScreenRendering(1)\n    renderWindow.AddRenderer(renderer)\n    renderWindow.SetSize(width, height)\n    renderWindow.Render()\n     \n    windowToImageFilter = vtk.vtkWindowToImageFilter()\n    windowToImageFilter.SetInput(renderWindow)\n    windowToImageFilter.Update()\n     \n    writer = vtk. vtkPNGWriter()\n    \n    if filename == \"\":\n        writer.SetWriteToMemory(1)\n        writer.SetInputConnection(windowToImageFilter.GetOutputPort())\n        writer.Write()    \n        return bytes(memoryview(writer.GetResult()))\n    else:\n        writer.SetFileName(filename+\".png\")\n        writer.SetInputConnection(windowToImageFilter.GetOutputPort())\n        writer.Write()    \n        return None\n    \ndef vtk_render_gif(renderer, N, name, Roll = False, Azimuth = False, Elevation = False, Actor = None, RotateX = False, RotateY = False, RotateZ = False, Zoom = 0, Dolly = 0, standard = True, width = default_width, height = default_height):    \n    if standard:\n        renderer.ResetCamera()\n        camera = renderer.MakeCamera()\n        renderer.ResetCameraClippingRange()\n        camera.SetPosition(0,0,0)\n    os.makedirs(name,exist_ok=True)\n    \n    if Zoom != 0:\n        renderer.GetActiveCamera().Zoom(Zoom)\n        \n    if Dolly != 0:\n        renderer.GetActiveCamera().Dolly(Dolly)\n        \n    #tmpN = 1\n    if N >0: # render gif\n        for fi in range(N):\n            if Roll:\n                renderer.GetActiveCamera().Roll(360//N) \n            if Azimuth:\n                renderer.GetActiveCamera().Azimuth(360//N) \n            if Elevation:\n                renderer.GetActiveCamera().Elevation(360//N)\n            if Actor is not None:\n                if RotateX:\n                    Actor.RotateX(360//N)\n                if RotateY:\n                    Actor.RotateY(360//N)\n                if RotateZ:\n                    Actor.RotateZ(360//N)                    \n            vtk_show(renderer,filename = name + \"/shot\"+str(fi), width = width, height = height)\n        # render gif and cleanup\n        img_list = []\n        for fi in range(N):\n            img_list.append(mpimg.imread(name + '/shot' + str(fi) + '.png'))\n        shutil.rmtree(name)\n        imageio.mimsave(name + \".gif\", img_list, duration=0.5)\n\n    #if N == 1: # render png\n       #vtk_show(renderer,filename = name + \".gif\")\n","execution_count":null,"outputs":[]},{"metadata":{"_kg_hide-input":true,"trusted":true,"id":"C_wjV1Q1GY6g"},"cell_type":"code","source":"%%time\nworkdir = '/kaggle/working/vtk_examples/'\nname = workdir + 'cylinder'\nos.makedirs(workdir,exist_ok=True)\ndisp = Display().start()\n\nren = vtk.vtkRenderer()\nren.AddActor(cylinderActor)\nren.SetBackground(colors.GetColor3d(\"BkgColor\"))\n\nvtk_render_gif(ren, N, name, Actor = cylinderActor, RotateX = True, Zoom = 1.5)\n\ndisp.stop()\n\nImage(filename=name + \".gif\", format='png')","execution_count":null,"outputs":[]},{"metadata":{"id":"b4trhSmoGY6k"},"cell_type":"markdown","source":"## Kitchen example\nThe example shows forty streamlines in a small kitchen. The room has two windows, a door(with air leakage), and a cooking area with a hot stove. The air leakage and temperature variation combine to produce air convection currents throughout the kitchen. The starting positions of the streamlines were defined by creating a rake, or curve (and its associated points). Here the rake was a straight line modeled with a vtkLineSource. These streamlines clearly show features of the flow field. By releasing many streamlines simultaneously we obtain even more information, as the eye tends to assemble nearby streamlines into a “global” understanding of flow field features:\n\nSources: \n* https://lorensen.github.io/VTKExamples/site/Python/Visualization/Kitchen/\n* https://www.youtube.com/watch?v=Uly_PRBoLgM"},{"metadata":{"_kg_hide-input":true,"trusted":true,"id":"biEWVr08GY6l"},"cell_type":"code","source":"%%time\n! wget -q https://raw.githubusercontent.com/lorensen/VTKExamples/master/src/Testing/Data/kitchen.vtk -O /kaggle/working/vtk_examples/kitchen.vtk\ndisp = Display().start()\nfileName = workdir +'kitchen.vtk'\ncolors = vtk.vtkNamedColors()\n# Set the furniture colors.\ncolors.SetColor(\"Furniture\", [204, 204, 153, 255])\n\nscalarRange = [0.0, 0.0]\nmaxTime = 0\n\naren = vtk.vtkRenderer()\n\n#\n# Read the data.\n#\nreader = vtk.vtkStructuredGridReader()\nreader.SetFileName(fileName)\nreader.Update()  # Force a read to occur.\nreader.GetOutput().GetLength()\n\nif reader.GetOutput().GetPointData().GetScalars():\n    reader.GetOutput().GetPointData().GetScalars().GetRange(scalarRange)\n\nif reader.GetOutput().GetPointData().GetVectors():\n    maxVelocity = reader.GetOutput().GetPointData().GetVectors().GetMaxNorm()\n    maxTime = 4.0 * reader.GetOutput().GetLength() / maxVelocity\n\n#\n# Outline around the data.\n#\noutlineF = vtk.vtkStructuredGridOutlineFilter()\noutlineF.SetInputConnection(reader.GetOutputPort())\noutlineMapper = vtk.vtkPolyDataMapper()\noutlineMapper.SetInputConnection(outlineF.GetOutputPort())\noutline = vtk.vtkActor()\noutline.SetMapper(outlineMapper)\noutline.GetProperty().SetColor(colors.GetColor3d(\"LampBlack\"))\n\n#\n# Set up shaded surfaces (i.e., supporting geometry).\n#\ndoorGeom = vtk.vtkStructuredGridGeometryFilter()\ndoorGeom.SetInputConnection(reader.GetOutputPort())\ndoorGeom.SetExtent(27, 27, 14, 18, 0, 11)\nmapDoor = vtk.vtkPolyDataMapper()\nmapDoor.SetInputConnection(doorGeom.GetOutputPort())\nmapDoor.ScalarVisibilityOff()\ndoor = vtk.vtkActor()\ndoor.SetMapper(mapDoor)\ndoor.GetProperty().SetColor(colors.GetColor3d(\"Burlywood\"))\n\nwindow1Geom = vtk.vtkStructuredGridGeometryFilter()\nwindow1Geom.SetInputConnection(reader.GetOutputPort())\nwindow1Geom.SetExtent(0, 0, 9, 18, 6, 12)\nmapWindow1 = vtk.vtkPolyDataMapper()\nmapWindow1.SetInputConnection(window1Geom.GetOutputPort())\nmapWindow1.ScalarVisibilityOff()\nwindow1 = vtk.vtkActor()\nwindow1.SetMapper(mapWindow1)\nwindow1.GetProperty().SetColor(colors.GetColor3d(\"SkyBlue\"))\nwindow1.GetProperty().SetOpacity(.6)\n\nwindow2Geom = vtk.vtkStructuredGridGeometryFilter()\nwindow2Geom.SetInputConnection(reader.GetOutputPort())\nwindow2Geom.SetExtent(5, 12, 23, 23, 6, 12)\nmapWindow2 = vtk.vtkPolyDataMapper()\nmapWindow2.SetInputConnection(window2Geom.GetOutputPort())\nmapWindow2.ScalarVisibilityOff()\nwindow2 = vtk.vtkActor()\nwindow2.SetMapper(mapWindow2)\nwindow2.GetProperty().SetColor(colors.GetColor3d(\"SkyBlue\"))\nwindow2.GetProperty().SetOpacity(.6)\n\nklower1Geom = vtk.vtkStructuredGridGeometryFilter()\nklower1Geom.SetInputConnection(reader.GetOutputPort())\nklower1Geom.SetExtent(17, 17, 0, 11, 0, 6)\nmapKlower1 = vtk.vtkPolyDataMapper()\nmapKlower1.SetInputConnection(klower1Geom.GetOutputPort())\nmapKlower1.ScalarVisibilityOff()\nklower1 = vtk.vtkActor()\nklower1.SetMapper(mapKlower1)\nklower1.GetProperty().SetColor(colors.GetColor3d(\"EggShell\"))\n\nklower2Geom = vtk.vtkStructuredGridGeometryFilter()\nklower2Geom.SetInputConnection(reader.GetOutputPort())\nklower2Geom.SetExtent(19, 19, 0, 11, 0, 6)\nmapKlower2 = vtk.vtkPolyDataMapper()\nmapKlower2.SetInputConnection(klower2Geom.GetOutputPort())\nmapKlower2.ScalarVisibilityOff()\nklower2 = vtk.vtkActor()\nklower2.SetMapper(mapKlower2)\nklower2.GetProperty().SetColor(colors.GetColor3d(\"EggShell\"))\n\nklower3Geom = vtk.vtkStructuredGridGeometryFilter()\nklower3Geom.SetInputConnection(reader.GetOutputPort())\nklower3Geom.SetExtent(17, 19, 0, 0, 0, 6)\nmapKlower3 = vtk.vtkPolyDataMapper()\nmapKlower3.SetInputConnection(klower3Geom.GetOutputPort())\nmapKlower3.ScalarVisibilityOff()\nklower3 = vtk.vtkActor()\nklower3.SetMapper(mapKlower3)\nklower3.GetProperty().SetColor(colors.GetColor3d(\"EggShell\"))\n\nklower4Geom = vtk.vtkStructuredGridGeometryFilter()\nklower4Geom.SetInputConnection(reader.GetOutputPort())\nklower4Geom.SetExtent(17, 19, 11, 11, 0, 6)\nmapKlower4 = vtk.vtkPolyDataMapper()\nmapKlower4.SetInputConnection(klower4Geom.GetOutputPort())\nmapKlower4.ScalarVisibilityOff()\nklower4 = vtk.vtkActor()\nklower4.SetMapper(mapKlower4)\nklower4.GetProperty().SetColor(colors.GetColor3d(\"EggShell\"))\n\nklower5Geom = vtk.vtkStructuredGridGeometryFilter()\nklower5Geom.SetInputConnection(reader.GetOutputPort())\nklower5Geom.SetExtent(17, 19, 0, 11, 0, 0)\nmapKlower5 = vtk.vtkPolyDataMapper()\nmapKlower5.SetInputConnection(klower5Geom.GetOutputPort())\nmapKlower5.ScalarVisibilityOff()\nklower5 = vtk.vtkActor()\nklower5.SetMapper(mapKlower5)\nklower5.GetProperty().SetColor(colors.GetColor3d(\"EggShell\"))\n\nklower6Geom = vtk.vtkStructuredGridGeometryFilter()\nklower6Geom.SetInputConnection(reader.GetOutputPort())\nklower6Geom.SetExtent(17, 19, 0, 7, 6, 6)\nmapKlower6 = vtk.vtkPolyDataMapper()\nmapKlower6.SetInputConnection(klower6Geom.GetOutputPort())\nmapKlower6.ScalarVisibilityOff()\nklower6 = vtk.vtkActor()\nklower6.SetMapper(mapKlower6)\nklower6.GetProperty().SetColor(colors.GetColor3d(\"EggShell\"))\n\nklower7Geom = vtk.vtkStructuredGridGeometryFilter()\nklower7Geom.SetInputConnection(reader.GetOutputPort())\nklower7Geom.SetExtent(17, 19, 9, 11, 6, 6)\nmapKlower7 = vtk.vtkPolyDataMapper()\nmapKlower7.SetInputConnection(klower7Geom.GetOutputPort())\nmapKlower7.ScalarVisibilityOff()\nklower7 = vtk.vtkActor()\nklower7.SetMapper(mapKlower7)\nklower7.GetProperty().SetColor(colors.GetColor3d(\"EggShell\"))\n\nhood1Geom = vtk.vtkStructuredGridGeometryFilter()\nhood1Geom.SetInputConnection(reader.GetOutputPort())\nhood1Geom.SetExtent(17, 17, 0, 11, 11, 16)\nmapHood1 = vtk.vtkPolyDataMapper()\nmapHood1.SetInputConnection(hood1Geom.GetOutputPort())\nmapHood1.ScalarVisibilityOff()\nhood1 = vtk.vtkActor()\nhood1.SetMapper(mapHood1)\nhood1.GetProperty().SetColor(colors.GetColor3d(\"Silver\"))\n\nhood2Geom = vtk.vtkStructuredGridGeometryFilter()\nhood2Geom.SetInputConnection(reader.GetOutputPort())\nhood2Geom.SetExtent(19, 19, 0, 11, 11, 16)\nmapHood2 = vtk.vtkPolyDataMapper()\nmapHood2.SetInputConnection(hood2Geom.GetOutputPort())\nmapHood2.ScalarVisibilityOff()\nhood2 = vtk.vtkActor()\nhood2.SetMapper(mapHood2)\nhood2.GetProperty().SetColor(colors.GetColor3d(\"Furniture\"))\n\nhood3Geom = vtk.vtkStructuredGridGeometryFilter()\nhood3Geom.SetInputConnection(reader.GetOutputPort())\nhood3Geom.SetExtent(17, 19, 0, 0, 11, 16)\nmapHood3 = vtk.vtkPolyDataMapper()\nmapHood3.SetInputConnection(hood3Geom.GetOutputPort())\nmapHood3.ScalarVisibilityOff()\nhood3 = vtk.vtkActor()\nhood3.SetMapper(mapHood3)\nhood3.GetProperty().SetColor(colors.GetColor3d(\"Furniture\"))\n\nhood4Geom = vtk.vtkStructuredGridGeometryFilter()\nhood4Geom.SetInputConnection(reader.GetOutputPort())\nhood4Geom.SetExtent(17, 19, 11, 11, 11, 16)\nmapHood4 = vtk.vtkPolyDataMapper()\nmapHood4.SetInputConnection(hood4Geom.GetOutputPort())\nmapHood4.ScalarVisibilityOff()\nhood4 = vtk.vtkActor()\nhood4.SetMapper(mapHood4)\nhood4.GetProperty().SetColor(colors.GetColor3d(\"Furniture\"))\n\nhood6Geom = vtk.vtkStructuredGridGeometryFilter()\nhood6Geom.SetInputConnection(reader.GetOutputPort())\nhood6Geom.SetExtent(17, 19, 0, 11, 16, 16)\nmapHood6 = vtk.vtkPolyDataMapper()\nmapHood6.SetInputConnection(hood6Geom.GetOutputPort())\nmapHood6.ScalarVisibilityOff()\nhood6 = vtk.vtkActor()\nhood6.SetMapper(mapHood6)\nhood6.GetProperty().SetColor(colors.GetColor3d(\"Furniture\"))\n\ncookingPlateGeom = vtk.vtkStructuredGridGeometryFilter()\ncookingPlateGeom.SetInputConnection(reader.GetOutputPort())\ncookingPlateGeom.SetExtent(17, 19, 7, 9, 6, 6)\nmapCookingPlate = vtk.vtkPolyDataMapper()\nmapCookingPlate.SetInputConnection(cookingPlateGeom.GetOutputPort())\nmapCookingPlate.ScalarVisibilityOff()\ncookingPlate = vtk.vtkActor()\ncookingPlate.SetMapper(mapCookingPlate)\ncookingPlate.GetProperty().SetColor(colors.GetColor3d(\"Tomato\"))\n\nfilterGeom = vtk.vtkStructuredGridGeometryFilter()\nfilterGeom.SetInputConnection(reader.GetOutputPort())\nfilterGeom.SetExtent(17, 19, 7, 9, 11, 11)\nmapFilter = vtk.vtkPolyDataMapper()\nmapFilter.SetInputConnection(filterGeom.GetOutputPort())\nmapFilter.ScalarVisibilityOff()\nsgfilter = vtk.vtkActor()\nsgfilter.SetMapper(mapFilter)\nsgfilter.GetProperty().SetColor(colors.GetColor3d(\"Furniture\"))\n#\n# regular streamlines\n#\nline = vtk.vtkLineSource()\nline.SetResolution(39)\nline.SetPoint1(0.08, 2.50, 0.71)\nline.SetPoint2(0.08, 4.50, 0.71)\nrakeMapper = vtk.vtkPolyDataMapper()\nrakeMapper.SetInputConnection(line.GetOutputPort())\nrake = vtk.vtkActor()\nrake.SetMapper(rakeMapper)\n\nstreamers = vtk.vtkStreamTracer()\n# streamers.DebugOn()\nstreamers.SetInputConnection(reader.GetOutputPort())\nstreamers.SetSourceConnection(line.GetOutputPort())\nstreamers.SetMaximumPropagation(maxTime)\nstreamers.SetInitialIntegrationStep(.5)\nstreamers.SetMinimumIntegrationStep(.1)\nstreamers.SetIntegratorType(2)\nstreamers.Update()\n\nstreamersMapper = vtk.vtkPolyDataMapper()\nstreamersMapper.SetInputConnection(streamers.GetOutputPort())\nstreamersMapper.SetScalarRange(scalarRange)\n\nlines = vtk.vtkActor()\nlines.SetMapper(streamersMapper)\nlines.GetProperty().SetColor(colors.GetColor3d(\"Black\"))\n\naren.TwoSidedLightingOn()\n\naren.AddActor(outline)\naren.AddActor(door)\naren.AddActor(window1)\naren.AddActor(window2)\naren.AddActor(klower1)\naren.AddActor(klower2)\naren.AddActor(klower3)\naren.AddActor(klower4)\naren.AddActor(klower5)\naren.AddActor(klower6)\naren.AddActor(klower7)\naren.AddActor(hood1)\naren.AddActor(hood2)\naren.AddActor(hood3)\naren.AddActor(hood4)\naren.AddActor(hood6)\naren.AddActor(cookingPlate)\naren.AddActor(sgfilter)\naren.AddActor(lines)\naren.AddActor(rake)\n\naren.SetBackground(colors.GetColor3d(\"SlateGray\"))\n\naCamera = vtk.vtkCamera()\naren.SetActiveCamera(aCamera)\naren.ResetCamera()\n\naCamera.SetFocalPoint(3.505, 2.505, 1.255)\naCamera.SetPosition(3.505, 24.6196, 1.255)\naCamera.SetViewUp(0, 0, 1)\naCamera.Azimuth(60)\naCamera.Elevation(30)\naren.ResetCameraClippingRange()\n\nname = workdir + 'kitchen'\nvtk_render_gif(aren, N, name, Azimuth = True, Dolly = 1.5)\ndisp.stop()\n\nImage(filename=name + \".gif\", format='png')    ","execution_count":null,"outputs":[]},{"metadata":{"id":"Dz1Uc9n4GY6p"},"cell_type":"markdown","source":"## Frog example\nWe will leave this example for now without explanation - enjoy the show ;)\n\n* Source: https://lorensen.github.io/VTKExamples/site/Python/Visualization/ViewFrog/"},{"metadata":{"_kg_hide-input":true,"trusted":true,"id":"7dSWOztwGY6q"},"cell_type":"code","source":"! wget -q https://raw.githubusercontent.com/lorensen/VTKWikiExamples/master/Testing/Data/frogtissue.mhd -P /kaggle/working/vtk_examples \n! wget -q https://raw.githubusercontent.com/lorensen/VTKWikiExamples/master/Testing/Data/frogtissue.zraw  -P /kaggle/working/vtk_examples\ndef CreateFrogLut():\n    colors = vtk.vtkNamedColors()\n\n    colorLut = vtk.vtkLookupTable()\n    colorLut.SetNumberOfColors(17)\n    colorLut.SetTableRange(0, 16)\n    colorLut.Build()\n\n    colorLut.SetTableValue(0, 0, 0, 0, 0)\n    colorLut.SetTableValue(1, colors.GetColor4d(\"salmon\"))  # blood\n    colorLut.SetTableValue(2, colors.GetColor4d(\"beige\"))  # brain\n    colorLut.SetTableValue(3, colors.GetColor4d(\"orange\"))  # duodenum\n    colorLut.SetTableValue(4, colors.GetColor4d(\"misty_rose\"))  # eye_retina\n    colorLut.SetTableValue(5, colors.GetColor4d(\"white\"))  # eye_white\n    colorLut.SetTableValue(6, colors.GetColor4d(\"tomato\"))  # heart\n    colorLut.SetTableValue(7, colors.GetColor4d(\"raspberry\"))  # ileum\n    colorLut.SetTableValue(8, colors.GetColor4d(\"banana\"))  # kidney\n    colorLut.SetTableValue(9, colors.GetColor4d(\"peru\"))  # l_intestine\n    colorLut.SetTableValue(10, colors.GetColor4d(\"pink\"))  # liver\n    colorLut.SetTableValue(11, colors.GetColor4d(\"powder_blue\"))  # lung\n    colorLut.SetTableValue(12, colors.GetColor4d(\"carrot\"))  # nerve\n    colorLut.SetTableValue(13, colors.GetColor4d(\"wheat\"))  # skeleton\n    colorLut.SetTableValue(14, colors.GetColor4d(\"violet\"))  # spleen\n    colorLut.SetTableValue(15, colors.GetColor4d(\"plum\"))  # stomach\n\n    return colorLut\n\ndef CreateTissueMap():\n    tissueMap = dict()\n    tissueMap[\"blood\"] = 1\n    tissueMap[\"brain\"] = 2\n    tissueMap[\"duodenum\"] = 3\n    tissueMap[\"eyeRetina\"] = 4\n    tissueMap[\"eyeWhite\"] = 5\n    tissueMap[\"heart\"] = 6\n    tissueMap[\"ileum\"] = 7\n    tissueMap[\"kidney\"] = 8\n    tissueMap[\"intestine\"] = 9\n    tissueMap[\"liver\"] = 10\n    tissueMap[\"lung\"] = 11\n    tissueMap[\"nerve\"] = 12\n    tissueMap[\"skeleton\"] = 13\n    tissueMap[\"spleen\"] = 14\n    tissueMap[\"stomach\"] = 15\n\n    return tissueMap\n\ndef CreateFrogActor(fileName, tissue):\n    reader = vtk.vtkMetaImageReader()\n    reader.SetFileName(fileName)\n    reader.Update()\n\n    selectTissue = vtk.vtkImageThreshold()\n    selectTissue.ThresholdBetween(tissue, tissue)\n    selectTissue.SetInValue(255)\n    selectTissue.SetOutValue(0)\n    selectTissue.SetInputConnection(reader.GetOutputPort())\n\n    gaussianRadius = 1\n    gaussianStandardDeviation = 2.0\n    gaussian = vtk.vtkImageGaussianSmooth()\n    gaussian.SetStandardDeviations(gaussianStandardDeviation, gaussianStandardDeviation, gaussianStandardDeviation)\n    gaussian.SetRadiusFactors(gaussianRadius, gaussianRadius, gaussianRadius)\n    gaussian.SetInputConnection(selectTissue.GetOutputPort())\n\n    isoValue = 127.5\n    mcubes = vtk.vtkMarchingCubes()\n    mcubes.SetInputConnection(gaussian.GetOutputPort())\n    mcubes.ComputeScalarsOff()\n    mcubes.ComputeGradientsOff()\n    mcubes.ComputeNormalsOff()\n    mcubes.SetValue(0, isoValue)\n\n    smoothingIterations = 5\n    passBand = 0.001\n    featureAngle = 60.0\n    smoother = vtk.vtkWindowedSincPolyDataFilter()\n    smoother.SetInputConnection(mcubes.GetOutputPort())\n    smoother.SetNumberOfIterations(smoothingIterations)\n    smoother.BoundarySmoothingOff()\n    smoother.FeatureEdgeSmoothingOff()\n    smoother.SetFeatureAngle(featureAngle)\n    smoother.SetPassBand(passBand)\n    smoother.NonManifoldSmoothingOn()\n    smoother.NormalizeCoordinatesOn()\n    smoother.Update()\n\n    normals = vtk.vtkPolyDataNormals()\n    normals.SetInputConnection(smoother.GetOutputPort())\n    normals.SetFeatureAngle(featureAngle)\n\n    stripper = vtk.vtkStripper()\n    stripper.SetInputConnection(normals.GetOutputPort())\n\n    mapper = vtk.vtkPolyDataMapper()\n    mapper.SetInputConnection(stripper.GetOutputPort())\n\n    actor = vtk.vtkActor()\n    actor.SetMapper(mapper)\n\n    return actor","execution_count":null,"outputs":[]},{"metadata":{"_kg_hide-input":true,"trusted":true,"id":"FU6FYH3kGY6u"},"cell_type":"code","source":"%%time\ndisp = Display().start()\n\nfileName = workdir + 'frogtissue.mhd'\ntissueMap = CreateTissueMap()\n\ncolors = vtk.vtkNamedColors()\n\ncolorLut = CreateFrogLut()\nrenderer = vtk.vtkRenderer()\n\nfor tissue in [t for t in tissueMap]:\n    actor = CreateFrogActor(fileName, tissueMap[tissue])\n    actor.GetProperty().SetDiffuseColor( colorLut.GetTableValue(tissueMap[tissue])[:3])\n    actor.GetProperty().SetSpecular(.5)\n    actor.GetProperty().SetSpecularPower(10)\n    renderer.AddActor(actor)\n\nrenderer.GetActiveCamera().SetViewUp(0, 0, -1)\nrenderer.GetActiveCamera().SetPosition(0, -1, 0)\n\nrenderer.GetActiveCamera().Azimuth(210)\nrenderer.GetActiveCamera().Elevation(30)\nrenderer.ResetCamera()\n\nname = workdir + 'frog'\nvtk_render_gif(renderer, N, name, Azimuth = True, Dolly = 1.5)\ndisp.stop()\n\nImage(filename=name + \".gif\", format='png')    ","execution_count":null,"outputs":[]},{"metadata":{"id":"_V_-xGu0GY75"},"cell_type":"markdown","source":"# Appendix: Slicer3D in Kaggle notebooks\n* https://www.slicer.org/\n\n3D Slicer is an open source software platform for medical image informatics, image processing, and three-dimensional visualization. Built over two decades through support from the National Institutes of Health and a worldwide developer community, Slicer brings free, powerful cross-platform processing tools to physicians, researchers, and the general public. "},{"metadata":{"id":"xD3cuTVNGY76"},"cell_type":"markdown","source":"## Installation\n* https://www.slicer.org/wiki/Documentation/Nightly/ScriptRepository\n* https://discourse.slicer.org/t/running-slicer-without-gui/11720/4"},{"metadata":{"trusted":true,"_kg_hide-input":false,"_kg_hide-output":true,"id":"Te9DlHfSGY76"},"cell_type":"code","source":"! wget https://download.slicer.org/bitstream/1023242 -O Slicer-4.3.0-linux-amd64.tar.gz -q >>/dev/null\n! tar xzf Slicer-4.3.0-linux-amd64.tar.gz -C ~/ >>/dev/null\n! apt-get install libglu1 -qq >>/dev/null\n! apt-get install libpulse-mainloop-glib0 -qq >>/dev/null\n! apt-get install libegl-mesa0 -y libegl1 -qq >>/dev/null","execution_count":null,"outputs":[]},{"metadata":{"id":"Wv7deiH-GY79"},"cell_type":"markdown","source":"From a certain point of view, it may seem senseless to run Slicer from the command line without a real x-server. However, this was my starting point to use beautiful VTK visualizations. The code below shows the loading and visualization of several sample data sets. \n\n## MRHead Example\nGenerally, it comes down to saving the python code in a separate file and running Slicer with this script:\n\n* Source: https://www.slicer.org/wiki/Documentation/Nightly/Developers/Python_scripting#How_to_run_pip_.3F"},{"metadata":{"trusted":true,"id":"5OMbjzGBGY79"},"cell_type":"code","source":"import os\nfilepath = ''\ndef MakeFile(file_name):\n    temp_path = filepath + file_name\n    with open(file_name, 'w') as f:\n        f.write('''\\\n# use a slicer scripted module logic\nfrom SampleData import SampleDataLogic\nSampleDataLogic().downloadMRHead()\nhead = slicer.util.getNode(\"MRHead\")\n\n# use a vtk class\nthreshold = vtk.vtkImageThreshold()\nthreshold.SetInputData(head.GetImageData())\nthreshold.ThresholdBetween(100, 200)\nthreshold.SetInValue(255)\nthreshold.SetOutValue(0)\n\n#  use a slicer-specific C++ class\nerode = slicer.vtkImageErode()\nerode.SetInputConnection(threshold.GetOutputPort())\nerode.SetNeighborTo4()  \nerode.Update()          \n\nhead.SetAndObserveImageData(erode.GetOutputDataObject(0))\n\nslicer.util.saveNode(head, \"/kaggle/working/eroded.nrrd\")\n\nimport ScreenCapture\nl=ScreenCapture.ScreenCaptureLogic()\nl.captureImageFromView(l.viewFromNode(slicer.util.getNode('vtkMRMLSliceNodeRed')), '/kaggle/working/red.png')\nl.captureImageFromView(l.viewFromNode(slicer.util.getNode('vtkMRMLSliceNodeGreen')), '/kaggle/working/green.png')\nl.captureImageFromView(l.viewFromNode(slicer.util.getNode('vtkMRMLSliceNodeYellow')), '/kaggle/working/yellow.png')\nexit()      \n''')\nMakeFile('slicer_code.py')\n\n! xvfb-run -a ~/Slicer-4.10.2-linux-amd64/Slicer --no-splash --python-script slicer_code.py  > /dev/null","execution_count":null,"outputs":[]},{"metadata":{"_kg_hide-input":true,"trusted":true,"id":"2EPGXrCVGY8B"},"cell_type":"code","source":"plt.figure(figsize=(20,10))\nimgs = ['red','green','yellow']\nfor fi in range(len(imgs)):\n    plt.subplot(1,3,fi+1) \n    im = mpimg.imread(imgs[fi]+'.png')\n    plt.imshow(im) ","execution_count":null,"outputs":[]},{"metadata":{"id":"Z89BdkNhGY8E"},"cell_type":"markdown","source":"## CTChest Example"},{"metadata":{"trusted":true,"_kg_hide-input":true,"id":"Gvtl0vk6GY8F"},"cell_type":"code","source":"def MakeFile(file_name):\n    temp_path = filepath + file_name\n    with open(file_name, 'w') as f:\n        f.write('''\\\n# use a slicer scripted module logic\nfrom SampleData import SampleDataLogic\nSampleDataLogic().downloadCTChest()\nhead = slicer.util.getNode(\"CTChest\")\n\n# use a vtk class\nthreshold = vtk.vtkImageThreshold()\nthreshold.SetInputData(head.GetImageData())\nthreshold.ThresholdBetween(100, 800)\nthreshold.SetInValue(255)\nthreshold.SetOutValue(0)\n\n#  use a slicer-specific C++ class\nerode = slicer.vtkImageErode()\nerode.SetInputConnection(threshold.GetOutputPort())\nerode.SetNeighborTo4()  \nerode.Update()          \n\nhead.SetAndObserveImageData(erode.GetOutputDataObject(0))\n\nslicer.util.saveNode(head, \"/kaggle/working/eroded.nrrd\")\n\nimport ScreenCapture\nl=ScreenCapture.ScreenCaptureLogic()\nl.captureImageFromView(l.viewFromNode(slicer.util.getNode('vtkMRMLSliceNodeRed')), '/kaggle/working/red.png')\nl.captureImageFromView(l.viewFromNode(slicer.util.getNode('vtkMRMLSliceNodeGreen')), '/kaggle/working/green.png')\nl.captureImageFromView(l.viewFromNode(slicer.util.getNode('vtkMRMLSliceNodeYellow')), '/kaggle/working/yellow.png')\n\nrenderer = slicer.app.layoutManager().threeDWidget(0).threeDView()\n#centerViewport = [0.33, 0.0, .66, 1.0]\n#renderer.SetViewport(centerViewport)\n\n#width = 1000\n#height = 1000\nrenderWindow = renderer.renderWindow()\n#renderWindow.SetSize(width, height)\nrenderWindow.SetAlphaBitPlanes(1)\nwti = vtk.vtkWindowToImageFilter()\nwti.SetInputBufferTypeToRGBA()\nwti.SetInput(renderWindow)\nwriter = vtk.vtkPNGWriter()\nwriter.SetFileName(\"screenshot.png\")\nwriter.SetInputConnection(wti.GetOutputPort())\nwriter.Write()\n\nexit()      \n''')\nMakeFile('slicer_code.py')\n\n# https://discourse.slicer.org/t/running-slicer-without-gui/11720/4\n! xvfb-run -a ~/Slicer-4.10.2-linux-amd64/Slicer --no-splash --python-script slicer_code.py  > /dev/null","execution_count":null,"outputs":[]},{"metadata":{"_kg_hide-input":true,"trusted":true,"id":"OpNs5bo9GY8H"},"cell_type":"code","source":"plt.figure(figsize=(20,10))\n\nimgs = ['red','green','yellow']\nfor fi in range(len(imgs)):\n    plt.subplot(1,4,fi+1) \n    im = mpimg.imread(imgs[fi]+'.png')\n    plt.imshow(im) ","execution_count":null,"outputs":[]},{"metadata":{"id":"lfKdhOslGY8K"},"cell_type":"markdown","source":"## OSIC training data Example"},{"metadata":{"trusted":true,"id":"cm0puuUuGY8K"},"cell_type":"code","source":"plt.figure(figsize=(5,5))\ndir = '/kaggle/input/osic-pulmonary-fibrosis-progression/train/ID00009637202177434476278/'\nimport pydicom as dcm\nplt.imshow(dcm.dcmread(dir + os.listdir(dir)[0]).pixel_array)","execution_count":null,"outputs":[]},{"metadata":{"trusted":true,"_kg_hide-input":true,"id":"sSBOfCYnGY8N"},"cell_type":"code","source":"def MakeFile(file_name):\n    temp_path = filepath + file_name\n    with open(file_name, 'w') as f:\n        f.write('''\\\n# use a slicer scripted module logic\nfrom SampleData import SampleDataLogic\n#SampleDataLogic().downloadCTChest()\n#head = slicer.util.getNode(\"CTChest\")\n\nfrom DICOMLib import DICOMUtils\ndicomDataDir = '/kaggle/input/osic-pulmonary-fibrosis-progression/train/ID00009637202177434476278/' # input folder with DICOM files\nloadedNodeIDs = []  # this list will contain the list of all loaded node IDs\n\ndb=slicer.dicomDatabase\n\nwith DICOMUtils.TemporaryDICOMDatabase() as db:\n    DICOMUtils.importDicom(dicomDataDir, db)\n    patientUIDs = db.patients()\n    for patientUID in patientUIDs:\n        #loadedNodeIDs.extend(DICOMUtils.loadPatientByUID(patientUID))\n        loadedNodeIDs.append(DICOMUtils.loadPatientByUID(patientUID))\n\nhead = slicer.util.getNode()\n# use a vtk class\nthreshold = vtk.vtkImageThreshold()\n#threshold.SetInputData(head.GetImageData())\n#threshold.ThresholdBetween(0, 500)\n#threshold.SetInValue(255)\n#threshold.SetOutValue(0)\n\nimport ScreenCapture\nl=ScreenCapture.ScreenCaptureLogic()\nnode_red = slicer.util.getNode('vtkMRMLSliceNodeRed')\nl.captureImageFromView(l.viewFromNode(node_red), '/kaggle/working/red.png')\nl.captureImageFromView(l.viewFromNode(slicer.util.getNode('vtkMRMLSliceNodeGreen')), '/kaggle/working/green.png')\nl.captureImageFromView(l.viewFromNode(slicer.util.getNode('vtkMRMLSliceNodeYellow')), '/kaggle/working/yellow.png')\n\n\nslicer.util.resetSliceViews()\n\nlayoutManager = slicer.app.layoutManager()\nfor sliceViewName in layoutManager.sliceViewNames():\n  controller = layoutManager.sliceWidget(sliceViewName).sliceController()\n  controller.setSliceVisible(True)\n  \nthreeDWidget = layoutManager.threeDWidget(0)\nthreeDView = threeDWidget.threeDView()\nthreeDView.resetFocalPoint()\n\nrenderWindow = threeDView.renderWindow()\nrenderWindow.SetAlphaBitPlanes(1)\nwti = vtk.vtkWindowToImageFilter()\nwti.SetInputBufferTypeToRGBA()\nwti.SetInput(renderWindow)\nwriter = vtk.vtkPNGWriter()\nwriter.SetFileName(\"screenshot.png\")\nwriter.SetInputConnection(wti.GetOutputPort())\nwriter.Write()\n\n#layoutManager = slicer.app.layoutManager()\n#threeDWidget = layoutManager.threeDWidget(0)\n#threeDView = threeDWidget.threeDView()\n#threeDView.resetFocalPoint()\n\nexit()      \n''')\nMakeFile('slicer_code.py')\n\n# https://discourse.slicer.org/t/running-slicer-without-gui/11720/4\n! xvfb-run -a ~/Slicer-4.10.2-linux-amd64/Slicer --no-splash --python-script slicer_code.py > /dev/null","execution_count":null,"outputs":[]},{"metadata":{"_kg_hide-input":true,"trusted":true,"id":"7zpiPDSmGY8Q"},"cell_type":"code","source":"plt.figure(figsize=(20,10))\n\nimgs = ['red','green','yellow']\nfor fi in range(len(imgs)):\n    plt.subplot(1,4,fi+1) \n    im = mpimg.imread(imgs[fi]+'.png')\n    plt.imshow(im) \n! rm Slicer-4.3.0-linux-amd64.tar.gz\n! rm slicer_code.py","execution_count":null,"outputs":[]}],"metadata":{"kernelspec":{"language":"python","display_name":"Python 3","name":"python3"},"language_info":{"pygments_lexer":"ipython3","nbconvert_exporter":"python","version":"3.6.4","file_extension":".py","codemirror_mode":{"name":"ipython","version":3},"name":"python","mimetype":"text/x-python"}},"nbformat":4,"nbformat_minor":4}