Examples with RigidModel2D mpicbg.trakem2.transform.RigidModel2D used on opensource projects

Example 6 with RigidModel2D

use of mpicbg.trakem2.transform.RigidModel2D in project TrakEM2 by trakem2.

the class Align method alignLayersLinearly.

/**
 * Align a range of layers by accumulating pairwise alignments of contiguous layers.
 *
 * @param layers The range of layers to align pairwise.
 * @param numThreads The number of threads to use.
 * @param filter The {@link Filter} to decide which {@link Patch} instances to use in each {@link Layer}. Can be null.
 */
public static final void alignLayersLinearly(final List<Layer> layers, final int numThreads, final Filter<Patch> filter) {
    param.sift.maxOctaveSize = 1600;
    if (!param.setup("Align layers linearly"))
        return;
    final Rectangle box = layers.get(0).getParent().getMinimalBoundingBox(Patch.class);
    final double scale = Math.min(1.0, Math.min((double) param.sift.maxOctaveSize / box.width, (double) param.sift.maxOctaveSize / box.height));
    final Param p = param.clone();
    p.maxEpsilon *= scale;
    final FloatArray2DSIFT sift = new FloatArray2DSIFT(p.sift);
    final SIFT ijSIFT = new SIFT(sift);
    Rectangle box1 = null;
    Rectangle box2 = null;
    final Collection<Feature> features1 = new ArrayList<Feature>();
    final Collection<Feature> features2 = new ArrayList<Feature>();
    final List<PointMatch> candidates = new ArrayList<PointMatch>();
    final List<PointMatch> inliers = new ArrayList<PointMatch>();
    final AffineTransform a = new AffineTransform();
    int i = 0;
    for (final Layer l : layers) {
        long s = System.currentTimeMillis();
        features1.clear();
        features1.addAll(features2);
        features2.clear();
        final Rectangle box3 = l.getMinimalBoundingBox(Patch.class);
        if (box3 == null)
            continue;
        box1 = box2;
        box2 = box3;
        final List<Patch> patches = l.getAll(Patch.class);
        if (null != filter) {
            for (final Iterator<Patch> it = patches.iterator(); it.hasNext(); ) {
                if (!filter.accept(it.next()))
                    it.remove();
            }
        }
        ijSIFT.extractFeatures(l.getProject().getLoader().getFlatImage(l, box2, scale, 0xffffffff, ImagePlus.GRAY8, Patch.class, patches, true).getProcessor(), features2);
        Utils.log(features2.size() + " features extracted in layer \"" + l.getTitle() + "\" (took " + (System.currentTimeMillis() - s) + " ms).");
        if (features1.size() > 0) {
            s = System.currentTimeMillis();
            candidates.clear();
            FeatureTransform.matchFeatures(features2, features1, candidates, p.rod);
            final AbstractAffineModel2D<?> model;
            switch(p.expectedModelIndex) {
                case 0:
                    model = new TranslationModel2D();
                    break;
                case 1:
                    model = new RigidModel2D();
                    break;
                case 2:
                    model = new SimilarityModel2D();
                    break;
                case 3:
                    model = new AffineModel2D();
                    break;
                default:
                    return;
            }
            boolean modelFound;
            boolean again = false;
            try {
                do {
                    again = false;
                    modelFound = model.filterRansac(candidates, inliers, 1000, p.maxEpsilon, p.minInlierRatio, p.minNumInliers, 3);
                    if (modelFound && p.rejectIdentity) {
                        final ArrayList<Point> points = new ArrayList<Point>();
                        PointMatch.sourcePoints(inliers, points);
                        if (Transforms.isIdentity(model, points, p.identityTolerance)) {
                            Utils.log("Identity transform for " + inliers.size() + " matches rejected.");
                            candidates.removeAll(inliers);
                            inliers.clear();
                            again = true;
                        }
                    }
                } while (again);
            } catch (final NotEnoughDataPointsException e) {
                modelFound = false;
            }
            if (modelFound) {
                Utils.log("Model found for layer \"" + l.getTitle() + "\" and its predecessor:\n  correspondences  " + inliers.size() + " of " + candidates.size() + "\n  average residual error  " + (model.getCost() / scale) + " px\n  took " + (System.currentTimeMillis() - s) + " ms");
                final AffineTransform b = new AffineTransform();
                b.translate(box1.x, box1.y);
                b.scale(1.0f / scale, 1.0f / scale);
                b.concatenate(model.createAffine());
                b.scale(scale, scale);
                b.translate(-box2.x, -box2.y);
                a.concatenate(b);
                l.apply(Displayable.class, a);
                Display.repaint(l);
            } else {
                Utils.log("No model found for layer \"" + l.getTitle() + "\" and its predecessor:\n  correspondence candidates  " + candidates.size() + "\n  took " + (System.currentTimeMillis() - s) + " ms");
                a.setToIdentity();
            }
        }
        IJ.showProgress(++i, layers.size());
    }
}

Example 7 with RigidModel2D

use of mpicbg.trakem2.transform.RigidModel2D in project TrakEM2 by trakem2.

the class Align method alignTileCollections.

/**
 * Align two collections of tiles
 * @param p
 * @param a
 * @param b
 */
public static final void alignTileCollections(final Param p, final Collection<AbstractAffineTile2D<?>> a, final Collection<AbstractAffineTile2D<?>> b) {
    final ArrayList<Patch> pa = new ArrayList<Patch>();
    final ArrayList<Patch> pb = new ArrayList<Patch>();
    for (final AbstractAffineTile2D<?> t : a) pa.add(t.getPatch());
    for (final AbstractAffineTile2D<?> t : b) pb.add(t.getPatch());
    final Layer la = pa.iterator().next().getLayer();
    final Layer lb = pb.iterator().next().getLayer();
    final Rectangle boxA = Displayable.getBoundingBox(pa, null);
    final Rectangle boxB = Displayable.getBoundingBox(pb, null);
    final double scale = Math.min(1.0, Math.min(Math.min((double) p.sift.maxOctaveSize / boxA.width, (double) p.sift.maxOctaveSize / boxA.height), Math.min((double) p.sift.maxOctaveSize / boxB.width, (double) p.sift.maxOctaveSize / boxB.height)));
    final Param pp = p.clone();
    pp.maxEpsilon *= scale;
    final FloatArray2DSIFT sift = new FloatArray2DSIFT(pp.sift);
    final SIFT ijSIFT = new SIFT(sift);
    final Collection<Feature> featuresA = new ArrayList<Feature>();
    final Collection<Feature> featuresB = new ArrayList<Feature>();
    final List<PointMatch> candidates = new ArrayList<PointMatch>();
    final List<PointMatch> inliers = new ArrayList<PointMatch>();
    long s = System.currentTimeMillis();
    ijSIFT.extractFeatures(la.getProject().getLoader().getFlatImage(la, boxA, scale, 0xffffffff, ImagePlus.GRAY8, null, pa, true, Color.GRAY).getProcessor(), featuresA);
    Utils.log(featuresA.size() + " features extracted in graph A in layer \"" + la.getTitle() + "\" (took " + (System.currentTimeMillis() - s) + " ms).");
    s = System.currentTimeMillis();
    ijSIFT.extractFeatures(lb.getProject().getLoader().getFlatImage(lb, boxB, scale, 0xffffffff, ImagePlus.GRAY8, null, pb, true, Color.GRAY).getProcessor(), featuresB);
    Utils.log(featuresB.size() + " features extracted in graph B in layer \"" + lb.getTitle() + "\" (took " + (System.currentTimeMillis() - s) + " ms).");
    if (featuresA.size() > 0 && featuresB.size() > 0) {
        s = System.currentTimeMillis();
        FeatureTransform.matchFeatures(featuresA, featuresB, candidates, pp.rod);
        final AbstractAffineModel2D<?> model;
        switch(p.expectedModelIndex) {
            case 0:
                model = new TranslationModel2D();
                break;
            case 1:
                model = new RigidModel2D();
                break;
            case 2:
                model = new SimilarityModel2D();
                break;
            case 3:
                model = new AffineModel2D();
                break;
            default:
                return;
        }
        boolean modelFound;
        boolean again = false;
        try {
            do {
                again = false;
                modelFound = model.filterRansac(candidates, inliers, 1000, p.maxEpsilon, p.minInlierRatio, p.minNumInliers, 3);
                if (modelFound && p.rejectIdentity) {
                    final ArrayList<Point> points = new ArrayList<Point>();
                    PointMatch.sourcePoints(inliers, points);
                    if (Transforms.isIdentity(model, points, p.identityTolerance)) {
                        Utils.log("Identity transform for " + inliers.size() + " matches rejected.");
                        candidates.removeAll(inliers);
                        inliers.clear();
                        again = true;
                    }
                }
            } while (again);
        } catch (final NotEnoughDataPointsException e) {
            modelFound = false;
        }
        if (modelFound) {
            Utils.log("Model found for graph A and B in layers \"" + la.getTitle() + "\" and \"" + lb.getTitle() + "\":\n  correspondences  " + inliers.size() + " of " + candidates.size() + "\n  average residual error  " + (model.getCost() / scale) + " px\n  took " + (System.currentTimeMillis() - s) + " ms");
            final AffineTransform at = new AffineTransform();
            at.translate(boxA.x, boxA.y);
            at.scale(1.0f / scale, 1.0f / scale);
            at.concatenate(model.createAffine());
            at.scale(scale, scale);
            at.translate(-boxB.x, -boxB.y);
            for (final Patch t : pa) t.preTransform(at, false);
            Display.repaint(la);
        } else
            Utils.log("No model found for graph A and B in layers \"" + la.getTitle() + "\" and \"" + lb.getTitle() + "\":\n  correspondence candidates  " + candidates.size() + "\n  took " + (System.currentTimeMillis() - s) + " ms");
    }
}

Example 8 with RigidModel2D

use of mpicbg.trakem2.transform.RigidModel2D in project TrakEM2 by trakem2.

the class AlignLayersTask method alignLayersLinearlyJob.

public static final void alignLayersLinearlyJob(final LayerSet layerSet, final int first, final int last, final boolean propagateTransform, final Rectangle fov, final Filter<Patch> filter) {
    // will reverse order if necessary
    final List<Layer> layerRange = layerSet.getLayers(first, last);
    final Align.Param p = Align.param.clone();
    // find the first non-empty layer, and remove all empty layers
    Rectangle box = fov;
    for (final Iterator<Layer> it = layerRange.iterator(); it.hasNext(); ) {
        final Layer la = it.next();
        if (!la.contains(Patch.class, true)) {
            it.remove();
            continue;
        }
        if (null == box) {
            // The first layer:
            // Only for visible patches
            box = la.getMinimalBoundingBox(Patch.class, true);
        }
    }
    if (0 == layerRange.size()) {
        Utils.log("All layers in range are empty!");
        return;
    }
    /* do not work if there is only one layer selected */
    if (layerRange.size() < 2)
        return;
    final double scale = Math.min(1.0, Math.min((double) p.sift.maxOctaveSize / (double) box.width, (double) p.sift.maxOctaveSize / (double) box.height));
    p.maxEpsilon *= scale;
    p.identityTolerance *= scale;
    // Utils.log2("scale: " + scale + "  maxOctaveSize: " + p.sift.maxOctaveSize + "  box: " + box.width + "," + box.height);
    final FloatArray2DSIFT sift = new FloatArray2DSIFT(p.sift);
    final SIFT ijSIFT = new SIFT(sift);
    Rectangle box1 = fov;
    Rectangle box2 = fov;
    final Collection<Feature> features1 = new ArrayList<Feature>();
    final Collection<Feature> features2 = new ArrayList<Feature>();
    final List<PointMatch> candidates = new ArrayList<PointMatch>();
    final List<PointMatch> inliers = new ArrayList<PointMatch>();
    final AffineTransform a = new AffineTransform();
    int s = 0;
    for (final Layer layer : layerRange) {
        if (Thread.currentThread().isInterrupted())
            return;
        final long t0 = System.currentTimeMillis();
        features1.clear();
        features1.addAll(features2);
        features2.clear();
        final Rectangle box3 = layer.getMinimalBoundingBox(Patch.class, true);
        // skipping empty layer
        if (box3 == null || (box3.width == 0 && box3.height == 0))
            continue;
        box1 = null == fov ? box2 : fov;
        box2 = null == fov ? box3 : fov;
        final List<Patch> patches = layer.getAll(Patch.class);
        if (null != filter) {
            for (final Iterator<Patch> it = patches.iterator(); it.hasNext(); ) {
                if (!filter.accept(it.next()))
                    it.remove();
            }
        }
        final ImageProcessor flatImage = layer.getProject().getLoader().getFlatImage(layer, box2, scale, 0xffffffff, ImagePlus.GRAY8, Patch.class, patches, true).getProcessor();
        ijSIFT.extractFeatures(flatImage, features2);
        IJ.log(features2.size() + " features extracted in layer \"" + layer.getTitle() + "\" (took " + (System.currentTimeMillis() - t0) + " ms).");
        if (features1.size() > 0) {
            final long t1 = System.currentTimeMillis();
            candidates.clear();
            FeatureTransform.matchFeatures(features2, features1, candidates, p.rod);
            final AbstractAffineModel2D<?> model;
            switch(p.expectedModelIndex) {
                case 0:
                    model = new TranslationModel2D();
                    break;
                case 1:
                    model = new RigidModel2D();
                    break;
                case 2:
                    model = new SimilarityModel2D();
                    break;
                case 3:
                    model = new AffineModel2D();
                    break;
                default:
                    return;
            }
            final AbstractAffineModel2D<?> desiredModel;
            switch(p.desiredModelIndex) {
                case 0:
                    desiredModel = new TranslationModel2D();
                    break;
                case 1:
                    desiredModel = new RigidModel2D();
                    break;
                case 2:
                    desiredModel = new SimilarityModel2D();
                    break;
                case 3:
                    desiredModel = new AffineModel2D();
                    break;
                default:
                    return;
            }
            boolean modelFound;
            boolean again = false;
            try {
                do {
                    again = false;
                    modelFound = model.filterRansac(candidates, inliers, 1000, p.maxEpsilon, p.minInlierRatio, p.minNumInliers, 3);
                    if (modelFound && p.rejectIdentity) {
                        final ArrayList<Point> points = new ArrayList<Point>();
                        PointMatch.sourcePoints(inliers, points);
                        if (Transforms.isIdentity(model, points, p.identityTolerance)) {
                            IJ.log("Identity transform for " + inliers.size() + " matches rejected.");
                            candidates.removeAll(inliers);
                            inliers.clear();
                            again = true;
                        }
                    }
                } while (again);
                if (modelFound)
                    desiredModel.fit(inliers);
            } catch (final NotEnoughDataPointsException e) {
                modelFound = false;
            } catch (final IllDefinedDataPointsException e) {
                modelFound = false;
            }
            if (Thread.currentThread().isInterrupted())
                return;
            if (modelFound) {
                IJ.log("Model found for layer \"" + layer.getTitle() + "\" and its predecessor:\n  correspondences  " + inliers.size() + " of " + candidates.size() + "\n  average residual error  " + (model.getCost() / scale) + " px\n  took " + (System.currentTimeMillis() - t1) + " ms");
                final AffineTransform b = new AffineTransform();
                b.translate(box1.x, box1.y);
                b.scale(1.0f / scale, 1.0f / scale);
                b.concatenate(desiredModel.createAffine());
                b.scale(scale, scale);
                b.translate(-box2.x, -box2.y);
                a.concatenate(b);
                AlignTask.transformPatchesAndVectorData(patches, a);
                Display.repaint(layer);
            } else {
                IJ.log("No model found for layer \"" + layer.getTitle() + "\" and its predecessor:\n  correspondence candidates  " + candidates.size() + "\n  took " + (System.currentTimeMillis() - s) + " ms");
                a.setToIdentity();
            }
        }
        IJ.showProgress(++s, layerRange.size());
    }
    if (Thread.currentThread().isInterrupted())
        return;
    if (propagateTransform) {
        if (last > first && last < layerSet.size() - 2)
            for (final Layer la : layerSet.getLayers(last + 1, layerSet.size() - 1)) {
                if (Thread.currentThread().isInterrupted())
                    return;
                AlignTask.transformPatchesAndVectorData(la, a);
            }
        else if (first > last && last > 0)
            for (final Layer la : layerSet.getLayers(0, last - 1)) {
                if (Thread.currentThread().isInterrupted())
                    return;
                AlignTask.transformPatchesAndVectorData(la, a);
            }
    }
}

Example 9 with RigidModel2D

use of mpicbg.trakem2.transform.RigidModel2D in project TrakEM2 by trakem2.

the class ElasticLayerAlignment method preAlignStack.

private void preAlignStack(final Param param, final Project project, final List<Layer> layerRange, final Rectangle box, final Filter<Patch> filter, final ArrayList<Triple<Integer, Integer, AbstractModel<?>>> pairs) {
    final double scale = Math.min(1.0, Math.min((double) param.ppm.sift.maxOctaveSize / (double) box.width, (double) param.ppm.sift.maxOctaveSize / (double) box.height));
    /* extract and save features, overwrite cached files if requested */
    try {
        AlignmentUtils.extractAndSaveLayerFeatures(layerRange, box, scale, filter, param.ppm.sift, param.ppm.clearCache, param.ppm.maxNumThreadsSift);
    } catch (final Exception e) {
        return;
    }
    /* match and filter feature correspondences */
    int numFailures = 0;
    final double pointMatchScale = param.layerScale / scale;
    for (int i = 0; i < layerRange.size(); ++i) {
        final ArrayList<Thread> threads = new ArrayList<Thread>(param.maxNumThreads);
        final int sliceA = i;
        final Layer layerA = layerRange.get(i);
        final int range = Math.min(layerRange.size(), i + param.maxNumNeighbors + 1);
        final String layerNameA = layerName(layerA);
        for (int j = i + 1; j < range; ) J: {
            final int numThreads = Math.min(param.maxNumThreads, range - j);
            final ArrayList<Triple<Integer, Integer, AbstractModel<?>>> models = new ArrayList<Triple<Integer, Integer, AbstractModel<?>>>(numThreads);
            for (int k = 0; k < numThreads; ++k) models.add(null);
            for (int t = 0; t < numThreads && j < range; ++t, ++j) {
                final int ti = t;
                final int sliceB = j;
                final Layer layerB = layerRange.get(j);
                final String layerNameB = layerName(layerB);
                final Thread thread = new Thread() {

                    @Override
                    public void run() {
                        IJ.showProgress(sliceA, layerRange.size() - 1);
                        Utils.log("matching " + layerNameB + " -> " + layerNameA + "...");
                        ArrayList<PointMatch> candidates = null;
                        if (!param.ppm.clearCache)
                            candidates = mpicbg.trakem2.align.Util.deserializePointMatches(project, param.ppm, "layer", layerB.getId(), layerA.getId());
                        if (null == candidates) {
                            final ArrayList<Feature> fs1 = mpicbg.trakem2.align.Util.deserializeFeatures(project, param.ppm.sift, "layer", layerA.getId());
                            final ArrayList<Feature> fs2 = mpicbg.trakem2.align.Util.deserializeFeatures(project, param.ppm.sift, "layer", layerB.getId());
                            candidates = new ArrayList<PointMatch>(FloatArray2DSIFT.createMatches(fs2, fs1, param.ppm.rod));
                            /* scale the candidates */
                            for (final PointMatch pm : candidates) {
                                final Point p1 = pm.getP1();
                                final Point p2 = pm.getP2();
                                final double[] l1 = p1.getL();
                                final double[] w1 = p1.getW();
                                final double[] l2 = p2.getL();
                                final double[] w2 = p2.getW();
                                l1[0] *= pointMatchScale;
                                l1[1] *= pointMatchScale;
                                w1[0] *= pointMatchScale;
                                w1[1] *= pointMatchScale;
                                l2[0] *= pointMatchScale;
                                l2[1] *= pointMatchScale;
                                w2[0] *= pointMatchScale;
                                w2[1] *= pointMatchScale;
                            }
                            if (!mpicbg.trakem2.align.Util.serializePointMatches(project, param.ppm, "layer", layerB.getId(), layerA.getId(), candidates))
                                Utils.log("Could not store point match candidates for layers " + layerNameB + " and " + layerNameA + ".");
                        }
                        AbstractModel<?> model;
                        switch(param.expectedModelIndex) {
                            case 0:
                                model = new TranslationModel2D();
                                break;
                            case 1:
                                model = new RigidModel2D();
                                break;
                            case 2:
                                model = new SimilarityModel2D();
                                break;
                            case 3:
                                model = new AffineModel2D();
                                break;
                            case 4:
                                model = new HomographyModel2D();
                                break;
                            default:
                                return;
                        }
                        final ArrayList<PointMatch> inliers = new ArrayList<PointMatch>();
                        boolean modelFound;
                        boolean again = false;
                        try {
                            do {
                                again = false;
                                modelFound = model.filterRansac(candidates, inliers, 1000, param.maxEpsilon * param.layerScale, param.minInlierRatio, param.minNumInliers, 3);
                                if (modelFound && param.rejectIdentity) {
                                    final ArrayList<Point> points = new ArrayList<Point>();
                                    PointMatch.sourcePoints(inliers, points);
                                    if (Transforms.isIdentity(model, points, param.identityTolerance * param.layerScale)) {
                                        IJ.log("Identity transform for " + inliers.size() + " matches rejected.");
                                        candidates.removeAll(inliers);
                                        inliers.clear();
                                        again = true;
                                    }
                                }
                            } while (again);
                        } catch (final NotEnoughDataPointsException e) {
                            modelFound = false;
                        }
                        if (modelFound) {
                            Utils.log(layerNameB + " -> " + layerNameA + ": " + inliers.size() + " corresponding features with an average displacement of " + (PointMatch.meanDistance(inliers) / param.layerScale) + "px identified.");
                            Utils.log("Estimated transformation model: " + model);
                            models.set(ti, new Triple<Integer, Integer, AbstractModel<?>>(sliceA, sliceB, model));
                        } else {
                            Utils.log(layerNameB + " -> " + layerNameA + ": no correspondences found.");
                            return;
                        }
                    }
                };
                threads.add(thread);
                thread.start();
            }
            try {
                for (final Thread thread : threads) thread.join();
            } catch (final InterruptedException e) {
                Utils.log("Establishing feature correspondences interrupted.");
                for (final Thread thread : threads) thread.interrupt();
                try {
                    for (final Thread thread : threads) thread.join();
                } catch (final InterruptedException f) {
                }
                return;
            }
            threads.clear();
            /* collect successfully matches pairs and break the search on gaps */
            for (int t = 0; t < models.size(); ++t) {
                final Triple<Integer, Integer, AbstractModel<?>> pair = models.get(t);
                if (pair == null) {
                    if (++numFailures > param.maxNumFailures) {
                        break J;
                    }
                } else {
                    numFailures = 0;
                    pairs.add(pair);
                }
            }
        }
    }
}