Revision: 10627
http://jmol.svn.sourceforge.net/jmol/?rev=10627&view=rev
Author: hansonr
Date: 2009-02-08 14:54:57 +0000 (Sun, 08 Feb 2009)
Log Message:
-----------
version=11.7.25_dev # code: rewritten code for ASimpleJvxlWriter allows streaming data
Modified Paths:
--------------
trunk/Jmol/src/org/jmol/jvxl/api/VertexDataServer.java
trunk/Jmol/src/org/jmol/jvxl/calc/MarchingCubes.java
trunk/Jmol/src/org/jmol/jvxl/calc/MarchingSquares.java
trunk/Jmol/src/org/jmol/jvxl/readers/IsoFxyReader.java
trunk/Jmol/src/org/jmol/jvxl/readers/IsoPlaneReader.java
trunk/Jmol/src/org/jmol/jvxl/readers/IsoShapeReader.java
trunk/Jmol/src/org/jmol/jvxl/readers/JvxlReader.java
trunk/Jmol/src/org/jmol/jvxl/readers/SurfaceReader.java
trunk/Jmol/src/org/jmol/jvxl/readers/VolumeDataReader.java
trunk/Jmol/src/org/jmol/jvxl/readers/VolumeFileReader.java
trunk/Jmol/src/org/jmol/shapesurface/Isosurface.java
trunk/Jmol/src/org/jmol/viewer/Jmol.properties
trunk/Jmol/src/org/openscience/jvxl/simplewriter/ASimpleJvxlWriter.java
trunk/Jmol/src/org/openscience/jvxl/simplewriter/JvxlWrite.java
trunk/Jmol/src/org/openscience/jvxl/simplewriter/SimpleMarchingCubes.java
trunk/Jmol/src/org/openscience/jvxl/simplewriter/VolumeData.java
trunk/Jmol/src/org/openscience/jvxl/simplewriter/VoxelDataCreator.java
Modified: trunk/Jmol/src/org/jmol/jvxl/api/VertexDataServer.java
===================================================================
--- trunk/Jmol/src/org/jmol/jvxl/api/VertexDataServer.java 2009-02-07 05:32:40 UTC (rev 10626)
+++ trunk/Jmol/src/org/jmol/jvxl/api/VertexDataServer.java 2009-02-08 14:54:57 UTC (rev 10627)
@@ -59,15 +59,16 @@
* @param pointA
* @param edgeVector vector from A to B
* @param isContourType
+ * @param fReturn
* @return new vertex index or Integer.MAX_VALUE
*/
- public abstract int getSurfacePointIndex(float cutoff,
+ public abstract int getSurfacePointIndexAndFraction(float cutoff,
boolean isCutoffAbsolute, int x,
int y, int z, Point3i offset,
int vertexA, int vertexB,
float valueA, float valueB,
Point3f pointA, Vector3f edgeVector,
- boolean isContourType);
+ boolean isContourType, float[] fReturn);
/**
* addVertexCopy is used by the Marching Squares algorithm to
@@ -96,5 +97,16 @@
* @param isAbsolute
*/
public abstract void addTriangleCheck(int iA, int iB, int iC, int check,
- boolean isAbsolute);
+ boolean isAbsolute);
+
+ /**
+ *
+ * for readers only
+ *
+ * @param x
+ * @param y
+ * @param z
+ * @return value[x][y][z]
+ */
+ public float getValue(int x, int y, int z);
}
Modified: trunk/Jmol/src/org/jmol/jvxl/calc/MarchingCubes.java
===================================================================
--- trunk/Jmol/src/org/jmol/jvxl/calc/MarchingCubes.java 2009-02-07 05:32:40 UTC (rev 10626)
+++ trunk/Jmol/src/org/jmol/jvxl/calc/MarchingCubes.java 2009-02-08 14:54:57 UTC (rev 10627)
@@ -31,6 +31,8 @@
import org.jmol.jvxl.api.VertexDataServer;
import org.jmol.jvxl.data.VolumeData;
+import org.jmol.jvxl.readers.JvxlReader;
+
//import org.jmol.util.Base64;
public class MarchingCubes {
@@ -43,6 +45,26 @@
*
* Author: Bob Hanson, hansonr@...
*
+ * inputs: surfaceReader: interface to other methods
+ * volumeData, containing all information relating to origin, axes, etc.
+ *
+ * MODE_CUBE volumeData.voxelData possibly with 3D voxel data
+ * MODE_BITSET bsVoxels -- optional alternative, JVXL-type BitSet indicating which vertices are inside and which outside
+ * in order for(x 0 to nX){for(y 0 to nY){for(z 0 to nZ){}}}
+ * MODE_GETXYZ If BOTH voxelData and bsVoxels are null,
+ * then we assume we can get the voxelData on the fly
+ * using surfaceReader.getValue(x, y, z)
+ *
+ * surfaceReader.getSurfacePointIndex gets the actual Point3f vertex points and
+ * returns the fraction information. It is exported here, because the JVXL reader may use this
+ * as an opportunity to deliver the fractional information, thus providing the
+ * actual position of the vertex point on the isosurface.
+ *
+ * outputs: Point3f positions of isosurface intersections with grid via SurfaceReader.addVertex()
+ * Point4f triangle data via surfaceReader.addTriangleCheck()
+ * bsVoxels -- returned same (JVXL file data) or filled (otherwise)
+ * edgeData -- encoded fraction data as a string
+ *
*/
private VertexDataServer surfaceReader;
@@ -51,56 +73,69 @@
private boolean isContoured;
private float cutoff;
private boolean isCutoffAbsolute;
+ private boolean isSquared;
+ private boolean isXLowToHigh;
private int cubeCountX, cubeCountY, cubeCountZ;
private int nY, nZ;
private BitSet bsVoxels;
- /* this was an idea .... but it results in 3x the file size
+ private StringBuffer edgeData = new StringBuffer();
- private StringBuffer maskData = new StringBuffer();
+ public BitSet getBsVoxels() {
+ return bsVoxels;
+ }
- public String getMaskData() {
- return Base64.getBase64(maskData).toString();
- }
+ /* this was an idea .... but it results in 3x the file size
+
+ private StringBuffer maskData = new StringBuffer();
+
+ public String getMaskData() {
+ return Base64.getBase64(maskData).toString();
+ }
- */
-
+ */
+
public MarchingCubes(VertexDataServer surfaceReader, VolumeData volumeData,
BitSet bsVoxels, boolean isContoured, int contourType, float cutoff,
- boolean isCutoffAbsolute) {
-
- // when just creating a JVXL file all you really need are:
+ boolean isCutoffAbsolute, boolean isSquared, boolean isXLowToHigh) {
+
+ // If just creating a JVXL file, see org.openscience.jmol.jvxl.simplewriter.SimpleMarchingCubes.java
//
- // volumeData.voxelData[x][y][z]
- // cutoff
- //
-
+
this.surfaceReader = surfaceReader;
this.volumeData = volumeData;
this.isContoured = isContoured;
this.cutoff = cutoff;
this.isCutoffAbsolute = isCutoffAbsolute;
this.contourType = contourType;
- if (bsVoxels == null) {
- cubeCountX = volumeData.voxelData.length - 1;
- cubeCountY = (nY = volumeData.voxelData[0].length) - 1;
- cubeCountZ = (nZ = volumeData.voxelData[0][0].length) - 1;
+ this.isSquared = isSquared;
+ this.isXLowToHigh = isXLowToHigh; // generally ixXLowToHigh is FALSE
+
+ cubeCountX = volumeData.voxelCounts[0] - 1;
+ cubeCountY = (nY = volumeData.voxelCounts[1]) - 1;
+ cubeCountZ = (nZ = volumeData.voxelCounts[2]) - 1;
+ if (volumeData.voxelData != null) {
+ mode = MODE_CUBE;
} else {
- cubeCountX = volumeData.voxelCounts[0] - 1;
- cubeCountY = (nY = volumeData.voxelCounts[1]) - 1;
- cubeCountZ = (nZ = volumeData.voxelCounts[2]) - 1;
+ mode = (bsVoxels == null ? MODE_GETXYZ : MODE_BITSET);
this.bsVoxels = bsVoxels;
}
yzCount = nY * nZ;
-
+ if (this.bsVoxels == null)
+ this.bsVoxels = new BitSet();
// calcVoxelVertexVectors is unnecessary if just creating a JVXL file:
-
+
calcVoxelVertexVectors();
}
- public final Vector3f[] voxelVertexVectors = new Vector3f[8];
+ private int mode;
+ private final static int MODE_CUBE = 1;
+ private final static int MODE_BITSET = 2;
+ private final static int MODE_GETXYZ = 3;
+
+ private final Vector3f[] voxelVertexVectors = new Vector3f[8];
private final float[] vertexValues = new float[8];
private final Point3i[] vertexPoints = new Point3i[8];
private final Vector3f[] edgeVectors = new Vector3f[12];
@@ -111,12 +146,10 @@
vertexPoints[i] = new Point3i();
}
- boolean isXLowToHigh;
int edgeCount;
-
private void calcVoxelVertexVectors() {
- setLinearOffsets();
+ setLinearOffsets();
volumeData.setMatrix();
for (int i = 8; --i >= 0;)
volumeData.transform(cubeVertexVectors[i],
@@ -125,48 +158,43 @@
edgeVectors[i].sub(voxelVertexVectors[edgeVertexes[i + i + 1]],
voxelVertexVectors[edgeVertexes[i + i]]);
}
-
-
- /* Note to Jason from Bob:
+
+ /* see also org.openscience.jmol.jvxl.simplewriter.SimpleMarchingCubes.java
+ * for a streamlined version of this method that does simple writing
+ * of JVXL files from cube data.
*
- * To just create a JVXL file, you need these five methods.
- * Their output is the fractionData string buffer and the
- * number of surface points
- *
- * inputs required:
- *
- * 1) volumeData.voxelData[x][y][z]
- * 2) cutoff
- * 3) values created in MarchingCubes constructor
- *
- * The first four methods are in org.jmol.jvxl.calc.MarchingCubes.java
- *
- * generateSurfaceData -- isXLowToHigh false; isContoured false
- * -- triangle stuff at end not needed
- * propagateNeighborPointIndexes -- EXACTLY as is, no changes allowed
- * isInside -- EXACTLY as is -- defines what "inside" means
- * processOneCubical -- EXACTLY as is, no changes at all
- * SurfaceReader.getSurfacePointIndex -- your job
- * -- receives the point value data and positions
- * -- responsible for creating the fractionData character buffer
- * -- just return 0 since you are not creating triangles
- *
*/
- public int generateSurfaceData(boolean isXLowToHigh) {
+ private static int[] xyPlanePts = new int[] { 0, 1, 1, 0, 0, 1, 1, 0 };
- // generally ixXLowToHigh is FALSE
+ public String getEdgeData() {
- this.isXLowToHigh = isXLowToHigh;
-
+ /* The Marching Cubes code creates a (nY-1)(nZ-1)(12) array
+ * to track all 12 edges of a full slice of cubes. I don't think that is
+ * necessary -- really each vertex is associated with at most three
+ * edges. We should be able to do this with an (nY)(nZ)(3) array.
+ *
+ * For now, I've left it as is.
+ *
+ * The changes introduced here But the other aspect of this is that
+ * we really only need (nY)(nZ)(2) points ever. So we now allow for the
+ * method to query another method to
+ */
+
+
+ /* isoPointIndexes
// set up the set of edge points in the YZ plane
// these are indexes into an array of Point3f values
// They will be initialized as -1 whenever a vertex is needed.
// But if just creating a JVXL file, all you need to do
// is set them to 0, not an index into any actual array.
+ *
+ */
int[][] isoPointIndexes = new int[cubeCountY * cubeCountZ][12];
+ float[][] xyPlanes = (mode == MODE_GETXYZ ? new float[2][yzCount] : null);
+
int insideCount = 0, outsideCount = 0, surfaceCount = 0;
edgeCount = 0;
@@ -188,14 +216,21 @@
// we are starting at the top corner, in the next to last
// cell on the next to last row of the next to last plane(!)
}
- for (int x = x0; x != x1; x += xStep, ptX += ptStep, pt = ptX) {
+ int cellIndex0 = cubeCountY * cubeCountZ - 1;
+ int cellIndex = cellIndex0;
+ for (int x = x0; x != x1; x += xStep, ptX += ptStep, pt = ptX, cellIndex = cellIndex0) {
+ if (mode == MODE_GETXYZ) {
+ float[] plane = xyPlanes[0];
+ xyPlanes[0] = xyPlanes[1];
+ xyPlanes[1] = plane;
+ }
for (int y = cubeCountY; --y >= 0; pt--) {
- for (int z = cubeCountZ; --z >= 0; pt--) {
+ for (int z = cubeCountZ; --z >= 0; pt--, cellIndex--) {
// set up the list of indices that need checking
int[] voxelPointIndexes = propagateNeighborPointIndexes(x, y, z,
- isoPointIndexes);
+ isoPointIndexes, cellIndex);
// create the bitset mask indicating which vertices are inside.
// 0xFF here means "all inside"; 0x00 means "all outside"
@@ -207,17 +242,30 @@
// to our base x,y,z cube position
boolean isInside;
- if (bsVoxels == null) {
- Point3i offset = cubeVertexOffsets[i];
+ Point3i offset = cubeVertexOffsets[i];
+ int pti = pt + linearOffsets[i];
+ switch (mode) {
+ case MODE_GETXYZ:
+ vertexValues[i] = getValue(i, x + offset.x, y + offset.y, z
+ + offset.z, pti, xyPlanes[xyPlanePts[i]]);
+ isInside = bsVoxels.get(pti);
+ break;
+ case MODE_BITSET:
+ isInside = bsVoxels.get(pti);
+ vertexValues[i] = (isInside ? 1 : 0);
+ break;
+ default:
+ case MODE_CUBE:
vertexValues[i] = volumeData.voxelData[x + offset.x][y + offset.y][z
+ offset.z];
+ if (isSquared)
+ vertexValues[i] *= vertexValues[i];
isInside = isInside(vertexValues[i], cutoff, isCutoffAbsolute);
- } else {
- isInside = bsVoxels.get(pt + linearOffsets[i]);
- vertexValues[i] = (isInside ? 1 : 0);
+ bsVoxels.set(pti);
}
- if (isInside)
+ if (isInside) {
insideMask |= 1 << i;
+ }
}
//maskData.append((char) insideMask);
@@ -234,6 +282,7 @@
// This cube is straddling the cutoff. We must check all edges
+ //System.out.print( " pt=" + pt + " ");
if (!processOneCubical(insideMask, voxelPointIndexes, x, y, z)
|| isContoured)
continue;
@@ -251,7 +300,6 @@
}
}
-
//System.out.println("marchingCubes: "
// + maskData.length()
// + " "
@@ -259,10 +307,25 @@
// Base64.getBase64(maskData).toString()).length());
// this came out about 3 times the length of what we have now
// could be a hair faster, but the files size cost is high, I think
-
- return edgeCount;
+
+ return edgeData.toString();
}
-
+
+ private BitSet bsValues = new BitSet();
+
+ private float getValue(int i, int x, int y, int z, int pt, float[] tempValues) {
+ if (bsValues.get(pt))
+ return tempValues[pt % yzCount];
+ bsValues.set(pt);
+ float value = surfaceReader.getValue(x, y, z);
+ if (isSquared)
+ value *= value;
+ tempValues[pt % yzCount] = value;
+ if (isInside(value, cutoff, isCutoffAbsolute))
+ bsVoxels.set(pt);
+ return value;
+ }
+
public static boolean isInside(float voxelValue, float max, boolean isAbsolute) {
return ((max > 0 && (isAbsolute ? Math.abs(voxelValue) : voxelValue) >= max) || (max <= 0 && voxelValue <= max));
}
@@ -271,7 +334,8 @@
-1, -1 };
private int[] propagateNeighborPointIndexes(int x, int y, int z,
- int[][] isoPointIndexes) {
+ int[][] isoPointIndexes,
+ int cellIndex) {
/* Y
* 4 --------4--------- 5
* /| /|
@@ -329,169 +393,153 @@
*
*
*/
-
+
/* DO NOT EVER CHANGE THIS */
-
- int cellIndex = y * cubeCountZ + z;
+
+ //int cellIndex = y * cubeCountZ + z;
int[] voxelPointIndexes = isoPointIndexes[cellIndex];
- boolean noXNeighbor = (x == cubeCountX - 1);
+ boolean noYNeighbor = (y == cubeCountY - 1);
+ int[] yNeighbor = noYNeighbor ? nullNeighbor
+ : isoPointIndexes[cellIndex + cubeCountZ];
+ boolean noZNeighbor = (z == cubeCountZ - 1);
+ int[] zNeighbor = noZNeighbor ? nullNeighbor
+ : isoPointIndexes[cellIndex + 1];
+ voxelPointIndexes[0] = -1;
+ voxelPointIndexes[2] = zNeighbor[0];
+ voxelPointIndexes[4] = yNeighbor[0];
+ voxelPointIndexes[6] = (noYNeighbor ? zNeighbor[4] : yNeighbor[2]);
+
if (isXLowToHigh) {
// reading x from low to high
- if (noXNeighbor) {
+ if (x == 0) {
voxelPointIndexes[3] = -1;
voxelPointIndexes[8] = -1;
- voxelPointIndexes[7] = -1;
- voxelPointIndexes[11] = -1;
+ voxelPointIndexes[7] = yNeighbor[3];
+ voxelPointIndexes[11] = zNeighbor[8];
} else {
voxelPointIndexes[3] = voxelPointIndexes[1];
voxelPointIndexes[7] = voxelPointIndexes[5];
voxelPointIndexes[8] = voxelPointIndexes[9];
voxelPointIndexes[11] = voxelPointIndexes[10];
}
+ voxelPointIndexes[1] = -1;
+ voxelPointIndexes[5] = yNeighbor[1];
+ voxelPointIndexes[9] = -1;
+ voxelPointIndexes[10] = zNeighbor[9];
} else {
// reading x from high to low
- if (noXNeighbor) {
- // the x neighbor is myself from my last pass through here
+ if (x == cubeCountX - 1) {
voxelPointIndexes[1] = -1;
+ voxelPointIndexes[5] = yNeighbor[1];
voxelPointIndexes[9] = -1;
- voxelPointIndexes[5] = -1;
- voxelPointIndexes[10] = -1;
+ voxelPointIndexes[10] = zNeighbor[9];
} else {
voxelPointIndexes[1] = voxelPointIndexes[3];
voxelPointIndexes[5] = voxelPointIndexes[7];
voxelPointIndexes[9] = voxelPointIndexes[8];
voxelPointIndexes[10] = voxelPointIndexes[11];
}
- }
- //from the y neighbor pick up the top
- boolean noYNeighbor = (y == cubeCountY - 1);
- int[] yNeighbor = noYNeighbor ? nullNeighbor : isoPointIndexes[cellIndex
- + cubeCountZ];
-
- voxelPointIndexes[4] = yNeighbor[0];
- voxelPointIndexes[6] = yNeighbor[2];
-
- if (isXLowToHigh) {
- voxelPointIndexes[5] = yNeighbor[1];
- if (noXNeighbor)
- voxelPointIndexes[7] = yNeighbor[3];
- } else {
+ voxelPointIndexes[3] = -1;
voxelPointIndexes[7] = yNeighbor[3];
- if (noXNeighbor)
- voxelPointIndexes[5] = yNeighbor[1];
+ voxelPointIndexes[8] = -1;
+ voxelPointIndexes[11] = zNeighbor[8];
}
- // from my z neighbor
- boolean noZNeighbor = (z == cubeCountZ - 1);
- int[] zNeighbor = noZNeighbor ? nullNeighbor
- : isoPointIndexes[cellIndex + 1];
- voxelPointIndexes[2] = zNeighbor[0];
- if (noYNeighbor)
- voxelPointIndexes[6] = zNeighbor[4];
- if (isXLowToHigh) {
- if (noXNeighbor)
- voxelPointIndexes[11] = zNeighbor[8];
- voxelPointIndexes[10] = zNeighbor[9];
- } else {
- if (noXNeighbor)
- voxelPointIndexes[10] = zNeighbor[9];
- voxelPointIndexes[11] = zNeighbor[8];
- }
- // these must always be calculated
- voxelPointIndexes[0] = -1;
- if (isXLowToHigh) {
- voxelPointIndexes[1] = -1;
- voxelPointIndexes[9] = -1;
- } else {
- voxelPointIndexes[3] = -1;
- voxelPointIndexes[8] = -1;
- }
return voxelPointIndexes;
}
private final Point3f pt0 = new Point3f();
private final Point3f pointA = new Point3f();
-
+
+ private static final int[] Pwr2 = new int[] { 1, 2, 4, 8, 16, 32, 64, 128,
+ 256, 512, 1024, 2048 };
+
private boolean processOneCubical(int insideMask, int[] voxelPointIndexes,
int x, int y, int z) {
-
+
// the key to the algorithm is that we have a catalog that
// maps the inside-vertex mask to an edge mask.
-
+
int edgeMask = insideMaskTable[insideMask];
boolean isNaN = false;
for (int iEdge = 12; --iEdge >= 0;) {
-
+
// bit set to one means it's a relevant edge
-
- if ((edgeMask & (1 << iEdge)) == 0)
+
+ if ((edgeMask & Pwr2[iEdge]) == 0)
continue;
-
+
// if we have a point already, we don't need to check this edge.
// for triangles, this will be an index into an array;
// for just creating JVXL files, this can just be 0
-
+
if (voxelPointIndexes[iEdge] >= 0)
continue; // propagated from neighbor
-
+
// here's an edge that has to be checked.
-
+
++edgeCount;
-
+
// get the vertex numbers 0 - 7
-
- int vertexA = edgeVertexes[2 * iEdge];
- int vertexB = edgeVertexes[2 * iEdge + 1];
-
+
+ int vertexA = edgeVertexes[iEdge << 1];
+ int vertexB = edgeVertexes[(iEdge << 1) + 1];
+
// pick up the actual value at each vertex
// this array of 8 values is updated as we go.
-
+
float valueA = vertexValues[vertexA];
float valueB = vertexValues[vertexB];
-
+
// we allow for NaN values -- missing triangles
-
+
if (Float.isNaN(valueA) || Float.isNaN(valueB))
isNaN = true;
-
+
// the exact point position -- not important for just
// creating the JVXL file. In that case, all you
// need are the two values valueA and valueB and the cutoff.
// from those you can define the fractional offset
-
+
volumeData.voxelPtToXYZ(x, y, z, pt0);
pointA.add(pt0, voxelVertexVectors[vertexA]);
-
+
// here is where we get the value and assign the point for that edge
// it is where the JVXL surface data line is appended
+
+ //System.out.print("edge" + iEdge);
+ voxelPointIndexes[iEdge] = surfaceReader.getSurfacePointIndexAndFraction(cutoff,
+ isCutoffAbsolute, x, y, z, cubeVertexOffsets[vertexA], vertexA,
+ vertexB, valueA, valueB, pointA, edgeVectors[iEdge],
+ edgeTypeTable[iEdge] == contourType, fReturn);
- voxelPointIndexes[iEdge] = surfaceReader.getSurfacePointIndex(
- cutoff, isCutoffAbsolute, x, y, z, cubeVertexOffsets[vertexA],
- vertexA, vertexB, valueA, valueB, pointA, edgeVectors[iEdge],
- edgeTypeTable[iEdge] == contourType);
+ edgeData.append(JvxlReader.jvxlFractionAsCharacter(fReturn[0]));
}
+ //System.out.println("");
return !isNaN;
}
+ private float[] fReturn = new float[1];
+
public void calcVertexPoint(int x, int y, int z, int vertex, Point3f pt) {
volumeData.voxelPtToXYZ(x, y, z, pt0);
pt.add(pt0, voxelVertexVectors[vertex]);
}
- final static Point3i[] cubeVertexOffsets = {
- new Point3i(0, 0, 0), //0 pt
- new Point3i(1, 0, 0), //1 pt + yz
- new Point3i(1, 0, 1), //2 pt + yz + 1
- new Point3i(0, 0, 1), //3 pt + 1
- new Point3i(0, 1, 0), //4 pt + z
- new Point3i(1, 1, 0), //5 pt + yz + z
- new Point3i(1, 1, 1), //6 pt + yz + z + 1
- new Point3i(0, 1, 1) //7 pt + z + 1
- };
+ final static Point3i[] cubeVertexOffsets = { new Point3i(0, 0, 0), //0 pt
+ new Point3i(1, 0, 0), //1 pt + yz
+ new Point3i(1, 0, 1), //2 pt + yz + 1
+ new Point3i(0, 0, 1), //3 pt + 1
+ new Point3i(0, 1, 0), //4 pt + z
+ new Point3i(1, 1, 0), //5 pt + yz + z
+ new Point3i(1, 1, 1), //6 pt + yz + z + 1
+ new Point3i(0, 1, 1) //7 pt + z + 1
+ };
private final int[] linearOffsets = new int[8];
int yzCount;
+
/* set the linear offsets for unique cell ID
* and for pointing into the inside/outside BitSet.
* Add offset to 0: x * (nY * nZ) + y * nZ + z
@@ -506,17 +554,15 @@
linearOffsets[6] = yzCount + nZ + 1;
linearOffsets[7] = nZ + 1;
}
-
+
public int getLinearOffset(int x, int y, int z, int offset) {
return x * yzCount + y * nZ + z + linearOffsets[offset];
}
-
+
static Vector3f[] cubeVertexVectors = { new Vector3f(0, 0, 0),
new Vector3f(1, 0, 0), new Vector3f(1, 0, 1), new Vector3f(0, 0, 1),
new Vector3f(0, 1, 0), new Vector3f(1, 1, 0), new Vector3f(1, 1, 1),
new Vector3f(0, 1, 1) };
-
-
/* Y
* 4 --------4--------- 5 +z --------4--------- +yz+z
@@ -570,17 +616,14 @@
*
*
*/
-
-
- private final static int edgeTypeTable[] = {
- 0, 2, 0, 2, 0, 2, 0, 2, 1, 1, 1, 1 };
+ private final static int edgeTypeTable[] = { 0, 2, 0, 2, 0, 2, 0, 2, 1, 1, 1,
+ 1 };
// 0=along X, 1=along Y, 2=along Z
- private final static byte edgeVertexes[] = {
- 0, 1, 1, 2, 2, 3, 3, 0, 4, 5,
+ private final static byte edgeVertexes[] = { 0, 1, 1, 2, 2, 3, 3, 0, 4, 5,
/*0 1 2 3 4 */
- 5, 6, 6, 7, 7, 4, 0, 4, 1, 5, 2, 6, 3, 7 };
+ 5, 6, 6, 7, 7, 4, 0, 4, 1, 5, 2, 6, 3, 7 };
/*5 6 7 8 9 10 11 */
private final static short insideMaskTable[] = { 0x0000, 0x0109, 0x0203,
Modified: trunk/Jmol/src/org/jmol/jvxl/calc/MarchingSquares.java
===================================================================
--- trunk/Jmol/src/org/jmol/jvxl/calc/MarchingSquares.java 2009-02-07 05:32:40 UTC (rev 10626)
+++ trunk/Jmol/src/org/jmol/jvxl/calc/MarchingSquares.java 2009-02-08 14:54:57 UTC (rev 10627)
@@ -312,7 +312,9 @@
}
void setValue(float value, VolumeData volumeData) {
- this.value = volumeData.voxelData[voxelLocation.x][voxelLocation.y][voxelLocation.z] = value;
+ this.value = value;
+ if (volumeData.voxelData != null)
+ volumeData.voxelData[voxelLocation.x][voxelLocation.y][voxelLocation.z] = value;
}
void setPixelLocation(Point3i pt) {
Modified: trunk/Jmol/src/org/jmol/jvxl/readers/IsoFxyReader.java
===================================================================
--- trunk/Jmol/src/org/jmol/jvxl/readers/IsoFxyReader.java 2009-02-07 05:32:40 UTC (rev 10626)
+++ trunk/Jmol/src/org/jmol/jvxl/readers/IsoFxyReader.java 2009-02-08 14:54:57 UTC (rev 10627)
@@ -58,7 +58,7 @@
JvxlReader.jvxlCreateHeaderWithoutTitleOrAtoms(volumeData, jvxlFileHeaderBuffer);
}
- protected float getValue(int x, int y, int z) {
+ public float getValue(int x, int y, int z) {
//System.out.println("isofxyreader getValue x y z " + x + " " + y + " " + z + " zreal " + (volumeData.origin[2] + z * volumeData.volumetricVectors[2].z) + " data " + data[x][y] + " diff " + ((volumeData.origin[2] + z * volumeData.volumetricVectors[2].z)-data[x][y]));
return (isPlanarMapping ? data[x][y] : data[x][y] - (volumeData.origin[2] + z * volumeData.volumetricVectors[2].z));
}
Modified: trunk/Jmol/src/org/jmol/jvxl/readers/IsoPlaneReader.java
===================================================================
--- trunk/Jmol/src/org/jmol/jvxl/readers/IsoPlaneReader.java 2009-02-07 05:32:40 UTC (rev 10626)
+++ trunk/Jmol/src/org/jmol/jvxl/readers/IsoPlaneReader.java 2009-02-08 14:54:57 UTC (rev 10627)
@@ -39,7 +39,7 @@
params.cutoff = 0;
}
- protected float getValue(int x, int y, int z) {
+ public float getValue(int x, int y, int z) {
return volumeData.calcVoxelPlaneDistance(x, y, z);
}
Modified: trunk/Jmol/src/org/jmol/jvxl/readers/IsoShapeReader.java
===================================================================
--- trunk/Jmol/src/org/jmol/jvxl/readers/IsoShapeReader.java 2009-02-07 05:32:40 UTC (rev 10626)
+++ trunk/Jmol/src/org/jmol/jvxl/readers/IsoShapeReader.java 2009-02-08 14:54:57 UTC (rev 10627)
@@ -108,7 +108,7 @@
setHeader(type + "\n");
}
- protected float getValue(int x, int y, int z) {
+ public float getValue(int x, int y, int z) {
volumeData.voxelPtToXYZ(x, y, z, ptPsi);
ptPsi.sub(center);
if (isEccentric)
Modified: trunk/Jmol/src/org/jmol/jvxl/readers/JvxlReader.java
===================================================================
--- trunk/Jmol/src/org/jmol/jvxl/readers/JvxlReader.java 2009-02-07 05:32:40 UTC (rev 10626)
+++ trunk/Jmol/src/org/jmol/jvxl/readers/JvxlReader.java 2009-02-08 14:54:57 UTC (rev 10627)
@@ -468,7 +468,7 @@
return dataOut.toString();
}
- protected BitSet getVoxelBitSet(int nPoints, StringBuffer sb) throws Exception {
+ protected BitSet getVoxelBitSet(int nPoints) throws Exception {
BitSet bs = new BitSet();
int bsVoxelPtr = 0;
if (surfaceDataCount <= 0)
@@ -486,8 +486,6 @@
endOfData = true;
nThisValue = 10000;
//throw new NullPointerException();
- } else if (sb != null) {
- sb.append(line).append('\n');
}
}
thisInside = !thisInside;
@@ -526,6 +524,27 @@
return (thisInside ? 1f : 0f);
}
+ public static void setSurfaceInfoFromBitSet(JvxlData jvxlData, BitSet bs,
+ Point4f thePlane) {
+ boolean inside = false;
+ int dataCount = 0;
+ StringBuffer sb = new StringBuffer();
+ int nSurfaceInts = 0;
+ int nPoints = jvxlData.nPointsX * jvxlData.nPointsY * jvxlData.nPointsZ;
+ for (int i = 0; i < nPoints; ++i) {
+ if (inside == bs.get(i)) {
+ dataCount++;
+ } else {
+ sb.append(' ').append(dataCount);
+ nSurfaceInts++;
+ dataCount = 1;
+ inside = !inside;
+ }
+ }
+ sb.append(' ').append(dataCount).append('\n');
+ setSurfaceInfo(jvxlData, thePlane, nSurfaceInts, sb);
+ }
+
protected static void setSurfaceInfo(JvxlData jvxlData, Point4f thePlane, int nSurfaceInts, StringBuffer surfaceData) {
jvxlData.jvxlSurfaceData = surfaceData.toString();
if (jvxlData.jvxlSurfaceData.indexOf("--") == 0)
@@ -534,14 +553,16 @@
jvxlData.nSurfaceInts = nSurfaceInts;
}
- protected float readSurfacePoint(float cutoff, boolean isCutoffAbsolute, float valueA,
+ protected float getSurfacePointAndFraction(float cutoff, boolean isCutoffAbsolute, float valueA,
float valueB, Point3f pointA, Vector3f edgeVector,
float[] fReturn, Point3f ptReturn) {
if (edgeDataCount <= 0)
- return super.readSurfacePoint(cutoff, isCutoffAbsolute, valueA, valueB,
+ return super.getSurfacePointAndFraction(cutoff, isCutoffAbsolute, valueA, valueB,
pointA, edgeVector, fReturn, ptReturn);
ptReturn.scaleAdd(fReturn[0] = jvxlGetNextFraction(edgeFractionBase, edgeFractionRange, 0.5f),
edgeVector, pointA);
+ //System.out.println(" jvxl: f=" + fReturn[0]);
+
return fReturn[0];
}
@@ -551,12 +572,9 @@
private float jvxlGetNextFraction(int base, int range, float fracOffset) {
if (fractionPtr >= strFractionTemp.length()) {
if (!endOfData)
- Logger.error("end of file reading compressed fraction data at point "
- + fractionData.length());
+ Logger.error("end of file reading compressed fraction data");
endOfData = true;
strFractionTemp = "" + (char) base;
- fractionData.append(strFractionTemp);
- fractionData.append('\n');
fractionPtr = 0;
}
return jvxlFractionFromCharacter(strFractionTemp.charAt(fractionPtr++),
@@ -571,7 +589,6 @@
int vertexCount = jvxlData.vertexCount = meshData.vertexCount;
short[] colixes = meshData.vertexColixes;
float[] vertexValues = meshData.vertexValues;
- fractionData = new StringBuffer();
strFractionTemp = (isJvxl ? jvxlColorDataRead : "");
if (isJvxl && strFractionTemp.length() == 0) {
Logger
@@ -744,38 +761,6 @@
.append(pt.z).append(" //BOGUS He ATOM ADDED FOR JVXL FORMAT\n");
}
- public static int jvxlCreateSurfaceData(JvxlData jvxlData, float[][][] voxelData, float cutoff, boolean isCutoffAbsolute, int nX, int nY, int nZ) {
- StringBuffer sb = new StringBuffer();
- boolean inside = false;
- int dataCount = 0;
- int nDataPoints = 0;
- int nSurfaceInts = 0;
- BitSet bs = new BitSet(); //TESTING
- for (int x = 0; x < nX; ++x)
- for (int y = 0; y < nY; ++y)
- for (int z = 0; z < nZ; ++z) {
- boolean itest = isInside(voxelData[x][y][z], cutoff, isCutoffAbsolute);
- if (itest)bs.set(nDataPoints);
- ++nDataPoints;
- if (inside == itest) {
- dataCount++;
- } else {
- if (dataCount != 0) {
- sb.append(' ').append(dataCount);
- ++nSurfaceInts;
- }
- dataCount = 1;
- inside = !inside;
- }
- }
- sb.append(' ').append(dataCount).append('\n');
- ++nSurfaceInts;
- //System.out.println("JvxlReader: " + sb.length() + " " + nDataPoints);
-
- setSurfaceInfo(jvxlData,null, nSurfaceInts, sb);
- return nDataPoints;
- }
-
public static String jvxlGetDefinitionLine(JvxlData jvxlData, boolean isInfo) {
String definitionLine = jvxlData.cutoff + " ";
@@ -972,6 +957,7 @@
return 0.999999f;
//if (logCompression)
//Logger.info("ffc: " + fraction + " <-- " + ich + " " + (char) ich);
+ //System.out.println("ffc: " + fraction + " <-- " + ich + " " + (char) ich);
return fraction;
}
@@ -1002,6 +988,10 @@
return jvxlFractionAsCharacter(fraction, base, range);
}
+ public static char jvxlFractionAsCharacter(float fraction) {
+ return jvxlFractionAsCharacter(fraction, defaultEdgeFractionBase, defaultEdgeFractionRange);
+ }
+
protected static char jvxlFractionAsCharacter(float fraction, int base, int range) {
if (fraction > 0.9999f)
fraction = 0.9999f;
@@ -1475,5 +1465,5 @@
private static void setNext(String data, String what, int[] next, int offset) {
next[0] = data.indexOf(what, next[0]) + what.length() + offset;
}
-
+
}
Modified: trunk/Jmol/src/org/jmol/jvxl/readers/SurfaceReader.java
===================================================================
--- trunk/Jmol/src/org/jmol/jvxl/readers/SurfaceReader.java 2009-02-07 05:32:40 UTC (rev 10626)
+++ trunk/Jmol/src/org/jmol/jvxl/readers/SurfaceReader.java 2009-02-08 14:54:57 UTC (rev 10627)
@@ -183,6 +183,7 @@
protected MeshData meshData;
protected JvxlData jvxlData;
protected VolumeData volumeData;
+ private String edgeData;
protected boolean isProgressive = false;
protected boolean isXLowToHigh = false; //can be overridden in some readers by --progressive
@@ -226,6 +227,7 @@
voxelCounts = v.voxelCounts;
voxelData = v.voxelData;
volumeData = v;
+
/* if (mustCalcPoint)
v.setDataSource(this);
*/ }
@@ -292,7 +294,7 @@
jvxlData.pointsPerAngstrom = 1f/volumeData.volumetricVectorLengths[0];
jvxlData.jvxlColorData = "";
jvxlData.jvxlPlane = params.thePlane;
- jvxlData.jvxlEdgeData = fractionData.toString();
+ jvxlData.jvxlEdgeData = edgeData;
jvxlData.isBicolorMap = params.isBicolorMap;
jvxlData.isContoured = params.isContoured;
jvxlData.nContours = (params.contourFromZero
@@ -390,22 +392,26 @@
protected MarchingSquares marchingSquares;
private MarchingCubes marchingCubes;
+ public float getValue(int x, int y, int z) {
+ return volumeData.voxelData[x][y][z];
+ }
+
private void generateSurfaceData() {
- fractionData = new StringBuffer();
+ edgeData = "";
if (vertexDataOnly) {
try {
readSurfaceData(false);
} catch (Exception e) {
e.printStackTrace();
- Logger.error("Exception in SurfaceReader::readSurfaceData: " + e.getMessage());
+ Logger.error("Exception in SurfaceReader::readSurfaceData: "
+ + e.getMessage());
}
return;
}
contourVertexCount = 0;
int contourType = -1;
marchingSquares = null;
- if (params.isSquared)
- volumeData.filterData(params.isSquared, Float.NaN);
+
if (params.thePlane != null || params.isContoured) {
marchingSquares = new MarchingSquares(this, volumeData, params.thePlane,
params.nContours, params.thisContour, params.contourFromZero);
@@ -413,34 +419,27 @@
marchingSquares.setMinMax(params.valueMappedToRed,
params.valueMappedToBlue);
}
-
- marchingCubes = new MarchingCubes(this, volumeData, jvxlVoxelBitSet, params.isContoured,
- contourType, params.cutoff, params.isCutoffAbsolute);
-
- edgeCount = marchingCubes.generateSurfaceData(isXLowToHigh);
-
+ marchingCubes = new MarchingCubes(this, volumeData, jvxlVoxelBitSet,
+ params.isContoured, contourType, params.cutoff,
+ params.isCutoffAbsolute, params.isSquared, isXLowToHigh);
+ edgeData = marchingCubes.getEdgeData();
+ if (volumeData.voxelData == null)
+ JvxlReader.setSurfaceInfoFromBitSet(jvxlData,
+ marchingCubes.getBsVoxels(), params.thePlane);
if (isJvxl)
- fractionData = new StringBuffer(jvxlEdgeDataRead);
- fractionData.append('\n');
+ edgeData = jvxlEdgeDataRead;
}
- protected static boolean isInside(float voxelValue, float max,
- boolean isAbsolute) {
- return MarchingCubes.isInside(voxelValue, max, isAbsolute);
- }
-
///////////////// MarchingReader Interface Methods ///////////////////
protected final Point3f ptTemp = new Point3f();
- final float[] fReturn = new float[1];
-
- public int getSurfacePointIndex(float cutoff, boolean isCutoffAbsolute,
+ public int getSurfacePointIndexAndFraction(float cutoff, boolean isCutoffAbsolute,
int x, int y, int z, Point3i offset, int vA,
int vB, float valueA, float valueB,
Point3f pointA, Vector3f edgeVector,
- boolean isContourType) {
- float thisValue = readSurfacePoint(cutoff, isCutoffAbsolute, valueA,
+ boolean isContourType, float[] fReturn) {
+ float thisValue = getSurfacePointAndFraction(cutoff, isCutoffAbsolute, valueA,
valueB, pointA, edgeVector, fReturn, ptTemp);
/*
* from MarchingCubes
@@ -460,7 +459,7 @@
: MarchingSquares.CONTOUR_POINT);
if (assocVertex >= 0)
assocVertex = marchingCubes.getLinearOffset(x, y, z, assocVertex);
- int iV = addVertexCopy(ptTemp, thisValue, assocVertex);
+ int n = addVertexCopy(ptTemp, thisValue, assocVertex);
if (params.iAddGridPoints) {
marchingCubes.calcVertexPoint(x, y, z, vB, ptTemp);
addVertexCopy(valueA < valueB ? pointA : ptTemp, Float.NaN,
@@ -468,17 +467,16 @@
addVertexCopy(valueA < valueB ? ptTemp : pointA, Float.NaN,
MarchingSquares.EDGE_POINT);
}
- return iV;
+ return n;
}
- protected float readSurfacePoint(float cutoff, boolean isCutoffAbsolute,
+ protected float getSurfacePointAndFraction(float cutoff, boolean isCutoffAbsolute,
float valueA, float valueB, Point3f pointA,
- Vector3f edgeVector,
- float[] fReturn, Point3f ptReturn) {
+ Vector3f edgeVector, float[] fReturn,
+ Point3f ptReturn) {
//JvxlReader may or may not call this
-
float diff = valueB - valueA;
float fraction = (cutoff - valueA) / diff;
if (isCutoffAbsolute && (fraction < 0 || fraction > 1))
@@ -490,10 +488,6 @@
fraction = Float.NaN;
}
fReturn[0] = fraction;
- if (!isJvxl)
- fractionData.append(JvxlReader.jvxlFractionAsCharacter(fraction,
- edgeFractionBase, edgeFractionRange));
-
ptReturn.scaleAdd(fraction, edgeVector, pointA);
//System.out.println("SurfaceReader " + ptReturn + " " + (valueA + fraction * diff));
return valueA + fraction * diff;
@@ -793,7 +787,4 @@
meshDataServer.fillMeshData(meshData, MeshData.MODE_PUT_SETS);
}
- public void getCalcPoint(Point3f pt) {
- // for VertexDataServer - isoShapeReader only
- }
}
Modified: trunk/Jmol/src/org/jmol/jvxl/readers/VolumeDataReader.java
===================================================================
--- trunk/Jmol/src/org/jmol/jvxl/readers/VolumeDataReader.java 2009-02-07 05:32:40 UTC (rev 10626)
+++ trunk/Jmol/src/org/jmol/jvxl/readers/VolumeDataReader.java 2009-02-08 14:54:57 UTC (rev 10627)
@@ -104,46 +104,22 @@
protected void readVoxelDataIndividually(boolean isMapData) throws Exception {
if (isMapData && !allowMapData)
return; //not applicable
- boolean inside = false;
- int dataCount = 0;
- voxelData = new float[nPointsX][nPointsY][nPointsZ];
- nDataPoints = 0;
- int nSurfaceInts = 0;
- StringBuffer sb = new StringBuffer();
- float cutoff = params.cutoff;
- boolean isCutoffAbsolute = params.isCutoffAbsolute;
+ voxelData = (isMapData ? new float[nPointsX][nPointsY][nPointsZ] : null);
+ volumeData.setVoxelData(voxelData);
+ if (!isMapData)
+ return;
for (int x = 0; x < nPointsX; ++x) {
float[][] plane = new float[nPointsY][];
voxelData[x] = plane;
for (int y = 0; y < nPointsY; ++y) {
float[] strip = plane[y] = new float[nPointsZ];
for (int z = 0; z < nPointsZ; ++z) {
- float voxelValue = strip[z] = getValue(x, y, z);
- ++nDataPoints;
- if (inside == isInside(voxelValue, cutoff, isCutoffAbsolute)) {
- dataCount++;
- } else {
- if (!isMapData)
- sb.append(' ').append(dataCount);
- ++nSurfaceInts;
- dataCount = 1;
- inside = !inside;
- }
+ strip[z] = getValue(x, y, z);
}
}
}
- //Jvxl getNextVoxelValue records the data read on its own.
- if (!isMapData) {
- sb.append(' ').append(dataCount).append('\n');
- JvxlReader.setSurfaceInfo(jvxlData, params.thePlane, nSurfaceInts, sb);
- }
- volumeData.setVoxelData(voxelData);
}
- protected float getValue(int x, int y, int z) {
- return 0;
- }
-
protected int setVoxelRange(int index, float min, float max, float ptsPerAngstrom,
int gridMax) {
if (min >= max) {
@@ -196,14 +172,10 @@
protected void readSurfaceData(boolean isMapData) throws Exception {
//precalculated -- just creating the JVXL equivalent
- if (!precalculateVoxelData) {
+ if (precalculateVoxelData)
+ generateCube();
+ else
readVoxelDataIndividually(isMapData);
- return;
- }
- generateCube();
- if (isMapData)
- return;
- nDataPoints = JvxlReader.jvxlCreateSurfaceData(jvxlData, volumeData.voxelData, params.cutoff, params.isCutoffAbsolute, nPointsX, nPointsY, nPointsZ);
}
protected void generateCube() {
Modified: trunk/Jmol/src/org/jmol/jvxl/readers/VolumeFileReader.java
===================================================================
--- trunk/Jmol/src/org/jmol/jvxl/readers/VolumeFileReader.java 2009-02-07 05:32:40 UTC (rev 10626)
+++ trunk/Jmol/src/org/jmol/jvxl/readers/VolumeFileReader.java 2009-02-08 14:54:57 UTC (rev 10627)
@@ -176,10 +176,8 @@
*/
next[0] = 0;
- boolean inside = false;
int downsampleFactor = params.downsampleFactor;
boolean isDownsampled = canDownsample && (downsampleFactor > 0);
- int dataCount = 0;
if (params.thePlane != null) {
params.cutoff = 0f;
} else if (isJvxl) {
@@ -187,15 +185,13 @@
}
nDataPoints = 0;
line = "";
- StringBuffer sb = new StringBuffer();
jvxlNSurfaceInts = 0;
if (isJvxl) {
nDataPoints = volumeData.setVoxelCounts(nPointsX, nPointsY, nPointsZ);
- jvxlVoxelBitSet = getVoxelBitSet(nDataPoints, sb);
+ jvxlVoxelBitSet = getVoxelBitSet(nDataPoints);
voxelData = null;
} else {
voxelData = new float[nPointsX][][];
- boolean collectData = (!isJvxl && params.thePlane == null);
int nSkipX = 0;
int nSkipY = 0;
int nSkipZ = 0;
@@ -211,83 +207,40 @@
//Note downsampling not allowed for JVXL files
- if (isMapData) {
- for (int x = 0; x < nPointsX; ++x) {
- float[][] plane = new float[nPointsY][];
- voxelData[x] = plane;
- for (int y = 0; y < nPointsY; ++y) {
- float[] strip = new float[nPointsZ];
- plane[y] = strip;
- for (int z = 0; z < nPointsZ; ++z) {
- strip[z] = getNextVoxelValue(sb);
- ++nDataPoints;
- if (isDownsampled)
- skipVoxels(nSkipX);
- }
+ for (int x = 0; x < nPointsX; ++x) {
+ float[][] plane = new float[nPointsY][];
+ voxelData[x] = plane;
+ for (int y = 0; y < nPointsY; ++y) {
+ float[] strip = new float[nPointsZ];
+ plane[y] = strip;
+ for (int z = 0; z < nPointsZ; ++z) {
+ strip[z] = getNextVoxelValue();
if (isDownsampled)
- skipVoxels(nSkipY);
+ skipVoxels(nSkipX);
}
if (isDownsampled)
- skipVoxels(nSkipZ);
+ skipVoxels(nSkipY);
}
- } else {
- float cutoff = params.cutoff;
- boolean isCutoffAbsolute = params.isCutoffAbsolute;
- for (int x = 0; x < nPointsX; ++x) {
- float[][] plane;
- plane = new float[nPointsY][];
- voxelData[x] = plane;
- for (int y = 0; y < nPointsY; ++y) {
- float[] strip = new float[nPointsZ];
- plane[y] = strip;
- for (int z = 0; z < nPointsZ; ++z) {
- float voxelValue = getNextVoxelValue(sb);
- strip[z] = voxelValue;
- ++nDataPoints;
- if (inside == isInside(voxelValue, cutoff, isCutoffAbsolute)) {
- dataCount++;
- } else {
- if (collectData && dataCount != 0) {
- sb.append(' ').append(dataCount);
- ++jvxlNSurfaceInts;
- }
- dataCount = 1;
- inside = !inside;
- }
- if (isDownsampled)
- skipVoxels(nSkipX);
- }
- if (isDownsampled)
- skipVoxels(nSkipY);
- }
- if (isDownsampled)
- skipVoxels(nSkipZ);
- }
+ if (isDownsampled)
+ skipVoxels(nSkipZ);
}
//Jvxl getNextVoxelValue records the data read on its own.
- if (collectData) {
- sb.append(' ').append(dataCount).append('\n');
- ++jvxlNSurfaceInts;
- }
}
- if (!isMapData)
- JvxlReader
- .setSurfaceInfo(jvxlData, params.thePlane, jvxlNSurfaceInts, sb);
volumeData.setVoxelData(voxelData);
}
private void skipVoxels(int n) throws Exception {
// not allowed for JVXL data
for (int i = n; --i >= 0;)
- getNextVoxelValue(null);
+ getNextVoxelValue();
}
- protected BitSet getVoxelBitSet(int nPoints, StringBuffer sb) throws Exception {
+ protected BitSet getVoxelBitSet(int nPoints) throws Exception {
// jvxlReader will use this to read the surface voxel data
return null;
}
- protected float getNextVoxelValue(StringBuffer sb) throws Exception {
+ protected float getNextVoxelValue() throws Exception {
//overloaded in JvxlReader, where sb is appended to
float voxelValue = 0;
if (nSurfaces > 1 && !params.blockCubeData) {
Modified: trunk/Jmol/src/org/jmol/shapesurface/Isosurface.java
===================================================================
--- trunk/Jmol/src/org/jmol/shapesurface/Isosurface.java 2009-02-07 05:32:40 UTC (rev 10626)
+++ trunk/Jmol/src/org/jmol/shapesurface/Isosurface.java 2009-02-08 14:54:57 UTC (rev 10627)
@@ -733,11 +733,11 @@
///////////// VertexDataServer interface methods ////////////////
- public int getSurfacePointIndex(float cutoff, boolean isCutoffAbsolute,
+ public int getSurfacePointIndexAndFraction(float cutoff, boolean isCutoffAbsolute,
int x, int y, int z, Point3i offset, int vA,
int vB, float valueA, float valueB,
Point3f pointA, Vector3f edgeVector,
- boolean isContourType) {
+ boolean isContourType, float[] fReturn) {
return 0;
}
@@ -752,6 +752,7 @@
boolean isAbsolute) {
if (isAbsolute && !MeshData.checkCutoff(iA, iB, iC, thisMesh.vertexValues))
return;
+ //System.out.println(" isosurface triangle check : " + iA + " " + iB + " " + iC);
thisMesh.addTriangleCheck(iA, iB, iC, check);
}
@@ -844,4 +845,8 @@
String colors = viewer.getColorSchemeList(schemeName, false);
return "\"" + (colors.length() == 0 ? schemeName : colors) + "\"";
}
+
+ public float getValue(int x, int y, int z) {
+ return 0;
+ }
}
Modified: trunk/Jmol/src/org/jmol/viewer/Jmol.properties
===================================================================
--- trunk/Jmol/src/org/jmol/viewer/Jmol.properties 2009-02-07 05:32:40 UTC (rev 10626)
+++ trunk/Jmol/src/org/jmol/viewer/Jmol.properties 2009-02-08 14:54:57 UTC (rev 10627)
@@ -3,6 +3,8 @@
version=11.7.25_dev
+# code: rewritten code for ASimpleJvxlWriter allows streaming data
+# code: rewritten code for Marching Cubes no longer requires 3D array of data
# new feature: CTRL-K toggles keystrokes with display; ALT-K or CTRL_ALT_K same without display
# color of frank displays condition -- blue (unsigned) or orange (signed) indicates keystrokes enabled
# bug fix: PDB load filter broken for ANISOU
Modified: trunk/Jmol/src/org/openscience/jvxl/simplewriter/ASimpleJvxlWriter.java
===================================================================
--- trunk/Jmol/src/org/openscience/jvxl/simplewriter/ASimpleJvxlWriter.java 2009-02-07 05:32:40 UTC (rev 10626)
+++ trunk/Jmol/src/org/openscience/jvxl/simplewriter/ASimpleJvxlWriter.java 2009-02-08 14:54:57 UTC (rev 10627)
@@ -37,32 +37,61 @@
public static void main(String[] args) {
- VolumeData volumeData = new VolumeData();
// parameters that need setting:
String outputFile = "c:/temp/simple.jvxl";
- float cutoff = 0.0f;
+ float cutoff = 0.01f;
boolean isCutoffAbsolute = false;
- int nX = 10;
- int nY = 10;
- int nZ = 10;
+ int nX = 30;
+ int nY = 30;
+ int nZ = 30;
+ String title = "created by SimpleJvxlWriter "
+ + new SimpleDateFormat("yyyy-MM-dd', 'HH:mm").format(new Date())
+ + "\naddional comment line\n";
+
+ VolumeData volumeData;
+ VoxelDataCreator vdc;
+ JvxlData jvxlData;
+
+ volumeData = new VolumeData();
volumeData.setVolumetricOrigin(0, 0, 0);
volumeData.setVolumetricVector(0, 1f, 0f, 0f);
volumeData.setVolumetricVector(1, 0f, 1f, 0f);
volumeData.setVolumetricVector(2, 0f, 0f, 1f);
+ volumeData.setUnitVectors();
+
volumeData.setVoxelCounts(nX, nY, nZ);
volumeData.setVoxelData(new float[nX][nY][nZ]);
- VoxelDataCreator vdc = new VoxelDataCreator();
- vdc.createVoxelData(volumeData);
- JvxlData jvxlData = new JvxlData();
+ vdc = new VoxelDataCreator(volumeData);
+ vdc.createVoxelData();
+ jvxlData = new JvxlData();
jvxlData.cutoff = cutoff;
jvxlData.isCutoffAbsolute = isCutoffAbsolute;
- StringBuffer sb = new StringBuffer("created by SimpleJvxlWriter "
- + new SimpleDateFormat("yyyy-MM-dd', 'HH:mm").format(new Date())
- + "\naddional comment line\n");
- writeFile(outputFile, JvxlWrite.jvxlGetData(jvxlData, volumeData, sb));
+ jvxlData.isXLowToHigh = true;
+ writeFile(outputFile + jvxlData.isXLowToHigh , JvxlWrite.jvxlGetData(null, jvxlData, volumeData, title));
+
+/*
+ volumeData = new VolumeData();
+ volumeData.setVolumetricOrigin(2, 2, 2);
+ volumeData.setVolumetricVector(0, 1f, 0f, 0f);
+ volumeData.setVolumetricVector(1, 0f, 1f, 0f);
+ volumeData.setVolumetricVector(2, 0f, 0f, 1f);
+ volumeData.setUnitVectors();
+
+
+ volumeData.setVoxelCounts(nX, nY, nZ);
+ volumeData.setVoxelData(new float[nX][nY][nZ]);
+ vdc = new VoxelDataCreator(volumeData);
+ vdc.createVoxelData();
+ jvxlData = new JvxlData();
+ jvxlData.cutoff = cutoff;
+ jvxlData.isCutoffAbsolute = isCutoffAbsolute;
+ jvxlData.isXLowToHigh = true;
+ writeFile(outputFile+ "B", JvxlWrite.jvxlGetData(vdc, jvxlData, volumeData, title));
+
+ */
System.out.flush();
System.exit(0);
}
Modified: trunk/Jmol/src/org/openscience/jvxl/simplewriter/JvxlWrite.java
===================================================================
--- trunk/Jmol/src/org/openscience/jvxl/simplewriter/JvxlWrite.java 2009-02-07 05:32:40 UTC (rev 10626)
+++ trunk/Jmol/src/org/openscience/jvxl/simplewriter/JvxlWrite.java 2009-02-08 14:54:57 UTC (rev 10627)
@@ -23,7 +23,10 @@
*/
package org.openscience.jvxl.simplewriter;
+import java.util.BitSet;
+
import javax.vecmath.Point3f;
+import javax.vecmath.Point4f;
public class JvxlWrite {
@@ -40,32 +43,84 @@
public JvxlWrite() {
}
- public static String jvxlGetData(JvxlData jvxlData, VolumeData volumeData, StringBuffer sb) {
+ public static String jvxlGetData(VoxelDataCreator vdc, JvxlData jvxlData, VolumeData volumeData, String title) {
// if the StringBuffer is not empty, it should have two comment lines
// that do not start with # already present.
+ StringBuffer sb = new StringBuffer();
+ if (title != null)
+ sb.append(title);
Point3f[] atomXYZ = null;
int[] atomNo = null;
int nAtoms = Integer.MAX_VALUE;
- jvxlCreateHeader(volumeData, nAtoms, atomXYZ, atomNo, sb);
+ jvxlCreateHeader(volumeData, nAtoms, atomXYZ, atomNo, jvxlData.isXLowToHigh, sb);
jvxlData.jvxlFileHeader = sb.toString();
int[] counts = volumeData.getVoxelCounts();
- jvxlCreateSurfaceData(jvxlData, volumeData.getVoxelData(), jvxlData.cutoff,
- jvxlData.isCutoffAbsolute, counts);
- SimpleMarchingCubes mc = new SimpleMarchingCubes(volumeData,
- jvxlData.cutoff, jvxlData.isCutoffAbsolute);
+ jvxlData.nPointsX = counts[0];
+ jvxlData.nPointsY = counts[1];
+ jvxlData.nPointsZ = counts[2];
+ SimpleMarchingCubes mc = new SimpleMarchingCubes(vdc, volumeData,
+ jvxlData.cutoff, jvxlData.isCutoffAbsolute, jvxlData.isXLowToHigh);
jvxlData.jvxlEdgeData = mc.getEdgeData();
- String jvxlMaskData = mc.getMaskData();
- System.out.println(jvxlMaskData.length() + " "
- + (jvxlMaskData = jvxlCompressString(jvxlMaskData)).length() + " " + jvxlData.jvxlSurfaceData.length());
+ setSurfaceInfoFromBitSet(jvxlData, mc.getBsVoxels(), null);
jvxlData.jvxlDefinitionLine = jvxlGetDefinitionLine(jvxlData);
return jvxlGetFile(jvxlData);
}
+ public static void setSurfaceInfoFromBitSet(JvxlData jvxlData, BitSet bs,
+ Point4f thePlane) {
+ boolean inside = false;
+ int dataCount = 0;
+ StringBuffer sb = new StringBuffer();
+ int nSurfaceInts = 0;
+ int nPoints = jvxlData.nPointsX * jvxlData.nPointsY * jvxlData.nPointsZ;
+ for (int i = 0; i < nPoints; ++i) {
+ if (inside == bs.get(i)) {
+ dataCount++;
+ } else {
+ sb.append(' ').append(dataCount);
+ nSurfaceInts++;
+ dataCount = 1;
+ inside = !inside;
+ }
+ }
+ sb.append(' ').append(dataCount).append('\n');
+ setSurfaceInfo(jvxlData, thePlane, nSurfaceInts, sb);
+ }
+
static char jvxlFractionAsCharacter(float fraction) {
+ //char ch = jvxlFractionAsCharacter(fraction, defaultEdgeFractionBase,
+ // defaultEdgeFractionRange);
+ //System.out.println(fraction + " " + ch + " " + jvxlFractionFromCharacter((int) ch, defaultEdgeFractionBase,
+ //defaultEdgeFractionRange,0));
return jvxlFractionAsCharacter(fraction, defaultEdgeFractionBase,
defaultEdgeFractionRange);
}
+ protected static float jvxlFractionFromCharacter(int ich, int base, int range,
+ float fracOffset) {
+ if (ich == base + range)
+ return Float.NaN;
+ if (ich < base)
+ ich = 92; // ! --> \
+ float fraction = (ich - base + fracOffset) / range;
+ if (fraction < 0f)
+ return 0f;
+ if (fraction > 1f)
+ return 0.999999f;
+ //if (logCompression)
+ //Logger.info("ffc: " + fraction + " <-- " + ich + " " + (char) ich);
+ return fraction;
+ }
+
+
+ protected static void setSurfaceInfo(JvxlData jvxlData, Point4f thePlane, int nSurfaceInts, StringBuffer surfaceData) {
+ jvxlData.jvxlSurfaceData = surfaceData.toString();
+ if (jvxlData.jvxlSurfaceData.indexOf("--") == 0)
+ jvxlData.jvxlSurfaceData = jvxlData.jvxlSurfaceData.substring(2);
+ jvxlData.jvxlPlane = thePlane;
+ jvxlData.nSurfaceInts = nSurfaceInts;
+ }
+
private static String jvxlGetFile(JvxlData jvxlData) {
StringBuffer data = new StringBuffer();
String s = jvxlData.jvxlFileHeader + jvxlExtraLine(jvxlData, 1);
@@ -102,23 +157,20 @@
// ascii-encoded fractional color data
// # optional comments
- private static void setSurfaceInfo(JvxlData jvxlData, int nSurfaceInts, StringBuffer surfaceData) {
- jvxlData.jvxlSurfaceData = surfaceData.toString();
- if (jvxlData.jvxlSurfaceData.indexOf("--") == 0)
- jvxlData.jvxlSurfaceData = jvxlData.jvxlSurfaceData.substring(2);
- jvxlData.nSurfaceInts = nSurfaceInts;
- }
//// methods for creating the JVXL code
private static void jvxlCreateHeader(VolumeData v, int nAtoms,
Point3f[] atomXyz, int[] atomNo,
+ boolean isXLowToHigh,
StringBuffer sb) {
// if the StringBuffer comes in non-empty, it should have two lines
// that do not start with # already present.
if (sb.length() == 0)
sb.append("Line 1\nLine 2\n");
- sb.append(nAtoms == Integer.MAX_VALUE ? -2 : -nAtoms).append(' ').append(
+ sb.append(isXLowToHigh ? "+" : "-");
+ sb.append(nAtoms == Integer.MAX_VALUE ? 2 : Math.abs(nAtoms));
+ sb.append(' ').append(
v.volumetricOrigin.x).append(' ').append(v.volumetricOrigin.y).append(
' ').append(v.volumetricOrigin.z).append(" ANGSTROMS\n");
for (int i = 0; i < 3; i++)
@@ -145,44 +197,6 @@
.append(pt.z).append(" //BOGUS He ATOM ADDED FOR JVXL FORMAT\n");
}
- private static int jvxlCreateSurfaceData(JvxlData jvxlData,
- float[][][] voxelData, float cutoff,
- boolean isCutoffAbsolute,
- int counts[]) {
- StringBuffer sb = new StringBuffer();
- boolean inside = false;
- int dataCount = 0;
- int nDataPoints = 0;
- int nSurfaceInts = 0;
- int nX = counts[0];
- int nY = counts[1];
- int nZ = counts[2];
- for (int x = 0; x < nX; ++x)
- for (int y = 0; y < nY; ++y)
- for (int z = 0; z < nZ; ++z) {
- ++nDataPoints;
- if (inside == isInside(voxelData[x][y][z], cutoff, isCutoffAbsolute)) {
- dataCount++;
- } else {
- if (dataCount != 0) {
- sb.append(' ').append(dataCount);
- ++nSurfaceInts;
- }
- dataCount = 1;
- inside = !inside;
- }
- }
- sb.append(' ').append(dataCount).append('\n');
- ++nSurfaceInts;
- setSurfaceInfo(jvxlData, nSurfaceInts, sb);
- return nDataPoints;
- }
-
- private static boolean isInside(float voxelValue, float max,
- boolean isAbsolute) {
- return SimpleMarchingCubes.isInside(voxelValue, max, isAbsolute);
- }
-
private static String jvxlGetDefinitionLine(JvxlData jvxlData) {
String definitionLine = jvxlData.cutoff + " ";
Modified: trunk/Jmol/src/org/openscience/jvxl/simplewriter/SimpleMarchingCubes.java
===================================================================
--- trunk/Jmol/src/org/openscience/jvxl/simplewriter/SimpleMarchingCubes.java 2009-02-07 05:32:40 UTC (rev 10626)
+++ trunk/Jmol/src/org/openscience/jvxl/simplewriter/SimpleMarchingCubes.java 2009-02-08 14:54:57 UTC (rev 10627)
@@ -23,10 +23,10 @@
*/
package org.openscience.jvxl.simplewriter;
+import java.util.BitSet;
+
import javax.vecmath.Point3i;
-import org.jmol.util.Base64;
-
public class SimpleMarchingCubes {
/*
@@ -47,16 +47,25 @@
private boolean isCutoffAbsolute;
private boolean isXLowToHigh;
private StringBuffer fractionData = new StringBuffer();
- private StringBuffer maskData = new StringBuffer();
-
- public String getMaskData() {
- return Base64.getBase64(maskData).toString();
+
+ private int cubeCountX, cubeCountY, cubeCountZ;
+ private int nY, nZ;
+
+ private BitSet bsVoxels = new BitSet();
+
+ public BitSet getBsVoxels() {
+ return bsVoxels;
}
- private int cubeCountX, cubeCountY, cubeCountZ;
+ private int mode;
+ private final static int MODE_CUBE = 1;
+ private final static int MODE_BITSET = 2;
+ private final static int MODE_GETXYZ = 3;
- public SimpleMarchingCubes(VolumeData volumeData, float cutoff,
- boolean isCutoffAbsolute) {
+ private VoxelDataCreator vdc;
+
+ public SimpleMarchingCubes(VoxelDataCreator vdc, VolumeData volumeData, float cutoff,
+ boolean isCutoffAbsolute , boolean isXLowToHigh) {
// when just creating a JVXL file all you really need are:
//
@@ -64,15 +73,23 @@
// cutoff
//
+ this.vdc = vdc;
this.volumeData = volumeData;
this.cutoff = cutoff;
this.isCutoffAbsolute = isCutoffAbsolute;
+ this.isXLowToHigh = isXLowToHigh;
+
+ if (vdc == null) {
+ mode = MODE_CUBE;
+ } else {
+ mode = MODE_GETXYZ;
+ }
- cubeCountX = volumeData.voxelData.length - 1;
- cubeCountY = volumeData.voxelData[0].length - 1;
- cubeCountZ = volumeData.voxelData[0][0].length - 1;
- xyCount = (cubeCountX + 1) * (cubeCountY + 1);
-
+ cubeCountX = volumeData.voxelCounts[0] - 1;
+ cubeCountY = (nY = volumeData.voxelCounts[1]) - 1;
+ cubeCountZ = (nZ = volumeData.voxelCounts[2]) - 1;
+ yzCount = nY * nZ;
+ setLinearOffsets();
}
private final float[] vertexValues = new float[8];
@@ -81,7 +98,7 @@
for (int i = 8; --i >= 0;)
vertexPoints[i] = new Point3i();
}
- int xyCount;
+
int edgeCount;
/* Note to Jason from Bob:
@@ -110,36 +127,54 @@
*
*/
+ private static int[] xyPlanePts = new int[] { 0, 1, 1, 0, 0, 1, 1, 0 };
+
public String getEdgeData() {
// set up the set of edge points in the YZ plane
- // these are indexes into an array of Point3f values
+ // isoPointIndixes are indices into an array of Point3f values
// They will be initialized as -1 whenever a vertex is needed.
// But if just creating a JVXL file, all you need to do
// is set them to 0, not an index into any actual array.
int[][] isoPointIndexes = new int[cubeCountY * cubeCountZ][12];
- edgeCount = 0;
+ float[][] xyPlanes = (mode == MODE_GETXYZ ? new float[2][yzCount] : null);
- int x0, x1, xStep;
+ int x0, x1, xStep, ptStep, pt, ptX;
if (isXLowToHigh) {
x0 = 0;
x1 = cubeCountX;
xStep = 1;
+ ptStep = yzCount;
+ pt = ptX = (yzCount - 1) - nZ - 1;
+ // we are starting at the top corner, in the next to last
+ // cell on the next to last row of the first plane
} else {
x0 = cubeCountX - 1;
x1 = -1;
xStep = -1;
+ ptStep = -yzCount;
+ pt = ptX = (cubeCountX * yzCount - 1) - nZ - 1;
+ // we are starting at the top corner, in the next to last
+ // cell on the next to last row of the next to last plane(!)
}
- for (int x = x0; x != x1; x += xStep) {
- for (int y = cubeCountY; --y >= 0;) {
- for (int z = cubeCountZ; --z >= 0;) {
+ int cellIndex0 = cubeCountY * cubeCountZ - 1;
+ int cellIndex = cellIndex0;
+ for (int x = x0; x != x1; x += xStep, ptX += ptStep, pt = ptX, cellIndex = cellIndex0) {
+ if (mode == MODE_GETXYZ) {
+ float[] plane = xyPlanes[0];
+ xyPlanes[0] = xyPlanes[1];
+ xyPlanes[1] = plane;
+ }
+ for (int y = cubeCountY; --y >= 0; pt--) {
+ for (int z = cubeCountZ; --z >= 0; pt--, cellIndex--) {
+
// set up the list of indices that need checking
- int[] voxelPointIndexes = propagateNeighborPointIndexes(x, y, z,
- isoPointIndexes);
+ int[] voxelPointIndexes = propagateNeighborPointIndexes(x, y, z, pt,
+ isoPointIndexes, cellIndex);
// create the bitset mask indicating which vertices are inside.
// 0xFF here means "all inside"; 0x00 means "all outside"
@@ -150,14 +185,32 @@
// cubeVertexOffsets just gets us the specific grid point relative
// to our base x,y,z cube position
+ boolean isInside;
Point3i offset = cubeVertexOffsets[i];
- if (isInside(
- (vertexValues[i] = volumeData.voxelData[x + offset.x][y
- + offset.y][z + offset.z]), cutoff, isCutoffAbsolute))
+ int pti = pt + linearOffsets[i];
+ switch (mode) {
+ case MODE_GETXYZ:
+ vertexValues[i] = getValue(i, x + offset.x, y + offset.y, z
+ + offset.z, pti, xyPlanes[xyPlanePts[i]]);
+ isInside = bsVoxels.get(pti);
+ break;
+ case MODE_BITSET:
+ isInside = bsVoxels.get(pti);
+ vertexValues[i] = (isInside ? 1 : 0);
+ break;
+ default:
+ case MODE_CUBE:
+ vertexValues[i] = volumeData.voxelData[x + offset.x][y + offset.y][z
+ + offset.z];
+ isInside = isInside(vertexValues[i], cutoff, isCutoffAbsolute);
+ if (isInside)
+ bsVoxels.set(pti);
+ }
+ if (isInside) {
insideMask |= 1 << i;
+ }
}
- maskData.append((char)insideMask);
if (insideMask == 0) {
continue;
}
@@ -166,7 +219,7 @@
}
// This cube is straddling the cutoff. We must check all edges
- processOneCubical(insideMask, voxelPointIndexes, x, y, z);
+ processOneCubical(insideMask, voxelPointIndexes, x, y, z, pt);
}
}
}
@@ -177,31 +230,28 @@
return ((max > 0 && (isAbsolute ? Math.abs(voxelValue) : voxelValue) >= max) || (max <= 0 && voxelValue <= max));
}
+ BitSet bsValues = new BitSet();
+
+ private float getValue(int i, int x, int y, int z, int pt, float[] tempValues) {
+ //if (bsValues.get(pt))
+ //return tempValues[pt % yzCount];
+ bsValues.set(pt);
+ float value = vdc.getValue(x, y, z);
+ tempValues[pt % yzCount] = value;
+ //System.out.println("xyz " + x + " " + y + " " + z + " v=" + value);
+ if (isInside(value, cutoff, isCutoffAbsolute))
+ bsVoxels.set(pt);
+ return value;
+ }
+
private final int[] nullNeighbor = { -1, -1, -1, -1, -1, -1, -1, -1, -1, -1,
-1, -1 };
- private int[] propagateNeighborPointIndexes(int x, int y, int z,
- int[][] isoPointIndexes) {
- /* Y
- * 4 --------4--------- 5
- * /| /|
- * / | / |
- * / | / |
- * 7 8 5 |
- * / | / 9
- * / | / |
- * 7 --------6--------- 6 |
- * | | | |
- * | 0 ---------0--|----- 1 X
- * | / | /
- * 11 / 10 /
- * | 3 | 1
- * | / | /
- * | / | /
- * 3 ---------2-------- 2
- * Z
+ private int[] propagateNeighborPointIndexes(int x, int y, int z, int pt,
+ int[][] isoPointIndexes,
+ int cellIndex) {
+ /*
*
- *
* We are running through the grid points in yz planes from high x --> low x
* and within those planes along strips from high y to low y
* and within those strips, from high z to low z. The "leading vertex" is 0,
@@ -237,102 +287,100 @@
* All we are really talking about is the JVXL reader, because we can certainly
* switch to progressive mode in all the other readers.
*
+ * Y
+ * 4 --------4--------- 5
+ * /| /|
+ * / | / |
+ * / | / |
+ * 7 8 5 |
+ * / | / 9
+ * / | / |
+ * 7 --------6--------- 6 |
+ * | | | |
+ * | 0 ---------0--|----- 1 X
+ * | / | /
+ * 11 / 10 /
+ * | 3 | 1
+ * | / | /
+ * | / | /
+ * 3 ---------2-------- 2
+ * Z
+ *
*
*/
-
+
/* DO NOT EVER CHANGE THIS */
+
- int cellIndex = y * cubeCountZ + z;
int[] voxelPointIndexes = isoPointIndexes[cellIndex];
- boolean noXNeighbor = (x == cubeCountX - 1);
+ boolean noYNeighbor = (y == cubeCountY - 1);
+ int[] yNeighbor = noYNeighbor ? nullNeighbor
+ : isoPointIndexes[cellIndex + cubeCountZ];
+ boolean noZNeighbor = (z == cubeCountZ - 1);
+ int[] zNeighbor = noZNeighbor ? nullNeighbor
+ : isoPointIndexes[cellIndex + 1];
+ voxelPointIndexes[0] = -1;
+ voxelPointIndexes[2] = zNeighbor[0];
+ voxelPointIndexes[4] = yNeighbor[0];
+ voxelPointIndexes[6] = (noYNeighbor ? zNeighbor[4] : yNeighbor[2]);
+
if (isXLowToHigh) {
// reading x from low to high
- if (noXNeighbor) {
+ if (x == 0) {
voxelPointIndexes[3] = -1;
voxelPointIndexes[8] = -1;
- voxelPointIndexes[7] = -1;
- voxelPointIndexes[11] = -1;
+ voxelPointIndexes[7] = yNeighbor[3];
+ voxelPointIndexes[11] = zNeighbor[8];
} else {
voxelPointIndexes[3] = voxelPointIndexes[1];
voxelPointIndexes[7] = voxelPointIndexes[5];
voxelPointIndexes[8] = voxelPointIndexes[9];
voxelPointIndexes[11] = voxelPointIndexes[10];
}
+ voxelPointIndexes[1] = -1;
+ voxelPointIndexes[5] = yNeighbor[1];
+ voxelPointIndexes[9] = -1;
+ voxelPointIndexes[10] = zNeighbor[9];
} else {
// reading x from high to low
- if (noXNeighbor) {
- // the x neighbor is myself from my last pass through here
+ if (x == cubeCountX - 1) {
voxelPointIndexes[1] = -1;
+ voxelPointIndexes[5] = yNeighbor[1];
voxelPointIndexes[9] = -1;
- voxelPointIndexes[5] = -1;
- voxelPointIndexes[10] = -1;
+ voxelPointIndexes[10] = zNeighbor[9];
} else {
voxelPointIndexes[1] = voxelPointIndexes[3];
voxelPointIndexes[5] = voxelPointIndexes[7];
voxelPointIndexes[9] = voxelPointIndexes[8];
voxelPointIndexes[10] = voxelPointIndexes[11];
}
- }
- //from the y neighbor pick up the top
- boolean noYNeighbor = (y == cubeCountY - 1);
- int[] yNeighbor = noYNeighbor ? nullNeighbor : isoPointIndexes[cellIndex
- + cubeCountZ];
-
- voxelPointIndexes[4] = yNeighbor[0];
- voxelPointIndexes[6] = yNeighbor[2];
-
- if (isXLowToHigh) {
- voxelPointIndexes[5] = yNeighbor[1];
- if (noXNeighbor)
- voxelPointIndexes[7] = yNeighbor[3];
- } else {
+ voxelPointIndexes[3] = -1;
voxelPointIndexes[7] = yNeighbor[3];
- if (noXNeighbor)
- voxelPointIndexes[5] = yNeighbor[1];
+ voxelPointIndexes[8] = -1;
+ voxelPointIndexes[11] = zNeighbor[8];
}
- // from my z neighbor
- boolean noZNeighbor = (z == cubeCountZ - 1);
- int[] zNeighbor = noZNeighbor ? nullNeighbor
- : isoPointIndexes[cellIndex + 1];
- voxelPointIndexes[2] = zNeighbor[0];
- if (noYNeighbor)
- voxelPointIndexes[6] = zNeighbor[4];
- if (isXLowToHigh) {
- if (noXNeighbor)
- voxelPointIndexes[11] = zNeighbor[8];
- voxelPointIndexes[10] = zNeighbor[9];
- } else {
- if (noXNeighbor)
- voxelPointIndexes[10] = zNeighbor[9];
- voxelPointIndexes[11] = zNeighbor[8];
- }
- // these must always be calculated
- voxelPointIndexes[0] = -1;
- if (isXLowToHigh) {
- voxelPointIndexes[1] = -1;
- voxelPointIndexes[9] = -1;
- } else {
- voxelPointIndexes[3] = -1;
- voxelPointIndexes[8] = -1;
- }
return voxelPointIndexes;
}
+ private static final int[] Pwr2 = new int[] { 1, 2, 4, 8, 16, 32, 64, 128,
+ 256, 512, 1024, 2048 };
+
private boolean processOneCubical(int insideMask, int[] voxelPointIndexes,
- int x, int y, int z) {
+ int x, int y, int z, int pt) {
// the key to the algorithm is that we have a catalog that
// maps the inside-vertex mask to an edge mask.
int edgeMask = insideMaskTable[insideMask];
+ //for (int i =0; i < 8; i++) System.out.print("\nvpi for cell " + pt + ": vertex " + i + ": " + voxelPointIndexes[i] + " " + Integer.toBinaryString(edgeMask));
boolean isNaN = false;
for (int iEdge = 12; --iEdge >= 0;) {
// bit set to one means it's a relevant edge
- if ((edgeMask & (1 << iEdge)) == 0)
+ if ((edgeMask & Pwr2[iEdge]) == 0)
continue;
// if we have a point already, we don't need to check this edge.
@@ -344,12 +392,10 @@
// here's an edge that has to be checked.
- ++edgeCount;
-
// get the vertex numbers 0 - 7
- int vertexA = edgeVertexes[2 * iEdge];
- int vertexB = edgeVertexes[2 * iEdge + 1];
+ int vertexA = edgeVertexes[iEdge << 1];
+ int vertexB = edgeVertexes[(iEdge << 1) + 1];
// pick up the actual value at each vertex
// this array of 8 values is updated as we go.
@@ -370,19 +416,46 @@
// here is where we get the value and assign the point for that edge
// it is where the JVXL surface data line is appended
- voxelPointIndexes[iEdge] = 0;
+ voxelPointIndexes[iEdge] = edgeCount++;
+ //System.out.println(" pt=" + pt + " edge" + iEdge + " xyz " + x + " " + y + " " + z + " vertexAB=" + vertexA + " " + vertexB + " valueAB=" + valueA + " " + valueB + " f= " + (cutoff - valueA) / (valueB - valueA));
fractionData.append(JvxlWrite.jvxlFractionAsCharacter((cutoff - valueA) / (valueB - valueA)));
}
return !isNaN;
}
- final static Point3i[] cubeVertexOffsets = { new Point3i(0, 0, 0),
- new Point3i(1, 0, 0), new Point3i(1, 0, 1), new Point3i(0, 0, 1),
- new Point3i(0, 1, 0), new Point3i(1, 1, 0), new Point3i(1, 1, 1),
- new Point3i(0, 1, 1) };
-
-
+ final static Point3i[] cubeVertexOffsets = { new Point3i(0, 0, 0), //0 pt
+ new Point3i(1, 0, 0), //1 pt + yz
+ new Point3i(1, 0, 1), //2 pt + yz + 1
+ new Point3i(0, 0, 1), //3 pt + 1
+ new Point3i(0, 1, 0), //4 pt + z
+ new Point3i(1, 1, 0), //5 pt + yz + z
+ new Point3i(1, 1, 1), //6 pt + yz + z + 1
+ new Point3i(0, 1, 1) //7 pt + z + 1
+};
+private final int[] linearOffsets = new int[8];
+int yzCount;
+
+/* set the linear offsets for unique cell ID
+ * and for pointing into the inside/outside BitSet.
+ * Add offset to 0: x * (nY * nZ) + y * nZ + z
+ */
+void setLinearOffsets() {
+ linearOffsets[0] = 0;
+ linearOffsets[1] = yzCount;
+ linearOffsets[2] = yzCount + 1;
+ linearOffsets[3] = 1;
+ linearOffsets[4] = nZ;
+ linearOffsets[5] = yzCount + nZ;
+ linearOffsets[6] = yzCount + nZ + 1;
+ linearOffsets[7] = nZ + 1;
+}
+
+public int getLinearOffset(int x, int y, int z, int offset) {
+ return x * yzCount + y * nZ + z + linearOffsets[offset];
+}
+
+
/* Y
* 4 --------4--------- 5 +z --------4--------- +yz+z
* /| /| /| /|
Modified: trunk/Jmol/src/org/openscience/jvxl/simplewriter/VolumeData.java
===================================================================
--- trunk/Jmol/src/org/openscience/jvxl/simplewriter/VolumeData.java 2009-02-07 05:32:40 UTC (rev 10626)
+++ trunk/Jmol/src/org/openscience/jvxl/simplewriter/VolumeData.java 2009-02-08 14:54:57 UTC (rev 10627)
@@ -161,10 +161,11 @@
return voxelCounts;
}
- public void setVoxelCounts(int nPointsX, int nPointsY, int nPointsZ) {
+ public int setVoxelCounts(int nPointsX, int nPointsY, int nPointsZ) {
voxelCounts[0] = nPointsX;
voxelCounts[1] = nPointsY;
voxelCounts[2] = nPointsZ;
+ return nPointsX * nPointsY * nPointsZ;
}
public float[][][] voxelData;
Modified: trunk/Jmol/src/org/openscience/jvxl/simplewriter/VoxelDataCreator.java
===================================================================
--- trunk/Jmol/src/org/openscience/jvxl/simplewriter/VoxelDataCreator.java 2009-02-07 05:32:40 UTC (rev 10626)
+++ trunk/Jmol/src/org/openscience/jvxl/simplewriter/VoxelDataCreator.java 2009-02-08 14:54:57 UTC (rev 10627)
@@ -4,30 +4,43 @@
public class VoxelDataCreator {
- VoxelDataCreator() {
+ VolumeData volumeData;
+
+ VoxelDataCreator(VolumeData volumeData) {
+ this.volumeData = volumeData;
}
/**
* Developer must customize this method
*
- * @param volumeData
*/
- void createVoxelData(VolumeData volumeData) {
+ void createVoxelData() {
float[][][] voxelData = volumeData.getVoxelData();
int[] counts = volumeData.getVoxelCounts();
int nX = counts[0];
int nY = counts[1];
int nZ = counts[2];
// whatever method here that is desired;
- Point3f pt = new Point3f();
+ // it is not necessary to create a whole block.
+ // you can set volumeData.voxelData = null, in
+ // which case the Marching Cube method will
+ // query getValue(x, y, z) itself.
+
for (int x = 0; x < nX; ++x)
for (int y = 0; y < nY; ++y)
for (int z = 0; z < nZ; ++z) {
- volumeData.voxelPtToXYZ(x, y, z, pt);
- // for instance...
- voxelData[x][y][z] = pt.x * pt.x + pt.y * pt.y - pt.z * pt.z;
+ voxelData[x][y][z] = getValue(x, y, z);
+ //System.out.println("draw a" + x + "" + y + "" + z + " {" + x + " " + y + " " + z + "} \"" + voxelData[x][y][z] + "\"");
}
}
+ Point3f pt = new Point3f();
+
+ public float getValue(int x, int y, int z) {
+ volumeData.voxelPtToXYZ(x, y, z, pt);
+ return pt.x * pt.x + pt.y * pt.y - pt.z * pt.z; // for instance
+ }
+
+
}
This was sent by the SourceForge.net collaborative development platform, the world's largest Open Source development site.
|