WIP update after removing phantom generator and ICV estimation

Author: SebastianoF

CONTRIBUTE.md

@@ -1,12 +1,11 @@
-# Contributing to nilabel
+# Contributing to NiLabels
 
 Thank you for being there!
 
-Nilabel (ex LABelsToolkit) started as a python package containing a range of heterogeneous imaging tools to perform
+NiLabels (ex LABelsToolkit) started as a python package containing a range of heterogeneous imaging tools to perform
 quick manipulations and measurement on segmentations from ipython or jupyter notebook. 
 Initially planned to support several projects undertook by the initial author, after some development and refactoring 
- it is not intended to be part of the Nipy ecosystem, to provide the community with another hopefully helpful tool 
-  in neuroimaging research.
+ it is now intended to be part of the Nipy ecosystem, to provide the neuroimaging developer community with another tool.
 
 ## Code of Conduct
 

README.md

@@ -52,10 +52,6 @@ nis_app.manipulate_labels.relabel('my_input_folder_path/my_segm.nii.gz', 'my_out
 ### How to install (in development mode) 
 
 
-+ Install python requirements in requirements.txt with
-    `pip install -r requirements.txt`
-in a [virtualenv](http://docs.python-guide.org/en/latest/dev/virtualenvs/).
-
 + Install NiLabels: the current version is not (yet) pip installable. It can be installed in development mode.
 To proceede, initialise a virtual environment and execute the following instructions:
 ```
@@ -67,6 +63,11 @@ pip install -e .
 In development mode every change made to your local code will be directly affecting the libray installed in the python distribution
 without the need of reinstalling.
 
+
++ Install python requirements in requirements.txt with
+    `pip install -r requirements.txt`
+in a [virtualenv](http://docs.python-guide.org/en/latest/dev/virtualenvs/).
+
 + For advanced method `symmetrise_wit_registration`, extra examples and quick arrays visualisation with ITK-snap you must
     + Install [NiftySeg](https://github.com/KCL-BMEIS/NiftySeg)
     + Install [NiftyReg](https://github.com/KCL-BMEIS/niftyreg)

examples/a_generate_headlike_phantoms.py

@@ -5,24 +5,116 @@
 import numpy as np
 import scipy.ndimage.filters as fil
 
-from nilabels.tools.defs import root_dir
+from nilabels.defs import root_dir
 from nilabels.tools.detections.get_segmentation import intensity_segmentation
-from nilabels.tools.phantoms_generator.shapes_for_phantoms import ellipsoid_shape, o_shape, c_shape, \
-    cube_shape
 
 
-def generate_figures(creation_list, segmentation_levels=7, sigma_smoothing=6, foreground=10, noise=5000):
-    print(creation_list)
-    print('\n.\n.\n\nGenerating figures for examples, may take some seconds and approx 150MB for the whole '
-          'dataset.\n.\n.')
-    pfo_examples = jph(root_dir, 'data_examples')
+# ---------- Simple shapes generators ---------------
+
+
+def o_shape(omega=(250, 250), radius=50,
+            background_intensity=0, foreground_intensity=20, dtype=np.uint8):
+    m = background_intensity * np.ones(omega, dtype=dtype)
+
+    if len(omega) == 2:
+        c = [omega[j] / 2 for j in range(len(omega))]
+        for x in range(omega[0]):
+            for y in range(omega[1]):
+                if (x - c[0]) ** 2 + (y - c[1]) ** 2 < radius ** 2:
+                    m[x, y] = foreground_intensity
+    elif len(omega) == 3:
+        c = [omega[j] / 2 for j in range(len(omega))]
+        for x in range(omega[0]):
+            for y in range(omega[1]):
+                for z in range(omega[2]):
+                    if (x - c[0]) ** 2 + (y - c[1]) ** 2 + (z - c[2]) ** 2 < radius ** 2:
+                        m[x, y, z] = foreground_intensity
+    return m
+
+
+def c_shape(omega=(250, 250), internal_radius=40, external_radius=60, opening_height=50,
+            background_intensity=0, foreground_intensity=20, dtype=np.uint8, margin=None):
+
+    def get_a_2d_c(omega2, internal_radius2d, external_radius2d, opening_height2d, background_intensity2d,
+                   foreground_intensity2d, dtype2d):
+
+        m = background_intensity2d * np.ones(omega2[:2], dtype=dtype2d)
+
+        c = [omega2[j] / 2 for j in range(len(omega2))]
+        # create the crown
+        for x in range(omega2[0]):
+            for y in range(omega2[1]):
+                if internal_radius2d ** 2 < (x - c[0]) ** 2 + (y - c[1]) ** 2 < external_radius2d ** 2:
+                    m[x, y] = foreground_intensity2d
+
+        # open the c
+        low_lim = int(omega2[0] / 2) - int(opening_height2d / 2)
+        high_lim = int(omega2[0] / 2) + int(opening_height2d / 2)
+
+        for x in range(omega2[0]):
+            for y in range(int(omega2[1] / 2), omega2[1]):
+                if low_lim < x < high_lim and m[x, y] == foreground_intensity2d:
+                    m[x, y] = background_intensity2d
+
+        return m
+
+    c_2d = get_a_2d_c(omega2=omega[:2], internal_radius2d=internal_radius, external_radius2d=external_radius,
+                      opening_height2d=opening_height, background_intensity2d=background_intensity,
+                      foreground_intensity2d=foreground_intensity, dtype2d=dtype)
+
+    if len(omega) == 2:
+        return c_2d
+
+    elif len(omega) == 3:
+        if margin is None:
+            return np.repeat(c_2d, omega[2]).reshape(omega)
+        else:
+            res = np.zeros(omega, dtype=dtype)
+            for z in range(margin, omega[2] - 2 * margin):
+                res[..., z] = c_2d
+            return res
+
+
+def ellipsoid_shape(omega, focus_1, focus_2, distance, background_intensity=0, foreground_intensity=100,
+                    dtype=np.uint8):
+    sky = background_intensity * np.ones(omega, dtype=dtype)
+    for xi in range(omega[0]):
+        for yi in range(omega[1]):
+            for zi in range(omega[2]):
+                if np.sqrt((focus_1[0] - xi) ** 2 + (focus_1[1] - yi) ** 2 + (focus_1[2] - zi) ** 2) + \
+                        np.sqrt((focus_2[0] - xi) ** 2 + (focus_2[1] - yi) ** 2 + (focus_2[2] - zi) ** 2) <= distance:
+                    sky[xi, yi, zi] = foreground_intensity
+    return sky
 
-    if creation_list['Examples folder']:
 
+def cube_shape(omega, center, side_length, background_intensity=0, foreground_intensity=100, dtype=np.uint8):
+    sky = background_intensity * np.ones(omega, dtype=dtype)
+    half_side_length = int(np.ceil(int(side_length / 2)))
+
+    for lx in range(-half_side_length, half_side_length + 1):
+        for ly in range(-half_side_length, half_side_length + 1):
+            for lz in range(-half_side_length, half_side_length + 1):
+                sky[center[0] + lx, center[1] + ly, center[2] + lz] = foreground_intensity
+    return sky
+
+
+def sphere_shape(omega, centre, radius, foreground_intensity=100, dtype=np.uint8):
+    sky = np.zeros(omega, dtype=dtype)
+    for xi in range(omega[0]):
+        for yi in range(omega[1]):
+            for zi in range(omega[2]):
+                if np.sqrt((centre[0] - xi) ** 2 + (centre[1] - yi) ** 2 + (centre[2] - zi) ** 2) <= radius:
+                    sky[xi, yi, zi] = foreground_intensity
+    return sky
+
+
+def generate_figures(creation_list, pfo_examples, segmentation_levels=7, sigma_smoothing=6, foreground=10):
+    print('\n.\n.\n\nGenerate figures for the examples, may take some seconds, and will take approx 150MB.\n.\n.')
+
+    if creation_list['Examples folder']:
         os.system('mkdir -p ' + pfo_examples)
 
     if creation_list['Punt e mes']:
-
         data_o_punt = o_shape(omega=(256, 256, 256), foreground_intensity=foreground, dtype=np.float64)
         data_o_punt = fil.gaussian_filter(data_o_punt, sigma=sigma_smoothing)
 
@@ -54,7 +146,6 @@ def generate_figures(creation_list, segmentation_levels=7, sigma_smoothing=6, fo
         print('mes generated')
 
     if creation_list['C']:
-
         data_c_ = c_shape(omega=(256, 256, 256), foreground_intensity=foreground)
         data_c_ = fil.gaussian_filter(data_c_, sigma_smoothing)
 
@@ -76,7 +167,7 @@ def generate_figures(creation_list, segmentation_levels=7, sigma_smoothing=6, fo
     if creation_list['Planetaruim']:
 
         omega = (256, 256, 256)
-        centers = [(32, 36, 40),    (100, 61, 94),   (130, 140, 99),  (220, 110, 210),
+        centers = [(32, 36, 40), (100, 61, 94), (130, 140, 99), (220, 110, 210),
                    (224, 220, 216), (156, 195, 162), (126, 116, 157), (36, 146, 46)]
         radii = [21, 55, 34, 13, 21, 55, 34, 13]
         intensities = [100] * 8
@@ -90,7 +181,7 @@ def generate_figures(creation_list, segmentation_levels=7, sigma_smoothing=6, fo
                         if (center[0] - x) ** 2 + (center[1] - y) ** 2 + (center[2] - z) ** 2 <= radius ** 2:
                             sky[x, y, z] = intensity
 
-        planetarium = fil.gaussian_filter(sky, np.min(radii)-np.std(radii))
+        planetarium = fil.gaussian_filter(sky, np.min(radii) - np.std(radii))
 
         nib_planetarium = nib.Nifti1Image(planetarium, affine=np.eye(4))
         nib.save(nib_planetarium, filename=jph(pfo_examples, 'planetarium.nii.gz'))
@@ -103,12 +194,12 @@ def generate_figures(creation_list, segmentation_levels=7, sigma_smoothing=6, fo
         omega = (120, 140, 160)
         epsilon = 5  # perturbation set to 0 to restore the axial symmetry
         foci_ellipses_left = [np.array([35, 70, 40]), np.array([50, 70, 120])]
-        foci_ellipses_right = [np.array([85 + epsilon, 70 - epsilon, 40+epsilon]), np.array([70, 70, 120])]
+        foci_ellipses_right = [np.array([85 + epsilon, 70 - epsilon, 40 + epsilon]), np.array([70, 70, 120])]
         d_left = 95
         d_right = 95
 
-        ellipsoid_left  = ellipsoid_shape(omega, foci_ellipses_left[0], foci_ellipses_left[1], d_left,
-                                          foreground_intensity=foreground)
+        ellipsoid_left = ellipsoid_shape(omega, foci_ellipses_left[0], foci_ellipses_left[1], d_left,
+                                         foreground_intensity=foreground)
         ellipsoid_right = ellipsoid_shape(omega, foci_ellipses_right[0], foci_ellipses_right[1], d_right,
                                           foreground_intensity=foreground)
 
@@ -130,12 +221,6 @@ def generate_figures(creation_list, segmentation_levels=7, sigma_smoothing=6, fo
         nib_ellipsoids_seg_half = nib.Nifti1Image(two_ellipsoids_half_segmentation, affine=np.eye(4))
         nib.save(nib_ellipsoids_seg_half, filename=jph(pfo_examples, 'ellipsoids_seg_half.nii.gz'))
 
-        salt_and_pepper = np.zeros_like(two_ellipsoids_segmentation)
-        for x, y, z in zip(*[np.random.randint(1, k, size=noise) for k in salt_and_pepper.shape]):
-            salt_and_pepper[x, y, z] = np.random.uniform(1, 4)
-        nib_ellipsoids_seg_salt_and_pepper = nib.Nifti1Image(two_ellipsoids_segmentation + salt_and_pepper, affine=np.eye(4))
-        nib.save(nib_ellipsoids_seg_salt_and_pepper, filename=jph(pfo_examples, 'ellipsoids_seg_noisy.nii.gz'))
-
         print('Buckle ellipsoids half segmented and whole segmented generated')
 
     if creation_list['Ellipsoids family']:
@@ -159,7 +244,6 @@ def generate_figures(creation_list, segmentation_levels=7, sigma_smoothing=6, fo
 
         # --- Generate random ellipsoids in omega: --
         for k in range(1, num_ellipsoids + 1):
-
             mid_point = np.array([50, 50, 50])
             # Get the first focus in a sphere at (25, 25, 50) with radius 15
             center_first_focus = np.array([25, 50, 50])
@@ -168,7 +252,7 @@ def generate_figures(creation_list, segmentation_levels=7, sigma_smoothing=6, fo
             u, v = np.random.uniform(), np.random.uniform()
             # sphere point picking (Wolfram)
             theta = 2 * np.pi * u
-            phi = np.arccos(2*v - 1)
+            phi = np.arccos(2 * v - 1)
             first_focus = np.array([radius_first_focus * np.cos(theta) * np.sin(phi) + center_first_focus[0],
                                     radius_first_focus * np.sin(theta) * np.sin(phi) + center_first_focus[1],
                                     radius_first_focus * np.cos(phi) + center_first_focus[2]])
@@ -212,7 +296,6 @@ def generate_figures(creation_list, segmentation_levels=7, sigma_smoothing=6, fo
         os.system(cmd)
 
     if creation_list['Cubes in the sky']:
-
         omega = [80, 80, 80]
         cube_a = [[10, 60, 55], 11, 1]
         cube_b = [[50, 55, 42], 17, 2]
@@ -235,13 +318,12 @@ def generate_figures(creation_list, segmentation_levels=7, sigma_smoothing=6, fo
         nib.save(im2, filename=jph(pfo_examples, 'cubes_in_space.nii.gz'))
 
     if creation_list['Sandwich']:
-
         omega = [9, 9, 10]
         sandwich = np.zeros(omega)
 
-        sandwich[:, :2, :]  = 2 * np.ones([9, 2, 10])
+        sandwich[:, :2, :] = 2 * np.ones([9, 2, 10])
         sandwich[:, 2:5, :] = 3 * np.ones([9, 3, 10])
-        sandwich[:, 5:, :]  = 4 * np.ones([9, 4, 10])
+        sandwich[:, 5:, :] = 4 * np.ones([9, 4, 10])
         im_sandwich = nib.Nifti1Image(sandwich, affine=np.diag([0.1, 0.2, 0.3, 1]))
 
         nib.save(im_sandwich, filename=jph(pfo_examples, 'sandwich.nii.gz'))
@@ -298,7 +380,7 @@ def generate_figures(creation_list, segmentation_levels=7, sigma_smoothing=6, fo
         # f = 15
         # a = np.sqrt( (f ** 2 + np.sqrt(f**4 + 4*r**4)) / 2. )
         a = 22
-        f = np.sqrt((a**4 - r ** 4) / float(a**2))
+        f = np.sqrt((a ** 4 - r ** 4) / float(a ** 2))
         assert f < a, [f, a]
         first_focus_23 = [75 - f, 25, 25]  # ellipsoid
         second_focus_23 = [75 + f, 25, 25]
@@ -340,16 +422,16 @@ def generate_figures(creation_list, segmentation_levels=7, sigma_smoothing=6, fo
 
 
 if __name__ == '__main__':
-
-    creation_list = {'Examples folder'    : True,
+    creation_steps = {'Examples folder'   : True,
                       'Punt e mes'        : True,
                       'C'                 : True,
                       'Planetaruim'       : True,
                       'Buckle ellipsoids' : True,
                       'Ellipsoids family' : True,
                       'Cubes in the sky'  : True,
                       'Sandwich'          : True,
-                      'Four-folds'        : True
-                      }
+                      'Four-folds'        : True}
+
+    path_folder_examples = jph(root_dir, 'data_examples')
 
-    generate_figures(creation_list)
+    generate_figures(creation_steps, path_folder_examples)

examples/a_generate_phantoms_for_examples.py

@@ -1,9 +1,8 @@
-import os
 from os.path import join as jph
 
 from nilabels.agents.agents_controller import AgentsController as LT
+from nilabels.defs import root_dir
 from nilabels.tools.caliber.distances import dice_score, covariance_distance, hausdorff_distance
-from nilabels.tools.defs import root_dir
 
 if __name__ == '__main__':
 

examples/caliber_usage_distances_example.py

@@ -4,7 +4,7 @@
 import numpy as np
 
 from nilabels.agents.agents_controller import AgentsController as NiL
-from nilabels.tools.defs import root_dir
+from nilabels.defs import root_dir
 
 if __name__ == '__main__':
 

examples/caliber_usage_example.py

@@ -7,8 +7,8 @@
 import nilabels.tools.image_colors_manipulations.relabeller as rel
 import nilabels.tools.image_shape_manipulations.splitter as sp
 from nilabels.agents.agents_controller import AgentsController
+from nilabels.defs import root_dir
 from nilabels.tools.aux_methods.utils_nib import set_new_data
-from nilabels.tools.defs import root_dir
 
 if __name__ == '__main__':
 

examples/manipulator_example.py

@@ -2,10 +2,9 @@
 from os.path import join as jph
 
 import a_generate_phantoms_for_examples as gen
-
 from nilabels.agents.agents_controller import AgentsController as NiL
+from nilabels.defs import root_dir
 from nilabels.tools.aux_methods.label_descriptor_manager import generate_dummy_label_descriptor
-from nilabels.tools.defs import root_dir
 
 # ---- GENERATE DATA ----
 

examples/prototype_check_imperfections.py

@@ -1,12 +1,12 @@
 import os
 from os.path import join as jph
 
-import a_generate_phantoms_for_examples as gen
 import nibabel as nib
 import numpy as np
 
+import a_generate_phantoms_for_examples as gen
 from nilabels.agents.agents_controller import AgentsController as NiL
-from nilabels.tools.defs import root_dir
+from nilabels.defs import root_dir
 
 # ---- GENERATE DATA ----
 

examples/prototype_clean_a_segmentation.py

@@ -2,9 +2,8 @@
 from os.path import join as jph
 
 import a_generate_phantoms_for_examples as gen
-
 from nilabels.agents.agents_controller import AgentsController as NiL
-from nilabels.tools.defs import root_dir
+from nilabels.defs import root_dir
 
 # ---- GENERATE DATA ----