use of ij.ImagePlus in project TrakEM2 by trakem2.
the class LayerStack method getImagePlus.
public synchronized ImagePlus getImagePlus() {
if (null == layer_imp) {
layer_imp = new ImagePlus("LayerSet Stack", this);
Calibration cal = layers.get(0).getParent().getCalibrationCopy();
// adjust to scale
cal.pixelWidth /= scale;
cal.pixelHeight /= scale;
// Simulate depth: assume all layers have the same thickness
// not pixelDepth
cal.pixelDepth = (layers.get(0).getThickness() * cal.pixelWidth) / scale;
layer_imp.setCalibration(cal);
}
return layer_imp;
}
use of ij.ImagePlus in project TrakEM2 by trakem2.
the class AlignLayersTask method alignLayersNonLinearlyJob.
public static final void alignLayersNonLinearlyJob(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();
// Remove all empty layers
for (final Iterator<Layer> it = layerRange.iterator(); it.hasNext(); ) {
if (!it.next().contains(Patch.class, true)) {
it.remove();
}
}
if (0 == layerRange.size()) {
Utils.log("No layers in range show any images!");
return;
}
/* do not work if there is only one layer selected */
if (layerRange.size() < 2)
return;
final List<Patch> all = new ArrayList<Patch>();
for (final Layer la : layerRange) {
for (final Patch patch : la.getAll(Patch.class)) {
if (null != filter && !filter.accept(patch))
continue;
all.add(patch);
}
}
AlignTask.transformPatchesAndVectorData(all, new Runnable() {
@Override
public void run() {
// ///
final Loader loader = layerSet.getProject().getLoader();
// Not concurrent safe! So two copies, one per layer and Thread:
final SIFT ijSIFT1 = new SIFT(new FloatArray2DSIFT(p.sift));
final SIFT ijSIFT2 = new SIFT(new FloatArray2DSIFT(p.sift));
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 int n_proc = Runtime.getRuntime().availableProcessors() > 1 ? 2 : 1;
final ExecutorService exec = Utils.newFixedThreadPool(n_proc, "alignLayersNonLinearly");
List<Patch> previousPatches = null;
int s = 0;
for (int i = 1; i < layerRange.size(); ++i) {
if (Thread.currentThread().isInterrupted())
break;
final Layer layer1 = layerRange.get(i - 1);
final Layer layer2 = layerRange.get(i);
final long t0 = System.currentTimeMillis();
features1.clear();
features2.clear();
final Rectangle box1 = null == fov ? layer1.getMinimalBoundingBox(Patch.class, true) : fov;
final Rectangle box2 = null == fov ? layer2.getMinimalBoundingBox(Patch.class, true) : fov;
/* calculate the common scale factor for both flat images */
final double scale = Math.min(1.0f, (double) p.sift.maxOctaveSize / (double) Math.max(box1.width, Math.max(box1.height, Math.max(box2.width, box2.height))));
final List<Patch> patches1;
if (null == previousPatches) {
patches1 = layer1.getAll(Patch.class);
if (null != filter) {
for (final Iterator<Patch> it = patches1.iterator(); it.hasNext(); ) {
if (!filter.accept(it.next()))
it.remove();
}
}
} else {
patches1 = previousPatches;
}
final List<Patch> patches2 = layer2.getAll(Patch.class);
if (null != filter) {
for (final Iterator<Patch> it = patches2.iterator(); it.hasNext(); ) {
if (!filter.accept(it.next()))
it.remove();
}
}
final Future<ImageProcessor> fu1 = exec.submit(new Callable<ImageProcessor>() {
@Override
public ImageProcessor call() {
final ImageProcessor ip1 = loader.getFlatImage(layer1, box1, scale, 0xffffffff, ImagePlus.GRAY8, Patch.class, patches1, true).getProcessor();
ijSIFT1.extractFeatures(ip1, features1);
Utils.log(features1.size() + " features extracted in layer \"" + layer1.getTitle() + "\" (took " + (System.currentTimeMillis() - t0) + " ms).");
return ip1;
}
});
final Future<ImageProcessor> fu2 = exec.submit(new Callable<ImageProcessor>() {
@Override
public ImageProcessor call() {
final ImageProcessor ip2 = loader.getFlatImage(layer2, box2, scale, 0xffffffff, ImagePlus.GRAY8, Patch.class, patches2, true).getProcessor();
ijSIFT2.extractFeatures(ip2, features2);
Utils.log(features2.size() + " features extracted in layer \"" + layer2.getTitle() + "\" (took " + (System.currentTimeMillis() - t0) + " ms).");
return ip2;
}
});
final ImageProcessor ip1, ip2;
try {
ip1 = fu1.get();
ip2 = fu2.get();
} catch (final Exception e) {
IJError.print(e);
return;
}
if (features1.size() > 0 && features2.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;
}
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);
} catch (final NotEnoughDataPointsException e) {
modelFound = false;
}
if (modelFound) {
IJ.log("Model found for layer \"" + layer2.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 ImagePlus imp1 = new ImagePlus("target", ip1);
final ImagePlus imp2 = new ImagePlus("source", ip2);
final List<Point> sourcePoints = new ArrayList<Point>();
final List<Point> targetPoints = new ArrayList<Point>();
PointMatch.sourcePoints(inliers, sourcePoints);
PointMatch.targetPoints(inliers, targetPoints);
imp2.setRoi(Util.pointsToPointRoi(sourcePoints));
imp1.setRoi(Util.pointsToPointRoi(targetPoints));
final ImageProcessor mask1 = ip1.duplicate();
mask1.threshold(1);
final ImageProcessor mask2 = ip2.duplicate();
mask2.threshold(1);
final Transformation warp = bUnwarpJ_.computeTransformationBatch(imp2, imp1, mask2, mask1, elasticParam);
final CubicBSplineTransform transf = new CubicBSplineTransform();
transf.set(warp.getIntervals(), warp.getDirectDeformationCoefficientsX(), warp.getDirectDeformationCoefficientsY(), imp2.getWidth(), imp2.getHeight());
final ArrayList<Future<?>> fus = new ArrayList<Future<?>>();
// Transform desired patches only
for (final Patch patch : patches2) {
try {
final Rectangle pbox = patch.getCoordinateTransformBoundingBox();
final AffineTransform at = patch.getAffineTransform();
final AffineTransform pat = new AffineTransform();
pat.scale(scale, scale);
pat.translate(-box2.x, -box2.y);
pat.concatenate(at);
pat.translate(-pbox.x, -pbox.y);
final mpicbg.trakem2.transform.AffineModel2D toWorld = new mpicbg.trakem2.transform.AffineModel2D();
toWorld.set(pat);
final CoordinateTransformList<CoordinateTransform> ctl = new CoordinateTransformList<CoordinateTransform>();
// move the patch into the global space where bUnwarpJ calculated the transformation
ctl.add(toWorld);
// Apply non-linear transformation
ctl.add(transf);
// move it back
ctl.add(toWorld.createInverse());
patch.appendCoordinateTransform(ctl);
fus.add(patch.updateMipMaps());
// Compensate for offset between boxes
final AffineTransform offset = new AffineTransform();
offset.translate(box1.x - box2.x, box1.y - box2.y);
offset.concatenate(at);
patch.setAffineTransform(offset);
} catch (final Exception e) {
e.printStackTrace();
}
}
// await regeneration of all mipmaps
Utils.wait(fus);
Display.repaint(layer2);
} else
IJ.log("No model found for layer \"" + layer2.getTitle() + "\" and its predecessor:\n correspondence candidates " + candidates.size() + "\n took " + (System.currentTimeMillis() - s) + " ms");
}
IJ.showProgress(++s, layerRange.size());
// for next iteration
previousPatches = patches2;
}
exec.shutdown();
if (propagateTransform)
Utils.log("Propagation not implemented yet for non-linear layer alignment.");
/* // CANNOT be done (at least not trivially:
* //an appropriate "scale" cannot be computed, and the box2 is part of the spline computation.
if ( propagateTransform && null != lastTransform )
{
for (final Layer la : l.getParent().getLayers(last > first ? last +1 : first -1, last > first ? l.getParent().size() -1 : 0)) {
// Transform visible patches only
final Rectangle box2 = la.getMinimalBoundingBox( Patch.class, true );
for ( final Displayable disp : la.getDisplayables( Patch.class, true ) )
{
// ...
}
}
}
*/
}
});
// end of transformPatchesAndVectorData
}
use of ij.ImagePlus in project TrakEM2 by trakem2.
the class ExportUnsignedShort method exportTEST.
public static final void exportTEST(final Layer layer, final int tileWidth, final int tileHeight) {
/* calculate intensity transfer */
final ArrayList<Displayable> patches = layer.getDisplayables(Patch.class);
final ArrayList<PatchIntensityRange> patchIntensityRanges = new ArrayList<PatchIntensityRange>();
double min = Double.MAX_VALUE;
double max = -Double.MAX_VALUE;
for (final Displayable d : patches) {
final Patch patch = (Patch) d;
final PatchIntensityRange pir = new PatchIntensityRange(patch);
if (pir.min < min)
min = pir.min;
if (pir.max > max)
max = pir.max;
patchIntensityRanges.add(pir);
}
/* render tiles */
/* TODO Do not render them into a stack but save them as files */
final ImageStack stack = new ImageStack(tileWidth, tileHeight);
ImagePlus imp = null;
final double minI = -min * 65535.0 / (max - min);
final double maxI = (1.0 - min) * 65535.0 / (max - min);
// ij.IJ.log("min, max: " + min + ", " + max + ", minI, maxI: " + minI + ", " + maxI);
final int nc = (int) Math.ceil(layer.getLayerWidth() / tileWidth);
final int nr = (int) Math.ceil(layer.getLayerHeight() / tileHeight);
for (int r = 0; r < nr; ++r) {
final int y0 = r * tileHeight;
for (int c = 0; c < nc; ++c) {
final int x0 = c * tileWidth;
final Rectangle box = new Rectangle(x0, y0, tileWidth, tileHeight);
final ShortProcessor sp = new ShortProcessor(tileWidth, tileHeight);
sp.setMinAndMax(minI, maxI);
for (final PatchIntensityRange pir : patchIntensityRanges) {
if (pir.patch.getBoundingBox().intersects(box))
map(new PatchTransform(pir), x0, y0, mapIntensities(pir, min, max), sp);
}
stack.addSlice(r + ", " + c, sp);
if (null == imp && stack.getSize() > 1) {
imp = new ImagePlus("tiles", stack);
imp.show();
}
if (null != imp) {
imp.setSlice(stack.getSize());
imp.updateAndDraw();
}
}
}
if (null == imp) {
// single-slice, non-StackWindow
new ImagePlus("tiles", stack).show();
}
}
use of ij.ImagePlus in project TrakEM2 by trakem2.
the class MatchIntensities method run.
/**
* @param layers
* @param radius
* @param scale
* @param numCoefficients
* @param lambda1
* @param lambda2
* @param neighborWeight
* @param roi
*/
public <M extends Model<M> & Affine1D<M>> void run(final List<Layer> layers, final int radius, final double scale, final int numCoefficients, final double lambda1, final double lambda2, final double neighborWeight, final Rectangle roi) throws InterruptedException, ExecutionException {
final int firstLayerIndex = layerset.getLayerIndex(layers.get(0).getId());
final int lastLayerIndex = layerset.getLayerIndex(layers.get(layers.size() - 1).getId());
// final PointMatchFilter filter = new RansacRegressionFilter();
final PointMatchFilter filter = new RansacRegressionReduceFilter();
/* collect patches */
Utils.log("Collecting patches ... ");
final ArrayList<Patch> patches = new ArrayList<Patch>();
for (final Layer layer : layers) patches.addAll((Collection) layer.getDisplayables(Patch.class, roi));
/* delete existing intensity coefficients */
Utils.log("Clearing existing intensity maps ... ");
for (final Patch p : patches) p.clearIntensityMap();
/* generate coefficient tiles for all patches
* TODO consider offering alternative models */
final HashMap<Patch, ArrayList<Tile<? extends M>>> coefficientsTiles = (HashMap) generateCoefficientsTiles(patches, new InterpolatedAffineModel1D<InterpolatedAffineModel1D<AffineModel1D, TranslationModel1D>, IdentityModel>(new InterpolatedAffineModel1D<AffineModel1D, TranslationModel1D>(new AffineModel1D(), new TranslationModel1D(), lambda1), new IdentityModel(), lambda2), numCoefficients * numCoefficients);
/* completed patches */
final HashSet<Patch> completedPatches = new HashSet<Patch>();
/* collect patch pairs */
Utils.log("Collecting patch pairs ... ");
final ArrayList<ValuePair<Patch, Patch>> patchPairs = new ArrayList<ValuePair<Patch, Patch>>();
for (final Patch p1 : patches) {
completedPatches.add(p1);
final Rectangle box1 = p1.getBoundingBox().intersection(roi);
final ArrayList<Patch> p2s = new ArrayList<Patch>();
/* across adjacent layers */
final int layerIndex = layerset.getLayerIndex(p1.getLayer().getId());
for (int i = Math.max(firstLayerIndex, layerIndex - radius); i <= Math.min(lastLayerIndex, layerIndex + radius); ++i) {
final Layer layer = layerset.getLayer(i);
if (layer != null)
p2s.addAll((Collection) layer.getDisplayables(Patch.class, box1));
}
for (final Patch p2 : p2s) {
/*
* if this patch had been processed earlier, all matches are
* already in
*/
if (completedPatches.contains(p2))
continue;
patchPairs.add(new ValuePair<Patch, Patch>(p1, p2));
}
}
final int numThreads = Integer.parseInt(layerset.getProperty("n_mipmap_threads", Integer.toString(Runtime.getRuntime().availableProcessors())));
Utils.log("Matching intensities using " + numThreads + " threads ... ");
final ExecutorService exec = Executors.newFixedThreadPool(numThreads);
final ArrayList<Future<?>> futures = new ArrayList<Future<?>>();
for (final ValuePair<Patch, Patch> patchPair : patchPairs) {
futures.add(exec.submit(new Matcher(roi, patchPair, (HashMap) coefficientsTiles, filter, scale, numCoefficients)));
}
for (final Future<?> future : futures) future.get();
/* connect tiles within patches */
Utils.log("Connecting coefficient tiles in the same patch ... ");
for (final Patch p1 : completedPatches) {
/* get the coefficient tiles */
final ArrayList<Tile<? extends M>> p1CoefficientsTiles = coefficientsTiles.get(p1);
for (int y = 1; y < numCoefficients; ++y) {
final int yr = numCoefficients * y;
final int yr1 = yr - numCoefficients;
for (int x = 0; x < numCoefficients; ++x) {
identityConnect(p1CoefficientsTiles.get(yr1 + x), p1CoefficientsTiles.get(yr + x), neighborWeight);
}
}
for (int y = 0; y < numCoefficients; ++y) {
final int yr = numCoefficients * y;
for (int x = 1; x < numCoefficients; ++x) {
final int yrx = yr + x;
identityConnect(p1CoefficientsTiles.get(yrx), p1CoefficientsTiles.get(yrx - 1), neighborWeight);
}
}
}
/* optimize */
Utils.log("Optimizing ... ");
final TileConfiguration tc = new TileConfiguration();
for (final ArrayList<Tile<? extends M>> coefficients : coefficientsTiles.values()) {
// for ( final Tile< ? > t : coefficients )
// if ( t.getMatches().size() == 0 )
// IJ.log( "bang" );
tc.addTiles(coefficients);
}
try {
tc.optimize(0.01f, iterations, iterations, 0.75f);
} catch (final NotEnoughDataPointsException e) {
// TODO Auto-generated catch block
e.printStackTrace();
} catch (final IllDefinedDataPointsException e) {
// TODO Auto-generated catch block
e.printStackTrace();
}
/* save coefficients */
final double[] ab = new double[2];
final FSLoader loader = (FSLoader) layerset.getProject().getLoader();
final String itsDir = loader.getUNUIdFolder() + "trakem2.its/";
for (final Entry<Patch, ArrayList<Tile<? extends M>>> entry : coefficientsTiles.entrySet()) {
final FloatProcessor as = new FloatProcessor(numCoefficients, numCoefficients);
final FloatProcessor bs = new FloatProcessor(numCoefficients, numCoefficients);
final Patch p = entry.getKey();
final double min = p.getMin();
final double max = p.getMax();
final ArrayList<Tile<? extends M>> tiles = entry.getValue();
for (int i = 0; i < numCoefficients * numCoefficients; ++i) {
final Tile<? extends M> t = tiles.get(i);
final Affine1D<?> affine = t.getModel();
affine.toArray(ab);
/* coefficients mapping into existing [min, max] */
as.setf(i, (float) ab[0]);
bs.setf(i, (float) ((max - min) * ab[1] + min - ab[0] * min));
}
final ImageStack coefficientsStack = new ImageStack(numCoefficients, numCoefficients);
coefficientsStack.addSlice(as);
coefficientsStack.addSlice(bs);
final String itsPath = itsDir + FSLoader.createIdPath(Long.toString(p.getId()), "it", ".tif");
new File(itsPath).getParentFile().mkdirs();
IJ.saveAs(new ImagePlus("", coefficientsStack), "tif", itsPath);
}
/* update mipmaps */
for (final Patch p : patches) p.getProject().getLoader().decacheImagePlus(p.getId());
final ArrayList<Future<Boolean>> mipmapFutures = new ArrayList<Future<Boolean>>();
for (final Patch p : patches) mipmapFutures.add(p.updateMipMaps());
for (final Future<Boolean> f : mipmapFutures) f.get();
Utils.log("Matching intensities done.");
}
use of ij.ImagePlus in project TrakEM2 by trakem2.
the class Render method main.
public static final void main(final String... args) {
new ImageJ();
final double scale = 0.05;
// final double scale = 1;
final int numCoefficients = 4;
final Project project = Project.openFSProject("/home/saalfeld/tmp/bock-lens-correction/subproject.xml", false);
final Layer layer = project.getRootLayerSet().getLayer(0);
final ArrayList<Patch> patches = (ArrayList) layer.getDisplayables(Patch.class);
final Patch patch1 = patches.get(0);
final Patch patch2 = patches.get(2);
final Rectangle box = patch1.getBoundingBox().intersection(patch2.getBoundingBox());
// final Rectangle box = patch1.getBoundingBox();
final int w = (int) (box.width * scale + 0.5);
final int h = (int) (box.height * scale + 0.5);
final FloatProcessor pixels1 = new FloatProcessor(w, h);
final FloatProcessor weights1 = new FloatProcessor(w, h);
final ColorProcessor coefficients1 = new ColorProcessor(w, h);
final FloatProcessor pixels2 = new FloatProcessor(w, h);
final FloatProcessor weights2 = new FloatProcessor(w, h);
final ColorProcessor coefficients2 = new ColorProcessor(w, h);
render(patch1, numCoefficients, numCoefficients, pixels1, weights1, coefficients1, box.x, box.y, scale);
render(patch2, numCoefficients, numCoefficients, pixels2, weights2, coefficients2, box.x, box.y, scale);
final ImageStack stack = new ImageStack(w, h);
stack.addSlice(pixels1);
stack.addSlice(pixels2);
stack.addSlice(weights1);
stack.addSlice(weights2);
stack.addSlice(coefficients1.convertToFloatProcessor());
stack.addSlice(coefficients2.convertToFloatProcessor());
new ImagePlus("", stack).show();
}
Aggregations