Search in sources :

Example 1 with Tree

use of cbit.util.graph.Tree in project vcell by virtualcell.

the class RegionImage method calculateRegions.

/**
 * Insert the method's description here.
 * Creation date: (3/28/2002 10:30:20 AM)
 * @param image cbit.image.VCImage
 */
private void calculateRegions(VCImage vcImage) throws cbit.image.ImageException {
    long time1 = System.currentTimeMillis();
    RegionMask[][] regionMasks = new RegionMask[numZ][];
    for (int k = 0; k < vcImage.getNumZ(); k++) {
        regionMasks[k] = calculateRegions3D(vcImage.getPixels(), k * numXY, numX, numY);
    }
    long time2 = System.currentTimeMillis();
    // time2 = System.currentTimeMillis();
    // RegionMask regionMasksFast[][] = new RegionMask[numZ][];
    // for (int k = 0; k < vcImage.getNumZ(); k++){
    // regionMasksFast[k] = calculateRegions3Dfaster(vcImage.getPixels(),k*numXY,numX,numY);
    // }
    // long time3 = System.currentTimeMillis();
    // System.out.println("4way recursive took "+((time2-time1)/1000.0)+" s, nonrecursive line took "+((time3-time2)/1000.0)+" s");
    // 
    // consolidate "off-plane contiguous" region indexes
    // build a graph of all 2D regions that touch
    // 
    // the final 3D regions are the set of nodes in each of the spanning trees.
    // 
    // 
    // build graph with 1 node per 2d region
    // 
    Graph connectionGraph = new Graph();
    for (int k = 0; k < numZ; k++) {
        for (int i = 0; i < regionMasks[k].length; i++) {
            Node node = new Node(k + "," + i, new org.vcell.util.ObjectReferenceWrapper(regionMasks[k][i]));
            connectionGraph.addNode(node);
        }
    }
    // 
    // add edges for any slice-slice touching of regions with same pixel value
    // 
    BitSet zeroBitSet = new BitSet();
    BitSet intersection = new BitSet(numXY);
    for (int k = 0; k < numZ - 1; k++) {
        for (int i = 0; i < regionMasks[k].length; i++) {
            Node node_i = connectionGraph.getNode(k + "," + i);
            for (int j = 0; j < regionMasks[k + 1].length; j++) {
                if (regionMasks[k][i].pixelValue == regionMasks[k + 1][j].pixelValue) {
                    // clear mask
                    intersection.and(zeroBitSet);
                    intersection.or(regionMasks[k][i].mask);
                    intersection.and(regionMasks[k + 1][j].mask);
                    if (!intersection.equals(zeroBitSet)) {
                        Node node_j = connectionGraph.getNode((k + 1) + "," + j);
                        connectionGraph.addEdge(new Edge(node_i, node_j));
                    }
                }
            }
        }
    }
    // 
    // get spanning forest, and for each spanning tree, assign a single regionID to each contained 2d mask
    // 
    Tree[] spanningForest = connectionGraph.getSpanningForest();
    regionInfos = new RegionInfo[spanningForest.length];
    for (int i = 0; i < spanningForest.length; i++) {
        Node[] nodes = spanningForest[i].getNodes();
        int pixelValue = -1;
        int numPixels = 0;
        BitSet fullMask = new BitSet(numXY * numZ);
        for (int j = 0; j < nodes.length; j++) {
            RegionMask regionMask = (RegionMask) ((ObjectReferenceWrapper) nodes[j].getData()).getObject();
            pixelValue = regionMask.pixelValue;
            for (int k = 0; k < numXY; k++) {
                if (regionMask.mask.get(k)) {
                    fullMask.set(regionMask.offset + k);
                    numPixels++;
                }
            }
        }
        regionInfos[i] = new RegionInfo(pixelValue, numPixels, i, fullMask);
    }
    long time3 = System.currentTimeMillis();
    System.out.println("4way recursive on slices took " + ((time2 - time1) / 1000.0) + " s, total RegionImage time " + ((time3 - time1) / 1000.0) + " s");
}
Also used : Node(cbit.util.graph.Node) BitSet(java.util.BitSet) ObjectReferenceWrapper(org.vcell.util.ObjectReferenceWrapper) Graph(cbit.util.graph.Graph) Tree(cbit.util.graph.Tree) Edge(cbit.util.graph.Edge)

Example 2 with Tree

use of cbit.util.graph.Tree in project vcell by virtualcell.

the class RayCasterBitSet method createGeometry.

public static Geometry createGeometry(GeometryThumbnailImageFactory geometryThumbnailImageFactory, SurfaceCollection surfaceCollection, Origin origin, Extent extent, ISize sampleSize) throws ImageException, PropertyVetoException, GeometryException, ExpressionException {
    int numX = sampleSize.getX();
    int numY = sampleSize.getY();
    int numZ = sampleSize.getZ();
    VolumeSamplesBitSet volumeSamples = volumeSampleSurface(surfaceCollection, sampleSize, origin, extent, false);
    // 
    // create a graph with:
    // a node for each surface face (mask)
    // an edge for each pair of surface faces share by a volume node.
    // 
    Graph connectivityGraph = new Graph();
    ArrayList<Long> masks = new ArrayList<Long>(volumeSamples.getIncidentSurfaceBitSets().keySet());
    for (int i = 0; i < masks.size(); i++) {
        cbit.util.graph.Node node = new cbit.util.graph.Node(masks.get(i).toString(), volumeSamples.getIncidentSurfaceBitSets().get(masks.get(i)));
        connectivityGraph.addNode(node);
    }
    for (int i = 0; i < masks.size(); i++) {
        long mask1 = masks.get(i);
        cbit.util.graph.Node node1 = connectivityGraph.getNode(Long.toString(mask1));
        BitSet bitSet1 = volumeSamples.getIncidentSurfaceBitSets().get(mask1);
        for (int j = i + 1; j < masks.size(); j++) {
            long mask2 = masks.get(j);
            cbit.util.graph.Node node2 = connectivityGraph.getNode(Long.toString(mask2));
            BitSet bitSet2 = volumeSamples.getIncidentSurfaceBitSets().get(mask2);
            if (bitSet1.intersects(bitSet2)) {
                BitSet bitSetIntersection = new BitSet();
                bitSetIntersection.or(bitSet1);
                bitSetIntersection.and(bitSet2);
                Edge edge = new Edge(node1, node2, bitSetIntersection);
                connectivityGraph.addEdge(edge);
            }
        }
    }
    // 
    for (int i = 0; i < surfaceCollection.getSurfaceCount(); i++) {
        Surface surface = surfaceCollection.getSurfaces(i);
        Tree[] spanningForest = connectivityGraph.getSpanningForest();
        for (Tree tree : spanningForest) {
            cbit.util.graph.Node insideNode = tree.getNode(Long.toString(surface.getInteriorMask()));
            cbit.util.graph.Node outsideNode = tree.getNode(Long.toString(surface.getExteriorMask()));
            if (insideNode != null && outsideNode != null) {
                System.out.println("++++++++++++++++++++++++++++++++ separating " + Long.toBinaryString(Long.parseLong(insideNode.getName())) + " with " + Long.toBinaryString(Long.parseLong(outsideNode.getName())));
                // 
                // remove sub-region covered by more loosely connected (insideNode or outsideNode) from this ensemble region
                // this could be improved by considering a minimal cut-set (too hard for today).
                // 
                int connectionCardinalityInsideNode = 0;
                Edge[] insideNodeEdges = connectivityGraph.getAdjacentEdges(insideNode);
                for (Edge insideEdge : insideNodeEdges) {
                    connectionCardinalityInsideNode += ((BitSet) insideEdge.getData()).cardinality();
                }
                int connectionCardinalityOutsideNode = 0;
                Edge[] outsideNodeEdges = connectivityGraph.getAdjacentEdges(outsideNode);
                for (Edge outsideEdge : outsideNodeEdges) {
                    connectionCardinalityOutsideNode += ((BitSet) outsideEdge.getData()).cardinality();
                }
                if (connectionCardinalityInsideNode < connectionCardinalityOutsideNode) {
                    // remove inside Node (clear its bits from all others in this tree)
                    for (cbit.util.graph.Node node : tree.getNodes()) {
                        if (node != insideNode) {
                            BitSet otherBitSet = (BitSet) node.getData();
                            otherBitSet.andNot((BitSet) insideNode.getData());
                        }
                    }
                    for (Edge edge : insideNodeEdges) {
                        connectivityGraph.remove(edge);
                    }
                } else {
                    // remove inside Node (clear its bits from all others in this tree)
                    for (cbit.util.graph.Node node : tree.getNodes()) {
                        if (node != outsideNode) {
                            BitSet otherBitSet = (BitSet) node.getData();
                            otherBitSet.andNot((BitSet) outsideNode.getData());
                        }
                    }
                    for (Edge edge : outsideNodeEdges) {
                        connectivityGraph.remove(edge);
                    }
                }
            }
        }
    }
    Tree[] spanningForest = connectivityGraph.getSpanningForest();
    BitSet[] regionBitSets = new BitSet[spanningForest.length];
    long[] regionMasks = new long[spanningForest.length];
    for (int i = 0; i < spanningForest.length; i++) {
        Tree tree = spanningForest[i];
        regionBitSets[i] = new BitSet();
        for (cbit.util.graph.Node node : tree.getNodes()) {
            regionBitSets[i].or((BitSet) node.getData());
            regionMasks[i] = regionMasks[i] | Long.valueOf(node.getName());
        }
        System.out.println("+++++++++++++++++++++++++++++++++++++ final region mask : " + Long.toBinaryString(regionMasks[i]) + ", region size = " + regionBitSets[i].cardinality());
    }
    byte[] pixels = new byte[numX * numY * numZ];
    byte pixelValue = 1;
    for (int r = 0; r < regionMasks.length; r++) {
        BitSet regionBitSet = regionBitSets[r];
        for (int i = 0; i < pixels.length; i++) {
            if (regionBitSet.get(i)) {
                pixels[i] = pixelValue;
            }
        }
        pixelValue++;
    }
    VCImageUncompressed vcImage = new VCImageUncompressed(null, pixels, extent, numX, numY, numZ);
    Geometry geometry = new Geometry("newGeometry", vcImage);
    geometry.getGeometrySpec().setExtent(extent);
    geometry.getGeometrySpec().setOrigin(origin);
    geometry.precomputeAll(geometryThumbnailImageFactory, true, true);
    return geometry;
}
Also used : ArrayList(java.util.ArrayList) BitSet(java.util.BitSet) VCImageUncompressed(cbit.image.VCImageUncompressed) Geometry(cbit.vcell.geometry.Geometry) Graph(cbit.util.graph.Graph) Tree(cbit.util.graph.Tree) Edge(cbit.util.graph.Edge)

Example 3 with Tree

use of cbit.util.graph.Tree in project vcell by virtualcell.

the class StructureSorter method sortStructures.

public static Structure[] sortStructures(Model model) {
    Graph graph = new Graph();
    // add all structures to the graph
    for (Structure structure : model.getStructures()) {
        graph.addNode(new Node(structure.getName(), structure));
    }
    // add an edge
    for (Structure structure : model.getStructures()) {
        Structure parentStructure = null;
        if (structure.getSbmlParentStructure() != null) {
            parentStructure = structure.getSbmlParentStructure();
        } else if (model.getStructureTopology().getParentStructure(structure) != null) {
            parentStructure = model.getStructureTopology().getParentStructure(structure);
        }
        if (parentStructure != null) {
            int index1 = graph.getIndex(graph.getNode(structure.getName()));
            int index2 = graph.getIndex(graph.getNode(parentStructure.getName()));
            StructureRelationshipEdge edge = getStructureRelationshipEdge(graph, index1, index2);
            edge.bParentChild = true;
        }
    }
    for (ReactionStep rs : model.getReactionSteps()) {
        int index1 = graph.getIndex(graph.getNode(rs.getStructure().getName()));
        for (ReactionParticipant rp : rs.getReactionParticipants()) {
            if (rp.getStructure() != rs.getStructure()) {
                int index2 = graph.getIndex(graph.getNode(rp.getStructure().getName()));
                StructureRelationshipEdge edge = getStructureRelationshipEdge(graph, index2, index1);
                edge.numConnections++;
            }
        }
    }
    // remove edges such that all nodes have degree <= 2
    for (Node node : graph.getNodes()) {
        Edge[] adjacentEdges = graph.getAdjacentEdges(node);
        while (adjacentEdges.length > 2) {
            StructureRelationshipEdge lowestCostEdge = getLowestCostEdge(Arrays.asList(adjacentEdges));
            graph.remove(lowestCostEdge);
            adjacentEdges = graph.getAdjacentEdges(node);
        }
    }
    // for each dependency loop, remove edge with lowest cost
    Path[] fundamentalCycles = graph.getFundamentalCycles();
    while (fundamentalCycles.length > 0) {
        StructureRelationshipEdge lowestCostEdge = getLowestCostEdge(Arrays.asList(fundamentalCycles[0].getEdges()));
        graph.remove(lowestCostEdge);
        fundamentalCycles = graph.getFundamentalCycles();
    }
    // all graphs in the forest are linear now ... find path from first to last and line up all trees in list.
    Tree[] spanningForest = graph.getSpanningForest();
    ArrayList<Structure> structures = new ArrayList<Structure>();
    for (Tree tree : spanningForest) {
        if (tree.getNodes().length > 1) {
            Node start = null;
            Node end = null;
            for (Node node : tree.getNodes()) {
                if (tree.getDegree(node) == 1) {
                    if (start == null) {
                        start = node;
                    } else if (end == null) {
                        end = node;
                    } else {
                        break;
                    }
                }
            }
            Path path = tree.getTreePath(start, end);
            for (Node node : path.getNodesTraversed()) {
                structures.add((Structure) node.getData());
            }
        } else {
            structures.add((Structure) tree.getNodes()[0].getData());
        }
    }
    return structures.toArray(new Structure[structures.size()]);
}
Also used : Path(cbit.util.graph.Path) Node(cbit.util.graph.Node) ArrayList(java.util.ArrayList) Graph(cbit.util.graph.Graph) Tree(cbit.util.graph.Tree) Edge(cbit.util.graph.Edge)

Example 4 with Tree

use of cbit.util.graph.Tree in project vcell by virtualcell.

the class ModelOptimizationSpec method removeUncoupledParameters.

public void removeUncoupledParameters() {
    try {
        localIssueList.clear();
        MathMapping mathMapping = getSimulationContext().createNewMathMapping();
        MathDescription mathDesc = mathMapping.getMathDescription();
        MathSystemHash mathSystemHash = fromMath(mathDesc);
        Graph graph = mathSystemHash.getDependencyGraph(mathSystemHash.getSymbols());
        Tree[] spanningTrees = graph.getSpanningForest();
        // 
        for (int i = 0; i < spanningTrees.length; i++) {
            Node[] treeNodes = spanningTrees[i].getNodes();
            boolean bHasStateVariables = false;
            for (int j = 0; j < treeNodes.length; j++) {
                Node node = treeNodes[j];
                Variable var = mathDesc.getVariable(node.getName());
                if (var instanceof VolVariable || var instanceof MemVariable || var instanceof FilamentVariable || var instanceof VolumeRegionVariable || var instanceof MembraneRegionVariable || var instanceof FilamentRegionVariable) {
                    bHasStateVariables = true;
                    break;
                }
            }
            if (!bHasStateVariables) {
                spanningTrees = (Tree[]) BeanUtils.removeElement(spanningTrees, spanningTrees[i]);
                i--;
            }
        }
        // 
        // remove parameters not mapped to a surviving tree (not coupled to any state variables
        // 
        ArrayList<ParameterMappingSpec> paramMappingSpecsList = new ArrayList<ParameterMappingSpec>();
        paramMappingSpecsList.addAll(Arrays.asList(fieldParameterMappingSpecs));
        for (int i = 0; i < paramMappingSpecsList.size(); i++) {
            Parameter parameter = paramMappingSpecsList.get(i).getModelParameter();
            String mathName = mathMapping.getMathSymbolMapping().getVariable(parameter).getName();
            boolean bFoundInTree = false;
            for (int j = 0; j < spanningTrees.length; j++) {
                Node node = spanningTrees[j].getNode(mathName);
                if (node != null) {
                    bFoundInTree = true;
                }
            }
            if (!bFoundInTree) {
                paramMappingSpecsList.remove(i);
                i--;
            }
        }
        ParameterMappingSpec[] parameterMappingSpecs = new ParameterMappingSpec[paramMappingSpecsList.size()];
        paramMappingSpecsList.toArray(parameterMappingSpecs);
        setParameterMappingSpecs(parameterMappingSpecs);
    } catch (Exception e) {
        e.printStackTrace(System.out);
        localIssueList.add(new Issue(this, localIssueContext, IssueCategory.ParameterEstimationGeneralWarning, e.getMessage(), Issue.SEVERITY_WARNING));
    // throw new RuntimeException(e.getMessage());
    }
}
Also used : MembraneRegionVariable(cbit.vcell.math.MembraneRegionVariable) FilamentVariable(cbit.vcell.math.FilamentVariable) VolVariable(cbit.vcell.math.VolVariable) VolumeRegionVariable(cbit.vcell.math.VolumeRegionVariable) MemVariable(cbit.vcell.math.MemVariable) FilamentRegionVariable(cbit.vcell.math.FilamentRegionVariable) Variable(cbit.vcell.math.Variable) Issue(org.vcell.util.Issue) MembraneRegionVariable(cbit.vcell.math.MembraneRegionVariable) MathDescription(cbit.vcell.math.MathDescription) Node(cbit.util.graph.Node) ArrayList(java.util.ArrayList) VolumeRegionVariable(cbit.vcell.math.VolumeRegionVariable) MemVariable(cbit.vcell.math.MemVariable) Tree(cbit.util.graph.Tree) VolVariable(cbit.vcell.math.VolVariable) FilamentRegionVariable(cbit.vcell.math.FilamentRegionVariable) PropertyVetoException(java.beans.PropertyVetoException) ExpressionException(cbit.vcell.parser.ExpressionException) Graph(cbit.util.graph.Graph) FilamentVariable(cbit.vcell.math.FilamentVariable) MathMapping(cbit.vcell.mapping.MathMapping) Parameter(cbit.vcell.model.Parameter) KineticsParameter(cbit.vcell.model.Kinetics.KineticsParameter) ModelParameter(cbit.vcell.model.Model.ModelParameter) KineticsProxyParameter(cbit.vcell.model.Kinetics.KineticsProxyParameter) SpeciesContextSpecParameter(cbit.vcell.mapping.SpeciesContextSpec.SpeciesContextSpecParameter)

Aggregations

Graph (cbit.util.graph.Graph)4 Tree (cbit.util.graph.Tree)4 Edge (cbit.util.graph.Edge)3 Node (cbit.util.graph.Node)3 ArrayList (java.util.ArrayList)3 BitSet (java.util.BitSet)2 VCImageUncompressed (cbit.image.VCImageUncompressed)1 Path (cbit.util.graph.Path)1 Geometry (cbit.vcell.geometry.Geometry)1 MathMapping (cbit.vcell.mapping.MathMapping)1 SpeciesContextSpecParameter (cbit.vcell.mapping.SpeciesContextSpec.SpeciesContextSpecParameter)1 FilamentRegionVariable (cbit.vcell.math.FilamentRegionVariable)1 FilamentVariable (cbit.vcell.math.FilamentVariable)1 MathDescription (cbit.vcell.math.MathDescription)1 MemVariable (cbit.vcell.math.MemVariable)1 MembraneRegionVariable (cbit.vcell.math.MembraneRegionVariable)1 Variable (cbit.vcell.math.Variable)1 VolVariable (cbit.vcell.math.VolVariable)1 VolumeRegionVariable (cbit.vcell.math.VolumeRegionVariable)1 KineticsParameter (cbit.vcell.model.Kinetics.KineticsParameter)1