use of net.imglib2.roi.labeling.LabelingType in project vcell by virtualcell.
the class ConstructTIRFGeometry method run.
@Override
public void run() {
// Calculate constant d in TIRF exponential decay function
// Angle of incidence in radians
theta = theta * 2 * Math.PI / 360;
// Refractive index of glass
final double n1 = 1.52;
// Refractive index of cytosol
final double n2 = 1.38;
final double d = lambda * Math.pow((Math.pow(n1, 2) * Math.pow(Math.sin(theta), 2) - Math.pow(n2, 2)), -0.5) / (4 * Math.PI);
System.out.println("d: " + d);
final double fluorPerMolecule = 250;
// Get frame of interest to define geometry
long maxX = data.dimension(0) - 1;
long maxY = data.dimension(1) - 1;
Interval interval = Intervals.createMinMax(0, 0, sliceIndex, maxX, maxY, sliceIndex);
RandomAccessibleInterval<T> croppedRAI = ops.transform().crop(data, interval, true);
// Subtract lowest pixel value
IterableInterval<T> dataII = Views.iterable(croppedRAI);
double min = ops.stats().min(dataII).getRealDouble();
Cursor<T> dataCursor = dataII.cursor();
while (dataCursor.hasNext()) {
double val = dataCursor.next().getRealDouble();
dataCursor.get().setReal(val - min);
}
// Perform Gaussian blur
RandomAccessibleInterval<T> blurredRAI = ops.filter().gauss(croppedRAI, 2);
IterableInterval<T> blurredII = Views.iterable(blurredRAI);
// Segment slice by threshold and fill holes
IterableInterval<BitType> thresholded = ops.threshold().huang(blurredII);
Img<BitType> thresholdedImg = ops.convert().bit(thresholded);
RandomAccessibleInterval<BitType> thresholdedRAI = ops.morphology().fillHoles(thresholdedImg);
// Get the largest region
RandomAccessibleInterval<LabelingType<ByteType>> labeling = ops.labeling().cca(thresholdedRAI, ConnectedComponents.StructuringElement.EIGHT_CONNECTED);
LabelRegions<ByteType> labelRegions = new LabelRegions<>(labeling);
Iterator<LabelRegion<ByteType>> iterator = labelRegions.iterator();
LabelRegion<ByteType> maxRegion = iterator.next();
while (iterator.hasNext()) {
LabelRegion<ByteType> currRegion = iterator.next();
if (currRegion.size() > maxRegion.size()) {
maxRegion = currRegion;
}
}
// Generate z index map
double iMax = ops.stats().max(dataII).getRealDouble();
Img<UnsignedShortType> dataImg = ops.convert().uint16(dataII);
Img<UnsignedShortType> zMap = ops.convert().uint16(ops.create().img(dataII));
LabelRegionCursor cursor = maxRegion.localizingCursor();
RandomAccess<UnsignedShortType> zMapRA = zMap.randomAccess();
RandomAccess<UnsignedShortType> dataRA = dataImg.randomAccess();
while (cursor.hasNext()) {
cursor.fwd();
zMapRA.setPosition(cursor);
dataRA.setPosition(cursor);
double val = dataRA.get().getRealDouble();
// Log of 0 is undefined
if (val < 1) {
val = 1;
}
int z = (int) Math.round(-d * Math.log(val / iMax) / zRes);
zMapRA.get().set(z);
}
System.out.println("6");
// Use map to construct 3D geometry
// Add 5 slices of padding on top
int maxZ = (int) ops.stats().max(zMap).getRealDouble() + 5;
long[] resultDimensions = { maxX + 1, maxY + 1, maxZ };
Img<BitType> result = new ArrayImgFactory<BitType>().create(resultDimensions, new BitType());
RandomAccess<BitType> resultRA = result.randomAccess();
System.out.println(maxZ);
cursor.reset();
while (cursor.hasNext()) {
cursor.fwd();
zMapRA.setPosition(cursor);
int zIndex = zMapRA.get().get();
int[] position = { cursor.getIntPosition(0), cursor.getIntPosition(1), zIndex };
while (position[2] < maxZ) {
resultRA.setPosition(position);
resultRA.get().set(true);
position[2]++;
}
}
output = datasetService.create(result);
CalibratedAxis[] axes = new DefaultLinearAxis[] { new DefaultLinearAxis(Axes.X), new DefaultLinearAxis(Axes.Y), new DefaultLinearAxis(Axes.Z) };
output.setAxes(axes);
System.out.println("Done constructing geometry");
}
use of net.imglib2.roi.labeling.LabelingType in project imagej-ops by imagej.
the class Watershed method compute.
@Override
public void compute(final RandomAccessibleInterval<T> in, final ImgLabeling<Integer, IntType> out) {
final RandomAccess<T> raIn = in.randomAccess();
RandomAccess<B> raMask = null;
if (mask != null) {
raMask = mask.randomAccess();
}
// stores the size of each dimension
final long[] dimensSizes = new long[in.numDimensions()];
in.dimensions(dimensSizes);
// calculates the number of points in the n-d space
long numPixels = Intervals.numElements(in);
// the pixels indices are stored in an array, which is sorted depending
// on the pixel values
final List<Long> imiList = new ArrayList<>();
if (mask != null) {
final Cursor<Void> c = Regions.iterable(mask).localizingCursor();
while (c.hasNext()) {
c.next();
imiList.add(IntervalIndexer.positionToIndex(c, in));
}
} else {
for (long i = 0; i < numPixels; i++) {
imiList.add(i);
}
}
final Long[] imi = imiList.toArray(new Long[imiList.size()]);
/*
* Sort the pixels of imi in the increasing order of their grey value
* (only the pixel indices are stored)
*/
Arrays.sort(imi, new Comparator<Long>() {
@Override
public int compare(final Long o1, final Long o2) {
IntervalIndexer.indexToPosition(o1, in, raIn);
final T value = raIn.get().copy();
IntervalIndexer.indexToPosition(o2, in, raIn);
return value.compareTo(raIn.get());
}
});
// lab and dist store the values calculated after each phase
final RandomAccessibleInterval<IntType> lab = ops().create().img(in, new IntType());
// extend border to be able to do a quick check, if a voxel is inside
final ExtendedRandomAccessibleInterval<IntType, RandomAccessibleInterval<IntType>> labExt = Views.extendBorder(lab);
final OutOfBounds<IntType> raLab = labExt.randomAccess();
final RandomAccessibleInterval<IntType> dist = ops().create().img(in, new IntType());
final RandomAccess<IntType> raDist = dist.randomAccess();
// initial values
for (final IntType pixel : Views.flatIterable(lab)) {
pixel.set(INIT);
}
int current_label = 0;
int current_dist;
final ArrayList<Long> fifo = new ArrayList<>();
// RandomAccess for Neighborhoods
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>> raNeighbor = neighborhoods.randomAccess();
/*
* Start flooding
*/
for (int j = 0; j < imi.length; j++) {
IntervalIndexer.indexToPosition(imi[j], in, raIn);
final T actualH = raIn.get().copy();
int i = j;
while (actualH.compareTo(raIn.get()) == 0) {
final long p = imi[i];
IntervalIndexer.indexToPosition(p, in, raIn);
raLab.setPosition(raIn);
raLab.get().set(MASK);
raNeighbor.setPosition(raIn);
final Cursor<T> neighborHood = raNeighbor.get().cursor();
while (neighborHood.hasNext()) {
neighborHood.fwd();
raLab.setPosition(neighborHood);
if (!raLab.isOutOfBounds()) {
final int f = raLab.get().get();
if ((f > 0) || (f == WSHED)) {
raDist.setPosition(raIn);
raDist.get().set(1);
fifo.add(p);
break;
}
}
}
i++;
if (i == imi.length) {
break;
}
IntervalIndexer.indexToPosition(imi[i], in, raIn);
}
current_dist = 1;
// add fictitious pixel
fifo.add(-1l);
while (true) {
long p = fifo.remove(0);
if (p == -1) {
if (fifo.isEmpty()) {
break;
}
fifo.add(-1l);
current_dist++;
p = fifo.remove(0);
}
IntervalIndexer.indexToPosition(p, in, raNeighbor);
final Cursor<T> neighborHood = raNeighbor.get().cursor();
raLab.setPosition(raNeighbor);
int labp = raLab.get().get();
final long[] posNeighbor = new long[neighborHood.numDimensions()];
while (neighborHood.hasNext()) {
neighborHood.fwd();
neighborHood.localize(posNeighbor);
raLab.setPosition(posNeighbor);
if (!raLab.isOutOfBounds()) {
raDist.setPosition(posNeighbor);
final int labq = raLab.get().get();
final int distq = raDist.get().get();
if ((distq < current_dist) && ((labq > 0) || (labq == WSHED))) {
// the watersheds
if (labq > 0) {
if ((labp == MASK) || (labp == WSHED)) {
labp = labq;
} else {
if (labp != labq) {
labp = WSHED;
}
}
} else {
if (labp == MASK) {
labp = WSHED;
}
}
raLab.setPosition(raNeighbor);
raLab.get().set(labp);
} else {
if ((labq == MASK) && (distq == 0)) {
raDist.setPosition(posNeighbor);
raDist.get().set(current_dist + 1);
fifo.add(IntervalIndexer.positionToIndex(posNeighbor, dimensSizes));
}
}
}
}
}
// checks if new minima have been discovered
IntervalIndexer.indexToPosition(imi[j], in, raIn);
i = j;
while (actualH.compareTo(raIn.get()) == 0) {
final long p = imi[i];
IntervalIndexer.indexToPosition(p, dist, raDist);
// the distance associated with p is reseted to 0
raDist.get().set(0);
raLab.setPosition(raDist);
if (raLab.get().get() == MASK) {
current_label++;
fifo.add(p);
raLab.get().set(current_label);
while (!fifo.isEmpty()) {
final long q = fifo.remove(0);
IntervalIndexer.indexToPosition(q, in, raNeighbor);
final Cursor<T> neighborHood = raNeighbor.get().cursor();
final long[] posNeighbor = new long[neighborHood.numDimensions()];
while (neighborHood.hasNext()) {
neighborHood.fwd();
neighborHood.localize(posNeighbor);
raLab.setPosition(posNeighbor);
if (!raLab.isOutOfBounds()) {
final long r = IntervalIndexer.positionToIndex(posNeighbor, dimensSizes);
if (raLab.get().get() == MASK) {
fifo.add(r);
raLab.get().set(current_label);
}
}
}
}
}
i++;
if (i == imi.length) {
break;
}
IntervalIndexer.indexToPosition(imi[i], in, raIn);
}
j = i - 1;
}
/*
* Draw output and remove as the case may be the watersheds
*/
final Cursor<LabelingType<Integer>> cursorOut = out.cursor();
while (cursorOut.hasNext()) {
cursorOut.fwd();
boolean maskValue = true;
if (mask != null) {
raMask.setPosition(cursorOut);
if (!raMask.get().get()) {
maskValue = false;
}
}
raLab.setPosition(cursorOut);
if (!maskValue) {
cursorOut.get().clear();
} else {
if (!drawWatersheds && raLab.get().get() == WSHED) {
raNeighbor.setPosition(cursorOut);
final Cursor<T> neighborHood = raNeighbor.get().cursor();
int newLab = WSHED;
while (neighborHood.hasNext()) {
neighborHood.fwd();
raLab.setPosition(neighborHood);
if (!raLab.isOutOfBounds()) {
newLab = raLab.get().get();
if (newLab > WSHED) {
break;
}
}
}
if (newLab == WSHED) {
cursorOut.get().clear();
} else {
cursorOut.get().add(newLab);
}
} else {
cursorOut.get().add(raLab.get().get());
}
}
}
/*
* Merge already present labels before calculation of watershed
*/
if (out() != null) {
final Cursor<LabelingType<Integer>> cursor = out().cursor();
final RandomAccess<LabelingType<Integer>> raOut = out.randomAccess();
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 net.imglib2.roi.labeling.LabelingType 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);
final Integer labelNeigh = raOut.get().iterator().next();
if (labelNeigh != INQUEUE && labelNeigh != OUTSIDE && !raOut.isOutOfBounds() && raMask.get().get() && raSeeds.get().isEmpty()) {
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 net.imglib2.roi.labeling.LabelingType in project imagej-ops by imagej.
the class CopyImgLabelingTest method createData.
@SuppressWarnings("unchecked")
@Before
public void createData() {
input = (ImgLabeling<String, IntType>) ops.run(DefaultCreateImgLabeling.class, new long[] { 10, 10 }, new IntType());
final Cursor<LabelingType<String>> inc = input.cursor();
while (inc.hasNext()) {
inc.next().add(Math.random() > 0.5 ? "A" : "B");
}
// and another loop to construct some ABs
while (inc.hasNext()) {
inc.next().add(Math.random() > 0.5 ? "A" : "B");
}
}
Aggregations