?  search:
summary | log | graphiclog1 | graphiclog2 | commit | commitdiff | tree | refs | edit | fork
blame | history | raw | HEAD
Corrected several long-standing issues with SIFT-linear and
[trakem2.git] / TrakEM2_ / src / main / java / mpicbg / trakem2 / align / AlignLayersTask.java
blob8b5be796e794989f4c6daf690c72a492d37be079
1 /**
2 * License: GPL
3 *
4 * This program is free software; you can redistribute it and/or
5 * modify it under the terms of the GNU General Public License 2
6 * as published by the Free Software Foundation.
7 *
8 * This program is distributed in the hope that it will be useful,
9 * but WITHOUT ANY WARRANTY; without even the implied warranty of
10 * MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the
11 * GNU General Public License for more details.
12 *
13 * You should have received a copy of the GNU General Public License
14 * along with this program; if not, write to the Free Software
15 * Foundation, Inc., 59 Temple Place - Suite 330, Boston, MA 02111-1307, USA.
16 */
17 package mpicbg.trakem2.align;
18
19 import ij.IJ;
20 import ij.ImagePlus;
21 import ij.gui.GenericDialog;
22 import ij.process.ImageProcessor;
23 import ini.trakem2.display.Display;
24 import ini.trakem2.display.Layer;
25 import ini.trakem2.display.LayerSet;
26 import ini.trakem2.display.Patch;
27 import ini.trakem2.persistence.Loader;
28 import ini.trakem2.utils.Bureaucrat;
29 import ini.trakem2.utils.Filter;
30 import ini.trakem2.utils.IJError;
31 import ini.trakem2.utils.Utils;
32 import ini.trakem2.utils.Worker;
33
34 import java.awt.Choice;
35 import java.awt.Rectangle;
36 import java.awt.event.ItemEvent;
37 import java.awt.event.ItemListener;
38 import java.awt.geom.AffineTransform;
39 import java.util.ArrayList;
40 import java.util.Collection;
41 import java.util.Iterator;
42 import java.util.List;
43 import java.util.Vector;
44 import java.util.concurrent.Callable;
45 import java.util.concurrent.ExecutorService;
46 import java.util.concurrent.Future;
47
48 import mpicbg.ij.FeatureTransform;
49 import mpicbg.ij.SIFT;
50 import mpicbg.ij.util.Util;
51 import mpicbg.imagefeatures.Feature;
52 import mpicbg.imagefeatures.FloatArray2DSIFT;
53 import mpicbg.models.AbstractAffineModel2D;
54 import mpicbg.models.AffineModel2D;
55 import mpicbg.models.IllDefinedDataPointsException;
56 import mpicbg.models.NotEnoughDataPointsException;
57 import mpicbg.models.Point;
58 import mpicbg.models.PointMatch;
59 import mpicbg.models.SimilarityModel2D;
60 import mpicbg.models.Transforms;
61 import mpicbg.trakem2.transform.CoordinateTransform;
62 import mpicbg.trakem2.transform.CoordinateTransformList;
63 import mpicbg.trakem2.transform.RigidModel2D;
64 import mpicbg.trakem2.transform.TranslationModel2D;
65 import bunwarpj.Transformation;
66 import bunwarpj.bUnwarpJ_;
67 import bunwarpj.trakem2.transform.CubicBSplineTransform;
68
69 /**
70 * Register a range of layers using linear or non-linear transformations.
71 *
72 * @author Stephan Saalfeld <saalfeld@mpi-cbg.de> and Ignacio Arganda <ignacio.arganda@gmail.com>
73 * @version 0.1a
74 */
75 final public class AlignLayersTask
76 {
77 static protected int LINEAR = 0, BUNWARPJ = 1, ELASTIC = 2;
78 static protected int mode = LINEAR;
79 final static String[] modeStrings = new String[]{
80 "least squares (linear feature correspondences)",
81 "bUnwarpJ (non-linear cubic B-Splines)",
82 "elastic (non-linear block correspondences)" };
83
84 static protected boolean propagateTransformBefore = false;
85 static protected boolean propagateTransformAfter = false;
86 static protected bunwarpj.Param elasticParam = new bunwarpj.Param();
87
88 private AlignLayersTask(){}
89 final static public Bureaucrat alignLayersTask ( final Layer l )
90 {
91 return alignLayersTask( l, null );
92 }
93
94 final static public Bureaucrat alignLayersTask ( final Layer l, final Rectangle fov )
95 {
96 final Worker worker = new Worker("Aligning layers", false, true) {
97 @Override
98 public void run() {
99 startedWorking();
100 try {
101 alignLayers(l, fov);
102 } catch (final Throwable e) {
103 IJError.print(e);
104 } finally {
105 finishedWorking();
106 }
107 }
108 };
109 return Bureaucrat.createAndStart(worker, l.getProject());
110 }
111
112 final static public void alignLayers( final Layer l ) throws Exception
113 {
114 alignLayers(l, null);
115 }
116
117 final static public void alignLayers( final Layer l, final Rectangle fov ) throws Exception
118 {
119 final List< Layer > layers = l.getParent().getLayers();
120 final String[] layerTitles = new String[ layers.size() ];
121 final String[] noneAndLayerTitles = new String[ layers.size() + 1 ];
122 noneAndLayerTitles[ 0 ] = "*None*";
123 for ( int i = 0; i < layers.size(); ++i )
124 {
125 layerTitles[ i ] = l.getProject().findLayerThing(layers.get( i )).toString();
126 noneAndLayerTitles[ i + 1 ] = layerTitles[ i ];
127 }
128
129 //Param p = Align.param;
130 //Align.param.sift.maxOctaveSize = 1024;
131
132 final GenericDialog gd = new GenericDialog( "Align Layers" );
133
134 gd.addChoice( "mode :", modeStrings, modeStrings[ LINEAR ] );
135
136 gd.addMessage( "Layer Range:" );
137 final int sel = l.getParent().indexOf(l);
138 gd.addChoice( "first :", layerTitles, layerTitles[ sel ] );
139 gd.addChoice( "reference :", noneAndLayerTitles, layerTitles[ sel ] );
140 gd.addChoice( "last :", layerTitles, layerTitles[ sel ] );
141 gd.addStringField("Use only images whose title matches:", "", 30);
142 gd.addCheckbox("Use visible images only", true);
143 gd.addCheckbox( "propagate transform to first layer", propagateTransformBefore );
144 gd.addCheckbox( "propagate transform to last layer", propagateTransformAfter );
145 final Vector<?> v = gd.getChoices();
146 final Choice cstart = (Choice) v.get(v.size() -3);
147 final Choice cref = (Choice) v.get(v.size() -2);
148 final Choice cend = (Choice) v.get(v.size() -1);
149 final ItemListener startEndListener = new ItemListener() {
150 @Override
151 public void itemStateChanged(final ItemEvent ie) {
152 final int startIndex = cstart.getSelectedIndex();
153 final int endIndex = cend.getSelectedIndex();
154 final int refIndex = cref.getSelectedIndex();
155
156 final int min = Math.min(startIndex, endIndex);
157 final int max = Math.max(startIndex, endIndex);
158
159 if (cref.getSelectedIndex() > 0) {
160 cref.select( Math.min( max + 1, Math.max( min + 1, refIndex ) ) );
161 }
162 }
163 };
164 cstart.addItemListener(startEndListener);
165 cend.addItemListener(startEndListener);
166 cref.addItemListener(new ItemListener() {
167 @Override
168 public void itemStateChanged(final ItemEvent ie) {
169 final int startIndex = cstart.getSelectedIndex();
170 final int endIndex = cend.getSelectedIndex();
171 final int refIndex = cref.getSelectedIndex() - 1;
172
173 final int min = Math.min(startIndex, endIndex);
174 final int max = Math.max(startIndex, endIndex);
175
176 if (refIndex >= 0) {
177 if (refIndex < min) {
178 if ( startIndex <= endIndex )
179 cstart.select(refIndex);
180 else
181 cend.select(refIndex);
182 }
183 else if (refIndex > max) {
184 if ( startIndex <= endIndex )
185 cend.select(refIndex);
186 else
187 cstart.select(refIndex);
188 }
189 }
190 }
191 });
192
193 gd.showDialog();
194 if ( gd.wasCanceled() ) return;
195
196 mode = gd.getNextChoiceIndex();
197
198 final int first = gd.getNextChoiceIndex();
199 final int ref = gd.getNextChoiceIndex() - 1;
200 final int last = gd.getNextChoiceIndex();
201
202 final String toMatch1 = gd.getNextString().trim();
203 final String toMatch2 = 0 == toMatch1.length() ? null : ".*" + toMatch1 + ".*";
204 final boolean visibleOnly = gd.getNextBoolean();
205
206 propagateTransformBefore = gd.getNextBoolean();
207 propagateTransformAfter = gd.getNextBoolean();
208
209 final Filter<Patch> filter = new Filter<Patch>() {
210 @Override
211 public final boolean accept(final Patch patch) {
212 if (visibleOnly && !patch.isVisible()) return false;
213 if (null != toMatch2 && !patch.getTitle().matches(toMatch2)) return false;
214 return true;
215 }
216 };
217
218 if ( mode == ELASTIC )
219 new ElasticLayerAlignment().exec( l.getParent(), first, last, ref, propagateTransformBefore, propagateTransformAfter, fov, filter );
220 else if ( mode == LINEAR )
221 new RegularizedAffineLayerAlignment().exec( l.getParent(), first, last, ref, propagateTransformBefore, propagateTransformAfter, fov, filter );
222 else
223 {
224 final GenericDialog gd2 = new GenericDialog( "Align Layers" );
225
226 Align.param.addFields( gd2 );
227
228 gd2.showDialog();
229 if ( gd2.wasCanceled() ) return;
230
231 Align.param.readFields( gd2 );
232
233 if ( mode == BUNWARPJ && !elasticParam.showDialog() ) return;
234
235 // From ref to first:
236 if (ref - first > 0)
237 {
238 if ( mode == BUNWARPJ )
239 alignLayersNonLinearlyJob(l.getParent(), ref, first, propagateTransformBefore, fov, filter);
240 else
241 alignLayersLinearlyJob(l.getParent(), ref, first, propagateTransformBefore, fov, filter);
242 }
243 // From ref to last:
244 if (last - ref > 0)
245 {
246 if ( mode == BUNWARPJ )
247 alignLayersNonLinearlyJob(l.getParent(), ref, last, propagateTransformAfter, fov, filter);
248 else
249 alignLayersLinearlyJob(l.getParent(), ref, last, propagateTransformAfter, fov, filter);
250 }
251 }
252 }
253
254
255 final static public void alignLayersLinearlyJob( final LayerSet layerSet, final int first, final int last,
256 final boolean propagateTransform, final Rectangle fov, final Filter<Patch> filter )
257 {
258 final List< Layer > layerRange = layerSet.getLayers(first, last); // will reverse order if necessary
259
260 final Align.Param p = Align.param.clone();
261 // find the first non-empty layer, and remove all empty layers
262 Rectangle box = fov;
263 for (final Iterator<Layer> it = layerRange.iterator(); it.hasNext(); ) {
264 final Layer la = it.next();
265 if (!la.contains(Patch.class, true)) {
266 it.remove();
267 continue;
268 }
269 if (null == box) {
270 // The first layer:
271 box = la.getMinimalBoundingBox(Patch.class, true); // Only for visible patches
272 }
273 }
274 if (0 == layerRange.size()) {
275 Utils.log("All layers in range are empty!");
276 return;
277 }
278
279 /* do not work if there is only one layer selected */
280 if ( layerRange.size() < 2 ) return;
281
282 final float scale = Math.min( 1.0f, Math.min( ( float )p.sift.maxOctaveSize / ( float )box.width, ( float )p.sift.maxOctaveSize / ( float )box.height ) );
283 p.maxEpsilon *= scale;
284 p.identityTolerance *= scale;
285
286 //Utils.log2("scale: " + scale + " maxOctaveSize: " + p.sift.maxOctaveSize + " box: " + box.width + "," + box.height);
287
288 final FloatArray2DSIFT sift = new FloatArray2DSIFT( p.sift );
289 final SIFT ijSIFT = new SIFT( sift );
290
291 Rectangle box1 = fov;
292 Rectangle box2 = fov;
293 final Collection< Feature > features1 = new ArrayList< Feature >();
294 final Collection< Feature > features2 = new ArrayList< Feature >();
295 final List< PointMatch > candidates = new ArrayList< PointMatch >();
296 final List< PointMatch > inliers = new ArrayList< PointMatch >();
297
298 final AffineTransform a = new AffineTransform();
299
300 int s = 0;
301 for ( final Layer layer : layerRange )
302 {
303 if ( Thread.currentThread().isInterrupted() ) return;
304
305 final long t0 = System.currentTimeMillis();
306
307 features1.clear();
308 features1.addAll( features2 );
309 features2.clear();
310
311 final Rectangle box3 = layer.getMinimalBoundingBox( Patch.class, true );
312
313 if ( box3 == null || ( box3.width == 0 && box3.height == 0 ) ) continue; // skipping empty layer
314
315 box1 = null == fov ? box2 : fov;
316 box2 = null == fov ? box3 : fov;
317
318 final List<Patch> patches = layer.getAll(Patch.class);
319 if (null != filter) {
320 for (final Iterator<Patch> it = patches.iterator(); it.hasNext(); ) {
321 if (!filter.accept(it.next())) it.remove();
322 }
323 }
324
325 final ImageProcessor flatImage = layer.getProject().getLoader().getFlatImage( layer, box2, scale, 0xffffffff, ImagePlus.GRAY8, Patch.class, patches, true ).getProcessor();
326
327 ijSIFT.extractFeatures(
328 flatImage,
329 features2 );
330 IJ.log( features2.size() + " features extracted in layer \"" + layer.getTitle() + "\" (took " + ( System.currentTimeMillis() - t0 ) + " ms)." );
331
332 if ( features1.size() > 0 )
333 {
334 final long t1 = System.currentTimeMillis();
335
336 candidates.clear();
337
338 FeatureTransform.matchFeatures(
339 features2,
340 features1,
341 candidates,
342 p.rod );
343
344 final AbstractAffineModel2D< ? > model;
345 switch ( p.expectedModelIndex )
346 {
347 case 0:
348 model = new TranslationModel2D();
349 break;
350 case 1:
351 model = new RigidModel2D();
352 break;
353 case 2:
354 model = new SimilarityModel2D();
355 break;
356 case 3:
357 model = new AffineModel2D();
358 break;
359 default:
360 return;
361 }
362
363 final AbstractAffineModel2D< ? > desiredModel;
364 switch ( p.desiredModelIndex )
365 {
366 case 0:
367 desiredModel = new TranslationModel2D();
368 break;
369 case 1:
370 desiredModel = new RigidModel2D();
371 break;
372 case 2:
373 desiredModel = new SimilarityModel2D();
374 break;
375 case 3:
376 desiredModel = new AffineModel2D();
377 break;
378 default:
379 return;
380 }
381
382
383
384 boolean modelFound;
385 boolean again = false;
386 try
387 {
388 do
389 {
390 again = false;
391 modelFound = model.filterRansac(
392 candidates,
393 inliers,
394 1000,
395 p.maxEpsilon,
396 p.minInlierRatio,
397 p.minNumInliers,
398 3 );
399 if ( modelFound && p.rejectIdentity )
400 {
401 final ArrayList< Point > points = new ArrayList< Point >();
402 PointMatch.sourcePoints( inliers, points );
403 if ( Transforms.isIdentity( model, points, p.identityTolerance ) )
404 {
405 IJ.log( "Identity transform for " + inliers.size() + " matches rejected." );
406 candidates.removeAll( inliers );
407 inliers.clear();
408 again = true;
409 }
410 }
411 }
412 while ( again );
413
414 if ( modelFound )
415 desiredModel.fit( inliers );
416 }
417 catch ( final NotEnoughDataPointsException e )
418 {
419 modelFound = false;
420 }
421 catch ( final IllDefinedDataPointsException e )
422 {
423 modelFound = false;
424 }
425
426 if ( Thread.currentThread().isInterrupted() ) return;
427
428 if ( modelFound )
429 {
430 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" );
431 final AffineTransform b = new AffineTransform();
432 b.translate( box1.x, box1.y );
433 b.scale( 1.0f / scale, 1.0f / scale );
434 b.concatenate( desiredModel.createAffine() );
435 b.scale( scale, scale );
436 b.translate( -box2.x, -box2.y);
437
438 a.concatenate( b );
439 AlignTask.transformPatchesAndVectorData(patches, a);
440 Display.repaint( layer );
441 }
442 else
443 {
444 IJ.log( "No model found for layer \"" + layer.getTitle() + "\" and its predecessor:\n correspondence candidates " + candidates.size() + "\n took " + ( System.currentTimeMillis() - s ) + " ms" );
445 a.setToIdentity();
446 }
447 }
448 IJ.showProgress( ++s, layerRange.size() );
449 }
450
451 if ( Thread.currentThread().isInterrupted() ) return;
452
453 if ( propagateTransform )
454 {
455 if ( last > first && last < layerSet.size() - 2 )
456 for ( final Layer la : layerSet.getLayers( last + 1, layerSet.size() - 1 ) ) {
457 if (Thread.currentThread().isInterrupted()) return;
458 AlignTask.transformPatchesAndVectorData( la, a );
459 }
460 else if ( first > last && last > 0 )
461 for ( final Layer la : layerSet.getLayers( 0, last - 1 ) ) {
462 if (Thread.currentThread().isInterrupted()) return;
463 AlignTask.transformPatchesAndVectorData( la, a );
464 }
465 }
466 }
467
468 final static public void alignLayersNonLinearlyJob( final LayerSet layerSet, final int first, final int last,
469 final boolean propagateTransform, final Rectangle fov, final Filter<Patch> filter )
470 {
471 final List< Layer > layerRange = layerSet.getLayers(first, last); // will reverse order if necessary
472
473 final Align.Param p = Align.param.clone();
474
475 // Remove all empty layers
476 for (final Iterator<Layer> it = layerRange.iterator(); it.hasNext(); ) {
477 if (!it.next().contains(Patch.class, true)) {
478 it.remove();
479 }
480 }
481 if (0 == layerRange.size()) {
482 Utils.log("No layers in range show any images!");
483 return;
484 }
485
486 /* do not work if there is only one layer selected */
487 if ( layerRange.size() < 2 ) return;
488
489 final List<Patch> all = new ArrayList<Patch>();
490 for (final Layer la : layerRange) {
491 for (final Patch patch : la.getAll(Patch.class)) {
492 if (null != filter && !filter.accept(patch)) continue;
493 all.add(patch);
494 }
495 }
496
497 AlignTask.transformPatchesAndVectorData(all, new Runnable() { @Override
498 public void run() {
499
500 /////
501
502 final Loader loader = layerSet.getProject().getLoader();
503
504 // Not concurrent safe! So two copies, one per layer and Thread:
505 final SIFT ijSIFT1 = new SIFT( new FloatArray2DSIFT( p.sift ) );
506 final SIFT ijSIFT2 = new SIFT( new FloatArray2DSIFT( p.sift ) );
507
508 final Collection< Feature > features1 = new ArrayList< Feature >();
509 final Collection< Feature > features2 = new ArrayList< Feature >();
510 final List< PointMatch > candidates = new ArrayList< PointMatch >();
511 final List< PointMatch > inliers = new ArrayList< PointMatch >();
512
513 final int n_proc = Runtime.getRuntime().availableProcessors() > 1 ? 2 : 1;
514 final ExecutorService exec = Utils.newFixedThreadPool(n_proc, "alignLayersNonLinearly");
515
516 List<Patch> previousPatches = null;
517
518 int s = 0;
519 for ( int i = 1; i < layerRange.size(); ++i )
520 {
521 if ( Thread.currentThread().isInterrupted() ) break;
522
523 final Layer layer1 = layerRange.get( i - 1 );
524 final Layer layer2 = layerRange.get( i );
525
526 final long t0 = System.currentTimeMillis();
527
528 features1.clear();
529 features2.clear();
530
531 final Rectangle box1 = null == fov ? layer1.getMinimalBoundingBox( Patch.class, true ) : fov;
532 final Rectangle box2 = null == fov ? layer2.getMinimalBoundingBox( Patch.class, true ) : fov;
533
534 /* calculate the common scale factor for both flat images */
535 final float scale = Math.min( 1.0f, ( float )p.sift.maxOctaveSize / ( float )Math.max( box1.width, Math.max( box1.height, Math.max( box2.width, box2.height ) ) ) );
536
537 final List<Patch> patches1;
538 if (null == previousPatches) {
539 patches1 = layer1.getAll(Patch.class);
540 if (null != filter) {
541 for (final Iterator<Patch> it = patches1.iterator(); it.hasNext(); ) {
542 if (!filter.accept(it.next())) it.remove();
543 }
544 }
545 } else {
546 patches1 = previousPatches;
547 }
548
549 final List<Patch> patches2 = layer2.getAll(Patch.class);
550 if (null != filter) {
551 for (final Iterator<Patch> it = patches2.iterator(); it.hasNext(); ) {
552 if (!filter.accept(it.next())) it.remove();
553 }
554 }
555
556 final Future<ImageProcessor> fu1 = exec.submit(new Callable<ImageProcessor>() {
557 @Override
558 public ImageProcessor call() {
559 final ImageProcessor ip1 = loader.getFlatImage( layer1, box1, scale, 0xffffffff, ImagePlus.GRAY8, Patch.class, patches1, true ).getProcessor();
560 ijSIFT1.extractFeatures(
561 ip1,
562 features1 );
563 Utils.log( features1.size() + " features extracted in layer \"" + layer1.getTitle() + "\" (took " + ( System.currentTimeMillis() - t0 ) + " ms)." );
564 return ip1;
565 }});
566 final Future<ImageProcessor> fu2 = exec.submit(new Callable<ImageProcessor>() {
567 @Override
568 public ImageProcessor call() {
569 final ImageProcessor ip2 = loader.getFlatImage( layer2, box2, scale, 0xffffffff, ImagePlus.GRAY8, Patch.class, patches2, true ).getProcessor();
570 ijSIFT2.extractFeatures(
571 ip2,
572 features2 );
573 Utils.log( features2.size() + " features extracted in layer \"" + layer2.getTitle() + "\" (took " + ( System.currentTimeMillis() - t0 ) + " ms)." );
574 return ip2;
575 }});
576
577 final ImageProcessor ip1, ip2;
578 try {
579 ip1 = fu1.get();
580 ip2 = fu2.get();
581 } catch (final Exception e) {
582 IJError.print(e);
583 return;
584 }
585
586 if ( features1.size() > 0 && features2.size() > 0 )
587 {
588 final long t1 = System.currentTimeMillis();
589
590 candidates.clear();
591
592 FeatureTransform.matchFeatures(
593 features2,
594 features1,
595 candidates,
596 p.rod );
597
598 final AbstractAffineModel2D< ? > model;
599 switch ( p.expectedModelIndex )
600 {
601 case 0:
602 model = new TranslationModel2D();
603 break;
604 case 1:
605 model = new RigidModel2D();
606 break;
607 case 2:
608 model = new SimilarityModel2D();
609 break;
610 case 3:
611 model = new AffineModel2D();
612 break;
613 default:
614 return;
615 }
616
617 boolean modelFound;
618 boolean again = false;
619 try
620 {
621 do
622 {
623 again = false;
624 modelFound = model.filterRansac(
625 candidates,
626 inliers,
627 1000,
628 p.maxEpsilon,
629 p.minInlierRatio,
630 p.minNumInliers,
631 3 );
632 if ( modelFound && p.rejectIdentity )
633 {
634 final ArrayList< Point > points = new ArrayList< Point >();
635 PointMatch.sourcePoints( inliers, points );
636 if ( Transforms.isIdentity( model, points, p.identityTolerance ) )
637 {
638 IJ.log( "Identity transform for " + inliers.size() + " matches rejected." );
639 candidates.removeAll( inliers );
640 inliers.clear();
641 again = true;
642 }
643 }
644 }
645 while ( again );
646 }
647 catch ( final NotEnoughDataPointsException e )
648 {
649 modelFound = false;
650 }
651
652 if ( modelFound )
653 {
654 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" );
655
656 final ImagePlus imp1 = new ImagePlus("target", ip1);
657 final ImagePlus imp2 = new ImagePlus("source", ip2);
658
659 final List< Point > sourcePoints = new ArrayList<Point>();
660 final List< Point > targetPoints = new ArrayList<Point>();
661
662 PointMatch.sourcePoints( inliers, sourcePoints );
663 PointMatch.targetPoints( inliers, targetPoints );
664
665 imp2.setRoi( Util.pointsToPointRoi(sourcePoints) );
666 imp1.setRoi( Util.pointsToPointRoi(targetPoints) );
667
668 final ImageProcessor mask1 = ip1.duplicate();
669 mask1.threshold(1);
670 final ImageProcessor mask2 = ip2.duplicate();
671 mask2.threshold(1);
672
673 final Transformation warp = bUnwarpJ_.computeTransformationBatch(imp2, imp1, mask2, mask1, elasticParam);
674
675 final CubicBSplineTransform transf = new CubicBSplineTransform();
676 transf.set(warp.getIntervals(), warp.getDirectDeformationCoefficientsX(), warp.getDirectDeformationCoefficientsY(),
677 imp2.getWidth(), imp2.getHeight());
678
679 final ArrayList<Future<?>> fus = new ArrayList<Future<?>>();
680
681 // Transform desired patches only
682 for ( final Patch patch : patches2 )
683 {
684 try {
685 final Rectangle pbox = patch.getCoordinateTransformBoundingBox();
686 final AffineTransform at = patch.getAffineTransform();
687 final AffineTransform pat = new AffineTransform();
688 pat.scale( scale, scale );
689 pat.translate( -box2.x, -box2.y );
690 pat.concatenate( at );
691 pat.translate( -pbox.x, -pbox.y );
692
693
694 final mpicbg.trakem2.transform.AffineModel2D toWorld = new mpicbg.trakem2.transform.AffineModel2D();
695 toWorld.set( pat );
696
697 final CoordinateTransformList< CoordinateTransform > ctl = new CoordinateTransformList< CoordinateTransform >();
698 // move the patch into the global space where bUnwarpJ calculated the transformation
699 ctl.add( toWorld );
700 // Apply non-linear transformation
701 ctl.add( transf );
702 // move it back
703 ctl.add( toWorld.createInverse() );
704
705 patch.appendCoordinateTransform( ctl );
706
707 fus.add(patch.updateMipMaps());
708
709 // Compensate for offset between boxes
710 final AffineTransform offset = new AffineTransform();
711 offset.translate( box1.x - box2.x , box1.y - box2.y );
712 offset.concatenate( at );
713 patch.setAffineTransform( offset );
714 } catch ( final Exception e ) {
715 e.printStackTrace();
716 }
717 }
718
719 // await regeneration of all mipmaps
720 Utils.wait(fus);
721
722 Display.repaint( layer2 );
723 }
724 else
725 IJ.log( "No model found for layer \"" + layer2.getTitle() + "\" and its predecessor:\n correspondence candidates " + candidates.size() + "\n took " + ( System.currentTimeMillis() - s ) + " ms" );
726 }
727 IJ.showProgress( ++s, layerRange.size() );
728
729 previousPatches = patches2; // for next iteration
730 }
731
732 exec.shutdown();
733
734 if (propagateTransform) Utils.log("Propagation not implemented yet for non-linear layer alignment.");
735
736 /* // CANNOT be done (at least not trivially:
737 * //an appropriate "scale" cannot be computed, and the box2 is part of the spline computation.
738 if ( propagateTransform && null != lastTransform )
739 {
740 for (final Layer la : l.getParent().getLayers(last > first ? last +1 : first -1, last > first ? l.getParent().size() -1 : 0)) {
741 // Transform visible patches only
742 final Rectangle box2 = la.getMinimalBoundingBox( Patch.class, true );
743 for ( final Displayable disp : la.getDisplayables( Patch.class, true ) )
744 {
745 // ...
746 }
747 }
748 }
749 */
750
751 }}); // end of transformPatchesAndVectorData
752 }
753 }
register/edit/analyze/3D model TB of images