Skip to content

Commit

Permalink
Merge pull request #135 from arti-sukasem/human2DAnnotation
Browse files Browse the repository at this point in the history
created lung fissure annotation
  • Loading branch information
rchristie authored May 11, 2021
2 parents 019c411 + aee2738 commit 5532f96
Show file tree
Hide file tree
Showing 3 changed files with 234 additions and 3 deletions.
3 changes: 3 additions & 0 deletions src/scaffoldmaker/annotation/lung_terms.py
Original file line number Diff line number Diff line change
Expand Up @@ -4,8 +4,11 @@

# convention: preferred name, preferred id, followed by any other ids and alternative names
lung_terms = [
("horizontal fissure of right lung", "None"),
("lung", "UBERON:0002048", "ILX:0726937"),
("left lung", "UBERON:0002168", "ILX:0733450"),
("oblique fissure of left lung", "None"),
("oblique fissure of right lung", "None"),
("right lung", "UBERON:0002167", "ILX:0729582"),
("upper lobe of left lung", "UBERON:0008952", "ILX:0735339"),
("lower lobe of left lung", "UBERON:0008953", "ILX:0735534"),
Expand Down
54 changes: 51 additions & 3 deletions src/scaffoldmaker/meshtypes/meshtype_3d_lung1.py
Original file line number Diff line number Diff line change
Expand Up @@ -2,7 +2,7 @@
Generates 3D lung surface mesh.
'''

from scaffoldmaker.annotation.annotationgroup import AnnotationGroup
from scaffoldmaker.annotation.annotationgroup import AnnotationGroup, findOrCreateAnnotationGroupForTerm, getAnnotationGroupForTerm
from scaffoldmaker.annotation.lung_terms import get_lung_term
from scaffoldmaker.meshtypes.scaffold_base import Scaffold_base
from scaffoldmaker.utils.eft_utils import remapEftLocalNodes, remapEftNodeValueLabel, setEftScaleFactorIds
Expand Down Expand Up @@ -70,8 +70,8 @@ def generateBaseMesh(cls, region, options):
:return: annotationGroups
'''
parameterSetName = options['Base parameter set']
isMouse = 'Mouse' in parameterSetName
isHuman = 'Human' in parameterSetName
isMouse = 'Mouse 1' in parameterSetName
isHuman = 'Human 1' in parameterSetName

fm = region.getFieldmodule()
coordinates = findOrCreateFieldCoordinates(fm)
Expand Down Expand Up @@ -685,6 +685,54 @@ def refineMesh(cls, meshrefinement, options):
refineElementsCount = options['Refine number of elements']
meshrefinement.refineAllElementsCubeStandard3d(refineElementsCount, refineElementsCount, refineElementsCount)

@classmethod
def defineFaceAnnotations(cls, region, options, annotationGroups):
"""
Add face annotation groups from the highest dimension mesh.
Must have defined faces and added subelements for highest dimension groups.
:param region: Zinc region containing model.
:param options: Dict containing options. See getDefaultOptions().
:param annotationGroups: List of annotation groups for top-level elements.
New face annotation groups are appended to this list.
"""
parameterSetName = options['Base parameter set']
isMouse = 'Mouse 1' in parameterSetName
isHuman = 'Human 1' in parameterSetName

# create fissure groups
fm = region.getFieldmodule()
mesh2d = fm.findMeshByDimension(2)

if isHuman:
upperLeftGroup = getAnnotationGroupForTerm(annotationGroups, get_lung_term("upper lobe of left lung"))
lowerLeftGroup = getAnnotationGroupForTerm(annotationGroups, get_lung_term("lower lobe of left lung"))

is_upperLeftGroup = upperLeftGroup.getFieldElementGroup(mesh2d)
is_lowerLeftGroup = lowerLeftGroup.getFieldElementGroup(mesh2d)

is_obliqueLeftGroup = fm.createFieldAnd(is_upperLeftGroup, is_lowerLeftGroup)
obliqueLeftGroup = findOrCreateAnnotationGroupForTerm(annotationGroups, region, get_lung_term("oblique fissure of left lung"))
obliqueLeftGroup.getMeshGroup(mesh2d).addElementsConditional(is_obliqueLeftGroup)

if isHuman or isMouse:
upperRightGroup = getAnnotationGroupForTerm(annotationGroups, get_lung_term("upper lobe of right lung"))
middleRightGroup = getAnnotationGroupForTerm(annotationGroups, get_lung_term("middle lobe of right lung"))
lowerRightGroup = getAnnotationGroupForTerm(annotationGroups, get_lung_term("lower lobe of right lung"))

is_upperRightGroup = upperRightGroup.getFieldElementGroup(mesh2d)
is_middleRightGroup = middleRightGroup.getFieldElementGroup(mesh2d)
is_lowerRightGroup = lowerRightGroup.getFieldElementGroup(mesh2d)

is_obliqueRightGroup = fm.createFieldAnd(fm.createFieldOr(is_middleRightGroup, is_upperRightGroup),
is_lowerRightGroup)
is_horizontalRightGroup = fm.createFieldAnd(is_upperRightGroup, is_middleRightGroup)

obliqueRightGroup = findOrCreateAnnotationGroupForTerm(annotationGroups, region, get_lung_term("oblique fissure of right lung"))
obliqueRightGroup.getMeshGroup(mesh2d).addElementsConditional(is_obliqueRightGroup)
horizontalRightGroup = findOrCreateAnnotationGroupForTerm(annotationGroups, region, get_lung_term("horizontal fissure of right lung"))
horizontalRightGroup.getMeshGroup(mesh2d).addElementsConditional(is_horizontalRightGroup)


def getLungNodes(lungSide, cache, coordinates, generateParameters, nodes, nodetemplate, nodeFieldParameters,
lElementsCount1, lElementsCount2, lElementsCount3,
uElementsCount1, uElementsCount2, uElementsCount3,
Expand Down
180 changes: 180 additions & 0 deletions tests/test_lung.py
Original file line number Diff line number Diff line change
@@ -0,0 +1,180 @@
import unittest
import copy
from opencmiss.utils.zinc.finiteelement import evaluateFieldNodesetRange, findNodeWithName
from opencmiss.utils.zinc.general import ChangeManager
from opencmiss.zinc.context import Context
from opencmiss.zinc.field import Field
from opencmiss.zinc.result import RESULT_OK
from scaffoldmaker.annotation.annotationgroup import AnnotationGroup, getAnnotationGroupForTerm
from scaffoldmaker.annotation.lung_terms import get_lung_term
from scaffoldmaker.meshtypes.meshtype_3d_lung1 import MeshType_3d_lung1
from scaffoldmaker.utils.meshrefinement import MeshRefinement
from testutils import assertAlmostEqualList

class LungScaffoldTestCase(unittest.TestCase):

def test_lung1(self):
"""
Test creation of heart scaffold.
"""
scaffold = MeshType_3d_lung1
parameterSetNames = scaffold.getParameterSetNames()
self.assertEqual(parameterSetNames, [ "Default", "Human 1", "Mouse 1" ])
options = scaffold.getDefaultOptions(["Human"])
self.assertEqual(3, len(options))
self.assertFalse(scaffold.checkOptions(options))

context = Context("Test")
region = context.getDefaultRegion()
self.assertTrue(region.isValid())
annotationGroups = scaffold.generateMesh(region, options)
self.assertEqual(11, len(annotationGroups))
fieldmodule = region.getFieldmodule()
mesh3d = fieldmodule.findMeshByDimension(3)
self.assertEqual(88, mesh3d.getSize())
mesh2d = fieldmodule.findMeshByDimension(2)
self.assertEqual(292, mesh2d.getSize())
mesh1d = fieldmodule.findMeshByDimension(1)
self.assertEqual(334, mesh1d.getSize())
nodes = fieldmodule.findNodesetByFieldDomainType(Field.DOMAIN_TYPE_NODES)
self.assertEqual(138, nodes.getSize())
datapoints = fieldmodule.findNodesetByFieldDomainType(Field.DOMAIN_TYPE_DATAPOINTS)
self.assertEqual(0, datapoints.getSize())

# check coordinates range, outside surface area and volume
coordinates = fieldmodule.findFieldByName("coordinates").castFiniteElement()
self.assertTrue(coordinates.isValid())
minimums, maximums = evaluateFieldNodesetRange(coordinates, nodes)
assertAlmostEqualList(self, minimums, [ 14.153, 62.444, -354.889 ], 1.0E-6)
assertAlmostEqualList(self, maximums, [ 297.502, 234.911, -23.486 ], 1.0E-6)
with ChangeManager(fieldmodule):
one = fieldmodule.createFieldConstant(1.0)
upperRightFissureGroup = fieldmodule.findFieldByName('upper lobe of right lung').castGroup()
self.assertTrue(upperRightFissureGroup.isValid())
upperRightFissureMeshGroup = upperRightFissureGroup.getFieldElementGroup(mesh2d).getMeshGroup()
self.assertTrue(upperRightFissureMeshGroup.isValid())
surfaceAreaField = fieldmodule.createFieldMeshIntegral(one, coordinates, upperRightFissureMeshGroup)
surfaceAreaField.setNumbersOfPoints(4)
volumeField = fieldmodule.createFieldMeshIntegral(one, coordinates, mesh3d)
volumeField.setNumbersOfPoints(3)
fieldcache = fieldmodule.createFieldcache()
result, surfaceArea = surfaceAreaField.evaluateReal(fieldcache, 1)
self.assertEqual(result, RESULT_OK)
self.assertAlmostEqual(surfaceArea, 139195.03021771502, delta=1.0E-2)
result, volume = volumeField.evaluateReal(fieldcache, 1)
self.assertEqual(result, RESULT_OK)
self.assertAlmostEqual(volume, 7169261.158690422, delta=1.0E-2)

# check some annotationGroups:
expectedSizes3d = {
"lung": 88,
"left lung": 44,
"right lung": 44,
"upper lobe of left lung": 24,
"lower lobe of left lung": 20,
"upper lobe of right lung": 16,
"middle lobe of right lung": 8,
"lower lobe of right lung": 20,
# "right lung accessory lobe": 0
}
for name in expectedSizes3d:
group = getAnnotationGroupForTerm(annotationGroups, get_lung_term(name))
size = group.getMeshGroup(mesh3d).getSize()
self.assertEqual(expectedSizes3d[name], size, name)
expectedSizes2d = {
"horizontal fissure of right lung": 4,
"oblique fissure of left lung": 8,
"oblique fissure of right lung": 8
}
for name in expectedSizes2d:
group = getAnnotationGroupForTerm(annotationGroups, get_lung_term(name))
size = group.getMeshGroup(mesh2d).getSize()
self.assertEqual(expectedSizes2d[name], size, name)

# test finding a marker in scaffold
markerGroup = fieldmodule.findFieldByName("marker").castGroup()
markerNodes = markerGroup.getFieldNodeGroup(nodes).getNodesetGroup()
self.assertEqual(6, markerNodes.getSize())
markerName = fieldmodule.findFieldByName("marker_name")
self.assertTrue(markerName.isValid())
markerLocation = fieldmodule.findFieldByName("marker_location")
self.assertTrue(markerLocation.isValid())
# test apex marker point
cache = fieldmodule.createFieldcache()
node = findNodeWithName(markerNodes, markerName, "apex of left lung")
self.assertTrue(node.isValid())
cache.setNode(node)
element, xi = markerLocation.evaluateMeshLocation(cache, 3)
self.assertEqual(37, element.getIdentifier())
assertAlmostEqualList(self, xi, [ 0.0, 0.0, 1.0 ], 1.0E-10)

# refine 4 and check result
# first remove any face (but not point) annotation groups as they are re-added by defineFaceAnnotations
removeAnnotationGroups = []
for annotationGroup in annotationGroups:
if (not annotationGroup.hasMeshGroup(mesh3d)) and (annotationGroup.hasMeshGroup(mesh2d) or annotationGroup.hasMeshGroup(mesh1d)):
removeAnnotationGroups.append(annotationGroup)
for annotationGroup in removeAnnotationGroups:
annotationGroups.remove(annotationGroup)
self.assertEqual(8, len(annotationGroups))

refineRegion = region.createRegion()
refineFieldmodule = refineRegion.getFieldmodule()
options['Refine'] = True
options['Refine number of elements'] = 4
refineNumberOfElements = options['Refine number of elements']
meshrefinement = MeshRefinement(region, refineRegion, annotationGroups)
scaffold.refineMesh(meshrefinement, options)
annotationGroups = meshrefinement.getAnnotationGroups()

refineFieldmodule.defineAllFaces()
oldAnnotationGroups = copy.copy(annotationGroups)
for annotationGroup in annotationGroups:
annotationGroup.addSubelements()
scaffold.defineFaceAnnotations(refineRegion, options, annotationGroups)
for annotation in annotationGroups:
if annotation not in oldAnnotationGroups:
annotationGroup.addSubelements()
self.assertEqual(11, len(annotationGroups))

mesh3d = refineFieldmodule.findMeshByDimension(3)
self.assertEqual(5632, mesh3d.getSize())
mesh2d = refineFieldmodule.findMeshByDimension(2)
self.assertEqual(17344, mesh2d.getSize())
mesh1d = refineFieldmodule.findMeshByDimension(1)
self.assertEqual(17848, mesh1d.getSize())
nodes = refineFieldmodule.findNodesetByFieldDomainType(Field.DOMAIN_TYPE_NODES)
self.assertEqual(6144, nodes.getSize())
datapoints = refineFieldmodule.findNodesetByFieldDomainType(Field.DOMAIN_TYPE_DATAPOINTS)
self.assertEqual(0, datapoints.getSize())

# check some refined annotationGroups:
for name in expectedSizes3d:
group = getAnnotationGroupForTerm(annotationGroups, get_lung_term(name))
size = group.getMeshGroup(mesh3d).getSize()
self.assertEqual(expectedSizes3d[name]*(refineNumberOfElements**3), size, name)
for name in expectedSizes2d:
group = getAnnotationGroupForTerm(annotationGroups, get_lung_term(name))
size = group.getMeshGroup(mesh2d).getSize()
self.assertEqual(expectedSizes2d[name]*(refineNumberOfElements**2), size, name)

# test finding a marker in refined scaffold
markerGroup = refineFieldmodule.findFieldByName("marker").castGroup()
refinedNodes = refineFieldmodule.findNodesetByFieldDomainType(Field.DOMAIN_TYPE_NODES)
markerNodes = markerGroup.getFieldNodeGroup(refinedNodes).getNodesetGroup()
self.assertEqual(6, markerNodes.getSize())
markerName = refineFieldmodule.findFieldByName("marker_name")
self.assertTrue(markerName.isValid())
markerLocation = refineFieldmodule.findFieldByName("marker_location")
self.assertTrue(markerLocation.isValid())
# test apex marker point
cache = refineFieldmodule.createFieldcache()
node = findNodeWithName(markerNodes, markerName, "apex of left lung")
self.assertTrue(node.isValid())
cache.setNode(node)
element, xi = markerLocation.evaluateMeshLocation(cache, 3)
self.assertEqual(2353, element.getIdentifier())
assertAlmostEqualList(self, xi, [ 0.0, 0.0, 1.0 ], 1.0E-10)

if __name__ == "__main__":
unittest.main()

0 comments on commit 5532f96

Please sign in to comment.