diff --git a/src/scaffoldmaker/annotation/uterus_terms.py b/src/scaffoldmaker/annotation/uterus_terms.py index 16a7dc48..0060690e 100644 --- a/src/scaffoldmaker/annotation/uterus_terms.py +++ b/src/scaffoldmaker/annotation/uterus_terms.py @@ -9,6 +9,7 @@ ("dorsal cervix junction with vagina", "None"), ("dorsal top left horn", "None"), ("dorsal top right horn", "None"), + ("dorsal uterus", "None"), ("external cervical os", "UBERON:0013760", "FMA:76836", "ILX:0736534"), ("fundus of uterus", "None"), ("internal cervical os", "UBERON:0013759", "FMA:17747", "ILX:0729495"), @@ -18,6 +19,7 @@ ("left transverse cervical ligament", "None"), ("left uterine horn", "UBERON:0009020"), ("left uterine tube", "UBERON:0001303", "FMA:18484", "ILX:0734218"), + ("left uterus", "None"), ("lumen of body of uterus", "None"), ("lumen of fallopian tube", "None"), ("lumen of left horn", "None"), @@ -33,6 +35,7 @@ ("right transverse cervical ligament", "None"), ("right uterine horn", "UBERON:0009022"), ("right uterine tube", "UBERON:0001302", "FMA:18483", "ILX:0724908"), + ("right uterus", "None"), ("serosa of body of uterus", "None"), ("serosa of left uterine tube", "None"), ("serosa of left horn", "None"), @@ -41,7 +44,7 @@ ("serosa of uterine cervix", "None"), ("serosa of uterus", "UBERON:0001297"), ("serosa of vagina", "None"), - ("uterine cervix", "UBERON:0000002","FMA:17740", "ILX:0724162"), + ("uterine cervix", "UBERON:0000002", "FMA:17740", "ILX:0724162"), ("uterine horn", "UBERON:000224"), ("uterine lumen", "UBERON:0013769"), ("uterine wall", "UBERON:0000459", "FMA:17560", "ILX:0735839"), @@ -50,15 +53,17 @@ ("vagina orifice", "UBERON:0012317", "FMA:19984", "ILX:0729556"), ("ventral cervix junction with vagina", "None"), ("ventral top left horn", "None"), - ("ventral top right horn", "None")] + ("ventral top right horn", "None"), + ("ventral uterus", "None")] + def get_uterus_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 ) + :return: ( preferred name, preferred id ) """ for term in uterus_terms: if name in term: - return (term[0], term[1]) - raise NameError("Uterus annotation term '" + name + "' not found.") \ No newline at end of file + return term[0], term[1] + raise NameError("Uterus annotation term '" + name + "' not found.") diff --git a/src/scaffoldmaker/meshtypes/meshtype_2d_tubenetwork1.py b/src/scaffoldmaker/meshtypes/meshtype_2d_tubenetwork1.py index ca175c2a..48b89cc7 100644 --- a/src/scaffoldmaker/meshtypes/meshtype_2d_tubenetwork1.py +++ b/src/scaffoldmaker/meshtypes/meshtype_2d_tubenetwork1.py @@ -34,8 +34,8 @@ def getDefaultOptions(cls, parameterSetName="Default"): options["Target element density along longest segment"] = 12.0 return options - @staticmethod - def getOrderedOptionNames(): + @classmethod + def getOrderedOptionNames(cls): return [ "Network layout", "Number of elements around", @@ -119,11 +119,11 @@ def generateBaseMesh(cls, region, options): tubeNetworkMeshBuilder = TubeNetworkMeshBuilder( networkMesh, targetElementDensityAlongLongestSegment=options["Target element density along longest segment"], - defaultElementsCountAround=options["Number of elements around"], - elementsCountThroughWall=1, layoutAnnotationGroups=networkLayout.getAnnotationGroups(), annotationElementsCountsAlong=options["Annotation numbers of elements along"], - annotationElementsCountsAround=options["Annotation numbers of elements around"]) + defaultElementsCountAround=options["Number of elements around"], + annotationElementsCountsAround=options["Annotation numbers of elements around"], + elementsCountThroughShell=1) tubeNetworkMeshBuilder.build() generateData = TubeNetworkMeshGenerateData( region, 2, diff --git a/src/scaffoldmaker/meshtypes/meshtype_3d_tubenetwork1.py b/src/scaffoldmaker/meshtypes/meshtype_3d_tubenetwork1.py index 0aede4e2..c69e4184 100644 --- a/src/scaffoldmaker/meshtypes/meshtype_3d_tubenetwork1.py +++ b/src/scaffoldmaker/meshtypes/meshtype_3d_tubenetwork1.py @@ -53,8 +53,8 @@ def getDefaultOptions(cls, parameterSetName="Default"): options["Target element density along longest segment"] = 12.0 return options - @staticmethod - def getOrderedOptionNames(): + @classmethod + def getOrderedOptionNames(cls): return [ "Network layout", "Number of elements around", @@ -212,31 +212,21 @@ def generateBaseMesh(cls, region, options): layoutRegion = region.createRegion() networkLayout = options["Network layout"] networkLayout.generate(layoutRegion) # ask scaffold to generate to get user-edited parameters - layoutAnnotationGroups = networkLayout.getAnnotationGroups() networkMesh = networkLayout.getConstructionObject() - defaultAroundCount = options["Number of elements around"] - coreTransitionCount = options["Number of elements across core transition"] - defaultCoreBoxMinorCount = options["Number of elements across core box minor"] - # implementation currently uses major count including transition - defaultCoreMajorCount = defaultAroundCount // 2 - defaultCoreBoxMinorCount + 2 * coreTransitionCount - annotationAroundCounts = options["Annotation numbers of elements around"] - annotationCoreBoxMinorCounts = options["Annotation numbers of elements across core box minor"] - tubeNetworkMeshBuilder = TubeNetworkMeshBuilder( networkMesh, targetElementDensityAlongLongestSegment=options["Target element density along longest segment"], - defaultElementsCountAround=defaultAroundCount, - elementsCountThroughWall=options["Number of elements through shell"], - layoutAnnotationGroups=layoutAnnotationGroups, + layoutAnnotationGroups=networkLayout.getAnnotationGroups(), annotationElementsCountsAlong=options["Annotation numbers of elements along"], - annotationElementsCountsAround=annotationAroundCounts, - defaultElementsCountCoreBoxMinor=defaultCoreBoxMinorCount, - elementsCountTransition=coreTransitionCount, - annotationElementsCountsCoreBoxMinor=annotationCoreBoxMinorCounts, + defaultElementsCountAround=options["Number of elements around"], + annotationElementsCountsAround=options["Annotation numbers of elements around"], + elementsCountThroughShell=options["Number of elements through shell"], isCore=options["Core"], + elementsCountTransition=options["Number of elements across core transition"], + defaultElementsCountCoreBoxMinor=options["Number of elements across core box minor"], + annotationElementsCountsCoreBoxMinor=options["Annotation numbers of elements across core box minor"], useOuterTrimSurfaces=options["Use outer trim surfaces"]) - tubeNetworkMeshBuilder.build() generateData = TubeNetworkMeshGenerateData( region, 3, diff --git a/src/scaffoldmaker/meshtypes/meshtype_3d_uterus2.py b/src/scaffoldmaker/meshtypes/meshtype_3d_uterus2.py index b1c18601..13f712e9 100644 --- a/src/scaffoldmaker/meshtypes/meshtype_3d_uterus2.py +++ b/src/scaffoldmaker/meshtypes/meshtype_3d_uterus2.py @@ -18,46 +18,30 @@ class UterusTubeNetworkMeshGenerateData(TubeNetworkMeshGenerateData): - def __init__(self, region, meshDimension, isLinearThroughWall, isShowTrimSurfaces, - coordinateFieldName="coordinates", startNodeIdentifier=1, startElementIdentifier=1): + def __init__(self, region, meshDimension, coordinateFieldName="coordinates", + startNodeIdentifier=1, startElementIdentifier=1, isLinearThroughWall=False, isShowTrimSurfaces=False): """ :param isLinearThroughWall: Callers should only set if 3-D with no core. :param isShowTrimSurfaces: Tells junction generateMesh to make 2-D trim surfaces. """ super(UterusTubeNetworkMeshGenerateData, self).__init__( - region, meshDimension, isLinearThroughWall, isShowTrimSurfaces, - coordinateFieldName, startNodeIdentifier, startElementIdentifier) + region, meshDimension, coordinateFieldName, startNodeIdentifier, startElementIdentifier, + isLinearThroughWall, isShowTrimSurfaces) self._fundusGroup = self.getOrCreateAnnotationGroup(get_uterus_term("fundus of uterus")) - self._leftGroup = self.getOrCreateAnnotationGroup(("left uterus", "None")) - self._rightGroup = self.getOrCreateAnnotationGroup(("right uterus", "None")) - self._dorsalGroup = self.getOrCreateAnnotationGroup(("dorsal uterus", "None")) - self._ventralGroup = self.getOrCreateAnnotationGroup(("ventral uterus", "None")) + # force these annotation group names in base class + self._leftGroup = self.getOrCreateAnnotationGroup(get_uterus_term("left uterus")) + self._rightGroup = self.getOrCreateAnnotationGroup(get_uterus_term("right uterus")) + self._dorsalGroup = self.getOrCreateAnnotationGroup(get_uterus_term("dorsal uterus")) + self._ventralGroup = self.getOrCreateAnnotationGroup(get_uterus_term("ventral uterus")) def getFundusMeshGroup(self): return self._fundusGroup.getMeshGroup(self._mesh) - def getLeftMeshGroup(self): - return self._leftGroup.getMeshGroup(self._mesh) - - def getRightMeshGroup(self): - return self._rightGroup.getMeshGroup(self._mesh) - - def getDorsalMeshGroup(self): - return self._dorsalGroup.getMeshGroup(self._mesh) - - def getVentralMeshGroup(self): - return self._ventralGroup.getMeshGroup(self._mesh) - class UterusTubeNetworkMeshBuilder(TubeNetworkMeshBuilder): - - def __init__(self, networkMesh: NetworkMesh, targetElementDensityAlongLongestSegment: float, - defaultElementsCountAround: int, elementsCountThroughWall: int, - layoutAnnotationGroups: list = [], annotationElementsCountsAlong: list = [], - annotationElementsCountsAround: list = [], useOuterTrimSurfaces=True): - super(UterusTubeNetworkMeshBuilder, self).__init__( - networkMesh, targetElementDensityAlongLongestSegment, defaultElementsCountAround, - elementsCountThroughWall, layoutAnnotationGroups, - annotationElementsCountsAlong, annotationElementsCountsAround, useOuterTrimSurfaces=useOuterTrimSurfaces) + """ + Adds left, right, dorsal, ventral, fundus annotations. + Future: derive from BodyTubeNetworkMeshBuilder to get left/right/dorsal/ventral. + """ def generateMesh(self, generateData): super(UterusTubeNetworkMeshBuilder, self).generateMesh(generateData) @@ -490,11 +474,12 @@ def generateBaseMesh(cls, region, options): uterusTubeNetworkMeshBuilder = UterusTubeNetworkMeshBuilder( networkMesh, targetElementDensityAlongLongestSegment=options["Target element density along longest segment"], - defaultElementsCountAround=options['Number of elements around'], - elementsCountThroughWall=options["Number of elements through wall"], layoutAnnotationGroups=layoutAnnotationGroups, annotationElementsCountsAlong=annotationAlongCounts, - annotationElementsCountsAround=annotationElementsCountsAround) + defaultElementsCountAround=options['Number of elements around'], + annotationElementsCountsAround=annotationElementsCountsAround, + elementsCountThroughShell=options["Number of elements through wall"], + useOuterTrimSurfaces=True) uterusTubeNetworkMeshBuilder.build() generateData = UterusTubeNetworkMeshGenerateData( region, 3, @@ -633,10 +618,10 @@ def defineFaceAnnotations(cls, region, options, annotationGroups): get_uterus_term("lumen of vagina")) lumenOfVagina.getMeshGroup(mesh2d).addElementsConditional(is_vagina_inner) - leftGroup = getAnnotationGroupForTerm(annotationGroups, ("left uterus", "None")) - rightGroup = getAnnotationGroupForTerm(annotationGroups, ("right uterus", "None")) - dorsalGroup = getAnnotationGroupForTerm(annotationGroups, ("dorsal uterus", "None")) - ventralGroup = getAnnotationGroupForTerm(annotationGroups, ("ventral uterus", "None")) + leftGroup = getAnnotationGroupForTerm(annotationGroups, get_uterus_term("left uterus")) + rightGroup = getAnnotationGroupForTerm(annotationGroups, get_uterus_term("right uterus")) + dorsalGroup = getAnnotationGroupForTerm(annotationGroups, get_uterus_term("dorsal uterus")) + ventralGroup = getAnnotationGroupForTerm(annotationGroups, get_uterus_term("ventral uterus")) if isHuman: leftUterineTubeGroup = getAnnotationGroupForTerm(annotationGroups, get_uterus_term("left uterine tube")) diff --git a/src/scaffoldmaker/meshtypes/meshtype_3d_wholebody2.py b/src/scaffoldmaker/meshtypes/meshtype_3d_wholebody2.py index d05e75ed..ae7f5f5f 100644 --- a/src/scaffoldmaker/meshtypes/meshtype_3d_wholebody2.py +++ b/src/scaffoldmaker/meshtypes/meshtype_3d_wholebody2.py @@ -15,9 +15,10 @@ computeCubicHermiteEndDerivative, getCubicHermiteArcLength, interpolateLagrangeHermiteDerivative, sampleCubicHermiteCurvesSmooth, smoothCubicHermiteDerivativesLine) from scaffoldmaker.utils.networkmesh import NetworkMesh -from scaffoldmaker.utils.tubenetworkmesh import TubeNetworkMeshBuilder, TubeNetworkMeshGenerateData +from scaffoldmaker.utils.tubenetworkmesh import BodyTubeNetworkMeshBuilder, TubeNetworkMeshGenerateData import math + class MeshType_1d_human_body_network_layout1(MeshType_1d_network_layout1): """ Defines body network layout. @@ -169,7 +170,6 @@ def checkOptions(cls, options): options[key] = 60.0 return dependentChanges - @classmethod def generateBaseMesh(cls, region, options): """ @@ -292,7 +292,7 @@ def generateBaseMesh(cls, region, options): footElementsCount = 2 legMeshGroup = legGroup.getMeshGroup(mesh) legToFootMeshGroup = legToFootGroup.getMeshGroup(mesh) - footMeshGroup = footGroup.getMeshGroup(mesh) + footMeshGroup = footGroup.getMeshGroup(mesh) for side in (left, right): sideLegGroup = leftLegGroup if (side == left) else rightLegGroup meshGroups = [bodyMeshGroup, legMeshGroup, legToFootMeshGroup, sideLegGroup.getMeshGroup(mesh)] @@ -379,7 +379,6 @@ def generateBaseMesh(cls, region, options): id2 = mult(d2, innerProportionDefault) id3 = mult(d3, innerProportionDefault) abdomenStartX = thoraxStartX + thoraxLength - px = None for i in range(abdomenElementsCount + 1): node = nodes.findNodeByIdentifier(nodeIdentifier) fieldcache.setNode(node) @@ -487,7 +486,6 @@ def generateBaseMesh(cls, region, options): nodeIdentifier += 1 # legs - cos45 = math.cos(0.25 * math.pi) legStartX = abdomenStartX + abdomenLength + pelvisDrop nonFootLegLength = legLength - footHeight legScale = nonFootLegLength / (legToFootElementsCount - 1) @@ -578,83 +576,6 @@ def getInteractiveFunctions(cls): return interactiveFunctions -class WholeBodyTubeNetworkMeshGenerateData(TubeNetworkMeshGenerateData): - - def __init__(self, region, meshDimension, isLinearThroughWall, isShowTrimSurfaces, - coordinateFieldName="coordinates", startNodeIdentifier=1, startElementIdentifier=1): - """ - :param isLinearThroughWall: Callers should only set if 3-D with no core. - :param isShowTrimSurfaces: Tells junction generateMesh to make 2-D trim surfaces. - """ - super(WholeBodyTubeNetworkMeshGenerateData, self).__init__( - region, meshDimension, isLinearThroughWall, isShowTrimSurfaces, - coordinateFieldName, startNodeIdentifier, startElementIdentifier) - # annotation groups are created on demand: - self._leftGroup = None - self._rightGroup = None - self._dorsalGroup = None - self._ventralGroup = None - - def getLeftMeshGroup(self): - if not self._leftGroup: - self._leftGroup = self.getOrCreateAnnotationGroup(("left", "")) - return self._leftGroup.getMeshGroup(self._mesh) - - def getRightMeshGroup(self): - if not self._rightGroup: - self._rightGroup = self.getOrCreateAnnotationGroup(("right", "")) - return self._rightGroup.getMeshGroup(self._mesh) - - def getDorsalMeshGroup(self): - if not self._dorsalGroup: - self._dorsalGroup = self.getOrCreateAnnotationGroup(("dorsal", "")) - return self._dorsalGroup.getMeshGroup(self._mesh) - - def getVentralMeshGroup(self): - if not self._ventralGroup: - self._ventralGroup = self.getOrCreateAnnotationGroup(("ventral", "")) - return self._ventralGroup.getMeshGroup(self._mesh) - - -class WholeBodyNetworkMeshBuilder(TubeNetworkMeshBuilder): - - def __init__(self, networkMesh: NetworkMesh, targetElementDensityAlongLongestSegment: float, - defaultElementsCountAround: int, elementsCountThroughWall: int, layoutAnnotationGroups: list = [], - annotationElementsCountsAlong: list = [], annotationElementsCountsAround: list = [], - defaultElementsCountCoreBoxMinor: int = 2, elementsCountTransition: int = 1, - annotationElementsCountsCoreBoxMinor: list = [], isCore=False, useOuterTrimSurfaces=True): - super(WholeBodyNetworkMeshBuilder, self).__init__( - networkMesh, targetElementDensityAlongLongestSegment, defaultElementsCountAround, - elementsCountThroughWall, layoutAnnotationGroups, - annotationElementsCountsAlong, annotationElementsCountsAround, defaultElementsCountCoreBoxMinor, - elementsCountTransition, annotationElementsCountsCoreBoxMinor, isCore, useOuterTrimSurfaces) - - def generateMesh(self, generateData): - super(WholeBodyNetworkMeshBuilder, self).generateMesh(generateData) - # build left, right, dorsal, ventral annotation groups - leftMeshGroup = generateData.getLeftMeshGroup() - rightMeshGroup = generateData.getRightMeshGroup() - dorsalMeshGroup = generateData.getDorsalMeshGroup() - ventralMeshGroup = generateData.getVentralMeshGroup() - for networkSegment in self._networkMesh.getNetworkSegments(): - # print("Segment", networkSegment.getNodeIdentifiers()) - segment = self._segments[networkSegment] - annotationTerms = segment.getAnnotationTerms() - for annotationTerm in annotationTerms: - if "left" in annotationTerm[0]: - segment.addAllElementsToMeshGroup(leftMeshGroup) - break - if "right" in annotationTerm[0]: - segment.addAllElementsToMeshGroup(rightMeshGroup) - break - else: - # segment on main axis - segment.addSideD2ElementsToMeshGroup(False, leftMeshGroup) - segment.addSideD2ElementsToMeshGroup(True, rightMeshGroup) - segment.addSideD3ElementsToMeshGroup(False, ventralMeshGroup) - segment.addSideD3ElementsToMeshGroup(True, dorsalMeshGroup) - - class MeshType_3d_wholebody2(Scaffold_base): """ Generates a 3-D hermite bifurcating tube network with core representing the human body. @@ -774,7 +695,8 @@ def getOptionScaffoldPackage(cls, optionName, scaffoldType, parameterSetName=Non @classmethod def checkOptions(cls, options): dependentChanges = False - if not options["Body network layout"].getScaffoldType() in cls.getOptionValidScaffoldTypes("Body network layout"): + if (options["Body network layout"].getScaffoldType() not in + cls.getOptionValidScaffoldTypes("Body network layout")): options["Body network layout"] = ScaffoldPackage(MeshType_1d_human_body_network_layout1) for key in [ "Number of elements along head", @@ -785,7 +707,7 @@ def checkOptions(cls, options): "Number of elements along hand", "Number of elements along leg to foot", "Number of elements along foot" - ]: + ]: if options[key] < 1: options[key] = 1 minElementsCountAround = None @@ -844,8 +766,7 @@ def generateBaseMesh(cls, region, options): elementsCountAroundTorso = options["Number of elements around torso"] elementsCountAroundArm = options["Number of elements around arm"] elementsCountAroundLeg = options["Number of elements around leg"] - coreBoxMinorCount = options["Number of elements across core box minor"] - coreTransitionCount = options['Number of elements across core transition'] + isCore = options["Use Core"] layoutRegion = region.createRegion() networkLayout.generate(layoutRegion) # ask scaffold to generate to get user-edited parameters @@ -854,7 +775,6 @@ def generateBaseMesh(cls, region, options): annotationAlongCounts = [] annotationAroundCounts = [] - # implementation currently uses major count including transition for layoutAnnotationGroup in layoutAnnotationGroups: alongCount = 0 aroundCount = 0 @@ -885,24 +805,24 @@ def generateBaseMesh(cls, region, options): aroundCount = elementsCountAroundLeg annotationAlongCounts.append(alongCount) annotationAroundCounts.append(aroundCount) - isCore = options["Use Core"] - tubeNetworkMeshBuilder = WholeBodyNetworkMeshBuilder( + tubeNetworkMeshBuilder = BodyTubeNetworkMeshBuilder( networkMesh, targetElementDensityAlongLongestSegment=2.0, # not used for body - defaultElementsCountAround=options["Number of elements around head"], - elementsCountThroughWall=options["Number of elements through shell"], layoutAnnotationGroups=layoutAnnotationGroups, annotationElementsCountsAlong=annotationAlongCounts, + defaultElementsCountAround=options["Number of elements around head"], annotationElementsCountsAround=annotationAroundCounts, - defaultElementsCountCoreBoxMinor=coreBoxMinorCount, - elementsCountTransition=coreTransitionCount, + elementsCountThroughShell=options["Number of elements through shell"], + isCore=isCore, + elementsCountTransition=options['Number of elements across core transition'], + defaultElementsCountCoreBoxMinor=options["Number of elements across core box minor"], annotationElementsCountsCoreBoxMinor=[], - isCore=isCore) + useOuterTrimSurfaces=True) meshDimension = 3 tubeNetworkMeshBuilder.build() - generateData = WholeBodyTubeNetworkMeshGenerateData( + generateData = TubeNetworkMeshGenerateData( region, meshDimension, isLinearThroughWall=False, isShowTrimSurfaces=options["Show trim surfaces"]) diff --git a/src/scaffoldmaker/utils/boxnetworkmesh.py b/src/scaffoldmaker/utils/boxnetworkmesh.py index 6c477ffd..10cc9cb1 100644 --- a/src/scaffoldmaker/utils/boxnetworkmesh.py +++ b/src/scaffoldmaker/utils/boxnetworkmesh.py @@ -186,7 +186,7 @@ def getNodeIdentifier(self): def setNodeIdentifier(self, nodeIdentifier): """ - Store the node identifier so it can be reused by adjacent segments. + Store the node identifier to reference from adjacent segments. :param nodeIdentifier: Identifier of generated node at junction. """ self._nodeIdentifier = nodeIdentifier @@ -222,7 +222,6 @@ def sample(self, targetElementLength): for segment, nodeIndex in segmentIndexes: segment.setSampledD1(nodeIndex, d1Mean) - def generateMesh(self, generateData: BoxNetworkMeshGenerateData): pass # nothing to do for box network diff --git a/src/scaffoldmaker/utils/networkmesh.py b/src/scaffoldmaker/utils/networkmesh.py index 58adc7a1..3d4f4f10 100644 --- a/src/scaffoldmaker/utils/networkmesh.py +++ b/src/scaffoldmaker/utils/networkmesh.py @@ -714,6 +714,17 @@ class NetworkMeshBuilder(ABC): def __init__(self, networkMesh: NetworkMesh, targetElementDensityAlongLongestSegment: float, layoutAnnotationGroups, annotationElementsCountsAlong=[]): + """ + Abstract base class for building meshes from a NetworkMesh network layout. + :param networkMesh: Description of the topology of the network layout. + :param targetElementDensityAlongLongestSegment: Real value which longest segment path in network is divided by + to get target element length, which is used to determine numbers of elements along except when set for a segment + through annotationElementsCountsAlong. + :param layoutAnnotationGroups: Annotation groups defined on the layout to mirror on the final mesh. + :param annotationElementsCountsAlong: List in same order as layoutAnnotationGroups, specifying fixed number of + elements along segment with any elements in the annotation group. Client must ensure exclusive map from + segments. Groups with zero value or past end of this list use the targetElementDensityAlongLongestSegment. + """ self._networkMesh = networkMesh self._targetElementDensityAlongLongestSegment = targetElementDensityAlongLongestSegment self._layoutAnnotationGroups = layoutAnnotationGroups diff --git a/src/scaffoldmaker/utils/tubenetworkmesh.py b/src/scaffoldmaker/utils/tubenetworkmesh.py index 89df3cb2..ac58b047 100644 --- a/src/scaffoldmaker/utils/tubenetworkmesh.py +++ b/src/scaffoldmaker/utils/tubenetworkmesh.py @@ -27,8 +27,8 @@ class TubeNetworkMeshGenerateData(NetworkMeshGenerateData): Data for passing to TubeNetworkMesh generateMesh functions. """ - def __init__(self, region, meshDimension, isLinearThroughWall, isShowTrimSurfaces, - coordinateFieldName="coordinates", startNodeIdentifier=1, startElementIdentifier=1): + def __init__(self, region, meshDimension, coordinateFieldName="coordinates", + startNodeIdentifier=1, startElementIdentifier=1, isLinearThroughWall=False, isShowTrimSurfaces=False): """ :param isLinearThroughWall: Callers should only set if 3-D with no core. :param isShowTrimSurfaces: Tells junction generateMesh to make 2-D trim surfaces. @@ -70,9 +70,14 @@ def __init__(self, region, meshDimension, isLinearThroughWall, isShowTrimSurface d3Defined, limitDirections=[None, [[0.0, 1.0, 0.0], [0.0, -1.0, 0.0]], None]) self._nodeLayoutTransitionTriplePoint = None - # annotation groups are created if core: + # annotation groups created if core: self._coreGroup = None self._shellGroup = None + # annotation groups created on demand: + self._leftGroup = None + self._rightGroup = None + self._dorsalGroup = None + self._ventralGroup = None def getStandardElementtemplate(self): return self._standardElementtemplate, self._standardEft @@ -160,6 +165,26 @@ def getShellMeshGroup(self): self._shellGroup = self.getOrCreateAnnotationGroup(("shell", "")) return self._shellGroup.getMeshGroup(self._mesh) + def getLeftMeshGroup(self): + if not self._leftGroup: + self._leftGroup = self.getOrCreateAnnotationGroup(("left", "")) + return self._leftGroup.getMeshGroup(self._mesh) + + def getRightMeshGroup(self): + if not self._rightGroup: + self._rightGroup = self.getOrCreateAnnotationGroup(("right", "")) + return self._rightGroup.getMeshGroup(self._mesh) + + def getDorsalMeshGroup(self): + if not self._dorsalGroup: + self._dorsalGroup = self.getOrCreateAnnotationGroup(("dorsal", "")) + return self._dorsalGroup.getMeshGroup(self._mesh) + + def getVentralMeshGroup(self): + if not self._ventralGroup: + self._ventralGroup = self.getOrCreateAnnotationGroup(("ventral", "")) + return self._ventralGroup.getMeshGroup(self._mesh) + def getNewTrimAnnotationGroup(self): self._trimAnnotationGroupCount += 1 return self.getOrCreateAnnotationGroup(("trim surface " + "{:03d}".format(self._trimAnnotationGroupCount), "")) @@ -167,13 +192,13 @@ def getNewTrimAnnotationGroup(self): class TubeNetworkMeshSegment(NetworkMeshSegment): - def __init__(self, networkSegment, pathParametersList, elementsCountAround, elementsCountThroughWall, + def __init__(self, networkSegment, pathParametersList, elementsCountAround, elementsCountThroughShell, isCore=False, elementsCountCoreBoxMinor: int = 2, elementsCountTransition: int = 1): """ :param networkSegment: NetworkSegment this is built from. :param pathParametersList: [pathParameters] if 2-D or [outerPathParameters, innerPathParameters] if 3-D :param elementsCountAround: Number of elements around this segment. - :param elementsCountThroughWall: Number of elements between inner and outer tube if 3-D, 1 if 2-D. + :param elementsCountThroughShell: Number of elements between inner and outer tube if 3-D, 1 if 2-D. :param isCore: True for generating a solid core inside the tube, False for regular tube network. :param elementsCountCoreBoxMinor: Number of elements across core box minor axis. :param elementsCountTransition: Number of elements across transition zone between core box elements and @@ -186,8 +211,8 @@ def __init__(self, networkSegment, pathParametersList, elementsCountAround, elem self._elementsCountCoreBoxMinor = elementsCountCoreBoxMinor self._elementsCountTransition = elementsCountTransition - assert elementsCountThroughWall > 0 - self._elementsCountThroughWall = elementsCountThroughWall + assert elementsCountThroughShell > 0 + self._elementsCountThroughShell = elementsCountThroughShell self._rawTubeCoordinatesList = [] self._rawTrackSurfaceList = [] for pathParameters in pathParametersList: @@ -403,7 +428,7 @@ def sample(self, fixedElementsCountAlong, targetElementLength): [[ring] for ring in self._sampledTubeCoordinates[0][2]], None) else: - wallFactor = 1.0 / self._elementsCountThroughWall + wallFactor = 1.0 / self._elementsCountThroughShell ox, od1, od2 = self._sampledTubeCoordinates[0][0:3] ix, id1, id2 = self._sampledTubeCoordinates[1][0:3] rx, rd1, rd2, rd3 = [], [], [], [] @@ -415,15 +440,15 @@ def sample(self, fixedElementsCountAlong, targetElementLength): itx, itd1, itd2 = ix[n2], id1[n2], id2[n2] # wx, wd3 = self._determineWallCoordinates(otx, otd1, otd2, itx, itd1, itd2, coreCentre, arcCentre) wd3 = [mult(sub(otx[n1], itx[n1]), wallFactor) for n1 in range(self._elementsCountAround)] - for n3 in range(self._elementsCountThroughWall + 1): - oFactor = n3 / self._elementsCountThroughWall + for n3 in range(self._elementsCountThroughShell + 1): + oFactor = n3 / self._elementsCountThroughShell iFactor = 1.0 - oFactor for r in (rx, rd1, rd2, rd3): r[n2].append([]) for n1 in range(self._elementsCountAround): if n3 == 0: x, d1, d2 = itx[n1], itd1[n1], itd2[n1] - elif n3 == self._elementsCountThroughWall: + elif n3 == self._elementsCountThroughShell: x, d1, d2 = otx[n1], otd1[n1], otd2[n1] else: x = add(mult(itx[n1], iFactor), mult(otx[n1], oFactor)) @@ -883,22 +908,22 @@ def _determineWallCoordinates(self, ox, od1, od2, ix, id1, id2, coreCentre, arcC it = cross(ic, id1[n1]) else: ot, it = cross(od1[n1], od2[n1]), cross(id1[n1], id2[n1]) - scalefactor = magnitude(sub(ox[n1], ix[n1])) / self._elementsCountThroughWall + scalefactor = magnitude(sub(ox[n1], ix[n1])) / self._elementsCountThroughShell od3 = mult(normalize(ot), scalefactor) id3 = mult(normalize(it), scalefactor) else: - wallFactor = 1.0 / self._elementsCountThroughWall + wallFactor = 1.0 / self._elementsCountThroughShell od3 = id3 = mult(sub(ox[n1], ix[n1]), wallFactor) txm, td3m, pe, pxi, psf = sampleCubicHermiteCurves( - [ix[n1], ox[n1]], [id3, od3], self._elementsCountThroughWall, arcLengthDerivatives=True) + [ix[n1], ox[n1]], [id3, od3], self._elementsCountThroughShell, arcLengthDerivatives=True) td3m = smoothCubicHermiteDerivativesLine(txm, td3m, fixStartDirection=True, fixEndDirection=True) tx.append(txm) td3.append(td3m) - for n3 in range(self._elementsCountThroughWall + 1): + for n3 in range(self._elementsCountThroughShell + 1): wx.append([]) wd3.append([]) for n1 in range(self._elementsCountAround): @@ -1333,7 +1358,7 @@ def addShellElementsToMeshGroup(self, meshGroup): :param meshGroup: Zinc MeshGroup to add elements to. """ elementsCountRim = self.getElementsCountRim() - elementsCountShell = self._elementsCountThroughWall + elementsCountShell = self._elementsCountThroughShell e3ShellStart = elementsCountRim - elementsCountShell self._addRimElementsToMeshGroup(0, self._elementsCountAround, e3ShellStart, elementsCountRim, meshGroup) @@ -2975,6 +3000,7 @@ def generateMesh(self, generateData: TubeNetworkMeshGenerateData): nids = [] nodeParameters = [] nodeLayouts = [] + elementIdentifier = generateData.nextElementIdentifier() lastTransition = self._isCore and (elementsCountTransition == (rim_e3 + 1)) needParameters = (e3 == 0) or lastTransition for n3 in [e3, e3 + 1] if (meshDimension == 3) else [e3]: @@ -3016,7 +3042,6 @@ def generateMesh(self, generateData: TubeNetworkMeshGenerateData): eft, scalefactors, nodeParameters, [5, 6, 7, 8], Node.VALUE_LABEL_D_DS3) elementtemplate.defineField(coordinates, -1, eft) - elementIdentifier = generateData.nextElementIdentifier() element = mesh.createElement(elementIdentifier, elementtemplate) element.setNodesByIdentifier(eft, nids) if scalefactors: @@ -3027,46 +3052,57 @@ def generateMesh(self, generateData: TubeNetworkMeshGenerateData): class TubeNetworkMeshBuilder(NetworkMeshBuilder): + """ + Builds contiguous tube network meshes with smooth element size transitions at junctions, optionally with solid core. + """ def __init__(self, networkMesh: NetworkMesh, targetElementDensityAlongLongestSegment: float, - defaultElementsCountAround: int, elementsCountThroughWall: int, layoutAnnotationGroups: list = [], - annotationElementsCountsAlong: list = [], annotationElementsCountsAround: list = [], - defaultElementsCountCoreBoxMinor: int = 2, elementsCountTransition: int = 1, - annotationElementsCountsCoreBoxMinor: list = [], isCore=False, useOuterTrimSurfaces=False): - """ + layoutAnnotationGroups: list=[], annotationElementsCountsAlong: list=[], + defaultElementsCountAround: int=8, annotationElementsCountsAround: list=[], + elementsCountThroughShell: int=1, isCore=False, elementsCountTransition: int=1, + defaultElementsCountCoreBoxMinor: int=2, annotationElementsCountsCoreBoxMinor: list=[], + useOuterTrimSurfaces=True): + """ + Builds contiguous tube network meshes with smooth element size transitions at junctions, optionally with solid + core. :param networkMesh: Description of the topology of the network layout. - :param targetElementDensityAlongLongestSegment: - :param defaultElementsCountAround: - :param elementsCountThroughWall: - :param layoutAnnotationGroups: - :param annotationElementsCountsAlong: List in same order as layoutAnnotationGroups, specifying fixed - number along segment with any elements in the annotation group. Client must ensure exclusive map from segments. - Groups with zero value or past end of this list use the targetElementDensityAlongLongestSegment. - :param annotationElementsCountsAround: List in same order as layoutAnnotationGroups, specifying fixed - number around segment with any elements in the annotation group. Client must ensure exclusive map from segments. - Groups with zero value or past end of this list use the defaultElementsCountAround. + :param targetElementDensityAlongLongestSegment: Real value which longest segment path in network is divided by + to get target element length, which is used to determine numbers of elements along except when set for a segment + through annotationElementsCountsAlong. + :param layoutAnnotationGroups: Annotation groups defined on the layout to mirror on the final mesh. + :param annotationElementsCountsAlong: List in same order as layoutAnnotationGroups, specifying fixed number of + elements along segment with any elements in the annotation group. Client must ensure exclusive map from + segments. Groups with zero value or past end of this list use the targetElementDensityAlongLongestSegment. + :param defaultElementsCountAround: Number of elements around segments to use unless overridden by + annotationElementsCountsAround. + :param annotationElementsCountsAround: List in same order as layoutAnnotationGroups, specifying fixed number of + elements around segments with any elements in the annotation group. Client must ensure exclusive map from + segments. Groups with zero value or past end of this list use the defaultElementsCountAround. + :param elementsCountThroughShell: Number of elements through shell/wall >= 1. + :param isCore: Set to True to define solid core box and transition elements. + :param elementsCountTransition: Number of rows of elements transitioning between core box and shell >= 1. :param defaultElementsCountCoreBoxMinor: Number of elements across the core box in the minor/d3 direction. - :param elementsCountTransition: :param annotationElementsCountsCoreBoxMinor: List in same order as layoutAnnotationGroups, specifying numbers of - elements across core box minor/d3 direction. - :param isCore: Set to True to define solid core box and transition elements. - :param useOuterTrimSurfaces: Set to True to use common trim surfaces calculated from outer. + elements across core box minor/d3 direction for segments with any elements in the annotation group. Client must + ensure exclusive map from segments. Groups with zero value or past end of this list use the + defaultElementsCountCoreBoxMinor. + :param useOuterTrimSurfaces: Set to False to use separate trim surfaces on inner and outer tubes. Ignored if + no inner path. """ super(TubeNetworkMeshBuilder, self).__init__( networkMesh, targetElementDensityAlongLongestSegment, layoutAnnotationGroups, annotationElementsCountsAlong) self._defaultElementsCountAround = defaultElementsCountAround - self._elementsCountThroughWall = elementsCountThroughWall - self._layoutAnnotationGroups = layoutAnnotationGroups self._annotationElementsCountsAround = annotationElementsCountsAround + self._elementsCountThroughShell = elementsCountThroughShell + self._isCore = isCore + self._elementsCountTransition = elementsCountTransition + self._defaultElementsCountCoreBoxMinor = defaultElementsCountCoreBoxMinor + self._annotationElementsCountsCoreBoxMinor = annotationElementsCountsCoreBoxMinor layoutFieldmodule = self._layoutRegion.getFieldmodule() self._layoutInnerCoordinates = layoutFieldmodule.findFieldByName("inner coordinates").castFiniteElement() if not self._layoutInnerCoordinates.isValid(): self._layoutInnerCoordinates = None - self._isCore = isCore - self._useOuterTrimSurfaces = useOuterTrimSurfaces - self._defaultElementsCountCoreBoxMinor = defaultElementsCountCoreBoxMinor - self._elementsCountTransition = elementsCountTransition - self._annotationElementsCountsCoreBoxMinor = annotationElementsCountsCoreBoxMinor + self._useOuterTrimSurfaces = useOuterTrimSurfaces if self._layoutInnerCoordinates else False def createSegment(self, networkSegment): pathParametersList = [get_nodeset_path_ordered_field_parameters( @@ -3102,7 +3138,7 @@ def createSegment(self, networkSegment): i += 1 return TubeNetworkMeshSegment(networkSegment, pathParametersList, elementsCountAround, - self._elementsCountThroughWall, self._isCore, elementsCountCoreBoxMinor, + self._elementsCountThroughShell, self._isCore, elementsCountCoreBoxMinor, self._elementsCountTransition) def createJunction(self, inSegments, outSegments): @@ -3125,6 +3161,40 @@ def generateMesh(self, generateData): segment.addShellElementsToMeshGroup(shellMeshGroup) +class BodyTubeNetworkMeshBuilder(TubeNetworkMeshBuilder): + """ + Specialization of TubeNetworkMeshBuilder adding annotations for left, right, dorsal, ventral regions. + Requires network layout to follow these conventions: + - consistently annotates fully left or right features with names including "left" or "right", respectively. + - along central body, +d2 direction is left, -d2 direction is right. + - +d3 direction is ventral, -d3 is dorsal. + """ + + def generateMesh(self, generateData): + super(BodyTubeNetworkMeshBuilder, self).generateMesh(generateData) + # build left, right, dorsal, ventral annotation groups + leftMeshGroup = generateData.getLeftMeshGroup() + rightMeshGroup = generateData.getRightMeshGroup() + dorsalMeshGroup = generateData.getDorsalMeshGroup() + ventralMeshGroup = generateData.getVentralMeshGroup() + for networkSegment in self._networkMesh.getNetworkSegments(): + segment = self._segments[networkSegment] + annotationTerms = segment.getAnnotationTerms() + for annotationTerm in annotationTerms: + if "left" in annotationTerm[0]: + segment.addAllElementsToMeshGroup(leftMeshGroup) + break + if "right" in annotationTerm[0]: + segment.addAllElementsToMeshGroup(rightMeshGroup) + break + else: + # segment on main axis + segment.addSideD2ElementsToMeshGroup(False, leftMeshGroup) + segment.addSideD2ElementsToMeshGroup(True, rightMeshGroup) + segment.addSideD3ElementsToMeshGroup(False, ventralMeshGroup) + segment.addSideD3ElementsToMeshGroup(True, dorsalMeshGroup) + + class TubeEllipseGenerator: """ Generates tube ellipse curves with even-sized elements with specified radius, phase angle, diff --git a/tests/test_network.py b/tests/test_network.py index 6f9763c3..1679c02b 100644 --- a/tests/test_network.py +++ b/tests/test_network.py @@ -133,8 +133,6 @@ def test_2d_tube_network_snake(self): """ scaffoldPackage = ScaffoldPackage(MeshType_2d_tubenetwork1, defaultParameterSetName="Snake") settings = scaffoldPackage.getScaffoldSettings() - networkLayoutScaffoldPackage = settings["Network layout"] - networkLayoutSettings = networkLayoutScaffoldPackage.getScaffoldSettings() self.assertEqual(12.0, settings["Target element density along longest segment"]) MeshType_2d_tubenetwork1.checkOptions(settings) @@ -177,7 +175,7 @@ def test_2d_tube_network_snake(self): def test_2d_tube_network_sphere_cube(self): """ - Test 2D bifurcation is generated correctly. + Test 2D sphere cube is generated correctly. """ scaffoldPackage = ScaffoldPackage(MeshType_2d_tubenetwork1, defaultParameterSetName="Sphere cube") settings = scaffoldPackage.getScaffoldSettings() @@ -297,8 +295,6 @@ def test_2d_tube_network_vase(self): """ scaffoldPackage = ScaffoldPackage(MeshType_2d_tubenetwork1, defaultParameterSetName="Vase") settings = scaffoldPackage.getScaffoldSettings() - networkLayoutScaffoldPackage = settings["Network layout"] - networkLayoutSettings = networkLayoutScaffoldPackage.getScaffoldSettings() self.assertEqual(12.0, settings["Target element density along longest segment"]) MeshType_2d_tubenetwork1.checkOptions(settings) @@ -368,7 +364,6 @@ def test_3d_tube_network_bifurcation(self): scaffoldPackage.generate(region) fieldmodule = region.getFieldmodule() - mesh3d = fieldmodule.findMeshByDimension(3) mesh3d = fieldmodule.findMeshByDimension(3) self.assertEqual(8 * 4 * 3, mesh3d.getSize()) @@ -445,7 +440,7 @@ def test_3d_tube_network_bifurcation_core(self): self.assertEqual((8 * 4 * 3) * 2 + (4 * 4 * 3), mesh3d.getSize()) nodes = fieldmodule.findNodesetByFieldDomainType(Field.DOMAIN_TYPE_NODES) - self.assertEqual((8 * 4 * 3 + 3 * 3 + 2) * 2 + (9 * 4 * 3 + 3 * 4), nodes.getSize()) + self.assertEqual((8 * 4 * 3 + 3 * 3 + 2) * 2 + (9 * 4 * 3 + 3 * 4), nodes.getSize()) coordinates = fieldmodule.findFieldByName("coordinates").castFiniteElement() self.assertTrue(coordinates.isValid()) @@ -526,7 +521,6 @@ def test_3d_tube_network_converging_bifurcation_core(self): self.assertTrue(findAnnotationGroupByName(annotationGroups, "shell") is not None) fieldmodule = region.getFieldmodule() - mesh3d = fieldmodule.findMeshByDimension(3) mesh3d = fieldmodule.findMeshByDimension(3) self.assertEqual(336, mesh3d.getSize()) @@ -847,7 +841,6 @@ def test_3d_tube_network_sphere_cube_core(self): self.assertEqual(result, RESULT_OK) self.assertAlmostEqual(volume, expectedSizes3d[name][1], delta=X_TOL) - def test_3d_tube_network_trifurcation_cross(self): """ Test trifurcation cross 3-D tube network is generated correctly with variable elements count around. @@ -892,7 +885,6 @@ def test_3d_tube_network_trifurcation_cross(self): self.assertEqual(1, len(annotationGroups)) fieldmodule = region.getFieldmodule() - mesh3d = fieldmodule.findMeshByDimension(3) mesh3d = fieldmodule.findMeshByDimension(3) self.assertEqual(144, mesh3d.getSize()) @@ -1038,8 +1030,6 @@ def test_3d_box_network_bifurcation(self): """ scaffoldPackage = ScaffoldPackage(MeshType_3d_boxnetwork1, defaultParameterSetName="Bifurcation") settings = scaffoldPackage.getScaffoldSettings() - networkLayoutScaffoldPackage = settings["Network layout"] - networkLayoutSettings = networkLayoutScaffoldPackage.getScaffoldSettings() self.assertEqual(3, len(settings)) self.assertEqual(4.0, settings["Target element density along longest segment"]) self.assertEqual([0], settings["Annotation numbers of elements along"]) @@ -1221,7 +1211,6 @@ def test_3d_tube_network_loop_core(self): self.assertAlmostEqual(volume, 0.0982033864405135, delta=1.0E-6) self.assertAlmostEqual(surfaceArea, 1.9683574196198823, delta=1.0E-6) - def test_3d_tube_network_loop_two_segments(self): """ Test loop 3-D tube network is generated with 2 segments with fixed element boundary between them. diff --git a/tests/test_uterus.py b/tests/test_uterus.py index 50b85d81..12c97e3e 100644 --- a/tests/test_uterus.py +++ b/tests/test_uterus.py @@ -134,8 +134,8 @@ def test_uterus1(self): for annotationGroup in annotationGroups: annotationGroup.addSubelements() scaffold.defineFaceAnnotations(refineRegion, options, annotationGroups) - for annotation in annotationGroups: - if annotation not in oldAnnotationGroups: + for annotationGroup in annotationGroups: + if annotationGroup not in oldAnnotationGroups: annotationGroup.addSubelements() self.assertEqual(34, len(annotationGroups)) #