Java tutorial
/* * To the extent possible under law, the Fiji developers have waived * all copyright and related or neighboring rights to this tutorial code. * * See the CC0 1.0 Universal license for details: * http://creativecommons.org/publicdomain/zero/1.0/ */ import org.apache.commons.lang3.ArrayUtils; import ij.gui.GenericDialog; import ij.process.ImageProcessor; //import mpicbg.ij.Mapping; //import mpicbg.ij.InverseTransformMapping; import mpicbg.ij.SIFT; import mpicbg.imagefeatures.*; import mpicbg.models.*; import ij.plugin.*; import ij.gui.*; import ij.*; import java.awt.*; import ij.measure.*; import java.util.ArrayList; import java.util.Arrays; import java.util.List; import java.util.Vector; import java.awt.Color; import java.awt.TextField; import java.awt.event.KeyEvent; import java.awt.event.KeyListener; import stitching.ImageInformation; import static stitching.CommonFunctions.LIN_BLEND; import static stitching.CommonFunctions.methodListCollection; import mpicbg.models.TranslationModel3D; import mpicbg.models.AffineModel3D; /** * Align and stitch two stacks of images assumed to be contiguous using automatically extracted robust landmark * correspondences. * @author Chloe Murtin <chloe.murtinl@gmail.com> and Carole Frindel <carole.frindel@creatis.insa-lyon.fr> * @version 0.0 */ public class SIFT_Volume_Stitching implements PlugIn, KeyListener { public static void main(String[] args) { // set the plugins.dir property to make the plugin appear in the Plugins menu Class<?> clazz = SIFT_Volume_Stitching.class; String url = clazz.getResource("/" + clazz.getName().replace('.', '/') + ".class").toString(); String pluginsDir = url.substring("file:".length(), url.length() - clazz.getName().length() - ".class".length()); System.setProperty("plugins.dir", pluginsDir); // start ImageJ new ImageJ(); // open the Clown sample ImagePlus image1 = IJ.openImage("https://www.creatis.insa-lyon.fr/~frindel/BackStack115.tif"); ImagePlus image2 = IJ.openImage("https://www.creatis.insa-lyon.fr/~frindel/FrontStack.tif"); image1.show(); image2.show(); // run the plugin IJ.runPlugIn(clazz.getName(), ""); } private List<Feature> fsf = new ArrayList<Feature>(); private List<Feature> fsb = new ArrayList<Feature>(); private ImagePlus firstImage; // Original Front Stack private ImagePlus secondImage; // Original Back Stack private ImagePlus impAlignedZYX; // Aligned Back Stack private ImagePlus fTemplate; private ImagePlus bTemplate; private ImagePlus channel2f; private ImagePlus channel2b; private ImagePlus channel3f; private ImagePlus channel3b; private ImageStack tempimpf; private ImageStack new_impf; private ImageStack tempimpb; private ImageStack new_impb; private ImagePlus channel2front; private ImagePlus channel2back; private ImagePlus channel3front; private ImagePlus channel3back; private boolean modelFound; private AffineModel3D BestModel3D = new AffineModel3D(); public boolean Reg3D = true; boolean template_bool = false; final static public String[] stitchingModelStrings = new String[] { "Front - Back", "Back - Front", "Left - Right", "Right - Left", "Top - Bottom", "Bottom - Top" }; public String stitchingMethod = "Front - Back"; final static public String[] OVMethod = new String[] { "Slice-by-Slice", "Block-by-Block", "Determined" }; public String myOVMethod = "Slice-by-Slice"; /** Fusion method*/ public String fusionMethod = methodListCollection[LIN_BLEND]; public double alpha = 1.5; static private class Param { public FloatArray2DSIFT.Param sift = new FloatArray2DSIFT.Param(); /**Closest/next closest neighbour distance ratio*/ public float rod = 0.92f; /** Maximal allowed alignment error in px*/ public float maxEpsilon = 25.0f; /** Inlier/candidates ratio*/ public float minInlierRatio = 0.05f; /** Implemeted transformation models for choice*/ public int modelIndex = 1; } static Param p = new Param(); /** * Main method of the plugin * @param args unused */ final public void run(final String args) { fsf.clear(); fsb.clear(); if (IJ.versionLessThan("1.41n")) return; /* Estimated Overlap size */ int ov = 100; /* Estimated start and end of overlap */ int startSlice = 100; int endSlice = 200; /* Image Choice Dialog Box */ Font myfont = new Font("SansSerif", Font.BOLD, 12); final int[] ids = WindowManager.getIDList(); if (ids == null || ids.length < 2) { IJ.showMessage("You should have at least two images open."); return; } final String[] titles = new String[ids.length]; final int[] sSizes = new int[ids.length]; for (int i = 0; i < ids.length; ++i) { titles[i] = (WindowManager.getImage(ids[i])).getTitle(); sSizes[i] = (WindowManager.getImage(ids[i])).getStackSize(); } final GenericDialog gd = new GenericDialog("2D-SIFT in 3D-Space Volume Stitching"); final GenericDialog gdChannels = new GenericDialog("Additional Channels"); final GenericDialog gdFiltered = new GenericDialog("Filtered Images"); gd.addMessage("* Image Selection", myfont); final String current = WindowManager.getCurrentImage().getTitle(); gd.addChoice("Front_Image", titles, current); gd.addChoice("Back_Image", titles, current.equals(titles[0]) ? titles[1] : titles[0]); gd.addChoice("Stitching Orientation", stitchingModelStrings, stitchingMethod); gd.addMessage("* Overlap Detection Method", myfont); gd.addChoice("Method Selection", OVMethod, myOVMethod); WindowManager.getCurrentImage().getStackSize(); gd.addNumericField("Split (Block-by-Block)", 5, 0, 0, ""); gd.addNumericField("Overlap Size (Slice-by-Slice)", 200, 0, 0, ""); gd.addNumericField("Start of Overlap (Determined)", startSlice, 0, 0, ""); gd.addNumericField("End of Overlap (Determined)", endSlice, 0, 0, ""); boolean cAll = false; gd.addMessage("* Additional Channels", myfont); gd.addCheckbox("More channels?", cAll); int it = 2; int MIP = 50; gd.addMessage("* 3D Registration Parameters", myfont); gd.addCheckbox("3D Registration", Reg3D); gd.addNumericField("Number_Of_Iterations", it, 0, 0, ""); gd.addNumericField("MIP size", MIP, 0, 0, "Slices"); float minSize = 3.5f; float maxSize = 14.5f; gd.addMessage("* SIFT Parameters", myfont); gd.addNumericField("Minimum Size of Structure :", minSize, 1, 0, "pixels"); gd.addNumericField("Maximum Size of Structure :", maxSize, 1, 0, "pixels"); gd.addMessage("* Filtered Image Selection", myfont); gd.addCheckbox("Compare Filtered Images?", template_bool); gd.addMessage("* Fusion Method", myfont); gd.addChoice("Method", methodListCollection, fusionMethod); /* ----- Window Channels -----*/ boolean c2 = false; boolean c3 = false; gdChannels.addCheckbox("Channel 2", c2); gdChannels.addChoice("Channel_2_Front_Image", titles, current); gdChannels.addChoice("Channel_2_Back_Image", titles, current.equals(titles[0]) ? titles[1] : titles[0]); gdChannels.addCheckbox("Channel 3", c3); gdChannels.addChoice("Channel_3_Front_Image", titles, current); gdChannels.addChoice("Channel_3_Back_Image", titles, current.equals(titles[0]) ? titles[1] : titles[0]); /* ----- -------------- -----*/ /* ------ Win Filtered ------*/ gdFiltered.addChoice("Front_Filtered_Image", titles, current); gdFiltered.addChoice("Back_Filtered_Image", titles, current.equals(titles[0]) ? titles[1] : titles[0]); /* ----- -------------- -----*/ gd.showDialog(); if (gd.wasCanceled()) return; /* Variables */ long start_time = System.currentTimeMillis(); firstImage = WindowManager.getImage(ids[gd.getNextChoiceIndex()]); secondImage = WindowManager.getImage(ids[gd.getNextChoiceIndex()]); ImagePlus impf = new ImagePlus(firstImage.getTitle(), firstImage.getStack()); impf.setCalibration(firstImage.getCalibration()); ImagePlus impb = new ImagePlus(secondImage.getTitle(), secondImage.getStack()); impb.setCalibration(secondImage.getCalibration()); if (impf.getCalibration().pixelWidth != impb.getCalibration().pixelWidth && impf.getCalibration().pixelHeight != impb.getCalibration().pixelHeight && impf.getCalibration().pixelDepth != impb.getCalibration().pixelDepth && impf.getCalibration().getUnit() == impb.getCalibration().getUnit()) { IJ.log("Calibration is different"); return; } Calibration copyCalibration = new Calibration(); copyCalibration = impf.getCalibration().copy(); stitchingMethod = stitchingModelStrings[gd.getNextChoiceIndex()]; if (stitchingMethod != "Front - Back") { impf = stackOrientation(impf, stitchingMethod); impb = stackOrientation(impb, stitchingMethod); } myOVMethod = OVMethod[gd.getNextChoiceIndex()]; int split = (int) gd.getNextNumber(); ov = (int) gd.getNextNumber(); startSlice = (int) gd.getNextNumber(); endSlice = (int) gd.getNextNumber(); int indb = impb.getStackSize(); cAll = gd.getNextBoolean(); Reg3D = gd.getNextBoolean(); it = (int) gd.getNextNumber(); MIP = (int) gd.getNextNumber(); minSize = (float) gd.getNextNumber(); maxSize = (float) gd.getNextNumber(); template_bool = gd.getNextBoolean(); fusionMethod = gd.getNextChoice(); /* ------- SIFT parameters ------- */ /*------------------------*/ IJ.log("Structures sizes, min:" + minSize + ", max:" + maxSize); int[] ImgSizef = impf.getDimensions(); int[] ImgSizeb = impb.getDimensions(); int[] AllSizes = ArrayUtils.addAll(Arrays.copyOfRange(ImgSizef, 0, 2), Arrays.copyOfRange(ImgSizeb, 0, 2)); int minOSize = AllSizes[0]; //IJ.log("Taille: "+AllSizes.length+" Contenu:"+AllSizes[0]+AllSizes[1]+AllSizes[2]+AllSizes[3]); for (int i = 1; i < AllSizes.length; i++) { if (AllSizes[i] < minOSize) { minOSize = AllSizes[i]; } } /*------------------------*/ p.sift.initialSigma = (float) (minSize / (2 * Math.sqrt(2))); p.sift.steps = (int) ((maxSize / minSize) + 1); p.sift.minOctaveSize = 64; p.sift.maxOctaveSize = (int) minOSize; IJ.log("Sift parameters, initial sigma:" + p.sift.initialSigma + "\n steps:" + p.sift.steps + " \n MinOctaveSize:" + p.sift.minOctaveSize + " \n MaxOctaveSize:" + p.sift.maxOctaveSize); /* Computations */ bTemplate = null; fTemplate = null; /* If use of a template exchange between template and real channel */ fTemplate = new ImagePlus(impf.getTitle(), impf.getStack()); bTemplate = new ImagePlus(impb.getTitle(), impb.getStack()); channel2f = null; channel2b = null; channel3f = null; channel3b = null; channel2f = null; channel2b = null; channel3f = null; channel3b = null; if (cAll == true && template_bool == true) { gdChannels.showDialog(); gd.setVisible(false); /* Variables */ c2 = gdChannels.getNextBoolean(); channel2f = WindowManager.getImage(ids[gdChannels.getNextChoiceIndex()]); channel2b = WindowManager.getImage(ids[gdChannels.getNextChoiceIndex()]); c3 = gdChannels.getNextBoolean(); channel3f = WindowManager.getImage(ids[gdChannels.getNextChoiceIndex()]); channel3b = WindowManager.getImage(ids[gdChannels.getNextChoiceIndex()]); if ((gdChannels.wasOKed())) { gdFiltered.showDialog(); /* Variables */ firstImage = WindowManager.getImage(ids[gdFiltered.getNextChoiceIndex()]); secondImage = WindowManager.getImage(ids[gdFiltered.getNextChoiceIndex()]); if (gdFiltered.wasOKed()) gd.setVisible(true); } } else if (cAll == false && template_bool == true) { gdFiltered.showDialog(); gd.setVisible(false); /* Variables */ firstImage = WindowManager.getImage(ids[gdFiltered.getNextChoiceIndex()]); secondImage = WindowManager.getImage(ids[gdFiltered.getNextChoiceIndex()]); if (gdFiltered.wasOKed()) { gd.setVisible(true); } } else if (cAll == true && template_bool == false) { gdChannels.showDialog(); gd.setVisible(false); /* Variables */ c2 = gdChannels.getNextBoolean(); channel2f = WindowManager.getImage(ids[gdChannels.getNextChoiceIndex()]); channel2b = WindowManager.getImage(ids[gdChannels.getNextChoiceIndex()]); c3 = gdChannels.getNextBoolean(); channel3f = WindowManager.getImage(ids[gdChannels.getNextChoiceIndex()]); channel3b = WindowManager.getImage(ids[gdChannels.getNextChoiceIndex()]); if (gdChannels.wasOKed()) { gd.setVisible(true); } } if (template_bool == true) { impf = new ImagePlus(firstImage.getTitle(), firstImage.getStack()); impf.setCalibration(firstImage.getCalibration()); impb = new ImagePlus(secondImage.getTitle(), secondImage.getStack()); impb.setCalibration(secondImage.getCalibration()); } if (cAll == true) { channel2front = new ImagePlus(channel2f.getTitle(), channel2f.getStack()); channel2front.setCalibration(channel2f.getCalibration()); channel2back = new ImagePlus(channel2b.getTitle(), channel2b.getStack()); channel2back.setCalibration(channel2b.getCalibration()); channel3front = new ImagePlus(channel3f.getTitle(), channel3f.getStack()); channel3front.setCalibration(channel3f.getCalibration()); channel3back = new ImagePlus(channel3b.getTitle(), channel3b.getStack()); channel3back.setCalibration(channel3b.getCalibration()); if (stitchingMethod != "Front - Back") { channel2front = stackOrientation(channel2f, stitchingMethod); channel2back = stackOrientation(channel2b, stitchingMethod); channel3front = stackOrientation(channel3f, stitchingMethod); channel3back = stackOrientation(channel3b, stitchingMethod); } if (myOVMethod == "Determined") { tempimpf = channel2front.getStack(); new_impf = tempimpf.crop(0, 0, 0, tempimpf.getWidth(), tempimpf.getHeight(), endSlice); channel2front.setStack(new_impf); channel2front.setCalibration(copyCalibration); tempimpb = channel2back.getStack(); new_impb = tempimpb.crop(0, 0, startSlice, tempimpb.getWidth(), tempimpb.getHeight(), tempimpb.getSize() - startSlice); channel2back.setStack(new_impb); channel2back.setCalibration(copyCalibration); tempimpf = channel3front.getStack(); new_impf = tempimpf.crop(0, 0, 0, tempimpf.getWidth(), tempimpf.getHeight(), endSlice); channel3front.setStack(new_impf); channel3front.setCalibration(copyCalibration); tempimpb = channel3back.getStack(); new_impb = tempimpb.crop(0, 0, startSlice, tempimpb.getWidth(), tempimpb.getHeight(), tempimpb.getSize() - startSlice); channel3back.setStack(new_impb); channel3back.setCalibration(copyCalibration); } } ImagePlus frontTemplate = new ImagePlus(fTemplate.getTitle(), fTemplate.getStack()); frontTemplate.setCalibration(fTemplate.getCalibration()); ImagePlus backTemplate = new ImagePlus(bTemplate.getTitle(), bTemplate.getStack()); backTemplate.setCalibration(bTemplate.getCalibration()); if (stitchingMethod != "Front - Back") { impf = stackOrientation(impf, stitchingMethod); impb = stackOrientation(impb, stitchingMethod); frontTemplate = stackOrientation(fTemplate, stitchingMethod); backTemplate = stackOrientation(bTemplate, stitchingMethod); } if (myOVMethod == "Determined") { ov = (int) (endSlice - startSlice); tempimpf = impf.getStack(); new_impf = tempimpf.crop(0, 0, 0, tempimpf.getWidth(), tempimpf.getHeight(), endSlice); impf.setStack(new_impf); impf.setCalibration(copyCalibration); tempimpb = impb.getStack(); new_impb = tempimpb.crop(0, 0, startSlice, tempimpb.getWidth(), tempimpb.getHeight(), tempimpb.getSize() - startSlice); impb.setStack(new_impb); impb.setCalibration(copyCalibration); tempimpf = frontTemplate.getStack(); new_impf = tempimpf.crop(0, 0, 0, tempimpf.getWidth(), tempimpf.getHeight(), endSlice); frontTemplate.setStack(new_impf); frontTemplate.setCalibration(copyCalibration); tempimpb = backTemplate.getStack(); new_impb = tempimpb.crop(0, 0, startSlice, tempimpb.getWidth(), tempimpb.getHeight(), tempimpb.getSize() - startSlice); backTemplate.setStack(new_impb); backTemplate.setCalibration(copyCalibration); } /*----------------*/ /* Copy */ ImagePlus impf2 = new ImagePlus(impf.getTitle(), impf.getStack()); impf2.setCalibration(impf.getCalibration()); impAlignedZYX = new ImagePlus(impb.getTitle(), impb.getStack()); impAlignedZYX.setCalibration(impb.getCalibration()); if (impb.getStackSize() < ov) { IJ.showMessage("The expected overlap exceed the back stack size."); return; } /* SIFT parameters Setting */ p.modelIndex = 3; // 1: Rigid Model, 3 : Affine Model (include the scale) p.sift.fdSize = 4; p.sift.fdBins = 8; // Number of tested directions p.rod = 0.92f; p.maxEpsilon = 25.0f; p.minInlierRatio = 0.05f; //IJ.log("psift parameters: "+p.sift.initialSigma+" "+p.sift.steps+" "+p.sift.minOctaveSize+" "+p.sift.maxOctaveSize); IJ.log("Size of front stack:" + impf.getStack().getSize() + "\n Size of back stack:" + impAlignedZYX.getStack().getSize()); /* Divers Variables */ ImagePlus impAlignedZYX2; boolean showStep = false; /** -------------------------- STEP 1: PREPROCESSING -------------------------- */ /* Log Display */ IJ.log("* IMAGE STACK REGISTRATION *"); IJ.log(" "); IJ.log("Stitching Orientation " + stitchingMethod); IJ.log("MIP Size " + MIP); IJ.log(" "); /* SIFT Object */ FloatArray2DSIFT sift = new FloatArray2DSIFT(p.sift); SIFT ijSIFT = new SIFT(sift); /** -------------------------- STEP 2: OVERLAP EXTRACTION -------------------------- */ double[] data = new double[12]; /* Start of the 3D Registration */ if (Reg3D) { for (int l = 1; l <= it; ++l) // Iterations { ImageStack stackZf = impf.getStack(); ImageStack stackZb = impAlignedZYX.getStack(); fsf.clear(); fsb.clear(); int l2 = 3 * l - 2; IJ.log("STEP " + l2 + ": BEST Z ROTATION"); /** Find Back Image Overlap according to the selected method */ if (myOVMethod == "Slice-by-Slice" || myOVMethod == "Determined") { indb = OverlapFinder(stackZf, stackZb, ijSIFT, 1, ov); } else if (myOVMethod == "Block-by-Block") { indb = recursiveOverlapFinder(stackZf, stackZb, ijSIFT, stackZb.getSize(), split); } IJ.log(" Overlap Size " + indb); ImageStack subStackZf = stackZf.crop(0, 0, stackZf.getSize() - indb, stackZf.getWidth(), stackZf.getHeight(), indb); ImageStack subStackZb = stackZb.crop(0, 0, 0, stackZb.getWidth(), stackZb.getHeight(), indb); fsf.clear(); fsb.clear(); AbstractAffineModel2D<?> BestModelZ = CompareCrossSection(subStackZf, subStackZb, ijSIFT, MIP); IJ.log("STEP 2"); /* Back Substack Affine Registration */ if (modelFound) { Model3D(null, null, BestModelZ); impAlignedZYX = Rotation3D(impb, BestModel3D); } else { IJ.log("No model found for the data"); } /* Show intermediate transformation matrix */ BestModel3D.toArray(data); IJ.log("1| " + data[0] + "\t| " + data[3] + "\t| " + data[6] + "\t| " + data[9]); IJ.log("2| " + data[1] + "\t| " + data[4] + "\t| " + data[7] + "\t| " + data[10]); IJ.log("3| " + data[2] + "\t| " + data[5] + "\t| " + data[8] + "\t| " + data[11]); /** -------------------------- STEP 3: BEST X ROTATION -------------------------- */ l2 = 3 * l - 1; IJ.log(" "); IJ.log("STEP " + l2 + ": BEST X ROTATION"); /* Parameter setting */ fsf.clear(); fsb.clear(); /* Rotation of the front substack around y */ impf = StackRotation(impf, 0, 90, 0); ImageStack stackYf = impf.getStack(); /* Rotation of the back substack around y */ impAlignedZYX = StackRotation(impAlignedZYX, 0, 90, 0); ImageStack stackYb = impAlignedZYX.getStack(); /* Cropping the overlap */ ImageStack subStackYf = stackYf.crop(stackYf.getWidth() - indb, 0, 0, indb, stackYf.getHeight(), stackYf.getSize()); ImageStack subStackYb = stackYb.crop(0, 0, 0, indb, stackYb.getHeight(), stackYb.getSize()); /* Comparison */ AbstractAffineModel2D<?> BestModelX = CompareCrossSection(subStackYf, subStackYb, ijSIFT, MIP); /* Alignment */ if (modelFound) { Model3D(BestModelX, null, null); impAlignedZYX = Rotation3D(impb, BestModel3D); if (showStep) { impAlignedZYX2 = StackRotation(impAlignedZYX, 0, 90, 0); impAlignedZYX2.setTitle(String.valueOf(l2)); impAlignedZYX2.show(); } } else { IJ.log("No model found for the data"); impAlignedZYX = StackRotation(impAlignedZYX, 0, -90, 0); } /* Transformation matrix */ BestModel3D.toArray(data); IJ.log("1| " + data[0] + "\t| " + data[3] + "\t| " + data[6] + "\t| " + data[9]); IJ.log("2| " + data[1] + "\t| " + data[4] + "\t| " + data[7] + "\t| " + data[10]); IJ.log("3| " + data[2] + "\t| " + data[5] + "\t| " + data[8] + "\t| " + data[11]); /* Show intermediate steps */ //impAlignedZYX.setTitle(String.valueOf( l2 )); //impAlignedZYX.show(); /** -------------------------- STEP 4: BEST Y ROTATION -------------------------- */ IJ.log(" "); IJ.log("STEP " + 3 * l + ": BEST Y ROTATION"); /* Parameter setting */ fsf.clear(); fsb.clear(); /* Rotation of the front substack around x */ impf = StackRotation(impf2, -90, 0, 0); ImageStack stackXf = impf.getStack(); /* Rotation of the back substack around x */ impAlignedZYX = StackRotation(impAlignedZYX, -90, 0, 0); ImageStack stackXb = impAlignedZYX.getStack(); /* Cropping the overlap */ ImageStack subStackXf = stackXf.crop(0, stackXf.getHeight() - indb, 0, stackXf.getWidth(), indb, stackXf.getSize()); ImageStack subStackXb = stackXb.crop(0, 0, 0, stackXb.getWidth(), indb, stackXb.getSize()); /* Comparison */ AbstractAffineModel2D<?> BestModelY = CompareCrossSection(subStackXf, subStackXb, ijSIFT, MIP); if (modelFound) { Model3D(null, BestModelY, null); impAlignedZYX = Rotation3D(impb, BestModel3D); if (showStep) { impAlignedZYX2 = StackRotation(impAlignedZYX, -90, 0, 0); impAlignedZYX2.setTitle(String.valueOf(3 * l)); impAlignedZYX2.show(); } } else { IJ.log("No model found for the data"); impAlignedZYX = StackRotation(impAlignedZYX, 90, 0, 0); } /* Transformation matrix */ BestModel3D.toArray(data); IJ.log("1| " + data[0] + "\t| " + data[3] + "\t| " + data[6] + "\t| " + data[9]); IJ.log("2| " + data[1] + "\t| " + data[4] + "\t| " + data[7] + "\t| " + data[10]); IJ.log("3| " + data[2] + "\t| " + data[5] + "\t| " + data[8] + "\t| " + data[11]); IJ.log(" "); /** -------------------------- REPLACE L'IMAGE COMME AU DEBUT -------------------------- */ impf = impf2.duplicate(); } // End of the iteration loop } // End of the 3D Registration /** -------------------------- LAST OVERLAP EXTRACTION -------------------------- */ int l2 = 3 * it + 1; IJ.log("STEP " + l2 + ": BEST Z ROTATION"); fsf.clear(); fsb.clear(); ImageStack stackZf = impf.getStack(); //impAlignedZYX=Rotation3D(impb, BestModel3D); ImageStack stackZb = impAlignedZYX.getStack(); /* Find Back Image Overlap */ fsf.clear(); fsb.clear(); int start = 1; int end = stackZb.getSize(); if (Reg3D == true) { start = indb - 250; if (start < 10) { start = 1; } end = indb + 250; if (end > stackZb.getSize()) { end = Math.max(stackZb.getSize(), stackZf.getSize()); } } //IJ.log("start: "+ start+" end:"+ end); /* Find the overlap using the slice-by-slice method */ indb = OverlapFinder(stackZf, stackZb, ijSIFT, start, end); //IJ.log( " size front:"+stackZf.getSize() +" size back:"+ stackZb.getSize()+" overlap:"+ indb); /* Cropping the overlap */ ImageStack subStackZf = stackZf.crop(0, 0, stackZf.getSize() - indb, stackZf.getWidth(), stackZf.getHeight(), indb); ImageStack subStackZb = stackZb.crop(0, 0, 0, stackZb.getWidth(), stackZb.getHeight(), indb); /* Comparison */ fsf.clear(); fsb.clear(); if (indb < MIP) { MIP = indb; } AbstractAffineModel2D<?> BestModelZ = CompareCrossSection(subStackZf, subStackZb, ijSIFT, MIP); /* Alignment */ if (modelFound) { Model3D(null, null, BestModelZ); impAlignedZYX = Rotation3D(impb, BestModel3D); impAlignedZYX.setCalibration(copyCalibration); if (template_bool == true) { impAlignedZYX = Rotation3D(backTemplate, BestModel3D); impAlignedZYX.setCalibration(copyCalibration); impf = new ImagePlus(frontTemplate.getTitle(), frontTemplate.getStack()); } } else { IJ.log("No model found for the data"); } /* Transformation matrix */ BestModel3D.toArray(data); IJ.log("1| " + data[0] + "\t|" + data[3] + "\t|" + data[6] + "\t|" + data[9]); IJ.log("2| " + data[1] + "\t|" + data[4] + "\t|" + data[7] + "\t|" + data[10]); IJ.log("3| " + data[2] + "\t|" + data[5] + "\t|" + data[8] + "\t|" + data[11]); /* Show final back stack */ //impAlignedZYX.setTitle("Final Back Stack"); //impAlignedZYX.show(); /** -------------------------- IMAGE FUSION -------------------------- */ IJ.log(" "); IJ.log("STEP: IMAGE FUSION"); String title = "Fused Image"; if (c2 || c3) { title = title + " Channel1"; } /* Channel 1 fusion */ ImagePlus FinalImg = fuseImages(impf, impAlignedZYX, indb, fusionMethod, title, copyCalibration); if (stitchingMethod != "Front - Back") { FinalImg = reverseStackOrientation(FinalImg, stitchingMethod); } FinalImg.show(); FinalImg.draw(); /* Channel 2 fusion */ if (c2) { ImagePlus impAligned2 = Rotation3D(channel2back, BestModel3D); ImagePlus FinalImg2 = fuseImages(channel2front, impAligned2, indb, fusionMethod, "Fused Image Channel2", copyCalibration); if (stitchingMethod != "Front - Back") { FinalImg2 = reverseStackOrientation(FinalImg2, stitchingMethod); } FinalImg2.show(); FinalImg2.draw(); } /* Channel 3 fusion */ if (c3) { ImagePlus impAligned3 = Rotation3D(channel3back, BestModel3D); ImagePlus FinalImg3 = fuseImages(channel3front, impAligned3, indb, fusionMethod, "Fused Image Channel3", copyCalibration); if (stitchingMethod != "Front - Back") { FinalImg3 = reverseStackOrientation(FinalImg3, stitchingMethod); } FinalImg3.show(); FinalImg3.draw(); } fsf.clear(); fsb.clear(); IJ.log(" took " + (System.currentTimeMillis() - start_time) + "ms"); IJ.log("* Done *"); IJ.log(" "); Toolkit.getDefaultToolkit().beep(); } /** Functions */ /** Find the overlap between front and back stacks @param front stack @param back stack @param sift object with parameters set as detailled in run method @param overlap size initialization @return real overlap size */ public int OverlapFinder(ImageStack stackf, ImageStack stackb, SIFT ijSIFT, int start, int end) { ImageProcessor ipb; ImageProcessor ipf = stackf.getProcessor(stackf.getSize()); /** Extration des features de la dernire slice du substack avant */ //long start_time = System.currentTimeMillis(); ijSIFT.extractFeatures(ipf, fsf);// fsf stores the features float[] bestModelInliers = new float[end - start + 1]; //number of matches float[] sliceNumber = new float[end - start + 1]; int ind = 1; //maximum indice float max = 0; //maximum value /** Comparison with slices from the back stack */ for (int i = start; i <= end; ++i) { ipb = stackb.getProcessor(i); Vector<PointMatch> inliers = searchBestInliers(ipf, ipb, ijSIFT, false); bestModelInliers[i - start] = (float) inliers.size(); sliceNumber[i - start] = (float) i; /** Best model index computation */ if (inliers.size() >= max) { ind = i; max = inliers.size(); } } Plot myplot = new Plot("Correspondence", "Slice Number", "Correspondence", sliceNumber, bestModelInliers); myplot.show(); IJ.log("(Info) Image Overlap Size : " + ind + " pixels"); searchBestInliers(ipf, stackb.getProcessor(ind), ijSIFT, true); return ind; } /** Compute the best 2D model between two substacks @param front substack @param back substack @param sift object with parameters set as detailled in run method @param MIP size @return computed model */ public AbstractAffineModel2D<?> CompareCrossSection(ImageStack subStack1, ImageStack subStack2, SIFT ijSIFT, int MIP) { /* Sauvegarde des Substacks */ ImageStack frontCOPY = subStack1.duplicate(); ImageStack backCOPY = subStack2.duplicate(); /* Construction des MIP */ if (MIP > 1) { int step = MIP; if (MIP > subStack1.getSize()) { step = subStack1.getSize(); } subStack1 = createMIP(subStack1, step); subStack2 = createMIP(subStack2, step); } /* Parameters Initialisation */ ImageProcessor ip2; ImageProcessor ip1; fsf.clear(); fsb.clear(); int taille = subStack2.getSize(); IJ.log("taille " + taille + " and size of stacks, 1= " + subStack1.getSize() + " 2= " + subStack2.getSize()); /* Comparisons of cross-section simultaneously */ for (int i = 1; i <= taille; ++i) { ip1 = subStack1.getProcessor(i); ijSIFT.extractFeatures(ip1, fsf); ip2 = subStack2.getProcessor(i); ijSIFT.extractFeatures(ip2, fsb); } /* Candidate computation */ System.out.print("identifying correspondences using brute force ..."); Vector<PointMatch> candidates = FloatArray2DSIFT.createMatches(fsb, fsf, 1.5f, null, Float.MAX_VALUE, p.rod); Vector<PointMatch> inliers = new Vector<PointMatch>(); /* Model Computation */ AbstractAffineModel2D<?> BestModel; switch (p.modelIndex) { case 0: BestModel = new TranslationModel2D(); break; case 1: BestModel = new RigidModel2D(); break; case 2: BestModel = new SimilarityModel2D(); break; case 3: BestModel = new AffineModel2D(); break; default: BestModel = new RigidModel2D(); } try { modelFound = BestModel.filterRansac(candidates, inliers, 1000, p.maxEpsilon, p.minInlierRatio); } catch (Exception e) { modelFound = false; System.err.println(e.getMessage()); } if (modelFound) { double[][] data = new double[2][3]; BestModel.toMatrix(data); IJ.log("Rotation : " + Math.acos(data[0][0]) * (180 / Math.PI) + ""); IJ.log("Horizontal Translation : " + data[0][2] + "pixels"); IJ.log("Vertical Translation : " + data[1][1] + "pixels"); ip1 = createMIP(frontCOPY, frontCOPY.getSize()).getProcessor(1); ip2 = createMIP(backCOPY, backCOPY.getSize()).getProcessor(1); displayFeatures(ip1, ip2, candidates, inliers, modelFound); IJ.log("(Info) Number of Matching Features : " + inliers.size()); } return BestModel; } /** Compute the partial MIP stack from a complete stack for a given step (number of slices) @param stack @param step (number of slices in each partial MIP) @return the real overlap size */ public ImageStack createMIP(ImageStack myStack, int step) { ImagePlus myIM = new ImagePlus("myStack", myStack); ImageProcessor myIP = myIM.getProcessor(); ImageStack myMIP = new ImageStack(myIM.getWidth(), myIM.getHeight()); int Zmax = myIM.getStackSize() / step; int extra = myIM.getStackSize() - Zmax * step; for (int j = 0; j < Zmax; j++) { ImageStack myTemp = new ImageStack(myIM.getWidth(), myIM.getHeight()); ImagePlus imTemp; for (int i = 1; i <= step; i++) { myIM.setSlice(i + step * j); myTemp.addSlice(myIP); } if ((j == Zmax - 1) && (extra != 0)) { for (int i = 1; i <= extra; i++) { myIM.setSlice(i + step * Zmax); myTemp.addSlice(myIP); } } imTemp = new ImagePlus("myTemp", myTemp); IJ.run(imTemp, "Z Project...", "projection=[Max Intensity]"); imTemp = WindowManager.getCurrentImage(); imTemp.setSlice(1); myMIP.addSlice(imTemp.getProcessor()); imTemp.hide(); } return myMIP; } /** Compute the matching key points between two image processors @param image processor 1 (front stack) @param image processor 1 (back stack) @param sift object @param boolen (information shown or not) @return the matching key points */ public Vector<PointMatch> searchBestInliers(ImageProcessor ip1, ImageProcessor ip2, SIFT ijSIFT, boolean showInfoBoolean) { fsb.clear(); ijSIFT.extractFeatures(ip2, fsb); System.out.print("identifying correspondences using brute force ..."); Vector<PointMatch> candidates = FloatArray2DSIFT.createMatches(fsb, fsf, 1.5f, null, Float.MAX_VALUE, p.rod); Vector<PointMatch> inliers = new Vector<PointMatch>(); AbstractAffineModel2D<?> BestModel; switch (p.modelIndex) { case 0: BestModel = new TranslationModel2D(); break; case 1: BestModel = new RigidModel2D(); break; case 2: BestModel = new SimilarityModel2D(); break; case 3: BestModel = new AffineModel2D(); break; default: BestModel = new RigidModel2D(); } try { modelFound = BestModel.filterRansac(candidates, inliers, 1000, p.maxEpsilon, p.minInlierRatio); } catch (Exception e) { modelFound = false; System.err.println(e.getMessage()); } if (modelFound) { if (showInfoBoolean) { double[][] data = new double[2][3]; BestModel.toMatrix(data); displayFeatures(ip1, ip2, candidates, inliers, modelFound); } } return inliers; } /** Compute the best transformation model given a SIFT object @param image processor @param sift object @return best transformation model */ public AbstractAffineModel2D<?> searchBestModel(ImageProcessor ipb, SIFT ijSIFT) { fsb.clear(); ijSIFT.extractFeatures(ipb, fsb); System.out.print("identifying correspondences using brute force ..."); Vector<PointMatch> candidates = FloatArray2DSIFT.createMatches(fsb, fsf, 1.5f, null, Float.MAX_VALUE, p.rod); Vector<PointMatch> inliers = new Vector<PointMatch>(); AbstractAffineModel2D<?> BestModel; switch (p.modelIndex) { case 0: BestModel = new TranslationModel2D(); break; case 1: BestModel = new RigidModel2D(); break; case 2: BestModel = new SimilarityModel2D(); break; case 3: BestModel = new AffineModel2D(); break; default: BestModel = new RigidModel2D(); } try { modelFound = BestModel.filterRansac(candidates, inliers, 1000, p.maxEpsilon, p.minInlierRatio); } catch (Exception e) { modelFound = false; System.err.println(e.getMessage()); } return BestModel; } public void keyPressed(KeyEvent e) { if ((e.getKeyCode() == KeyEvent.VK_F1) && (e.getSource() instanceof TextField)) { } } public void keyReleased(KeyEvent e) { } public void keyTyped(KeyEvent e) { } /** Fuses the front stack with the registered back stack @param image 1 @param image 2 @param overlap size @param fusion method @return Fused image */ private ImagePlus fuseImages(ImagePlus imp1, ImagePlus imp2, int ov, String fusionMethod, String name, Calibration finalcal) { /* Parameter setting */ final int dim = 3; ImageInformation i1 = new ImageInformation(dim, 1, null); // i1.closeAtEnd = false; i1.imp = imp1; i1.size[0] = imp1.getWidth(); i1.size[1] = imp1.getHeight(); i1.size[2] = imp1.getStack().getSize(); i1.position = new float[dim]; i1.position[0] = i1.position[1] = i1.position[2] = 0; i1.imageName = imp1.getTitle(); i1.invalid = false; i1.imageType = imp1.getType(); ImageInformation i2 = new ImageInformation(dim, 2, null); i2.closeAtEnd = false; i2.imp = imp2; i2.size[0] = imp2.getWidth(); i2.size[1] = imp2.getHeight(); i2.size[2] = imp2.getStack().getSize(); i2.position = new float[dim]; i2.position[0] = 0; i2.position[1] = 0; i2.position[2] = imp1.getStack().getSize() - ov; i2.imageName = imp2.getTitle(); i2.invalid = false; i2.imageType = imp2.getType(); ArrayList<ImageInformation> imageInformationList = new ArrayList<ImageInformation>(); imageInformationList.add(i1); imageInformationList.add(i2); /* Fusion */ float[] maximum = getAndApplyMinMax(imageInformationList, dim); ImagePlus image = Stitch_Image_Collection.fuseImages(imageInformationList, maximum, name, fusionMethod, "rgb", dim, alpha, true); image.setCalibration(finalcal); return image; } float[] getAndApplyMinMax(final ArrayList<ImageInformation> imageInformationList, final int dim) { float min[] = new float[dim]; float max[] = new float[dim]; for (int i = 0; i < min.length; i++) { min[i] = Float.MAX_VALUE; max[i] = Float.MIN_VALUE; } for (ImageInformation iI : imageInformationList) for (int i = 0; i < min.length; i++) { if (iI.position[i] < min[i]) min[i] = iI.position[i]; if (iI.position[i] + iI.size[i] > max[i]) max[i] = iI.position[i] + iI.size[i]; } for (ImageInformation iI : imageInformationList) for (int i = 0; i < min.length; i++) iI.position[i] -= min[i]; for (int i = 0; i < min.length; i++) { max[i] -= min[i]; min[i] = 0; } return max; } /** Rotation at 90, 180 and 270 only of the image stack @param image @param rotation around x @param rotation around y @param rotation around z @return rotated image */ public ImagePlus StackRotation(ImagePlus imp, int Rx, int Ry, int Rz) { /* Parameters */ int Xm = imp.getWidth(); int Ym = imp.getHeight(); int Zm = imp.getStackSize(); ImageStack impStack = imp.getStack(); int x = Xm; int y = Ym; int z = Zm; int cosRx = (int) Math.round(Math.cos(Math.toRadians(Rx))); int sinRx = (int) Math.round(Math.sin(Math.toRadians(Rx))); int cosRy = (int) Math.round(Math.cos(Math.toRadians(Ry))); int sinRy = (int) Math.round(Math.sin(Math.toRadians(Ry))); int cosRz = (int) Math.round(Math.cos(Math.toRadians(Rz))); int sinRz = (int) Math.round(Math.sin(Math.toRadians(Rz))); /* New Stack Size */ int newZm = ((x * sinRz + y * cosRz) * sinRx + z * cosRx) * cosRy - (x * cosRz - y * sinRz) * sinRy; int newXm = ((x * sinRz + y * cosRz) * sinRx + z * cosRx) * sinRy + (x * cosRz - y * sinRz) * cosRy; int newYm = (x * sinRz + y * cosRz) * cosRx - z * sinRx; /* Image Creation */ ImagePlus impR = NewImage.createImage("Rotation", Math.abs(newXm), Math.abs(newYm), Math.abs(newZm), imp.getBitDepth(), 1); ImageStack imgR = impR.getStack(); /* Calibration */ Calibration cal = new Calibration(); double Cx = imp.getCalibration().pixelWidth; double Cy = imp.getCalibration().pixelHeight; double Cz = imp.getCalibration().pixelDepth; double newCz = Math .abs(((Cx * sinRz + Cy * cosRz) * sinRx + Cz * cosRx) * cosRy - (Cx * cosRz - Cy * sinRz) * sinRy); double newCx = Math .abs(((Cx * sinRz + Cy * cosRz) * sinRx + Cz * cosRx) * sinRy + (Cx * cosRz - Cy * sinRz) * cosRy); double newCy = Math.abs((Cx * sinRz + Cy * cosRz) * cosRx - Cz * sinRx); cal.pixelWidth = newCx; cal.pixelHeight = newCy; cal.pixelDepth = newCz; cal.setUnit(imp.getCalibration().getUnit()); /* Rotation */ int max = 0; int newX = 0; int newY = 0; int newZ = 0; for (z = 0; z < Zm; z++) { for (x = 0; x < Xm; x++) { for (y = 0; y < Ym; y++) { newZ = ((x * sinRz + y * cosRz) * sinRx + z * cosRx) * cosRy - (x * cosRz - y * sinRz) * sinRy; if (newZm < 0) { newZ = newZ - newZm - 1; } newX = ((x * sinRz + y * cosRz) * sinRx + z * cosRx) * sinRy + (x * cosRz - y * sinRz) * cosRy; if (newXm < 0) { newX = newX - newXm - 1; } newY = (x * sinRz + y * cosRz) * cosRx - z * sinRx; if (newYm < 0) { newY = newY - newYm - 1; } imgR.setVoxel(newX, newY, newZ, impStack.getVoxel(x, y, z)); if (impStack.getVoxel(x, y, z) > max) { max = (int) impStack.getVoxel(x, y, z); } } } } ImagePlus imp2 = new ImagePlus("imp", imgR); imp2.setCalibration(cal); /** 12 bit Images */ if (255 < max && max <= 4095) { ImagePlus.setDefault16bitRange(12); } return imp2; } /** Rotation in 3D of an image stack using a transformation model @param image to be rotated @param transformation model @return rotated image */ public ImagePlus Rotation3D(ImagePlus imp, AffineModel3D model) { Calibration c = imp.getCalibration(); float zFactor = (float) (c.pixelDepth / c.pixelWidth); AffineModel3D unScale = new AffineModel3D(); unScale.set(1.0f, 0.0f, 0.0f, 0.0f, 0.0f, 1.0f, 0.0f, 0.0f, 0.0f, 0.0f, zFactor, 0.0f); /* center shift + un-shift */ TranslationModel3D centerShift = new TranslationModel3D(); centerShift.set(-imp.getWidth() / 2, -imp.getHeight() / 2, -imp.getStack().getSize() / 2 * zFactor); TranslationModel3D centerUnShift = new TranslationModel3D(); centerUnShift.set(imp.getWidth() / 2, imp.getHeight() / 2, imp.getStack().getSize() / 2 * zFactor); AffineModel3D transform = new AffineModel3D(); transform.preConcatenate(model); /* bounding volume */ int w = imp.getWidth(); int h = imp.getHeight(); int d = imp.getStackSize(); /* render target stack */ mpicbg.ij.stack.InverseTransformMapping<AffineModel3D> mapping = new mpicbg.ij.stack.InverseTransformMapping<AffineModel3D>( transform); ImageProcessor ip = imp.getStack().getProcessor(1).createProcessor(imp.getWidth(), imp.getHeight()); ImageStack targetStack = new ImageStack(w, h); for (int s = 0; s < d; ++s) { ip = ip.createProcessor(w, h); mapping.setSlice(s); try { mapping.mapInterpolated(imp.getStack(), ip); } catch (Exception e) { e.printStackTrace(); } targetStack.addSlice("", ip); } /* set proper calibration (it's isotropic at the former x,y-scale now) */ ImagePlus impTarget = new ImagePlus("target", targetStack); return impTarget; } /** Concatenation of 2D transformation models @param transformation model perpendicular to x @param transformation model perpendicular to y @param transformation model perpendicular to z @return 3D transformation model */ public AffineModel3D Model3D(AbstractAffineModel2D<?> modelX, AbstractAffineModel2D<?> modelY, AbstractAffineModel2D<?> modelZ) { AffineModel3D modelTemp = new AffineModel3D(); double[] data = new double[6]; /* Transformation model (Rotation around Z) */ if (modelZ != null) { modelZ.toArray(data); modelTemp.set(data[0], data[2], 0.0f, data[4], data[1], data[3], 0.0f, data[5], 0.0f, 0.0f, 1.0f, 0.0f); BestModel3D.concatenate(modelTemp); } /* Transformation model (Rotation around Y) */ if (modelX != null) { modelX.toArray(data); modelTemp.set(1.0f, 0.0f, 0.0f, 0.0f, 0.0f, data[3], data[1], data[5], //0.0f, data[2], data[0], data[4]); 0.0f, data[2], data[0], 0.0f); BestModel3D.concatenate(modelTemp); } /* Transformation model (Rotation around Y) */ if (modelY != null) { /* center shift + un-shift*/ // ATTENTION Change the origin of the rotation TranslationModel3D Shift = new TranslationModel3D(); TranslationModel3D UnShift = new TranslationModel3D(); modelY.toArray(data); modelTemp.set(data[0], 0.0f, data[2], data[4], 0.0f, 1.0f, 0.0f, 0.0f, //data[1], 0.0f, data[3], data[5]); data[1], 0.0f, data[3], 0.0f); BestModel3D.concatenate(UnShift); BestModel3D.concatenate(modelTemp); BestModel3D.concatenate(Shift); } return BestModel3D; } /** Display SIFT position and matching between two images @param image 1 @param image 2 @param candidate vectors (red) @param matching vectors (green) @param boolean @return features (display) */ public void displayFeatures(ImageProcessor ip1, ImageProcessor ip2, Vector<PointMatch> candidates, Vector<PointMatch> inliers, boolean modelfound) { /* Features vizualisation */ ImageProcessor ip1wf = null; ImageProcessor ip2wf = null; float vis_scale = 1; ip1wf = ip1.convertToRGB().duplicate(); ip2wf = ip2.convertToRGB().duplicate(); ip1wf.setColor(Color.red); ip2wf.setColor(Color.red); ip1wf.setLineWidth(2); ip2wf.setLineWidth(2); for (PointMatch m : candidates) { double[] m_p1 = m.getP1().getL(); double[] m_p2 = m.getP2().getL(); //float[] m_p1 = m.getP1().getL(); //float[] m_p2 = m.getP2().getL(); ip1wf.drawDot((int) Math.round(vis_scale * m_p2[0]), (int) Math.round(vis_scale * m_p2[1])); ip2wf.drawDot((int) Math.round(vis_scale * m_p1[0]), (int) Math.round(vis_scale * m_p1[1])); } if (modelfound) { ip1wf.setColor(Color.green); ip2wf.setColor(Color.green); ip1wf.setLineWidth(2); ip2wf.setLineWidth(2); for (PointMatch m : inliers) { //float[] m_p1 = m.getP1().getL(); //float[] m_p2 = m.getP2().getL(); double[] m_p1 = m.getP1().getL(); double[] m_p2 = m.getP2().getL(); ip1wf.drawDot((int) Math.round(vis_scale * m_p2[0]), (int) Math.round(vis_scale * m_p2[1])); ip2wf.drawDot((int) Math.round(vis_scale * m_p1[0]), (int) Math.round(vis_scale * m_p1[1])); } ImageStack stackInfo = new ImageStack(Math.round(ip1.getWidth()), Math.round(ip1.getHeight())); ImageProcessor tmp; tmp = ip1wf.createProcessor(stackInfo.getWidth(), stackInfo.getHeight()); tmp.insert(ip1wf, 0, 0); stackInfo.addSlice(null, tmp); tmp = ip2wf.createProcessor(stackInfo.getWidth(), stackInfo.getHeight()); tmp.insert(ip2wf, 0, 0); stackInfo.addSlice(null, tmp); ImagePlus impInfo = new ImagePlus("Alignment info", stackInfo); impInfo.show(); } } /** Orientate the stack 1 and 2 according to the user selection (right-left/top-bottom...) @param image @param orientation selection on the interface @return correctly oriented image */ public ImagePlus stackOrientation(ImagePlus imp, String orientation) { if (orientation == "Back - Front") { imp.revert(); } else if (orientation == "Left - Right") { imp = StackRotation(imp, 0, -90, 0); } else if (orientation == "Right - Left") { imp = StackRotation(imp, 0, 90, 0); } else if (orientation == "Top - Bottom") { imp = StackRotation(imp, 90, 0, 0); } else if (orientation == "Bottom - Top") { imp = StackRotation(imp, -90, 0, 0); } return imp; } /** Replace the stack in its original orientation @param image @param orientation selection on the interface @return correctly oriented image */ public ImagePlus reverseStackOrientation(ImagePlus imp, String orientation) { if (orientation == "Back - Front") { imp.revert(); } else if (orientation == "Left - Right") { imp = StackRotation(imp, 0, 90, 0); } else if (orientation == "Right - Left") { imp = StackRotation(imp, 0, -90, 0); } else if (orientation == "Top - Bottom") { imp = StackRotation(imp, -90, 0, 0); } else if (orientation == "Bottom - Top") { imp = StackRotation(imp, 90, 0, 0); } return imp; } /** Find the overlap position using the block-by-block technique @param image 1 @param image 2 @param SIFT @param estimated overlap position @param image spliting parameter @return overlap position */ public int recursiveOverlapFinder(ImageStack stackf, ImageStack stackb, SIFT ijSIFT, int Cov, int split) { int sb = stackb.getSize(); int sf = stackf.getSize(); ImageStack subStackf = new ImageStack(); ImageStack subStackb = new ImageStack(); int myMIP = (int) sb / split; if (myMIP > 1) { subStackb = createMIP(stackb, myMIP); if (sf >= (int) sb / split) { subStackf = makeSubstack(stackf, (int) sf - sb / split + 1, sf); subStackf = createMIP(subStackf, subStackf.getSize()); } else { subStackf = createMIP(stackf, stackf.getSize()); } } else if (myMIP <= 1) { subStackb = stackb; subStackf = stackf; } /* Parameters Initialisation */ ImageProcessor ipf; ImageProcessor ipb; fsf.clear(); fsb.clear(); /* Comparisons */ int taille = subStackb.getSize(); ipf = subStackf.getProcessor(subStackf.getSize()); // Last slice of the front stack ijSIFT.extractFeatures(ipf, fsf); int ind = 1; int max = 0; for (int i = 1; i <= taille; ++i) // Find the best match { ipb = subStackb.getProcessor(i); ijSIFT.extractFeatures(ipb, fsb); Vector<PointMatch> inliers = searchBestInliers(ipf, ipb, ijSIFT, false); if (inliers.size() >= max) { ind = i; max = inliers.size(); } } if (myMIP > 1) { Cov = (int) (Cov - sb + ind * sb / split); stackb = makeSubstack(stackb, (int) ((ind - 1) * (myMIP) + 1), (int) (ind * myMIP)); stackf = makeSubstack(stackf, (int) sf - myMIP + 1, sf); } else { Cov = (int) (Cov - sb + ind); stackb = makeSubstack(stackb, (int) ((ind - 1) * (sb / taille) + 1), (int) (ind * sb / taille)); stackf = makeSubstack(stackf, (int) sf - sb / taille, sf); } if (sb <= split || split == 1) { IJ.log("Return OV: " + Cov); return Cov; } else { Cov = recursiveOverlapFinder(stackf, stackb, ijSIFT, Cov, split); return Cov; } } /** Make a substack between two given positions in a 3D image stack @param stack @param position 1 @param position 2 @return substack */ public ImageStack makeSubstack(ImageStack stack, int z1, int z2) { ImageStack substack = new ImageStack(stack.getWidth(), stack.getHeight()); for (int i = z1; i <= z2; ++i) { substack.addSlice(stack.getProcessor(i)); } return substack; } }