Commit 0348adb5 authored by Benjy Marks's avatar Benjy Marks
Browse files

Merge branch 'master' of gricad-gitlab.univ-grenoble-alpes.fr:ttk/radioSphere

parents 32ecee84 ed210da4
[report]
exclude_lines =
# no cover
if debug:
if graphShow:
if GRAPH:
[run]
omit =
*/_*
source =
radioSphere
concurrency = multiprocessing
[html]
directory = coverage
......@@ -7,6 +7,7 @@ test:
- DEBIAN_FRONTEND=noninteractive apt-get install -y build-essential gcc cmake python3-pip python3 python3-dev libeigen3-dev
- pip3 install -U -r requirements.txt
- python3 setup.py install
#- python setup.py test
- sphinx-build -b html docs/source public
only:
- branches
......@@ -23,6 +24,7 @@ pages:
- sphinx-build -b html docs/source public
# - conda build .
# - tar --exclude='.[^/]*' -czvf public/radioSphere.tar.gz .
# -pytest --cov-report html --cov=radioSphere tests/
artifacts:
paths:
......
......@@ -9,7 +9,7 @@ from scipy.spatial import distance
import numpy
import tifffile
plt.style.use('./tools/radioSphere.mplstyle')
#plt.style.use('./tools/radioSphere.mplstyle')
zoomLevel = 10
sourceObjectDistMM = 20.0
......@@ -127,4 +127,4 @@ dumpfolder = '~/code/NDDEM/Samples/radioSphere/'
# write_dumpfile(positionsXYZmmClean,radiiMM,2,dumpfolder=dumpfolder)
plt.show()
#plt.show()
......@@ -7,18 +7,65 @@ If a new function is added it should be called at the bottom of the function
"""
import radioSphere
import numpy
import unittest
import matplotlib.pyplot as plt
import radioSphere
import shutil
class TestDetection(unittest.TestCase):
def tearDown(self):
try:
shutil.rmtree("./output")
pass
except OSError:
pass
def test_fxToPeaks(self):
# 2D fx
fx = numpy.zeros((200, 100))
fx[40, 70] = 0.5
peaks = radioSphere.detectSpheres.indicatorFunctionToDetectorPositions(fx, scanFixedNumber=1)
self.assertTrue(numpy.allclose([0, 40, 70], peaks))
# 3D fx
fx = numpy.zeros((30, 200, 100))
fx[5, 40, 70] = 0.5
peaksa = radioSphere.detectSpheres.indicatorFunctionToDetectorPositions(fx, scanFixedNumber=1)
peaksb = radioSphere.detectSpheres.indicatorFunctionToDetectorPositions(fx, massThreshold=0.4)
self.assertTrue(numpy.allclose([5, 40, 70], peaksa))
self.assertTrue(numpy.allclose([5, 40, 70], peaksb))
def test_peaksTo3DPos(self):
# peaks at middle of detector so real pos should be at 0 of coord system
peaks = numpy.array([[2, 100, 50]])
CORxPos = numpy.arange(3)
posXYZ = radioSphere.detectSpheres.detectorPeaksTo3DCoordinates(peaks, CORxPos, detectorResolution=[200, 100])
self.assertTrue(numpy.allclose(posXYZ, [2, 0, 0]))
# backwards compatibility with projectSphere
detectorResolution=[200, 100]
sourceDetectorDistMM = 50
pixelSizeMM = 0.5
peaks = [2, 20, 80]
CORxPos = numpy.arange(3)
posXYZ = radioSphere.detectSpheres.detectorPeaksTo3DCoordinates(numpy.array([peaks]), CORxPos,
pixelSizeMM=pixelSizeMM,
detectorResolution=detectorResolution,
sourceDetectorDistMM=sourceDetectorDistMM)
p = radioSphere.projectSphere.projectSphereMM(posXYZ,
numpy.ones(posXYZ.shape[0]),
detectorResolution=[200, 100],
sourceDetectorDistMM=sourceDetectorDistMM,
pixelSizeMM=pixelSizeMM)
sphereCentreDetector = numpy.argwhere(p == p.max())[0]
# 1px tol
self.assertTrue(numpy.allclose(peaks[1:], sphereCentreDetector, atol=1))
def test_tomopackOneSphere(self):
COR = 50
r = 1
......@@ -35,8 +82,8 @@ class TestDetection(unittest.TestCase):
detectorResolution = detectorResolution,
sourceDetectorDistMM = numpy.inf)
f_x1 = radioSphere.detectSpheres.tomopack(p1, psi)
peaksJI1 = radioSphere.detectSpheres.indicatorFunctionToDetectorPositions(f_x1, scanFixedNumber=1)
fx1 = radioSphere.detectSpheres.tomopack(p1, psi)
peaksJI1 = radioSphere.detectSpheres.indicatorFunctionToDetectorPositions(fx1, scanFixedNumber=1)[:, 1:]
self.assertEqual(peaksJI1[0][0], detectorResolution[0]//2)
self.assertEqual(peaksJI1[0][1], detectorResolution[1]//2)
......@@ -47,17 +94,17 @@ class TestDetection(unittest.TestCase):
detectorResolution = detectorResolution,
sourceDetectorDistMM = numpy.inf)
f_x2 = radioSphere.detectSpheres.tomopack(p2, psi)
peaksJI2 = radioSphere.detectSpheres.indicatorFunctionToDetectorPositions(f_x2, massThreshold=0.1)
fx2 = radioSphere.detectSpheres.tomopack(p2, psi)
peaksJI2 = radioSphere.detectSpheres.indicatorFunctionToDetectorPositions(fx2, massThreshold=0.1)[:, 1:]
spherePosDetector2 = numpy.where(p2 == p2.max())
estimatedSpherePosDetector2 = numpy.where(f_x2 == f_x2.max())
estimatedSpherePosDetector2 = numpy.where(fx2 == fx2.max())
# 1px max diff
self.assertTrue(abs(spherePosDetector2[0]-peaksJI2[0][0]) <= 1)
self.assertTrue(abs(spherePosDetector2[1]-peaksJI2[0][1]) <= 1)
# Case 3: divergent beam, simplest case with 1 sphere in the middle of the detector
# Case 3: A bit of zoom, simplest case with 1 sphere in the middle of the detector
p3 = radioSphere.projectSphere.projectSphereMM(numpy.array([[COR, 0, 0]]),
numpy.array([r]),
detectorResolution = detectorResolution,
......@@ -68,27 +115,27 @@ class TestDetection(unittest.TestCase):
detectorResolution = detectorResolution,
sourceDetectorDistMM = 100)
f_x3 = radioSphere.detectSpheres.tomopack(p3, psi3)
peaksJI3 = radioSphere.detectSpheres.indicatorFunctionToDetectorPositions(f_x3, scanFixedNumber=1)
fx3 = radioSphere.detectSpheres.tomopack(p3, psi3)
peaksJI3 = radioSphere.detectSpheres.indicatorFunctionToDetectorPositions(fx3, scanFixedNumber=1)[:, 1:]
self.assertEqual(peaksJI3[0][0], detectorResolution[0]//2)
self.assertEqual(peaksJI3[0][1], detectorResolution[1]//2)
# Case 4: divergent beam, 1 sphere offset from the middle
# Case 4: A bit of zoom, 1 sphere offset from the middle
p4 = radioSphere.projectSphere.projectSphereMM(numpy.array([[COR, -1.2, -2.6]]),
numpy.array([r]),
detectorResolution = detectorResolution,
sourceDetectorDistMM = 100)
f_x4 = radioSphere.detectSpheres.tomopack(p4, psi3)
peaksJI4 = radioSphere.detectSpheres.indicatorFunctionToDetectorPositions(f_x4, scanFixedNumber=1)
fx4 = radioSphere.detectSpheres.tomopack(p4, psi3)
peaksJI4 = radioSphere.detectSpheres.indicatorFunctionToDetectorPositions(fx4, scanFixedNumber=1)[:, 1:]
spherePosDetector4 = numpy.where(p4 == p4.max())
# 1px max diff
self.assertTrue(abs(spherePosDetector4[0][0]-peaksJI4[0][0]) <= 1)
self.assertTrue(abs(spherePosDetector4[1][0]-peaksJI4[0][1]) <= 1)
# tomopack 10 random spheres with divergent beam
# Tomopack 10 random spheres
def test_tomopack10Spheres(self):
COR = 50
r = 1
......@@ -118,25 +165,94 @@ class TestDetection(unittest.TestCase):
sourceDetectorDistMM = sourceDetectorDistMM)
f_x = radioSphere.detectSpheres.tomopack(p, psi)
peaksJI = radioSphere.detectSpheres.indicatorFunctionToDetectorPositions(f_x, scanFixedNumber=pos.shape[0])
posTomopackXYZ = numpy.zeros([peaksJI.shape[0], 3])
posTomopackXYZ[:, 0] += COR
fx = radioSphere.detectSpheres.tomopack(p, psi)
# TODO: this should be a separate helper function for both a divergent and a parallel beam
for i in range(posTomopackXYZ.shape[0]):
# detector I gives real position Y in mm and scale by zoom
posTomopackXYZ[i, 1] = -1*(peaksJI[i,1] - p.shape[1]/2.0)*pixelSizeMM
posTomopackXYZ[i, 1] /= sourceDetectorDistMM/COR
# call + test helper which transforms fx guess to detector peaks
peaksJI = radioSphere.detectSpheres.indicatorFunctionToDetectorPositions(fx, scanFixedNumber=pos.shape[0])
# detector J gives real position Z in mm and scale by zoom
posTomopackXYZ[i, 2] = -1*(peaksJI[i,0] - p.shape[0]/2.0)*pixelSizeMM
posTomopackXYZ[i, 2] /= sourceDetectorDistMM/COR
# call + test helper which transforms det peaks to 3d coord
posTomopackXYZ = radioSphere.detectSpheres.detectorPeaksTo3DCoordinates(peaksJI,
numpy.ones(peaksJI.shape[0])*COR,
pixelSizeMM=pixelSizeMM,
detectorResolution=p.shape,
sourceDetectorDistMM=sourceDetectorDistMM
)
# sort tomopack guess to calculate residual directly on positions
posTomopackXYZ = posTomopackXYZ[numpy.argsort(posTomopackXYZ[:, 1])]
self.assertTrue(numpy.linalg.norm(pos - posTomopackXYZ) < r)
self.assertTrue(numpy.mean(numpy.linalg.norm(pos - posTomopackXYZ, axis=1) < 0.5*r))
def test_tomopackDivergentBeamOneSphere(self):
r = 1
COR = 50
detectorResolution = [400, 600]
# Case 1: simplest case with 1 sphere in the middle of the detector
pos1 = numpy.array([[COR, 0, 0]])
p1 = radioSphere.projectSphere.projectSphereMM(pos1,
numpy.array([r]),
detectorResolution = detectorResolution,
sourceDetectorDistMM = 100)
posXYZ1a = radioSphere.detectSpheres.tomopackDivergentScanTo3DPositions(p1, r, CORxRef=COR, CORxNumber=5, scanFixedNumber=1, saveSeries=False, verbose=False)
posXYZ1b = radioSphere.detectSpheres.tomopackDivergentScanTo3DPositions(p1, r, CORxRef=COR, CORxNumber=5, saveSeries=False, massThreshold=0.2)
self.assertTrue(numpy.allclose(pos1, posXYZ1a[0], atol=0.1))
self.assertTrue(numpy.allclose(pos1, posXYZ1b[0], atol=0.1))
# Case 2: 1 sphere shifted from COR and centre of detector
pos2 = numpy.array([[numpy.random.uniform(46, 48), numpy.random.uniform(5, -5), numpy.random.uniform(-3, 3)]])
#pos2 = numpy.array([[47.20396842, 0, 0]])
#pos2 = numpy.array([[46.22945869, -1.63507697, 0.96052506]])
p2 = radioSphere.projectSphere.projectSphereMM(pos2,
numpy.array([r]),
detectorResolution= detectorResolution,
sourceDetectorDistMM=100)
posXYZ2a = radioSphere.detectSpheres.tomopackDivergentScanTo3DPositions(p2, r, CORxMin=COR-4*r, CORxMax=COR+4*r, CORxNumber=20, scanFixedNumber=1, saveSeries=False, verbose=False)
posXYZ2b = radioSphere.detectSpheres.tomopackDivergentScanTo3DPositions(p2, r, CORxMin=COR-4*r, CORxMax=COR+4*r, CORxNumber=20, saveSeries=False, verbose=False, massThreshold=0.1)
self.assertTrue(numpy.allclose(pos2, posXYZ2a[0], atol=0.5))
self.assertTrue(numpy.allclose(pos2, posXYZ2b[0], atol=0.5))
def test_tomopackDivergentBeam10Spheres(self):
r = 1
COR = 50
detectorResolution = [400, 600]
pixelSizeMM = 0.1
xMin, xMax = 40, 60
yMin, yMax = -5, 5
zMin, zMax = -7, 7
pos = numpy.zeros((10, 3))
pos[:, 0] = numpy.arange(xMin, xMax, 2*r)
pos[:, 1] = yMin + r + (yMax - yMin - 2*r)*numpy.random.rand(10)
pos[:, 2] = zMin + r + (zMax - zMin - 2*r)*numpy.random.rand(10)
radii = numpy.ones_like(pos[:, 0])*r
p = radioSphere.projectSphere.projectSphereMM(pos,
numpy.ones_like(pos[:, 0])*r,
pixelSizeMM=pixelSizeMM,
detectorResolution= detectorResolution,
sourceDetectorDistMM=100)
posXYZa = radioSphere.detectSpheres.tomopackDivergentScanTo3DPositions(p, r,
CORxMin=COR-11*r, CORxMax=COR+11*r, CORxNumber=50,
scanFixedNumber=pos.shape[0],
saveSeries=True, saveSeriesDirectory="./output")
posXYZb = radioSphere.detectSpheres.tomopackDivergentScanTo3DPositions(p, r, CORxMin=COR-11*r, CORxMax=COR+11*r, CORxNumber=50,
fXseries = "./output/fXseries.tif", psiXseries = "./output/psiXseries.tif",
massThreshold=0.1,
saveSeries=False)
# sort guess based on X axis to calculate residual directly on position
posXYZa = posXYZa[numpy.argsort(posXYZa[:, 0])]
posXYZb = posXYZb[numpy.argsort(posXYZb[:, 0])]
self.assertTrue(numpy.mean(numpy.linalg.norm(pos - posXYZa, axis=1)) < 1.5*r)
self.assertTrue(numpy.mean(numpy.linalg.norm(pos - posXYZb, axis=1)) < 1.5*r)
if __name__ == "__main__":
unittest.main()
import os
os.environ['OPENBLAS_NUM_THREADS'] = '1'
os.environ['OMP_NUM_THREADS'] = '1'
__all__ = ["optimisePositions", "projectSphere", "DEM", "detectSpheres", "calibrateAttenuation"]
#from .blender_viewer import *
......
This diff is collapsed.
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