{
"cells": [
{
"cell_type": "code",
"execution_count": 0,
"metadata": {
"colab": {
"autoexec": {
"startup": false,
"wait_interval": 0
}
},
"colab_type": "code",
"id": "mmKBfDbTa_4M"
},
"outputs": [],
"source": [
"import numpy as np\n",
"import scipy.ndimage "
]
},
{
"cell_type": "code",
"execution_count": 0,
"metadata": {
"colab": {
"autoexec": {
"startup": false,
"wait_interval": 0
}
},
"colab_type": "code",
"id": "-r22-U0ZnrRz"
},
"outputs": [],
"source": [
"# neighbour_code_to_normals is a lookup table.\n",
"# For every binary neighbour code \n",
"# (2x2x2 neighbourhood = 8 neighbours = 8 bits = 256 codes) \n",
"# it contains the surface normals of the triangles (called \"surfel\" for \n",
"# \"surface element\" in the following). The length of the normal \n",
"# vector encodes the surfel area.\n",
"#\n",
"# created by compute_surface_area_lookup_table.ipynb using the \n",
"# marching_cube algorithm, see e.g. https://en.wikipedia.org/wiki/Marching_cubes\n",
"#\n",
"neighbour_code_to_normals = [\n",
" [[0,0,0]],\n",
" [[0.125,0.125,0.125]],\n",
" [[-0.125,-0.125,0.125]],\n",
" [[-0.25,-0.25,0.0],[0.25,0.25,-0.0]],\n",
" [[0.125,-0.125,0.125]],\n",
" [[-0.25,-0.0,-0.25],[0.25,0.0,0.25]],\n",
" [[0.125,-0.125,0.125],[-0.125,-0.125,0.125]],\n",
" [[0.5,0.0,-0.0],[0.25,0.25,0.25],[0.125,0.125,0.125]],\n",
" [[-0.125,0.125,0.125]],\n",
" [[0.125,0.125,0.125],[-0.125,0.125,0.125]],\n",
" [[-0.25,0.0,0.25],[-0.25,0.0,0.25]],\n",
" [[0.5,0.0,0.0],[-0.25,-0.25,0.25],[-0.125,-0.125,0.125]],\n",
" [[0.25,-0.25,0.0],[0.25,-0.25,0.0]],\n",
" [[0.5,0.0,0.0],[0.25,-0.25,0.25],[-0.125,0.125,-0.125]],\n",
" [[-0.5,0.0,0.0],[-0.25,0.25,0.25],[-0.125,0.125,0.125]],\n",
" [[0.5,0.0,0.0],[0.5,0.0,0.0]],\n",
" [[0.125,-0.125,-0.125]],\n",
" [[0.0,-0.25,-0.25],[0.0,0.25,0.25]],\n",
" [[-0.125,-0.125,0.125],[0.125,-0.125,-0.125]],\n",
" [[0.0,-0.5,0.0],[0.25,0.25,0.25],[0.125,0.125,0.125]],\n",
" [[0.125,-0.125,0.125],[0.125,-0.125,-0.125]],\n",
" [[0.0,0.0,-0.5],[0.25,0.25,0.25],[-0.125,-0.125,-0.125]],\n",
" [[-0.125,-0.125,0.125],[0.125,-0.125,0.125],[0.125,-0.125,-0.125]],\n",
" [[-0.125,-0.125,-0.125],[-0.25,-0.25,-0.25],[0.25,0.25,0.25],[0.125,0.125,0.125]],\n",
" [[-0.125,0.125,0.125],[0.125,-0.125,-0.125]],\n",
" [[0.0,-0.25,-0.25],[0.0,0.25,0.25],[-0.125,0.125,0.125]],\n",
" [[-0.25,0.0,0.25],[-0.25,0.0,0.25],[0.125,-0.125,-0.125]],\n",
" [[0.125,0.125,0.125],[0.375,0.375,0.375],[0.0,-0.25,0.25],[-0.25,0.0,0.25]],\n",
" [[0.125,-0.125,-0.125],[0.25,-0.25,0.0],[0.25,-0.25,0.0]],\n",
" [[0.375,0.375,0.375],[0.0,0.25,-0.25],[-0.125,-0.125,-0.125],[-0.25,0.25,0.0]],\n",
" [[-0.5,0.0,0.0],[-0.125,-0.125,-0.125],[-0.25,-0.25,-0.25],[0.125,0.125,0.125]],\n",
" [[-0.5,0.0,0.0],[-0.125,-0.125,-0.125],[-0.25,-0.25,-0.25]],\n",
" [[0.125,-0.125,0.125]],\n",
" [[0.125,0.125,0.125],[0.125,-0.125,0.125]],\n",
" [[0.0,-0.25,0.25],[0.0,0.25,-0.25]],\n",
" [[0.0,-0.5,0.0],[0.125,0.125,-0.125],[0.25,0.25,-0.25]],\n",
" [[0.125,-0.125,0.125],[0.125,-0.125,0.125]],\n",
" [[0.125,-0.125,0.125],[-0.25,-0.0,-0.25],[0.25,0.0,0.25]],\n",
" [[0.0,-0.25,0.25],[0.0,0.25,-0.25],[0.125,-0.125,0.125]],\n",
" [[-0.375,-0.375,0.375],[-0.0,0.25,0.25],[0.125,0.125,-0.125],[-0.25,-0.0,-0.25]],\n",
" [[-0.125,0.125,0.125],[0.125,-0.125,0.125]],\n",
" [[0.125,0.125,0.125],[0.125,-0.125,0.125],[-0.125,0.125,0.125]],\n",
" [[-0.0,0.0,0.5],[-0.25,-0.25,0.25],[-0.125,-0.125,0.125]],\n",
" [[0.25,0.25,-0.25],[0.25,0.25,-0.25],[0.125,0.125,-0.125],[-0.125,-0.125,0.125]],\n",
" [[0.125,-0.125,0.125],[0.25,-0.25,0.0],[0.25,-0.25,0.0]],\n",
" [[0.5,0.0,0.0],[0.25,-0.25,0.25],[-0.125,0.125,-0.125],[0.125,-0.125,0.125]],\n",
" [[0.0,0.25,-0.25],[0.375,-0.375,-0.375],[-0.125,0.125,0.125],[0.25,0.25,0.0]],\n",
" [[-0.5,0.0,0.0],[-0.25,-0.25,0.25],[-0.125,-0.125,0.125]],\n",
" [[0.25,-0.25,0.0],[-0.25,0.25,0.0]],\n",
" [[0.0,0.5,0.0],[-0.25,0.25,0.25],[0.125,-0.125,-0.125]],\n",
" [[0.0,0.5,0.0],[0.125,-0.125,0.125],[-0.25,0.25,-0.25]],\n",
" [[0.0,0.5,0.0],[0.0,-0.5,0.0]],\n",
" [[0.25,-0.25,0.0],[-0.25,0.25,0.0],[0.125,-0.125,0.125]],\n",
" [[-0.375,-0.375,-0.375],[-0.25,0.0,0.25],[-0.125,-0.125,-0.125],[-0.25,0.25,0.0]],\n",
" [[0.125,0.125,0.125],[0.0,-0.5,0.0],[-0.25,-0.25,-0.25],[-0.125,-0.125,-0.125]],\n",
" [[0.0,-0.5,0.0],[-0.25,-0.25,-0.25],[-0.125,-0.125,-0.125]],\n",
" [[-0.125,0.125,0.125],[0.25,-0.25,0.0],[-0.25,0.25,0.0]],\n",
" [[0.0,0.5,0.0],[0.25,0.25,-0.25],[-0.125,-0.125,0.125],[-0.125,-0.125,0.125]],\n",
" [[-0.375,0.375,-0.375],[-0.25,-0.25,0.0],[-0.125,0.125,-0.125],[-0.25,0.0,0.25]],\n",
" [[0.0,0.5,0.0],[0.25,0.25,-0.25],[-0.125,-0.125,0.125]],\n",
" [[0.25,-0.25,0.0],[-0.25,0.25,0.0],[0.25,-0.25,0.0],[0.25,-0.25,0.0]],\n",
" [[-0.25,-0.25,0.0],[-0.25,-0.25,0.0],[-0.125,-0.125,0.125]],\n",
" [[0.125,0.125,0.125],[-0.25,-0.25,0.0],[-0.25,-0.25,0.0]],\n",
" [[-0.25,-0.25,0.0],[-0.25,-0.25,0.0]],\n",
" [[-0.125,-0.125,0.125]],\n",
" [[0.125,0.125,0.125],[-0.125,-0.125,0.125]],\n",
" [[-0.125,-0.125,0.125],[-0.125,-0.125,0.125]],\n",
" [[-0.125,-0.125,0.125],[-0.25,-0.25,0.0],[0.25,0.25,-0.0]],\n",
" [[0.0,-0.25,0.25],[0.0,-0.25,0.25]],\n",
" [[0.0,0.0,0.5],[0.25,-0.25,0.25],[0.125,-0.125,0.125]],\n",
" [[0.0,-0.25,0.25],[0.0,-0.25,0.25],[-0.125,-0.125,0.125]],\n",
" [[0.375,-0.375,0.375],[0.0,-0.25,-0.25],[-0.125,0.125,-0.125],[0.25,0.25,0.0]],\n",
" [[-0.125,-0.125,0.125],[-0.125,0.125,0.125]],\n",
" [[0.125,0.125,0.125],[-0.125,-0.125,0.125],[-0.125,0.125,0.125]],\n",
" [[-0.125,-0.125,0.125],[-0.25,0.0,0.25],[-0.25,0.0,0.25]],\n",
" [[0.5,0.0,0.0],[-0.25,-0.25,0.25],[-0.125,-0.125,0.125],[-0.125,-0.125,0.125]],\n",
" [[-0.0,0.5,0.0],[-0.25,0.25,-0.25],[0.125,-0.125,0.125]],\n",
" [[-0.25,0.25,-0.25],[-0.25,0.25,-0.25],[-0.125,0.125,-0.125],[-0.125,0.125,-0.125]],\n",
" [[-0.25,0.0,-0.25],[0.375,-0.375,-0.375],[0.0,0.25,-0.25],[-0.125,0.125,0.125]],\n",
" [[0.5,0.0,0.0],[-0.25,0.25,-0.25],[0.125,-0.125,0.125]],\n",
" [[-0.25,0.0,0.25],[0.25,0.0,-0.25]],\n",
" [[-0.0,0.0,0.5],[-0.25,0.25,0.25],[-0.125,0.125,0.125]],\n",
" [[-0.125,-0.125,0.125],[-0.25,0.0,0.25],[0.25,0.0,-0.25]],\n",
" [[-0.25,-0.0,-0.25],[-0.375,0.375,0.375],[-0.25,-0.25,0.0],[-0.125,0.125,0.125]],\n",
" [[0.0,0.0,-0.5],[0.25,0.25,-0.25],[-0.125,-0.125,0.125]],\n",
" [[-0.0,0.0,0.5],[0.0,0.0,0.5]],\n",
" [[0.125,0.125,0.125],[0.125,0.125,0.125],[0.25,0.25,0.25],[0.0,0.0,0.5]],\n",
" [[0.125,0.125,0.125],[0.25,0.25,0.25],[0.0,0.0,0.5]],\n",
" [[-0.25,0.0,0.25],[0.25,0.0,-0.25],[-0.125,0.125,0.125]],\n",
" [[-0.0,0.0,0.5],[0.25,-0.25,0.25],[0.125,-0.125,0.125],[0.125,-0.125,0.125]],\n",
" [[-0.25,0.0,0.25],[-0.25,0.0,0.25],[-0.25,0.0,0.25],[0.25,0.0,-0.25]],\n",
" [[0.125,-0.125,0.125],[0.25,0.0,0.25],[0.25,0.0,0.25]],\n",
" [[0.25,0.0,0.25],[-0.375,-0.375,0.375],[-0.25,0.25,0.0],[-0.125,-0.125,0.125]],\n",
" [[-0.0,0.0,0.5],[0.25,-0.25,0.25],[0.125,-0.125,0.125]],\n",
" [[0.125,0.125,0.125],[0.25,0.0,0.25],[0.25,0.0,0.25]],\n",
" [[0.25,0.0,0.25],[0.25,0.0,0.25]],\n",
" [[-0.125,-0.125,0.125],[0.125,-0.125,0.125]],\n",
" [[0.125,0.125,0.125],[-0.125,-0.125,0.125],[0.125,-0.125,0.125]],\n",
" [[-0.125,-0.125,0.125],[0.0,-0.25,0.25],[0.0,0.25,-0.25]],\n",
" [[0.0,-0.5,0.0],[0.125,0.125,-0.125],[0.25,0.25,-0.25],[-0.125,-0.125,0.125]],\n",
" [[0.0,-0.25,0.25],[0.0,-0.25,0.25],[0.125,-0.125,0.125]],\n",
" [[0.0,0.0,0.5],[0.25,-0.25,0.25],[0.125,-0.125,0.125],[0.125,-0.125,0.125]],\n",
" [[0.0,-0.25,0.25],[0.0,-0.25,0.25],[0.0,-0.25,0.25],[0.0,0.25,-0.25]],\n",
" [[0.0,0.25,0.25],[0.0,0.25,0.25],[0.125,-0.125,-0.125]],\n",
" [[-0.125,0.125,0.125],[0.125,-0.125,0.125],[-0.125,-0.125,0.125]],\n",
" [[-0.125,0.125,0.125],[0.125,-0.125,0.125],[-0.125,-0.125,0.125],[0.125,0.125,0.125]],\n",
" [[-0.0,0.0,0.5],[-0.25,-0.25,0.25],[-0.125,-0.125,0.125],[-0.125,-0.125,0.125]],\n",
" [[0.125,0.125,0.125],[0.125,-0.125,0.125],[0.125,-0.125,-0.125]],\n",
" [[-0.0,0.5,0.0],[-0.25,0.25,-0.25],[0.125,-0.125,0.125],[0.125,-0.125,0.125]],\n",
" [[0.125,0.125,0.125],[-0.125,-0.125,0.125],[0.125,-0.125,-0.125]],\n",
" [[0.0,-0.25,-0.25],[0.0,0.25,0.25],[0.125,0.125,0.125]],\n",
" [[0.125,0.125,0.125],[0.125,-0.125,-0.125]],\n",
" [[0.5,0.0,-0.0],[0.25,-0.25,-0.25],[0.125,-0.125,-0.125]],\n",
" [[-0.25,0.25,0.25],[-0.125,0.125,0.125],[-0.25,0.25,0.25],[0.125,-0.125,-0.125]],\n",
" [[0.375,-0.375,0.375],[0.0,0.25,0.25],[-0.125,0.125,-0.125],[-0.25,0.0,0.25]],\n",
" [[0.0,-0.5,0.0],[-0.25,0.25,0.25],[-0.125,0.125,0.125]],\n",
" [[-0.375,-0.375,0.375],[0.25,-0.25,0.0],[0.0,0.25,0.25],[-0.125,-0.125,0.125]],\n",
" [[-0.125,0.125,0.125],[-0.25,0.25,0.25],[0.0,0.0,0.5]],\n",
" [[0.125,0.125,0.125],[0.0,0.25,0.25],[0.0,0.25,0.25]],\n",
" [[0.0,0.25,0.25],[0.0,0.25,0.25]],\n",
" [[0.5,0.0,-0.0],[0.25,0.25,0.25],[0.125,0.125,0.125],[0.125,0.125,0.125]],\n",
" [[0.125,-0.125,0.125],[-0.125,-0.125,0.125],[0.125,0.125,0.125]],\n",
" [[-0.25,-0.0,-0.25],[0.25,0.0,0.25],[0.125,0.125,0.125]],\n",
" [[0.125,0.125,0.125],[0.125,-0.125,0.125]],\n",
" [[-0.25,-0.25,0.0],[0.25,0.25,-0.0],[0.125,0.125,0.125]],\n",
" [[0.125,0.125,0.125],[-0.125,-0.125,0.125]],\n",
" [[0.125,0.125,0.125],[0.125,0.125,0.125]],\n",
" [[0.125,0.125,0.125]],\n",
" [[0.125,0.125,0.125]],\n",
" [[0.125,0.125,0.125],[0.125,0.125,0.125]],\n",
" [[0.125,0.125,0.125],[-0.125,-0.125,0.125]],\n",
" [[-0.25,-0.25,0.0],[0.25,0.25,-0.0],[0.125,0.125,0.125]],\n",
" [[0.125,0.125,0.125],[0.125,-0.125,0.125]],\n",
" [[-0.25,-0.0,-0.25],[0.25,0.0,0.25],[0.125,0.125,0.125]],\n",
" [[0.125,-0.125,0.125],[-0.125,-0.125,0.125],[0.125,0.125,0.125]],\n",
" [[0.5,0.0,-0.0],[0.25,0.25,0.25],[0.125,0.125,0.125],[0.125,0.125,0.125]],\n",
" [[0.0,0.25,0.25],[0.0,0.25,0.25]],\n",
" [[0.125,0.125,0.125],[0.0,0.25,0.25],[0.0,0.25,0.25]],\n",
" [[-0.125,0.125,0.125],[-0.25,0.25,0.25],[0.0,0.0,0.5]],\n",
" [[-0.375,-0.375,0.375],[0.25,-0.25,0.0],[0.0,0.25,0.25],[-0.125,-0.125,0.125]],\n",
" [[0.0,-0.5,0.0],[-0.25,0.25,0.25],[-0.125,0.125,0.125]],\n",
" [[0.375,-0.375,0.375],[0.0,0.25,0.25],[-0.125,0.125,-0.125],[-0.25,0.0,0.25]],\n",
" [[-0.25,0.25,0.25],[-0.125,0.125,0.125],[-0.25,0.25,0.25],[0.125,-0.125,-0.125]],\n",
" [[0.5,0.0,-0.0],[0.25,-0.25,-0.25],[0.125,-0.125,-0.125]],\n",
" [[0.125,0.125,0.125],[0.125,-0.125,-0.125]],\n",
" [[0.0,-0.25,-0.25],[0.0,0.25,0.25],[0.125,0.125,0.125]],\n",
" [[0.125,0.125,0.125],[-0.125,-0.125,0.125],[0.125,-0.125,-0.125]],\n",
" [[-0.0,0.5,0.0],[-0.25,0.25,-0.25],[0.125,-0.125,0.125],[0.125,-0.125,0.125]],\n",
" [[0.125,0.125,0.125],[0.125,-0.125,0.125],[0.125,-0.125,-0.125]],\n",
" [[-0.0,0.0,0.5],[-0.25,-0.25,0.25],[-0.125,-0.125,0.125],[-0.125,-0.125,0.125]],\n",
" [[-0.125,0.125,0.125],[0.125,-0.125,0.125],[-0.125,-0.125,0.125],[0.125,0.125,0.125]],\n",
" [[-0.125,0.125,0.125],[0.125,-0.125,0.125],[-0.125,-0.125,0.125]],\n",
" [[0.0,0.25,0.25],[0.0,0.25,0.25],[0.125,-0.125,-0.125]],\n",
" [[0.0,-0.25,-0.25],[0.0,0.25,0.25],[0.0,0.25,0.25],[0.0,0.25,0.25]],\n",
" [[0.0,0.0,0.5],[0.25,-0.25,0.25],[0.125,-0.125,0.125],[0.125,-0.125,0.125]],\n",
" [[0.0,-0.25,0.25],[0.0,-0.25,0.25],[0.125,-0.125,0.125]],\n",
" [[0.0,-0.5,0.0],[0.125,0.125,-0.125],[0.25,0.25,-0.25],[-0.125,-0.125,0.125]],\n",
" [[-0.125,-0.125,0.125],[0.0,-0.25,0.25],[0.0,0.25,-0.25]],\n",
" [[0.125,0.125,0.125],[-0.125,-0.125,0.125],[0.125,-0.125,0.125]],\n",
" [[-0.125,-0.125,0.125],[0.125,-0.125,0.125]],\n",
" [[0.25,0.0,0.25],[0.25,0.0,0.25]],\n",
" [[0.125,0.125,0.125],[0.25,0.0,0.25],[0.25,0.0,0.25]],\n",
" [[-0.0,0.0,0.5],[0.25,-0.25,0.25],[0.125,-0.125,0.125]],\n",
" [[0.25,0.0,0.25],[-0.375,-0.375,0.375],[-0.25,0.25,0.0],[-0.125,-0.125,0.125]],\n",
" [[0.125,-0.125,0.125],[0.25,0.0,0.25],[0.25,0.0,0.25]],\n",
" [[-0.25,-0.0,-0.25],[0.25,0.0,0.25],[0.25,0.0,0.25],[0.25,0.0,0.25]],\n",
" [[-0.0,0.0,0.5],[0.25,-0.25,0.25],[0.125,-0.125,0.125],[0.125,-0.125,0.125]],\n",
" [[-0.25,0.0,0.25],[0.25,0.0,-0.25],[-0.125,0.125,0.125]],\n",
" [[0.125,0.125,0.125],[0.25,0.25,0.25],[0.0,0.0,0.5]],\n",
" [[0.125,0.125,0.125],[0.125,0.125,0.125],[0.25,0.25,0.25],[0.0,0.0,0.5]],\n",
" [[-0.0,0.0,0.5],[0.0,0.0,0.5]],\n",
" [[0.0,0.0,-0.5],[0.25,0.25,-0.25],[-0.125,-0.125,0.125]],\n",
" [[-0.25,-0.0,-0.25],[-0.375,0.375,0.375],[-0.25,-0.25,0.0],[-0.125,0.125,0.125]],\n",
" [[-0.125,-0.125,0.125],[-0.25,0.0,0.25],[0.25,0.0,-0.25]],\n",
" [[-0.0,0.0,0.5],[-0.25,0.25,0.25],[-0.125,0.125,0.125]],\n",
" [[-0.25,0.0,0.25],[0.25,0.0,-0.25]],\n",
" [[0.5,0.0,0.0],[-0.25,0.25,-0.25],[0.125,-0.125,0.125]],\n",
" [[-0.25,0.0,-0.25],[0.375,-0.375,-0.375],[0.0,0.25,-0.25],[-0.125,0.125,0.125]],\n",
" [[-0.25,0.25,-0.25],[-0.25,0.25,-0.25],[-0.125,0.125,-0.125],[-0.125,0.125,-0.125]],\n",
" [[-0.0,0.5,0.0],[-0.25,0.25,-0.25],[0.125,-0.125,0.125]],\n",
" [[0.5,0.0,0.0],[-0.25,-0.25,0.25],[-0.125,-0.125,0.125],[-0.125,-0.125,0.125]],\n",
" [[-0.125,-0.125,0.125],[-0.25,0.0,0.25],[-0.25,0.0,0.25]],\n",
" [[0.125,0.125,0.125],[-0.125,-0.125,0.125],[-0.125,0.125,0.125]],\n",
" [[-0.125,-0.125,0.125],[-0.125,0.125,0.125]],\n",
" [[0.375,-0.375,0.375],[0.0,-0.25,-0.25],[-0.125,0.125,-0.125],[0.25,0.25,0.0]],\n",
" [[0.0,-0.25,0.25],[0.0,-0.25,0.25],[-0.125,-0.125,0.125]],\n",
" [[0.0,0.0,0.5],[0.25,-0.25,0.25],[0.125,-0.125,0.125]],\n",
" [[0.0,-0.25,0.25],[0.0,-0.25,0.25]],\n",
" [[-0.125,-0.125,0.125],[-0.25,-0.25,0.0],[0.25,0.25,-0.0]],\n",
" [[-0.125,-0.125,0.125],[-0.125,-0.125,0.125]],\n",
" [[0.125,0.125,0.125],[-0.125,-0.125,0.125]],\n",
" [[-0.125,-0.125,0.125]],\n",
" [[-0.25,-0.25,0.0],[-0.25,-0.25,0.0]],\n",
" [[0.125,0.125,0.125],[-0.25,-0.25,0.0],[-0.25,-0.25,0.0]],\n",
" [[-0.25,-0.25,0.0],[-0.25,-0.25,0.0],[-0.125,-0.125,0.125]],\n",
" [[-0.25,-0.25,0.0],[-0.25,-0.25,0.0],[-0.25,-0.25,0.0],[0.25,0.25,-0.0]],\n",
" [[0.0,0.5,0.0],[0.25,0.25,-0.25],[-0.125,-0.125,0.125]],\n",
" [[-0.375,0.375,-0.375],[-0.25,-0.25,0.0],[-0.125,0.125,-0.125],[-0.25,0.0,0.25]],\n",
" [[0.0,0.5,0.0],[0.25,0.25,-0.25],[-0.125,-0.125,0.125],[-0.125,-0.125,0.125]],\n",
" [[-0.125,0.125,0.125],[0.25,-0.25,0.0],[-0.25,0.25,0.0]],\n",
" [[0.0,-0.5,0.0],[-0.25,-0.25,-0.25],[-0.125,-0.125,-0.125]],\n",
" [[0.125,0.125,0.125],[0.0,-0.5,0.0],[-0.25,-0.25,-0.25],[-0.125,-0.125,-0.125]],\n",
" [[-0.375,-0.375,-0.375],[-0.25,0.0,0.25],[-0.125,-0.125,-0.125],[-0.25,0.25,0.0]],\n",
" [[0.25,-0.25,0.0],[-0.25,0.25,0.0],[0.125,-0.125,0.125]],\n",
" [[0.0,0.5,0.0],[0.0,-0.5,0.0]],\n",
" [[0.0,0.5,0.0],[0.125,-0.125,0.125],[-0.25,0.25,-0.25]],\n",
" [[0.0,0.5,0.0],[-0.25,0.25,0.25],[0.125,-0.125,-0.125]],\n",
" [[0.25,-0.25,0.0],[-0.25,0.25,0.0]],\n",
" [[-0.5,0.0,0.0],[-0.25,-0.25,0.25],[-0.125,-0.125,0.125]],\n",
" [[0.0,0.25,-0.25],[0.375,-0.375,-0.375],[-0.125,0.125,0.125],[0.25,0.25,0.0]],\n",
" [[0.5,0.0,0.0],[0.25,-0.25,0.25],[-0.125,0.125,-0.125],[0.125,-0.125,0.125]],\n",
" [[0.125,-0.125,0.125],[0.25,-0.25,0.0],[0.25,-0.25,0.0]],\n",
" [[0.25,0.25,-0.25],[0.25,0.25,-0.25],[0.125,0.125,-0.125],[-0.125,-0.125,0.125]],\n",
" [[-0.0,0.0,0.5],[-0.25,-0.25,0.25],[-0.125,-0.125,0.125]],\n",
" [[0.125,0.125,0.125],[0.125,-0.125,0.125],[-0.125,0.125,0.125]],\n",
" [[-0.125,0.125,0.125],[0.125,-0.125,0.125]],\n",
" [[-0.375,-0.375,0.375],[-0.0,0.25,0.25],[0.125,0.125,-0.125],[-0.25,-0.0,-0.25]],\n",
" [[0.0,-0.25,0.25],[0.0,0.25,-0.25],[0.125,-0.125,0.125]],\n",
" [[0.125,-0.125,0.125],[-0.25,-0.0,-0.25],[0.25,0.0,0.25]],\n",
" [[0.125,-0.125,0.125],[0.125,-0.125,0.125]],\n",
" [[0.0,-0.5,0.0],[0.125,0.125,-0.125],[0.25,0.25,-0.25]],\n",
" [[0.0,-0.25,0.25],[0.0,0.25,-0.25]],\n",
" [[0.125,0.125,0.125],[0.125,-0.125,0.125]],\n",
" [[0.125,-0.125,0.125]],\n",
" [[-0.5,0.0,0.0],[-0.125,-0.125,-0.125],[-0.25,-0.25,-0.25]],\n",
" [[-0.5,0.0,0.0],[-0.125,-0.125,-0.125],[-0.25,-0.25,-0.25],[0.125,0.125,0.125]],\n",
" [[0.375,0.375,0.375],[0.0,0.25,-0.25],[-0.125,-0.125,-0.125],[-0.25,0.25,0.0]],\n",
" [[0.125,-0.125,-0.125],[0.25,-0.25,0.0],[0.25,-0.25,0.0]],\n",
" [[0.125,0.125,0.125],[0.375,0.375,0.375],[0.0,-0.25,0.25],[-0.25,0.0,0.25]],\n",
" [[-0.25,0.0,0.25],[-0.25,0.0,0.25],[0.125,-0.125,-0.125]],\n",
" [[0.0,-0.25,-0.25],[0.0,0.25,0.25],[-0.125,0.125,0.125]],\n",
" [[-0.125,0.125,0.125],[0.125,-0.125,-0.125]],\n",
" [[-0.125,-0.125,-0.125],[-0.25,-0.25,-0.25],[0.25,0.25,0.25],[0.125,0.125,0.125]],\n",
" [[-0.125,-0.125,0.125],[0.125,-0.125,0.125],[0.125,-0.125,-0.125]],\n",
" [[0.0,0.0,-0.5],[0.25,0.25,0.25],[-0.125,-0.125,-0.125]],\n",
" [[0.125,-0.125,0.125],[0.125,-0.125,-0.125]],\n",
" [[0.0,-0.5,0.0],[0.25,0.25,0.25],[0.125,0.125,0.125]],\n",
" [[-0.125,-0.125,0.125],[0.125,-0.125,-0.125]],\n",
" [[0.0,-0.25,-0.25],[0.0,0.25,0.25]],\n",
" [[0.125,-0.125,-0.125]],\n",
" [[0.5,0.0,0.0],[0.5,0.0,0.0]],\n",
" [[-0.5,0.0,0.0],[-0.25,0.25,0.25],[-0.125,0.125,0.125]],\n",
" [[0.5,0.0,0.0],[0.25,-0.25,0.25],[-0.125,0.125,-0.125]],\n",
" [[0.25,-0.25,0.0],[0.25,-0.25,0.0]],\n",
" [[0.5,0.0,0.0],[-0.25,-0.25,0.25],[-0.125,-0.125,0.125]],\n",
" [[-0.25,0.0,0.25],[-0.25,0.0,0.25]],\n",
" [[0.125,0.125,0.125],[-0.125,0.125,0.125]],\n",
" [[-0.125,0.125,0.125]],\n",
" [[0.5,0.0,-0.0],[0.25,0.25,0.25],[0.125,0.125,0.125]],\n",
" [[0.125,-0.125,0.125],[-0.125,-0.125,0.125]],\n",
" [[-0.25,-0.0,-0.25],[0.25,0.0,0.25]],\n",
" [[0.125,-0.125,0.125]],\n",
" [[-0.25,-0.25,0.0],[0.25,0.25,-0.0]],\n",
" [[-0.125,-0.125,0.125]],\n",
" [[0.125,0.125,0.125]],\n",
" [[0,0,0]]]"
]
},
{
"cell_type": "code",
"execution_count": 0,
"metadata": {
"colab": {
"autoexec": {
"startup": false,
"wait_interval": 0
}
},
"colab_type": "code",
"id": "-kruspp5Y1Ip"
},
"outputs": [],
"source": [
"def compute_surface_distances(mask_gt, mask_pred, spacing_mm):\n",
" \"\"\"Compute closest distances from all surface points to the other surface.\n",
"\n",
" Finds all surface elements \"surfels\" in the ground truth mask `mask_gt` and\n",
" the predicted mask `mask_pred`, computes their area in mm^2 and the distance\n",
" to the closest point on the other surface. It returns two sorted lists of\n",
" distances together with the corresponding surfel areas. If one of the masks\n",
" is empty, the corresponding lists are empty and all distances in the other\n",
" list are `inf` \n",
" \n",
" Args:\n",
" mask_gt: 3-dim Numpy array of type bool. The ground truth mask.\n",
" mask_pred: 3-dim Numpy array of type bool. The predicted mask.\n",
" spacing_mm: 3-element list-like structure. Voxel spacing in x0, x1 and x2\n",
" direction \n",
"\n",
" Returns:\n",
" A dict with \n",
" \"distances_gt_to_pred\": 1-dim numpy array of type float. The distances in mm\n",
" from all ground truth surface elements to the predicted surface, \n",
" sorted from smallest to largest\n",
" \"distances_pred_to_gt\": 1-dim numpy array of type float. The distances in mm\n",
" from all predicted surface elements to the ground truth surface, \n",
" sorted from smallest to largest \n",
" \"surfel_areas_gt\": 1-dim numpy array of type float. The area in mm^2 of \n",
" the ground truth surface elements in the same order as \n",
" distances_gt_to_pred\n",
" \"surfel_areas_pred\": 1-dim numpy array of type float. The area in mm^2 of \n",
" the predicted surface elements in the same order as \n",
" distances_pred_to_gt\n",
" \n",
" \"\"\"\n",
" \n",
" # compute the area for all 256 possible surface elements \n",
" # (given a 2x2x2 neighbourhood) according to the spacing_mm\n",
" neighbour_code_to_surface_area = np.zeros([256])\n",
" for code in range(256):\n",
" normals = np.array(neighbour_code_to_normals[code])\n",
" sum_area = 0\n",
" for normal_idx in range(normals.shape[0]):\n",
" # normal vector\n",
" n = np.zeros([3])\n",
" n[0] = normals[normal_idx,0] * spacing_mm[1] * spacing_mm[2]\n",
" n[1] = normals[normal_idx,1] * spacing_mm[0] * spacing_mm[2]\n",
" n[2] = normals[normal_idx,2] * spacing_mm[0] * spacing_mm[1]\n",
" area = np.linalg.norm(n)\n",
" sum_area += area\n",
" neighbour_code_to_surface_area[code] = sum_area\n",
"\n",
" # compute the bounding box of the masks to trim\n",
" # the volume to the smallest possible processing subvolume\n",
" mask_all = mask_gt | mask_pred\n",
" bbox_min = np.zeros(3, np.int64)\n",
" bbox_max = np.zeros(3, np.int64)\n",
"\n",
" # max projection to the x0-axis\n",
" proj_0 = np.max(np.max(mask_all, axis=2), axis=1)\n",
" idx_nonzero_0 = np.nonzero(proj_0)[0]\n",
" if len(idx_nonzero_0) == 0:\n",
" return {\"distances_gt_to_pred\": np.array([]), \n",
" \"distances_pred_to_gt\": np.array([]), \n",
" \"surfel_areas_gt\": np.array([]), \n",
" \"surfel_areas_pred\": np.array([])}\n",
" \n",
" bbox_min[0] = np.min(idx_nonzero_0)\n",
" bbox_max[0] = np.max(idx_nonzero_0)\n",
"\n",
" # max projection to the x1-axis\n",
" proj_1 = np.max(np.max(mask_all, axis=2), axis=0)\n",
" idx_nonzero_1 = np.nonzero(proj_1)[0]\n",
" bbox_min[1] = np.min(idx_nonzero_1)\n",
" bbox_max[1] = np.max(idx_nonzero_1)\n",
"\n",
" # max projection to the x2-axis\n",
" proj_2 = np.max(np.max(mask_all, axis=1), axis=0)\n",
" idx_nonzero_2 = np.nonzero(proj_2)[0]\n",
" bbox_min[2] = np.min(idx_nonzero_2)\n",
" bbox_max[2] = np.max(idx_nonzero_2)\n",
"\n",
" print(\"bounding box min = {}\".format(bbox_min))\n",
" print(\"bounding box max = {}\".format(bbox_max))\n",
"\n",
" # crop the processing subvolume.\n",
" # we need to zeropad the cropped region with 1 voxel at the lower, \n",
" # the right and the back side. This is required to obtain the \"full\" \n",
" # convolution result with the 2x2x2 kernel\n",
" cropmask_gt = np.zeros((bbox_max - bbox_min)+2, np.uint8)\n",
" cropmask_pred = np.zeros((bbox_max - bbox_min)+2, np.uint8)\n",
"\n",
" cropmask_gt[0:-1, 0:-1, 0:-1] = mask_gt[bbox_min[0]:bbox_max[0]+1,\n",
" bbox_min[1]:bbox_max[1]+1,\n",
" bbox_min[2]:bbox_max[2]+1]\n",
"\n",
" cropmask_pred[0:-1, 0:-1, 0:-1] = mask_pred[bbox_min[0]:bbox_max[0]+1,\n",
" bbox_min[1]:bbox_max[1]+1,\n",
" bbox_min[2]:bbox_max[2]+1]\n",
"\n",
" # compute the neighbour code (local binary pattern) for each voxel\n",
" # the resultsing arrays are spacially shifted by minus half a voxel in each axis.\n",
" # i.e. the points are located at the corners of the original voxels\n",
" kernel = np.array([[[128,64],\n",
" [32,16]],\n",
" [[8,4],\n",
" [2,1]]])\n",
" neighbour_code_map_gt = scipy.ndimage.filters.correlate(cropmask_gt.astype(np.uint8), kernel, mode=\"constant\", cval=0) \n",
" neighbour_code_map_pred = scipy.ndimage.filters.correlate(cropmask_pred.astype(np.uint8), kernel, mode=\"constant\", cval=0) \n",
"\n",
" # create masks with the surface voxels\n",
" borders_gt = ((neighbour_code_map_gt != 0) & (neighbour_code_map_gt != 255))\n",
" borders_pred = ((neighbour_code_map_pred != 0) & (neighbour_code_map_pred != 255))\n",
"\n",
" # compute the distance transform (closest distance of each voxel to the surface voxels)\n",
" if borders_gt.any():\n",
" distmap_gt = scipy.ndimage.morphology.distance_transform_edt(~borders_gt, sampling=spacing_mm)\n",
" else:\n",
" distmap_gt = np.Inf * np.ones(borders_gt.shape)\n",
"\n",
" if borders_pred.any(): \n",
" distmap_pred = scipy.ndimage.morphology.distance_transform_edt(~borders_pred, sampling=spacing_mm)\n",
" else:\n",
" distmap_pred = np.Inf * np.ones(borders_pred.shape)\n",
"\n",
" # compute the area of each surface element\n",
" surface_area_map_gt = neighbour_code_to_surface_area[neighbour_code_map_gt]\n",
" surface_area_map_pred = neighbour_code_to_surface_area[neighbour_code_map_pred]\n",
"\n",
" # create a list of all surface elements with distance and area\n",
" distances_gt_to_pred = distmap_pred[borders_gt]\n",
" distances_pred_to_gt = distmap_gt[borders_pred]\n",
" surfel_areas_gt = surface_area_map_gt[borders_gt]\n",
" surfel_areas_pred = surface_area_map_pred[borders_pred]\n",
"\n",
" # sort them by distance\n",
" if distances_gt_to_pred.shape != (0,):\n",
" sorted_surfels_gt = np.array(sorted(zip(distances_gt_to_pred, surfel_areas_gt)))\n",
" distances_gt_to_pred = sorted_surfels_gt[:,0]\n",
" surfel_areas_gt = sorted_surfels_gt[:,1]\n",
"\n",
" if distances_pred_to_gt.shape != (0,):\n",
" sorted_surfels_pred = np.array(sorted(zip(distances_pred_to_gt, surfel_areas_pred)))\n",
" distances_pred_to_gt = sorted_surfels_pred[:,0]\n",
" surfel_areas_pred = sorted_surfels_pred[:,1]\n",
"\n",
"\n",
" return {\"distances_gt_to_pred\": distances_gt_to_pred, \n",
" \"distances_pred_to_gt\": distances_pred_to_gt, \n",
" \"surfel_areas_gt\": surfel_areas_gt, \n",
" \"surfel_areas_pred\": surfel_areas_pred}\n"
]
},
{
"cell_type": "code",
"execution_count": 0,
"metadata": {
"colab": {
"autoexec": {
"startup": false,
"wait_interval": 0
}
},
"colab_type": "code",
"id": "y8shBcxlrUtJ"
},
"outputs": [],
"source": [
"def compute_average_surface_distance(surface_distances):\n",
" distances_gt_to_pred = surface_distances[\"distances_gt_to_pred\"]\n",
" distances_pred_to_gt = surface_distances[\"distances_pred_to_gt\"]\n",
" surfel_areas_gt = surface_distances[\"surfel_areas_gt\"]\n",
" surfel_areas_pred = surface_distances[\"surfel_areas_pred\"]\n",
" average_distance_gt_to_pred = np.sum( distances_gt_to_pred * surfel_areas_gt) / np.sum(surfel_areas_gt)\n",
" average_distance_pred_to_gt = np.sum( distances_pred_to_gt * surfel_areas_pred) / np.sum(surfel_areas_pred)\n",
" return (average_distance_gt_to_pred, average_distance_pred_to_gt)\n",
"\n",
"def compute_robust_hausdorff(surface_distances, percent):\n",
" distances_gt_to_pred = surface_distances[\"distances_gt_to_pred\"]\n",
" distances_pred_to_gt = surface_distances[\"distances_pred_to_gt\"]\n",
" surfel_areas_gt = surface_distances[\"surfel_areas_gt\"]\n",
" surfel_areas_pred = surface_distances[\"surfel_areas_pred\"]\n",
" if len(distances_gt_to_pred) > 0:\n",
" surfel_areas_cum_gt = np.cumsum(surfel_areas_gt) / np.sum(surfel_areas_gt)\n",
" idx = np.searchsorted(surfel_areas_cum_gt, percent/100.0)\n",
" perc_distance_gt_to_pred = distances_gt_to_pred[min(idx, len(distances_gt_to_pred)-1)]\n",
" else:\n",
" perc_distance_gt_to_pred = np.Inf\n",
" \n",
" if len(distances_pred_to_gt) > 0:\n",
" surfel_areas_cum_pred = np.cumsum(surfel_areas_pred) / np.sum(surfel_areas_pred)\n",
" idx = np.searchsorted(surfel_areas_cum_pred, percent/100.0)\n",
" perc_distance_pred_to_gt = distances_pred_to_gt[min(idx, len(distances_pred_to_gt)-1)]\n",
" else:\n",
" perc_distance_pred_to_gt = np.Inf\n",
" \n",
" return max( perc_distance_gt_to_pred, perc_distance_pred_to_gt)\n",
"\n",
"def compute_surface_overlap_at_tolerance(surface_distances, tolerance_mm):\n",
" distances_gt_to_pred = surface_distances[\"distances_gt_to_pred\"]\n",
" distances_pred_to_gt = surface_distances[\"distances_pred_to_gt\"]\n",
" surfel_areas_gt = surface_distances[\"surfel_areas_gt\"]\n",
" surfel_areas_pred = surface_distances[\"surfel_areas_pred\"]\n",
" rel_overlap_gt = np.sum(surfel_areas_gt[distances_gt_to_pred <= tolerance_mm]) / np.sum(surfel_areas_gt)\n",
" rel_overlap_pred = np.sum(surfel_areas_pred[distances_pred_to_gt <= tolerance_mm]) / np.sum(surfel_areas_pred)\n",
" return (rel_overlap_gt, rel_overlap_pred)\n",
"\n",
"def compute_surface_dice_at_tolerance(surface_distances, tolerance_mm):\n",
" distances_gt_to_pred = surface_distances[\"distances_gt_to_pred\"]\n",
" distances_pred_to_gt = surface_distances[\"distances_pred_to_gt\"]\n",
" surfel_areas_gt = surface_distances[\"surfel_areas_gt\"]\n",
" surfel_areas_pred = surface_distances[\"surfel_areas_pred\"]\n",
" overlap_gt = np.sum(surfel_areas_gt[distances_gt_to_pred <= tolerance_mm])\n",
" overlap_pred = np.sum(surfel_areas_pred[distances_pred_to_gt <= tolerance_mm])\n",
" surface_dice = (overlap_gt + overlap_pred) / (\n",
" np.sum(surfel_areas_gt) + np.sum(surfel_areas_pred))\n",
" return surface_dice"
]
},
{
"cell_type": "code",
"execution_count": 0,
"metadata": {
"colab": {
"autoexec": {
"startup": false,
"wait_interval": 0
}
},
"colab_type": "code",
"id": "8XgG0aHkHIej"
},
"outputs": [],
"source": [
"def compute_dice_coefficient(mask_gt, mask_pred):\n",
" \"\"\"Compute soerensen-dice coefficient.\n",
"\n",
" compute the soerensen-dice coefficient between the ground truth mask `mask_gt`\n",
" and the predicted mask `mask_pred`. \n",
" \n",
" Args:\n",
" mask_gt: 3-dim Numpy array of type bool. The ground truth mask.\n",
" mask_pred: 3-dim Numpy array of type bool. The predicted mask.\n",
"\n",
" Returns:\n",
" the dice coeffcient as float. If both masks are empty, the result is NaN\n",
" \"\"\"\n",
" volume_sum = mask_gt.sum() + mask_pred.sum()\n",
" if volume_sum == 0:\n",
" return np.NaN\n",
" volume_intersect = (mask_gt & mask_pred).sum()\n",
" return 2*volume_intersect / volume_sum\n",
" "
]
},
{
"cell_type": "markdown",
"metadata": {
"colab_type": "text",
"id": "4VLjvdslp1ZW"
},
"source": [
"# Some Simple Tests"
]
},
{
"cell_type": "code",
"execution_count": 0,
"metadata": {
"colab": {
"autoexec": {
"startup": false,
"wait_interval": 0
},
"base_uri": "https://localhost:8080/",
"height": 153
},
"colab_type": "code",
"executionInfo": {
"elapsed": 520,
"status": "ok",
"timestamp": 1528999072261,
"user": {
"displayName": "Olaf Ronneberger",
"photoUrl": "//lh5.googleusercontent.com/-Hp2EnPvUOmU/AAAAAAAAAAI/AAAAAAAAAFc/QR20zxrgwZ0/s50-c-k-no/photo.jpg",
"userId": "110531354013919487725"
},
"user_tz": -60
},
"id": "Vl25cjRakoSb",
"outputId": "8b9bec94-7fa4-41cf-d8bf-f1a92db4dda8"
},
"outputs": [
{
"name": "stdout",
"output_type": "stream",
"text": [
"bounding box min = [50 60 70]\n",
"bounding box max = [50 60 72]\n",
"average surface distance: (1.5, 1.5) mm\n",
"hausdorff (100%): 2.0 mm\n",
"hausdorff (95%): 2.0 mm\n",
"surface overlap at 1mm: (0.5, 0.5)\n",
"surface dice at 1mm: 0.5\n",
"volumetric dice: 0.0\n"
]
}
],
"source": [
"# single pixels, 2mm away\n",
"mask_gt = np.zeros((128,128,128), np.uint8)\n",
"mask_pred = np.zeros((128,128,128), np.uint8)\n",
"mask_gt[50,60,70] = 1\n",
"mask_pred[50,60,72] = 1\n",
"surface_distances = compute_surface_distances(mask_gt, mask_pred, spacing_mm=(3,2,1))\n",
"print(\"surface dice at 1mm: {}\".format(compute_surface_dice_at_tolerance(surface_distances, 1)))\n",
"print(\"volumetric dice: {}\".format(compute_dice_coefficient(mask_gt, mask_pred)))"
]
},
{
"cell_type": "code",
"execution_count": 0,
"metadata": {
"colab": {
"autoexec": {
"startup": false,
"wait_interval": 0
},
"base_uri": "https://localhost:8080/",
"height": 204
},
"colab_type": "code",
"executionInfo": {
"elapsed": 763,
"status": "ok",
"timestamp": 1528999073141,
"user": {
"displayName": "Olaf Ronneberger",
"photoUrl": "//lh5.googleusercontent.com/-Hp2EnPvUOmU/AAAAAAAAAAI/AAAAAAAAAFc/QR20zxrgwZ0/s50-c-k-no/photo.jpg",
"userId": "110531354013919487725"
},
"user_tz": -60
},
"id": "HKp9I52Ik26R",
"outputId": "f1ecf88c-8114-4139-f5a9-885330249a97"
},
"outputs": [
{
"name": "stdout",
"output_type": "stream",
"text": [
"bounding box min = [0 0 0]\n",
"bounding box max = [50 99 99]\n",
"average surface distance: (0.3224470726766679, 0.3393962930605405) mm\n",
"hausdorff (100%): 2.0 mm\n",
"hausdorff (95%): 2.0 mm\n",
"surface overlap at 1mm: (0.8420663973547752, 0.8303018534697297)\n",
"surface dice at 1mm: 0.836145008498\n",
"volumetric dice: 0.990099009901\n",
"\n",
"expected average_distance_gt_to_pred = 1./6 * 2mm = 0.333333333333mm\n",
"expected volumetric dice: 0.990099009901\n"
]
}
],
"source": [
"# two cubes. cube 1 is 100x100x100 mm^3 and cube 2 is 102x100x100 mm^3\n",
"mask_gt = np.zeros((100,100,100), np.uint8)\n",
"mask_pred = np.zeros((100,100,100), np.uint8)\n",
"spacing_mm=(2,1,1)\n",
"mask_gt[0:50, :, :] = 1\n",
"mask_pred[0:51, :, :] = 1\n",
"surface_distances = compute_surface_distances(mask_gt, mask_pred, spacing_mm)\n",
"print(\"surface dice at 1mm: {}\".format(compute_surface_dice_at_tolerance(surface_distances, 1)))\n",
"print(\"volumetric dice: {}\".format(compute_dice_coefficient(mask_gt, mask_pred)))\n",
"print(\"\")\n",
"print(\"expected average_distance_gt_to_pred = 1./6 * 2mm = {}mm\".format(1./6 * 2)) \n",
"print(\"expected volumetric dice: {}\".format(2.*100*100*100 / (100*100*100 + 102*100*100) ))"
]
},
{
"cell_type": "code",
"execution_count": 0,
"metadata": {
"colab": {
"autoexec": {
"startup": false,
"wait_interval": 0
},
"base_uri": "https://localhost:8080/",
"height": 204
},
"colab_type": "code",
"executionInfo": {
"elapsed": 463,
"status": "ok",
"timestamp": 1528999073635,
"user": {
"displayName": "Olaf Ronneberger",
"photoUrl": "//lh5.googleusercontent.com/-Hp2EnPvUOmU/AAAAAAAAAAI/AAAAAAAAAFc/QR20zxrgwZ0/s50-c-k-no/photo.jpg",
"userId": "110531354013919487725"
},
"user_tz": -60
},
"id": "1EgXUCqNqI3X",
"outputId": "6f7946b5-1bc3-4cc4-b489-82801e7aad2c"
},
"outputs": [
{
"name": "stdout",
"output_type": "stream",
"text": [
"bounding box min = [50 60 70]\n",
"bounding box max = [50 60 70]\n",
"average surface distance: (inf, nan) mm\n",
"hausdorff (100%): inf mm\n",
"hausdorff (95%): inf mm\n",
"surface overlap at 1mm: (0.0, nan)\n",
"surface dice at 1mm: 0.0\n",
"volumetric dice: 0.0\n"
]
},
{
"name": "stderr",
"output_type": "stream",
"text": [
"/usr/local/lib/python2.7/dist-packages/ipykernel_launcher.py:7: RuntimeWarning: invalid value encountered in double_scalars\n",
" import sys\n",
"/usr/local/lib/python2.7/dist-packages/ipykernel_launcher.py:37: RuntimeWarning: invalid value encountered in double_scalars\n"
]
}
],
"source": [
"# test empty mask in prediction\n",
"mask_gt = np.zeros((128,128,128), np.uint8)\n",
"mask_pred = np.zeros((128,128,128), np.uint8)\n",
"mask_gt[50,60,70] = 1\n",
"#mask_pred[50,60,72] = 1\n",
"surface_distances = compute_surface_distances(mask_gt, mask_pred, spacing_mm=(3,2,1))\n",
"print(\"average surface distance: {} mm\".format(compute_average_surface_distance(surface_distances)))\n",
"print(\"hausdorff (100%): {} mm\".format(compute_robust_hausdorff(surface_distances, 100)))\n",
"print(\"hausdorff (95%): {} mm\".format(compute_robust_hausdorff(surface_distances, 95)))\n",
"print(\"surface overlap at 1mm: {}\".format(compute_surface_overlap_at_tolerance(surface_distances, 1)))\n",
"print(\"surface dice at 1mm: {}\".format(compute_surface_dice_at_tolerance(surface_distances, 1)))\n",
"print(\"volumetric dice: {}\".format(compute_dice_coefficient(mask_gt, mask_pred)))"
]
},
{
"cell_type": "code",
"execution_count": 0,
"metadata": {
"colab": {
"autoexec": {
"startup": false,
"wait_interval": 0
},
"base_uri": "https://localhost:8080/",
"height": 204
},
"colab_type": "code",
"executionInfo": {
"elapsed": 468,
"status": "ok",
"timestamp": 1528999074187,
"user": {
"displayName": "Olaf Ronneberger",
"photoUrl": "//lh5.googleusercontent.com/-Hp2EnPvUOmU/AAAAAAAAAAI/AAAAAAAAAFc/QR20zxrgwZ0/s50-c-k-no/photo.jpg",
"userId": "110531354013919487725"
},
"user_tz": -60
},
"id": "TcBi7tcNaK3w",
"outputId": "a0f488d1-d11a-4a71-f6be-38f8dbaca57b"
},
"outputs": [
{
"name": "stdout",
"output_type": "stream",
"text": [
"bounding box min = [50 60 72]\n",
"bounding box max = [50 60 72]\n",
"average surface distance: (nan, inf) mm\n",
"hausdorff (100%): inf mm\n",
"hausdorff (95%): inf mm\n",
"surface overlap at 1mm: (nan, 0.0)\n",
"surface dice at 1mm: 0.0\n",
"volumetric dice: 0.0\n"
]
},
{
"name": "stderr",
"output_type": "stream",
"text": [
"/usr/local/lib/python2.7/dist-packages/ipykernel_launcher.py:6: RuntimeWarning: invalid value encountered in double_scalars\n",
" \n",
"/usr/local/lib/python2.7/dist-packages/ipykernel_launcher.py:36: RuntimeWarning: invalid value encountered in double_scalars\n"
]
}
],
"source": [
"# test empty mask in ground truth\n",
"mask_gt = np.zeros((128,128,128), np.uint8)\n",
"mask_pred = np.zeros((128,128,128), np.uint8)\n",
"#mask_gt[50,60,70] = 1\n",
"mask_pred[50,60,72] = 1\n",
"surface_distances = compute_surface_distances(mask_gt, mask_pred, spacing_mm=(3,2,1))\n",
"print(\"average surface distance: {} mm\".format(compute_average_surface_distance(surface_distances)))\n",
"print(\"hausdorff (100%): {} mm\".format(compute_robust_hausdorff(surface_distances, 100)))\n",
"print(\"hausdorff (95%): {} mm\".format(compute_robust_hausdorff(surface_distances, 95)))\n",
"print(\"surface overlap at 1mm: {}\".format(compute_surface_overlap_at_tolerance(surface_distances, 1)))\n",
"print(\"surface dice at 1mm: {}\".format(compute_surface_dice_at_tolerance(surface_distances, 1)))\n",
"print(\"volumetric dice: {}\".format(compute_dice_coefficient(mask_gt, mask_pred)))"
]
},
{
"cell_type": "code",
"execution_count": 0,
"metadata": {
"colab": {
"autoexec": {
"startup": false,
"wait_interval": 0
},
"base_uri": "https://localhost:8080/",
"height": 238
},
"colab_type": "code",
"executionInfo": {
"elapsed": 476,
"status": "ok",
"timestamp": 1528999074739,
"user": {
"displayName": "Olaf Ronneberger",
"photoUrl": "//lh5.googleusercontent.com/-Hp2EnPvUOmU/AAAAAAAAAAI/AAAAAAAAAFc/QR20zxrgwZ0/s50-c-k-no/photo.jpg",
"userId": "110531354013919487725"
},
"user_tz": -60
},
"id": "bcgpOy5EqMLv",
"outputId": "f929a69f-b1e9-4a51-f424-c02f7340aff7"
},
"outputs": [
{
"name": "stdout",
"output_type": "stream",
"text": [
"average surface distance: (nan, nan) mm\n",
"hausdorff (100%): inf mm\n",
"hausdorff (95%): inf mm\n",
"surface overlap at 1mm: (nan, nan)\n",
"surface dice at 1mm: nan\n",
"volumetric dice: nan\n"
]
},
{
"name": "stderr",
"output_type": "stream",
"text": [
"/usr/local/lib/python2.7/dist-packages/ipykernel_launcher.py:6: RuntimeWarning: invalid value encountered in double_scalars\n",
" \n",
"/usr/local/lib/python2.7/dist-packages/ipykernel_launcher.py:7: RuntimeWarning: invalid value encountered in double_scalars\n",
" import sys\n",
"/usr/local/lib/python2.7/dist-packages/ipykernel_launcher.py:36: RuntimeWarning: invalid value encountered in double_scalars\n",
"/usr/local/lib/python2.7/dist-packages/ipykernel_launcher.py:37: RuntimeWarning: invalid value encountered in double_scalars\n",
"/usr/local/lib/python2.7/dist-packages/ipykernel_launcher.py:48: RuntimeWarning: invalid value encountered in double_scalars\n"
]
}
],
"source": [
"# test both masks empty\n",
"mask_gt = np.zeros((128,128,128), np.uint8)\n",
"mask_pred = np.zeros((128,128,128), np.uint8)\n",
"#mask_gt[50,60,70] = 1\n",
"#mask_pred[50,60,72] = 1\n",
"surface_distances = compute_surface_distances(mask_gt, mask_pred, spacing_mm=(3,2,1))\n",
"print(\"average surface distance: {} mm\".format(compute_average_surface_distance(surface_distances)))\n",
"print(\"hausdorff (100%): {} mm\".format(compute_robust_hausdorff(surface_distances, 100)))\n",
"print(\"hausdorff (95%): {} mm\".format(compute_robust_hausdorff(surface_distances, 95)))\n",
"print(\"surface overlap at 1mm: {}\".format(compute_surface_overlap_at_tolerance(surface_distances, 1)))\n",
"print(\"surface dice at 1mm: {}\".format(compute_surface_dice_at_tolerance(surface_distances, 1)))\n",
"print(\"volumetric dice: {}\".format(compute_dice_coefficient(mask_gt, mask_pred)))"
]
},
{
"cell_type": "code",
"execution_count": 0,
"metadata": {
"colab": {
"autoexec": {
"startup": false,
"wait_interval": 0
}
},
"colab_type": "code",
"id": "sVyRElyxYzzV"
},
"outputs": [],
"source": []
}
],
"metadata": {
"colab": {
"collapsed_sections": [],
"default_view": {},
"name": "Surface_distance_based_measures.ipynb",
"provenance": [],
"version": "0.3.2",
"views": {}
},
"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.0"
}
},
"nbformat": 4,
"nbformat_minor": 1
}