diff --git a/bob/bio/vein/preprocessor/utils.py b/bob/bio/vein/preprocessor/utils.py
index 84cce9a92bcd05361ed24b1d38af42172602ee66..b5ac3ff571a297d38122a69a0bc359d0341c9e6f 100644
--- a/bob/bio/vein/preprocessor/utils.py
+++ b/bob/bio/vein/preprocessor/utils.py
@@ -179,7 +179,7 @@ def show_image(image):
   img.show()
 
 
-def show_mask_over_image(image, mask, color='red'):
+def draw_mask_over_image(image, mask, color='red'):
   """Plots the mask over the image of a finger, for debugging purposes
 
   Parameters:
@@ -190,6 +190,11 @@ def show_mask_over_image(image, mask, color='red'):
     mask (numpy.ndarray): A 2D numpy.ndarray compose of boolean values
       containing the calculated mask
 
+
+  Returns:
+
+    PIL.Image: An image in PIL format
+
   """
 
   from PIL import Image
@@ -198,7 +203,24 @@ def show_mask_over_image(image, mask, color='red'):
   msk = Image.fromarray((~mask).astype('uint8')*80)
   red = Image.new('RGBA', img.size, color=color)
   img.paste(red, mask=msk)
-  img.show()
+
+  return img
+
+
+def show_mask_over_image(image, mask, color='red'):
+  """Plots the mask over the image of a finger, for debugging purposes
+
+  Parameters:
+
+    image (numpy.ndarray): A 2D numpy.ndarray compose of 8-bit unsigned
+      integers containing the original image
+
+    mask (numpy.ndarray): A 2D numpy.ndarray compose of boolean values
+      containing the calculated mask
+
+  """
+
+  draw_mask_over_image(image, mask, color).show()
 
 
 def jaccard_index(a, b):
diff --git a/bob/bio/vein/script/__init__.py b/bob/bio/vein/script/__init__.py
new file mode 100644
index 0000000000000000000000000000000000000000..e69de29bb2d1d6434b8b29ae775ad8c2e48c5391
diff --git a/bob/bio/vein/script/compare_rois.py b/bob/bio/vein/script/compare_rois.py
new file mode 100644
index 0000000000000000000000000000000000000000..f8578619dfa4d21df458088a42d1432f3a36c2a4
--- /dev/null
+++ b/bob/bio/vein/script/compare_rois.py
@@ -0,0 +1,214 @@
+#!/usr/bin/env python
+# vim: set fileencoding=utf-8 :
+
+"""Compares two set of masks and prints some error metrics
+
+This program requires that the masks for both databases (one representing the
+ground-truth and a second the database with an automated method) are
+represented as HDF5 files containing a ``mask`` object, which should be boolean
+in nature.
+
+
+Usage: %(prog)s [-v...] [-n X] <ground-truth> <database>
+       %(prog)s --help
+       %(prog)s --version
+
+
+Arguments:
+  <ground-truth>  Path to a set of files that contain ground truth annotations
+                  for the ROIs you wish to compare.
+  <database>      Path to a similar set of files as in `<ground-truth>`, but
+                  with ROIs calculated automatically. Every HDF5 in this
+                  directory will be matched to an equivalent file in the
+                  `<ground-truth>` database and their masks will be compared
+
+
+Options:
+  -h, --help          Shows this help message and exits
+  -V, --version       Prints the version and exits
+  -v, --verbose       Increases the output verbosity level
+  -n N, --annotate=N  Print out the N worst cases available in the database,
+                      taking into consideration the various metrics analyzed
+
+
+Example:
+
+  1. Just run for basic statistics:
+
+     $ %(prog)s -vvv gt/ automatic/
+
+  2. Identify worst 5 samples in the database according to a certain criterion:
+
+     $ %(prog)s -vv -n 5 gt/ automatic/
+
+"""
+
+import os
+import sys
+import fnmatch
+import operator
+
+import numpy
+
+import bob.core
+logger = bob.core.log.setup("bob.measure")
+
+import bob.io.base
+
+
+def make_catalog(d):
+  """Returns a catalog dictionary containing the file stems available in ``d``
+
+  Parameters:
+
+    d (str): A path representing a directory that will be scanned for .hdf5
+      files
+
+
+  Returns
+
+    list: A list of stems, from the directory ``d``, that represent files of
+    type HDF5 in that directory. Each file should contain two objects:
+    ``image`` and ``mask``.
+
+  """
+
+  logger.info("Scanning directory `%s'..." % d)
+  retval = []
+  for path, dirs, files in os.walk(d):
+    basedir = os.path.relpath(path, d)
+    logger.debug("Scanning sub-directory `%s'..." % basedir)
+    candidates = fnmatch.filter(files, '*.hdf5')
+    if not candidates: continue
+    logger.debug("Found %d files" % len(candidates))
+    retval += [os.path.join(basedir, k) for k in candidates]
+  logger.info("Found a total of %d files at `%s'" % (len(retval), d))
+  return sorted(retval)
+
+
+def sort_table(table, cols):
+  """Sorts a table by multiple columns
+
+
+  Parameters:
+
+    table (:py:class:`list` of :py:class:`list`): Or tuple of tuples, where
+      each inner list represents a row
+
+    cols (list, tuple): Specifies the column numbers to sort by e.g. (1,0)
+      would sort by column 1, then by column 0
+
+
+  Returns:
+
+    list: of lists, with the table re-ordered as you see fit.
+
+  """
+
+  for col in reversed(cols):
+      table = sorted(table, key=operator.itemgetter(col))
+  return table
+
+
+def mean_std_for_column(table, column):
+  """Calculates the mean and standard deviation for the column in question
+
+
+  Parameters:
+
+    table (:py:class:`list` of :py:class:`list`): Or tuple of tuples, where
+      each inner list represents a row
+
+    col (int): The number of the column from where to extract the data for
+      calculating the mean and the standard-deviation.
+
+
+  Returns:
+
+    float: mean
+
+    float: (unbiased) standard deviation
+
+  """
+
+  z = numpy.array([k[column] for k in table])
+  return z.mean(), z.std(ddof=1)
+
+
+def main(user_input=None):
+
+  if user_input is not None:
+    argv = user_input
+  else:
+    argv = sys.argv[1:]
+
+  import docopt
+  import pkg_resources
+
+  completions = dict(
+      prog=os.path.basename(sys.argv[0]),
+      version=pkg_resources.require('bob.bio.vein')[0].version
+      )
+
+  args = docopt.docopt(
+      __doc__ % completions,
+      argv=argv,
+      version=completions['version'],
+      )
+
+  # Sets-up logging
+  verbosity = int(args['--verbose'])
+  bob.core.log.set_verbosity_level(logger, verbosity)
+
+  # Catalogs
+  gt = make_catalog(args['<ground-truth>'])
+  db = make_catalog(args['<database>'])
+
+  if gt != db:
+    raise RuntimeError("Ground-truth and database have different files!")
+
+  # Calculate all metrics required
+  from ..preprocessor import utils
+  metrics = []
+  for k in gt:
+    logger.info("Evaluating metrics for `%s'..." % k)
+    gt_file = os.path.join(args['<ground-truth>'], k)
+    db_file = os.path.join(args['<database>'], k)
+    gt_roi = bob.io.base.HDF5File(gt_file).read('mask')
+    db_roi = bob.io.base.HDF5File(db_file).read('mask')
+    metrics.append((
+      k,
+      utils.jaccard_index(gt_roi, db_roi),
+      utils.intersect_ratio(gt_roi, db_roi),
+      utils.intersect_ratio_of_complement(gt_roi, db_roi),
+      ))
+
+  # Print statistics
+  names = (
+      (1, 'Jaccard index'),
+      (2, 'Intersection ratio (m1)'),
+      (3, 'Intersection ratio of complement (m2)'),
+      )
+  print("Statistics:")
+  for k, name in names:
+    mean, std = mean_std_for_column(metrics, k)
+    print(name + ': ' + '%.2e +- %.2e' % (mean, std))
+
+  # Print worst cases, if the user asked so
+  if args['--annotate'] is not None:
+    N = int(args['--annotate'])
+    if N <= 0:
+      raise docopt.DocoptExit("Argument to --annotate should be >0")
+
+    print("Worst cases by metric:")
+    for k, name in names:
+      print(name + ':')
+
+      if k in (1,2):
+        worst = sort_table(metrics, (k,))[:N]
+      else:
+        worst = reversed(sort_table(metrics, (k,))[-N:])
+
+      for n, l in enumerate(worst):
+        fname = os.path.join(args['<database>'], l[0])
+        print('  %d. [%.2e] %s' % (n, l[k], fname))
diff --git a/bob/bio/vein/script/view_mask.py b/bob/bio/vein/script/view_mask.py
new file mode 100644
index 0000000000000000000000000000000000000000..2156039f7184cbc3c6284c53dda765d383065803
--- /dev/null
+++ b/bob/bio/vein/script/view_mask.py
@@ -0,0 +1,82 @@
+#!/usr/bin/env python
+# vim: set fileencoding=utf-8 :
+# Mon 07 Nov 2016 15:20:26 CET
+
+
+"""Visualizes masks applied to vein imagery
+
+Usage: %(prog)s [-v...] [options] <file> [<file>...]
+       %(prog)s --help
+       %(prog)s --version
+
+
+Arguments:
+  <file>  The HDF5 file to load image and mask from
+
+
+Options:
+  -h, --help             Shows this help message and exits
+  -V, --version          Prints the version and exits
+  -v, --verbose          Increases the output verbosity level
+  -s path, --save=path   If set, saves image into a file instead of displaying
+                         it
+
+
+Examples:
+
+  Visualize the mask on a single image:
+
+     $ %(prog)s data.hdf5
+
+  Visualize multiple masks (like in a proof-sheet):
+
+     $ %(prog)s *.hdf5
+
+"""
+
+
+import os
+import sys
+
+import bob.core
+logger = bob.core.log.setup("bob.bio.vein")
+
+from ..preprocessor import utils
+
+
+def main(user_input=None):
+
+  if user_input is not None:
+    argv = user_input
+  else:
+    argv = sys.argv[1:]
+
+  import docopt
+  import pkg_resources
+
+  completions = dict(
+      prog=os.path.basename(sys.argv[0]),
+      version=pkg_resources.require('bob.bio.vein')[0].version
+      )
+
+  args = docopt.docopt(
+      __doc__ % completions,
+      argv=argv,
+      version=completions['version'],
+      )
+
+  # Sets-up logging
+  verbosity = int(args['--verbose'])
+  bob.core.log.set_verbosity_level(logger, verbosity)
+
+  # Loads the image, the mask and save it to a PNG file
+  from ..preprocessor import utils
+  for filename in args['<file>']:
+    f = bob.io.base.HDF5File(filename)
+    image = f.read('image')
+    mask  = f.read('mask')
+    img = utils.draw_mask_over_image(image, mask)
+    if args['--save']:
+      img.save(args['--save'])
+    else:
+      img.show()
diff --git a/setup.py b/setup.py
index a91004f45240d9ef6d8e06874800d19130eb7ef8..78a4761d545589f9b4db2ec974b8ac990587dff6 100644
--- a/setup.py
+++ b/setup.py
@@ -45,6 +45,10 @@ setup(
         'parallel = bob.bio.vein.configurations.parallel',
         ],
 
+      'console_scripts': [
+        'compare_rois.py = bob.bio.vein.script.compare_rois:main',
+        'view_mask.py = bob.bio.vein.script.view_mask:main',
+        ]
       },
 
     classifiers = [