generators

Generate Artificial Images

This module contains a variety of functions for generating artificial images of porous materials, generally for testing, validation, debugging, and illustration purposes.

porespy.generators.blobs(shape, porosity, …) Generates an image containing amorphous blobs
porespy.generators.bundle_of_tubes(shape, …) Create a 3D image of a bundle of tubes, in the form of a rectangular plate with randomly sized holes through it.
porespy.generators.cylinders(shape, radius, …) Generates a binary image of overlapping cylinders.
porespy.generators.generate_noise(shape[, …]) Generate a field of spatially correlated random noise using the Perlin noise algorithm, or the updated Simplex noise algorithm.
porespy.generators.insert_shape(im, element) Inserts sub-image into a larger image at the specified location.
porespy.generators.lattice_spheres(shape, …) Generates a cubic packing of spheres in a specified lattice arrangement
porespy.generators.line_segment(X0, X1) Calculate the voxel coordinates of a straight line between the two given end points
porespy.generators.overlapping_spheres(…) Generate a packing of overlapping mono-disperse spheres
porespy.generators.polydisperse_spheres(…) Create an image of randomly place, overlapping spheres with a distribution of radii.
porespy.generators.RSA(im, radius, …) Generates a sphere or disk packing using Random Sequential Addition
porespy.generators.voronoi_edges(shape, …) Create an image of the edges in a Voronoi tessellation
porespy.generators.blobs(shape: List[int], porosity: float = 0.5, blobiness: int = 1)[source]

Generates an image containing amorphous blobs

Parameters:
  • shape (list) – The size of the image to generate in [Nx, Ny, Nz] where N is the number of voxels
  • porosity (float) – If specified, this will threshold the image to the specified value prior to returning. If None is specified, then the scalar noise field is converted to a uniform distribution and returned without thresholding.
  • blobiness (int or list of ints(default = 1)) – Controls the morphology of the blobs. A higher number results in a larger number of small blobs. If a list is supplied then the blobs are anisotropic.
Returns:

image – A boolean array with True values denoting the pore space

Return type:

ND-array

See also

norm_to_uniform()

porespy.generators.bundle_of_tubes(shape: List[int], spacing: int)[source]

Create a 3D image of a bundle of tubes, in the form of a rectangular plate with randomly sized holes through it.

Parameters:
  • shape (list) – The size the image, with the 3rd dimension indicating the plate thickness. If the 3rd dimension is not given then a thickness of 1 voxel is assumed.
  • spacing (scalar) – The center to center distance of the holes. The hole sizes will be randomly distributed between this values down to 3 voxels.
Returns:

image – A boolean array with True values denoting the pore space

Return type:

ND-array

porespy.generators.cylinders(shape: List[int], radius: int, ncylinders: int, phi_max: float = 0, theta_max: float = 90)[source]

Generates a binary image of overlapping cylinders. This is a good approximation of a fibrous mat.

Parameters:
  • shape (list) – The size of the image to generate in [Nx, Ny, Nz] where N is the number of voxels. 2D images are not permitted.
  • radius (scalar) – The radius of the cylinders in voxels
  • ncylinders (scalar) – The number of cylinders to add to the domain. Adjust this value to control the final porosity, which is not easily specified since cylinders overlap and intersect different fractions of the domain.
  • theta_max (scalar) – A value between 0 and 90 that controls the amount of rotation in the XY plane, with 0 meaning all fibers point in the X-direction, and 90 meaning they are randomly rotated about the Z axis by as much as +/- 90 degrees.
  • phi_max (scalar) – A value between 0 and 90 that controls the amount that the fibers lie out of the XY plane, with 0 meaning all fibers lie in the XY plane, and 90 meaning that fibers are randomly oriented out of the plane by as much as +/- 90 degrees.
Returns:

image – A boolean array with True values denoting the pore space

Return type:

ND-array

porespy.generators.generate_noise(shape: List[int], porosity=None, octaves: int = 3, frequency: int = 32, mode: str = 'simplex')[source]

Generate a field of spatially correlated random noise using the Perlin noise algorithm, or the updated Simplex noise algorithm.

Parameters:
  • shape (array_like) – The size of the image to generate in [Nx, Ny, Nz] where N is the number of voxels.
  • porosity (float) – If specified, this will threshold the image to the specified value prior to returning. If no value is given (the default), then the scalar noise field is returned.
  • octaves (int) – Controls the texture of the noise, with higher octaves giving more complex features over larger length scales.
  • frequency (array_like) – Controls the relative sizes of the features, with higher frequencies giving larger features. A scalar value will apply the same frequency in all directions, given an isotropic field; a vector value will apply the specified values along each axis to create anisotropy.
  • mode (string) – Which noise algorithm to use, either 'simplex' (default) or 'perlin'.
Returns:

image – If porosity is given, then a boolean array with True values denoting the pore space is returned. If not, then normally distributed and spatially correlated randomly noise is returned.

Return type:

ND-array

Notes

This method depends the a package called ‘noise’ which must be compiled. It is included in the Anaconda distribution, or a platform specific binary can be downloaded.

porespy.generators.insert_shape(im, element, center=None, corner=None, value=1, mode='overwrite')[source]

Inserts sub-image into a larger image at the specified location.

If the inserted image extends beyond the boundaries of the image it will be cropped accordingly.

Parameters:
  • im (ND-array) – The image into which the sub-image will be inserted
  • element (ND-array) – The sub-image to insert
  • center (tuple) – Coordinates indicating the position in the main image where the inserted imaged will be centered. If center is given then corner cannot be specified. Note that center can only be used if all dimensions of element are odd, otherwise the meaning of center is not defined.
  • corner (tuple) – Coordinates indicating the position in the main image where the lower corner (i.e. [0, 0, 0]) of the inserted image should be anchored. If corner is given then corner cannot be specified.
  • value (scalar) – A scalar value to apply to the sub-image. The default is 1.
  • mode (string) – If ‘overwrite’ (default) the inserted image replaces the values in the main image. If ‘overlay’ the inserted image is added to the main image. In both cases the inserted image is multiplied by value first.
Returns:

im – A copy of im with the supplied element inserted.

Return type:

ND-array

porespy.generators.lattice_spheres(shape: List[int], radius: int, offset: int = 0, lattice: str = 'sc')[source]

Generates a cubic packing of spheres in a specified lattice arrangement

Parameters:
  • shape (list) – The size of the image to generate in [Nx, Ny, Nz] where N is the number of voxels in each direction. For a 2D image, use [Nx, Ny].
  • radius (scalar) – The radius of spheres (circles) in the packing
  • offset (scalar) – The amount offset (+ or -) to add between sphere centers.
  • lattice (string) –

    Specifies the type of lattice to create. Options are:

    ’sc’ - Simple Cubic (default)

    ’fcc’ - Face Centered Cubic

    ’bcc’ - Body Centered Cubic

    For 2D images, ‘sc’ gives a square lattice and both ‘fcc’ and ‘bcc’ give a triangular lattice.

Returns:

image – A boolean array with True values denoting the pore space

Return type:

ND-array

porespy.generators.line_segment(X0, X1)[source]

Calculate the voxel coordinates of a straight line between the two given end points

Parameters:and X1 (X0) – The [x, y] or [x, y, z] coordinates of the start and end points of the line.
Returns:coords – A list of lists containing the X, Y, and Z coordinates of all voxels that should be drawn between the start and end points to create a solid line.
Return type:list of lists
porespy.generators.overlapping_spheres(shape: List[int], radius: int, porosity: float, iter_max: int = 10, tol: float = 0.01)[source]

Generate a packing of overlapping mono-disperse spheres

Parameters:
  • shape (list) – The size of the image to generate in [Nx, Ny, Nz] where Ni is the number of voxels in the i-th direction.
  • radius (scalar) – The radius of spheres in the packing.
  • porosity (scalar) – The porosity of the final image, accurate to the given tolerance.
  • iter_max (int) – Maximum number of iterations for the iterative algorithm that improves the porosity of the final image to match the given value.
  • tol (float) – Tolerance for porosity of the final image compared to the given value.
Returns:

image – A boolean array with True values denoting the pore space

Return type:

ND-array

Notes

This method can also be used to generate a dispersion of hollows by treating porosity as solid volume fraction and inverting the returned image.

porespy.generators.polydisperse_spheres(shape: List[int], porosity: float, dist, nbins: int = 5, r_min: int = 5)[source]

Create an image of randomly place, overlapping spheres with a distribution of radii.

Parameters:
  • shape (list) – The size of the image to generate in [Nx, Ny, Nz] where Ni is the number of voxels in each direction. If shape is only 2D, then an image of polydisperse disks is returns
  • porosity (scalar) – The porosity of the image, defined as the number of void voxels divided by the number of voxels in the image. The specified value is only matched approximately, so it’s suggested to check this value after the image is generated.
  • dist (scipy.stats distribution object) – This should be an initialized distribution chosen from the large number of options in the scipy.stats submodule. For instance, a normal distribution with a mean of 20 and a standard deviation of 10 can be obtained with dist = scipy.stats.norm(loc=20, scale=10)
  • nbins (scalar) – The number of discrete sphere sizes that will be used to generate the image. This function generates nbins images of monodisperse spheres that span 0.05 and 0.95 of the possible values produced by the provided distribution, then overlays them to get polydispersivity.
Returns:

image – A boolean array with True values denoting the pore space

Return type:

ND-array

porespy.generators.RSA(im: numpy.array, radius: int, volume_fraction: int = 1, mode: str = 'extended')[source]

Generates a sphere or disk packing using Random Sequential Addition

This which ensures that spheres do not overlap but does not guarantee they are tightly packed.

Parameters:
  • im (ND-array) – The image into which the spheres should be inserted. By accepting an image rather than a shape, it allows users to insert spheres into an already existing image. To begin the process, start with an array of zero such as im = np.zeros([200, 200], dtype=bool).
  • radius (int) – The radius of the disk or sphere to insert.
  • volume_fraction (scalar) – The fraction of the image that should be filled with spheres. The spheres are addeds 1’s, so each sphere addition increases the volume_fraction until the specified limit is reach.
  • mode (string) –

    Controls how the edges of the image are handled. Options are:

    ’extended’ - Spheres are allowed to extend beyond the edge of the image

    ’contained’ - Spheres are all completely within the image

    ’periodic’ - The portion of a sphere that extends beyond the image is inserted into the opposite edge of the image (Not Implemented Yet!)

Returns:

image – A copy of im with spheres of specified radius added to the background.

Return type:

ND-array

Notes

Each sphere is filled with 1’s, but the center is marked with a 2. This allows easy boolean masking to extract only the centers, which can be converted to coordinates using scipy.where and used for other purposes. The obtain only the spheres, use``im = im == 1``.

This function adds spheres to the background of the received im, which allows iteratively adding spheres of different radii to the unfilled space.

References

[1] Random Heterogeneous Materials, S. Torquato (2001)

porespy.generators.voronoi_edges(shape: List[int], radius: int, ncells: int, flat_faces: bool = True)[source]

Create an image of the edges in a Voronoi tessellation

Parameters:
  • shape (array_like) – The size of the image to generate in [Nx, Ny, Nz] where Ni is the number of voxels in each direction.
  • radius (scalar) – The radius to which Voronoi edges should be dilated in the final image.
  • ncells (scalar) – The number of Voronoi cells to include in the tesselation.
  • flat_faces (Boolean) – Whether the Voronoi edges should lie on the boundary of the image (True), or if edges outside the image should be removed (False).
Returns:

image – A boolean array with True values denoting the pore space

Return type:

ND-array