use of com.forstudy.efjava.ch03.item10.Point in project imagej-ops by imagej.
the class DefaultDetectRidges method getNextPoint.
/**
* Recursively determines the next line point and adds it to the running list
* of line points.
*
* @param gradientRA - the {@link RandomAccess} of the gradient image.
* @param pRA - the {@link RandomAccess} of the eigenvector image.
* @param nRA - the {@link RandomAccess} of the subpixel line location image.
* @param points - the {@link ArrayList} containing the line points.
* @param octant - integer denoting the octant of the last gradient vector,
* oriented with 1 being 0 degrees and increasing in the
* counterclockwise direction.
* @param lastnx - the x component of the gradient vector of the last line
* point.
* @param lastny - the y component of the gradient vector of the last line
* point.
* @param lastpx - the x component of the subpixel line location of the last
* line point.
* @param lastpy - the y component of the subpixel line location of the last
* line point.
*/
private void getNextPoint(RandomAccess<DoubleType> gradientRA, RandomAccess<DoubleType> pRA, RandomAccess<DoubleType> nRA, List<RealPoint> points, int octant, double lastnx, double lastny, double lastpx, double lastpy) {
Point currentPos = new Point(gradientRA);
// variables for the best line point of the three.
Point salientPoint = new Point(gradientRA);
double salientnx = 0;
double salientny = 0;
double salientpx = 0;
double salientpy = 0;
double bestSalience = Double.MAX_VALUE;
boolean lastPointInLine = true;
// check the three possible points that could continue the line, starting at
// the octant after the given octant and rotating clockwise around the
// current pixel.
double lastAngle = RidgeDetectionUtils.getAngle(lastnx, lastny);
for (int i = 1; i < 4; i++) {
int[] modifier = RidgeDetectionUtils.getOctantCoords(octant + i);
gradientRA.move(modifier[0], 0);
gradientRA.move(modifier[1], 1);
// there.
if (gradientRA.get().get() > lowerThreshold) /*&& isMaxRA.get().get() > 0*/
{
long[] vectorArr = { gradientRA.getLongPosition(0), gradientRA.getLongPosition(1), 0 };
nRA.setPosition(vectorArr);
double nx = nRA.get().get();
nRA.fwd(2);
double ny = nRA.get().get();
pRA.setPosition(vectorArr);
double px = pRA.get().get();
pRA.fwd(2);
double py = pRA.get().get();
double currentAngle = RidgeDetectionUtils.getAngle(nx, ny);
double subpixelDiff = Math.sqrt(Math.pow(px - lastpx, 2) + Math.pow(py - lastpy, 2));
double angleDiff = Math.abs(currentAngle - lastAngle);
lastPointInLine = false;
// numbers relative to other potential line points.
if (subpixelDiff + angleDiff < bestSalience) {
// record the values of the new most salient pixel
salientPoint = new Point(gradientRA);
salientnx = nx;
salientny = ny;
salientpx = px;
salientpy = py;
bestSalience = subpixelDiff + angleDiff;
}
// set the values to zero so that they are not added to another line.
gradientRA.get().set(0);
}
// reset our randomAccess for the next check
gradientRA.setPosition(currentPos);
}
// set the current pixel to 0 in the first slice of eigenRA!
gradientRA.get().setReal(0);
// find the next line point as long as there is one to find
if (!lastPointInLine) {
// take the most salient point
gradientRA.setPosition(salientPoint);
points.add(RidgeDetectionUtils.get2DRealPoint(gradientRA.getDoublePosition(0) + salientpx, gradientRA.getDoublePosition(1) + salientpy));
// the gradient vector itself refers to the greatest change in intensity,
// and for a pixel on a line this vector will be perpendicular to the
// direction of the line. But this vector can point to either the left or
// the right of the line from the perspective of the detector, and there
// is no guarantee that the vectors at line point will point off the same
// side of the line. So if they point off different sides, set the current
// vector by 180 degrees for the purposes of this detector. We set the
// threshold for angle fixing just above 90 degrees since any lower would
// prevent ridges curving.
double potentialGradient = RidgeDetectionUtils.getAngle(salientnx, salientny);
// even though they are close enough to satisfy.
if (lastAngle < angleThreshold)
lastAngle += 360;
if (potentialGradient < angleThreshold)
potentialGradient += 360;
if (Math.abs(potentialGradient - lastAngle) > angleThreshold) {
salientnx = -salientnx;
salientny = -salientny;
}
// perform the operation again on the new end of the line being formed.
getNextPoint(gradientRA, pRA, nRA, points, RidgeDetectionUtils.getOctant(salientnx, salientny), salientnx, salientny, salientpx, salientpy);
}
}
use of com.forstudy.efjava.ch03.item10.Point in project imagej-ops by imagej.
the class WatershedSeeded method compute.
@SuppressWarnings("unchecked")
@Override
public void compute(final RandomAccessibleInterval<T> in, final ImgLabeling<Integer, IntType> out) {
// extend border to be able to do a quick check, if a voxel is inside
final LabelingType<Integer> oustide = out.firstElement().copy();
oustide.clear();
oustide.add(OUTSIDE);
final ExtendedRandomAccessibleInterval<LabelingType<Integer>, ImgLabeling<Integer, IntType>> outExt = Views.extendValue(out, oustide);
final OutOfBounds<LabelingType<Integer>> raOut = outExt.randomAccess();
// if no mask provided, set the mask to the whole image
if (mask == null) {
mask = (RandomAccessibleInterval<B>) ops().create().img(in, new BitType());
for (B b : Views.flatIterable(mask)) {
b.set(true);
}
}
// initialize output labels
final Cursor<B> maskCursor = Views.flatIterable(mask).cursor();
while (maskCursor.hasNext()) {
maskCursor.fwd();
if (maskCursor.get().get()) {
raOut.setPosition(maskCursor);
raOut.get().clear();
raOut.get().add(INIT);
}
}
// RandomAccess for Mask, Seeds and Neighborhoods
final RandomAccess<B> raMask = mask.randomAccess();
final RandomAccess<LabelingType<Integer>> raSeeds = seeds.randomAccess();
final Shape shape;
if (useEightConnectivity) {
shape = new RectangleShape(1, true);
} else {
shape = new DiamondShape(1);
}
final RandomAccessible<Neighborhood<T>> neighborhoods = shape.neighborhoodsRandomAccessible(in);
final RandomAccess<Neighborhood<T>> raNeigh = neighborhoods.randomAccess();
/*
* Carry over the seeding points to the new label and adds them to a
* voxel priority queue
*/
final PriorityQueue<WatershedVoxel> pq = new PriorityQueue<>();
// Only iterate seeds that are not excluded by the mask
final IterableRegion<B> maskRegions = Regions.iterable(mask);
final IterableInterval<LabelingType<Integer>> seedsMasked = Regions.sample(maskRegions, seeds);
final Cursor<LabelingType<Integer>> cursorSeeds = seedsMasked.localizingCursor();
while (cursorSeeds.hasNext()) {
final Set<Integer> l = cursorSeeds.next();
if (l.isEmpty()) {
continue;
}
if (l.size() > 1) {
throw new IllegalArgumentException("Seeds must have exactly one label!");
}
final Integer label = l.iterator().next();
if (label < 0) {
throw new IllegalArgumentException("Seeds must have positive integers as labels!");
}
raNeigh.setPosition(cursorSeeds);
final Cursor<T> neighborhood = raNeigh.get().cursor();
// Add unlabeled neighbors to priority queue
while (neighborhood.hasNext()) {
neighborhood.fwd();
raSeeds.setPosition(neighborhood);
raMask.setPosition(neighborhood);
raOut.setPosition(neighborhood);
if (raOut.isOutOfBounds() || !raMask.get().get() || !raSeeds.get().isEmpty())
continue;
final Integer labelNeigh = raOut.get().iterator().next();
if (labelNeigh != INQUEUE && labelNeigh != OUTSIDE) {
raOut.setPosition(neighborhood);
pq.add(new WatershedVoxel(IntervalIndexer.positionToIndex(neighborhood, in), neighborhood.get().getRealDouble()));
raOut.get().clear();
raOut.get().add(INQUEUE);
}
}
// Overwrite label in output with the seed label
raOut.setPosition(cursorSeeds);
raOut.get().clear();
raOut.get().add(label);
}
/*
* Pop the head of the priority queue, label and push all unlabeled
* neighbored pixels.
*/
// list to store neighbor labels
final ArrayList<Integer> neighborLabels = new ArrayList<>();
// list to store neighbor voxels
final ArrayList<WatershedVoxel> neighborVoxels = new ArrayList<>();
// iterate the queue
final Point pos = new Point(in.numDimensions());
while (!pq.isEmpty()) {
IntervalIndexer.indexToPosition(pq.poll().getPos(), out, pos);
// reset list of neighbor labels
neighborLabels.clear();
// reset list of neighbor voxels
neighborVoxels.clear();
// iterate the neighborhood of the pixel
raNeigh.setPosition(pos);
final Cursor<T> neighborhood = raNeigh.get().cursor();
while (neighborhood.hasNext()) {
neighborhood.fwd();
// Unlabeled neighbors go into the queue if they are not there
// yet
raOut.setPosition(neighborhood);
raMask.setPosition(raOut);
if (!raOut.get().isEmpty()) {
final Integer label = raOut.get().iterator().next();
if (label == INIT && raMask.get().get()) {
neighborVoxels.add(new WatershedVoxel(IntervalIndexer.positionToIndex(neighborhood, out), neighborhood.get().getRealDouble()));
} else {
if (label > WSHED && (!drawWatersheds || !neighborLabels.contains(label))) {
// store labels of neighbors in a list
neighborLabels.add(label);
}
}
}
}
if (drawWatersheds) {
// if the neighbors of the extracted voxel that have already
// been labeled
// all have the same label, then the voxel is labeled with their
// label.
raOut.setPosition(pos);
raOut.get().clear();
if (neighborLabels.size() == 1) {
raOut.get().add(neighborLabels.get(0));
// list
for (final WatershedVoxel v : neighborVoxels) {
IntervalIndexer.indexToPosition(v.getPos(), out, raOut);
raOut.get().clear();
raOut.get().add(INQUEUE);
pq.add(v);
}
} else if (neighborLabels.size() > 1)
raOut.get().add(WSHED);
} else {
if (neighborLabels.size() > 0) {
raOut.setPosition(pos);
raOut.get().clear();
// take the label which most of the neighbors have
if (neighborLabels.size() > 2) {
final Map<Integer, Long> countLabels = neighborLabels.stream().collect(Collectors.groupingBy(e -> e, Collectors.counting()));
final Integer keyMax = Collections.max(countLabels.entrySet(), Comparator.comparingLong(Map.Entry::getValue)).getKey();
raOut.get().add(keyMax);
} else {
raOut.get().add(neighborLabels.get(0));
}
// list
for (final WatershedVoxel v : neighborVoxels) {
IntervalIndexer.indexToPosition(v.getPos(), out, raOut);
raOut.get().clear();
raOut.get().add(INQUEUE);
pq.add(v);
}
}
}
}
/*
* Merge already present labels before calculation of watershed
*/
if (out() != null) {
final Cursor<LabelingType<Integer>> cursor = out().cursor();
while (cursor.hasNext()) {
cursor.fwd();
raOut.setPosition(cursor);
final List<Integer> labels = new ArrayList<>();
cursor.get().iterator().forEachRemaining(labels::add);
raOut.get().addAll(labels);
}
}
}
use of com.forstudy.efjava.ch03.item10.Point in project imagej-ops by imagej.
the class ConvolveTest method testCreateAndConvolvePoints.
/**
* tests fft based convolve
*/
@Test
public void testCreateAndConvolvePoints() {
final int xSize = 128;
final int ySize = 128;
final int zSize = 128;
int[] size = new int[] { xSize, ySize, zSize };
Img<DoubleType> phantom = ops.create().img(size);
RandomAccess<DoubleType> randomAccess = phantom.randomAccess();
randomAccess.setPosition(new long[] { xSize / 2, ySize / 2, zSize / 2 });
randomAccess.get().setReal(255.0);
randomAccess.setPosition(new long[] { xSize / 4, ySize / 4, zSize / 4 });
randomAccess.get().setReal(255.0);
Point location = new Point(phantom.numDimensions());
location.setPosition(new long[] { 3 * xSize / 4, 3 * ySize / 4, 3 * zSize / 4 });
HyperSphere<DoubleType> hyperSphere = new HyperSphere<>(phantom, location, 5);
for (DoubleType value : hyperSphere) {
value.setReal(16);
}
// create psf using the gaussian kernel op (alternatively PSF could be an
// input to the script)
RandomAccessibleInterval<DoubleType> psf = ops.create().kernelGauss(new double[] { 5, 5, 5 }, new DoubleType());
RandomAccessibleInterval<DoubleType> convolved = ops.create().img(size);
// convolve psf with phantom
convolved = ops.filter().convolve(convolved, phantom, psf);
DoubleType sum = new DoubleType();
DoubleType max = new DoubleType();
DoubleType min = new DoubleType();
ops.stats().sum(sum, Views.iterable(convolved));
ops.stats().max(max, Views.iterable(convolved));
ops.stats().min(min, Views.iterable(convolved));
assertEquals(sum.getRealDouble(), 8750.000184601617, 0.0);
assertEquals(max.getRealDouble(), 3.154534101486206, 0.0);
assertEquals(min.getRealDouble(), -2.9776862220387557E-7, 0.0);
}
use of com.forstudy.efjava.ch03.item10.Point in project imagej-ops by imagej.
the class FFTTest method placeSphereInCenter.
/**
* utility that places a sphere in the center of the image
*
* @param img
*/
private void placeSphereInCenter(final Img<FloatType> img) {
final Point center = new Point(img.numDimensions());
for (int d = 0; d < img.numDimensions(); d++) center.setPosition(img.dimension(d) / 2, d);
final HyperSphere<FloatType> hyperSphere = new HyperSphere<>(img, center, 2);
for (final FloatType value : hyperSphere) {
value.setReal(1);
}
}
use of com.forstudy.efjava.ch03.item10.Point in project imagej-ops by imagej.
the class DeconvolveTest method placeSphereInCenter.
// utility to place a small sphere at the center of the image
private void placeSphereInCenter(Img<FloatType> img) {
final Point center = new Point(img.numDimensions());
for (int d = 0; d < img.numDimensions(); d++) center.setPosition(img.dimension(d) / 2, d);
HyperSphere<FloatType> hyperSphere = new HyperSphere<>(img, center, 30);
for (final FloatType value : hyperSphere) {
value.setReal(1);
}
}
Aggregations