Commit cc817be4 authored by gonciarz's avatar gonciarz
Browse files

Merge branch 'particleTrackerWithClij2' into develop

parents ed021cae 965ad9ee
......@@ -16,7 +16,7 @@ import mosaic.core.imageUtils.convolution.Kernel1D;
import mosaic.core.imageUtils.convolution.Kernel2D;
import mosaic.core.imageUtils.convolution.Kernel3D;
import mosaic.core.utils.DilateImage;
import mosaic.core.utils.DilateImageClij;
/**
* FeaturePointDetector detects the "real" particles in provided frames.
......@@ -194,8 +194,13 @@ public class FeaturePointDetector {
private void pointLocationsEstimation(ImageStack ips) {
float threshold = findThreshold(ips, iPercentile, iAbsIntensityThreshold);
/* do a grayscale dilation */
final ImageStack dilated_ips = DilateImage.dilate(ips, iRadius, 4, iUseCLIJ);
final ImageStack dilated_ips = iUseCLIJ ?
DilateImageClij.dilate(ips, iRadius) :
DilateImage.dilate(ips, iRadius, 4);
// new StackWindow(new ImagePlus("dilated ", dilated_ips));
iParticles = new Vector<Particle>();
/* loop over all pixels */
final int height = ips.getHeight();
......@@ -618,7 +623,7 @@ public class FeaturePointDetector {
*
* @see #generateDilationMasks(int)
*/
public boolean setDetectionParameters(double cutoff, float percentile, int radius, float Threshold, boolean absolute, boolean aUseCLIJ) {
public boolean setDetectionParameters(double cutoff, float percentile, int radius, float Threshold, boolean absolute, boolean aUseClij) {
final boolean changed = (radius != iRadius || cutoff != iCutoff || (percentile != iPercentile));// && intThreshold != absIntensityThreshold || mode != getThresholdMode() || thsmode != getThresholdMode();
iCutoff = cutoff;
......@@ -632,9 +637,10 @@ public class FeaturePointDetector {
iThresholdMode = Mode.PERCENTILE_MODE;
}
iUseCLIJ = aUseCLIJ;
iUseCLIJ = aUseClij;
logger.info("Detection options: radius=" + iRadius + " cutoff=" + iCutoff + " percentile=" + iPercentile + " threshold=" + iAbsIntensityThreshold + " mode=" + (absolute ? "THRESHOLD" : "PERCENTILE") + " useCLIJ=" + iUseCLIJ);
System.out.println("Detection options: radius=" + iRadius + " cutoff=" + iCutoff + " percentile=" + iPercentile + " threshold=" + iAbsIntensityThreshold + " mode=" + (absolute ? "THRESHOLD" : "PERCENTILE") + " useCLIJ=" + iUseCLIJ);
// create Mask for Dilation with the user defined radius
generateDilationMasks(iRadius);
......
......@@ -23,6 +23,16 @@ import ij.gui.StackWindow;
* TODO: Obviously it must be refactored by moving into much more proper places.
*/
public class GUIhelper {
private static boolean isClijAvailable() {
try {
net.haesleinhuepf.clij2.CLIJ2.getInstance();
return true;
}
catch (final NoClassDefFoundError err) {
return false;
}
}
/**
* gd has to be shown with showDialog and handles the fields added to the dialog
* with addUserDefinedParamtersDialog(gd).
......@@ -30,13 +40,13 @@ public class GUIhelper {
* @param gd <code>GenericDialog</code> at which the UserDefinedParameter fields where added.
* @return true if user changed the parameters and false if the user didn't changed them.
*/
public static Boolean getUserDefinedParameters(GenericDialog gd, FeaturePointDetector fpd) {
public static boolean getUserDefinedParameters(GenericDialog gd, FeaturePointDetector fpd) {
final int rad = (int) gd.getNextNumber();
final double cut = gd.getNextNumber();
final float per = ((float) gd.getNextNumber()) / 100;
final float intThreshold = per * 100;
final boolean absolute = gd.getNextBoolean();
final boolean useCLIJ = gd.getNextBoolean();
final boolean useCLIJ = isClijAvailable() && gd.getNextBoolean();
return fpd.setDetectionParameters(cut, per, rad, intThreshold, absolute, useCLIJ);
}
......@@ -57,9 +67,9 @@ public class GUIhelper {
final float per = (Float.parseFloat((vec.elementAt(2)).getText())) / 100;
final float intThreshold = per * 100;
final boolean absolute = vecb.elementAt(0).getState();
final boolean useCLIJ = vecb.elementAt(1).getState();
final boolean useClij = isClijAvailable() && vecb.elementAt(1).getState();
return fpd.setDetectionParameters(cut, per, rad, intThreshold, absolute, useCLIJ);
return fpd.setDetectionParameters(cut, per, rad, intThreshold, absolute, useClij);
}
public static void addUserDefinedParametersDialog(GenericDialog gd, FeaturePointDetector fpd) {
......@@ -70,7 +80,9 @@ public class GUIhelper {
gd.addNumericField("Per/Abs", fpd.getPercentile() * 100, 3, 7, null);
gd.addCheckbox("Absolute", fpd.getThresholdMode() == FeaturePointDetector.Mode.ABS_THRESHOLD_MODE);
gd.addCheckbox("Accelerate_with_CLIJ2 (experimental)", false);
if (isClijAvailable()) {
gd.addCheckbox("Accelerate_with_CLIJ2 (experimental)", false);
}
}
......
......@@ -2,12 +2,9 @@ package mosaic.core.utils;
import ij.IJ;
import ij.ImagePlus;
import ij.ImageStack;
import ij.process.FloatProcessor;
import ij.process.ImageProcessor;
import net.haesleinhuepf.clij2.CLIJ2;
import net.haesleinhuepf.clij.clearcl.ClearCLBuffer;
import org.apache.log4j.Logger;
import java.util.Vector;
......@@ -25,52 +22,33 @@ public class DilateImage {
* @param ips ImageProcessor to do the dilation with
* @return the dilated copy of the given <code>ImageProcessor</code>
*/
public static ImageStack dilate(ImageStack ips, int radius, int number_of_threads, boolean aUseCLIJ) {
if (!aUseCLIJ) {
logger.debug("dilate start");
final FloatProcessor[] dilated_procs = new FloatProcessor[ips.getSize()];
final AtomicInteger z = new AtomicInteger(-1);
final Vector<Thread> threadsVector = new Vector<Thread>(number_of_threads);
for (int thread_counter = 0; thread_counter < number_of_threads; thread_counter++) {
threadsVector.add(new DilateThread(ips, radius, dilated_procs, z));
}
for (final Thread t : threadsVector) {
t.start();
}
logger.debug("Threads started");
for (final Thread t : threadsVector) {
try {
t.join();
} catch (final InterruptedException ie) {
IJ.showMessage("Calculation interrupted. An error occured in parallel dilation:\n" + ie.getMessage());
}
}
logger.debug("dilate threads done");
final ImageStack dilated_ips = new ImageStack(ips.getWidth(), ips.getHeight());
for (int s = 0; s < ips.getSize(); s++) {
dilated_ips.addSlice(null, dilated_procs[s]);
public static ImageStack dilate(ImageStack ips, int radius, int number_of_threads) {
logger.debug("dilate start");
final FloatProcessor[] dilated_procs = new FloatProcessor[ips.getSize()];
final AtomicInteger z = new AtomicInteger(-1);
final Vector<Thread> threadsVector = new Vector<Thread>(number_of_threads);
for (int thread_counter = 0; thread_counter < number_of_threads; thread_counter++) {
threadsVector.add(new DilateThread(ips, radius, dilated_procs, z));
}
for (final Thread t : threadsVector) {
t.start();
}
logger.debug("Threads started");
for (final Thread t : threadsVector) {
try {
t.join();
} catch (final InterruptedException ie) {
IJ.showMessage("Calculation interrupted. An error occured in parallel dilation:\n" + ie.getMessage());
}
logger.debug("dilate stop");
return dilated_ips;
}
logger.debug("dilate threads done");
final ImageStack dilated_ips = new ImageStack(ips.getWidth(), ips.getHeight());
for (int s = 0; s < ips.getSize(); s++) {
dilated_ips.addSlice(null, dilated_procs[s]);
}
logger.debug("dilate stop");
// get CLIJ2 instance with default device
CLIJ2 clij2 = CLIJ2.getInstance();
// push image and create empty one for result (with same parameters as input image)
ClearCLBuffer imgBuffer = clij2.push(new ImagePlus("", ips));
ClearCLBuffer outBuffer = clij2.create(imgBuffer);
// Run kernel
logger.debug("maximum2dSphere before");
clij2.maximum2DSphere(imgBuffer, outBuffer, radius, radius);
logger.debug("maximum2dSphere after");
// Get output and clear CLIJ2 buffers
ImagePlus imgRes = clij2.pull(outBuffer);
clij2.clear();
return imgRes.getImageStack();
return dilated_ips;
}
/**
......
package mosaic.core.utils;
import ij.ImagePlus;
import ij.ImageStack;
import net.haesleinhuepf.clij.clearcl.ClearCLBuffer;
import net.haesleinhuepf.clij2.CLIJ2;
import org.apache.log4j.Logger;
public class DilateImageClij {
private static final Logger logger = Logger.getLogger(DilateImageClij.class);
/**
* Dilates all values and returns a copy of the input image.
* A spherical structuring element of radius <code>radius</code> is used.
* Adapted as is from Ingo Oppermann implementation
*
* @param ips ImageProcessor to do the dilation with
* @return the dilated copy of the given <code>ImageProcessor</code>
*/
public static ImageStack dilate(ImageStack ips, int radius) {
// get CLIJ2 instance with default device
CLIJ2 clij2 = CLIJ2.getInstance();
// push image and create empty one for result (with same parameters as input image)
ClearCLBuffer imgBuffer = clij2.push(new ImagePlus("", ips));
ClearCLBuffer outBuffer = clij2.create(imgBuffer);
// Run kernel
if (ips.getSize() == 1) {
logger.debug("maximum2dSphere before");
clij2.maximum2DSphere(imgBuffer, outBuffer, radius, radius);
logger.debug("maximum2dSphere after");
}
else {
logger.debug("maximum3dSphere before");
clij2.maximum3DSphere(imgBuffer, outBuffer, radius, radius, radius);
logger.debug("maximum3dSphere after");
}
// Get output and clear CLIJ2 buffers
ImagePlus imgRes = clij2.pull(outBuffer);
clij2.clear();
return imgRes.getImageStack();
}
}
......@@ -666,8 +666,6 @@ public class ParticleTracker3DModular_ implements PlugInFilter, PreviewInterface
* Detects particles in the current displayed frame according to the parameters currently set
* Draws dots on the positions of the detected partciles on the frame and circles them
*
* @see #getUserDefinedPreviewParams()
* @see MyFrame#featurePointDetection()
* @see PreviewCanvas
*/
private synchronized void preview(int zDepth) {
......@@ -802,7 +800,6 @@ public class ParticleTracker3DModular_ implements PlugInFilter, PreviewInterface
* The trajectories are drawn on <code>TrajectoryCanvas</code> when it's constructed and not on the image
*
* @param duplicated_imp the image upon which the view will be updated - can be null
* @see TrajectoryStackWindow
*/
public void generateView(ImagePlus duplicated_imp, Img<ARGBType> out) {
if (iTrajectories.size() == 0) {
......@@ -900,11 +897,9 @@ public class ParticleTracker3DModular_ implements PlugInFilter, PreviewInterface
* The new Stack will be made of RGB ImageProcessors upon which the trajectory will be drawn
*
* @param trajectory_index the trajectory index in the <code>trajectories</code> Vector (starts with 0)
* @param magnification the scale factor to use for rescaling the original image
* @see IJ#run(java.lang.String, java.lang.String)
* @see ImagePlus#getRoi()
* @see StackConverter#convertToRGB()
* @see Trajectory#animate(int)
*/
public void generateTrajFocusView(int trajectory_index) {
int magnification = results_window.magnification_factor;
......@@ -959,7 +954,6 @@ public class ParticleTracker3DModular_ implements PlugInFilter, PreviewInterface
* @see IJ#run(java.lang.String, java.lang.String)
* @see ImagePlus#getRoi()
* @see StackConverter#convertToRGB()
* @see Trajectory#animate(int)
*/
public void generateAreaFocusView(int magnification) {
......@@ -1464,8 +1458,7 @@ public class ParticleTracker3DModular_ implements PlugInFilter, PreviewInterface
/**
* Rescale and filter particles
* @param aCalibration
* @param aParticles
* @param aParticles
*/
private void filterParticles(Vector<Particle> aParticles) {
int size = 0;
......
package mosaic.core.utils;
import static org.junit.Assert.*;
import ij.ImagePlus;
import ij.ImageStack;
import ij.gui.StackWindow;
import mosaic.test.framework.CommonBase;
import net.haesleinhuepf.clij.clearcl.ClearCLBuffer;
import org.apache.log4j.Logger;
import org.junit.Test;
import net.haesleinhuepf.clij2.CLIJ2;
import net.haesleinhuepf.clij.CLIJ;
import java.util.ArrayList;
public class DilateImageTest extends CommonBase {
private static final Logger logger = Logger.getLogger(DilateImageTest.class);
@Test
public void testDilate() {
int radius = 5;
ImagePlus img = loadImagePlus("/Users/gonciarz/Desktop/pt/singleFrame.tif");
logger.debug("============> CPU start");
// ImageStack dimg = DilateImage.dilate(img.getImageStack().convertToFloat(), radius, 4);
logger.debug("<============ CPU stop");
// new StackWindow(new ImagePlus("dilated ", dimg));
logger.debug("============> CLIJ2 start");
CLIJ.debug = true;
// CLIJ2 clij2 = CLIJ2.getInstance("GeForce GT 750M");
ArrayList<String> clDevs = CLIJ.getAvailableDeviceNames();
System.out.println(clDevs);
CLIJ2 clij2 = CLIJ2.getInstance(clDevs.get(1));
ClearCLBuffer imgBuffer = clij2.push(img);
ClearCLBuffer outBuffer = clij2.create(imgBuffer);
logger.debug("maximum2dSphere before");
clij2.maximum2DSphere(imgBuffer, outBuffer, radius, radius);
// clij2.maximum3DSphere(imgBuffer, outBuffer, radius, radius, radius);
logger.debug("maximum2dSphere after");
ImagePlus imgRes = clij2.pull(outBuffer);
imgBuffer.close();
outBuffer.close();
logger.debug("<============ CLIJ2 stop");
// imgRes.show();
// Maximum2DSphere.maximum2DSphere(clij2, imgBuffer, outBuffer, radius, radius);
// sleep(8000);
}
}
package mosaic.particleTracker;
import static org.junit.Assert.assertEquals;
import org.junit.FixMethodOrder;
import org.junit.Test;
import org.junit.runners.MethodSorters;
import mosaic.core.detection.Particle;
import mosaic.test.framework.CommonBase;
/**
* This class is responsible for testing {@link TrajectoryAnalysis} class.
* @author Krzysztof Gonciarz <gonciarz@mpi-cbg.de>
*/
@FixMethodOrder(MethodSorters.NAME_ASCENDING)
public class TrajectoryAnalysisTest extends CommonBase {
/**
* Tests simple trajectory with particles moving in equal length intervals
* along x axis (x,y) = {(1, 0), (2,0), ...}
*/
@Test
public void testLinearMovement() {
// Create trajectory
final int trajectoryLen = 6;
final Particle[] particles = new Particle[trajectoryLen];
for (int i = 0; i < trajectoryLen; ++i) {
particles[i] = new Particle(i + 1, 0, 0, i);
}
// Prepare Trajectory Analysis for calculations
final TrajectoryAnalysis ta = new TrajectoryAnalysis(particles);
ta.setTimeInterval(1.0);
ta.setLengthOfAPixel(1.0);
// Set some tolerance on double numbers comparisons
assertEquals("Calculation should be successful", TrajectoryAnalysis.SUCCESS, ta.calculateAll());
final double epsilon = 0.000001;
assertEquals("MSS slope", 1.0, ta.getMSSlinear(), epsilon);
assertEquals("MSS y-axis intercept", 0.0, ta.getMSSlinearY0(), epsilon);
assertEquals("D2 diffusion coefficient", 0.25, ta.getDiffusionCoefficients()[1], epsilon);
assertEquals("MSD slope", 2.0, ta.getGammasLogarithmic()[1], epsilon);
assertEquals("MSD y-axis intercept", 0.0, ta.getGammasLogarithmicY0()[1], epsilon);
assertEquals("Track lenght", 5.0, ta.getDistance(), epsilon);
assertEquals("Avg step (per frame) lenght", 1.0, ta.getAvgDistance(), epsilon);
assertEquals("Straightness", 1.0, ta.getStraightness(), epsilon);
assertEquals("Bending", 0.0, ta.getBending(), epsilon);
assertEquals("Bending (linear)", 0.0, ta.getBendingLinear(), epsilon);
assertEquals("Efficiency", 1.0, ta.getEfficiency(), epsilon);
}
/**
* Tests trajectory features calculated for zigzag movement
*/
@Test
public void testZigZagMovement() {
// Create trajectory
final int trajectoryLen = 6;
final Particle[] particles = new Particle[trajectoryLen];
particles[0] = new Particle(0, 0, 0, 1);
particles[1] = new Particle(1, 1, 0, 2);
particles[2] = new Particle(0, 2, 0, 3);
particles[3] = new Particle(1, 3, 0, 4);
particles[4] = new Particle(0, 4, 0, 5);
particles[5] = new Particle(1, 5, 0, 6);
// Prepare Trajectory Analysis for calculations
final TrajectoryAnalysis ta = new TrajectoryAnalysis(particles);
ta.setTimeInterval(1.0);
ta.setLengthOfAPixel(1.0);
// Set some tolerance on double numbers comparisons
assertEquals("Calculation should be successful", TrajectoryAnalysis.SUCCESS, ta.calculateAll());
final double epsilon = 0.000001;
assertEquals("Track lenght", 5 * Math.sqrt(2), ta.getDistance(), epsilon);
assertEquals("Avg step (per frame) lenght", Math.sqrt(2), ta.getAvgDistance(), epsilon);
assertEquals("Straightness", 6.123e-17, ta.getStraightness(), epsilon);
assertEquals("Bending", 0.0, ta.getBending(), epsilon);
assertEquals("Bending (linear)", 0.0, ta.getBendingLinear(), epsilon);
assertEquals("Efficiency", 0.52, ta.getEfficiency(), epsilon);
}
/**
* Test Trajectory with enough number of frames but with every second frame missing
* (containing only frames 0, 2, 4... ) In such case it is not possible to perform calculations
* since there are no points enough to calculate slopes (at least two needed)
*/
@Test
public void testLinearMovementWithManySkippedFrames() {
// Create trajectory
final int trajectoryLen = 6;
final Particle[] particles = new Particle[trajectoryLen];
for (int i = 0; i < trajectoryLen; ++i) {
particles[i] = new Particle(i + 1, 0, 0, i*2);
}
// Prepare Trajectory Analysis for calculations
final TrajectoryAnalysis ta = new TrajectoryAnalysis(particles);
ta.setTimeInterval(1.0);
ta.setLengthOfAPixel(1.0);
// Set some tolerance on double numbers comparisons
assertEquals("Calculation should fail", TrajectoryAnalysis.FAILURE, ta.calculateAll());
}
/**
* With trajectory shorter that 6 points it is impossible to calculate MSS/MSD because
* deltas are in range {1 ... length/3}
*/
@Test
public void testTooShortTrajectory() {
// Create too short trajectory
final int trajectoryLen = 5;
final Particle[] particles = new Particle[trajectoryLen];
for (int i = 0; i < trajectoryLen; ++i) {
particles[i] = new Particle(i + 1, 0, 0, i);
}
// Prepare Trajectory Analysis for calculations
final TrajectoryAnalysis ta = new TrajectoryAnalysis(particles);
ta.setTimeInterval(1.0);
ta.setLengthOfAPixel(1.0);
// Set some tolerance on double numbers comparisons
assertEquals("Calculation should fail", TrajectoryAnalysis.FAILURE, ta.calculateAll());
}
/**
* With trajectory shorter that 6 points it is impossible to calculate MSS/MSD because
* deltas are in range {1 ... length/3}
*/
@Test
public void testTooShortTrajectoryWithZeroElements() {
// Create simple trajectory with particles moving in equal length intervals
// along x axis
final int trajectoryLen = 0;
final Particle[] particles = new Particle[trajectoryLen];
for (int i = 0; i < trajectoryLen; ++i) {
particles[i] = new Particle(i + 1, 0, 0, i);
}
// Prepare Trajectory Analysis for calculations
final TrajectoryAnalysis ta = new TrajectoryAnalysis(particles);
// Set some tolerance on double numbers comparisons
assertEquals("Calculation should fail", TrajectoryAnalysis.FAILURE, ta.calculateAll());
}
/**
* Test extreme case with null trajectory. (test both constructors)
*/
@Test
public void testNullTrajectory() {
final Particle[] particles = null;
final TrajectoryAnalysis ta1 = new TrajectoryAnalysis(particles);
// Set some tolerance on double numbers comparisons
assertEquals("Calculation should fail", TrajectoryAnalysis.FAILURE, ta1.calculateAll());
final Trajectory trajectory = null;
final TrajectoryAnalysis ta2 = new TrajectoryAnalysis(trajectory);
// Set some tolerance on double numbers comparisons
assertEquals("Calculation should fail", TrajectoryAnalysis.FAILURE, ta2.calculateAll());
}
}
Markdown is supported
0% or .
You are about to add 0 people to the discussion. Proceed with caution.
Finish editing this message first!
Please register or to comment