Commit f40784da authored by Edward Andò's avatar Edward Andò
Browse files

nice function, probably ready for merging back to master

parent 252fd0ce
......@@ -522,7 +522,7 @@ def slicePadded(im, startStop, createMask=False, padValue=0):
return imSliced
def splitImage(im, zDiv, yDiv, xDiv, margin, verbose=False):
def splitImage(im, divisions, margin, verbose=False):
"""
Divides the image in zDiv x yDiv x xDiv blocks, each block is padded
with a margin.
......@@ -531,21 +531,15 @@ def splitImage(im, zDiv, yDiv, xDiv, margin, verbose=False):
----------
im : 3D numpy array
The image to be splitted
zDiv : int
Number of block divisions along Z axis
yDiv : int
Number of block divisions along Y axis
xDiv : int
Number of block partitions along X axis
divisions : 3-component list of ints
Desired number of blocks along Z, Y, X axes
margin : int
Overlapping margin between each block.
It is recommended that the margin is at least 1.5 times
bigger than the particles largest axis
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)
......@@ -554,13 +548,14 @@ def splitImage(im, zDiv, yDiv, xDiv, margin, verbose=False):
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')
......@@ -574,7 +569,7 @@ def splitImage(im, zDiv, yDiv, xDiv, margin, verbose=False):
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')
......@@ -592,23 +587,23 @@ def splitImage(im, zDiv, yDiv, xDiv, margin, verbose=False):
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
# Get the origin of each block
blockOrigin = numpy.array([zBlock*zSize,
yBlock*ySize,
xBlock*xSize])
......@@ -620,11 +615,11 @@ def splitImage(im, zDiv, yDiv, xDiv, margin, verbose=False):
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]))
(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).zfill(2) + str(yBlock).zfill(2) + str(xBlock).zfill(2): [blockOrigin, imBlock]})
# Save the margin
......@@ -634,27 +629,27 @@ def splitImage(im, zDiv, yDiv, xDiv, margin, verbose=False):
return output
def rebuildImage(listBlocks, listCoordinates, margin, mode, keepLabels = False):
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
......@@ -663,45 +658,45 @@ def rebuildImage(listBlocks, listCoordinates, margin, mode, keepLabels = False):
-------
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,
......@@ -710,14 +705,14 @@ def rebuildImage(listBlocks, listCoordinates, margin, mode, keepLabels = False):
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
......@@ -735,12 +730,12 @@ def rebuildImage(listBlocks, listCoordinates, margin, mode, keepLabels = False):
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
......@@ -758,7 +753,7 @@ def rebuildImage(listBlocks, listCoordinates, margin, mode, keepLabels = False):
else:
# Both are outside
outside.append(j)
# Create true block array, keep only particles fully inside
trueBlock = spam.label.removeLabels(block, partial+outside)
......@@ -769,7 +764,7 @@ def rebuildImage(listBlocks, listCoordinates, margin, mode, keepLabels = False):
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
# 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]
......@@ -778,7 +773,7 @@ def rebuildImage(listBlocks, listCoordinates, margin, mode, keepLabels = False):
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
......@@ -806,16 +801,16 @@ def rebuildImage(listBlocks, listCoordinates, margin, mode, keepLabels = False):
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')
print('spam.helpers.imageManipulation.splitImage(): Incorrect mode, check your input')
return -1
# private functions
......
......@@ -158,30 +158,30 @@ class testAll(unittest.TestCase):
self.assertEqual(list(imSliced.shape), [20, 20, 20])
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)
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)
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)
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])
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])
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])
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)
split = spam.helpers.imageManipulation.splitImage(im, (3, 3, 3), 10)
# Extract
listCoordinates = []
listBlocks = []
......@@ -195,35 +195,35 @@ class testAll(unittest.TestCase):
# 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')
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')
mode='grey')
self.assertEqual(res4, -1)
# case 5: MODE = Grey: is the reuslt of the rebuilding equal to the original image?
# 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')
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
......@@ -238,15 +238,24 @@ class testAll(unittest.TestCase):
# 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 = numpy.zeros((boxSize, boxSize, boxSize), dtype="<f8")
spam.kalisphera.makeSphere(Box, centres, radii)
#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)
res = spam.helpers.imageManipulation.splitImage(Box, (2, 1, 1), 15)
# Create list
listBlocks = []
listCoordinates = []
......@@ -263,9 +272,9 @@ class testAll(unittest.TestCase):
margin = res['margin']
# Rebuild the image
imTest = spam.helpers.imageManipulation.rebuildImage(listBlocks,
listCoordinates,
margin,
mode = 'label')
listCoordinates,
margin,
mode = 'label')
# Get the volumes & sort
volOr = spam.label.volumes(imLab)
volMod = spam.label.volumes(imTest)
......@@ -274,24 +283,23 @@ class testAll(unittest.TestCase):
# 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)
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')
listCoordinates,
margin,
mode = 'notMode')
self.assertEqual(res7, -1)
if __name__ == '__main__':
unittest.main()
Markdown is supported
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