Skip to content
GitLab
Commit 7a2c07da authored by Edward Andò's avatar Edward Andò
Browse files

Merge branch 'splitImage'

parents 2480c2de da8f5cf4
Pipeline #64837 passed with stages
in 26 minutes and 4 seconds
Hide whitespace changes
Inline Side-by-side
tools/helpers/imageManipulation.py
View file @ 7a2c07da
......@@ -20,7 +20,7 @@ from __future__ import print_function
import numpy
import tifffile
import spam.label
def stackToArray(prefix, sufix='.tif', stack=range(10), digits='05d', erosion=False, verbose=False):
"""
......@@ -521,6 +521,298 @@ def slicePadded(im, startStop, createMask=False, padValue=0):
return imSliced
def splitImage(im, divisions, margin, verbose=False):
"""
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)
Returns
-------
output: 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==True:
print('spam.helpers.imageManipulation.splitImage: Working with margin of ', margin)
print('spam.helpers.imageManipulation.splitImage: 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
def rebuildImage(listBlocks, listCoordinates, margin, mode, keepLabels=False):
"""
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
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 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]
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]
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
block = listBlocks[i]
# 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 not keepLabels:
imBuild = spam.label.makeLabelsSequential(imBuild)
return imBuild
else:
# The mode is not correct
print('spam.helpers.imageManipulation.splitImage(): Incorrect mode, check your input')
return -1
# private functions
#def _mask2D(im, erosion=False, structure=None, ):
#"""
......
tools/tests/test_imageManipulation.py
View file @ 7a2c07da
......@@ -8,7 +8,8 @@ import numpy
import spam.helpers
import spam.datasets
import spam.kalisphera
import spam.label
class testAll(unittest.TestCase):
......@@ -158,5 +159,147 @@ class testAll(unittest.TestCase):
self.assertEqual(numpy.isfinite(imSliced).sum(), 0)
self.assertEqual(mask.sum(), 0)
def test_splitRebuildImage(self):
# Importing snow data as im
im = spam.datasets.loadSnow()
# case 1: can be the slices be made?
res1a = spam.helpers.imageManipulation.splitImage(im, (im.shape[0], 1, 1), 0)
self.assertEqual(res1a, -1)
res1b = spam.helpers.imageManipulation.splitImage(im, (1, im.shape[1], 1), 0)
self.assertEqual(res1b, -1)
res1c = spam.helpers.imageManipulation.splitImage(im, (1, 1, im.shape[2]), 0)
self.assertEqual(res1c, -1)
# case 2: margin not greater than the slice
res2a = spam.helpers.imageManipulation.splitImage(im, (1, 1, 1), im.shape[0])
self.assertEqual(res2a, -1)
res2b = spam.helpers.imageManipulation.splitImage(im, (1, 1, 1), im.shape[1])
self.assertEqual(res2b, -1)
res2c = spam.helpers.imageManipulation.splitImage(im, (1, 1, 1), im.shape[2])
self.assertEqual(res2c, -1)
# case 3 : check if coordinates are repeating themselves
# Split
split = spam.helpers.imageManipulation.splitImage(im, (3, 3, 3), 10)
# Extract
listCoordinates = []
listBlocks = []
for key in split:
if key != 'margin':
coord, block = split[key]
listCoordinates.append(coord)
listBlocks.append(block)
# Get the margin
margin = split['margin']
# Appending an existing vector to the list of coordinates
listCoordinatesCopy = listCoordinates[:]
listCoordinatesCopy.append(listCoordinates[0])
res3 = spam.helpers.imageManipulation.rebuildImage(listBlocks,
listCoordinatesCopy,
margin,
mode='grey')
self.assertEqual(res3, -1)
# case 4: check if listCoordinates and listBlocks have same length
# Copy of listCoordinates
listCoordinatesCopy = listCoordinates[:]
listCoordinatesCopy.pop()
res4 = spam.helpers.imageManipulation.rebuildImage(listBlocks,
listCoordinatesCopy,
margin,
mode='grey')
self.assertEqual(res4, -1)
# case 5: MODE = Grey: is the result of the rebuilding equal to the original image?
rebuild = spam.helpers.imageManipulation.rebuildImage(listBlocks,
listCoordinates,
margin,
mode='grey')
diff = im - rebuild
res5 = numpy.sum(diff)
self.assertEqual(res5, 0)
# case 6: MODE = Label:
# as per Kalisphera example
# m/pixel
pixelSize = 40.e-6
# The standard deviation of the image blurring to be applied first
blurSTD = 0
# The standard deviation of the random noise to be added to each voxel
noiseSTD = 0
boxSizeDEM, centres, radii = spam.datasets.loadDEMboxsizeCentreRadius()
# get maximum radius to pad our image (periodic boundaries...)
rMax = numpy.amax(radii)
boxSize = boxSizeDEM + 3 * rMax
# move the positions to the new center of the image
centres[:, :] = centres[:, :] + 1.5 * rMax
# turn the mm measures into pixels
boxSize = int(numpy.ceil(numpy.max(boxSize[:]) / pixelSize))
#centres = centres / pixelSize
#radii = radii / pixelSize
Box = spam.kalisphera.makeBlurryNoisySphere((boxSize, boxSize, boxSize),
centres/pixelSize,
radii/pixelSize,
blur=0,
noise=0,
flatten=True,
background=0.0,
foreground=1.0)
#Box = numpy.zeros((boxSize, boxSize, boxSize), dtype="<f8")
#spam.kalisphera.makeSphere(Box, centres, radii)
# Create Gold-standard segmentation
imLab = spam.label.watershed((Box>0.5).astype(int))
# Split Greylevels
res = spam.helpers.imageManipulation.splitImage(Box, (2, 1, 1), 15)
# Create list
listBlocks = []
listCoordinates = []
# Label the blocks
for key in res:
if key != 'margin':
coord, block = res[key]
# Run the watershed
imLabBlock = spam.label.watershed((block>0.5).astype(int))
# Add them to the list
listBlocks.append(imLabBlock)
listCoordinates.append(coord)
# Get the margin
margin = res['margin']
# Rebuild the image
imTest = spam.helpers.imageManipulation.rebuildImage(listBlocks,
listCoordinates,
margin,
mode = 'label')
# Get the volumes & sort
volOr = spam.label.volumes(imLab)
volMod = spam.label.volumes(imTest)
volOr.sort()
volMod.sort()
# Get the centres of mass & sort
comOr = spam.label.centresOfMass(imLab)
comMod = spam.label.centresOfMass(imTest)
comOr = numpy.sort(comOr, axis=0)
comMod = numpy.sort(comMod, axis=0)
# Case 6A: Check the volumes
self.assertEqual(numpy.sum(volOr - volMod), 0)
# Case 6B: Check the centres of mass
self.assertEqual(numpy.sum(comOr - comMod), 0)
# Case 6C: Check the size of the final image
self.assertTrue(imLab.shape == imTest.shape)
# Case 7: Check for a weird mode
res7 = spam.helpers.imageManipulation.rebuildImage(listBlocks,
listCoordinates,
margin,
mode = 'notMode')
self.assertEqual(res7, -1)
if __name__ == '__main__':
unittest.main()
Supports Markdown
0% or .
You are about to add 0 people to the discussion. Proceed with caution.
Finish editing this message first!
Please register or to comment