# Library of SPAM functions for manipulating images
# Copyright (C) 2020 SPAM Contributors
#
# This program is free software: you can redistribute it and/or modify it
# under the terms of the GNU General Public License as published by the Free
# Software Foundation, either version 3 of the License, or (at your option)
# any later version.
#
# This program is distributed in the hope that it will be useful, but WITHOUT
# ANY WARRANTY; without even the implied warranty of MERCHANTABILITY or
# FITNESS FOR A PARTICULAR PURPOSE. See the GNU General Public License for
# more details.
#
# You should have received a copy of the GNU General Public License along with
# this program. If not, see <http://www.gnu.org/licenses/>.
import numpy
import progressbar
import spam.label
import tifffile
[docs]
def stackToArray(prefix, suffix=".tif", stack=range(10), digits="05d"):
"""
Convert of stack of 2D sequential tif images into a 3D array.
Parameters
----------
prefix : string
The common name of the 2D images files before the sequential number
suffix : string, default='.tif'
The common name and extension of the 2D images after the sequential number
stack : sequence, default=range(10)
The numbers of the slices with no formating (with no leading zeros)
digits : string, default='05d'
The format (number of digits) of the numbers (add leading zeros).
Returns
-------
numpy array
The 3D image
"""
def _slice_name(pr, st, di, su):
return f"{pr}{st:{di}}{su}"
# Step 1 If nBytes is not defined: we open the first slice just for the dimensions
slice_name = _slice_name(prefix, stack[0], digits, suffix)
slice_im = tifffile.imread(slice_name)
# Step 2 create empty array of good size and type
ny, nx = slice_im.shape
nz = len(stack)
im = numpy.zeros([nz, ny, nx], dtype=slice_im.dtype)
# Step 3 stack all the slices
for i, s in enumerate(stack):
slice_name = _slice_name(prefix, s, digits, suffix)
im[i, :, :] = tifffile.imread(slice_name)
return im
[docs]
def crop(im, boxSize, boxOrigin=None):
"""
This function crops an image using slicing.
Parameters
----------
im: array
The image to crop.
boxSize: int
The size in voxels of the crop box (from boxOrigin). If a int, the size is the same for each axis. If a sequence, ``boxSize`` should contain one value for each axis.
boxOrigin : int, default=None
The coordinates in voxels of the origin of the crop box. If a int, the coordinates are the same for each axis.
If a tuple, ``boxOrigin`` should contain one value for each axis. If ``None`` the image is cropped from its centre.
Returns
-------
array
The cropped image.
"""
# get center of the image
imCentre = [int(s) / 2 for s in im.shape]
# get sizes
(sz, sy, sx) = (boxSize, boxSize, boxSize) if isinstance(boxSize, int) else boxSize
# get box origin
if boxOrigin is not None:
(cz, cy, cx) = (boxOrigin, boxOrigin, boxOrigin) if isinstance(boxOrigin, int) else boxOrigin
else:
(cz, cy, cx) = (
imCentre[0] - sz // 2,
imCentre[1] - sy // 2,
imCentre[2] - sx // 2,
)
# test sizes
if cz + sz > im.shape[0] or cy + sy > im.shape[1] or cx + sx > im.shape[2]:
print("spam.helpers.imageManipulation.crop: box bigger than image.")
print("exit function.")
return -1
return im[int(cz) : int(cz + sz), int(cy) : int(cy + sy), int(cx) : int(cx + sx)]
[docs]
def rescale(im, scale=(0, 1)):
"""
This function **rescales** the values of an image according to a scale
and save it to as 4 bytes floats (float32).
Parameters
----------
im: array
The image to rescale
scale : (float, float), default=(0 1)
The min and max of the rescaled image
Returns
-------
array, float
The rescaled image.
Examples
--------
>>> im = numpy.random.randn( 100, 100, 100 ).astype( '<f4' )
produce float32 array of positive and negative numbers
>>> imRescaled = rescale( im, scale=[-1, 1] )
produce float32 array of numbers between -1 and 1
"""
im_max = float(im.max())
im_min = float(im.min())
return (min(scale) + (max(scale) - min(scale)) * ((im.astype("<f4") - im_min) / (im_max - im_min))).astype("<f4")
[docs]
def rescaleToInteger(im, nBytes=1, scale=None):
"""
This function **rescales** a 4 bytes float image values to a unsigned integers of ``nBytes``.
Parameters
----------
im: float32 numpy array
The image to rescale
nBytes : int, default=1
The number of bytes of the unsigned interger output.
Possible values are power of 2
.. code-block:: text
reminder
1 byte = 8 bits -> ouput from 0 to 255
2 bytes = 16 bits -> ouput from 0 to 65 535
4 bytes = 32 bits -> ouput from 0 to 4 294 967 295
scale : (float, float), default=None
If None, the maximum and minimum use for the rescaling is the maximum and the minimum of the image
Returns
-------
numpy array, uint
The rescaled image
Examples
--------
>>> im = numpy.random.randn( 100, 100, 100 ).astype( '<f4' )
produce float32 array of positive and negative numbers
>>> imRescaled = rescaleToInteger( im, nBytes=4 )
produce uint32 array of numbers between 0 and 4 294 967 295
"""
nBytes = int(nBytes)
if (nBytes & (nBytes - 1)) != 0:
raise ValueError(f"nBytes should be a power of 2 ({nBytes} given)")
# check if float32 given
if im.dtype is not numpy.dtype("float32"):
raise ValueError(f"image should be encode in float32 ({im.dtype} given)")
if scale is None:
# if no scale is given: it takes the max and min of the image
im_max = im.max()
im_min = im.min()
else:
# if a scale is given take it if larger (smaller) than max (min) of image
# im_max = max(scale) if max(scale) > im.max() else im.max()
# im_min = min(scale) if min(scale) < im.min() else im.min()
im_max = max(scale)
im_min = min(scale)
im[im > im_max] = im_max
im[im < im_min] = im_min
im_min = float(im_min)
im_max = float(im_max)
return ((2 ** (8 * nBytes) - 1) * ((im.astype("<f4") - im_min) / (im_max - im_min))).astype(f"<u{nBytes}")
[docs]
def convertUnsignedIntegers(im, nBytes=1):
"""
This function **converts** an images of unsigned integers.
Note: this function does not rescale.
Parameters
----------
im: array, uint
The image to convert.
nBytes : int, default=1
The number of bytes of the unsigned interger output.
Possible values are power of 2.
.. code-block:: text
reminder
1 byte = 8 bits -> ouput from 0 to 255
2 bytes = 16 bits -> ouput from 0 to 65 535
4 bytes = 32 bits -> ouput from 0 to 4 294 967 295
Returns
-------
array, uint
The converted image.
Examples
--------
>>> im = numpy.random.randint( 12, high=210, size=(100, 100, 100) ).astype( '<u1' )
produce an uint8 array of numbers between 12 and 210
>>> imRescaled = rescaleToInteger( im, nBytes=2 )
produce an uint16 array 3084 and 53970
"""
nBytes = int(nBytes)
if (nBytes & (nBytes - 1)) != 0:
raise ValueError(f"nBytes should be a power of 2 ({nBytes} given)")
# number of bits of the output
nbo = 8 * nBytes
# number of bits of the input
inputType = im.dtype
if inputType == numpy.uint8:
nbi = 8
elif inputType == numpy.uint16:
nbi = 16
elif inputType == numpy.uint32:
nbi = 32
elif inputType == numpy.uint64:
nbi = 64
else:
raise ValueError(f"input image type should be unisgned integers ({inputType} given)")
return (float(2**nbo - 1) * (im) / float(2**nbi - 1)).astype("<u{}".format(nBytes))
[docs]
def singleShift(im, shift, axis, sub=0):
"""
This function shift the image and replace the border by an substitution value.
It uses ``numpy.roll``.
Parameters
-----------
im : array
The input to shift.
shift : int
The number of places by which elements are shifted (from numpy.roll).
Default: 1
axis : int
The axis along which elements are shifted (from numpy.rool).
sub : foat, default=0
The substitution value of the border
Returns
-------
array :
The shifted image.
"""
# Step 1: Cyclic permutation on im
im = numpy.roll(im, shift, axis=axis)
# Step 2: get image dimension
dim = len(im.shape)
# Step 3: modify the boundary with replacement value
if dim == 2: # if 2D image
if shift == 1 and axis == 0:
im[0, :] = sub
elif shift == -1 and axis == 0:
im[-1, :] = sub
elif shift == 1 and axis == 1:
im[:, 0] = sub
elif shift == -1 and axis == 1:
im[:, -1] = sub
elif dim == 3: # if 3D image
if shift >= 1 and axis == 0:
im[0:shift, :, :] = sub
elif shift <= -1 and axis == 0:
im[shift:, :, :] = sub
elif shift >= 1 and axis == 1:
im[:, 0:shift, :] = sub
elif shift <= -1 and axis == 1:
im[:, shift:, :] = sub
elif shift >= 1 and axis == 2:
im[:, :, 0:shift] = sub
elif shift <= -1 and axis == 2:
im[:, :, shift:] = sub
else:
print("spam.helpers.imageManipulation.singleShift: dim={}. Should be 2 or 3.".format(dim))
print("exit function.")
return -1
return im
[docs]
def multipleShifts(im, shifts, sub=0):
"""
This function call ``singleShift`` multiple times.
Parameters
----------
im : array
The input to shift.
shifts : [int, int, int]
Defines the number of shifts to apply in every axis.
.. code-block:: text
shift = [s_x, s_y, s_z] applies a shift of:
. s_x on axis 0
. s_y on axis 1
. s_z on axis 2
sub : float, default=0
The substitution value of the border
Returns
-------
array :
The shifted image.
"""
# loop over the n axis
for i in range(len(shifts)):
# if the value of the shift is not 0 on axis i
if shifts[i]:
# we call singleShift (only once)
im = singleShift(im, shift=shifts[i], axis=i)
return im
def _binarisation(im, threshold=0.0, boolean=False, op=">", mask=None):
"""
This function binarise an input image according to a given threshold
It has an option to apply a mask to the binarized image to ignore the
outside of the sample/specimen
Parameters
-----------
im: array
The image to binarise.
threshold : float
the input limit value for binarization
boolean : bool
Changes the output format and phase distribution (see output)
op : string, default='>'
defines the thresholding operation
mask : array, default=None
The mask of the input image: is 0 outside the boundary(specimen) and 1 inside
Returns
--------
phases : array
The repartition of phases resulting the binarisation.
For operator '>' it gives, if ``boolean=True``:
.. code-block:: text
0 - masked parts (where mask equals 0)
0 - below threshold
1 - above threshold
and if ``boolean=False``
.. code-block:: text
0 - masked parts (where mask equals 0)
1 - below threshold
2 - above threshold
"""
import operator
# Step 1: Get operator
operation = {
">": operator.gt,
"<": operator.lt,
">=": operator.ge,
"<=": operator.le,
"=": operator.eq,
}.get(op)
# Step 2: binarisation
phases = operation(im, threshold).astype("<u1")
# Step 3: rescaling if bool
if not boolean:
phases += 1
# Step 4: apply mask
if mask is not None:
phases = phases * mask
return phases
[docs]
def slicePadded(im, startStop, createMask=False, padValue=0, verbose=False):
"""
Extract slice from im, padded with zeros, which is always of the dimensions asked (given from startStop)
Parameters
----------
im : 3D numpy array
The image to be sliced
startStop : 6 component list of ints
This array contains:
[Zmin, Zmax, Ymin, Ymax, Xmin, Xmax]
createMask : bool, optional
If True, return a padded slice, which is False when the slice falls outside im
Default = False
Returns
-------
imSliced : 3D numpy array
The sliced image
mask : 3D numpy array of bools
The 3D mask, only returned if createMask is True
"""
startStop = numpy.array(startStop).astype(int)
assert startStop[1] > startStop[0], "spam.helpers.slicePadded(): Zmax should be bigger than Zmin"
assert startStop[3] > startStop[2], "spam.helpers.slicePadded(): Ymax should be bigger than Ymin"
assert startStop[5] > startStop[4], "spam.helpers.slicePadded(): Xmax should be bigger than Xmin"
imSliced = (
numpy.zeros(
(
startStop[1] - startStop[0],
startStop[3] - startStop[2],
startStop[5] - startStop[4],
),
dtype=im.dtype,
)
+ padValue
)
start = numpy.array([startStop[0], startStop[2], startStop[4]])
stop = numpy.array([startStop[1], startStop[3], startStop[5]])
startOffset = numpy.array([max(0, start[0]), max(0, start[1]), max(0, start[2])])
stopOffset = numpy.array(
[
min(im.shape[0], stop[0]),
min(im.shape[1], stop[1]),
min(im.shape[2], stop[2]),
]
)
startLocal = startOffset - start
stopLocal = startLocal + stopOffset - startOffset
# Check condition that we're asking for a slice of data wholly outside im
# This means either that the stop values are all smaller than 0
# OR the start are all bigger than im.shape
if numpy.any(stop < numpy.array([0, 0, 0])) or numpy.any(start >= numpy.array(im.shape)):
if verbose:
print("spam.helpers.slicePadded(): The extracted padded slice doesn't not touch the image!")
imSliced = imSliced.astype("<f4")
imSliced *= numpy.nan
if createMask:
return imSliced, numpy.zeros_like(imSliced, dtype=bool)
else:
imSliced[startLocal[0] : stopLocal[0], startLocal[1] : stopLocal[1], startLocal[2] : stopLocal[2],] = im[
startOffset[0] : stopOffset[0],
startOffset[1] : stopOffset[1],
startOffset[2] : stopOffset[2],
]
if createMask:
mask = numpy.zeros_like(imSliced, dtype=bool)
mask[
startLocal[0] : stopLocal[0],
startLocal[1] : stopLocal[1],
startLocal[2] : stopLocal[2],
] = 1
return imSliced, mask
return imSliced
[docs]
def splitImage(im, divisions, margin, verbose=True):
"""
Divides the image in zDiv x yDiv x xDiv blocks, each block is padded
with a margin.
Parameters
----------
im : 3D numpy array
The image to be splitted
divisions : 3-component list of ints
Desired number of blocks along Z, Y, X axes
margin : int
Overlapping margin between each block.
For applying a filter on subvolumes, it is recommended to use a margin of 1.5 times the filter diameter.
For labelled data it is recommended that the margin is at least 1.5 times bigger than the particles largest axis
verbose : bool
Print the parameters of the operations (number of blocks and margin)
Default = True
Returns
-------
Dictionary
Dictionary with keys labelled acoording to the position of the block along each axis (e.g., 000, 001, 002,...)
Each element (e.g., 001) within the dictionary carries the block origin and the resulting block, in that order
Note
----
This function should be used along `spam.helpers.imageManipulation.rebuildImage()`
"""
zDiv, yDiv, xDiv = divisions
# Check if the slices can be made
if zDiv >= im.shape[0]:
print("spam.helpers.imageManipulation.splitImage: Incorrect number of slices for axis z")
print("exit function.")
return -1
if yDiv >= im.shape[1]:
print("spam.helpers.imageManipulation.splitImage: Incorrect number of slices for axis y")
print("exit function.")
return -1
if xDiv >= im.shape[2]:
print("spam.helpers.imageManipulation.splitImage: Incorrect number of slices for axis x")
print("exit function.")
return -1
# Check that margin is not greater than the slice
if margin >= im.shape[0] / zDiv:
print("spam.helpers.imageManipulation.splitImage: Margin is too big for z axis")
print("exit function.")
return -1
if margin >= im.shape[1] / yDiv:
print("spam.helpers.imageManipulation.splitImage: Margin is too big for y axis")
print("exit function.")
return -1
if margin >= im.shape[2] / xDiv:
print("spam.helpers.imageManipulation.splitImage: Margin is too big for x axis")
print("exit function.")
return -1
# Print parameters if needed
if verbose:
print(
"spam.helpers.imageManipulation.splitImage: Working with margin of ",
margin,
". The total number of blocks is ",
zDiv * yDiv * xDiv,
)
# Pad initial image with zeros on the edge
imPad = numpy.pad(im, margin, mode="edge")
# Compute size of blocks
zSize = int(im.shape[0] / zDiv)
ySize = int(im.shape[1] / yDiv)
xSize = int(im.shape[2] / xDiv)
# Create return dictionary
output = {}
# Iterate through each block
for zBlock in range(zDiv):
for yBlock in range(yDiv):
for xBlock in range(xDiv):
# Get the origin of each block
blockOrigin = numpy.array([zBlock * zSize, yBlock * ySize, xBlock * xSize])
blockSize = [zSize, ySize, xSize]
# Check if the size needs to be changed to fit the image
if zBlock == zDiv - 1 and im.shape[0] % zDiv != 0:
blockSize[0] = zSize + (im.shape[0] % zDiv)
if yBlock == yDiv - 1 and im.shape[1] % yDiv != 0:
blockSize[1] = ySize + (im.shape[1] % yDiv)
if xBlock == xDiv - 1 and im.shape[2] % xDiv != 0:
blockSize[2] = xSize + (im.shape[2] % xDiv)
# Generate block with the margin on all sides
imBlock = crop(
imPad,
(
blockSize[0] + 2 * margin,
blockSize[1] + 2 * margin,
blockSize[2] + 2 * margin,
),
boxOrigin=(blockOrigin[0], blockOrigin[1], blockOrigin[2]),
)
# Save the results
output.update({str(zBlock) + str(yBlock) + str(xBlock): [blockOrigin, imBlock]})
# Save the margin
output.update({"margin": margin})
# Return
return output
[docs]
def rebuildImage(listBlocks, listCoordinates, margin, mode, keepLabels=False, verbose=True):
"""
Rebuilds splitted image from `spam.helpers.imageManipulation.splitImage()`.
Parameters
----------
listBlocks : list
List of the 3D blocks that will form the re-built the image.
Note: The order of listBlocks should be equivalent to the order of listCoordinates
listCoordinates : list
List of the origin coordinates of each block. (Usually taken from `spam.helpers.imageManipulation.splitImage()`)
Note: The order of listCoordinates should be equivalent to the order of listBlocks
margin : integer
Value of the margin used for the images. (Usually taken from `spam.helpers.imageManipulation.splitImage()`)
mode : string
'grey' : re-builds 3D greyscale arrays
'label' : re-builds 3D labelled arrays
keepLabels : bool
Do we need to want to keep the current labels from the blocks, or create a new one?
Default = False
verbose : bool
Print the evolution of the operation
Default = True
Returns
-------
imBuild : 3D numpy array
Re-built image without the margins
Note
----
This function should be used along with `spam.helpers.imageManipulation.splitImage()`
"""
# Checking if listBlocks and listCoordinates have the same length
if len(listBlocks) != len(listCoordinates):
print("spam.helpers.imageManipulation.splitImage: listBlocks and listCoordinates must have the same length")
return -1
# Transform listCoordinates into array
arrayCoord = numpy.asarray(listCoordinates)
# Checking if all the origin coordinates are different
_, counts = numpy.unique(arrayCoord, axis=0, return_counts=True)
if len(counts) != len(arrayCoord):
print("spam.helpers.imageManipulation.splitImage: coordinates in listCoordinates must be all different")
return -1
if verbose:
# Create progressbar
widgets = [
progressbar.FormatLabel(""),
" ",
progressbar.Bar(),
" ",
progressbar.AdaptiveETA(),
]
pbar = progressbar.ProgressBar(widgets=widgets, maxval=len(listCoordinates))
pbar.start()
finishedBlocks = 0
if mode == "grey":
# Shape of the block opposite to the origine
shapeLast = listBlocks[numpy.argmax(numpy.sum(arrayCoord, axis=1))].shape
# Tentative size of the final image
zSize = numpy.amax(arrayCoord[:, 0]) + shapeLast[0] - 2 * margin
ySize = numpy.amax(arrayCoord[:, 1]) + shapeLast[1] - 2 * margin
xSize = numpy.amax(arrayCoord[:, 2]) + shapeLast[2] - 2 * margin
# Initialising rebuild image
imBuild = numpy.zeros((zSize, ySize, xSize))
# Loop on the length to lists, so to replace zeros in imBuild with the actual values at the right position
for i in range(len(listCoordinates)):
origin = listCoordinates[i]
# GP 01/03/2022: Changing to ease the charge on the memory
blockPad = listBlocks.pop(0)
# blockPad = listBlocks[i]
if margin == 0:
imBuild[
origin[0] : origin[0] + blockPad.shape[0] - 2 * margin,
origin[1] : origin[1] + blockPad.shape[1] - 2 * margin,
origin[2] : origin[2] + blockPad.shape[2] - 2 * margin,
] = blockPad
else:
imBuild[
origin[0] : origin[0] + blockPad.shape[0] - 2 * margin,
origin[1] : origin[1] + blockPad.shape[1] - 2 * margin,
origin[2] : origin[2] + blockPad.shape[2] - 2 * margin,
] = blockPad[margin:-margin, margin:-margin, margin:-margin]
if verbose:
# Update the progressbar
finishedBlocks += 1
widgets[0] = progressbar.FormatLabel("{}/{} ".format(finishedBlocks, len(listCoordinates)))
pbar.update(finishedBlocks)
print("\n\t")
return imBuild
if mode == "label":
# Shape of the block opposite to the origine
shapeLast = listBlocks[numpy.argmax(numpy.sum(arrayCoord, axis=1))].shape
# Size of the final image
zSize = numpy.amax(arrayCoord[:, 0]) + shapeLast[0] - 2 * margin
ySize = numpy.amax(arrayCoord[:, 1]) + shapeLast[1] - 2 * margin
xSize = numpy.amax(arrayCoord[:, 2]) + shapeLast[2] - 2 * margin
# Initialising rebuild image
imBuild = numpy.zeros((zSize, ySize, xSize))
# Loop over the lists
for i in range(len(listCoordinates)):
# Get the origin of the block
origin = listCoordinates[i]
# Get the block
# GP 01/03/2022: Changing to ease the charge on the memory
# block = listBlocks[i]
block = listBlocks.pop(0)
# Compute the bounding boxes
boundingBoxes = spam.label.boundingBoxes(block)
# Compute centre of mass
centresOfMass = spam.label.centresOfMass(block, boundingBoxes=boundingBoxes)
# List for classifying the labels
inside = []
outside = []
partial = []
# Check if each label is inside the true block - i.e., it is inside the block without the margin
for j in range(1, len(boundingBoxes), 1):
# Get the box
box = boundingBoxes[j]
# Check if the origin is inside the true block
checkOrigin = margin < box[0] < block.shape[0] - margin and margin < box[2] < block.shape[1] - margin and margin < box[4] < block.shape[2] - margin
# Check if the coordinate opposite to the origin is inside the true block
checkOpp = margin < box[1] < block.shape[0] - margin and margin < box[3] < block.shape[1] - margin and margin < box[5] < block.shape[2] - margin
if checkOrigin and checkOpp:
# Both points are inside
inside.append(j)
elif checkOrigin or checkOpp:
# Only one is inside
partial.append(j)
else:
# Both are outside
outside.append(j)
# Create true block array, keep only particles fully inside
trueBlock = spam.label.removeLabels(block, partial + outside)
if not keepLabels:
# Check that we have labels inside the block
if len(numpy.unique(trueBlock)) > 1:
# Make labels consecutive
trueBlock = spam.label.makeLabelsSequential(trueBlock)
trueBlock = numpy.where(trueBlock != 0, trueBlock + numpy.max(imBuild), trueBlock)
# IV (04-03-21): Info needed to avoid chopping grains that would be considered as outside particles
imBuildSubSet = imBuild[
origin[0] : origin[0] + trueBlock.shape[0] - 2 * margin,
origin[1] : origin[1] + trueBlock.shape[1] - 2 * margin,
origin[2] : origin[2] + trueBlock.shape[2] - 2 * margin,
]
# Add the true block preserving previous info
imBuild[origin[0] : origin[0] + trueBlock.shape[0] - 2 * margin, origin[1] : origin[1] + trueBlock.shape[1] - 2 * margin, origin[2] : origin[2] + trueBlock.shape[2] - 2 * margin,] = (
trueBlock[margin:-margin, margin:-margin, margin:-margin] + imBuildSubSet
)
# Get current maximum label
tempMax = numpy.max(imBuild)
# Label counter
labCounter = 1
# Solve the labels inside partial list
if len(partial) > 0:
# Iterate trough each label
for label in partial:
# Get the bounding box
box = boundingBoxes[label]
# Get subset
labelSubset = spam.label.getLabel(
block,
label,
boundingBoxes=boundingBoxes,
centresOfMass=centresOfMass,
)
labelSubvol = labelSubset["subvol"].astype(int)
labelSubvol = numpy.where(labelSubvol != 0, labCounter + tempMax, labelSubvol)
# IV+GP implementation - Just put the box there, overwrite whatever you find there
# Get the subset from imBuild
imBuildSubset = imBuild[
origin[0] + box[0] - margin : origin[0] + box[1] + 1 - margin,
origin[1] + box[2] - margin : origin[1] + box[3] + 1 - margin,
origin[2] + box[4] - margin : origin[2] + box[5] + 1 - margin,
]
# Change
imBuildSubset = numpy.where(labelSubvol != 0, labelSubvol, imBuildSubset)
# Put back
imBuild[
origin[0] + box[0] - margin : origin[0] + box[1] + 1 - margin,
origin[1] + box[2] - margin : origin[1] + box[3] + 1 - margin,
origin[2] + box[4] - margin : origin[2] + box[5] + 1 - margin,
] = imBuildSubset
# Update label counter
labCounter += 1
if verbose:
# Update the progressbar
finishedBlocks += 1
widgets[0] = progressbar.FormatLabel("{}/{} ".format(finishedBlocks, len(listCoordinates)))
pbar.update(finishedBlocks)
if not keepLabels:
imBuild = spam.label.makeLabelsSequential(imBuild)
print("\n\t")
return imBuild
else:
# The mode is not correct
print("spam.helpers.imageManipulation.splitImage(): Incorrect mode, check your input")
return -1
[docs]
def checkerBoard(im1, im2, n=5, inv=False, rescale=True):
"""
This function generates a "checkerboard" mix of two 2D images of the same size.
This is useful to see if they have been properly aligned, especially if the two images are
quantitatively different (i.e., one is a neutron tomography and the other is an x-ray tomography).
Parameters
----------
im1 : 2D numpy array
This is the first image
im2 : 2D/3D numpy array
This is the second image, should be same shape as first image
n : integer, optional
The number of divisions of the checkerboard to aim for.
Default = 5
inv : bool, optional
Whether im2 should be -im2 in the checkerboard.
Default = False
rescale : bool, optional
Whether greylevels should be rescaled with spam.helpers.rescale.
Default = True
Returns
-------
im1G : checkerBoard mix of im1 and im2
"""
if inv:
c = -1.0
else:
c = 1.0
if rescale:
import spam.helpers
im1 = spam.helpers.rescale(im1)
im2 = spam.helpers.rescale(im2)
# 2D version
if len(im1.shape) == 2:
# initialize
im1G = im1.copy()
# get number of pixel / square based on min size
nP = int(min(im1.shape) / n)
for x in range(im1.shape[0]):
for y in range(im1.shape[1]):
if int((x % (2 * nP)) / nP) + int((y % (2 * nP)) / nP) - 1:
im1G[x, y] = c * im2[x, y]
else:
print("checkerBoard works only with dim2 images")
return 0
return im1G
# private functions
# def _mask2D(im, erosion=False, structure=None, ):
# """
# get contour of 2D image.
# """
# import cv2
# from scipy import ndimage
# step 2: convert into uint8 if not the case
# if im.dtype != 'uint8':
# actually it rescales the image but it doesn't really amtter
# im = rescaleToInteger(im, nBytes=1)
# Step 3: ...
# blur = cv2.GaussianBlur(im, (5, 5), 0)
# _, thresh = cv2.threshold(blur, 0, 255, cv2.THRESH_BINARY + cv2.THRESH_OTSU)
# _, contours, _ = cv2.findContours(thresh, cv2.RETR_TREE, cv2.CHAIN_APPROX_NONE)
# largest = 0
# biggest = []
# for contour in contours:
# area = cv2.contourArea(contour)
# if largest < area:
# largest = area
# biggest = contour
# mask = numpy.zeros(im.shape, dtype='<u1')
# cv2.drawContours(mask, [biggest], 0, 1, -1)
# Step 4: apply erosion of the mask (which corresponds to an erosion of the specimen)
# if erosion:
# mask = ndimage.morphology.binary_erosion(
# mask, structure=structure).astype(mask.dtype)
# return mask