Examples with Vector3d org.scijava.vecmath.Vector3d used on opensource projects

Example 11 with Vector3d

use of org.scijava.vecmath.Vector3d in project TrakEM2 by trakem2.

the class Polyline method simplify.

/**
 * Resample the curve to, first, a number of points as resulting from resampling to a point interdistance of delta, and second, as adjustment by random walk of those points towards the original points.
 */
public static Object[] simplify(final double[] px, final double[] py, final long[] p_layer_ids, /*final double delta, final double allowed_error_per_point,*/
final int max_iterations, final LayerSet layer_set) throws Exception {
    if (px.length != py.length || py.length != p_layer_ids.length)
        return null;
    final double[] pz = new double[px.length];
    for (int i = 0; i < pz.length; i++) {
        pz[i] = layer_set.getLayer(p_layer_ids[i]).getZ();
    }
    /*
		// Resample:
		VectorString3D vs = new VectorString3D(px, py, pz, false);
		Calibration cal = layer_set.getCalibrationCopy();
		vs.calibrate(cal);
		vs.resample(delta);
		cal.pixelWidth = 1 / cal.pixelWidth;
		cal.pixelHeight = 1 / cal.pixelHeight;
		vs.calibrate(cal); // undo it, since source points are in pixels

		Pth path = new Pth(vs.getPoints(0), vs.getPoints(1), vs.getPoints(2));
		vs = null;

		// The above fails with strangely jagged lines.
		*/
    // Instead, just a line:
    final Calibration cal = layer_set.getCalibrationCopy();
    // in pixels
    final double one_unit = 1 / cal.pixelWidth;
    final double traced_length = M.distance(px[0], py[0], pz[0], px[px.length - 1], py[py.length - 1], pz[pz.length - 1]);
    final double segment_length = (one_unit * 5);
    final int n_new_points = (int) (traced_length / segment_length) + 1;
    double[] rx = new double[n_new_points];
    double[] ry = new double[n_new_points];
    double[] rz = new double[n_new_points];
    rx[0] = px[0];
    rx[rx.length - 1] = px[px.length - 1];
    ry[0] = py[0];
    ry[ry.length - 1] = py[py.length - 1];
    rz[0] = pz[0];
    rz[rz.length - 1] = pz[pz.length - 1];
    final Vector3d v = new Vector3d(rx[n_new_points - 1] - rx[0], ry[n_new_points - 1] - ry[0], rz[n_new_points - 1] - rz[0]);
    v.normalize();
    v.scale(segment_length);
    for (int i = 1; i < n_new_points - 1; i++) {
        rx[i] = rx[0] + v.x * i;
        ry[i] = ry[0] + v.y * i;
        rz[i] = rz[0] + v.z * i;
    }
    Pth path = new Pth(rx, ry, rz);
    rx = ry = rz = null;
    // final double lowest_error = px.length * allowed_error_per_point;
    final double d = 1;
    final Random rand = new Random(System.currentTimeMillis());
    double current_error = Double.MAX_VALUE;
    int i = 0;
    // double[] er = new double[max_iterations];
    // double[] index = new double[max_iterations];
    int missed_in_a_row = 0;
    for (; i < max_iterations; i++) {
        final Pth copy = path.copy().shakeUpOne(d, rand);
        final double error = copy.measureErrorSq(px, py, pz);
        // If less error, keep the copy
        if (error < current_error) {
            current_error = error;
            path = copy;
            // er[i] = error;
            // index[i] = i;
            missed_in_a_row = 0;
        } else {
            // er[i] = current_error;
            // index[i] = i;
            missed_in_a_row++;
            if (missed_in_a_row > 10 * path.px.length) {
                Utils.log2("Stopped random walk at iteration " + i);
                break;
            }
            continue;
        }
    /*
			// If below lowest_error, quit searching
			if (current_error < lowest_error) {
				Utils.log2("Stopped at iteration " + i);
				break;
			}
			*/
    }
    if (max_iterations == i) {
        Utils.log2("Reached max iterations -- current error: " + current_error);
    }
    // Approximate new Z to a layer id:
    final long[] plids = new long[path.px.length];
    // first point untouched
    plids[0] = p_layer_ids[0];
    for (int k = 1; k < path.pz.length; k++) {
        plids[k] = layer_set.getNearestLayer(path.pz[k]).getId();
    }
    return new Object[] { path.px, path.py, plids };
}