Skip to content

GitLab

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

[skip-ci] spam-pixelSearch with initial guess again

parent 0876e502
Pipeline #54344 skipped
Hide whitespace changes
Inline Side-by-side
scripts/spam-pixelSearch
View file @ b5989395
......@@ -82,6 +82,7 @@ if args.LAB1 is not None:
boundingBoxes = spam.label.boundingBoxes(lab1)
nodePositions = spam.label.centresOfMass(lab1, boundingBoxes=boundingBoxes)
numberOfNodes = nodePositions.shape[0]
im1mask = None
### Otherwise we are in node spacing and half-window size mode
else:
......@@ -108,7 +109,6 @@ PhiField = numpy.zeros((nodePositions.shape[0], 4, 4))
for node in range(nodePositions.shape[0]):
PhiField[node] = numpy.eye(4)
if args.PHIFILE is not None:
PhiFromFile = spam.helpers.readCorrelationTSV(args.PHIFILE.name, fieldBinRatio=args.PHIFILE_BIN_RATIO)
......@@ -241,8 +241,8 @@ Returns
# WARNING ------------------------- VVVVVVVVVVV--- gets easily overwritten, pass a .copy()!
def getImagettes(nodePosition, Phi, searchRange, boundingBox, im1, im2, im1mask, minMaskCoverage, greyThreshold):
returnStatus = 2
returnStatus = 1
imagette1mask = None
initialDisplacement = Phi[0:3, 3].astype(int)
# point in im2 that we are searching around
......@@ -251,91 +251,100 @@ def getImagettes(nodePosition, Phi, searchRange, boundingBox, im1, im2, im1mask,
#searchCentre = halfWindowSize - searchRange[0::2]
#Add initial displacement to searchRange
searchRange[0::2] += initialDisplacement
searchRange[1::2] += initialDisplacement
#searchRange[0::2] += initialDisplacement
#searchRange[1::2] += initialDisplacement
# 2020-09-25 OS and EA: Prepare startStop array for imagette 1 to be extracted with new slicePadded
PhiNoDisp = Phi
PhiNoDisp[0:3,-1] = 0.0
# Catch bad bounding boxes:
if numpy.all((boundingBox[1::2] - boundingBox[0::2]) > [0, 0, 0]) and numpy.all(boundingBox[1::2] - boundingBox[0::2] != [0, 0, 0]):
# If F is not the identity, create a pad to be able to apply F to imagette 1
if numpy.allclose(PhiNoDisp, numpy.eye(4)):
applyPhiPad = 0
else:
# 2020-10-06 OS and EA: Add a pad to each dimension of 25% of max(halfWindowSize) to allow space to apply F (no displacement) to imagette1
applyPhiPad = int(0.5*numpy.ceil(max(boundingBox[1::2]-boundingBox[0::2])))
# 2020-09-25 OS and EA: Prepare startStop array for imagette 1 to be extracted with new slicePadded
PhiNoDisp = Phi.copy()
PhiNoDisp[0:3,-1] = 0.0
startStopIm1 = [int(boundingBox[0] - applyPhiPad), int(boundingBox[1] + applyPhiPad),
int(boundingBox[2] - applyPhiPad), int(boundingBox[3] + applyPhiPad),
int(boundingBox[4] - applyPhiPad), int(boundingBox[5] + applyPhiPad)]
# If F is not the identity, create a pad to be able to apply F to imagette 1
if numpy.allclose(PhiNoDisp, numpy.eye(4)):
applyPhiPad = 0
else:
# 2020-10-06 OS and EA: Add a pad to each dimension of 25% of max(halfWindowSize) to allow space to apply F (no displacement) to imagette1
applyPhiPad = int(0.5*numpy.ceil(max(boundingBox[1::2]-boundingBox[0::2])))
# In either case, extract imagette1, now guaranteed to be the right size
imagette1padded = spam.helpers.slicePadded(im1, startStopIm1)
startStopIm1 = [int(boundingBox[0] - applyPhiPad), int(boundingBox[1] + applyPhiPad + 1),
int(boundingBox[2] - applyPhiPad), int(boundingBox[3] + applyPhiPad + 1),
int(boundingBox[4] - applyPhiPad), int(boundingBox[5] + applyPhiPad + 1)]
# If F is not the identity, apply F and undo crop
if numpy.allclose(PhiNoDisp, numpy.eye(4)):
# In this case there is is no padding (despite the name) and we can just keep going
imagette1def = imagette1padded
else:
# apply F to imagette 1 padded
imagette1paddedDef = spam.DIC.applyPhi(imagette1padded, PhiNoDisp)
# undo padding
imagette1def = imagette1paddedDef[applyPhiPad:-applyPhiPad,
applyPhiPad:-applyPhiPad,
applyPhiPad:-applyPhiPad]
### Check mask
if im1mask is None:
# no mask1 --> always pas this test (e.g., labelled image)
maskVolumeCondition = True
imagette1mask = None
else:
imagette1mask = spam.helpers.slicePadded(im1mask, boundingBox) != 0
maskVolumeCondition = imagette1mask.mean() > minMaskCoverage
# Make sure imagette is not 0-dimensional in any dimension
# Check minMaskVolume
if numpy.all(numpy.array(imagette1def.shape) > 0):
# ------------ Grey threshold low --------------- and -------------- Grey threshold high -----------
if numpy.nanmean(imagette1def) > greyThreshold[0] and numpy.nanmean(imagette1def) < greyThreshold[1]:
if maskVolumeCondition:
# Slice for image 2
## 2020-09-25 OS and EA: Prepare startStop array for imagette 2 to be extracted with new slicePadded
## Extract it...
startStopIm2 = [int(boundingBox[0] + searchRange[0] ),
int(boundingBox[1] + searchRange[1] + 1),
int(boundingBox[2] + searchRange[2] ),
int(boundingBox[3] + searchRange[3] + 1),
int(boundingBox[4] + searchRange[4] ),
int(boundingBox[5] + searchRange[5] + 1)]
imagette2 = spam.helpers.slicePadded(im2, startStopIm2)
# Failed minMaskVolume condition
# In either case, extract imagette1, now guaranteed to be the right size
imagette1padded = spam.helpers.slicePadded(im1, startStopIm1)
# If F is not the identity, apply F and undo crop
if numpy.allclose(PhiNoDisp, numpy.eye(4)):
# In this case there is is no padding (despite the name) and we can just keep going
imagette1def = imagette1padded
else:
# apply F to imagette 1 padded
imagette1paddedDef = spam.DIC.applyPhi(imagette1padded, PhiNoDisp)
# undo padding
imagette1def = imagette1paddedDef[applyPhiPad:-applyPhiPad,
applyPhiPad:-applyPhiPad,
applyPhiPad:-applyPhiPad]
### Check mask
if im1mask is None:
# no mask1 --> always pas this test (e.g., labelled image)
maskVolumeCondition = True
imagette1mask = None
else:
imagette1mask = spam.helpers.slicePadded(im1mask, boundingBox + numpy.array([0, 1, 0, 1, 0, 1])) != 0
maskVolumeCondition = imagette1mask.mean() > minMaskCoverage
# Make sure imagette is not 0-dimensional in any dimension
# Check minMaskVolume
if numpy.all(numpy.array(imagette1def.shape) > 0):
# ------------ Grey threshold low --------------- and -------------- Grey threshold high -----------
if numpy.nanmean(imagette1def) > greyThreshold[0] and numpy.nanmean(imagette1def) < greyThreshold[1]:
if maskVolumeCondition:
# Slice for image 2
## 2020-09-25 OS and EA: Prepare startStop array for imagette 2 to be extracted with new slicePadded
## Extract it...
startStopIm2 = [int(boundingBox[0] + initialDisplacement[0] + searchRange[0] ),
int(boundingBox[1] + initialDisplacement[0] + searchRange[1] + 1),
int(boundingBox[2] + initialDisplacement[1] + searchRange[2] ),
int(boundingBox[3] + initialDisplacement[1] + searchRange[3] + 1),
int(boundingBox[4] + initialDisplacement[2] + searchRange[4] ),
int(boundingBox[5] + initialDisplacement[2] + searchRange[5] + 1)]
imagette2 = spam.helpers.slicePadded(im2, startStopIm2)
# Failed minMaskVolume condition
else:
returnStatus = -5
imagette1def = None
imagette2 = None
# Failed greylevel condition
else:
returnStatus = -5
imagette1def = None
imagette2 = None
# Failed greylevel condition
# Failed 0-dimensional imagette test
else:
returnStatus = -5
imagette1def = None
imagette2 = None
# Failed 0-dimensional imagette test
# Failed bounding box test
else:
returnStatus = -5
returnStatus = -7
imagette1def = None
imagette2 = None
return {'imagette1': imagette1def,
'imagette1mask': imagette1mask,
'imagette2': imagette2,
'returnStatus': returnStatus,
'initialDisplacement': initialDisplacement,
'searchRange': searchRange,
'searchCentre': searchCentre}
return {'imagette1': imagette1def,
'imagette1mask': imagette1mask,
'imagette2': imagette2,
'returnStatus': returnStatus,
'initialDisplacement': initialDisplacement,
'searchRange': searchRange,
'searchCentre': searchCentre}
mpi = False
......@@ -358,10 +367,12 @@ else:
# Create pixelSearchCC vector
pixelSearchCC = numpy.zeros((numberOfNodes))
pixelSearchCC = numpy.zeros((numberOfNodes), dtype=float)
returnStatus = numpy.ones((numberOfNodes), dtype=int)
# Add nodes to a queue -- mostly useful for MPI
q = queue.Queue()
#for node in range(3):
for node in range(numberOfNodes):
q.put(node)
finishedNodes = 0
......@@ -380,9 +391,9 @@ while finishedNodes != numberOfNodes:
# Get the next node off the queue
nodeNumber = q.get()
imagetteReturns = getImagettes(nodePositions[nodeNumber], PhiField[nodeNumber], searchRange.copy(), boundingBoxes[nodeNumber], im1, im2, im1mask, args.MASK_COVERAGE, greyThreshold)
imagetteReturns = getImagettes(nodePositions[nodeNumber], PhiField[nodeNumber].copy(), searchRange.copy(), boundingBoxes[nodeNumber], im1, im2, im1mask, args.MASK_COVERAGE, greyThreshold)
if imagetteReturns['returnStatus'] == 2:
if imagetteReturns['returnStatus'] == 1:
if mpi:
# build message for pixel search worker
m = {'nodeNumber': nodeNumber,
......@@ -391,7 +402,7 @@ while finishedNodes != numberOfNodes:
'im1mask': imagetteReturns['imagette1mask'],
'searchRange': imagetteReturns['searchRange'],
'searchCentre': imagetteReturns['searchCentre'],
'initialDisplacement': imagetteReturns['initialDisplacement']
'initialDisplacement': imagetteReturns['initialDisplacement'],
}
# print "\tBoss: sending node {} to worker {}".format( nodeNumber, worker )
......@@ -416,7 +427,8 @@ while finishedNodes != numberOfNodes:
pixelSearchCC[nodeNumber] = 0.0
finishedNodes += 1
PhiField[nodeNumber, 0:3, 0:3] = numpy.eye(3)
PhiField[nodeNumber, 0:3, 3] = numpy.nan
PhiField[nodeNumber, 0:3, 3] = numpy.nan
returnStatus[nodeNumber] = imagetteReturns['returnStatus']
# Otherwise spend time looking waiting for replies from workers
elif mpi:
......@@ -438,7 +450,7 @@ while finishedNodes != numberOfNodes:
finishedNodes += 1
writeReturns = False
# set translation for this node
PhiField[nodeNumber, 0:3, 3] = numpy.array(returns['t'])
PhiField[nodeNumber, 0:3, 3] += numpy.array(returns['t'])
pixelSearchCC[nodeNumber] = returns['cc']
widgets[0] = progressbar.FormatLabel(" CC={:0>7.5f} ".format(pixelSearchCC[nodeNumber]))
......@@ -454,13 +466,13 @@ if args.TSV:
# 3 node positions,
# F[0:3,0:2]
# Pixel-search CC
TSVheader = "NodeNumber\tZpos\tYpos\tXpos\tFzz\tFzy\tFzx\tZdisp\tFyz\tFyy\tFyx\tYdisp\tFxz\tFxy\tFxx\tXdisp\tpixelSearchCC"
TSVheader = "NodeNumber\tZpos\tYpos\tXpos\tFzz\tFzy\tFzx\tZdisp\tFyz\tFyy\tFyx\tYdisp\tFxz\tFxy\tFxx\tXdisp\tpixelSearchCC\treturnStatus"
outMatrix = numpy.array([numpy.array(range(nodePositions.shape[0])),
nodePositions[:, 0], nodePositions[:, 1], nodePositions[:, 2],
PhiField[:, 0, 0], PhiField[:, 0, 1], PhiField[:, 0, 2], PhiField[:, 0, 3],
PhiField[:, 1, 0], PhiField[:, 1, 1], PhiField[:, 1, 2], PhiField[:, 1, 3],
PhiField[:, 2, 0], PhiField[:, 2, 1], PhiField[:, 2, 2], PhiField[:, 2, 3],
pixelSearchCC]).T
pixelSearchCC, returnStatus]).T
numpy.savetxt(args.OUT_DIR+"/"+args.PREFIX+".tsv",
outMatrix,
......@@ -500,7 +512,8 @@ if args.TIFF:
#tifffile.imsave(args.OUT_DIR+"/"+args.PREFIX+"-Zdisp-regsub.tif", PhiFieldMinusRigid[:, 0, -1].astype('<f4').reshape(nodesDim))
tifffile.imsave(args.OUT_DIR+"/"+args.PREFIX+"-Ydisp.tif", PhiField[:, 1, -1].astype('<f4').reshape(nodesDim))
tifffile.imsave(args.OUT_DIR+"/"+args.PREFIX+"-Xdisp.tif", PhiField[:, 2, -1].astype('<f4').reshape(nodesDim))
tifffile.imsave(args.OUT_DIR+"/"+args.PREFIX+"-pixelSearchCC.tif", pixelSearchCC.astype('<f4').reshape(nodesDim))
tifffile.imsave(args.OUT_DIR+"/"+args.PREFIX+"-CC.tif", pixelSearchCC.astype('<f4').reshape(nodesDim))
tifffile.imsave(args.OUT_DIR+"/"+args.PREFIX+"-returnStatus.tif", returnStatus.astype('<f4').reshape(nodesDim))
#if args.REGSUB:
#tifffile.imsave(args.OUT_DIR+"/"+args.PREFIX+"-Ydisp-regsub.tif", PhiFieldMinusRigid[:, 1, -1].astype('<f4').reshape(nodesDim))
#tifffile.imsave(args.OUT_DIR+"/"+args.PREFIX+"-Xdisp-regsub.tif", PhiFieldMinusRigid[:, 2, -1].astype('<f4').reshape(nodesDim))
......
tools/DIC/correlate.py
View file @ b5989395
......@@ -821,7 +821,7 @@ def pixelSearch(imagette1, imagette2, imagette1mask=None, searchCentre=None, sea
else: searchCentre = numpy.array(searchCentre, dtype='<f4')
if imagette1mask is not None:
assert imagette1.shape == imagette1mask.shape, "spam.DIC.correlate.pixelSearch: imagette1mask should have same size as imagette1 {} {} ".format( imagette1.shape, imagette1mask.shape)
assert imagette1.shape == imagette1mask.shape, "spam.DIC.correlate.pixelSearch: imagette1mask ({}) should have same size as imagette1 ({})".format( imagette1mask.shape, imagette1.shape)
imagette1 = imagette1.astype('<f4')
imagette1[imagette1mask == 0] = numpy.nan
......@@ -831,14 +831,12 @@ def pixelSearch(imagette1, imagette2, imagette1mask=None, searchCentre=None, sea
botDiff = numpy.array(imagette2.shape) - searchCentre.astype(int) - imagette1halfWindowSize - searchRange[1::2]
# Check that the search range doesn't go outside imagette2, if it does, redact it
#print("SR in", searchRange)
if topDiff[0] < 0: searchRange[0] -= topDiff[0]
if topDiff[1] < 0: searchRange[2] -= topDiff[1]
if topDiff[2] < 0: searchRange[4] -= topDiff[2]
if botDiff[0] < 0: searchRange[1] += botDiff[0]
if botDiff[1] < 0: searchRange[3] += botDiff[1]
if botDiff[2] < 0: searchRange[5] += botDiff[2]
#print("SR out", searchRange)
# Run the actual pixel search
# print imagette1.shape, imagette2.shape, searchCentre, searchRange
......
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