Vous avez reçu un message "Your GitLab account has been locked ..." ? Pas d'inquiétude : lisez cet article https://docs.gricad-pages.univ-grenoble-alpes.fr/help/unlock/

Commit 463815a9 authored by Edward Andò's avatar Edward Andò
Browse files

Closes #203 with a parallel implementation

parent ce027491
Pipeline #66793 passed with stages
in 27 minutes and 5 seconds
......@@ -44,6 +44,7 @@ parser = argparse.ArgumentParser(description="spam-ddic "+spam.helpers.optionsPa
# Parse arguments with external helper function
args = spam.helpers.optionsParser.ddicParser(parser)
if args.PROCESSES is None: args.PROCESSES = multiprocessing.cpu_count()
print("spam-ddic -- Current Settings:")
argsDict = vars(args)
......@@ -129,14 +130,6 @@ if args.PHIFILE is not None:
PhiField[label] = PhiInit.copy()
PhiField[label][0:3, -1] = spam.deformation.decomposePhi(PhiInit.copy(), PhiCentre=PhiFromFile["fieldCoords"][0], PhiPoint=centresOfMass[label])["t"]
## Now recompute F to be only rigid, and calculate rigid-body translations for each point
#if args.REGSUB:
#PhiInitRigid = spam.deformation.computeRigidPhi(PhiInit.copy())
#for label in range(labelPositions.shape[0]):
#rigidDisp[label] = spam.deformation.decomposePhi(PhiInitRigid.copy(),
#PhiCentre=regCentre,
#PhiPoint=labelPositions[label])["t"]
# If the read Phi-file contains multiple lines it's an F field!
else:
#print("spam-ddic: Assuming loaded PhiFile is coherent with the current run (i.e., labels are the same).")
......@@ -175,9 +168,9 @@ def correlateOneLabel(label):
im1, im2,
[0, 0, 0, 0, 0, 0], # Search range, don't worry about it
boundingBoxes, centresOfMass,
margin=args.LABEL_CORRELATE_MARGIN,
margin=args.MARGIN,
labelDilate=labelDilateCurrent,
maskOtherLabels=args.LABEL_CORRELATE_MASK_OTHERS,
maskOtherLabels=args.MASK_OTHERS,
applyF='no',
volumeThreshold=args.VOLUME_THRESHOLD)
......@@ -197,15 +190,15 @@ def correlateOneLabel(label):
registerReturns = spam.DIC.correlate.registerMultiscale(imagetteReturns['imagette1'],
imagetteReturns['imagette2'],
args.LABEL_CORRELATE_MULTISCALE_BINNING,
args.MULTISCALE_BINNING,
im1mask=imagetteReturns['imagette1mask'],
margin=1,
PhiInit=PhiTemp,
PhiRigid=args.LABEL_CORRELATE_RIGID,
updateGradient=args.LABEL_CORRELATE_UPDATE_GRADIENT,
maxIterations=args.LABEL_CORRELATE_MAX_ITERATIONS,
deltaPhiMin=args.LABEL_CORRELATE_MIN_PHI_CHANGE,
interpolationOrder=args.LABEL_CORRELATE_INTERPOLATION_ORDER,
PhiRigid=args.CORRELATE_RIGID,
updateGradient=args.UPDATE_GRADIENT,
maxIterations=args.MAX_ITERATIONS,
deltaPhiMin=args.MIN_PHI_CHANGE,
interpolationOrder=args.INTERPOLATION_ORDER,
verbose=args.DEBUG,
imShowProgress=args.DEBUG)
goodPhi = registerReturns['Phi']
......@@ -223,8 +216,6 @@ if args.SKIP_PARTICLES:
else:
labelsToCorrelate = numpy.delete(labelsToCorrelate, 0)
if args.PROCESSES is None: args.PROCESSES = multiprocessing.cpu_count()
print("\n\tStarting Discrete DIC (with {} process{})".format(args.PROCESSES, 'es' if args.PROCESSES > 1 else ''))
widgets = [progressbar.FormatLabel(''), ' ', progressbar.Bar(), ' ', progressbar.AdaptiveETA()]
pbar = progressbar.ProgressBar(widgets=widgets, maxval=len(labelsToCorrelate))
......
......@@ -52,6 +52,8 @@ if len(args.inFiles) < 2:
print("\nldic: Did not receive enough input images... you need (at least) two to tango...")
exit()
if args.PROCESSES is None: args.PROCESSES = multiprocessing.cpu_count()
print("spam-ldic -- Current Settings:")
argsDict = vars(args)
for key in sorted(argsDict):
......@@ -117,10 +119,6 @@ def correlateOneNode(nodeNumber):
initialDisplacement = numpy.round(PhiInit[0:3, 3]).astype(int)
PhiInit[0:3,-1] -= initialDisplacement
#print(imagetteReturns['imagette1'].shape)
#print(imagetteReturns['imagette2'].shape)
#print()
registerReturns = spam.DIC.register(imagetteReturns['imagette1'],
imagetteReturns['imagette2'],
im1mask=imagetteReturns['imagette1mask'],
......@@ -208,21 +206,12 @@ for im2number in range(1, len(args.inFiles)):
nodePositions, nodesDim = spam.DIC.makeGrid(im1.shape, args.NS)
numberOfNodes = nodePositions.shape[0]
# We have a registration to apply to all points.
# This is done in 2 steps:
# 1. by copying the registration's little F to the Fs of all points
# 2. by calling the decomposePhi function to compute the translation of each point
# Now that we know how many points we want to correlate, initalise PhiField
PhiField = numpy.zeros((numberOfNodes, 4, 4))
for node in range(numberOfNodes):
if args.APPLY_F == "all":
PhiField[node] = PhiInit.copy()
elif args.APPLY_F == "rigid":
PhiField[node] = spam.deformation.computeRigidPhi(PhiInit.copy())
else:
PhiField[node] = numpy.eye(4)
# trasnlation application
PhiField[node][0:3, -1] = spam.deformation.decomposePhi(PhiInit.copy(), PhiCentre=PhiFromFile["fieldCoords"][0], PhiPoint=nodePositions[node])["t"]
PhiField = spam.DIC.applyRegistrationToPoints(PhiInit.copy(),
PhiFromFile["fieldCoords"][0],
nodePositions,
applyF = args.APPLY_F,
nProcesses = args.PROCESSES,
verbose = False)
else: # we have a Phi field and not a registration
nodePositions = PhiFromFile["fieldCoords"]
......@@ -263,8 +252,6 @@ for im2number in range(1, len(args.inFiles)):
returnStatus = numpy.zeros((numberOfNodes))
deltaPhiNorm = numpy.zeros((numberOfNodes))
if args.PROCESSES is None: args.PROCESSES = multiprocessing.cpu_count()
print("\n\tStarting local dic (with {} process{})".format(args.PROCESSES, 'es' if args.PROCESSES > 1 else ''))
widgets = [progressbar.FormatLabel(''), ' ', progressbar.Bar(), ' ', progressbar.AdaptiveETA()]
......
......@@ -105,7 +105,7 @@ if PhiFromFile is None:
if len(args.PHIFILE2) > 0:
print(f"\n\nspam-passPhiField: I see {len(args.PHIFILE2)} -pf2 file{'s' if len(args.PHIFILE2) > 1 else ''}, so will merge grid + discrete -> grid.")
# check that -pf file is a grid
assert(PhiFromFile['numberOfLabels'] == 0), "in merge mode, -pf1 should be a grid file from spam-ldic or grid pixelSearch"
......@@ -133,7 +133,7 @@ if len(args.PHIFILE2) > 0:
binningLabelled=1,
alwaysLabel=args.MERGE_PREFER_LABEL)
# merge
print(f"\n\ndone. Now saving (without 'mergeSource' :( )...")
print(f"\n\ndone. Now saving (without 'mergeSource' field :( )...")
outputNumberOfNodes = PhiFromFile['fieldCoords'].shape[0]
outputNodePositions = PhiFromFile['fieldCoords']
......@@ -176,9 +176,6 @@ else:
# start setting up
outputNumberOfNodes = outputNodePositions.shape[0]
# Now that we know how many points we want to correlate, initalise PhiField
outputPhiField = numpy.zeros((outputNumberOfNodes, 4, 4))
# If the read Phi-file has only one line -- it's a single point registration!
# We can either apply it to a grid or to labels
......@@ -202,31 +199,19 @@ else:
print("\t\t", decomposedPhiInit['z'])
del decomposedPhiInit
# We have a registration to apply to all points.
# This is done in 2 steps:
# 1. by copying the registration's little F to the Fs of all points
# 2. by calling the decomposePhi function to compute the translation of each point
if args.MODE is "all":
for node in range(outputNumberOfNodes):
outputPhiField[node] = PhiInit
elif args.MODE is "rigid":
PhiInitRigid = spam.deformation.computeRigidPhi(PhiInit)
for node in range(outputNumberOfNodes):
outputPhiField[node] = PhiInitRigid
else: # args.MODE is "no"
for node in range(outputNumberOfNodes):
outputPhiField[node] = numpy.eye(4)
# In any case, apply the displacements to the points
for node in range(outputNumberOfNodes):
outputPhiField[node][0:3, -1] = spam.deformation.decomposePhi(PhiInit,
PhiCentre=PhiFromFile["fieldCoords"][0],
PhiPoint=outputNodePositions[node])["t"]
outputPhiField = spam.DIC.applyRegistrationToPoints(PhiInit,
PhiFromFile["fieldCoords"][0],
outputNodePositions,
applyF=args.APPLY_F,
nProcesses=args.PROCESSES,
verbose=False)
###############################################################
### Input Phi file is a Phi FIELD
###############################################################
else:
outputPhiField = numpy.zeros((outputNumberOfNodes, 4, 4))
# Interpolate these?
pixelSearchCC = numpy.zeros((outputNumberOfNodes), dtype=float)
error = numpy.zeros((outputNumberOfNodes), dtype=float)
......@@ -244,9 +229,6 @@ else:
goodInputNodePositions = inputNodePositions[goodPointsMask]
goodInputPhiField = inputPhiField[goodPointsMask]
print(goodInputNodePositions.shape)
print(goodInputPhiField.shape)
# Check neighbour inputs, either args.NEIGHBOUR_RADIUS or args.NUMBER_OF_NEIGHBOURS should be set.
if args.NEIGHBOUR_RADIUS is not None and args.NUMBER_OF_NEIGHBOURS is not None:
print("Both number of neighbours and neighbour radius are set, I'm taking the radius and ignoring the number of neighbours")
......@@ -268,7 +250,7 @@ else:
outputNodePositions,
nNeighbours=args.NUMBER_OF_NEIGHBOURS,
neighbourRadius=args.NEIGHBOUR_RADIUS,
mode=args.MODE,
interpolateF=args.APPLY_F,
nProcesses=args.PROCESSES,
verbose=True)
......
......@@ -49,6 +49,8 @@ parser = argparse.ArgumentParser(description="spam-pixelSearch "+spam.helpers.op
# Parse arguments with external helper function
args = spam.helpers.optionsParser.pixelSearch(parser)
if args.PROCESSES is None: args.PROCESSES = multiprocessing.cpu_count()
print("spam-pixelSearch -- Current Settings:")
argsDict = vars(args)
for key in sorted(argsDict):
......@@ -131,21 +133,12 @@ if args.PHIFILE is not None:
nodePositions, nodesDim = spam.DIC.makeGrid(im1.shape, args.NS)
numberOfNodes = nodePositions.shape[0]
# We have a registration to apply to all points.
# This is done in 2 steps:
# 1. by copying the registration's little F to the Fs of all points
# 2. by calling the decomposePhi function to compute the translation of each point
# Now that we know how many points we want to correlate, initalise PhiField
PhiField = numpy.zeros((numberOfNodes, 4, 4))
for node in range(numberOfNodes):
if args.APPLY_F == "all":
PhiField[node] = PhiInit.copy()
elif args.APPLY_F == "rigid":
PhiField[node] = spam.deformation.computeRigidPhi(PhiInit.copy())
else:
PhiField[node] = numpy.eye(4)
# trasnlation application
PhiField[node][0:3, -1] = spam.deformation.decomposePhi(PhiInit.copy(), PhiCentre=PhiFromFile["fieldCoords"][0], PhiPoint=nodePositions[node])["t"]
PhiField = spam.DIC.applyRegistrationToPoints(PhiInit,
PhiFromFile["fieldCoords"][0],
nodePositions,
applyF=args.APPLY_F,
nProcesses=args.PROCESSES,
verbose=False)
# If the read Phi-file contains multiple lines it's an F field!
else:
......@@ -181,7 +174,7 @@ else: # no Phi file
print("spam-pixelSearch: You're in regular grid mode, but no -ns is set and no Phi Field TSV has been passed, exiting.")
exit()
nodePositions, nodesDim = spam.DIC.makeGrid(im1.shape, args.NS)
numberOfNodes = nodePositions.shape[0]
numberOfNodes = nodePositions.shape[0]
PhiField = numpy.zeros((numberOfNodes, 4, 4))
for node in range(numberOfNodes):
......@@ -262,8 +255,6 @@ if args.LAB1 is not None: firstNode = 1; finishedNodes = 1; returnStatus[0] = 0
else: firstNode = 0; finishedNodes = 0
if args.PROCESSES is None: args.PROCESSES = multiprocessing.cpu_count()
print("\n\tStarting Pixel search (with {} process{})".format(args.PROCESSES, 'es' if args.PROCESSES > 1 else ''))
widgets = [progressbar.FormatLabel(''), ' ', progressbar.Bar(), ' ', progressbar.AdaptiveETA()]
......
......@@ -822,7 +822,7 @@ class testAll(unittest.TestCase):
"-lab1", testFolder + "Lab0.tif",
"-od", testFolder + "",
"-pre", testFolder + "balls-2c-reg",
"-m", "rigid"
"-F", "rigid"
])
self.assertEqual(exitCode, 0)
DDICresult2c = spam.helpers.readCorrelationTSV(testFolder + "balls-2c-reg-passed-labelled.tsv")
......@@ -856,8 +856,8 @@ class testAll(unittest.TestCase):
testFolder + "Step0.tif",
testFolder + "Lab0.tif",
testFolder + "Step1.tif",
"-lcmsb", "2",
"-lcnr",
"-msb", "2",
"-nr",
"-ld", "2",
"-od", testFolder + "",
"-pre", "balls-3a"])
......@@ -1074,7 +1074,7 @@ class testAll(unittest.TestCase):
testFolder + "extreme-im1.tif",
testFolder + "extreme-im1lab.tif",
testFolder + "extreme-im2.tif",
"-lcm", "10", "-ld", "2"])
"-m", "10", "-ld", "2"])
self.assertEqual(exitCode, 0)
TSVextreme = spam.helpers.readCorrelationTSV(testFolder + "extreme-im1-extreme-im2-ddic.tsv")
#self.assertAlmostEqual(TSV['PhiField'][i][u,-1], TSV2['PhiField'][i][u,-1], places=1)
......
......@@ -351,7 +351,7 @@ def estimateDisplacementFromQuadraticFit(fieldCoords, displacements, pointsToEst
return estimatedDisplacements
def interpolatePhiField(fieldCoords, PhiField, pointsToInterpolate, nNeighbours=None, neighbourRadius=None, mode="all", neighbourDistanceWeight='inverse', nProcesses=nProcessesDefault, verbose=False):
def interpolatePhiField(fieldCoords, PhiField, pointsToInterpolate, nNeighbours=None, neighbourRadius=None, interpolateF="all", neighbourDistanceWeight='inverse', nProcesses=nProcessesDefault, verbose=False):
"""
This function interpolates components of a Phi field at a given number of points, using scipy's KD-tree to find neighbours.
......@@ -378,9 +378,9 @@ def interpolatePhiField(fieldCoords, PhiField, pointsToInterpolate, nNeighbours=
This OR nNeighbours must be set.
Default = None
mode : string, optional
interpolateF : string, optional
Interpolate the whole Phi, just the rigid part, or just the displacement?
Corresponding options are 'all', 'rigid', 'disp'
Corresponding options are 'all', 'rigid', 'no'
Default = "all"
neighbourDistanceWeight : string, optional
......@@ -412,7 +412,7 @@ def interpolatePhiField(fieldCoords, PhiField, pointsToInterpolate, nNeighbours=
# extract the Phi matrices of the bad points
outputPhiField = numpy.zeros((numberOfPointsToInterpolate, 4, 4), dtype=PhiField.dtype)
assert mode in ["all", "rigid", "disp"], "spam.DIC.interpolatePhiField(): mode argument should either be 'all', 'rigid', or 'disp'"
assert interpolateF in ["all", "rigid", "no"], "spam.DIC.interpolatePhiField(): interpolateF argument should either be 'all', 'rigid', or 'no'"
assert neighbourDistanceWeight in ["inverse", "gaussian", "mean", "median"], "spam.DIC.interpolatePhiField(): neighbourDistanceWeight argument should be 'inverse', 'gaussian', 'mean', or 'median'"
# check if neighbours are selected based on radius or based on number, default is radius
assert (nNeighbours is None) != (neighbourRadius is None), "spam.DIC.interpolatePhiField(): One and only one of nNeighbours and neighbourRadius must be passed"
......@@ -427,7 +427,7 @@ def interpolatePhiField(fieldCoords, PhiField, pointsToInterpolate, nNeighbours=
pointToInterpolate = pointsToInterpolate[pointNumber]
outputPhi = numpy.zeros((4, 4), dtype=PhiField.dtype)
outputPhi[-1] = [0, 0, 0, 1]
if mode == 'disp':
if interpolateF == 'no':
outputPhi[0:3, 0:3] = numpy.eye(3)
#######################################################
......@@ -451,9 +451,9 @@ def interpolatePhiField(fieldCoords, PhiField, pointsToInterpolate, nNeighbours=
### Check if there is only one neighbour
#######################################################
if indices.size == 1:
if mode in ["all", "rigid"]: # We need to interpolate all 12 components of Phi
if interpolateF in ["all", "rigid"]: # We need to interpolate all 12 components of Phi
outputPhi = PhiField[indices].copy()
if mode == "rigid":
if interpolateF == "rigid":
outputPhi = spam.deformation.computeRigidPhi(outputPhi)
else: # interpolate only displacements
outputPhi[0:3, -1] = PhiField[indices, 0:3, -1].copy()
......@@ -475,12 +475,12 @@ def interpolatePhiField(fieldCoords, PhiField, pointsToInterpolate, nNeighbours=
weights = numpy.zeros_like(distances)
weights[len(distances)//2] = 1
if mode == 'disp':
if interpolateF == 'no':
outputPhi[0:3,-1] = numpy.sum(PhiField[indices, 0:3,-1] * weights[:, numpy.newaxis], axis=0) / weights.sum()
else:
outputPhi[:-1] = numpy.sum(PhiField[indices, :-1] * weights[:, numpy.newaxis, numpy.newaxis], axis=0) / weights.sum()
if mode == "rigid":
if interpolateF == "rigid":
outputPhi = spam.deformation.computeRigidPhi(outputPhi)
return pointNumber, outputPhi
......@@ -794,6 +794,85 @@ def getDisplacementFromNeighbours(labIm, DVC, fileName, method='getLabel', centr
print('spam.deformation.deformationField.getDisplacementFromNeighbours(): Aborting')
def applyRegistrationToPoints(Phi, PhiCentre, points, applyF='all', nProcesses=nProcessesDefault, verbose=False):
"""
This function takes a whole-image registration and applies it to a set of points
Parameters
----------
Phi : 4x4 numpy array of floats
Measured Phi function to apply to points
PhiCentre : 3-component list of floats
Origin where the Phi is measured (normally the middle of the image unless masked)
points : nx3 numpy array of floats
Points to apply the Phi to
applyF : string, optional
The whole Phi is *always* applied to the positions of the points to get their displacement.
This mode *only* controls what is copied into the F for each point, everything, only rigid, or only displacements?
Corresponding options are 'all', 'rigid', 'no'
Default = "all"
nProcesses : integer, optional
Number of processes for multiprocessing
Default = number of CPUs in the system
verbose : boolean, optional
Print progress?
Default = True
Returns
-------
PhiField : nx4x4 numpy array of floats
Output Phi field
"""
assert applyF in ["all", "rigid", "no"], "spam.DIC.applyRegistrationToPoints(): applyF should be 'all', 'rigid', or 'no'"
numberOfPoints = points.shape[0]
PhiField = numpy.zeros((numberOfPoints, 4, 4), dtype=float)
if applyF is "rigid":
PhiRigid = spam.deformation.computeRigidPhi(Phi)
global applyPhiToPoint
def applyPhiToPoint(n):
# We have a registration to apply to all points.
# This is done in 2 steps:
# 1. by copying the registration's little F to the Fs of all points (depending on mode)
# 2. by calling the decomposePhi function to compute the translation of each point
if applyF is "all":
outputPhi = Phi.copy()
if applyF is "rigid":
outputPhi = PhiRigid.copy()
else: # applyF is displacement only
outputPhi = numpy.eye(4)
outputPhi[0:3, -1] = spam.deformation.decomposePhi(Phi,
PhiCentre=PhiCentre,
PhiPoint=points[n])["t"]
return n, outputPhi
if verbose:
pbar = progressbar.ProgressBar(maxval=numberOfPoints).start()
finishedPoints = 0
with multiprocessing.Pool(processes=nProcesses) as pool:
for returns in pool.imap_unordered(applyPhiToPoint, range(numberOfPoints)):
if verbose:
finishedPoints += 1
pbar.update(finishedPoints)
PhiField[returns[0]] = returns[1]
pool.close()
pool.join()
if verbose: pbar.finish()
return PhiField
def mergeRegistrationAndDiscreteFields(regTSV, discreteTSV, fileName, regTSVBinRatio=1):
"""
This function merges a registration TSV with a discrete TSV.
......
......@@ -454,7 +454,7 @@ def _binarisation(im, threshold=0.0, boolean=False, op='>', mask=None):
return phases
def slicePadded(im, startStop, createMask=False, padValue=0):
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)
......@@ -500,7 +500,7 @@ def slicePadded(im, startStop, createMask=False, padValue=0):
# 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)):
print("spam.helpers.slicePadded(): The extracted padded slice doesn't not touch the image!")
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:
......
This diff is collapsed.
......@@ -578,7 +578,116 @@ class TestFunctionDVC(unittest.TestCase):
#self.assertEqual(mergeSource[0], 0)
#self.assertEqual(mergeSource[1], 1)
def test_interpolatePhiField(self):
# Fake Phi field -- this will be defined at 0,0,0 and +100 in all directions
# z = 0 --> transformation below:
transformation = {'t': [5., 0., 0.],
'r': [0., 0., 0.]}
Phi = spam.deformation.computePhi(transformation)
#print(Phi)
centre = [50, 50, 50]
# z = 100 --> nothing
fieldCoords = numpy.array([[ 0, 0, 0],
[ 0, 0,100],
[ 0,100, 0],
[ 0,100,100],
[100, 0, 0],
[100, 0,100],
[100,100, 0],
[100,100,100]])
fieldValues = numpy.zeros([fieldCoords.shape[0],4,4])
#fieldValues[0:4] = spam.deformation.computePhi(transformation)
#fieldValues[4:8] = numpy.eye(4)
for n in range(fieldCoords.shape[0]):
fieldValues[n] = numpy.eye(4)
fieldValues[n, 0:3, 0:3] = Phi[0:3, 0:3]
fieldValues[n, 0:3, -1] = spam.deformation.decomposePhi(Phi, PhiPoint=fieldCoords[n], PhiCentre=centre)['t']
interpCoords = numpy.array([[50.,50.,50.]])
interpolatedPhi = spam.DIC.interpolatePhiField(fieldCoords, fieldValues, interpCoords, nNeighbours=1, verbose=True)
self.assertTrue(numpy.allclose(interpolatedPhi[0, 0:3, -1], transformation['t'], atol=0.01))
self.assertTrue(numpy.allclose(interpolatedPhi[0, 0:3, 0:3], Phi[0:3, 0:3], atol=0.01))
interpolatedPhi = spam.DIC.interpolatePhiField(fieldCoords,
fieldValues,
interpCoords,
nNeighbours=1,
interpolateF='rigid')
self.assertTrue(numpy.allclose(interpolatedPhi[0, 0:3, -1], transformation['t'], atol=0.01))
self.assertTrue(numpy.allclose(interpolatedPhi[0, 0:3, 0:3], Phi[0:3, 0:3], atol=0.01))
interpolatedPhi = spam.DIC.interpolatePhiField(fieldCoords,
fieldValues,
interpCoords,
nNeighbours=1,
interpolateF='no')
self.assertTrue(numpy.allclose(interpolatedPhi[0, 0:3, -1], transformation['t'], atol=0.01))
self.assertTrue(numpy.allclose(interpolatedPhi[0, 0:3, 0:3], numpy.eye(3), atol=0.01))
interpolatedPhi = spam.DIC.interpolatePhiField(fieldCoords,
fieldValues,
interpCoords,
nNeighbours=2)
self.assertTrue(numpy.allclose(interpolatedPhi[0, 0:3, -1], transformation['t'], atol=0.01))
self.assertTrue(numpy.allclose(interpolatedPhi[0, 0:3, 0:3], Phi[0:3, 0:3], atol=0.01))
interpolatedPhi = spam.DIC.interpolatePhiField(fieldCoords,
fieldValues,
interpCoords,
nNeighbours=2,
interpolateF='rigid')
self.assertTrue(numpy.allclose(interpolatedPhi[0, 0:3, -1], transformation['t'], atol=0.01))
self.assertTrue(numpy.allclose(interpolatedPhi[0, 0:3, 0:3], Phi[0:3, 0:3], atol=0.01))
interpolatedPhi = spam.DIC.interpolatePhiField(fieldCoords,
fieldValues,
interpCoords,
nNeighbours=2,
interpolateF='no')
self.assertTrue(numpy.allclose(interpolatedPhi[0, 0:3, -1], transformation['t'], atol=0.01))
self.assertTrue(numpy.allclose(interpolatedPhi[0, 0:3, 0:3], numpy.eye(3), atol=0.01))
interpolatedPhi = spam.DIC.interpolatePhiField(fieldCoords,
fieldValues,
interpCoords,
neighbourRadius=10)
self.assertTrue(numpy.isnan(interpolatedPhi[0]).sum() == 16)
interpolatedPhi = spam.DIC.interpolatePhiField(fieldCoords,
fieldValues,
interpCoords,
neighbourRadius=100)
self.assertTrue(numpy.allclose(interpolatedPhi[0, 0:3, -1], transformation['t'], atol=0.01))
self.assertTrue(numpy.allclose(interpolatedPhi[0, 0:3, 0:3], numpy.eye(3), atol=0.01))
# Test all new interplation options
interpolatedPhi = spam.DIC.interpolatePhiField(fieldCoords,
fieldValues,
interpCoords,
nNeighbours=8,
neighbourDistanceWeight='gaussian')
self.assertTrue(numpy.allclose(interpolatedPhi[0, 0:3, -1], transformation['t'], atol=0.01))
self.assertTrue(numpy.allclose(interpolatedPhi[0, 0:3, 0:3], Phi[0:3, 0:3], atol=0.01))
interpolatedPhi = spam.DIC.interpolatePhiField(fieldCoords,
fieldValues,
interpCoords,
nNeighbours=8,
neighbourDistanceWeight='mean')
self.assertTrue(numpy.allclose(interpolatedPhi[0, 0:3, -1], transformation['t'], atol=0.01))
self.assertTrue(numpy.allclose(interpolatedPhi[0, 0:3, 0:3], Phi[0:3, 0:3], atol=0.01))
interpolatedPhi = spam.DIC.interpolatePhiField(fieldCoords,
fieldValues,
interpCoords,
nNeighbours=8,
neighbourDistanceWeight='median')
self.assertTrue(numpy.allclose(interpolatedPhi[0, 0:3, -1], transformation['t'], atol=0.01))
self.assertTrue(numpy.allclose(interpolatedPhi[0, 0:3, 0:3], Phi[0:3, 0:3], atol=0.01))
def test_LQC(self):
points, _ = spam.DIC.grid.makeGrid((10, 10, 10), 2)
......
......@@ -724,81 +724,6 @@ class testAll(unittest.TestCase):
self.assertEqual(decomposedF['dev'][otherTets].tolist(), numpy.zeros(otherTets.sum()).tolist())
self.assertTrue(numpy.sum(decomposedF['dev'][tetsWithPointToMove] != 0) == len(tetsWithPointToMove))
def test_interpolatePhiField(self):
# Fake Phi field -- this will be defined at 0,0,0 and +100 in all directions
# z = 0 --> transformation below:
transformation = {'t': [5., 0., 0.],
'r': [0., 0., 0.]}
Phi = spam.deformation.computePhi(transformation)
#print(Phi)
centre = [50, 50, 50]
# z = 100 --> nothing