Stomach

src/scaffoldmaker/annotation/stomach_terms.py

@@ -0,0 +1,30 @@
+"""
+Common resource for stomach annotation terms.
+"""
+
+# convention: preferred name, preferred id, followed by any other ids and alternative names
+stomach_terms = [
+    ( "body of stomach", "UBERON:0001161", " FMA:14560", "ILX:0724929"),
+    ( "cardia of stomach", "UBERON:0001162", " FMA:14561", "ILX:0729096"),
+    ( "duodenum", "UBERON:0002114", " FMA:7206", "ILX:0726125"),
+    ( "esophagus", "UBERON:0001043", "FMA: 7131", "ILX:0735017"),
+    ( "esophagogastric junction", "UBERON:0007650", "FMA: 9434", "ILX:0733910"),
+    ( "forestomach-glandular stomach junction", "UBERON:0012270", "ILX:0729974"),
+    ( "forestomach-glandular stomach junction on inner wall", None),
+    ( "forestomach-glandular stomach junction on outer wall", None),
+    ( "fundus of stomach", "UBERON:0001160", " FMA:14559", "ILX:0724443"),
+    ( "pyloric antrum", "UBERON:0001165", " FMA:14579", "ILX:0728672"),
+    ( "pylorus", "UBERON:0001166", " FMA:14581", "ILX:0734150"),
+    ( "stomach", "UBERON:0000945", "FMA:7148", "ILX:0736697")
+    ]
+
+def get_stomach_term(name : str):
+    """
+    Find term by matching name to any identifier held for a term.
+    Raise exception if name not found.
+    :return ( preferred name, preferred id )
+    """
+    for term in stomach_terms:
+        if name in term:
+            return ( term[0], term[1] )
+    raise NameError("Stomach annotation term '" + name + "' not found.")

src/scaffoldmaker/scaffolds.py

@@ -38,6 +38,7 @@
 from scaffoldmaker.meshtypes.meshtype_3d_sphereshell1 import MeshType_3d_sphereshell1
 from scaffoldmaker.meshtypes.meshtype_3d_sphereshellseptum1 import MeshType_3d_sphereshellseptum1
 from scaffoldmaker.meshtypes.meshtype_3d_stellate1 import MeshType_3d_stellate1
+from scaffoldmaker.meshtypes.meshtype_3d_stomach1 import MeshType_3d_stomach1
 from scaffoldmaker.meshtypes.meshtype_3d_stomachhuman1 import MeshType_3d_stomachhuman1
 from scaffoldmaker.meshtypes.meshtype_3d_tube1 import MeshType_3d_tube1
 from scaffoldmaker.meshtypes.meshtype_3d_tubeseptum1 import MeshType_3d_tubeseptum1
@@ -84,6 +85,7 @@ def __init__(self):
             MeshType_3d_sphereshell1,
             MeshType_3d_sphereshellseptum1,
             MeshType_3d_stellate1,
+            MeshType_3d_stomach1,
             MeshType_3d_stomachhuman1,
             MeshType_3d_tube1,
             MeshType_3d_tubeseptum1,

src/scaffoldmaker/utils/eftfactory_bicubichermitelinear.py

@@ -249,6 +249,70 @@ def createEftWedgeXi1Zero(self):
         assert eft.validate(), 'eftfactory_tricubichermite.createEftWedgeXi1Zero:  Failed to validate eft'
         return eft
 
+    def createEftWedgeCollapseXi1Quadrant(self, collapseNodes):
+        '''
+        Create a bicubic hermite linear element field for a wedge element, where xi1 collapsed on xi3 = 0 or xi3 = 1.
+        :return: Element field template
+        '''
+        eft = self.createEftBasic()
+        setEftScaleFactorIds(eft, [1], [])
+
+        valid = True
+        if collapseNodes in [[1, 3], [2, 4]]:
+            nodes = [1, 2, 3, 4]
+            if collapseNodes == [1, 3]:
+                remapEftNodeValueLabel(eft, nodes, Node.VALUE_LABEL_D_DS1, [])
+            elif collapseNodes == [2, 4]:
+                remapEftNodeValueLabel(eft, nodes, Node.VALUE_LABEL_D_DS1, [])
+            else:
+                valid = False
+            ln_map = [1, 1, 2, 2, 3, 4, 5, 6]
+        elif collapseNodes in [[5, 7], [6, 8]]:
+            nodes = [5, 6, 7, 8]
+            if collapseNodes == [5, 7]:
+                remapEftNodeValueLabel(eft, nodes, Node.VALUE_LABEL_D_DS1, [])
+            elif collapseNodes == [6, 8]:
+                remapEftNodeValueLabel(eft, nodes, Node.VALUE_LABEL_D_DS1, [])
+            else:
+                valid = False
+            ln_map = [1, 2, 3, 4, 5, 5, 6, 6]
+        elif collapseNodes in [[1, 5], [2, 6]]:
+            nodes = [1, 2, 5, 6]
+            # remap parameters on xi2 = 0 before collapsing nodes
+            if collapseNodes == [1, 5]:
+                remapEftNodeValueLabel(eft, collapseNodes, Node.VALUE_LABEL_D_DS2, [(Node.VALUE_LABEL_D_DS1, [1])])
+                remapEftNodeValueLabel(eft, nodes, Node.VALUE_LABEL_D_DS1, [])
+            elif collapseNodes == [2, 6]:
+                remapEftNodeValueLabel(eft, nodes, Node.VALUE_LABEL_D_DS1, [])
+                remapEftNodeValueLabel(eft, collapseNodes, Node.VALUE_LABEL_D_DS2, [(Node.VALUE_LABEL_D_DS1, [])])
+            else:
+                valid = False
+            ln_map = [1, 1, 2, 3, 4, 4, 5, 6]
+        elif collapseNodes in [[3, 7], [4, 8]]:
+            nodes = [3, 4, 7, 8]
+            # remap parameters on xi2 = 1 before collapsing nodes
+            if collapseNodes == [3, 7]:
+                remapEftNodeValueLabel(eft, nodes, Node.VALUE_LABEL_D_DS1, [])
+                remapEftNodeValueLabel(eft, collapseNodes, Node.VALUE_LABEL_D_DS2, [(Node.VALUE_LABEL_D_DS1, [])])
+            elif collapseNodes == [4, 8]:
+                remapEftNodeValueLabel(eft, collapseNodes, Node.VALUE_LABEL_D_DS2, [(Node.VALUE_LABEL_D_DS1, [1])])
+                remapEftNodeValueLabel(eft, nodes, Node.VALUE_LABEL_D_DS1, [])
+            else:
+                valid = False
+            ln_map = [1, 2, 3, 3, 4, 5, 6, 6]
+        else:
+            valid = False
+
+        if not valid:
+            assert False, "createEftWedgeCollapseXi1Quadrant.  Not implemented for collapse nodes " + str(collapseNodes)
+
+        # zero cross derivative parameters
+        remapEftNodeValueLabel(eft, nodes, Node.VALUE_LABEL_D2_DS1DS2, [])
+
+        remapEftLocalNodes(eft, 6, ln_map)
+        assert eft.validate(), 'eftfactory_bicubichermitelinear.createEftWedgeCollapseXi1Quadrant:  Failed to validate eft'
+        return eft
+
     def createEftWedgeXi1ZeroOpenTube(self):
         '''
         Create a basic bicubic hermite linear element template for elements
@@ -431,3 +495,30 @@ def createEftPyramidBottomSimple(self, nodeScaleFactorOffset0, nodeScaleFactorOf
 
         assert eft.validate(), 'eftfactory_bicubichermitelinear.createEftPyramidBottomSimple:  Failed to validate eft'
         return eft
+
+    def createEftWedgeCollapseXi2(self, collapseNodes):
+        '''
+        Create a bicubic hermite linear element field for a wedge element, where xi2 collapsed on xi1 = 1.
+        :return: Element field template
+        '''
+        eft = self.createEftBasic()
+
+        if collapseNodes in [[4, 8]]:
+            setEftScaleFactorIds(eft, [1], [])
+            nodes = [2, 4, 6, 8]
+            remapEftNodeValueLabel(eft, [4, 8], Node.VALUE_LABEL_D_DS1, [(Node.VALUE_LABEL_D_DS2, [1])])
+            remapEftNodeValueLabel(eft, nodes, Node.VALUE_LABEL_D_DS2, [])
+            ln_map = [1, 2, 3, 2, 4, 5, 6, 5]
+
+        elif collapseNodes in [[3, 7]]:
+            nodes = [1, 3, 5, 7]
+            remapEftNodeValueLabel(eft, nodes, Node.VALUE_LABEL_D_DS2, [])
+            remapEftNodeValueLabel(eft, [3, 7], Node.VALUE_LABEL_D_DS1, [(Node.VALUE_LABEL_D_DS2, [])])
+            ln_map = [1, 2, 1, 3, 4, 5, 4, 6]
+
+        remapEftNodeValueLabel(eft, nodes, Node.VALUE_LABEL_D2_DS1DS2, [])
+
+        remapEftLocalNodes(eft, 6, ln_map)
+        assert eft.validate(), 'eftfactory_bicubichermitelinear.createEftWedgeCollapseXi2:  Failed to validate eft'
+
+        return eft

src/scaffoldmaker/utils/eftfactory_tricubichermite.py

@@ -1663,3 +1663,31 @@ def replaceTwoElementWithInlet6(self, origElement1, origElement2, startElementId
         self._mesh.destroyElement(origElement1)
         self._mesh.destroyElement(origElement2)
         fm.endChange()
+
+    def createEftWedgeCollapseXi2(self, collapseNodes):
+        '''
+        Create a tricubic hermite element field for a wedge element, where xi2 collapsed on xi1 = 1.
+        :return: Element field template
+        '''
+        eft = self.createEftBasic()
+
+        if collapseNodes in [[4, 8]]:
+            setEftScaleFactorIds(eft, [1], [])
+            nodes = [2, 4, 6, 8]
+            remapEftNodeValueLabel(eft, [4, 8], Node.VALUE_LABEL_D_DS1, [(Node.VALUE_LABEL_D_DS2, [1])])
+            remapEftNodeValueLabel(eft, nodes, Node.VALUE_LABEL_D_DS2, [])
+            ln_map = [1, 2, 3, 2, 4, 5, 6, 5]
+
+        elif collapseNodes in [[3, 7]]:
+            nodes = [1, 3, 5, 7]
+            remapEftNodeValueLabel(eft, nodes, Node.VALUE_LABEL_D_DS2, [])
+            remapEftNodeValueLabel(eft, [3, 7], Node.VALUE_LABEL_D_DS1, [(Node.VALUE_LABEL_D_DS2, [])])
+            ln_map = [1, 2, 1, 3, 4, 5, 4, 6]
+
+        remapEftNodeValueLabel(eft, nodes, Node.VALUE_LABEL_D2_DS1DS2, [])
+        remapEftNodeValueLabel(eft, nodes, Node.VALUE_LABEL_D2_DS2DS3, [])
+        remapEftNodeValueLabel(eft, nodes, Node.VALUE_LABEL_D3_DS1DS2DS3, [])
+
+        remapEftLocalNodes(eft, 6, ln_map)
+        assert eft.validate(), 'eftfactory_tricubichermite.createEftWedgeCollapseXi2:  Failed to validate eft'
+        return eft

tests/test_colon.py

@@ -112,14 +112,14 @@ def test_colon1(self):
         coordinates = fieldmodule.findFieldByName("coordinates").castFiniteElement()
         self.assertTrue(coordinates.isValid())
         minimums, maximums = evaluateFieldNodesetRange(coordinates, nodes)
-        assertAlmostEqualList(self, minimums, [ 108.03959866945482, -36.876103983560014, -25.93217903595996 ], 1.0E-6)
-        assertAlmostEqualList(self, maximums, [ 185.43446761757468, 48.07433517752282, 34.995316052158934 ], 1.0E-6)
+        assertAlmostEqualList(self, minimums, [ 108.02506479907721, -36.876103983560014, -25.89741158484918 ], 1.0E-6)
+        assertAlmostEqualList(self, maximums, [ 185.46457506076914, 48.1011574894518, 34.995316052158934 ], 1.0E-6)
 
         flatCoordinates = fieldmodule.findFieldByName("flat coordinates").castFiniteElement()
         self.assertTrue(flatCoordinates.isValid())
         minimums, maximums = evaluateFieldNodesetRange(flatCoordinates, nodes)
         assertAlmostEqualList(self, minimums, [ 0.0, 0.0, 0.0 ], 1.0E-6)
-        assertAlmostEqualList(self, maximums, [ 186.72988844629867, 75.79769125771575, 3.2000000000000006 ], 1.0E-6)
+        assertAlmostEqualList(self, maximums, [ 186.72988844629867, 77.41781871321301, 3.2000000000000006 ], 1.0E-6)
 
         textureCoordinates = fieldmodule.findFieldByName("texture coordinates").castFiniteElement()
         minimums, maximums = evaluateFieldNodesetRange(textureCoordinates, nodes)
@@ -136,10 +136,10 @@ def test_colon1(self):
         fieldcache = fieldmodule.createFieldcache()
         result, surfaceArea = surfaceAreaField.evaluateReal(fieldcache, 1)
         self.assertEqual(result, RESULT_OK)
-        self.assertAlmostEqual(surfaceArea, 14623.340352853866, delta=1.0E-6)
+        self.assertAlmostEqual(surfaceArea, 14612.416789097502, delta=1.0E-6)
         result, volume = volumeField.evaluateReal(fieldcache, 1)
         self.assertEqual(result, RESULT_OK)
-        self.assertAlmostEqual(volume, 26869.74823158621, delta=1.0E-6)
+        self.assertAlmostEqual(volume, 26826.06954301569, delta=1.0E-6)
 
     def test_mousecolon1(self):
         """

tests/test_smallintestine.py

@@ -20,7 +20,7 @@ def test_smallintestine1(self):
         Test creation of small intestine scaffold.
         """
         parameterSetNames = MeshType_3d_smallintestine1.getParameterSetNames()
-        self.assertEqual(parameterSetNames, ["Default", "Mouse 1"])
+        self.assertEqual(parameterSetNames, ["Default", "Cattle 1", "Mouse 1"])
         centralPathDefaultScaffoldPackages = {
             'Test line': ScaffoldPackage(MeshType_1d_path1, {
                 'scaffoldSettings': {

tests/test_stomach.py

@@ -0,0 +1,84 @@
+import unittest
+from opencmiss.utils.zinc.finiteelement import evaluateFieldNodesetRange
+from opencmiss.utils.zinc.general import ChangeManager
+from opencmiss.zinc.context import Context
+from opencmiss.zinc.element import Element
+from opencmiss.zinc.field import Field
+from opencmiss.zinc.result import RESULT_OK
+from scaffoldmaker.meshtypes.meshtype_3d_stomach1 import MeshType_3d_stomach1
+from scaffoldmaker.utils.zinc_utils import createFaceMeshGroupExteriorOnFace
+from testutils import assertAlmostEqualList
+
+class StomachScaffoldTestCase(unittest.TestCase):
+
+    def test_stomach1(self):
+        """
+        Test creation of stomach scaffold.
+        """
+        parameterSetNames = MeshType_3d_stomach1.getParameterSetNames()
+        self.assertEqual(parameterSetNames, [ "Default", "Human 1", "Rat 1" ])
+        options = MeshType_3d_stomach1.getDefaultOptions("Rat 1")
+        self.assertEqual(17, len(options))
+        self.assertEqual(12, options.get("Number of elements around esophagus"))
+        self.assertEqual(14, options.get("Number of elements around duodenum"))
+        self.assertEqual(2, options.get("Number of elements between annulus and duodenum"))
+        self.assertEqual(1, options.get("Number of elements through wall"))
+        self.assertEqual(1, options.get("Number of radial elements in annulus"))
+        self.assertEqual(0.5, options.get("Wall thickness"))
+        self.assertEqual(True, options.get("Limiting ridge"))
+        ostiumOptions = options['Gastro-esophagal junction']
+        ostiumSettings = ostiumOptions.getScaffoldSettings()
+        self.assertEqual(1, ostiumSettings.get("Number of vessels"))
+        self.assertEqual(8, ostiumSettings.get("Number of elements around ostium"))
+        self.assertEqual(1, ostiumSettings.get("Number of elements through wall"))
+        self.assertEqual(4.0, ostiumSettings.get("Ostium diameter"))
+        self.assertEqual(3.5, ostiumSettings.get("Ostium length"))
+        self.assertEqual(0.5, ostiumSettings.get("Ostium wall thickness"))
+        self.assertEqual(1.25, ostiumSettings.get("Vessel inner diameter"))
+        self.assertEqual(0.5, ostiumSettings.get("Vessel wall thickness"))
+        self.assertEqual(0.0, ostiumSettings.get("Vessel angle 1 degrees"))
+        self.assertEqual(0.55, options.get("Gastro-esophagal junction position along factor"))
+        self.assertEqual(0.2, options.get("Annulus derivative factor"))
+
+        context = Context("Test")
+        region = context.getDefaultRegion()
+        self.assertTrue(region.isValid())
+        annotationGroups = MeshType_3d_stomach1.generateBaseMesh(region, options)
+        self.assertEqual(8, len(annotationGroups))
+
+        fieldmodule = region.getFieldmodule()
+        self.assertEqual(RESULT_OK, fieldmodule.defineAllFaces())
+        mesh3d = fieldmodule.findMeshByDimension(3)
+        self.assertEqual(158, mesh3d.getSize())
+        mesh2d = fieldmodule.findMeshByDimension(2)
+        self.assertEqual(643, mesh2d.getSize())
+        mesh1d = fieldmodule.findMeshByDimension(1)
+        self.assertEqual(823, mesh1d.getSize())
+        nodes = fieldmodule.findNodesetByFieldDomainType(Field.DOMAIN_TYPE_NODES)
+        self.assertEqual(341, nodes.getSize())
+        datapoints = fieldmodule.findNodesetByFieldDomainType(Field.DOMAIN_TYPE_DATAPOINTS)
+        self.assertEqual(0, datapoints.getSize())
+
+        coordinates = fieldmodule.findFieldByName("coordinates").castFiniteElement()
+        self.assertTrue(coordinates.isValid())
+        minimums, maximums = evaluateFieldNodesetRange(coordinates, nodes)
+        assertAlmostEqualList(self, minimums, [-17.977095466754566, -15.437578891036395, -8.706694306455596], 1.0E-6)
+        assertAlmostEqualList(self, maximums, [17.943222678234513, 15.191836205539767, 8.725549026319936], 1.0E-6)
+
+        with ChangeManager(fieldmodule):
+            one = fieldmodule.createFieldConstant(1.0)
+            faceMeshGroup = createFaceMeshGroupExteriorOnFace(fieldmodule, Element.FACE_TYPE_XI3_1)
+            surfaceAreaField = fieldmodule.createFieldMeshIntegral(one, coordinates, faceMeshGroup)
+            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, 2436.4955183926895, delta=1.0E-6)
+        result, volume = volumeField.evaluateReal(fieldcache, 1)
+        self.assertEqual(result, RESULT_OK)
+        self.assertAlmostEqual(volume, 1168.9418891341106, delta=1.0E-6)
+
+if __name__ == "__main__":
+    unittest.main()