SoftwareX 16 (2021) 100823
Contents lists available at ScienceDirect
SoftwareX
journal homepage: [Link]/locate/softx
Original software publication
QuantImPy: Minkowski functionals and functions with Python
∗
Arnout M.P. Boelens, Hamdi A. Tchelepi
Department of Energy Resources Engineering, Stanford University, Stanford, CA, 94305, USA
article info a b s t r a c t
Article history: The Minkowski functionals and functions are a family of morphological measures and can be used
Received 19 May 2021 to describe both the morphology (shape) and topology (connectedness) of a system. This paper
Received in revised form 10 September 2021 presents the QuantImPy Python package which can compute both the Minkowski functionals and
Accepted 13 September 2021
functions. In addition, this package can efficiently perform basic morphological operations and compute
Keywords: their distance maps. QuantImPy is easy to install, well-documented, integrated with existing Python
Minkowski functionals packages, and open source.
Minkowski functions © 2021 The Authors. Published by Elsevier B.V. This is an open access article under the CC BY-NC-ND
Morphological operations license ([Link]
Python
Code metadata
Current code version 0.2.1
Permanent link to code/repository used for this code version [Link]
Code Ocean compute capsule –
Legal Code License GPL v3+
Code versioning system used git
Software code languages, tools, and services used C, Cython, Python
Compilation requirements, operating environments & dependencies Linux, MacOS, Windows, Python (v3xx), edt, numpy
If available Link to developer documentation/manual [Link]
Support email for questions boelens@[Link]
1. Motivation and significance lack of well-documented and easy to install libraries to compute
the Minkowski functionals and functions in both 2D and 3D.
In nature, many examples can be found where there is an In this paper the QuantImPy Python package is presented.
interaction between the geometry and thermodynamics of a sys- This package is inspired by and partly based on the QuantIm
tem, or where the structure of a material affects its properties [1– C/C++ library for scientific image processing [10]. The QuantImPy
4]. One way of studying these kinds of problems is by using the package can compute the Minkowski functionals and functions
Minkowski functionals [5]. These functionals are a family of mor- of both 2D and 3D images [11,12]. In addition, it can perform a
phological measures and can be used to describe both the mor- number of fast morphological operations based on the Euclidean
phology (shape) and topology (connectedness) of a system. They distance transform [13,14]. By releasing this package the authors
can either be applied to a geometry directly, or as function of cho- aim to provide a Python package that is (i) easy to install, (ii) well-
sen threshold value, in which case they are called the Minkowski documented, (iii) integrated with existing Python packages, and
functions or measures. Examples of a threshold value include the (iv) open source. This will allows anyone to easily compute the
distance from an interface or a density [1]. The many applications Minkowski functionals and functions.
of Minkowski functionals and functions include medical imaging
[6], astronomy [7–9], statistical physics [1], and porous media [2]. 2. Software description
However, despite these many applications there is an apparent
2.1. Design choices
∗ Corresponding author.
E-mail addresses: boelens@[Link] (Arnout M.P. Boelens), The following design choices where made in the development
tchelepi@[Link] (Hamdi A. Tchelepi). of QuantImPy.
[Link]
2352-7110/© 2021 The Authors. Published by Elsevier B.V. This is an open access article under the CC BY-NC-ND license ([Link]
nc-nd/4.0/).
Arnout M.P. Boelens and Hamdi A. Tchelepi SoftwareX 16 (2021) 100823
Fig. 1. Absolute relative error of the Minkowski functionals, δ M, as function of the maximum size of the domain, L, in pixels in both 2D (a) and 3D (b). The absolute
relative error for M2 (2D) and M3 (3D) is machine precision or zero for the used geometry and has therefore been omitted from this plot.
Easy installation The QuantImPy package has been made avail- • Closing (A • B = (A ⊕ B) ⊖ B)
able on the Python Package Index (PyPI) repository (https:
where B is the disk or sphere shaped structuring element, and
//[Link]/project/quantimpy/). As a consequence it can be
A is the set on which the morphological operation is performed
easily installed on various operating systems.
(phase 1). The dilation of A by B is defined as the set of all the
High-quality documentation Extensive documentation, includ- points such that B hits A. The erosion of A by B is defined as the
ing various examples, is available on [Link] (https:// set of all the points such that B is included in A [22]. In addition,
[Link]/quantimpy/). This should allow anyone the same morphological operations are also available as distance
to get up and running with the code quickly. maps. The erosion and dilation maps are distance maps from the
interface into phase 1 and phase 0, respectively. The opening and
Full integration To avoid duplicating functionality, QuantImPy closing maps give the size of the largest disk or sphere that locally
focuses on fast morphological operations and the computa- fits into phase 1 or phase 0, respectively.
tion of the Minkowski functionals and functions. Because it The morphology module contains the following functions:
operates on 2D and 3D Numpy arrays [15], there is seam-
less integration with many existing Python packages, in-
erode(image, dist, res=None) This function performs the
morphological erosion on the binary Numpy array image.
cluding Scikit-image [16], Pillow [17], and Matplotlib [18].
Both 2D and 3D arrays are supported. dist is the distance
Open source The source code is released under the GNU Gen- away from the interface to which an array is dilated. Op-
eral Public License v3.0 on Github ([Link] tionally, the (anisotropic) resolution of the array can be
boeleman/quantimpy). This allows anyone to download, provided using the Numpy array res. When a resolution
review, and modify the source code. array is provided it needs to be of the same dimension as
the image array.
2.2. Modules dilate(image, dist, res=None) This function performs the
morphological dilation operation on the binary Numpy
The QuantImPy python package consists of two modules: array image. Both 2D and 3D arrays are supported. dist
morphology, and minkowski. This section explains their func- is the distance away from the interface to which an array
tionality. is dilated. Optionally, the (anisotropic) resolution of the
array can be provided using the Numpy array res. When
2.2.1. Morphology a resolution array is provided it needs to be of the same
The morphology module performs basic morphological op- dimension as the image array.
erations and computes various distance maps on both 2D and
3D images with an isotropic or anisotropic resolution. For com- open(erosion, dist, res=None) This function is an alias for
putational performance reasons, both the morphological opera- the function dilate(). Together with the erode() func-
tions and the distance maps are computed using the Euclidean tion this function performs the morphological opening op-
distance transform [13]. For the computation of the Euclidean eration on the binary Numpy array erosion. Both 2D and
distance transform the Multi-Label Anisotropic 3D Euclidean Dis- 3D arrays are supported. dist is the distance away from
tance Transform (MLAEDT-3D) library is used [14], which is based the interface to which an array is opened. Optionally, the
on work found in Refs. [19,20] and [21]. This library computes (anisotropic) resolution of the array can be provided using
exact Euclidean distance maps. The module currently only accepts the Numpy array res. When a resolution array is provided
binary images and the following basic morphological operations it needs to be of the same dimension as the dilation
are supported: array.
• Erosion (A ⊖ B) close(dilation, dist, res=None) This function is an alias
• Dilation (A ⊕ B) for the function erode(). Together with the dilate()
• Opening (A ◦ B = (A ⊖ B) ⊕ B) function, this function performs the morphological closing
2
Arnout M.P. Boelens and Hamdi A. Tchelepi SoftwareX 16 (2021) 100823
Fig. 2. Cross-section of a Bentheimer sandstone sample (a) and two different distance maps, an erosion distance map (b) and opening distance map (c). The size of
this sample is 10243 pixels and the uniform resolution is 1.53 µm per pixel.
operation on the binary Numpy array dilation. Both 2D X , with a smooth boundary, δ X , the following functionals are
and 3D arrays are supported. dist is the distance away computed:
from the interface to which an array is closed. Optionally, ∫
the (anisotropic) resolution of the array can be provided M0 (X ) = ds, (1)
using the Numpy array res. When a resolution array is X
∫
provided it needs to be of the same dimension as the 1
M1 (X ) = dc , and (2)
dilation array. 2π δ X
∫ [ ]
1 1
erode_map(image, res=None) This function computes a mor- M2 (X ) = dc , (3)
2π 2 δ X R
phological erosion map of the binary Numpy array im-
age. Both 2D and 3D arrays are supported. Optionally, the where ds is a surface element and dc is a circumference element.
(anisotropic) resolution of the array can be provided using R is the radius of the local curvature. This results in the follow-
the Numpy array res. When a resolution array is provided ing definitions for the surface area, S = M0 (X ), circumference,
it needs to be of the same dimension as the image array. C = 2π M1 (X ), and the 2D Euler characteristic, χ (X ) = π M2 (X ).
Considering a 3D body, X , with a smooth boundary surface, δ X ,
dilate_map(image, res=None) This function computes a the following functionals are computed:
morphological dilation map of the binary Numpy array ∫
image. Both 2D and 3D arrays are supported. Optionally, M0 (X ) = V = dv, (4)
X
the (anisotropic) resolution of the array can be provided ∫
1
using the Numpy array res. When a resolution array is M1 (X ) = ds, (5)
provided it needs to be of the same dimension as the 8 δX
∫ [ ]
dilation array. 1 1 1 1
M2 (X ) = + ds, and (6)
2π 2 δX 2 R2 R1
open_map(erosion_map, res=None) Together with the ∫ [ ]
3 1
erode_map() function this function computes a morpho- M3 (X ) = ds, (7)
logical opening map of the Numpy array erosion_map. (4π )2 δX R1 R2
Both 2D and 3D arrays are supported. Optionally, the where dv is a volume element and ds is a surface element. R1
(anisotropic) resolution of the array can be provided using and R2 are the principal radii of curvature of surface element
the Numpy array res. When a resolution array is provided ds. This results in the following definitions for the volume, V =
it needs to be of the same dimension as the dilation M0 (X ), surface area, S = 8M1 (X ), integral mean curvature, H =
array. 2π 2 M2 (X ), and the 3D Euler characteristic, χ (X ) = 4π /3M3 (X ).
In addition to the functionals, this module can also compute
close_map(dilation_map, res=None) Together with the the Minkowski functions. In this case the gray-scale value of the
dilate_map() function this function computes a mor- input image is used as the threshold value. When using an erosion
phological closing map of the Numpy array dilation_map. or dilation distance map, this results in the computation of the
Both 2D and 3D arrays are supported. Optionally, the Minkowski functions of a parallel surface [1]. The code used to
(anisotropic) resolution of the array can be provided using compute the Minkowski functionals is based on work found in
the Numpy array res. When a resolution array is provided Refs. [10,12] and [11]. For computational performance the code
it needs to be of the same dimension as the dilation is written in C and wrapped by Cython [23] for a native Python
array. interface.
The minkowski module contains the following functions:
2.2.2. Minkowski
The minkowski module computes the Minkowski functionals functionals(image, res) This function computes the
from binary images and the Minkowski functions from distance Minkowski functionals for the Numpy array image. Both
maps or other gray-scale images. These images can be either 2D 2D and 3D arrays are supported. Optionally, the (anisotropic)
or 3D and both isotropic and anisotropic resolutions are sup- resolution of the array can be provided using the Numpy
ported. The definition of the Minkowski functionals follows the array res. When a resolution array is provided it needs to
convention in the physics literature [1]. Considering a 2D body, be of the same dimension as the image array.
3
Arnout M.P. Boelens and Hamdi A. Tchelepi SoftwareX 16 (2021) 100823
functions_close(closing, res) This function computes 2.3. Illustrative examples
the Minkowski functionals as function of the gray-scale
values in the Numpy array closing. Both 2D and 3D ar-
rays are supported. Optionally, the (anisotropic) resolution This section covers a number of examples to illustrate how
of the array can be provided using the Numpy array res. to use the QuantImPy package. First, a small selection of the
When a resolution array is provided it needs to be of the many examples that can be found in the online documentation is
same dimension as the closing array. presented. The last example is a use case where the Minkowski
functions of a sandstone sample are computed.
The algorithm iterates over all gray-scale values present in
The first example shows how to use the dilate and close
the array, starting at the largest value (white). For every
functions on a 3D image of a cube.
gray-scale value the array is converted into a binary im-
age where values larger than the gray-scale value become
one (white) and all other values become zero (black). For
each of these binary images the minkowski functionals are
computed according to the function functionals().
This function can be used in combination with the mor-
phology() module to compute the Minkowski functions
of different morphological distance maps.
functions_open(opening, res) This function computes the
Minkowski functionals as function of the gray-scale values
in the Numpy array opening. Both 2D and 3D arrays are
supported. Optionally, the (anisotropic) resolution of the
array can be provided using the Numpy array res. When
a resolution array is provided it needs to be of the same
dimension as the opening array.
The algorithm iterates over all gray-scale values present in
the array, starting at the smallest value (black). For every
gray-scale value the array is converted into a binary im-
age where values larger than the gray-scale value become
one (white) and all other values become zero (black). For
each of these binary images the minkowski functionals are
computed according to the function functionals().
This function can be used in combination with the mor-
phology() module to compute the Minkowski functions On a slow laptop, both the dilate and close functions take
of different morphological distance maps. about 0.1 s. The next example shows how to compute the
Minkowski functionals for an image of a ball with both an
To check the convergence of the algorithm that computes the isotropic and an anisotropic resolution.
Minkowski functionals two different geometries have been used.
A square of 512 µm × 512 µm in which a disk is placed with
a radius of 192 µm for the 2D Minkowski functionals, and a
box with a size of 512 µm × 512 µm × 512 µm in which a
sphere is placed with a radius of 192 µm for the 3D Minkowski
functionals. Fig. 1 shows the convergence of the relative error of
the Minkowski functionals, δ M, as function of the maximum size
of the domain in pixels, L. The domain size in µm is constant
which implies that the resolution changes as the domain size in
pixels changes. The absolute relative error is computed for both
an isotropic resolution and an anisotropic resolution, where the
resolution in the x direction is twice as large as in the other direc-
tion(s). The relative errors for M2 (2D) and M3 (3D) are machine
precision or zero for the used geometry and have therefore been
omitted from the plot.
Both Fig. 1(a) and Fig. 1(b) show that the absolute relative
error is generally slightly larger when the resolution is anisotropic
instead of isotropic. However, for all the Minkowski functionals
the absolute relative error decreases as the domain size in pixels
increases. While it is known that M1 is biased as the resolution,
∆ → 0 [11], this Minkowski functional shows convergence
as well. For the computation of the Minkowski functions the The functionals function needs about 0.02 s to execute. The
functionals are used and they thus show the same convergence last example from the online documentation shows how to com-
behavior. pute the Minkowski functions of a ball.
4
Arnout M.P. Boelens and Hamdi A. Tchelepi SoftwareX 16 (2021) 100823
Fig. 3. Minkowski functions of the Bentheimer sandstone sample shown in Fig. 2(a). The functions are based on the erosion map and the opening map show in
Figs. 2 (b) and 2(c), respectively.
This example took several hours to compute. Fig. 2 shows a cross-
section of the original sample, the erosion distance map, and the
opening distance map. The erosion distance map in Fig. 2(b) is
the same as an Euclidean distance transform and a lighter gray-
scale value indicates a larger distance from the interface. In the
opening distance map in Fig. 2(c), the gray scale value indicates
The erode_map function takes about 0.1 s to complete and the the radius of the largest sphere that fits inside a pore space.
In Fig. 3 the different Minkowski functions are shown for the
functions_open function needs about 2 s to compute. erosion and opening maps. The definitions for M0 through M3 are
The final example shows the use case of computing the shown in Eqs. (4) through (7). The two curves show the difference
between gradually eroding all pores at the same time (erosion
Minkowski functions for two different distance maps of a Ben-
map) versus eroding away the pore space pore by pore, starting
theimer sandstone sample. The size of this sample is 10243 with the smallest pores (opening map).
pixels and the uniform resolution is 1.53 µm per pixel. The code 3. Impact
snippet below shows how to load a data file [Link],
The main impact of the QuantImPy package is to lower the
compute an erosion and opening distance map, and determine
barriers to using the Minkowski functionals and functions in sci-
the Minkowski functions. entific research. This software is relevant to any research question
5
Arnout M.P. Boelens and Hamdi A. Tchelepi SoftwareX 16 (2021) 100823
where one has to focus on spatial patterns in order to under- [6] Yamagishi K, Naruto N, Mizukami T, Saito J, et al. Histological classification
stand the physical properties of a system. Examples of possible of non-small cell lung cancer using gabor filtering minkowski functionals
and neighborhood component analysis on ct images. Research Square
applications include, but are not limited to, transport proper-
2021. [Link]
ties in porous media [2], micro-emulsions and colloidal suspen- [7] Schmalzing J, Górski KM. Minkowski functionals used in the morphological
sions, dissipative structures in hydrodynamics, mesophases of liq- analysis of cosmic microwave background anisotropy maps. Mon Not R
uid crystals, Turing patterns occurring in chemical reactions [1], Astron Soc 1998;297:355–65. [Link]
phase behavior (e.g. spinodal decomposition [24] and capillary 01467.x.
[8] Kerscher M, Schmalzing J, Retzlaff J, Borgani S, Buchert T, Gottlöber S,
condensation [25]), nano coatings [26], polymer films, and block- et al. Minkowski functionals of abell/aco clusters. Mon Not R Astron Soc
copolymers phases [27]. In addition, this package can be useful 1997;284:73–84. [Link]
for image analysis applications, potentially in combination with [9] Parroni C, Cardone VF, Maoli R, Scaramella R. Going deep with minkowski
existing machine learning libraries. Overall, there is a very wide functionals of convergence maps. Astron Astrophys 2020;633:A71. http:
//[Link]/10.1051/0004-6361/201935988.
range of applications for QuantImPy in many different fields of
[10] Vogel H-J, Weller U, Schlüter S. Quantification of soil structure based on
science and engineering. minkowski functions. Comput Geosci 2010;36:1236–45. [Link]
10.1016/[Link].2010.03.007.
4. Conclusions [11] Ohser J, Mücklich F. Statistical analysis of microstructures in materials
science. Wiley and Sons; 2000.
[12] Mantz H, Jacobs K, Mecke K. Utilizing minkowski functionals for im-
To summarize, the aim of the QuantImPy package is to provide age analysis: a marching square algorithm. J Stat Mech Theory Exp
a Python package that can compute the Minkowski functionals 2008;2008:P12015. [Link]
and functions which is (i) easy to install, (ii) well-documented, [13] Ragnemalm I. Fast erosion and dilation by contour processing and thresh-
(iii) integrated with standard Python libraries, and (iv) open olding of distance maps. Pattern Recognit Lett 1992;13:161–6. [Link]
[Link]/10.1016/0167-8655(92)90055-5.
source. Installation instructions can be found on [Link]: https://
[14] Seung Lab. Multi-label anisotropic 3d euclidean distance transform
[Link]/project/quantimpy/, documentation on [Link]: https: (mlaedt-3d). 2021, URL [Link]
//[Link]/quantimpy/, and the source code on: https: transform-3d.
//[Link]/boeleman/quantimpy. [15] Harris CR, Millman KJ, van der Walt SJ, Gommers R, Virtanen P, Courna-
peau D, et al. Array programming with numpy. Nature 2020;585:357–62.
[Link]
Declaration of competing interest [16] van der Walt S, Schönberger JL, Nunez-Iglesias J, Boulogne F, Warner JD,
Yager N, et al. Scikit-image: image processing in Python. PeerJ
The authors declare that they have no known competing finan- 2014;2:e453. [Link]
cial interests or personal relationships that could have appeared [17] Clark A, et al. Pillow (pil fork) documentation. 2021, URL [Link]
[Link].
to influence the work reported in this paper.
[18] Hunter JD. Matplotlib: A 2D graphics environment. IEEE Ann Hist Comput
2007;9:90–5. [Link]
Acknowledgments [19] Rosenfeld A, Pfaltz JL. Sequential operations in digital picture processing.
J ACM 1966;13:471–94. [Link]
This work was supported as part of the Center for Mechanistic [20] Meijster A, Roerdink JB, Hesselink WH. A general algorithm for computing
distance transforms in linear time. In: Mathematical morphology and its
Control of Water-Hydrocarbon-Rock Interactions in Unconven-
applications to image and signal processing. Springer; 2002, p. 331–40.
tional and Tight Oil Formations (CMC-UF), an Energy Frontier [Link]
Research Center funded by the U.S. Department of Energy, Office [21] Felzenszwalb PF, Huttenlocher DP. Distance transforms of sampled func-
of Science under DOE (BES) Award DE-SC0019165. tions. Theory Comput 2012;8:415–28. [Link]
v008a019.
In addition, the authors thank Kelly Guan from the Kovscek re-
[22] Serra J. Introduction to mathematical morphology. Comput Vis Graph
search group at Stanford University for providing the Bentheimer Image Process 1986;35:283–305.
sandstone sample shown in the examples. [23] Behnel S, Bradshaw R, Citro C, Dalcin L, Seljebotn DS, Smith K. Cython: The
best of both worlds. Computing in Science & Engineering 2010;13:31–9.
References [Link]
[24] Sofonea V, Mecke K. Morphological characterization of spinodal decom-
position kinetics. Eur Phys J B 1999;8:99–112. [Link]
[1] Mecke KR. Additivity, convexity, and beyond: applications of minkowski
s100510050672.
functionals in statistical physics. In: Statistical physics and spatial statistics.
[25] Boelens AM, Tchelepi HA. The effect of topology on phase behavior under
Springer; 2000, p. 111–84. [Link]
confinement. Processes 2021;9:1220. [Link]
[2] Armstrong RT, McClure JE, Robins V, Liu Z, Arns CH, Schlüter S, et al. Porous
[26] Korpi AG, Ţău Ş, Bramowicz M, Arman A, Kulesza S, Pszczolkowski B, et
media characterization using minkowski functionals: Theories, applications
al. Minkowski functional characterization and fractal analysis of surfaces
and future directions. Transp Porous Media 2018;1–31. [Link]
of titanium nitride films. Mater Res Express 2019;6:086463. [Link]
10.1007/s11242-018-1201-4.
org/10.1088/2053-1591/ab26be.
[3] Noshay A, McGrath JE. Block copolymers: overview and critical survey.
[27] Pinna M, Zvelindovsky AV, Todd S, Goldbeck-Wood G. Cubic phases
Elsevier; 2013.
of block copolymers under shear and electric fields by cell dynamics
[4] Weiner S, Wagner HD. The material bone: structure-mechanical function
simulation. i. spherical phase. J Chem Phys 2006;125:154905. [Link]
relations. Annu Rev Mater Sci 1998;28:271–98. [Link]
[Link]/10.1063/1.2356468.
[Link].28.1.271.
[5] Hadwiger H. Vorlesungenüber inhalt, oberfläche und isoperimetrie, vol. 93.
Springer-Verlag; 1957, [Link]