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

Example 1 with SurfaceCollection

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

the class PDEDataViewer method createDataValueSurfaceViewer.

/**
 * Insert the method's description here.
 * Creation date: (9/25/2005 1:53:00 PM)
 */
private DataValueSurfaceViewer createDataValueSurfaceViewer(ClientTaskStatusSupport clientTaskStatusSupport) throws ImageException, UserCancelException {
    // try{
    // if(fieldDataValueSurfaceViewer == null){
    // Surfaces
    CartesianMesh cartesianMesh = getPdeDataContext().getCartesianMesh();
    if (cartesianMesh.getMembraneElements() == null || cartesianMesh.getMembraneElements().length == 0 || cartesianMesh.isChomboMesh()) {
        // return fieldDataValueSurfaceViewer;
        return getDataValueSurfaceViewer();
    }
    meshRegionSurfaces = new MeshDisplayAdapter(cartesianMesh).generateMeshRegionSurfaces(clientTaskStatusSupport);
    SurfaceCollection surfaceCollection = meshRegionSurfaces.getSurfaceCollection();
    // SurfaceNames
    final String[] surfaceNames = new String[meshRegionSurfaces.getSurfaceCollection().getSurfaceCount()];
    for (int i = 0; i < meshRegionSurfaces.getSurfaceCollection().getSurfaceCount(); i++) {
        // Get the first element, any will do, all have same inside/outside volumeIndex
        MembraneElement me = cartesianMesh.getMembraneElements()[meshRegionSurfaces.getMembraneIndexForPolygon(i, 0)];
        if (getSimulationModelInfo() != null) {
            surfaceNames[i] = getSimulationModelInfo().getMembraneName(cartesianMesh.getSubVolumeFromVolumeIndex(me.getInsideVolumeIndex()), cartesianMesh.getSubVolumeFromVolumeIndex(me.getOutsideVolumeIndex()), false);
        } else {
            surfaceNames[i] = i + "";
        }
    }
    // SurfaceAreas
    final Double[] surfaceAreas = new Double[meshRegionSurfaces.getSurfaceCollection().getSurfaceCount()];
    for (int i = 0; i < meshRegionSurfaces.getSurfaceCollection().getSurfaceCount(); i++) {
        surfaceAreas[i] = new Double(cartesianMesh.getRegionMembraneSurfaceAreaFromMembraneIndex(meshRegionSurfaces.getMembraneIndexForPolygon(i, 0)));
    }
    // DataValueSurfaceViewer fieldDataValueSurfaceViewer0 = new DataValueSurfaceViewer();
    TaubinSmoothing taubinSmoothing = new TaubinSmoothingWrong();
    TaubinSmoothingSpecification taubinSpec = TaubinSmoothingSpecification.getInstance(.3);
    taubinSmoothing.smooth(surfaceCollection, taubinSpec);
    getDataValueSurfaceViewer().init(meshRegionSurfaces.getSurfaceCollection(), cartesianMesh.getOrigin(), cartesianMesh.getExtent(), surfaceNames, surfaceAreas, cartesianMesh.getGeometryDimension());
    return getDataValueSurfaceViewer();
// }
// }catch(UserCancelException e){
// throw e;
// }catch(Exception e){
// PopupGenerator.showErrorDialog(PDEDataViewer.this, e.getMessage(), e);
// }
// return fieldDataValueSurfaceViewer;
}

Example 2 with SurfaceCollection

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

the class PDEDataViewer method updateDataValueSurfaceViewer0.

// private AsynchClientTask[] getDataVlaueSurfaceViewerTasks(){
// AsynchClientTask createDataValueSurfaceViewerTask = new AsynchClientTask("Create surface viewer...",AsynchClientTask.TASKTYPE_NONSWING_BLOCKING) {
// @Override
// public void run(Hashtable<String, Object> hashTable) throws Exception {
// if(getDataValueSurfaceViewer().getSurfaceCollectionDataInfoProvider() == null){
// createDataValueSurfaceViewer(getClientTaskStatusSupport());
// }
// }
// };
// 
// AsynchClientTask updateDataValueSurfaceViewerTask = new AsynchClientTask("Update surface viewer...",AsynchClientTask.TASKTYPE_SWING_BLOCKING) {
// @Override
// public void run(Hashtable<String, Object> hashTable) throws Exception {
// updateDataValueSurfaceViewer0();
// }
// };
// 
// AsynchClientTask resetDataValueSurfaceViewerTask = new AsynchClientTask("Reset tab...",AsynchClientTask.TASKTYPE_SWING_NONBLOCKING,false,false) {
// @Override
// public void run(Hashtable<String, Object> hashTable) throws Exception {
// if(getDataValueSurfaceViewer().getSurfaceCollectionDataInfoProvider() == null){
// viewDataTabbedPane.setSelectedIndex(0);
// }
// }
// };
// return new AsynchClientTask[] {createDataValueSurfaceViewerTask,updateDataValueSurfaceViewerTask,resetDataValueSurfaceViewerTask};
// }
// private Timer dataValueSurfaceTimer;
// //private boolean bPdeIsParamScan=false;
// private void updateDataValueSurfaceViewer(){
// //	if((dataValueSurfaceTimer = ClientTaskDispatcher.getBlockingTimer(this,getPdeDataContext(),null,dataValueSurfaceTimer,new ActionListener() {@Override public void actionPerformed(ActionEvent e2) {updateDataValueSurfaceViewer();}}))!=null){
// //		return;
// //	}
// if(bSkipSurfaceCalc){
// return;
// }
// if(getDataValueSurfaceViewer().getSurfaceCollectionDataInfoProvider() == null){
// if((dataValueSurfaceTimer = ClientTaskDispatcher.getBlockingTimer(this,getPdeDataContext(),null,dataValueSurfaceTimer,new ActionListener() {@Override public void actionPerformed(ActionEvent e2) {updateDataValueSurfaceViewer();}}))!=null){
// return;
// }
// ClientTaskDispatcher.dispatch(this, new Hashtable<String, Object>(), getDataVlaueSurfaceViewerTasks(),true,true,null);
// }else{
// try{
// updateDataValueSurfaceViewer0();
// }catch(Exception e){
// e.printStackTrace();
// DialogUtils.showErrorDialog(this, e.getMessage());
// }
// }
// }
/**
 * Insert the method's description here.
 * Creation date: (9/25/2005 2:00:05 PM)
 */
private void updateDataValueSurfaceViewer0() {
    // viewDataTabbedPane.addTab(CurrentView.SURFACE_VIEW.title, getDataValueSurfaceViewer());
    if (viewDataTabbedPane.getSelectedIndex() != CurrentView.SURFACE_VIEW.ordinal()) {
        return;
    }
    // SurfaceColors and DataValues
    if (getDataValueSurfaceViewer().getSurfaceCollectionDataInfo() == null) {
        // happens with PostProcessingImageData version of PDEDataViewer
        return;
    }
    SurfaceCollection surfaceCollection = getDataValueSurfaceViewer().getSurfaceCollectionDataInfo().getSurfaceCollection();
    DisplayAdapterService das = getPDEDataContextPanel1().getdisplayAdapterService1();
    final int[][] surfaceColors = new int[surfaceCollection.getSurfaceCount()][];
    final double[][] surfaceDataValues = new double[surfaceCollection.getSurfaceCount()][];
    boolean bMembraneVariable = getPdeDataContext().getDataIdentifier().getVariableType().equals(VariableType.MEMBRANE);
    RecodeDataForDomainInfo recodeDataForDomainInfo = getPDEDataContextPanel1().getRecodeDataForDomainInfo();
    double[] membraneValues = (recodeDataForDomainInfo.isRecoded() ? recodeDataForDomainInfo.getRecodedDataForDomain() : getPdeDataContext().getDataValues());
    for (int i = 0; i < surfaceCollection.getSurfaceCount(); i += 1) {
        Surface surface = surfaceCollection.getSurfaces(i);
        surfaceColors[i] = new int[surface.getPolygonCount()];
        surfaceDataValues[i] = new double[surface.getPolygonCount()];
        for (int j = 0; j < surface.getPolygonCount(); j += 1) {
            int membraneIndexForPolygon = meshRegionSurfaces.getMembraneIndexForPolygon(i, j);
            if (bMembraneVariable) {
                surfaceDataValues[i][j] = membraneValues[membraneIndexForPolygon];
            } else {
                // get membrane region index from membrane index
                surfaceDataValues[i][j] = membraneValues[getPdeDataContext().getCartesianMesh().getMembraneRegionIndex(membraneIndexForPolygon)];
            }
            surfaceColors[i][j] = das.getColorFromValue(surfaceDataValues[i][j]);
        }
    }
    DataValueSurfaceViewer.SurfaceCollectionDataInfoProvider svdp = new DataValueSurfaceViewer.SurfaceCollectionDataInfoProvider() {

        private DisplayAdapterService updatedDAS = new DisplayAdapterService(getPDEDataContextPanel1().getdisplayAdapterService1());

        private String updatedVariableName = getPdeDataContext().getVariableName();

        private double updatedTimePoint = getPdeDataContext().getTimePoint();

        private double[] updatedVariableValues = getPdeDataContext().getDataValues();

        private VCDataIdentifier updatedVCDataIdentifier = getPdeDataContext().getVCDataIdentifier();

        public void makeMovie(SurfaceCanvas surfaceCanvas) {
            makeSurfaceMovie(surfaceCanvas, updatedVariableValues.length, updatedVariableName, updatedDAS, updatedVCDataIdentifier);
        }

        public double getValue(int surfaceIndex, int polygonIndex) {
            return updatedVariableValues[meshRegionSurfaces.getMembraneIndexForPolygon(surfaceIndex, polygonIndex)];
        }

        public String getValueDescription(int surfaceIndex, int polygonIndex) {
            return updatedVariableName;
        }

        public int[][] getSurfacePolygonColors() {
            return surfaceColors;
        }

        public Coordinate getCentroid(int surfaceIndex, int polygonIndex) {
            return getPdeDataContext().getCartesianMesh().getMembraneElements()[meshRegionSurfaces.getMembraneIndexForPolygon(surfaceIndex, polygonIndex)].getCentroid();
        }

        public float getArea(int surfaceIndex, int polygonIndex) {
            return getPdeDataContext().getCartesianMesh().getMembraneElements()[meshRegionSurfaces.getMembraneIndexForPolygon(surfaceIndex, polygonIndex)].getArea();
        }

        public Vect3d getNormal(int surfaceIndex, int polygonIndex) {
            return getPdeDataContext().getCartesianMesh().getMembraneElements()[meshRegionSurfaces.getMembraneIndexForPolygon(surfaceIndex, polygonIndex)].getNormal();
        }

        public int getMembraneIndex(int surfaceIndex, int polygonIndex) {
            return meshRegionSurfaces.getMembraneIndexForPolygon(surfaceIndex, polygonIndex);
        }

        public Color getROIHighlightColor() {
            return new Color(getPDEDataContextPanel1().getdisplayAdapterService1().getSpecialColors()[cbit.image.DisplayAdapterService.FOREGROUND_HIGHLIGHT_COLOR_OFFSET]);
        }

        @Override
        public boolean isMembrIndexInVarDomain(int membrIndex) {
            return (getPDEDataContextPanel1().getRecodeDataForDomainInfo() != null ? getPDEDataContextPanel1().getRecodeDataForDomainInfo().isIndexInDomain(membrIndex) : true);
        }

        // public void showComponentInFrame(Component comp,String title){
        // PDEDataViewer.this.showComponentInFrame(comp,title);
        // }
        public void plotTimeSeriesData(int[][] indices, boolean bAllTimes, boolean bTimeStats, boolean bSpaceStats) throws DataAccessException {
            double[] timePoints = getPdeDataContext().getTimePoints();
            double beginTime = (bAllTimes ? timePoints[0] : updatedTimePoint);
            double endTime = (bAllTimes ? timePoints[timePoints.length - 1] : beginTime);
            String[] varNames = new String[indices.length];
            for (int i = 0; i < varNames.length; i += 1) {
                varNames[i] = updatedVariableName;
            }
            TimeSeriesJobSpec timeSeriesJobSpec = new TimeSeriesJobSpec(varNames, indices, beginTime, 1, endTime, bSpaceStats, bTimeStats, VCDataJobID.createVCDataJobID(getDataViewerManager().getUser(), true));
            Hashtable<String, Object> hash = new Hashtable<String, Object>();
            hash.put(StringKey_timeSeriesJobSpec, timeSeriesJobSpec);
            AsynchClientTask task1 = new TimeSeriesDataRetrievalTask("Retrieve data", PDEDataViewer.this, getPdeDataContext());
            AsynchClientTask task2 = new AsynchClientTask("Showing surface", AsynchClientTask.TASKTYPE_SWING_BLOCKING) {

                @Override
                public void run(Hashtable<String, Object> hashTable) throws Exception {
                    TSJobResultsSpaceStats tsJobResultsSpaceStats = (TSJobResultsSpaceStats) hashTable.get(StringKey_timeSeriesJobResults);
                    plotSpaceStats(tsJobResultsSpaceStats);
                }
            };
            ClientTaskDispatcher.dispatch(PDEDataViewer.this, hash, new AsynchClientTask[] { task1, task2 }, true, true, null);
        }
    };
    getDataValueSurfaceViewer().setSurfaceCollectionDataInfoProvider(svdp);
}

Example 3 with SurfaceCollection

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

the class SmoldynFileWriter method writeSurfaces.

private void writeSurfaces() throws SolverException, ImageException, PropertyVetoException, GeometryException, ExpressionException {
    GeometrySurfaceDescription geometrySurfaceDescription = resampledGeometry.getGeometrySurfaceDescription();
    SurfaceClass[] surfaceClasses = geometrySurfaceDescription.getSurfaceClasses();
    GeometrySpec geometrySpec = resampledGeometry.getGeometrySpec();
    SubVolume[] surfaceGeometrySubVolumes = geometrySpec.getSubVolumes();
    GeometricRegion[] AllGeometricRegions = resampledGeometry.getGeometrySurfaceDescription().getGeometricRegions();
    ArrayList<SurfaceGeometricRegion> surfaceRegionList = new ArrayList<SurfaceGeometricRegion>();
    ArrayList<VolumeGeometricRegion> volumeRegionList = new ArrayList<VolumeGeometricRegion>();
    for (GeometricRegion geometricRegion : AllGeometricRegions) {
        if (geometricRegion instanceof SurfaceGeometricRegion) {
            surfaceRegionList.add((SurfaceGeometricRegion) geometricRegion);
        } else if (geometricRegion instanceof VolumeGeometricRegion) {
            volumeRegionList.add((VolumeGeometricRegion) geometricRegion);
        } else {
            throw new SolverException("unsupported geometric region type " + geometricRegion.getClass());
        }
    }
    printWriter.println("# geometry");
    printWriter.println(SmoldynVCellMapper.SmoldynKeyword.dim + " " + dimension);
    if (bHasNoSurface) {
        printWriter.println(SmoldynVCellMapper.SmoldynKeyword.max_compartment + " " + surfaceGeometrySubVolumes.length);
    } else {
        printWriter.println(SmoldynVCellMapper.SmoldynKeyword.max_compartment + " " + (surfaceGeometrySubVolumes.length + 1));
        // plus the surface which are bounding walls
        printWriter.println(SmoldynVCellMapper.SmoldynKeyword.max_surface + " " + (surfaceClasses.length + dimension));
    }
    printWriter.println();
    // write boundaries and wall surfaces
    writeWallSurfaces();
    // for 3D ... smoldyn normal convension is triangle right-hand-rule normal points to the outside compartment subdomain.
    if (!bHasNoSurface) {
        membraneSubdomainTriangleMap = new HashMap<MembraneSubDomain, ArrayList<TrianglePanel>>();
        // write surfaces
        printWriter.println("# surfaces");
        int triangleGlobalCount = 0;
        int membraneIndex = -1;
        SurfaceCollection surfaceCollection = geometrySurfaceDescription.getSurfaceCollection();
        // pre-allocate collections used repeatedly in following loops; clear before reusing
        HashMap<Node, Set<String>> nodeTriMap = new HashMap<>();
        ArrayList<TrianglePanel> triList = new ArrayList<TrianglePanel>();
        // use a sorted set to ensure neighbors written out is same order for reproducibility
        SortedSet<String> neighborsForCurrentNode = new TreeSet<String>();
        for (int sci = 0; sci < surfaceClasses.length; sci++) {
            nodeTriMap.clear();
            triList.clear();
            int triLocalCount = 0;
            SurfaceClass surfaceClass = surfaceClasses[sci];
            GeometricRegion[] geometricRegions = geometrySurfaceDescription.getGeometricRegions(surfaceClass);
            for (GeometricRegion gr : geometricRegions) {
                SurfaceGeometricRegion sgr = (SurfaceGeometricRegion) gr;
                VolumeGeometricRegion volRegion0 = (VolumeGeometricRegion) sgr.getAdjacentGeometricRegions()[0];
                VolumeGeometricRegion volRegion1 = (VolumeGeometricRegion) sgr.getAdjacentGeometricRegions()[1];
                SubVolume subVolume0 = volRegion0.getSubVolume();
                SubVolume subVolume1 = volRegion1.getSubVolume();
                CompartmentSubDomain compart0 = mathDesc.getCompartmentSubDomain(subVolume0.getName());
                CompartmentSubDomain compart1 = mathDesc.getCompartmentSubDomain(subVolume1.getName());
                MembraneSubDomain membraneSubDomain = mathDesc.getMembraneSubDomain(compart0, compart1);
                if (membraneSubDomain == null) {
                    throw new SolverException(VCellErrorMessages.getSmoldynUnexpectedSurface(compart0, compart1));
                }
                int exteriorRegionID = volRegion0.getRegionID();
                int interiorRegionID = volRegion1.getRegionID();
                if (membraneSubDomain.getInsideCompartment() == compart0) {
                    exteriorRegionID = volRegion1.getRegionID();
                    interiorRegionID = volRegion0.getRegionID();
                }
                for (int j = 0; j < surfaceCollection.getSurfaceCount(); j++) {
                    Surface surface = surfaceCollection.getSurfaces(j);
                    if ((surface.getInteriorRegionIndex() == exteriorRegionID && surface.getExteriorRegionIndex() == interiorRegionID) || (surface.getInteriorRegionIndex() == interiorRegionID && surface.getExteriorRegionIndex() == exteriorRegionID)) {
                        // Polygon polygon = surface.getPolygons(k);
                        for (Polygon polygon : surface) {
                            if (polygonMembaneElementMap != null) {
                                membraneIndex = polygonMembaneElementMap.get(polygon).getMembraneIndex();
                            }
                            Node[] nodes = polygon.getNodes();
                            if (dimension == 2) {
                                // ignore z
                                Vect3d unitNormal = new Vect3d();
                                polygon.getUnitNormal(unitNormal);
                                unitNormal.set(unitNormal.getX(), unitNormal.getY(), 0);
                                int point0 = 0;
                                Vect3d v0 = new Vect3d(nodes[point0].getX(), nodes[point0].getY(), 0);
                                int point1 = 1;
                                Vect3d v1 = null;
                                for (point1 = 1; point1 < nodes.length; point1++) {
                                    if (v0.getX() != nodes[point1].getX() || v0.getY() != nodes[point1].getY()) {
                                        v1 = new Vect3d(nodes[point1].getX(), nodes[point1].getY(), 0);
                                        break;
                                    }
                                }
                                if (v1 == null) {
                                    throw new RuntimeException("failed to generate surface");
                                }
                                Vect3d v01 = Vect3d.sub(v1, v0);
                                Vect3d unit01n = v01.cross(unitNormal);
                                unit01n.unit();
                                if (Math.abs(unit01n.getZ() - 1.0) < 1e-6) {
                                    // v0 to v1 opposes vcell surface normal. it's already flipped.
                                    Triangle triangle;
                                    if (surface.getInteriorRegionIndex() == interiorRegionID) {
                                        // we have to flipped it back
                                        triangle = new Triangle(nodes[point1], nodes[point0], null);
                                    } else {
                                        triangle = new Triangle(nodes[point0], nodes[point1], null);
                                    }
                                    triList.add(new TrianglePanel(triLocalCount++, triangleGlobalCount++, membraneIndex, triangle));
                                } else if (Math.abs(unit01n.getZ() + 1.0) < 1e-6) {
                                    // v0 to v1 is in direction of vcell surface normal.
                                    Triangle triangle;
                                    if (surface.getInteriorRegionIndex() == interiorRegionID) {
                                        triangle = new Triangle(nodes[point0], nodes[point1], null);
                                    } else {
                                        triangle = new Triangle(nodes[point1], nodes[point0], null);
                                    }
                                    triList.add(new TrianglePanel(triLocalCount++, triangleGlobalCount++, membraneIndex, triangle));
                                } else {
                                    throw new RuntimeException("failed to generate surface");
                                }
                            } else if (dimension == 3) {
                                Triangle triangle1;
                                Triangle triangle2;
                                if (surface.getInteriorRegionIndex() == interiorRegionID) {
                                    // interior
                                    triangle1 = new Triangle(nodes[0], nodes[1], nodes[2]);
                                    triangle2 = new Triangle(nodes[0], nodes[2], nodes[3]);
                                } else {
                                    triangle1 = new Triangle(nodes[2], nodes[1], nodes[0]);
                                    triangle2 = new Triangle(nodes[3], nodes[2], nodes[0]);
                                }
                                triList.add(new TrianglePanel(triLocalCount++, triangleGlobalCount++, membraneIndex, triangle1));
                                triList.add(new TrianglePanel(triLocalCount++, triangleGlobalCount++, membraneIndex, triangle2));
                            }
                        }
                    }
                }
            }
            // add triangles to node hash
            for (TrianglePanel triPanel : triList) {
                for (Node node : triPanel.triangle.getNodes()) {
                    if (node == null) {
                        continue;
                    }
                    Set<String> triNameSet = nodeTriMap.get(node);
                    if (triNameSet == null) {
                        triNameSet = new HashSet<String>();
                        nodeTriMap.put(node, triNameSet);
                    }
                    triNameSet.add(triPanel.name);
                }
            }
            SubVolume[] adjacentSubvolums = surfaceClass.getAdjacentSubvolumes().toArray(new SubVolume[0]);
            CompartmentSubDomain csd0 = simulation.getMathDescription().getCompartmentSubDomain(adjacentSubvolums[0].getName());
            CompartmentSubDomain csd1 = simulation.getMathDescription().getCompartmentSubDomain(adjacentSubvolums[1].getName());
            MembraneSubDomain membraneSubDomain = simulation.getMathDescription().getMembraneSubDomain(csd0, csd1);
            membraneSubdomainTriangleMap.put(membraneSubDomain, triList);
            final boolean initialMoleculesOnMembrane = (closestTriangles != null);
            if (initialMoleculesOnMembrane) {
                findClosestTriangles(membraneSubDomain, triList);
            }
            printWriter.println(SmoldynVCellMapper.SmoldynKeyword.start_surface + " " + surfaceClass.getName());
            printWriter.println(SmoldynVCellMapper.SmoldynKeyword.action + " " + SmoldynVCellMapper.SmoldynKeyword.all + "(" + SmoldynVCellMapper.SmoldynKeyword.all + ") " + SmoldynVCellMapper.SmoldynKeyword.both + " " + SmoldynVCellMapper.SmoldynKeyword.reflect);
            // printWriter.println(SmoldynKeyword.action + " " + SmoldynKeyword.all + "(" + SmoldynKeyword.up + ") " + SmoldynKeyword.both + " " + SmoldynKeyword.reflect);
            // get color after species
            Color c = colors[sci + particleVariableList.size()];
            printWriter.println(SmoldynVCellMapper.SmoldynKeyword.color + " " + SmoldynVCellMapper.SmoldynKeyword.both + " " + c.getRed() / 255.0 + " " + c.getGreen() / 255.0 + " " + c.getBlue() / 255.0 + " 0.1");
            printWriter.println(SmoldynVCellMapper.SmoldynKeyword.polygon + " " + SmoldynVCellMapper.SmoldynKeyword.front + " " + SmoldynVCellMapper.SmoldynKeyword.edge);
            printWriter.println(SmoldynVCellMapper.SmoldynKeyword.polygon + " " + SmoldynVCellMapper.SmoldynKeyword.back + " " + SmoldynVCellMapper.SmoldynKeyword.edge);
            printWriter.println(SmoldynVCellMapper.SmoldynKeyword.max_panels + " " + SmoldynVCellMapper.SmoldynKeyword.tri + " " + triList.size());
            for (TrianglePanel trianglePanel : triList) {
                Triangle triangle = trianglePanel.triangle;
                printWriter.print(SmoldynVCellMapper.SmoldynKeyword.panel + " " + SmoldynVCellMapper.SmoldynKeyword.tri);
                switch(dimension) {
                    case 1:
                        printWriter.print(" " + triangle.getNodes(0).getX());
                        break;
                    case 2:
                        printWriter.print(" " + triangle.getNodes(0).getX() + " " + triangle.getNodes(0).getY());
                        printWriter.print(" " + triangle.getNodes(1).getX() + " " + triangle.getNodes(1).getY());
                        break;
                    case 3:
                        for (Node node : triangle.getNodes()) {
                            printWriter.print(" " + node.getX() + " " + node.getY() + " " + node.getZ());
                        }
                        break;
                }
                printWriter.println(" " + trianglePanel.name);
            }
            for (TrianglePanel triPanel : triList) {
                neighborsForCurrentNode.clear();
                for (Node node : triPanel.triangle.getNodes()) {
                    if (node == null) {
                        continue;
                    }
                    neighborsForCurrentNode.addAll(nodeTriMap.get(node));
                }
                neighborsForCurrentNode.remove(triPanel.name);
                // printWriter.print(SmoldynKeyword.neighbors + " " +triPanel.name);
                // to allow smoldyn read line length as 256, chop the neighbors to multiple lines
                int maxNeighborCount = 4;
                // 
                int count = 0;
                for (String neigh : neighborsForCurrentNode) {
                    if (count % maxNeighborCount == 0) {
                        printWriter.println();
                        printWriter.print(SmoldynVCellMapper.SmoldynKeyword.neighbors + " " + triPanel.name);
                    }
                    printWriter.print(" " + neigh);
                    count++;
                }
            }
            printWriter.println();
            printWriter.println(SmoldynVCellMapper.SmoldynKeyword.end_surface);
            printWriter.println();
        }
    // write compartment
    // printWriter.println("# bounding wall compartment");
    // printWriter.println(SmoldynKeyword.start_compartment + " " + VCellSmoldynKeyword.bounding_wall_compartment);
    // printWriter.println(SmoldynKeyword.surface + " " + VCellSmoldynKeyword.bounding_wall_surface_X);
    // if (dimension > 1) {
    // printWriter.println(SmoldynKeyword.surface + " " + VCellSmoldynKeyword.bounding_wall_surface_Y);
    // if (dimension > 2) {
    // printWriter.println(SmoldynKeyword.surface + " " + VCellSmoldynKeyword.bounding_wall_surface_Z);
    // }
    // }
    // printWriter.println(SmoldynKeyword.end_compartment);
    // printWriter.println();
    }
}

Example 4 with SurfaceCollection

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

the class SmoldynSurfaceTessellator method writeSurfaces.

protected void writeSurfaces() throws SolverException, ImageException, PropertyVetoException, GeometryException, ExpressionException {
    GeometrySurfaceDescription geometrySurfaceDescription = resampledGeometry.getGeometrySurfaceDescription();
    SurfaceClass[] surfaceClasses = geometrySurfaceDescription.getSurfaceClasses();
    GeometrySpec geometrySpec = resampledGeometry.getGeometrySpec();
    SubVolume[] surfaceGeometrySubVolumes = geometrySpec.getSubVolumes();
    GeometricRegion[] AllGeometricRegions = resampledGeometry.getGeometrySurfaceDescription().getGeometricRegions();
    ArrayList<SurfaceGeometricRegion> surfaceRegionList = new ArrayList<SurfaceGeometricRegion>();
    ArrayList<VolumeGeometricRegion> volumeRegionList = new ArrayList<VolumeGeometricRegion>();
    for (GeometricRegion geometricRegion : AllGeometricRegions) {
        if (geometricRegion instanceof SurfaceGeometricRegion) {
            surfaceRegionList.add((SurfaceGeometricRegion) geometricRegion);
        } else if (geometricRegion instanceof VolumeGeometricRegion) {
            volumeRegionList.add((VolumeGeometricRegion) geometricRegion);
        } else {
            throw new SolverException("unsupported geometric region type " + geometricRegion.getClass());
        }
    }
    printWriter.println("# geometry");
    printWriter.println(SmoldynVCellMapper.SmoldynKeyword.dim + " " + dimension);
    if (bHasNoSurface) {
        printWriter.println(SmoldynVCellMapper.SmoldynKeyword.max_compartment + " " + surfaceGeometrySubVolumes.length);
    } else {
        printWriter.println(SmoldynVCellMapper.SmoldynKeyword.max_compartment + " " + (surfaceGeometrySubVolumes.length + 1));
        // plus the surface which are bounding walls
        printWriter.println(SmoldynVCellMapper.SmoldynKeyword.max_surface + " " + (surfaceClasses.length + dimension));
    }
    printWriter.println();
    // write boundaries and wall surfaces
    writeWallSurfaces();
    // for 3D ... smoldyn normal convension is triangle right-hand-rule normal points to the outside compartment subdomain.
    if (!bHasNoSurface) {
        membraneSubdomainTriangleMap = new HashMap<MembraneSubDomain, ArrayList<TrianglePanel>>();
        // write surfaces
        printWriter.println("# surfaces");
        int triangleGlobalCount = 0;
        int membraneIndex = -1;
        SurfaceCollection surfaceCollection = geometrySurfaceDescription.getSurfaceCollection();
        // pre-allocate collections used repeatedly in following loops; clear before reusing
        HashMap<Node, Set<String>> nodeTriMap = new HashMap<>();
        ArrayList<TrianglePanel> triList = new ArrayList<TrianglePanel>();
        // use a sorted set to ensure neighbors written out is same order for reproducibility
        SortedSet<String> neighborsForCurrentNode = new TreeSet<String>();
        for (int sci = 0; sci < surfaceClasses.length; sci++) {
            nodeTriMap.clear();
            triList.clear();
            int triLocalCount = 0;
            SurfaceClass surfaceClass = surfaceClasses[sci];
            GeometricRegion[] geometricRegions = geometrySurfaceDescription.getGeometricRegions(surfaceClass);
            for (GeometricRegion gr : geometricRegions) {
                SurfaceGeometricRegion sgr = (SurfaceGeometricRegion) gr;
                VolumeGeometricRegion volRegion0 = (VolumeGeometricRegion) sgr.getAdjacentGeometricRegions()[0];
                VolumeGeometricRegion volRegion1 = (VolumeGeometricRegion) sgr.getAdjacentGeometricRegions()[1];
                SubVolume subVolume0 = volRegion0.getSubVolume();
                SubVolume subVolume1 = volRegion1.getSubVolume();
                CompartmentSubDomain compart0 = mathDesc.getCompartmentSubDomain(subVolume0.getName());
                CompartmentSubDomain compart1 = mathDesc.getCompartmentSubDomain(subVolume1.getName());
                MembraneSubDomain membraneSubDomain = mathDesc.getMembraneSubDomain(compart0, compart1);
                if (membraneSubDomain == null) {
                    throw new SolverException(VCellErrorMessages.getSmoldynUnexpectedSurface(compart0, compart1));
                }
                int exteriorRegionID = volRegion0.getRegionID();
                int interiorRegionID = volRegion1.getRegionID();
                if (membraneSubDomain.getInsideCompartment() == compart0) {
                    exteriorRegionID = volRegion1.getRegionID();
                    interiorRegionID = volRegion0.getRegionID();
                }
                for (int j = 0; j < surfaceCollection.getSurfaceCount(); j++) {
                    Surface surface = surfaceCollection.getSurfaces(j);
                    if ((surface.getInteriorRegionIndex() == exteriorRegionID && surface.getExteriorRegionIndex() == interiorRegionID) || (surface.getInteriorRegionIndex() == interiorRegionID && surface.getExteriorRegionIndex() == exteriorRegionID)) {
                        // Polygon polygon = surface.getPolygons(k);
                        for (Polygon polygon : surface) {
                            if (polygonMembaneElementMap != null) {
                                membraneIndex = polygonMembaneElementMap.get(polygon).getMembraneIndex();
                            }
                            Node[] nodes = polygon.getNodes();
                            if (dimension == 2) {
                                // ignore z
                                Vect3d unitNormal = new Vect3d();
                                polygon.getUnitNormal(unitNormal);
                                unitNormal.set(unitNormal.getX(), unitNormal.getY(), 0);
                                int point0 = 0;
                                Vect3d v0 = new Vect3d(nodes[point0].getX(), nodes[point0].getY(), 0);
                                int point1 = 1;
                                Vect3d v1 = null;
                                for (point1 = 1; point1 < nodes.length; point1++) {
                                    if (v0.getX() != nodes[point1].getX() || v0.getY() != nodes[point1].getY()) {
                                        v1 = new Vect3d(nodes[point1].getX(), nodes[point1].getY(), 0);
                                        break;
                                    }
                                }
                                if (v1 == null) {
                                    throw new RuntimeException("failed to generate surface");
                                }
                                Vect3d v01 = Vect3d.sub(v1, v0);
                                Vect3d unit01n = v01.cross(unitNormal);
                                unit01n.unit();
                                if (Math.abs(unit01n.getZ() - 1.0) < 1e-6) {
                                    // v0 to v1 opposes vcell surface normal. it's already flipped.
                                    Triangle triangle;
                                    if (surface.getInteriorRegionIndex() == interiorRegionID) {
                                        // we have to flipped it back
                                        triangle = new Triangle(nodes[point1], nodes[point0], null);
                                    } else {
                                        triangle = new Triangle(nodes[point0], nodes[point1], null);
                                    }
                                    triList.add(new TrianglePanel(triLocalCount++, triangleGlobalCount++, membraneIndex, triangle));
                                } else if (Math.abs(unit01n.getZ() + 1.0) < 1e-6) {
                                    // v0 to v1 is in direction of vcell surface normal.
                                    Triangle triangle;
                                    if (surface.getInteriorRegionIndex() == interiorRegionID) {
                                        triangle = new Triangle(nodes[point0], nodes[point1], null);
                                    } else {
                                        triangle = new Triangle(nodes[point1], nodes[point0], null);
                                    }
                                    triList.add(new TrianglePanel(triLocalCount++, triangleGlobalCount++, membraneIndex, triangle));
                                } else {
                                    throw new RuntimeException("failed to generate surface");
                                }
                            } else if (dimension == 3) {
                                Triangle triangle1;
                                Triangle triangle2;
                                if (surface.getInteriorRegionIndex() == interiorRegionID) {
                                    // interior
                                    triangle1 = new Triangle(nodes[0], nodes[1], nodes[2]);
                                    triangle2 = new Triangle(nodes[0], nodes[2], nodes[3]);
                                } else {
                                    triangle1 = new Triangle(nodes[2], nodes[1], nodes[0]);
                                    triangle2 = new Triangle(nodes[3], nodes[2], nodes[0]);
                                }
                                triList.add(new TrianglePanel(triLocalCount++, triangleGlobalCount++, membraneIndex, triangle1));
                                triList.add(new TrianglePanel(triLocalCount++, triangleGlobalCount++, membraneIndex, triangle2));
                            }
                        }
                    }
                }
            }
            // add triangles to node hash
            for (TrianglePanel triPanel : triList) {
                for (Node node : triPanel.triangle.getNodes()) {
                    if (node == null) {
                        continue;
                    }
                    Set<String> triNameSet = nodeTriMap.get(node);
                    if (triNameSet == null) {
                        triNameSet = new HashSet<String>();
                        nodeTriMap.put(node, triNameSet);
                    }
                    triNameSet.add(triPanel.name);
                }
            }
            SubVolume[] adjacentSubvolums = surfaceClass.getAdjacentSubvolumes().toArray(new SubVolume[0]);
            CompartmentSubDomain csd0 = simulation.getMathDescription().getCompartmentSubDomain(adjacentSubvolums[0].getName());
            CompartmentSubDomain csd1 = simulation.getMathDescription().getCompartmentSubDomain(adjacentSubvolums[1].getName());
            MembraneSubDomain membraneSubDomain = simulation.getMathDescription().getMembraneSubDomain(csd0, csd1);
            membraneSubdomainTriangleMap.put(membraneSubDomain, triList);
            final boolean initialMoleculesOnMembrane = (closestTriangles != null);
            if (initialMoleculesOnMembrane) {
                findClosestTriangles(membraneSubDomain, triList);
            }
            printWriter.println(SmoldynVCellMapper.SmoldynKeyword.start_surface + " " + surfaceClass.getName());
            printWriter.println(SmoldynVCellMapper.SmoldynKeyword.action + " " + SmoldynVCellMapper.SmoldynKeyword.all + "(" + SmoldynVCellMapper.SmoldynKeyword.all + ") " + SmoldynVCellMapper.SmoldynKeyword.both + " " + SmoldynVCellMapper.SmoldynKeyword.reflect);
            // printWriter.println(SmoldynKeyword.action + " " + SmoldynKeyword.all + "(" + SmoldynKeyword.up + ") " + SmoldynKeyword.both + " " + SmoldynKeyword.reflect);
            Color c = colorForSurface(sci);
            printWriter.println(SmoldynVCellMapper.SmoldynKeyword.color + " " + SmoldynVCellMapper.SmoldynKeyword.both + " " + c.getRed() / 255.0 + " " + c.getGreen() / 255.0 + " " + c.getBlue() / 255.0 + " 0.1");
            printWriter.println(SmoldynVCellMapper.SmoldynKeyword.polygon + " " + SmoldynVCellMapper.SmoldynKeyword.front + " " + SmoldynVCellMapper.SmoldynKeyword.edge);
            printWriter.println(SmoldynVCellMapper.SmoldynKeyword.polygon + " " + SmoldynVCellMapper.SmoldynKeyword.back + " " + SmoldynVCellMapper.SmoldynKeyword.edge);
            printWriter.println(SmoldynVCellMapper.SmoldynKeyword.max_panels + " " + SmoldynVCellMapper.SmoldynKeyword.tri + " " + triList.size());
            for (TrianglePanel trianglePanel : triList) {
                Triangle triangle = trianglePanel.triangle;
                printWriter.print(SmoldynVCellMapper.SmoldynKeyword.panel + " " + SmoldynVCellMapper.SmoldynKeyword.tri);
                switch(dimension) {
                    case 1:
                        printWriter.print(" " + triangle.getNodes(0).getX());
                        break;
                    case 2:
                        printWriter.print(" " + triangle.getNodes(0).getX() + " " + triangle.getNodes(0).getY());
                        printWriter.print(" " + triangle.getNodes(1).getX() + " " + triangle.getNodes(1).getY());
                        break;
                    case 3:
                        for (Node node : triangle.getNodes()) {
                            printWriter.print(" " + node.getX() + " " + node.getY() + " " + node.getZ());
                        }
                        break;
                }
                printWriter.println(" " + trianglePanel.name);
            }
            for (TrianglePanel triPanel : triList) {
                neighborsForCurrentNode.clear();
                for (Node node : triPanel.triangle.getNodes()) {
                    if (node == null) {
                        continue;
                    }
                    neighborsForCurrentNode.addAll(nodeTriMap.get(node));
                }
                neighborsForCurrentNode.remove(triPanel.name);
                // printWriter.print(SmoldynKeyword.neighbors + " " +triPanel.name);
                // to allow smoldyn read line length as 256, chop the neighbors to multiple lines
                int maxNeighborCount = 4;
                // 
                int count = 0;
                for (String neigh : neighborsForCurrentNode) {
                    if (count % maxNeighborCount == 0) {
                        printWriter.println();
                        printWriter.print(SmoldynVCellMapper.SmoldynKeyword.neighbors + " " + triPanel.name);
                    }
                    printWriter.print(" " + neigh);
                    count++;
                }
            }
            printWriter.println();
            printWriter.println(SmoldynVCellMapper.SmoldynKeyword.end_surface);
            printWriter.println();
        }
    }
}

Example 5 with SurfaceCollection

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

the class RegionImage method generateSurfaceCollection.

private void generateSurfaceCollection(int numRegions, // int[] mapImageIndexToLinkRegion,int[] mapLinkRegionToDistinctRegion,
VCImage vcImage, BitSet xSurfElements, BitSet ySurfElements, BitSet zSurfElements, int dimension, Extent extent, Origin origin) {
    bMembraneNeighborCalculationFailed = false;
    int masterIndex = 0;
    double dX = extent.getX() / (vcImage.getNumX() - 1);
    double dY = extent.getY() / (vcImage.getNumY() - 1);
    double dZ = extent.getZ() / (vcImage.getNumZ() - 1);
    // 
    double loX = origin.getX() - 0.5 * dX;
    double loY = origin.getY() - 0.5 * dY;
    double loZ = origin.getZ() - 0.5 * dZ;
    Vector<Vector<Quadrilateral>> surfQuadsV = new Vector<Vector<Quadrilateral>>();
    int[][] mapInsideOutsideToSurface = new int[numRegions][];
    Vector<cbit.vcell.geometry.surface.Node> nodeListV = new Vector<cbit.vcell.geometry.surface.Node>();
    cbit.vcell.geometry.surface.Node[][][] mapImageIndexToNode = new cbit.vcell.geometry.surface.Node[2][vcImage.getNumY() + 1][vcImage.getNumX() + 1];
    double zComp, zpComp;
    double yComp, ypComp;
    double xComp, xpComp;
    final int yStep = vcImage.getNumX();
    final int zStep = vcImage.getNumX() * vcImage.getNumY();
    boolean bMvXm, bMvYm, bMvZm;
    boolean bMvXp, bMvYp, bMvZp;
    double xmEdge = origin.getX();
    double xpEdge = origin.getX() + extent.getX();
    double ymEdge = origin.getY();
    double ypEdge = origin.getY() + extent.getY();
    double zmEdge = origin.getZ();
    double zpEdge = origin.getZ() + extent.getZ();
    for (int zIndex = 0; zIndex < vcImage.getNumZ(); zIndex++) {
        mapImageIndexToNode[0] = mapImageIndexToNode[1];
        mapImageIndexToNode[1] = new cbit.vcell.geometry.surface.Node[vcImage.getNumY() + 1][vcImage.getNumX() + 1];
        zComp = loZ + zIndex * dZ;
        zpComp = zComp + dZ;
        bMvZm = (zIndex != 0);
        bMvZp = (zIndex != (vcImage.getNumZ() - 1));
        for (int yIndex = 0; yIndex < vcImage.getNumY(); yIndex++) {
            yComp = loY + yIndex * dY;
            ypComp = yComp + dY;
            bMvYm = (yIndex != 0);
            bMvYp = (yIndex != (vcImage.getNumY() - 1));
            for (int xIndex = 0; xIndex < vcImage.getNumX(); xIndex++) {
                if (xSurfElements.get(masterIndex)) {
                    bMvXp = (xIndex != (vcImage.getNumX() - 1));
                    xpComp = loX + (xIndex + 1) * dX;
                    cbit.vcell.geometry.surface.Node[] nodeArr = new cbit.vcell.geometry.surface.Node[4];
                    nodeArr[0] = mapImageIndexToNode[0][yIndex][xIndex + 1];
                    if (nodeArr[0] == null) {
                        // nodeArr[0] = new cbit.vcell.geometry.surface.Node(xpComp,yComp,zComp);
                        nodeArr[0] = new cbit.vcell.geometry.surface.Node((bMvXp ? xpComp : xpEdge), (bMvYm ? yComp : ymEdge), (bMvZm ? zComp : zmEdge));
                        nodeArr[0].setGlobalIndex(nodeListV.size());
                        nodeListV.add(nodeArr[0]);
                        nodeArr[0].setMoveX(bMvXp);
                        nodeArr[0].setMoveY(bMvYm);
                        nodeArr[0].setMoveZ(bMvZm);
                        mapImageIndexToNode[0][yIndex][xIndex + 1] = nodeArr[0];
                    }
                    nodeArr[1] = mapImageIndexToNode[0][yIndex + 1][xIndex + 1];
                    if (nodeArr[1] == null) {
                        // nodeArr[1] = new cbit.vcell.geometry.surface.Node(xpComp,ypComp,zComp);
                        nodeArr[1] = new cbit.vcell.geometry.surface.Node((bMvXp ? xpComp : xpEdge), (bMvYp ? ypComp : ypEdge), (bMvZm ? zComp : zmEdge));
                        nodeArr[1].setGlobalIndex(nodeListV.size());
                        nodeListV.add(nodeArr[1]);
                        nodeArr[1].setMoveX(bMvXp);
                        nodeArr[1].setMoveY(bMvYp);
                        nodeArr[1].setMoveZ(bMvZm);
                        mapImageIndexToNode[0][yIndex + 1][xIndex + 1] = nodeArr[1];
                    }
                    nodeArr[2] = mapImageIndexToNode[1][yIndex + 1][xIndex + 1];
                    if (nodeArr[2] == null) {
                        // nodeArr[2] = new cbit.vcell.geometry.surface.Node(xpComp,ypComp,zpComp);
                        nodeArr[2] = new cbit.vcell.geometry.surface.Node((bMvXp ? xpComp : xpEdge), (bMvYp ? ypComp : ypEdge), (bMvZp ? zpComp : zpEdge));
                        nodeArr[2].setGlobalIndex(nodeListV.size());
                        nodeListV.add(nodeArr[2]);
                        nodeArr[2].setMoveX(bMvXp);
                        nodeArr[2].setMoveY(bMvYp);
                        nodeArr[2].setMoveZ(bMvZp);
                        mapImageIndexToNode[1][yIndex + 1][xIndex + 1] = nodeArr[2];
                    }
                    nodeArr[3] = mapImageIndexToNode[1][yIndex][xIndex + 1];
                    if (nodeArr[3] == null) {
                        // nodeArr[3] = new cbit.vcell.geometry.surface.Node(xpComp,yComp,zpComp);
                        nodeArr[3] = new cbit.vcell.geometry.surface.Node((bMvXp ? xpComp : xpEdge), (bMvYm ? yComp : ymEdge), (bMvZp ? zpComp : zpEdge));
                        nodeArr[3].setGlobalIndex(nodeListV.size());
                        nodeListV.add(nodeArr[3]);
                        nodeArr[3].setMoveX(bMvXp);
                        nodeArr[3].setMoveY(bMvYm);
                        nodeArr[3].setMoveZ(bMvZp);
                        mapImageIndexToNode[1][yIndex][xIndex + 1] = nodeArr[3];
                    }
                    addQuadToSurface(surfQuadsV, mapInsideOutsideToSurface, // masterIndex+1,mapImageIndexToRegionIndex.getValue(masterIndex+1),
                    masterIndex, mapLinkRegionToDistinctRegion[mapImageIndexToLinkRegion[masterIndex]], masterIndex + 1, mapLinkRegionToDistinctRegion[mapImageIndexToLinkRegion[masterIndex + 1]], nodeArr, numRegions, MembraneElementIdentifier.PerpendicularTo.X);
                // surfQuad =
                // new Quadrilateral(nodeArr,
                // mapLinkRegionToDistinctRegion[mapImageIndexToLinkRegion[masterIndex]],
                // mapLinkRegionToDistinctRegion[mapImageIndexToLinkRegion[masterIndex+1]]);
                // if(surfQuad.getVolIndexNeighbor1() > surfQuad.getVolIndexNeighbor2()){
                // surfQuad.reverseDirection();
                // }
                // if(mapInsideOutsideToSurface[surfQuad.getVolIndexNeighbor1()] == null){
                // mapInsideOutsideToSurface[surfQuad.getVolIndexNeighbor1()] = new int[numRegions];
                // Arrays.fill(mapInsideOutsideToSurface[surfQuad.getVolIndexNeighbor1()], -1);
                // }
                // if(mapInsideOutsideToSurface[surfQuad.getVolIndexNeighbor1()][surfQuad.getVolIndexNeighbor2()] == -1){
                // mapInsideOutsideToSurface[surfQuad.getVolIndexNeighbor1()][surfQuad.getVolIndexNeighbor2()] = surfQuadsV.size();
                // surfQuadsV.add(new Vector<Quadrilateral>());
                // }
                // surfQuadsV.elementAt(mapInsideOutsideToSurface[surfQuad.getVolIndexNeighbor1()][surfQuad.getVolIndexNeighbor2()]).add(surfQuad);
                }
                if (ySurfElements.get(masterIndex)) {
                    bMvXm = (xIndex != 0);
                    bMvXp = (xIndex != (vcImage.getNumX() - 1));
                    xComp = loX + xIndex * dX;
                    xpComp = loX + (xIndex + 1) * dX;
                    cbit.vcell.geometry.surface.Node[] nodeArr = new cbit.vcell.geometry.surface.Node[4];
                    nodeArr[0] = mapImageIndexToNode[0][yIndex + 1][xIndex + 1];
                    if (nodeArr[0] == null) {
                        // nodeArr[0] = new cbit.vcell.geometry.surface.Node(xpComp,ypComp,zComp);
                        nodeArr[0] = new cbit.vcell.geometry.surface.Node((bMvXp ? xpComp : xpEdge), (bMvYp ? ypComp : ypEdge), (bMvZm ? zComp : zmEdge));
                        nodeArr[0].setGlobalIndex(nodeListV.size());
                        nodeListV.add(nodeArr[0]);
                        nodeArr[0].setMoveX(bMvXp);
                        nodeArr[0].setMoveY(bMvYp);
                        nodeArr[0].setMoveZ(bMvZm);
                        mapImageIndexToNode[0][yIndex + 1][xIndex + 1] = nodeArr[0];
                    }
                    nodeArr[1] = mapImageIndexToNode[0][yIndex + 1][xIndex];
                    if (nodeArr[1] == null) {
                        // nodeArr[1] = new cbit.vcell.geometry.surface.Node(xComp,ypComp,zComp);
                        nodeArr[1] = new cbit.vcell.geometry.surface.Node((bMvXm ? xComp : xmEdge), (bMvYp ? ypComp : ypEdge), (bMvZm ? zComp : zmEdge));
                        nodeArr[1].setGlobalIndex(nodeListV.size());
                        nodeListV.add(nodeArr[1]);
                        nodeArr[1].setMoveX(bMvXm);
                        nodeArr[1].setMoveY(bMvYp);
                        nodeArr[1].setMoveZ(bMvZm);
                        mapImageIndexToNode[0][yIndex + 1][xIndex] = nodeArr[1];
                    }
                    nodeArr[2] = mapImageIndexToNode[1][yIndex + 1][xIndex];
                    if (nodeArr[2] == null) {
                        // nodeArr[2] = new cbit.vcell.geometry.surface.Node(xComp,ypComp,zpComp);
                        nodeArr[2] = new cbit.vcell.geometry.surface.Node((bMvXm ? xComp : xmEdge), (bMvYp ? ypComp : ypEdge), (bMvZp ? zpComp : zpEdge));
                        nodeArr[2].setGlobalIndex(nodeListV.size());
                        nodeListV.add(nodeArr[2]);
                        nodeArr[2].setMoveX(bMvXm);
                        nodeArr[2].setMoveY(bMvYp);
                        nodeArr[2].setMoveZ(bMvZp);
                        mapImageIndexToNode[1][yIndex + 1][xIndex] = nodeArr[2];
                    }
                    nodeArr[3] = mapImageIndexToNode[1][yIndex + 1][xIndex + 1];
                    if (nodeArr[3] == null) {
                        // nodeArr[3] = new cbit.vcell.geometry.surface.Node(xpComp,ypComp,zpComp);
                        nodeArr[3] = new cbit.vcell.geometry.surface.Node((bMvXp ? xpComp : xpEdge), (bMvYp ? ypComp : ypEdge), (bMvZp ? zpComp : zpEdge));
                        nodeArr[3].setGlobalIndex(nodeListV.size());
                        nodeListV.add(nodeArr[3]);
                        nodeArr[3].setMoveX(bMvXp);
                        nodeArr[3].setMoveY(bMvYp);
                        nodeArr[3].setMoveZ(bMvZp);
                        mapImageIndexToNode[1][yIndex + 1][xIndex + 1] = nodeArr[3];
                    }
                    addQuadToSurface(surfQuadsV, mapInsideOutsideToSurface, // masterIndex+yStep,mapImageIndexToRegionIndex.getValue(masterIndex+yStep),
                    masterIndex, mapLinkRegionToDistinctRegion[mapImageIndexToLinkRegion[masterIndex]], masterIndex + yStep, mapLinkRegionToDistinctRegion[mapImageIndexToLinkRegion[masterIndex + yStep]], nodeArr, numRegions, MembraneElementIdentifier.PerpendicularTo.Y);
                // surfQuad =
                // new Quadrilateral(nodeArr,
                // mapLinkRegionToDistinctRegion[mapImageIndexToLinkRegion[masterIndex]],
                // mapLinkRegionToDistinctRegion[mapImageIndexToLinkRegion[masterIndex+yStep]]);
                // if(surfQuad.getVolIndexNeighbor1() > surfQuad.getVolIndexNeighbor2()){
                // surfQuad.reverseDirection();
                // }
                // surfQuadsV[surfQuad.getVolIndexNeighbor1()].add(surfQuad);
                }
                if (zSurfElements.get(masterIndex)) {
                    bMvXm = (xIndex != 0);
                    bMvXp = (xIndex != (vcImage.getNumX() - 1));
                    xComp = loX + xIndex * dX;
                    xpComp = loX + (xIndex + 1) * dX;
                    cbit.vcell.geometry.surface.Node[] nodeArr = new cbit.vcell.geometry.surface.Node[4];
                    nodeArr[0] = mapImageIndexToNode[1][yIndex][xIndex + 1];
                    if (nodeArr[0] == null) {
                        // nodeArr[0] = new cbit.vcell.geometry.surface.Node(xpComp,yComp,zpComp);
                        nodeArr[0] = new cbit.vcell.geometry.surface.Node((bMvXp ? xpComp : xpEdge), (bMvYm ? yComp : ymEdge), (bMvZp ? zpComp : zpEdge));
                        nodeArr[0].setGlobalIndex(nodeListV.size());
                        nodeListV.add(nodeArr[0]);
                        nodeArr[0].setMoveX(bMvXp);
                        nodeArr[0].setMoveY(bMvYm);
                        nodeArr[0].setMoveZ(bMvZp);
                        mapImageIndexToNode[1][yIndex][xIndex + 1] = nodeArr[0];
                    }
                    nodeArr[1] = mapImageIndexToNode[1][yIndex + 1][xIndex + 1];
                    if (nodeArr[1] == null) {
                        // nodeArr[1] = new cbit.vcell.geometry.surface.Node(xpComp,ypComp,zpComp);
                        nodeArr[1] = new cbit.vcell.geometry.surface.Node((bMvXp ? xpComp : xpEdge), (bMvYp ? ypComp : ypEdge), (bMvZp ? zpComp : zpEdge));
                        nodeArr[1].setGlobalIndex(nodeListV.size());
                        nodeListV.add(nodeArr[1]);
                        nodeArr[1].setMoveX(bMvXp);
                        nodeArr[1].setMoveY(bMvYp);
                        nodeArr[1].setMoveZ(bMvZp);
                        mapImageIndexToNode[1][yIndex + 1][xIndex + 1] = nodeArr[1];
                    }
                    nodeArr[2] = mapImageIndexToNode[1][yIndex + 1][xIndex];
                    if (nodeArr[2] == null) {
                        // nodeArr[2] = new cbit.vcell.geometry.surface.Node(xComp,ypComp,zpComp);
                        nodeArr[2] = new cbit.vcell.geometry.surface.Node((bMvXm ? xComp : xmEdge), (bMvYp ? ypComp : ypEdge), (bMvZp ? zpComp : zpEdge));
                        nodeArr[2].setGlobalIndex(nodeListV.size());
                        nodeListV.add(nodeArr[2]);
                        nodeArr[2].setMoveX(bMvXm);
                        nodeArr[2].setMoveY(bMvYp);
                        nodeArr[2].setMoveZ(bMvZp);
                        mapImageIndexToNode[1][yIndex + 1][xIndex] = nodeArr[2];
                    }
                    nodeArr[3] = mapImageIndexToNode[1][yIndex][xIndex];
                    if (nodeArr[3] == null) {
                        // nodeArr[3] = new cbit.vcell.geometry.surface.Node(xComp,yComp,zpComp);
                        nodeArr[3] = new cbit.vcell.geometry.surface.Node((bMvXm ? xComp : xmEdge), (bMvYm ? yComp : ymEdge), (bMvZp ? zpComp : zpEdge));
                        nodeArr[3].setGlobalIndex(nodeListV.size());
                        nodeListV.add(nodeArr[3]);
                        nodeArr[3].setMoveX(bMvXm);
                        nodeArr[3].setMoveY(bMvYm);
                        nodeArr[3].setMoveZ(bMvZp);
                        mapImageIndexToNode[1][yIndex][xIndex] = nodeArr[3];
                    }
                    addQuadToSurface(surfQuadsV, mapInsideOutsideToSurface, // masterIndex+zStep,mapImageIndexToRegionIndex.getValue(masterIndex+zStep),
                    masterIndex, mapLinkRegionToDistinctRegion[mapImageIndexToLinkRegion[masterIndex]], masterIndex + zStep, mapLinkRegionToDistinctRegion[mapImageIndexToLinkRegion[masterIndex + zStep]], nodeArr, numRegions, MembraneElementIdentifier.PerpendicularTo.Z);
                // surfQuad =
                // new Quadrilateral(nodeArr,
                // mapLinkRegionToDistinctRegion[mapImageIndexToLinkRegion[masterIndex]],
                // mapLinkRegionToDistinctRegion[mapImageIndexToLinkRegion[masterIndex+zStep]]);
                // if(surfQuad.getVolIndexNeighbor1() > surfQuad.getVolIndexNeighbor2()){
                // surfQuad.reverseDirection();
                // }
                // surfQuadsV[surfQuad.getVolIndexNeighbor1()].add(surfQuad);
                }
                masterIndex += 1;
            }
        }
    }
    remapQuadIndexes = new int[surfQuadsV.size()][];
    int quadCounter = 0;
    surfaceCollection = new SurfaceCollection();
    cbit.vcell.geometry.surface.Node[] allNodes = new cbit.vcell.geometry.surface.Node[nodeListV.size()];
    nodeListV.copyInto(allNodes);
    surfaceCollection.setNodes(allNodes);
    for (int surfaceIndex = 0; surfaceIndex < surfQuadsV.size(); surfaceIndex++) {
        Vector<Quadrilateral> surfV = surfQuadsV.elementAt(surfaceIndex);
        remapQuadIndexes[surfaceIndex] = new int[surfV.size()];
        OrigSurface surface = new OrigSurface(// surfV.elementAt(0).getVolIndexNeighbor1(),
        mapLinkRegionToDistinctRegion[mapImageIndexToLinkRegion[surfV.elementAt(0).getVolIndexNeighbor1()]], // surfV.elementAt(0).getVolIndexNeighbor2()
        mapLinkRegionToDistinctRegion[mapImageIndexToLinkRegion[surfV.elementAt(0).getVolIndexNeighbor2()]]);
        for (int faceIndex = 0; faceIndex < surfV.size(); faceIndex++) {
            surface.addPolygon(surfV.elementAt(faceIndex));
            remapQuadIndexes[surfaceIndex][faceIndex] = quadCounter;
            quadIndexToSurfAndFace.add(new SurfAndFace(surfaceIndex, faceIndex));
            quadCounter++;
        }
        surfaceCollection.addSurface(surface);
    }
    try {
        if (!bMembraneNeighborCalculationFailed) {
            calculateNeighbors();
            surfaceCollection.setMembraneEdgeNeighbors(membraneEdgeNeighbors);
        }
    } catch (Exception e) {
        // If MembraneNeighbors fails somewhere let Surfacecollection generate without membraneneighbors (original behavior just in case)
        e.printStackTrace();
        bMembraneNeighborCalculationFailed = true;
    } finally {
        // these aren't needed after fail or surfacecollection.set(...)
        remapQuadIndexes = null;
        edgeMap = null;
    }
// RegionImage.sortSurfaceCollection(surfaceCollection);
// //check
// for (int i = 0; i < surfaceCollection.getSurfaceCount(); i++) {
// Surface surface = surfaceCollection.getSurfaces(i);
// for (int j = 0; j < surface.getPolygonCount(); j++) {
// Polygon polygon = surface.getPolygons(j);
// for (int k = 0; k < polygon.getNodes().length; k++) {
// if(surfaceCollection.getNodes()[polygon.getNodes()[k].getGlobalIndex()] != polygon.getNodes()[k]){
// throw new RuntimeException("Nodes not match");
// }
// }
// }
// }
}