method: DROP (#39)

+ method: DROP #8
* convert MHD images
* warping regist. results
* loading mhd deform
* update docs
* update params

Author: Borda

README.md

@@ -215,6 +215,16 @@ For each benchmark experiment, the explanation about how to install and use a pa
         -script ./scripts/Python/run_ANTsPy.py \
         --visual --unique
     ```
+* **[DROP](https://www.mrf-registration.net)** is image registration and motion estimation based on Markov Random Fields.
+    ```bash
+    python bm_experiments/bm_DROP.py \
+        -c ./data_images/pairs-imgs-lnds_histol.csv \
+        -d ./data_images \
+        -o ./results \
+        -DROP ~/Applications/DROP/dropreg2d \
+        --path_config ./configs/drop.txt \
+        --visual --unique
+    ```
 * ...
 
 Some more image registration methods integrated in ImageJ are listed in [Registration](https://imagej.net/Registration).

birl/cls_benchmark.py

@@ -130,8 +130,11 @@ class ImRegBenchmark(Experiment):
     >>> from birl.utilities.data_io import create_folder, update_path
     >>> path_out = create_folder('temp_results')
     >>> path_csv = os.path.join(update_path('data_images'), 'pairs-imgs-lnds_mix.csv')
-    >>> params = {'nb_workers': 1, 'unique': False, 'visual': True,
-    ...           'path_out': path_out, 'path_cover': path_csv}
+    >>> params = {'path_cover': path_csv,
+    ...           'path_out': path_out,
+    ...           'nb_workers': 1,
+    ...           'unique': False,
+    ...           'visual': True}
     >>> benchmark = ImRegBenchmark(params)
     >>> benchmark.run()
     True
@@ -142,8 +145,11 @@ class ImRegBenchmark(Experiment):
     >>> from birl.utilities.data_io import create_folder, update_path
     >>> path_out = create_folder('temp_results')
     >>> path_csv = os.path.join(update_path('data_images'), 'pairs-imgs-lnds_mix.csv')
-    >>> params = {'nb_workers': 2, 'unique': False, 'visual': True,
-    ...           'path_out': path_out, 'path_cover': path_csv}
+    >>> params = {'path_cover': path_csv,
+    ...           'path_out': path_out,
+    ...           'nb_workers': 2,
+    ...           'unique': False,
+    ...           'visual': True}
     >>> benchmark = ImRegBenchmark(params)
     >>> benchmark.run()
     True
@@ -371,7 +377,7 @@ def __convert_gray(path_img_col):
                                (path_img_move, COL_IMAGE_MOVE)]
                 # TODO: find a way how to convert images in parallel inside mproc pool
                 for path_img, col in iterate_mproc_map(__convert_gray, argv_params,
-                                                       nb_workers=1):
+                                                       nb_workers=1, desc=None):
                     record[col + COL_IMAGE_EXT_TEMP] = path_img
             else:
                 logging.warning('unrecognized pre-processing: %s', pproc)

birl/utilities/data_io.py

@@ -11,8 +11,10 @@
 
 import numpy as np
 import pandas as pd
+import SimpleITK as sitk
 from PIL import Image
-from skimage.color import gray2rgb, rgb2hsv, hsv2rgb, rgb2lab, lab2rgb, lch2lab, lab2lch
+from skimage.color import (
+    rgb2gray, gray2rgb, rgb2hsv, hsv2rgb, rgb2lab, lab2rgb, lch2lab, lab2lch)
 
 #: landmarks coordinates, loading from CSV file
 LANDMARK_COORDS = ['X', 'Y']
@@ -344,6 +346,83 @@ def save_image(path_image, image):
     image.save(path_image)
 
 
+@io_image_decorate
+def convert_from_mhd(path_image, path_out_dir=None, img_ext='.png'):
+    """ convert standard image to MHD format
+
+    .. ref:: https://www.programcreek.com/python/example/96382/SimpleITK.WriteImage
+
+    :param str path_image: path to the input image
+    :param str path_out_dir: path to output directory, if None use the input dir
+    :param str img_ext: image extension like PNG or JPEG
+    :return str: path to exported image
+
+    >>> path_img = os.path.join(update_path('data_images'), 'images',
+    ...                         'artificial_reference.jpg')
+    >>> path_img = convert_to_mhd(path_img)
+    >>> convert_from_mhd(path_img)  # doctest: +ELLIPSIS
+    '...artificial_reference.png'
+    """
+    path_image = update_path(path_image)
+    assert os.path.isfile(path_image), 'missing image: %s' % path_image
+    # Reads the image using SimpleITK
+    itk_image = sitk.ReadImage(path_image)
+
+    # Convert the image to a numpy array first and then shuffle the dimensions
+    # to get axis in the order z,y,x
+    img = sitk.GetArrayFromImage(itk_image)
+
+    # define output/destination path
+    img_name = os.path.splitext(os.path.basename(path_image))[0]
+    if not path_out_dir:
+        path_out_dir = os.path.dirname(path_image)
+    path_image = os.path.join(path_out_dir, img_name + img_ext)
+    save_image(path_image, img)
+    return path_image
+
+
+@io_image_decorate
+def convert_to_mhd(path_image, path_out_dir=None, to_gray=True, overwrite=True):
+    """ converting standard image to MHD (Nifty format)
+
+    .. ref:: https://stackoverflow.com/questions/37290631
+
+    :param str path_image: path to the input image
+    :param str path_out_dir: path to output directory, if None use the input dir
+    :param bool overwrite: allow overwrite existing image
+    :return str: path to exported image
+
+    >>> path_img = os.path.join(update_path('data_images'), 'images',
+    ...                         'artificial_moving-affine.jpg')
+    >>> convert_to_mhd(path_img)  # doctest: +ELLIPSIS
+    '...artificial_moving-affine.mhd'
+    """
+    path_image = update_path(path_image)
+    # define output/destination path
+    img_name = os.path.splitext(os.path.basename(path_image))[0]
+    if not path_out_dir:
+        path_out_dir = os.path.dirname(path_image)
+    path_image_new = os.path.join(path_out_dir, img_name + '.mhd')
+    # in case the image exists and you are not allowed to overwrite it
+    if os.path.isfile(path_image_new) and not overwrite:
+        return path_image_new
+
+    img = load_image(path_image)
+    # if required and RGB on input convert to gray-scale
+    if to_gray and img.ndim == 3 and img.shape[2] in (3, 4):
+        img = rgb2gray(img)
+    # the MHD usually require pixel value range (0, 255)
+    if np.max(img) <= 1.5:
+        img = np.round(img * 255)
+    np.clip(img, a_min=0, a_max=255, out=img)
+
+    image = sitk.GetImageFromArray(img.astype(np.uint8), isVector=False)
+
+    # do not use text in MHD, othwerwise it crash DROP method
+    sitk.WriteImage(image, path_image_new, False)
+    return path_image_new
+
+
 def image_histogram_matching(source, reference, use_color='hsv'):
     """ adjust image histogram between two images
 

birl/utilities/experiments.py

@@ -229,8 +229,10 @@ def parse_arg_params(parser, upper_dirs=None):
 def exec_commands(commands, path_logger=None, timeout=None):
     """ run the given commands in system Command Line
 
-    SEE: https://stackoverflow.com/questions/1996518
-    https://www.quora.com/Whats-the-difference-between-os-system-and-subprocess-call-in-Python
+    See refs:
+
+    * https://stackoverflow.com/questions/1996518
+    * https://www.quora.com/Whats-the-difference-between-os-system-and-subprocess-call-in-Python
 
     :param [str] commands: commands to be executed
     :param str path_logger: path to the logger
@@ -258,9 +260,11 @@ def exec_commands(commands, path_logger=None, timeout=None):
     for cmd in commands:
         outputs += [cmd.encode('utf-8')]
         try:
+            # os.system(cmd)
+            # NOTE: for some reason " in the command makes it crash, e.g with DROP
             outputs += [subprocess.check_output(cmd.split(), **options)]
         except subprocess.CalledProcessError as e:
-            logging.exception(commands)
+            logging.exception(cmd)
             outputs += [e.output]
             success = False
     # export the output if path exists

bm_ANHIR/evaluate_submission.py

@@ -228,7 +228,7 @@ def parse_landmarks(idx_row):
 def compute_scores(df_experiments, min_landmarks=1.):
     """ compute all main metrics
 
-    SEE: https://anhir.grand-challenge.org/Evaluation/
+    .. ref:: https://anhir.grand-challenge.org/Evaluation/
 
     :param DF df_experiments: complete experiments
     :param float min_landmarks: required number of submitted landmarks in range (0, 1),

bm_dataset/rescale_tissue_landmarks.py

@@ -51,9 +51,9 @@
 def arg_parse_params():
     """ argument parser from cmd
 
-    SEE: https://docs.python.org/3/library/argparse.html
     :return {str: ...}:
     """
+    # SEE: https://docs.python.org/3/library/argparse.html
     parser = argparse.ArgumentParser()
     parser.add_argument('-a', '--path_annots', type=str, required=False,
                         help='path to folder with annotations')

bm_experiments/bm_ANTsPy.py

@@ -73,11 +73,13 @@ class BmANTsPy(ImRegBenchmark):
     >>> from birl.utilities.data_io import create_folder, update_path
     >>> path_out = create_folder('temp_results')
     >>> fn_path_conf = lambda n: os.path.join(update_path('configs'), n)
-    >>> params = {'nb_workers': 1, 'unique': False,
+    >>> path_csv = os.path.join(update_path('data_images'), 'pairs-imgs-lnds_mix.csv')
+    >>> params = {'path_cover': path_csv,
     ...           'path_out': path_out,
-    ...           'path_cover': os.path.join(update_path('data_images'),
-    ...                                      'pairs-imgs-lnds_mix.csv'),
-    ...           'exec_Python': 'python', 'path_script': '.'}
+    ...           'nb_workers': 2,
+    ...           'unique': False,
+    ...           'exec_Python': 'python',
+    ...           'path_script': '.'}
     >>> benchmark = BmANTsPy(params)
     >>> benchmark.run()  # doctest: +SKIP
     >>> del benchmark