Examples with VolumeGeometricRegion cbit.vcell.geometry.surface.VolumeGeometricRegion used on opensource projects

Example 16 with VolumeGeometricRegion

use of cbit.vcell.geometry.surface.VolumeGeometricRegion in project vcell by virtualcell.

the class SurfaceRegionObject method gatherIssues.

@Override
public void gatherIssues(IssueContext issueContext, List<Issue> issueList) {
    if (simulationContext.getGeometry().getGeometrySurfaceDescription() != null) {
        GeometricRegion[] regions = simulationContext.getGeometry().getGeometrySurfaceDescription().getGeometricRegions();
        boolean bFound = false;
        if (regions != null) {
            for (GeometricRegion region : regions) {
                if (region instanceof SurfaceGeometricRegion) {
                    SurfaceGeometricRegion sr = (SurfaceGeometricRegion) region;
                    if (sr.getAdjacentGeometricRegions() != null && sr.getAdjacentGeometricRegions().length == 2) {
                        VolumeGeometricRegion vr1 = (VolumeGeometricRegion) sr.getAdjacentGeometricRegions()[0];
                        VolumeGeometricRegion vr2 = (VolumeGeometricRegion) sr.getAdjacentGeometricRegions()[1];
                        if (vr1.getSubVolume() == insideSubVolume && vr1.getRegionID() == insideRegionID && vr2.getSubVolume() == outsideSubVolume && vr2.getRegionID() == outsideRegionID) {
                            bFound = true;
                        }
                        if (vr1.getSubVolume() == outsideSubVolume && vr1.getRegionID() == outsideRegionID && vr2.getSubVolume() == insideSubVolume && vr2.getRegionID() == insideRegionID) {
                            bFound = true;
                        }
                    }
                }
            }
        }
        if (!bFound) {
            issueList.add(new Issue(this, issueContext, IssueCategory.Identifiers, "could not find corresponding surface region in geometry", Issue.Severity.ERROR));
        }
    }
}

Example 17 with VolumeGeometricRegion

use of cbit.vcell.geometry.surface.VolumeGeometricRegion in project vcell by virtualcell.

the class XmlReader method getGeometrySurfaceDescription.

private GeometrySurfaceDescription getGeometrySurfaceDescription(Element param, Geometry geom) throws XmlParseException {
    GeometrySurfaceDescription gsd = geom.getGeometrySurfaceDescription();
    String cutoffStr = param.getAttributeValue(XMLTags.CutoffFrequencyAttrTag);
    String xDim = param.getAttributeValue(XMLTags.NumSamplesXAttrTag);
    String yDim = param.getAttributeValue(XMLTags.NumSamplesYAttrTag);
    String zDim = param.getAttributeValue(XMLTags.NumSamplesZAttrTag);
    if (cutoffStr == null || xDim == null || yDim == null || zDim == null) {
        throw new XmlParseException("Attributes for element Surface Description not properly set, under geometry: " + ((Element) param.getParent()).getAttributeValue(XMLTags.NameAttrTag));
    }
    try {
        ISize isize = new ISize(Integer.parseInt(xDim), Integer.parseInt(yDim), Integer.parseInt(zDim));
        gsd.setVolumeSampleSize(isize);
        gsd.setFilterCutoffFrequency(new Double(cutoffStr));
        // these lists are allowed to be empty.
        ArrayList<Element> memRegions = new ArrayList<Element>(param.getChildren(XMLTags.MembraneRegionTag, vcNamespace));
        ArrayList<Element> volRegions = new ArrayList<Element>(param.getChildren(XMLTags.VolumeRegionTag, vcNamespace));
        ArrayList<GeometricRegion> regions = new ArrayList<GeometricRegion>();
        GeometryUnitSystem geometryUnitSystem = geom.getUnitSystem();
        for (Element temp : volRegions) {
            String regionID = temp.getAttributeValue(XMLTags.RegionIDAttrTag);
            String name = temp.getAttributeValue(XMLTags.NameAttrTag);
            String subvolumeRef = temp.getAttributeValue(XMLTags.SubVolumeAttrTag);
            if (regionID == null || name == null || subvolumeRef == null) {
                throw new XmlParseException("Attributes for element Volume Region not properly set, under geometry: " + ((Element) param.getParent()).getAttributeValue(XMLTags.NameAttrTag));
            }
            SubVolume subvolume = geom.getGeometrySpec().getSubVolume(subvolumeRef);
            if (subvolume == null) {
                throw new XmlParseException("The subvolume " + subvolumeRef + " could not be resolved.");
            }
            double size = -1;
            VCUnitDefinition unit = null;
            String sizeStr = temp.getAttributeValue(XMLTags.SizeAttrTag);
            if (sizeStr != null) {
                size = Double.parseDouble(sizeStr);
                String unitSymbol = temp.getAttributeValue(XMLTags.VCUnitDefinitionAttrTag);
                if (unitSymbol != null) {
                    unit = geometryUnitSystem.getInstance(unitSymbol);
                }
            }
            VolumeGeometricRegion vgr = new VolumeGeometricRegion(name, size, unit, subvolume, Integer.parseInt(regionID));
            regions.add(vgr);
        }
        for (Element temp : memRegions) {
            String volRegion_1 = temp.getAttributeValue(XMLTags.VolumeRegion_1AttrTag);
            String volRegion_2 = temp.getAttributeValue(XMLTags.VolumeRegion_2AttrTag);
            String name = temp.getAttributeValue(XMLTags.NameAttrTag);
            if (volRegion_1 == null || volRegion_2 == null || name == null) {
                throw new XmlParseException("Attributes for element Membrane Region not properly set, under geometry: " + ((Element) param.getParent()).getAttributeValue(XMLTags.NameAttrTag));
            }
            VolumeGeometricRegion region1 = getAdjacentVolumeRegion(regions, volRegion_1);
            VolumeGeometricRegion region2 = getAdjacentVolumeRegion(regions, volRegion_2);
            if (region1 == null || region2 == null) {
                throw new XmlParseException("Element Membrane Region refernces invalid volume regions, under geometry: " + ((Element) param.getParent()).getAttributeValue(XMLTags.NameAttrTag));
            }
            double size = -1;
            VCUnitDefinition unit = null;
            String sizeStr = temp.getAttributeValue(XMLTags.SizeAttrTag);
            if (sizeStr != null) {
                size = Double.parseDouble(sizeStr);
                String unitSymbol = temp.getAttributeValue(XMLTags.VCUnitDefinitionAttrTag);
                if (unitSymbol != null) {
                    unit = geometryUnitSystem.getInstance(unitSymbol);
                }
            }
            SurfaceGeometricRegion rsl = new SurfaceGeometricRegion(name, size, unit);
            rsl.addAdjacentGeometricRegion(region1);
            region1.addAdjacentGeometricRegion(rsl);
            rsl.addAdjacentGeometricRegion(region2);
            region2.addAdjacentGeometricRegion(rsl);
            regions.add(rsl);
        }
        if (regions.size() > 0) {
            gsd.setGeometricRegions((GeometricRegion[]) regions.toArray(new GeometricRegion[regions.size()]));
        }
    } catch (Exception e) {
        System.err.println("Unable to read geometry surface description from XML, for geometry: " + ((Element) param.getParent()).getAttributeValue(XMLTags.NameAttrTag));
        e.printStackTrace();
    }
    return gsd;
}

Example 18 with VolumeGeometricRegion

use of cbit.vcell.geometry.surface.VolumeGeometricRegion in project vcell by virtualcell.

the class GeometryFileWriter method write.

/**
 * Insert the method's description here.
 * Creation date: (7/19/2004 10:54:30 AM)
 * @param geometrySurfaceDescription cbit.vcell.geometry.surface.GeometrySurfaceDescription
 * @throws IOException
 */
public static void write(Writer writer, Geometry resampledGeometry) throws IOException {
    // 
    // "name" name
    // "dimension" dimension
    // "extent" extentx extenty extentz
    // "origin" originx originy originz
    // "volumeRegions" num
    // name totalVolume featureHandle
    // "membraneRegions" num
    // name totalArea volumeRegionIndex1 volumeRegionIndex2
    // "volumeSamples" numX, numY, numZ
    // uncompressed regionIndexs for each volume element
    // compressed regionIndexs for each volume element
    // "nodes" num
    // nodeIndex x y z
    // "cells" num
    // cellIndex patchIndex node1 node2 node3 node4
    // "celldata"
    // insideVolumeIndex outsideVolumeIndex area normalx normaly normalz
    // 
    // 
    // When we are writing volume regions, we sort regions so that ID is equal to index
    // 
    writer.write("name " + resampledGeometry.getName() + "\n");
    writer.write("dimension " + resampledGeometry.getDimension() + "\n");
    org.vcell.util.Extent extent = resampledGeometry.getExtent();
    org.vcell.util.Origin origin = resampledGeometry.getOrigin();
    switch(resampledGeometry.getDimension()) {
        case 1:
            writer.write("size " + extent.getX() + "\n");
            writer.write("origin " + origin.getX() + "\n");
            break;
        case 2:
            writer.write("size " + extent.getX() + " " + extent.getY() + "\n");
            writer.write("origin " + origin.getX() + " " + origin.getY() + "\n");
            break;
        case 3:
            writer.write("size " + extent.getX() + " " + extent.getY() + " " + extent.getZ() + "\n");
            writer.write("origin " + origin.getX() + " " + origin.getY() + " " + origin.getZ() + "\n");
            break;
    }
    GeometrySurfaceDescription geoSurfaceDesc = resampledGeometry.getGeometrySurfaceDescription();
    RegionImage regionImage = geoSurfaceDesc.getRegionImage();
    SurfaceCollection surfaceCollection = geoSurfaceDesc.getSurfaceCollection();
    GeometricRegion[] geometricRegions = geoSurfaceDesc.getGeometricRegions();
    int numVolumeRegions = 0;
    int numMembraneRegions = 0;
    Vector<VolumeGeometricRegion> volRegionList = new Vector<VolumeGeometricRegion>();
    if (geometricRegions != null) {
        for (int i = 0; i < geometricRegions.length; i++) {
            if (geometricRegions[i] instanceof VolumeGeometricRegion) {
                numVolumeRegions++;
                volRegionList.add((VolumeGeometricRegion) geometricRegions[i]);
            } else if (geometricRegions[i] instanceof SurfaceGeometricRegion) {
                numMembraneRegions++;
            }
        }
    }
    // 
    // get ordered array of volume regions (where "id" == index into array)... fail if impossible
    // 
    java.util.Collections.sort(volRegionList, new Comparator<VolumeGeometricRegion>() {

        public int compare(VolumeGeometricRegion reg1, VolumeGeometricRegion reg2) {
            if (reg1.getRegionID() < reg2.getRegionID()) {
                return -1;
            } else if (reg1.getRegionID() > reg2.getRegionID()) {
                return 1;
            } else {
                return 0;
            }
        }

        public boolean equals(Object obj) {
            return this == obj;
        }
    });
    VolumeGeometricRegion[] volRegions = (VolumeGeometricRegion[]) org.vcell.util.BeanUtils.getArray(volRegionList, VolumeGeometricRegion.class);
    writer.write("volumeRegions " + numVolumeRegions + "\n");
    for (int i = 0; i < volRegions.length; i++) {
        if (volRegions[i].getRegionID() != i) {
            throw new RuntimeException("Region ID != Region Index, they must be the same!");
        }
        writer.write(volRegions[i].getName() + " " + volRegions[i].getSize() + " " + volRegions[i].getSubVolume().getHandle() + "\n");
    }
    writer.write("membraneRegions " + numMembraneRegions + "\n");
    if (geometricRegions != null) {
        for (int i = 0; i < geometricRegions.length; i++) {
            if (geometricRegions[i] instanceof SurfaceGeometricRegion) {
                SurfaceGeometricRegion surfaceRegion = (SurfaceGeometricRegion) geometricRegions[i];
                GeometricRegion[] neighbors = surfaceRegion.getAdjacentGeometricRegions();
                VolumeGeometricRegion insideRegion = (VolumeGeometricRegion) neighbors[0];
                VolumeGeometricRegion outsideRegion = (VolumeGeometricRegion) neighbors[1];
                writer.write(surfaceRegion.getName() + " " + surfaceRegion.getSize() + " " + insideRegion.getRegionID() + " " + outsideRegion.getRegionID() + "\n");
            }
        }
    }
    // 
    // write volume samples
    // 
    ISize volumeSampleSize = geoSurfaceDesc.getVolumeSampleSize();
    switch(resampledGeometry.getDimension()) {
        case 1:
            writer.write("volumeSamples " + volumeSampleSize.getX() + "\n");
            break;
        case 2:
            writer.write("volumeSamples " + volumeSampleSize.getX() + " " + volumeSampleSize.getY() + "\n");
            break;
        case 3:
            writer.write("volumeSamples " + volumeSampleSize.getX() + " " + volumeSampleSize.getY() + " " + volumeSampleSize.getZ() + "\n");
            break;
    }
    // regionImage
    if (regionImage != null) {
        if (regionImage.getNumRegions() > 65536) {
            throw new RuntimeException("cannot process a geometry with more than 65536 volume regions");
        }
        byte[] uncompressedRegionIDs = new byte[2 * regionImage.getNumX() * regionImage.getNumY() * regionImage.getNumZ()];
        for (int i = 0, j = 0; i < uncompressedRegionIDs.length; i += 2, j++) {
            int regindex = regionImage.getRegionInfoFromOffset(j).getRegionIndex();
            uncompressedRegionIDs[i] = (byte) (regindex & 0x000000ff);
            uncompressedRegionIDs[i + 1] = (byte) ((regindex & 0x0000ff00) >> 8);
        }
        ByteArrayOutputStream bos = new ByteArrayOutputStream();
        DeflaterOutputStream dos = new DeflaterOutputStream(bos);
        dos.write(uncompressedRegionIDs, 0, uncompressedRegionIDs.length);
        dos.close();
        byte[] compressedRegionIDs = bos.toByteArray();
        writer.write(org.vcell.util.Hex.toString(compressedRegionIDs) + "\n");
    } else {
        writer.write("\n");
    }
    // 
    if (surfaceCollection == null) {
        throw new RuntimeException("geometry is not updated");
    }
    int numCells = surfaceCollection.getTotalPolygonCount();
    writer.write("cells " + numCells + "\n");
    // "celldata"
    // insideVolumeIndex outsideVolumeIndex area normalx normaly normalz
    // 
    int cellID = 0;
    int dimension = resampledGeometry.getDimension();
    double correctCoeff = 1;
    if (dimension == 1) {
        correctCoeff = extent.getY() * extent.getZ();
    } else if (dimension == 2) {
        correctCoeff = extent.getZ();
    }
    if (surfaceCollection != null) {
        for (int i = 0; i < surfaceCollection.getSurfaceCount(); i++) {
            Surface surface = surfaceCollection.getSurfaces(i);
            int region1Outside = 0;
            int region1Inside = 0;
            for (int j = 0; j < surface.getPolygonCount(); j++) {
                Quadrilateral polygon = (Quadrilateral) surface.getPolygons(j);
                Node[] node = polygon.getNodes();
                cbit.vcell.render.Vect3d elementCoord = new cbit.vcell.render.Vect3d();
                int nodesOnBoundary = 0;
                for (int k = 0; k < node.length; k++) {
                    if (!node[k].getMoveX() || (dimension > 1 && !node[k].getMoveY()) || (dimension == 3 && !node[k].getMoveZ())) {
                        nodesOnBoundary++;
                    }
                }
                if (nodesOnBoundary == 0) {
                    for (int k = 0; k < node.length; k++) {
                        elementCoord.add(new cbit.vcell.render.Vect3d(node[k].getX(), node[k].getY(), node[k].getZ()));
                    }
                    elementCoord.scale(0.25);
                } else if (nodesOnBoundary == 2) {
                    for (int k = 0; k < node.length; k++) {
                        if (!node[k].getMoveX() || !node[k].getMoveY() || !node[k].getMoveZ()) {
                            elementCoord.add(new cbit.vcell.render.Vect3d(node[k].getX(), node[k].getY(), node[k].getZ()));
                        }
                    }
                    elementCoord.scale(0.5);
                } else if (nodesOnBoundary == 3) {
                    for (int k = 0; k < node.length; k++) {
                        if (!node[k].getMoveX() && !node[k].getMoveY() || !node[k].getMoveY() && !node[k].getMoveZ() || !node[k].getMoveX() && !node[k].getMoveZ()) {
                            elementCoord.set(node[k].getX(), node[k].getY(), node[k].getZ());
                        }
                    }
                } else {
                    throw new RuntimeException("Unexcepted number of nodes on boundary for a polygon: " + nodesOnBoundary);
                }
                cbit.vcell.render.Vect3d unitNormal = new cbit.vcell.render.Vect3d();
                polygon.getUnitNormal(unitNormal);
                int volNeighbor1Region = regionImage.getRegionInfoFromOffset(polygon.getVolIndexNeighbor1()).getRegionIndex();
                int volNeighbor2Region = regionImage.getRegionInfoFromOffset(polygon.getVolIndexNeighbor2()).getRegionIndex();
                if (surface.getExteriorRegionIndex() == volNeighbor1Region && surface.getInteriorRegionIndex() == volNeighbor2Region) {
                    region1Outside++;
                }
                if (surface.getExteriorRegionIndex() == volNeighbor2Region && surface.getInteriorRegionIndex() == volNeighbor1Region) {
                    region1Inside++;
                }
                writer.write(cellID + " " + polygon.getVolIndexNeighbor1() + " " + polygon.getVolIndexNeighbor2() + " " + polygon.getArea() / correctCoeff + " " + elementCoord.getX() + " " + elementCoord.getY() + " " + elementCoord.getZ() + " " + unitNormal.getX() + " " + unitNormal.getY() + " " + unitNormal.getZ() + "\n");
                cellID++;
            }
            if (region1Inside != surface.getPolygonCount() && region1Outside != surface.getPolygonCount()) {
                throw new RuntimeException("Volume neighbor regions not consistent: [total, inside, outside]=" + surface.getPolygonCount() + "," + region1Inside + "," + region1Outside + "]");
            }
        }
    }
}

Example 19 with VolumeGeometricRegion

use of cbit.vcell.geometry.surface.VolumeGeometricRegion in project vcell by virtualcell.

the class SBMLExporter method addGeometry.

private void addGeometry() throws SbmlException {
    SpatialModelPlugin mplugin = (SpatialModelPlugin) sbmlModel.getPlugin(SBMLUtils.SBML_SPATIAL_NS_PREFIX);
    // Creates a geometry object via SpatialModelPlugin object.
    org.sbml.jsbml.ext.spatial.Geometry sbmlGeometry = mplugin.createGeometry();
    sbmlGeometry.setCoordinateSystem(GeometryKind.cartesian);
    sbmlGeometry.setSpatialId("vcell");
    Geometry vcGeometry = getSelectedSimContext().getGeometry();
    Model vcModel = getSelectedSimContext().getModel();
    // 
    // list of CoordinateComponents : 1 if geometry is 1-d, 2 if geometry is 2-d, 3 if geometry is 3-d
    // 
    int dimension = vcGeometry.getDimension();
    Extent vcExtent = vcGeometry.getExtent();
    Origin vcOrigin = vcGeometry.getOrigin();
    // add x coordinate component
    CoordinateComponent xComp = sbmlGeometry.createCoordinateComponent();
    xComp.setSpatialId(vcModel.getX().getName());
    xComp.setType(CoordinateKind.cartesianX);
    final UnitDefinition sbmlUnitDef_length = getOrCreateSBMLUnit(vcModel.getUnitSystem().getLengthUnit());
    xComp.setUnits(sbmlUnitDef_length);
    Boundary minX = new Boundary();
    xComp.setBoundaryMinimum(minX);
    minX.setSpatialId("Xmin");
    minX.setValue(vcOrigin.getX());
    Boundary maxX = new Boundary();
    xComp.setBoundaryMaximum(maxX);
    maxX.setSpatialId("Xmax");
    maxX.setValue(vcOrigin.getX() + (vcExtent.getX()));
    org.sbml.jsbml.Parameter pX = sbmlModel.createParameter();
    pX.setId(vcModel.getX().getName());
    pX.setValue(0.0);
    pX.setConstant(false);
    pX.setUnits(sbmlUnitDef_length);
    SpatialParameterPlugin spPluginPx = (SpatialParameterPlugin) pX.getPlugin(SBMLUtils.SBML_SPATIAL_NS_PREFIX);
    SpatialSymbolReference spSymRefPx = new SpatialSymbolReference();
    spPluginPx.setParamType(spSymRefPx);
    spSymRefPx.setSpatialRef(xComp.getSpatialId());
    // add y coordinate component
    if (dimension == 2 || dimension == 3) {
        CoordinateComponent yComp = sbmlGeometry.createCoordinateComponent();
        yComp.setSpatialId(vcModel.getY().getName());
        yComp.setType(CoordinateKind.cartesianY);
        yComp.setUnits(sbmlUnitDef_length);
        Boundary minY = new Boundary();
        yComp.setBoundaryMinimum(minY);
        minY.setSpatialId("Ymin");
        minY.setValue(vcOrigin.getY());
        Boundary maxY = new Boundary();
        yComp.setBoundaryMaximum(maxY);
        maxY.setSpatialId("Ymax");
        maxY.setValue(vcOrigin.getY() + (vcExtent.getY()));
        org.sbml.jsbml.Parameter pY = sbmlModel.createParameter();
        pY.setId(vcModel.getY().getName());
        pY.setValue(0.0);
        pY.setConstant(false);
        pY.setUnits(sbmlUnitDef_length);
        SpatialParameterPlugin spPluginPy = (SpatialParameterPlugin) pY.getPlugin(SBMLUtils.SBML_SPATIAL_NS_PREFIX);
        SpatialSymbolReference spSymRefPy = new SpatialSymbolReference();
        spPluginPy.setParamType(spSymRefPy);
        spSymRefPy.setSpatialRef(yComp.getSpatialId());
    }
    // add z coordinate component
    if (dimension == 3) {
        CoordinateComponent zComp = sbmlGeometry.createCoordinateComponent();
        zComp.setSpatialId(vcModel.getZ().getName());
        zComp.setType(CoordinateKind.cartesianZ);
        zComp.setUnits(sbmlUnitDef_length);
        Boundary minZ = new Boundary();
        zComp.setBoundaryMinimum(minZ);
        minZ.setSpatialId("Zmin");
        minZ.setValue(vcOrigin.getZ());
        Boundary maxZ = new Boundary();
        zComp.setBoundaryMaximum(maxZ);
        maxZ.setSpatialId("Zmax");
        maxZ.setValue(vcOrigin.getZ() + (vcExtent.getZ()));
        org.sbml.jsbml.Parameter pZ = sbmlModel.createParameter();
        pZ.setId(vcModel.getZ().getName());
        pZ.setValue(0.0);
        pZ.setConstant(false);
        pZ.setUnits(sbmlUnitDef_length);
        SpatialParameterPlugin spPluginPz = (SpatialParameterPlugin) pZ.getPlugin(SBMLUtils.SBML_SPATIAL_NS_PREFIX);
        SpatialSymbolReference spSymRefPz = new SpatialSymbolReference();
        spPluginPz.setParamType(spSymRefPz);
        spSymRefPz.setSpatialRef(zComp.getSpatialId());
    }
    // 
    // list of compartmentMappings : VC structureMappings
    // 
    GeometryContext vcGeoContext = getSelectedSimContext().getGeometryContext();
    StructureMapping[] vcStrucMappings = vcGeoContext.getStructureMappings();
    for (int i = 0; i < vcStrucMappings.length; i++) {
        StructureMapping vcStructMapping = vcStrucMappings[i];
        String structName = vcStructMapping.getStructure().getName();
        Compartment comp = sbmlModel.getCompartment(TokenMangler.mangleToSName(structName));
        SpatialCompartmentPlugin cplugin = (SpatialCompartmentPlugin) comp.getPlugin(SBMLUtils.SBML_SPATIAL_NS_PREFIX);
        GeometryClass gc = vcStructMapping.getGeometryClass();
        if (!goodPointer(gc, GeometryClass.class, structName)) {
            continue;
        }
        CompartmentMapping compMapping = new CompartmentMapping();
        cplugin.setCompartmentMapping(compMapping);
        String geomClassName = gc.getName();
        String id = TokenMangler.mangleToSName(geomClassName + structName);
        compMapping.setSpatialId(id);
        compMapping.setDomainType(TokenMangler.mangleToSName(DOMAIN_TYPE_PREFIX + geomClassName));
        try {
            StructureMappingParameter usp = vcStructMapping.getUnitSizeParameter();
            Expression e = usp.getExpression();
            if (goodPointer(e, Expression.class, id)) {
                compMapping.setUnitSize(e.evaluateConstant());
            }
        } catch (ExpressionException e) {
            e.printStackTrace(System.out);
            throw new RuntimeException("Unable to create compartment mapping for structureMapping '" + compMapping.getId() + "' : " + e.getMessage());
        }
    }
    // 
    // list of domain types : subvolumes and surface classes from VC
    // 
    boolean bAnyAnalyticSubvolumes = false;
    boolean bAnyImageSubvolumes = false;
    boolean bAnyCSGSubvolumes = false;
    GeometryClass[] vcGeomClasses = vcGeometry.getGeometryClasses();
    int numSubVols = 0;
    for (int i = 0; i < vcGeomClasses.length; i++) {
        DomainType domainType = sbmlGeometry.createDomainType();
        domainType.setSpatialId(DOMAIN_TYPE_PREFIX + vcGeomClasses[i].getName());
        if (vcGeomClasses[i] instanceof SubVolume) {
            if (((SubVolume) vcGeomClasses[i]) instanceof AnalyticSubVolume) {
                bAnyAnalyticSubvolumes = true;
            } else if (((SubVolume) vcGeomClasses[i]) instanceof ImageSubVolume) {
                bAnyImageSubvolumes = true;
            } else if (((SubVolume) vcGeomClasses[i]) instanceof CSGObject) {
                bAnyCSGSubvolumes = true;
            }
            domainType.setSpatialDimensions(3);
            numSubVols++;
        } else if (vcGeomClasses[i] instanceof SurfaceClass) {
            domainType.setSpatialDimensions(2);
        }
    }
    // 
    // list of domains, adjacent domains : from VC geometricRegions
    // 
    GeometrySurfaceDescription vcGSD = vcGeometry.getGeometrySurfaceDescription();
    if (vcGSD.getRegionImage() == null) {
        try {
            vcGSD.updateAll();
        } catch (Exception e) {
            e.printStackTrace(System.out);
            throw new RuntimeException("Unable to generate region images for geometry");
        }
    }
    GeometricRegion[] vcGeometricRegions = vcGSD.getGeometricRegions();
    ISize sampleSize = vcGSD.getVolumeSampleSize();
    int numX = sampleSize.getX();
    int numY = sampleSize.getY();
    int numZ = sampleSize.getZ();
    double ox = vcOrigin.getX();
    double oy = vcOrigin.getY();
    double oz = vcOrigin.getZ();
    RegionInfo[] regionInfos = vcGSD.getRegionImage().getRegionInfos();
    for (int i = 0; i < vcGeometricRegions.length; i++) {
        // domains
        Domain domain = sbmlGeometry.createDomain();
        domain.setSpatialId(vcGeometricRegions[i].getName());
        if (vcGeometricRegions[i] instanceof VolumeGeometricRegion) {
            domain.setDomainType(DOMAIN_TYPE_PREFIX + ((VolumeGeometricRegion) vcGeometricRegions[i]).getSubVolume().getName());
            // 
            // get a list of interior points ... should probably use the distance map to find a point
            // furthest inside (or several points associated with the morphological skeleton).
            // 
            InteriorPoint interiorPt = domain.createInteriorPoint();
            int regionID = ((VolumeGeometricRegion) vcGeometricRegions[i]).getRegionID();
            boolean bFound = false;
            int regInfoIndx = 0;
            for (int j = 0; j < regionInfos.length; j++) {
                regInfoIndx = j;
                if (regionInfos[j].getRegionIndex() == regionID) {
                    int volIndx = 0;
                    for (int z = 0; z < numZ && !bFound; z++) {
                        for (int y = 0; y < numY && !bFound; y++) {
                            for (int x = 0; x < numX && !bFound; x++) {
                                if (regionInfos[j].isIndexInRegion(volIndx)) {
                                    bFound = true;
                                    double unit_z = (numZ > 1) ? ((double) z) / (numZ - 1) : 0.5;
                                    double coordZ = oz + vcExtent.getZ() * unit_z;
                                    double unit_y = (numY > 1) ? ((double) y) / (numY - 1) : 0.5;
                                    double coordY = oy + vcExtent.getY() * unit_y;
                                    double unit_x = (numX > 1) ? ((double) x) / (numX - 1) : 0.5;
                                    double coordX = ox + vcExtent.getX() * unit_x;
                                    interiorPt.setCoord1(coordX);
                                    interiorPt.setCoord2(coordY);
                                    interiorPt.setCoord3(coordZ);
                                }
                                volIndx++;
                            }
                        // end - for x
                        }
                    // end - for y
                    }
                // end - for z
                }
            // end if
            }
            // end for regionInfos
            if (!bFound) {
                throw new RuntimeException("Unable to find interior point for region '" + regionInfos[regInfoIndx].toString());
            }
        } else if (vcGeometricRegions[i] instanceof SurfaceGeometricRegion) {
            SurfaceGeometricRegion vcSurfaceGeomReg = (SurfaceGeometricRegion) vcGeometricRegions[i];
            GeometricRegion geomRegion0 = vcSurfaceGeomReg.getAdjacentGeometricRegions()[0];
            GeometricRegion geomRegion1 = vcSurfaceGeomReg.getAdjacentGeometricRegions()[1];
            SurfaceClass surfaceClass = vcGSD.getSurfaceClass(((VolumeGeometricRegion) geomRegion0).getSubVolume(), ((VolumeGeometricRegion) geomRegion1).getSubVolume());
            domain.setDomainType(DOMAIN_TYPE_PREFIX + surfaceClass.getName());
            // adjacent domains : 2 adjacent domain objects for each surfaceClass in VC.
            // adjacent domain 1
            GeometricRegion adjGeomRegion0 = vcSurfaceGeomReg.getAdjacentGeometricRegions()[0];
            GeometricRegion adjGeomRegion1 = vcSurfaceGeomReg.getAdjacentGeometricRegions()[1];
            AdjacentDomains adjDomain = new AdjacentDomains();
            adjDomain.setSpatialId(TokenMangler.mangleToSName(vcSurfaceGeomReg.getName() + "_" + adjGeomRegion0.getName()));
            adjDomain.setDomain1(vcSurfaceGeomReg.getName());
            adjDomain.setDomain2(adjGeomRegion0.getName());
            sbmlGeometry.addAdjacentDomain(adjDomain);
            // adj domain 2
            adjDomain = new AdjacentDomains();
            adjDomain.setSpatialId(TokenMangler.mangleToSName(vcSurfaceGeomReg.getName() + "_" + adjGeomRegion1.getName()));
            adjDomain.setDomain1(vcSurfaceGeomReg.getName());
            adjDomain.setDomain2(adjGeomRegion1.getName());
            sbmlGeometry.addAdjacentDomain(adjDomain);
        }
    }
    // 
    if (bAnyAnalyticSubvolumes && !bAnyImageSubvolumes && !bAnyCSGSubvolumes) {
        AnalyticGeometry sbmlAnalyticGeomDefinition = sbmlGeometry.createAnalyticGeometry();
        sbmlAnalyticGeomDefinition.setSpatialId(TokenMangler.mangleToSName("Analytic_" + vcGeometry.getName()));
        sbmlAnalyticGeomDefinition.setIsActive(true);
        for (int i = 0; i < vcGeomClasses.length; i++) {
            if (vcGeomClasses[i] instanceof AnalyticSubVolume) {
                AnalyticVolume analyticVol = sbmlAnalyticGeomDefinition.createAnalyticVolume();
                analyticVol.setSpatialId(vcGeomClasses[i].getName());
                analyticVol.setDomainType(DOMAIN_TYPE_PREFIX + vcGeomClasses[i].getName());
                analyticVol.setFunctionType(FunctionKind.layered);
                analyticVol.setOrdinal(numSubVols - (i + 1));
                Expression expr = ((AnalyticSubVolume) vcGeomClasses[i]).getExpression();
                try {
                    String mathMLStr = ExpressionMathMLPrinter.getMathML(expr, true, MathType.BOOLEAN);
                    ASTNode mathMLNode = ASTNode.readMathMLFromString(mathMLStr);
                    analyticVol.setMath(mathMLNode);
                } catch (Exception e) {
                    e.printStackTrace(System.out);
                    throw new RuntimeException("Error converting VC subvolume expression to mathML" + e.getMessage());
                }
            }
        }
    }
    // 
    if (!bAnyAnalyticSubvolumes && !bAnyImageSubvolumes && bAnyCSGSubvolumes) {
        CSGeometry sbmlCSGeomDefinition = new CSGeometry();
        sbmlGeometry.addGeometryDefinition(sbmlCSGeomDefinition);
        sbmlCSGeomDefinition.setSpatialId(TokenMangler.mangleToSName("CSG_" + vcGeometry.getName()));
        for (int i = 0; i < vcGeomClasses.length; i++) {
            if (vcGeomClasses[i] instanceof CSGObject) {
                CSGObject vcellCSGObject = (CSGObject) vcGeomClasses[i];
                org.sbml.jsbml.ext.spatial.CSGObject sbmlCSGObject = new org.sbml.jsbml.ext.spatial.CSGObject();
                sbmlCSGeomDefinition.addCSGObject(sbmlCSGObject);
                sbmlCSGObject.setSpatialId(vcellCSGObject.getName());
                sbmlCSGObject.setDomainType(DOMAIN_TYPE_PREFIX + vcellCSGObject.getName());
                // the ordinal should the the least for the default/background subVolume
                sbmlCSGObject.setOrdinal(numSubVols - (i + 1));
                org.sbml.jsbml.ext.spatial.CSGNode sbmlcsgNode = getSBMLCSGNode(vcellCSGObject.getRoot());
                sbmlCSGObject.setCSGNode(sbmlcsgNode);
            }
        }
    }
    // 
    // add "Segmented" and "DistanceMap" SampledField Geometries
    // 
    final boolean bVCGeometryIsImage = bAnyImageSubvolumes && !bAnyAnalyticSubvolumes && !bAnyCSGSubvolumes;
    // 55if (bAnyAnalyticSubvolumes || bAnyImageSubvolumes || bAnyCSGSubvolumes){
    if (bVCGeometryIsImage) {
        // 
        // add "Segmented" SampledFieldGeometry
        // 
        SampledFieldGeometry segmentedImageSampledFieldGeometry = sbmlGeometry.createSampledFieldGeometry();
        segmentedImageSampledFieldGeometry.setSpatialId(TokenMangler.mangleToSName("SegmentedImage_" + vcGeometry.getName()));
        segmentedImageSampledFieldGeometry.setIsActive(true);
        // 55boolean bVCGeometryIsImage = bAnyImageSubvolumes && !bAnyAnalyticSubvolumes && !bAnyCSGSubvolumes;
        Geometry vcImageGeometry = null;
        {
            if (bVCGeometryIsImage) {
                // make a resampled image;
                if (dimension == 3) {
                    try {
                        ISize imageSize = vcGeometry.getGeometrySpec().getDefaultSampledImageSize();
                        vcGeometry.precomputeAll(new GeometryThumbnailImageFactoryAWT());
                        vcImageGeometry = RayCaster.resampleGeometry(new GeometryThumbnailImageFactoryAWT(), vcGeometry, imageSize);
                    } catch (Throwable e) {
                        e.printStackTrace(System.out);
                        throw new RuntimeException("Unable to convert the original analytic or constructed solid geometry to image-based geometry : " + e.getMessage());
                    }
                } else {
                    try {
                        vcGeometry.precomputeAll(new GeometryThumbnailImageFactoryAWT(), true, false);
                        GeometrySpec origGeometrySpec = vcGeometry.getGeometrySpec();
                        VCImage newVCImage = origGeometrySpec.getSampledImage().getCurrentValue();
                        // 
                        // construct the new geometry with the sampled VCImage.
                        // 
                        vcImageGeometry = new Geometry(vcGeometry.getName() + "_asImage", newVCImage);
                        vcImageGeometry.getGeometrySpec().setExtent(vcGeometry.getExtent());
                        vcImageGeometry.getGeometrySpec().setOrigin(vcGeometry.getOrigin());
                        vcImageGeometry.setDescription(vcGeometry.getDescription());
                        vcImageGeometry.getGeometrySurfaceDescription().setFilterCutoffFrequency(vcGeometry.getGeometrySurfaceDescription().getFilterCutoffFrequency());
                        vcImageGeometry.precomputeAll(new GeometryThumbnailImageFactoryAWT(), true, true);
                    } catch (Exception e) {
                        e.printStackTrace(System.out);
                        throw new RuntimeException("Unable to convert the original analytic or constructed solid geometry to image-based geometry : " + e.getMessage());
                    }
                }
                GeometryClass[] vcImageGeomClasses = vcImageGeometry.getGeometryClasses();
                for (int j = 0; j < vcImageGeomClasses.length; j++) {
                    if (vcImageGeomClasses[j] instanceof ImageSubVolume) {
                        SampledVolume sampledVol = segmentedImageSampledFieldGeometry.createSampledVolume();
                        sampledVol.setSpatialId(vcGeomClasses[j].getName());
                        sampledVol.setDomainType(DOMAIN_TYPE_PREFIX + vcGeomClasses[j].getName());
                        sampledVol.setSampledValue(((ImageSubVolume) vcImageGeomClasses[j]).getPixelValue());
                    }
                }
                // add sampledField to sampledFieldGeometry
                SampledField segmentedImageSampledField = sbmlGeometry.createSampledField();
                VCImage vcImage = vcImageGeometry.getGeometrySpec().getImage();
                segmentedImageSampledField.setSpatialId("SegmentedImageSampledField");
                segmentedImageSampledField.setNumSamples1(vcImage.getNumX());
                segmentedImageSampledField.setNumSamples2(vcImage.getNumY());
                segmentedImageSampledField.setNumSamples3(vcImage.getNumZ());
                segmentedImageSampledField.setInterpolationType(InterpolationKind.nearestneighbor);
                segmentedImageSampledField.setCompression(CompressionKind.uncompressed);
                segmentedImageSampledField.setDataType(DataKind.UINT8);
                segmentedImageSampledFieldGeometry.setSampledField(segmentedImageSampledField.getId());
                try {
                    byte[] vcImagePixelsBytes = vcImage.getPixels();
                    // imageData.setCompression("");
                    StringBuffer sb = new StringBuffer();
                    for (int i = 0; i < vcImagePixelsBytes.length; i++) {
                        int uint8_sample = ((int) vcImagePixelsBytes[i]) & 0xff;
                        sb.append(uint8_sample + " ");
                    }
                    segmentedImageSampledField.setSamplesLength(vcImage.getNumXYZ());
                    segmentedImageSampledField.setSamples(sb.toString().trim());
                } catch (ImageException e) {
                    e.printStackTrace(System.out);
                    throw new RuntimeException("Unable to export image from VCell to SBML : " + e.getMessage());
                }
            }
        }
    /*		
		//
		// add "DistanceMap" SampledFieldGeometry if there are exactly two subvolumes (else need more fields) and geometry is 3d.
		//
		if (numSubVols==2 && dimension == 3){
			SignedDistanceMap[] distanceMaps = null;
			try {
				distanceMaps = DistanceMapGenerator.computeDistanceMaps(vcImageGeometry, vcImageGeometry.getGeometrySpec().getImage(), false, false);
			} catch (ImageException e) {
				e.printStackTrace(System.out);
				System.err.println("Unable to export distance map sampled field from VCell to SBML : " + e.getMessage());
				// throw new RuntimeException("Unable to export distance map sampled field from VCell to SBML : " + e.getMessage());
				
				// don't want to throw an exception and stop export because distance map geometry couldn't be exported. 
				// just 'return' from method (since this is the last thing that is being done in this method).
				return;
			}
			//
			// the two distanceMaps should be redundant (one is negation of the other) ... so choose first one for field.
			//
			double[] signedDistances = distanceMaps[0].getSignedDistances();
			SampledFieldGeometry distanceMapSampledFieldGeometry = sbmlGeometry.createSampledFieldGeometry();
			distanceMapSampledFieldGeometry.setSpatialId(TokenMangler.mangleToSName("DistanceMap_"+vcGeometry.getName()));
			SampledField distanceMapSampledField = distanceMapSampledFieldGeometry.createSampledField();
			distanceMapSampledField.setSpatialId("DistanceMapSampledField");
			distanceMapSampledField.setNumSamples1(distanceMaps[0].getSamplesX().length);
			distanceMapSampledField.setNumSamples2(distanceMaps[0].getSamplesY().length);
			distanceMapSampledField.setNumSamples3(distanceMaps[0].getSamplesZ().length);
			distanceMapSampledField.setDataType("real");
System.err.println("do we need distanceMapSampleField.setDataType()?");
			distanceMapSampledField.setInterpolationType("linear");
			ImageData distanceMapImageData = distanceMapSampledField.createImageData();
			distanceMapImageData.setDataType("int16");
System.err.println("should be:\n  distanceMapImageData.setDataType(\"float32\")");
//					distanceMapImageData.setCompression("");

			double maxAbsValue = 0;
			for (int i = 0; i < signedDistances.length; i++) {
				maxAbsValue = Math.max(maxAbsValue,Math.abs(signedDistances[i]));
			}
			if (maxAbsValue==0.0){
				throw new RuntimeException("computed distance map all zeros");
			}
			double scale = (Short.MAX_VALUE-1)/maxAbsValue;
			int[] scaledIntegerDistanceMap = new int[signedDistances.length];
			for (int i = 0; i < signedDistances.length; i++) {
				scaledIntegerDistanceMap[i] = (int)(scale * signedDistances[i]);
			}
			distanceMapImageData.setSamples(scaledIntegerDistanceMap, signedDistances.length);
System.err.println("should be:\n  distanceMapImageData.setSamples((float[])signedDistances,signedDistances.length)");
			SampledVolume sampledVol = distanceMapSampledFieldGeometry.createSampledVolume();
			sampledVol.setSpatialId(distanceMaps[0].getInsideSubvolumeName());
			sampledVol.setDomainType(DOMAIN_TYPE_PREFIX+distanceMaps[0].getInsideSubvolumeName());
			sampledVol.setSampledValue(255);
			sampledVol = distanceMapSampledFieldGeometry.createSampledVolume();
			sampledVol.setSpatialId(distanceMaps[1].getInsideSubvolumeName());
			sampledVol.setDomainType(DOMAIN_TYPE_PREFIX+distanceMaps[1].getInsideSubvolumeName());
			sampledVol.setSampledValue(1);
		}
*/
    }
// 
// add "SurfaceMesh" ParametricGeometry
// 
// if (bAnyAnalyticSubvolumes || bAnyImageSubvolumes || bAnyCSGSubvolumes){
// ParametricGeometry sbmlParametricGeomDefinition = sbmlGeometry.createParametricGeometry();
// sbmlParametricGeomDefinition.setSpatialId(TokenMangler.mangleToSName("SurfaceMesh_"+vcGeometry.getName()));
// xxxx
// }
}

Example 20 with VolumeGeometricRegion

use of cbit.vcell.geometry.surface.VolumeGeometricRegion in project vcell by virtualcell.

the class SBMLImporter method getAssociatedSurfaceGeometricRegion.

private SurfaceGeometricRegion getAssociatedSurfaceGeometricRegion(GeometrySurfaceDescription vcGsd, Vector<VolumeGeometricRegion> volGeomRegionsVector) {
    GeometricRegion[] geomeRegions = vcGsd.getGeometricRegions();
    // adjVolGeomRegionsVector should have only 2 elements - the 2 adj
    // volGeomRegions for any surfaceRegion.
    VolumeGeometricRegion[] volGeomRegionsArray = volGeomRegionsVector.toArray(new VolumeGeometricRegion[0]);
    for (int i = 0; i < geomeRegions.length; i++) {
        if (geomeRegions[i] instanceof SurfaceGeometricRegion) {
            SurfaceGeometricRegion surfaceRegion = (SurfaceGeometricRegion) geomeRegions[i];
            GeometricRegion[] adjVolGeomRegs = surfaceRegion.getAdjacentGeometricRegions();
            // if the 2 arrays do not have 2 elements each, throw exception
            if (volGeomRegionsArray.length != 2 && adjVolGeomRegs.length != 2) {
                throw new SBMLImportException("There should be 2 adjacent geometric regions for surfaceRegion '" + surfaceRegion.getName() + "'");
            }
            // winner! - return surfaceRegion
            if ((adjVolGeomRegs[0].compareEqual(volGeomRegionsArray[0]) && adjVolGeomRegs[1].compareEqual(volGeomRegionsArray[1])) || (adjVolGeomRegs[0].compareEqual(volGeomRegionsArray[1]) && adjVolGeomRegs[1].compareEqual(volGeomRegionsArray[0]))) {
                return surfaceRegion;
            }
        }
    }
    return null;
}