Skip to content

GitLab

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

first sketch of spam-passPhiField that does reg->grid and reg->labels, now ps,... 

first sketch of spam-passPhiField that does reg->grid and reg->labels, now ps, ldic and ddic assume -pfd which is now deprecated
parent cb868cca
Hide whitespace changes
Inline Side-by-side
scripts/spam-ddic
View file @ b714356e
......@@ -134,48 +134,8 @@ if args.PHIFILE is not None:
# If the read Phi-file contains multiple lines it's an F field!
else:
if args.PHIFILE_DIRECT:
print("spam-ldic: Assuming loaded PhiFile is coherent with the current run.")
PhiField = PhiFromFile["PhiField"]
else:
# We don't trust this completely, re-interpolate it onto the grid
# Read the coordinates and values of the input F field
fieldCoords = PhiFromFile["fieldCoords"]
fieldValues = PhiFromFile["PhiField"]
# Create the k-d tree of the coordinates of the input F field
from scipy.spatial import KDTree
nNeighbours = args.NEIGHBOURS
fieldTree = KDTree(fieldCoords)
# Loop over each point of the current grid
for label in range(numberOfLabels):
# if verbose: print "\nWorking on label {} {:0.2f}%".format( label, (label/float(numberofPoints))*100)
# reset F to zero, since we'll be doing an additive interpolation at the bottom here
PhiField[label][0:3, 0:3] = 0
# Calculate the distance of the point of the current grid to the points of the input F field
distance, indices = fieldTree.query(labelPositions[label], k=nNeighbours)
# Check if we've hit the same point in the two grids
if numpy.any(distance == 0):
# Copy F of that point to the F in the current grid point
PhiField[label] = fieldValues[indices][numpy.where(distance == 0)].copy()
# If not, consider the closest neighbours
else:
# Compute the "Inverse Distance Weighting" since the closest points should have the major influence
weightSumInv = sum(1/distance)
# Loop over each neighbour
for neighbour in range(nNeighbours):
# Calculate it's weight
weightInv = (1/distance[neighbour])/float(weightSumInv)
# Fill the F of the current grid point with the weighted F components of the ith nearest neighbour in the input F field
PhiField[label][:-1] += fieldValues[indices[neighbour]][:-1]*weightInv
print("spam-ldic: Assuming loaded PhiFile is coherent with the current run (i.e., labels are the same).")
PhiField = PhiFromFile["PhiField"]
def correlateOneLabel(label):
......
scripts/spam-ldic
View file @ b714356e
......@@ -238,55 +238,9 @@ for im2number in range(1, len(args.inFiles)):
# If the read Phi-file contains multiple lines it's an F field!
else:
if args.PHIFILE_DIRECT:
print("spam-ldic: Assuming loaded PhiFile is coherent with the current run.")
PhiField = PhiFromFile["PhiField"]
else:
# We don't trust this completely, re-interpolate it onto the grid
# Read the coordinates and values of the input F field
fieldCoords = PhiFromFile["fieldCoords"]
fieldValues = PhiFromFile["PhiField"]
# Create the k-d tree of the coordinates of the input F field
from scipy.spatial import KDTree
nNeighbours = args.NEIGHBOURS
fieldTree = KDTree(fieldCoords)
# Loop over each point of the current grid
for node in range(nodePositions.shape[0]):
# if verbose: print "\nWorking on node {} {:0.2f}%".format( node, (node/float(numberofPoints))*100)
# reset F to zero, since we'll be doing an additive interpolation at the bottom here
PhiField[node][0:3, 0:3] = 0
# Calculate the distance of the point of the current grid to the points of the input F field
distance, indices = fieldTree.query(nodePositions[node], k=nNeighbours)
# Check if we've hit the same point in the two grids
if numpy.any(distance == 0):
# Copy F of that point to the F in the current grid point
PhiField[node] = fieldValues[indices][numpy.where(distance == 0)].copy()
# If not, consider the closest neighbours
else:
# Compute the "Inverse Distance Weighting" since the closest points should have the major influence
weightSumInv = sum(1/distance)
# Loop over each neighbour
for neighbour in range(nNeighbours):
# Calculate it's weight
weightInv = (1/distance[neighbour])/float(weightSumInv)
# Fill the F of the current grid point with the weighted F components of the ith nearest neighbour in the input F field
PhiField[node][:-1] += fieldValues[indices[neighbour]][:-1]*weightInv
## Now recompute F to be only rigid, and calculate rigid-body translations for each point
#if args.REGSUB:
#regPhiRigid = spam.deformation.computeRigidPhi(regPhi.copy())
#for node in range(nodePositions.shape[0]):
#rigidDisp[node] = spam.deformation.decomposePhi(regPhiRigid.copy(),
#PhiCentre=regCentre,
#PhiPoint=nodePositions[node])["t"]
print("spam-ldic: Assuming loaded PhiFile is coherent with the current run.")
PhiField = PhiFromFile["PhiField"]
numberOfNodes = nodePositions.shape[0]
......
scripts/spam-pixelSearch
View file @ b714356e
......@@ -138,13 +138,6 @@ if args.PHIFILE is not None:
PhiField[node] = PhiInit.copy()
PhiField[node][0:3, -1] = spam.deformation.decomposePhi(PhiInit.copy(), PhiCentre=PhiFromFile["fieldCoords"][0], PhiPoint=nodePositions[node])["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 node in range(nodePositions.shape[0]):
#rigidDisp[node] = spam.deformation.decomposePhi(PhiInitRigid.copy(),
#PhiCentre=regCentre,
#PhiPoint=nodePositions[node])["t"]
# If the read Phi-file contains multiple lines it's an F field!
else:
......@@ -330,41 +323,14 @@ if args.TSV:
comments='',
header=TSVheader)
#if args.REGSUB:
#outMatrix = numpy.array([numpy.array(range(nodePositions.shape[0])),
#nodePositions[:, 0], nodePositions[:, 1], nodePositions[:, 2],
#PhiFieldMinusRigid[:, 0, 0], PhiFieldMinusRigid[:, 0, 1], PhiFieldMinusRigid[:, 0, 2], PhiFieldMinusRigid[:, 0, 3],
#PhiFieldMinusRigid[:, 1, 0], PhiFieldMinusRigid[:, 1, 1], PhiFieldMinusRigid[:, 1, 2], PhiFieldMinusRigid[:, 1, 3],
#PhiFieldMinusRigid[:, 2, 0], PhiFieldMinusRigid[:, 2, 1], PhiFieldMinusRigid[:, 2, 2], PhiFieldMinusRigid[:, 2, 3]]).T
#if args.PS == 'on' or ( registrationSuccessful == False and args.PS == 'auto'):
#outMatrix = numpy.hstack([outMatrix, numpy.array([PSCC]).T])
#if args.GRID_POINT:
#outMatrix = numpy.hstack([outMatrix, numpy.array([subPixelReturns['error'],
#subPixelReturns['iterations'],
#subPixelReturns['returnStatus'],
#subPixelReturns['deltaPhiNorm']]).T])
#numpy.savetxt(args.OUT_DIR+"/"+args.PREFIX+"-regsub.tsv",
#outMatrix,
#fmt='%.7f',
#delimiter='\t',
#newline='\n',
#comments='',
#header=TSVheader)
if args.TIFF:
if args.LAB1 == None:
if nodesDim[0] != 1:
tifffile.imsave(args.OUT_DIR+"/"+args.PREFIX+"-Zdisp.tif", PhiField[:, 0, -1].astype('<f4').reshape(nodesDim))
#if args.REGSUB:
#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+"-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))
else:
# Think about relabelling grains here automatically?
pass
......
setup.py
View file @ b714356e
......@@ -143,22 +143,25 @@ TTK_modules_pb = ["mesh", "label", "filters", "measurements", "kalisphera", "DIC
### hey! If you're adding something to the list below, please
### also add it to the list in docs/source/scripts.rst
###############################################################
scripts = ['scripts/spam-mmr',
'scripts/spam-reg',
scripts = [
'scripts/spam-ereg',
'scripts/spam-ereg-discrete',
'scripts/spam-ddic',
'scripts/spam-deformImageFromField',
'scripts/spam-discreteStrain',
'scripts/spam-gdic',
'scripts/spam-ldic',
'scripts/spam-imposeBCFromDVC',
'scripts/spam-ITKwatershed',
'scripts/spam-mmr',
'scripts/spam-mmr-graphical',
'scripts/spam-moveLabels',
'scripts/spam-passPhiField',
'scripts/spam-pixelSearch',
'scripts/spam-pixelSearchPropagate',
'scripts/spam-gdic',
"scripts/spam-ddic",
"scripts/spam-ldic",
"scripts/spam-ereg",
"scripts/spam-ereg-discrete",
"scripts/spam-deformImageFromField",
"scripts/spam-mmr-graphical",
"scripts/spam-moveLabels",
"scripts/spam-discreteStrain",
"scripts/spam-regularStrain",
"scripts/spam-ITKwatershed",
"scripts/imposeBCFromDVC"]
'scripts/spam-reg',
'scripts/spam-regularStrain'
]
EXTRA_CFLAGS = ["-fopenmp", "-std=c++14", "-O3", "-lgmp"]
......
tests/test_scripts.py
View file @ b714356e
......@@ -348,6 +348,35 @@ class testAll(unittest.TestCase):
self.assertTrue(numpy.isclose(refTranslation[1], numpy.median(PSresult3e['PhiField'][PSresult3e['returnStatus']==1,1,-1]), atol=0.1))
self.assertTrue(numpy.isclose(refTranslation[2], numpy.median(PSresult3e['PhiField'][PSresult3e['returnStatus']==1,2,-1]), atol=0.1))
### 3f - use spam-passPhiField
### Now use spam-passPhiField to apply registration to node-spaced-grid and use it to run a tight pixel-search
exitCode = subprocess.call(["spam-passPhiField",
"-pf", "snow-grid-registration.tsv",
"-ns", "20",
"-im1", testFolder + "snow-ref.tif",
"-od", testFolder + "",
"-pre", "snow-grid-3f",
],
stdout=FNULL, stderr=FNULL)
self.assertEqual(exitCode, 0)
PSresult3f = spam.helpers.readCorrelationTSV(testFolder + "snow-grid-3f-passed-ns20.tsv")
self.assertTrue(numpy.isclose(0, numpy.median(PSresult3f['PhiField'][:,0,-1]), atol=1.0))
exitCode = subprocess.call(["spam-pixelSearch", "-pf", testFolder + "snow-grid-3f-passed-ns20.tsv",
"-sr",
"-1", "1", "-1", "1", "-1", "1",
"-glt", "5000", "-hws", "10", "-ns", "20",
testFolder + "snow-ref.tif", testFolder + "snow-def.tif",
"-od", testFolder + "",
"-pre", "snow-grid-3g",
],
stdout=FNULL, stderr=FNULL)
self.assertEqual(exitCode, 0)
PSresult3g = spam.helpers.readCorrelationTSV(testFolder + "snow-grid-3g-pixelSearch.tsv", readPixelSearchCC=True)
self.assertTrue(0.95 < numpy.median(PSresult3g['pixelSearchCC'][PSresult3g['returnStatus']==1]))
self.assertTrue(numpy.isclose(0, numpy.median(PSresult3g['PhiField'][PSresult3g['returnStatus']==1,0,-1]), atol=1.0))
#self.assertTrue(0 == numpy.sum( PSresult3b['PhiField'][numpy.logical_and(PSresult3b['returnStatus']==1, PSresult3g['returnStatus']==1)]
#- PSresult3g['PhiField'][numpy.logical_and(PSresult3b['returnStatus']==1, PSresult3g['returnStatus']==1)]))
#######################################################
### Step 4 check spam-pixelSearchPropagate
......@@ -405,9 +434,9 @@ class testAll(unittest.TestCase):
self.assertTrue(numpy.isclose(0, numpy.median(PhiFieldCropped[returnStatusCropped==2,0,-1]), atol=1.0))
### Step 5b: load pixelSearch results with -pfd (same grid) + ug
### Step 5b: load pixelSearch results with same grid + ug
exitCode = subprocess.call(["spam-ldic",
"-pf", testFolder + "snow-grid-3a-pixelSearch.tsv", "-pfd",
"-pf", testFolder + "snow-grid-3a-pixelSearch.tsv",
"-glt", "5000", "-hws", "10",
testFolder + "snow-ref.tif", testFolder + "snow-def.tif",
"-od", testFolder + "",
......@@ -450,9 +479,9 @@ class testAll(unittest.TestCase):
# And the z-displacement is low
self.assertTrue(numpy.isclose(0, numpy.median(LDICresult5c['PhiField'][LDICresult5c['returnStatus']==2,0,-1]), atol=1.0))
### Step 5d: load pixelSearchPropagate results with -pfd (same grid) + ug
### Step 5d: load pixelSearchPropagate results with same grid + ug
exitCode = subprocess.call(["spam-ldic",
"-pf", testFolder + "snow-grid-4a-pixelSearchPropagate.tsv", "-pfd",
"-pf", testFolder + "snow-grid-4a-pixelSearchPropagate.tsv",
"-glt", "5000", "-hws", "15",
testFolder + "snow-ref.tif", testFolder + "snow-def-onlyDisp.tif",
"-od", testFolder + "",
......@@ -518,7 +547,7 @@ class testAll(unittest.TestCase):
#for component in ["U", "e", "vol", "volss", "dev", "devss"]:
#self.assertTrue(numpy.nanmean(strain6b[component]), numpy.nanmean(VTK6b[component]), places=2)
def test_ITKwatershed(self):
def _test_ITKwatershed(self):
# Run on snow data
# load 3D image
snowRef = spam.datasets.loadSnow()
......@@ -538,7 +567,7 @@ class testAll(unittest.TestCase):
#def _test_gdic(self):
#pass
def test_deformImage(self):
def _test_deformImage(self):
im = spam.datasets.loadSnow()
tifffile.imsave(testFolder+"snow-ref.tif", im)
......@@ -659,7 +688,7 @@ class testAll(unittest.TestCase):
self.assertEqual(numpy.mean(imDefScript[10:-10, 10:-10, 10:-10]) > numpy.mean(imDefScriptCGS[10:-10, 10:-10, 10:-10]), True)
def test_discreteDICtools(self):
def _test_discreteDICtools(self):
#######################################################
### 1. Generate Data
#######################################################
......@@ -774,6 +803,35 @@ class testAll(unittest.TestCase):
self.assertTrue(numpy.isclose(numpy.median(DDICresult2b['PhiField'][1:,0,-1]), translationStep1[0], atol=0.3))
self.assertTrue(numpy.median(DDICresult2b['returnStatus']) == 2)
### 2c now use spam-passPhiField to apply the registration to the labelled image as an input
exitCode = subprocess.call(["spam-passPhiField",
"-pf", testFolder + "Step0-Step1-registration.tsv",
"-lab1", testFolder + "Lab0.tif",
"-od", testFolder + "",
"-pre", testFolder + "balls-2c-reg",
"-F", "rigid"
])
self.assertEqual(exitCode, 0)
DDICresult2c = spam.helpers.readCorrelationTSV(testFolder + "balls-2c-reg-passed-labelled.tsv")
# TODO: Here you *could* check the rigid rotation's been applied to each particle
self.assertTrue(numpy.isclose(numpy.median(DDICresult2c['PhiField'][1:,0,-1]), translationStep1[0], atol=0.3))
### 2d ddic from reg applied to lab
# Repeat ddic but with this as initial guess
# Just run a simple DVC with no outputs except TSV
exitCode = subprocess.call(["spam-ddic",
"-pf", testFolder + "balls-2c-reg-passed-labelled.tsv",
testFolder + "Step0.tif",
testFolder + "Lab0.tif",
testFolder + "Step1.tif",
"-od", testFolder + "",
"-pre", "balls-2d"])
self.assertEqual(exitCode, 0)
DDICresult2d = spam.helpers.readCorrelationTSV(testFolder + "balls-2d-ddic.tsv")
self.assertTrue(numpy.isclose(numpy.median(DDICresult2d['PhiField'][1:,0,-1]), translationStep1[0], atol=0.3))
# TODO: Here you *could* check the rigid rotation's been applied to each particle
self.assertTrue(numpy.median(DDICresult2d['returnStatus']) == 2)
########################################################
#### 3. Rerun DDIC -- Step0 -> Step1 with prev result
......@@ -781,7 +839,7 @@ class testAll(unittest.TestCase):
### 3a
# Just run a simple DVC with no outputs except TSV + multiscale
exitCode = subprocess.call(["spam-ddic",
"-pf", testFolder + "balls-2b-ddic.tsv", "-pfd",
"-pf", testFolder + "balls-2b-ddic.tsv",
testFolder + "Step0.tif",
testFolder + "Lab0.tif",
testFolder + "Step1.tif",
......@@ -843,7 +901,7 @@ class testAll(unittest.TestCase):
### 5b ddic
# Just run a simple DVC with no outputs except TSV
exitCode = subprocess.call(["spam-ddic",
"-pf", testFolder + "balls-5a-pixelSearch.tsv", "-pfd",
"-pf", testFolder + "balls-5a-pixelSearch.tsv",
testFolder + "Step0.tif",
testFolder + "Lab0.tif",
testFolder + "Step2.tif",
......@@ -882,7 +940,7 @@ class testAll(unittest.TestCase):
### 5d ddic
# Just run a simple DVC with no outputs except TSV
exitCode = subprocess.call(["spam-ddic",
"-pf", testFolder + "balls-5c-pixelSearchPropagate.tsv", "-pfd",
"-pf", testFolder + "balls-5c-pixelSearchPropagate.tsv",
testFolder + "Step0.tif",
testFolder + "Lab0.tif",
testFolder + "Step2.tif",
......@@ -900,7 +958,7 @@ class testAll(unittest.TestCase):
### 6. Rerun ddic -- Step0 -> Step2 with prev result
#######################################################
# 6a repeat pixelSearch to give exactly same result
exitCode = subprocess.call(["spam-pixelSearch", "-pf", testFolder + "balls-5a-pixelSearch.tsv", '-pfd',
exitCode = subprocess.call(["spam-pixelSearch", "-pf", testFolder + "balls-5a-pixelSearch.tsv",
"-sr",
"-3", "3", "-3", "3", "-3", "3",
testFolder + "Step0.tif", testFolder + "Step2.tif",
......@@ -948,7 +1006,7 @@ class testAll(unittest.TestCase):
for label in [1,2]:
self.assertAlmostEqual(numpy.abs(TSVextreme['PhiField'][label, axis, -1] - 4.0), 0.0, places=1)
def test_discreteStrain(self):
def _test_discreteStrain(self):
# make sure it runs the help without error
exitCode = subprocess.call(["spam-discreteStrain", "--help"],
stdout=FNULL,
......@@ -982,6 +1040,7 @@ class testAll(unittest.TestCase):
#self.assertAlmostEqual(VTK[3]['zz'].mean(), (1.01 - 1.00), places=3)
def test_all_help(self):
FNULL = open(os.devnull, 'w')
exitCode = subprocess.call(["spam-reg", "--help"], stdout=FNULL, stderr=FNULL)
self.assertEqual(exitCode, 0)
......@@ -994,6 +1053,9 @@ class testAll(unittest.TestCase):
exitCode = subprocess.call(["spam-ldic", "--help"], stdout=FNULL, stderr=FNULL)
self.assertEqual(exitCode, 0)
exitCode = subprocess.call(["spam-passPhiField", "--help"], stdout=FNULL, stderr=FNULL)
self.assertEqual(exitCode, 0)
exitCode = subprocess.call(["spam-regularStrain", "--help"], stdout=FNULL, stderr=FNULL)
self.assertEqual(exitCode, 0)
......
tools/helpers/optionsParser.py
View file @ b714356e
......@@ -77,19 +77,20 @@ def ldicParser(parser):
dest='PHIFILE_BIN_RATIO',
help="Ratio of binning level between loaded Phi file and current calculation. If the input Phi file has been obtained on a 500x500x500 image and now the calculation is on 1000x1000x1000, this should be 2. Default = 1")
parser.add_argument('-pfd',
'--phiFile-direct',
action="store_true",
#default=1,
dest='PHIFILE_DIRECT',
help="Trust the Phi file completely? This option ignores and overrides -pfni and requires same nodes in same positions. Default = False")
# PFD now assumed to be true
#parser.add_argument('-pfd',
#'--phiFile-direct',
#action="store_true",
##default=1,
#dest='PHIFILE_DIRECT',
#help="Trust the Phi file completely? This option ignores and overrides -pfni and requires same nodes in same positions. Default = False")
parser.add_argument('-pfni',
'--neighbours-for-phi-field-interpolation',
type=int,
default=6,
dest='NEIGHBOURS',
help="Number of neighbours for field interpolation. Default = 6")
#parser.add_argument('-pfni',
#'--neighbours-for-phi-field-interpolation',
#type=int,
#default=6,
#dest='NEIGHBOURS',
#help="Number of neighbours for field interpolation. Default = 6")
parser.add_argument('-glt',
'--grey-low-threshold',
......@@ -400,12 +401,13 @@ def ddicParser(parser):
dest='PHIFILE_BIN_RATIO',
help="Ratio of binning level between loaded Phi file and current calculation. If the input Phi file has been obtained on a 500x500x500 image and now the calculation is on 1000x1000x1000, this should be 2. Default = 1")
parser.add_argument('-pfd',
'--phiFile-direct',
action="store_true",
#default=1,
dest='PHIFILE_DIRECT',
help="Trust the Phi file completely? This option ignores and overrides -cif and requires labels to be aligned between Phi file and labelled image. Default = False")
# PFD now strictly assumed
#parser.add_argument('-pfd',
#'--phiFile-direct',
#action="store_true",
##default=1,
#dest='PHIFILE_DIRECT',
#help="Trust the Phi file completely? This option ignores and overrides -cif and requires labels to be aligned between Phi file and labelled image. Default = False")
parser.add_argument('-cif',
'--correct-input-field',
......@@ -2468,3 +2470,185 @@ def pixelSearchPropagate(parser):
args.HWS = [args.HWS[0], args.HWS[0], args.HWS[0]]
return args
def passPhiField(parser):
parser.add_argument('-pf',
'-phiFile',
dest='PHIFILE',
default=None,
type=argparse.FileType('r'),
help="Path to TSV file containing initial F guess, can be single-point registration or multiple point correlation. Default = None")
parser.add_argument('-pfb',
'--phiFile-bin-ratio',
type=int,
default=1,
dest='PHIFILE_BIN_RATIO',
help="Ratio of binning level between loaded Phi file and current calculation. If the input Phi file has been obtained on a 500x500x500 image and now the calculation is on 1000x1000x1000, this should be 2. Default = 1")
parser.add_argument('-F',
'--apply-F',
type=str,
default='all',
dest='APPLY_F',
help="Apply the F part of Phi guess? Accepted values are:\n\t\"all\": apply all of F" +
"\n\t\"rigid\": apply rigid part (mostly rotation) \n\t\"no\": don't apply it \"all\" is default")
#parser.add_argument('-pfd',
#'--phiFile-direct',
#action="store_true",
##default=1,
#dest='PHIFILE_DIRECT',
#help="Trust the Phi file completely? This option ignores and overrides -pfni and requires same nodes in same positions. Default = False")
parser.add_argument('-lab1',
'--labelledFile1',
dest='LAB1',
default=None,
type=argparse.FileType('r'),
help="Path to tiff file containing a labelled image 1 that defines zones to correlate. Disactivates -hws and -ns options")
# Default: node spacing equal in all three directions
parser.add_argument('-ns',
'--node-spacing',
nargs=1,
type=int,
default=None,
dest='NS',
help="Node spacing in pixels (assumed equal in all 3 directions -- see -ns3 for different setting). Default = 10px")
parser.add_argument('-ns3',
'--node-spacing-3',
nargs=3,
type=int,
default=None,
dest='NS',
help="Node spacing in pixels (different in 3 directions). Default = 10, 10, 10px")
parser.add_argument('-im1',
'--image1',
dest='im1',
default=None,
type=argparse.FileType('r'),
help="Path to tiff file containing refence image, just to know the image size for the node spacing")
parser.add_argument('-im1shape',
'--image1-shape',
nargs=3,
type=int,
default=None,
dest='im1shape',
help="Size of im1 in pixels Z Y X")
parser.add_argument('-nr',
'--neighbourhood-radius',
type=float,
default=None,
dest='RADIUS',
help="Radius (in pixels) inside which to select neighbours. Default = mean(hws)+mean(sr)")
parser.add_argument('-od',
'--out-dir',
type=str,
default=None,
dest='OUT_DIR',
help="Output directory, default is the dirname of im1 file")
parser.add_argument('-pre',
'--prefix',
type=str,
default=None,
dest='PREFIX',
help='Prefix for output files (without extension). Default is basename of im1 and im2 files')
parser.add_argument('-vtk',
'--VTKout',
action="store_true",
dest='VTK',
help='Activate VTK output format. Default = False')
parser.add_argument('-tif',
'-tiff',
'--TIFFout',
'--TIFout',
action="store_true",
dest='TIFF',
help='Activate TIFFoutput format. Default = False')
args = parser.parse_args()
if args.PHIFILE is None:
print("This function definitely needs a TSV Phi file input")
exit()
# If we have no out dir specified, deliver on our default promise -- this can't be done inline before since parser.parse_args() has not been run at that stage.
if args.OUT_DIR is None:
args.OUT_DIR = os.path.dirname(args.PHIFILE.name)
# However if we have no dir, notice this and make it the current directory.
if args.OUT_DIR == "":
args.OUT_DIR = "./"
else:
# Check existence of output directory
try:
if args.OUT_DIR:
os.makedirs(args.OUT_DIR)
else:
args.DIR_out = os.path.dirname(args.PHIFILE.name)
except OSError:
if not os.path.isdir(args.OUT_DIR):
raise
# Output file name prefix
if args.PREFIX is None:
args.PREFIX = os.path.splitext(os.path.basename(args.PHIFILE.name))[0] + "-passed"
else:
args.PREFIX += "-passed"
if args.LAB1 is not None:
# We have a labelled image and so no nodeSpacing or halfWindowSize
print("I have been passed a labelled image and so I am disactivating:")
print("\t- node spacing")
args.NS = None
args.im1 = None
args.im1shape = None
# Output file name prefix
args.PREFIX += "-labelled"
else:
print("No labelled image so I'm in grid mode")
# We are in grid, with a nodeSpacing and halfWindowSize
# Catch interdependent node spacing and correlation window sizes
if args.NS is None:
print("...actually no node spacing either so, output basis not defined!")
exit()
else:
# Catch 3D options
if len(args.NS) == 1:
args.PREFIX += f"-ns{args.NS[0]}"
args.NS = [args.NS[0], args.NS[0], args.NS[0]]
else:
# 3 NSs are passed
args.PREFIX += f"-ns