use of mpicbg.imglib.type.numeric.integer.ByteType in project TrakEM2 by trakem2.
the class AreaUtils method generateTriangles.
/**
* Expects areas in local coordinates to the Displayable @param d.
* @param d
* @param scale The scaling of the entire universe, to limit the overall box
* @param resample_ The optimization parameter for marching cubes (i.e. a value of 2 will scale down to half, then apply marching cubes, then scale up by 2 the vertices coordinates).
* @param areas
* @return The List of triangles involved, specified as three consecutive vertices. A list of Point3f vertices.
*/
public static List<Point3f> generateTriangles(final Displayable d, final double scale, final int resample_, final Map<Layer, Area> areas) {
// in the LayerSet, layers are ordered by Z already.
try {
int n = areas.size();
if (0 == n)
return null;
final int resample;
if (resample_ <= 0) {
resample = 1;
Utils.log2("Fixing zero or negative resampling value to 1.");
} else
resample = resample_;
final LayerSet layer_set = d.getLayerSet();
final AffineTransform aff = d.getAffineTransformCopy();
final Rectangle r = d.getBoundingBox(null);
// remove translation from a copy of the Displayable's AffineTransform
final AffineTransform at_translate = new AffineTransform();
at_translate.translate(-r.x, -r.y);
aff.preConcatenate(at_translate);
// incorporate resampling scaling into the transform
final AffineTransform atK = new AffineTransform();
// Utils.log("resample: " + resample + " scale: " + scale);
// 'scale' is there to limit gigantic universes
final double K = (1.0 / resample) * scale;
atK.scale(K, K);
aff.preConcatenate(atK);
final Calibration cal = layer_set.getCalibrationCopy();
// Find first layer, compute depth, and fill in the depth vs area map
Layer first_layer = null, last_layer = null;
final int w = (int) Math.ceil(r.width * K);
final int h = (int) Math.ceil(r.height * K);
int depth = 0;
final Map<Integer, Area> ma = new HashMap<Integer, Area>();
for (final Layer la : layer_set.getLayers()) {
// layers sorted by Z ASC
final Area area = areas.get(la);
if (null != area) {
ma.put(depth, area);
if (null == first_layer) {
first_layer = la;
}
// Utils.log("area at depth " + depth + " for layer " + la);
depth++;
n--;
} else if (0 != depth) {
// Utils.log("Empty area at depth " + depth);
// an empty layer
depth++;
}
if (0 == n) {
last_layer = la;
// no more areas to paint
break;
}
}
if (0 == depth) {
Utils.log("ERROR could not find any areas for " + d);
return null;
}
if (0 != n) {
Utils.log("WARNING could not find all areas for " + d);
}
// No zero-padding: Marching Cubes now can handle edges
final ShapeList<ByteType> shapeList = new ShapeListCached<ByteType>(new int[] { w, h, depth }, new ByteType(), 32);
final Image<ByteType> shapeListImage = new Image<ByteType>(shapeList, shapeList.getBackground(), "ShapeListContainer");
// 255 or -1 don't work !? So, giving the highest value (127) that is both a byte and an int.
final ByteType intensity = new ByteType((byte) 127);
for (final Map.Entry<Integer, Area> e : ma.entrySet()) {
Area a = e.getValue();
if (!aff.isIdentity()) {
a = M.areaInIntsByRounding(a.createTransformedArea(aff));
}
shapeList.addShape(a, intensity, new int[] { e.getKey() });
}
// debug:
// ImagePlus imp = ImageJFunctions.displayAsVirtualStack(shapeListImage);
// imp.getProcessor().setMinAndMax( 0, 255 );
// imp.show();
// Utils.log2("Using imglib Shape List Image Container");
// Now marching cubes
// origins at 0,0,0: uncalibrated
final List<Point3f> list = new MCTriangulator().getTriangles(shapeListImage, 1, new float[3]);
// The list of triangles has coordinates:
// - in x,y: in pixels, scaled by K = (1 / resample) * scale,
// translated by r.x, r.y (the top-left coordinate of this AreaList bounding box)
// - in z: in stack slice indices
// So all x,y,z must be corrected in x,y and z of the proper layer
// final double offset = first_layer.getZ();
final int i_first_layer = layer_set.indexOf(first_layer);
// The x,y translation to correct each point by:
final float dx = (float) (r.x * scale * cal.pixelWidth);
final float dy = (float) (r.y * scale * cal.pixelHeight);
// Correct x,y by resampling and calibration, but not scale
// scale is already in the pixel coordinates
final float rsw = (float) (resample * cal.pixelWidth);
final float rsh = (float) (resample * cal.pixelHeight);
// no resampling in Z. and Uses pixelWidth, not pixelDepth.
final double sz = scale * cal.pixelWidth;
// debug:
/*
// which p.z types exist?
final TreeSet<Float> ts = new TreeSet<Float>();
for (final Iterator it = list.iterator(); it.hasNext(); ) {
ts.add(((Point3f)it.next()).z);
}
for (final Float pz : ts) Utils.log2("A z: " + pz);
*/
// debug: How many different Z?
/*
HashSet<Float> zs = new HashSet<Float>();
for (Point3f p : list) {
zs.add(p.z);
}
ArrayList<Float> a = new ArrayList<Float>(zs);
java.util.Collections.sort(a);
for (Float f : a) {
Utils.log("f: " + f);
}
*/
// Utils.log2("Number of slices: " + imp.getNSlices());
// Fix all points:
// Read from list, modify and put into verts
// and don't modify it if the verts already has it (it's just coincident)
final Point3f[] verts = new Point3f[list.size()];
// Utils.log("number of verts: " + verts.length + " mod 3: " + (verts.length % 3));
final TreeMap<Integer, Point3f> output = new TreeMap<Integer, Point3f>();
// The first section generates vertices at -1 and 0
// The last section generates them at last_section_index and last_section_index +1
// Capture from -1 to 0
fix3DPoints(list, output, verts, first_layer.getZ(), 0, -1, dx, dy, rsw, rsh, sz, 1);
int slice_index = 0;
for (final Layer la : layer_set.getLayers().subList(i_first_layer, i_first_layer + depth)) {
// If layer is empty, continue
/* // YEAH don't! At least the immediate next layer would have points, like the extra Z level after last layer, to account for the thickness of the layer!
if (empty_layers.contains(la)) {
slice_index++;
continue;
}
*/
fix3DPoints(list, output, verts, la.getZ(), la.getThickness(), slice_index, dx, dy, rsw, rsh, sz, 1);
slice_index++;
}
// Do the last layer again. The last layer has two Z planes in which it has pixels:
try {
// Capture from last_section_index to last_section_index+1, inclusive
fix3DPoints(list, output, verts, last_layer.getZ() + last_layer.getThickness(), 0, slice_index, dx, dy, rsw, rsh, sz, 2);
} catch (final Exception ee) {
IJError.print(ee);
}
// Handle potential errors:
if (0 != list.size() - output.size()) {
Utils.log2("Unprocessed/unused points: " + (list.size() - output.size()));
for (int i = 0; i < verts.length; i++) {
if (null == verts[i]) {
final Point3f p = (Point3f) list.get(i);
Utils.log2("verts[" + i + "] = " + p.x + ", " + p.y + ", " + p.z + " p.z as int: " + ((int) (p.z + 0.05f)));
}
}
return new ArrayList<Point3f>(output.values());
} else {
return java.util.Arrays.asList(verts);
}
} catch (final Exception e) {
e.printStackTrace();
}
return null;
}
Aggregations