Skip to content

Commit

Permalink
Merge pull request #184 from rchristie/user_markers
Browse files Browse the repository at this point in the history
Add user marker point annotation groups
  • Loading branch information
rchristie authored Mar 25, 2022
2 parents 54932e5 + eaab998 commit 5e3b453
Show file tree
Hide file tree
Showing 10 changed files with 734 additions and 161 deletions.
354 changes: 338 additions & 16 deletions src/scaffoldmaker/annotation/annotationgroup.py

Large diffs are not rendered by default.

152 changes: 55 additions & 97 deletions src/scaffoldmaker/meshtypes/meshtype_3d_brainstem.py
Original file line number Diff line number Diff line change
Expand Up @@ -11,6 +11,7 @@
from opencmiss.zinc.field import FieldFindMeshLocation
from opencmiss.utils.zinc.field import Field, findOrCreateFieldCoordinates, findOrCreateFieldGroup, findOrCreateFieldNodeGroup, findOrCreateFieldStoredMeshLocation, \
findOrCreateFieldStoredString
from opencmiss.utils.zinc.general import ChangeManager
from opencmiss.utils.zinc.finiteelement import getMaximumNodeIdentifier
from scaffoldmaker.annotation.annotationgroup import AnnotationGroup, findOrCreateAnnotationGroupForTerm
from scaffoldmaker.annotation.brainstem_terms import get_brainstem_term
Expand All @@ -37,7 +38,7 @@ class MeshType_3d_brainstem1(Scaffold_base):
'D2 derivatives': True,
'D3 derivatives': True,
'Length': 3.0,
'Number of elements': 6
'Number of elements': 8
},
'meshEdits': exnodeStringFromNodeValues( # dimensional.
[Node.VALUE_LABEL_VALUE, Node.VALUE_LABEL_D_DS1, Node.VALUE_LABEL_D_DS2, Node.VALUE_LABEL_D2_DS1DS2,
Expand All @@ -56,21 +57,21 @@ class MeshType_3d_brainstem1(Scaffold_base):
{
'_AnnotationGroup': True,
'dimension': 1,
'identifierRanges': '1',
'identifierRanges': '1-3',
'name': get_brainstem_term('medulla oblongata')[0],
'ontId': get_brainstem_term('medulla oblongata')[1]
},
{
'_AnnotationGroup': True,
'dimension': 1,
'identifierRanges': '2',
'identifierRanges': '4-6',
'name': get_brainstem_term('pons')[0],
'ontId': get_brainstem_term('pons')[1]
},
{
'_AnnotationGroup': True,
'dimension': 1,
'identifierRanges': '3',
'identifierRanges': '7-8',
'name': get_brainstem_term('midbrain')[0],
'ontId': get_brainstem_term('midbrain')[1]
}]
Expand Down Expand Up @@ -103,21 +104,21 @@ class MeshType_3d_brainstem1(Scaffold_base):
{
'_AnnotationGroup': True,
'dimension': 1,
'identifierRanges': '1',
'identifierRanges': '1-2',
'name': get_brainstem_term('medulla oblongata')[0],
'ontId': get_brainstem_term('medulla oblongata')[1]
},
{
'_AnnotationGroup': True,
'dimension': 1,
'identifierRanges': '2',
'identifierRanges': '3-4',
'name': get_brainstem_term('pons')[0],
'ontId': get_brainstem_term('pons')[1]
},
{
'_AnnotationGroup': True,
'dimension': 1,
'identifierRanges': '3',
'identifierRanges': '5-6',
'name': get_brainstem_term('midbrain')[0],
'ontId': get_brainstem_term('midbrain')[1]
}]
Expand Down Expand Up @@ -501,37 +502,18 @@ def generateBaseMesh(cls, region, options):
annotationGroupAlong = [[brainstemGroup, midbrainGroup],
[brainstemGroup, ponsGroup],
[brainstemGroup, medullaGroup]]

# point markers
# centralCanal = findOrCreateAnnotationGroupForTerm(annotationGroups, region,
# get_brainstem_term('central canal of spinal cord'))
# cerebralAqueduct = findOrCreateAnnotationGroupForTerm(annotationGroups, region,
# get_brainstem_term('cerebral aqueduct'))
# foramenCaecum = findOrCreateAnnotationGroupForTerm(annotationGroups, region,
# get_brainstem_term('foramen caecum of medulla oblongata'))
dorsalMidCaudalGroup = findOrCreateAnnotationGroupForTerm(annotationGroups, region,
get_brainstem_term('brainstem dorsal midline caudal point'))
ventralMidCaudalGroup = findOrCreateAnnotationGroupForTerm(annotationGroups, region,
get_brainstem_term('brainstem ventral midline caudal point'))
dorsalMidCranGroup = findOrCreateAnnotationGroupForTerm(annotationGroups, region,
get_brainstem_term('brainstem dorsal midline cranial point'))
ventralMidCranGroup = findOrCreateAnnotationGroupForTerm(annotationGroups, region,
get_brainstem_term('brainstem ventral midline cranial point'))
dorsalMidMedullaPonsJunction = findOrCreateAnnotationGroupForTerm(annotationGroups, region,
get_brainstem_term('brainstem dorsal midline pons-medulla junction'))
ventralMidMedullaPonsJunction = findOrCreateAnnotationGroupForTerm(annotationGroups, region,
get_brainstem_term('brainstem ventral midline pons-medulla junction'))
dorsalMidMidbrainPonsJunction = findOrCreateAnnotationGroupForTerm(annotationGroups, region,
get_brainstem_term('brainstem dorsal midline midbrain-pons junction'))
ventralMidMidbrainPonsJunction = findOrCreateAnnotationGroupForTerm(annotationGroups, region,
get_brainstem_term('brainstem ventral midline midbrain-pons junction'))

#######################
# CREATE MAIN BODY MESH
#######################
cylinderShape = CylinderShape.CYLINDER_SHAPE_FULL if not halfBrainStem else CylinderShape.CYLINDER_SHAPE_LOWER_HALF

# Body coordinates
# brainstem coordinates
cylinderCentralPath = CylinderCentralPath(region, centralPath, elementsCountAlong)
base = CylinderEnds(elementsCountAcrossMajor, elementsCountAcrossMinor,
centre=[0.0, 0.0, 0.0],
Expand All @@ -544,40 +526,43 @@ def generateBaseMesh(cls, region, options):
cylinderCentralPath=cylinderCentralPath,
useCrossDerivatives=False)

# workaround for old Zinc field wrapper bug: must create brainstem coordinates field before reading file
brainstem_coordinates = findOrCreateFieldCoordinates(fm, name="brainstem coordinates")

# Brain coordinates
# generate brainstem coordinates field in temporary region
tmp_region = region.createRegion()
tmp_fm = tmp_region.getFieldmodule()
tmp_brainstem_coordinates = findOrCreateFieldCoordinates(tmp_fm, name="brainstem coordinates")

cylinderCentralPath1 = CylinderCentralPath(tmp_region, brainstemPath, elementsCountAlong)

base1 = CylinderEnds(elementsCountAcrossMajor, elementsCountAcrossMinor,
centre=[0.0, 0.0, 0.0],
alongAxis=cylinderCentralPath1.alongAxis[0],
majorAxis=cylinderCentralPath1.majorAxis[0],
minorRadius=cylinderCentralPath1.minorRadii[0])

cylinder2 = CylinderMesh(tmp_fm, tmp_brainstem_coordinates, elementsCountAlong, base1,
cylinderShape=cylinderShape,
cylinderCentralPath=cylinderCentralPath1,
useCrossDerivatives=False)

# Write two coordinates
sir = tmp_region.createStreaminformationRegion()
srm = sir.createStreamresourceMemory()
tmp_region.write(sir)
result, buffer = srm.getBuffer()

sir = region.createStreaminformationRegion()
srm = sir.createStreamresourceMemoryBuffer(buffer)
region.read(sir)

del srm
del sir
with ChangeManager(tmp_fm):
tmp_brainstem_coordinates = findOrCreateFieldCoordinates(tmp_fm, name="brainstem coordinates")

cylinderCentralPath1 = CylinderCentralPath(tmp_region, brainstemPath, elementsCountAlong)

base1 = CylinderEnds(elementsCountAcrossMajor, elementsCountAcrossMinor,
centre=[0.0, 0.0, 0.0],
alongAxis=cylinderCentralPath1.alongAxis[0],
majorAxis=cylinderCentralPath1.majorAxis[0],
minorRadius=cylinderCentralPath1.minorRadii[0])

cylinder2 = CylinderMesh(tmp_fm, tmp_brainstem_coordinates, elementsCountAlong, base1,
cylinderShape=cylinderShape,
cylinderCentralPath=cylinderCentralPath1,
useCrossDerivatives=False)

# write to memory buffer
sir = tmp_region.createStreaminformationRegion()
srm = sir.createStreamresourceMemory()
tmp_region.write(sir)
result, buffer = srm.getBuffer()

# read into main region
sir = region.createStreaminformationRegion()
srm = sir.createStreamresourceMemoryBuffer(buffer)
region.read(sir)

del srm
del sir
del tmp_brainstem_coordinates
del tmp_fm
del tmp_brainstem_coordinates
del tmp_region

# Annotating groups
Expand All @@ -596,49 +581,22 @@ def generateBaseMesh(cls, region, options):
################
# point markers
################
pointMarkers = [
{"group": dorsalMidCaudalGroup, "marker_brainstem_coordinates": [0.0, 1.0, 0.0]},
{"group": ventralMidCaudalGroup, "marker_brainstem_coordinates": [0.0, -1.0, 0.0]},
{"group": dorsalMidCranGroup, "marker_brainstem_coordinates": [0.0, 1.0, 8.0]},
{"group": ventralMidCranGroup, "marker_brainstem_coordinates": [0.0, -1.0, 8.0]},
{"group": dorsalMidMedullaPonsJunction, "marker_brainstem_coordinates": [0.0, 1.0, 3.0]},
{"group": ventralMidMedullaPonsJunction, "marker_brainstem_coordinates": [0.0, -1.0, 3.0]},
{"group": dorsalMidMidbrainPonsJunction, "marker_brainstem_coordinates": [0.0, 1.0, 6.0]},
{"group": ventralMidMidbrainPonsJunction, "marker_brainstem_coordinates": [0.0, -1.0, 6.0]},
]

markerGroup = findOrCreateFieldGroup(fm, "marker")
markerName = findOrCreateFieldStoredString(fm, name="marker_name")
markerLocation = findOrCreateFieldStoredMeshLocation(fm, mesh, name="marker_location")

markerTermNameBrainstemCoordinatesMap = {
'brainstem dorsal midline caudal point': [0.0, 1.0, 0.0],
'brainstem ventral midline caudal point': [0.0, -1.0, 0.0],
'brainstem dorsal midline cranial point': [0.0, 1.0, 8.0],
'brainstem ventral midline cranial point': [0.0, -1.0, 8.0],
'brainstem dorsal midline pons-medulla junction': [0.0, 1.0, 3.0],
'brainstem ventral midline pons-medulla junction': [0.0, -1.0, 3.0],
'brainstem dorsal midline midbrain-pons junction': [0.0, 1.0, 6.0],
'brainstem ventral midline midbrain-pons junction': [0.0, -1.0, 6.0]
}
nodes = fm.findNodesetByFieldDomainType(Field.DOMAIN_TYPE_NODES)
markerPoints = findOrCreateFieldNodeGroup(markerGroup, nodes).getNodesetGroup()
markerBrainstemCoordinates = findOrCreateFieldCoordinates(fm, name="marker_body_coordinates")
markerTemplateInternal = nodes.createNodetemplate()
markerTemplateInternal.defineField(markerName)
markerTemplateInternal.defineField(markerLocation)
markerTemplateInternal.defineField(markerBrainstemCoordinates)

cache = fm.createFieldcache()

brainstemNodesetGroup = brainstemGroup.getNodesetGroup(nodes)

nodeIdentifier = max(1, getMaximumNodeIdentifier(nodes) + 1)
findMarkerLocation = fm.createFieldFindMeshLocation(markerBrainstemCoordinates, brainstem_coordinates, mesh)
findMarkerLocation.setSearchMode(FieldFindMeshLocation.SEARCH_MODE_EXACT)

for pointMarker in pointMarkers:
group = pointMarker["group"]
markerPoint = markerPoints.createNode(nodeIdentifier, markerTemplateInternal)
cache.setNode(markerPoint)

markerBrainstemCoordinates.assignReal(cache, pointMarker["marker_brainstem_coordinates"])
markerName.assignString(cache, group.getName())

element, xi = findMarkerLocation.evaluateMeshLocation(cache, 3)
markerLocation.assignMeshLocation(cache, element, xi)
group.getNodesetGroup(nodes).addNode(markerPoint)
brainstemNodesetGroup.addNode(markerPoint)
for termName, brainstemCoordinatesValues in markerTermNameBrainstemCoordinatesMap.items():
annotationGroup = findOrCreateAnnotationGroupForTerm(annotationGroups, region, get_brainstem_term(termName))
annotationGroup.createMarkerNode(nodeIdentifier, brainstem_coordinates, brainstemCoordinatesValues)
nodeIdentifier += 1

return annotationGroups
Expand Down Expand Up @@ -701,7 +659,7 @@ def createCranialNerveEmergentMarkers(region, mesh, coordinatesName):
# return element xi
# use findMeshLocation to find the elementxi in an arbitrary mesh of given number of elements.

if coordinatesName == "bodyCoordinates":
if coordinatesName == "brainstem coordinates":
# brainstem_coordinates: the left-side nerves
nerveDict = {'OCULOMOTOR_left': [-0.13912257342955267, -0.5345161733750351, -0.7374762051676923],
'TROCHLEAR_left': [-0.13148279719950992, 0.4218745504359067, -0.7375838988856348],
Expand Down
39 changes: 36 additions & 3 deletions src/scaffoldmaker/meshtypes/meshtype_3d_heartatria1.py
Original file line number Diff line number Diff line change
Expand Up @@ -12,7 +12,7 @@
findOrCreateFieldNodeGroup, findOrCreateFieldStoredMeshLocation, findOrCreateFieldStoredString
from opencmiss.utils.zinc.finiteelement import getMaximumElementIdentifier, getMaximumNodeIdentifier
from opencmiss.zinc.element import Element, Elementbasis
from opencmiss.zinc.field import Field
from opencmiss.zinc.field import Field, FieldGroup
from opencmiss.zinc.node import Node
from scaffoldmaker.annotation.annotationgroup import AnnotationGroup, findOrCreateAnnotationGroupForTerm, getAnnotationGroupForTerm
from scaffoldmaker.annotation.heart_terms import get_heart_term
Expand Down Expand Up @@ -837,10 +837,11 @@ def generateBaseMesh(cls, region, options):
svcGroup = AnnotationGroup(region, get_heart_term("superior vena cava"))
laeVenousMidpointGroup = AnnotationGroup(region, get_heart_term("left atrium epicardium venous midpoint"))
raeVenousMidpointGroup = AnnotationGroup(region, get_heart_term("right atrium epicardium venous midpoint"))
annotationGroups += [ laeVenousMidpointGroup, ivcGroup, ivcInletGroup, raeVenousMidpointGroup, svcGroup, svcInletGroup ]
# av boundary nodes are put in left and right fibrous ring groups only so they can be found by heart1
# av boundary nodes are put in left and right fibrous ring groups so they can be found by heart1
lFibrousRingGroup = AnnotationGroup(region, get_heart_term("left fibrous ring"))
rFibrousRingGroup = AnnotationGroup(region, get_heart_term("right fibrous ring"))
annotationGroups += [laeVenousMidpointGroup, ivcGroup, ivcInletGroup, raeVenousMidpointGroup,
svcGroup, svcInletGroup, lFibrousRingGroup, rFibrousRingGroup]

# annotation fiducial points
markerGroup = findOrCreateFieldGroup(fm, "marker")
Expand Down Expand Up @@ -3312,6 +3313,38 @@ def defineFaceAnnotations(cls, region, options, annotationGroups):
raaEpiGroup = findOrCreateAnnotationGroupForTerm(annotationGroups, region, get_heart_term("epicardium of right auricle"))
raaEpiGroup.getMeshGroup(mesh2d).addElementsConditional(is_raa_epi)

lFibrousRingGroup = getAnnotationGroupForTerm(annotationGroups, get_heart_term("left fibrous ring"))
rFibrousRingGroup = getAnnotationGroupForTerm(annotationGroups, get_heart_term("right fibrous ring"))
if (lFibrousRingGroup.getDimension() == 0) or (lFibrousRingGroup.getDimension() == 0):
# not already added by full heart scaffold
nodes = fm.findNodesetByFieldDomainType(Field.DOMAIN_TYPE_NODES)
lFibrousRingNodeGroup = lFibrousRingGroup.getNodesetGroup(nodes)
rFibrousRingNodeGroup = rFibrousRingGroup.getNodesetGroup(nodes)
# make temp group containing all elements, faces etc. if have any nodes in fibrous ring groups
tmpGroup = fm.createFieldGroup()
tmpGroup.setSubelementHandlingMode(FieldGroup.SUBELEMENT_HANDLING_MODE_FULL)
mesh3d = fm.findMeshByDimension(3)
tmp3dMeshGroup = tmpGroup.createFieldElementGroup(mesh3d).getMeshGroup()
coordinates = fm.findFieldByName("coordinates").castFiniteElement()
elementIter = mesh3d.createElementiterator()
element = elementIter.next()
while element.isValid():
eft = element.getElementfieldtemplate(coordinates, -1)
if eft.isValid():
for n in range(eft.getNumberOfLocalNodes()):
node = element.getNode(eft, n + 1)
if lFibrousRingNodeGroup.containsNode(node) or rFibrousRingNodeGroup.containsNode(node):
tmp3dMeshGroup.addElement(element)
break
element = elementIter.next()
tmp2dElementGroup = tmpGroup.getFieldElementGroup(mesh2d)
if tmp2dElementGroup.isValid():
is_exterior_face_xi2_0 = fm.createFieldAnd(is_exterior, fm.createFieldIsOnFace(Element.FACE_TYPE_XI2_0))
is_fibrous_ring = fm.createFieldAnd(tmp2dElementGroup, is_exterior_face_xi2_0)
is_left_fibrous_ring = fm.createFieldAnd(is_la, is_fibrous_ring)
lFibrousRingGroup.getMeshGroup(mesh2d).addElementsConditional(is_left_fibrous_ring)
is_right_fibrous_ring = fm.createFieldAnd(is_ra, is_fibrous_ring)
rFibrousRingGroup.getMeshGroup(mesh2d).addElementsConditional(is_right_fibrous_ring)

def getLeftAtriumPulmonaryVeinOstiaElementsCounts(elementsCountAroundLeftAtriumFreeWall, elementsCountOverAtria, commonLeftRightPvOstium):
'''
Expand Down
Loading

0 comments on commit 5e3b453

Please sign in to comment.