Java Code Examples for net.imglib2.RandomAccessibleInterval#numDimensions()

@Override
public void compute(final RandomAccessibleInterval<T> input,
	final IterableInterval<V> output)
{
	final Cursor<V> cursor = output.localizingCursor();
	final RandomAccess<T> access = input.randomAccess();

	while (cursor.hasNext()) {
		cursor.fwd();
		for (int d = 0; d < input.numDimensions(); d++) {
			if (d != dim) {
				access.setPosition(cursor.getIntPosition(d - (d > dim ? -1 : 0)), d);
			}
		}

		method.compute(new DimensionIterable(input.dimension(dim), access),
			cursor.get());
	}
}
 

protected static long[] getMaxCoords(
	RandomAccessibleInterval<DoubleType> input, boolean useAbsoluteValue)
{
	long[] dims = new long[input.numDimensions()];
	double max = Double.MIN_VALUE;
	Cursor<DoubleType> cursor = Views.iterable(input).localizingCursor();
	while (cursor.hasNext()) {
		cursor.fwd();
		double current = useAbsoluteValue ? Math.abs(cursor.get().get()) : cursor
			.get().get();
		if (current > max) {
			max = current;
			for (int d = 0; d < input.numDimensions(); d++) {
				dims[d] = cursor.getLongPosition(d);
			}
		}
	}
	return dims;
}
 

public Align(final RandomAccessibleInterval< T > template, final ImgFactory< FloatType > factory, WarpFunction model)
{
	this.template = template;

	n = template.numDimensions();
	warpFunction = model;
	numParameters = warpFunction.numParameters();
	
	currentTransform = new AffineTransform( n );
	
	final long[] dim = new long[n + 1];
	for ( int d = 0; d < n; ++d )
		dim[d] = template.dimension( d );
	dim[n] = n;
	final Img< FloatType > gradients = factory.create( dim, new FloatType() );
	gradients( Views.extendBorder( template ), gradients );

	dim[n] = numParameters;
	descent = factory.create( dim, new FloatType() );
	computeSteepestDescents( gradients, warpFunction, descent );

	Hinv = computeInverseHessian( descent );

	error = factory.create( template, new FloatType() );
}
 

/**
 * @param source
 *            {@link RandomAccessibleInterval} which will be virtually
 *            cropped
 * @param axesOfInterest
 *            axes which define a plane, cube, hypercube, ...! All other
 *            axes will be iterated.
 * @param dropSingletonDimensions
 *            if true, dimensions of size one will be discarded in the
 *            sliced images
 */
public SlicesII(final RandomAccessibleInterval<T> source, final int[] axesOfInterest,
		final boolean dropSingletonDimensions) {
	super(initIntervals(source, axesOfInterest));

	final long[] sliceMin = new long[source.numDimensions()];
	final long[] sliceMax = new long[source.numDimensions()];

	for (int d = 0; d < source.numDimensions(); d++) {
		if (dimension(d) == 1) {
			sliceMin[d] = source.min(d);
			sliceMax[d] = source.max(d);
		}
	}

	this.dropSingltonDimensions = dropSingletonDimensions;
	this.slice = new FinalInterval(sliceMin, sliceMax);
	this.source = source;
}
 

/**
 * Gaussian Smooth of the input image using intermediate float format.
 * 
 * @param <T>
 * @param img
 * @param sigma
 * @return
 */
public static <T extends RealType<T> & NativeType<T>>
	Img<T> gaussianSmooth(RandomAccessibleInterval<T> img,
		double[] sigma)
{
	Interval interval = Views.iterable(img);

	ImgFactory<T> outputFactory = new ArrayImgFactory<>(Util.getTypeFromInterval(img));
	final long[] dim = new long[img.numDimensions()];
	img.dimensions(dim);
	Img<T> output = outputFactory.create(dim);

	final long[] pos = new long[img.numDimensions()];
	Arrays.fill(pos, 0);
	Localizable origin = new Point(pos);

	ImgFactory<FloatType> tempFactory = new ArrayImgFactory<>(new FloatType());
	RandomAccessible<T> input = Views.extendMirrorSingle(img);
	Gauss.inFloat(sigma, input, interval, output, origin, tempFactory);

	return output;
}
 

/**
 * Removes leading 0s from integral image after composite creation.
 *
 * @param input Input RAI (can be a RAI of Composite)
 * @return An extended and cropped version of input
 */
private <T> RandomAccessibleInterval<T> removeLeadingZeros(
	final RandomAccessibleInterval<T> input)
{
	// Remove 0s from integralImg by shifting its interval by +1
	final long[] min = Intervals.minAsLongArray(input);
	final long[] max = Intervals.maxAsLongArray(input);

	for (int d = 0; d < input.numDimensions(); ++d) {
		int correctedSpan = getShape().getSpan() - 1;
		min[d] += (1 + correctedSpan);
		max[d] -= correctedSpan;
	}

	// Define the Interval on the infinite random accessibles
	final FinalInterval interval = new FinalInterval(min, max);

	final RandomAccessibleInterval<T> extendedImg = Views.offsetInterval(Views
		.extendBorder(input), interval);
	return extendedImg;
}
 

/**
 * Converts the given image to a form renderable by scientific notebooks.
 *
 * @param <T>
 * @param source The image to render.
 * @return an object that the notebook knows how to draw onscreen.
 */
default <T extends RealType<T>> Object display(
        final RandomAccessibleInterval<T> source) {
    // NB: Assume <=3 samples in the 3rd dimension means channels. Of course,
    // we have no metadata with a vanilla RAI, but this is a best guess;
    // 3rd dimensions with >3 samples are probably something like Z or time.
    final int cAxis
            = //
            source.numDimensions() > 2 && source.dimension(2) <= 3 ? 2 : -1;

    return RAIToPNG(source, 0, 1, cAxis, ValueScaling.AUTO);
}
 

/**
 * Copies the outlines of the objects in the input interval into the output
 *
 * @param input an N-dimensional binary interval
 * @param excludeEdges are elements on stack edges outline or not
 *          <p>
 *          For example, a 2D square:<br>
 *          0 0 0 0<br>
 *          1 1 1 0<br>
 *          E 1 1 0<br>
 *          1 1 1 0<br>
 *          0 0 0 0<br>
 *          Element E is removed if parameter true, kept if false
 *          </p>
 * @param output outlines of the objects in interval
 */
@Override
public void compute(final RandomAccessibleInterval<B> input,
	final Boolean excludeEdges,
	final RandomAccessibleInterval<BitType> output)
{
	final ExtendedRandomAccessibleInterval<B, RandomAccessibleInterval<B>> extendedInput =
			extendInterval(input);
	final int nThreads = Runtime.getRuntime().availableProcessors();
	final ExecutorService pool = Executors.newFixedThreadPool(nThreads);
	final List<Future<?>> futures = new ArrayList<>();
	final long[] dimensions = new long[input.numDimensions()];
	input.dimensions(dimensions);
	final long size = stream(dimensions).reduce((a, b) -> a * b).orElse(0);
	final long perThread = size / nThreads;
	final long remainder = size % perThread;
	final IterableInterval<B> iterable = Views.iterable(input);

	for (int i = 0; i < nThreads; i++) {
		// Advancing a Cursor once sets it at [0, 0... 0]
		final long start = i * perThread + 1;
		final long share = i == 0 ? perThread + remainder : perThread;
		final Cursor<B> cursor = iterable.cursor();
		final OutOfBounds<B> inputAccess = extendedInput.randomAccess();
		final RandomAccess<BitType> outputAccess = output.randomAccess();
		final Runnable runnable = createTask(cursor, inputAccess,
				outputAccess, start, share);
		futures.add(pool.submit(runnable));
	}

	for (Future<?> future : futures) {
		try {
			future.get();
		} catch (InterruptedException | ExecutionException e) {
			logService.error(e);
		}
	}
}
 

public static <T extends RealType<T>> RandomAccessible<T> extendImageToSize(RandomAccessibleInterval<T> img, Dimensions extDims)
{
	int[] extEachSide = getSizeDifference(img, extDims);
	for (int i = 0; i< img.numDimensions(); i++){
		extEachSide[i] /= 2;
	}
	return new BlendedExtendedMirroredRandomAccesible2<T>(img, extEachSide);
}
 

public static <T extends RealType<T>> RandomAccessible<T> extendImageByFactor(RandomAccessibleInterval<T> img, int [] extension)
{
	int[] extEachSide = new int[img.numDimensions()];
	for (int i = 0; i <img.numDimensions(); i++){
		extEachSide[i] = (int) (img.dimension(i) < extension[i] ? img.dimension(i) : extension[i]);	
	}
	return new BlendedExtendedMirroredRandomAccesible2<T>(img, extEachSide);
}
 

@Override
public void compute(RandomAccessibleInterval<T> input, RandomAccessibleInterval<T> output) {
	RandomAccessibleInterval<T> in = input;
	for (int i = input.numDimensions() - 1; i >= 0; i--) {
		RandomAccessibleInterval<T> derivative = createRAI.calculate(input);
		if (dimension == i) {
			kernelBConvolveOp.compute(Views.interval(Views.extendMirrorDouble(in), input), derivative);
		} else {
			kernelAConvolveOps[i].compute(Views.interval(Views.extendMirrorDouble(in), input), derivative);
		}
		in = derivative;
	}
	addOp.compute(output, in, output);
}
 

/**
 * Creates a new {@link DataContainer} for a specific image channel
 * combination. We create default thresholds here that are the max and min of
 * the data type of the source image channels.
 *
 * @param src1 The channel one image source
 * @param src2 The channel two image source
 * @param ch1 The channel one image channel
 * @param ch2 The channel two image channel
 */
public DataContainer(RandomAccessibleInterval<T> src1,
		RandomAccessibleInterval<T> src2, int ch1, int ch2,
		String name1, String name2) {
	sourceImage1 = src1;
	sourceImage2 = src2;
	sourceImage1Name = name1;
	sourceImage2Name = name2;

	// create a mask that is true at all pixels.
	final long[] dims = new long[src1.numDimensions()];
	src1.dimensions(dims);
	mask = MaskFactory.createMask(dims, true);
	this.ch1 = ch1;
	this.ch2 = ch2;
	// fill mask dimension information, here the whole image
	maskBBOffset = new long[mask.numDimensions()];
	Arrays.fill(maskBBOffset, 0);
	maskBBSize = new long[mask.numDimensions()];
	mask.dimensions(maskBBSize);
	// indicated that there is actually no mask
	maskType = MaskType.None;

	maskHash = mask.hashCode();
	// create a jobName so ResultHandler instances can all use the same object
	// for the job name.
	jobName = "Colocalization_of_" + sourceImage1Name + "_versus_" + sourceImage2Name + "_" + maskHash;

	calculateStatistics();
}
 

/**
 * A method to combine a foreground image and a background image.
 * If data on the foreground image is above zero, it will be
 * placed on the background. While doing that, the image data from
 * the foreground is scaled to be in range of the background.
 */
public static <T extends RealType<T>> void combineImages(RandomAccessibleInterval<T> background,
		RandomAccessibleInterval<T> foreground) {
	final long[] dim = new long[ background.numDimensions() ];
	background.dimensions(dim);
	RandomAccessibleInterval<BitType> alwaysTrueMask = MaskFactory.createMask(dim, true);
	TwinCursor<T> cursor = new TwinCursor<T>(
			background.randomAccess(),
			foreground.randomAccess(),
			Views.iterable(alwaysTrueMask).localizingCursor());
	// find a scaling factor for scale forground range into background
	double bgMin = ImageStatistics.getImageMin(background).getRealDouble();
	double bgMax = ImageStatistics.getImageMax(background).getRealDouble();
	double fgMin = ImageStatistics.getImageMin(foreground).getRealDouble();
	double fgMax = ImageStatistics.getImageMax(foreground).getRealDouble();

	double scaling = (bgMax - bgMin ) / (fgMax - fgMin);
	// iterate over both images
	while (cursor.hasNext()) {
		cursor.fwd();
		T bgData = cursor.getFirst();
		double fgData = cursor.getSecond().getRealDouble() * scaling;
		if (fgData > 0.01) {
			/* if the foreground data is above zero, copy
			 * it to the background.
			 */
			bgData.setReal(fgData);
		}
	}
}
 

/**
 * Compute an ND inverse FFT
 */
@Override
public void mutate(final RandomAccessibleInterval<C> inout) {
	for (int d = inout.numDimensions() - 1; d >= 0; d--)
		FFTMethods.complexToComplex(inout, d, false, true, ts
			.getExecutorService());
}
 

@SuppressWarnings("unchecked")
@Override
public void compute(RandomAccessibleInterval<T> in, RandomAccessibleInterval<T> out) {
	final RandomAccessible<FloatType> convertedIn = Converters.convert(Views.extendMirrorDouble(in),
			converterToFloat, new FloatType());

	// compute partial derivative for each dimension
	RandomAccessibleInterval<FloatType> derivative0 = createImgOp.calculate();
	RandomAccessibleInterval<FloatType> derivative1 = createImgOp.calculate();

	// case of grayscale image
	if (in.numDimensions() == 2) {
		PartialDerivative.gradientCentralDifference(convertedIn, derivative0, 0);
		PartialDerivative.gradientCentralDifference(convertedIn, derivative1, 1);
	}
	// case of color image
	else {
		List<RandomAccessibleInterval<FloatType>> listDerivs0 = new ArrayList<>();
		List<RandomAccessibleInterval<FloatType>> listDerivs1 = new ArrayList<>();
		for (int i = 0; i < in.dimension(2); i++) {
			final RandomAccessibleInterval<FloatType> deriv0 = createImgOp.calculate();
			final RandomAccessibleInterval<FloatType> deriv1 = createImgOp.calculate();
			PartialDerivative.gradientCentralDifference(
					Views.interval(convertedIn, new long[] { 0, 0, i }, new long[] { in.max(0), in.max(1), i }),
					deriv0, 0);
			PartialDerivative.gradientCentralDifference(
					Views.interval(convertedIn, new long[] { 0, 0, i }, new long[] { in.max(0), in.max(1), i }),
					deriv1, 1);
			listDerivs0.add(deriv0);
			listDerivs1.add(deriv1);
		}
		derivative0 = Converters.convert(Views.collapse(Views.stack(listDerivs0)), converterGetMax,
				new FloatType());
		derivative1 = Converters.convert(Views.collapse(Views.stack(listDerivs1)), converterGetMax,
				new FloatType());
	}
	final RandomAccessibleInterval<FloatType> finalderivative0 = derivative0;
	final RandomAccessibleInterval<FloatType> finalderivative1 = derivative1;

	// compute angles and magnitudes
	final RandomAccessibleInterval<FloatType> angles = createImgOp.calculate();
	final RandomAccessibleInterval<FloatType> magnitudes = createImgOp.calculate();

	final CursorBasedChunk chunkable = new CursorBasedChunk() {

		@Override
		public void execute(long startIndex, long stepSize, long numSteps) {
			final Cursor<FloatType> cursorAngles = Views.flatIterable(angles).localizingCursor();
			final Cursor<FloatType> cursorMagnitudes = Views.flatIterable(magnitudes).localizingCursor();
			final Cursor<FloatType> cursorDerivative0 = Views.flatIterable(finalderivative0).localizingCursor();
			final Cursor<FloatType> cursorDerivative1 = Views.flatIterable(finalderivative1).localizingCursor();

			setToStart(cursorAngles, startIndex);
			setToStart(cursorMagnitudes, startIndex);
			setToStart(cursorDerivative0, startIndex);
			setToStart(cursorDerivative1, startIndex);

			for (long i = 0; i < numSteps; i++) {
				final float x = cursorDerivative0.get().getRealFloat();
				final float y = cursorDerivative1.get().getRealFloat();
				cursorAngles.get().setReal(getAngle(x, y));
				cursorMagnitudes.get().setReal(getMagnitude(x, y));

				cursorAngles.jumpFwd(stepSize);
				cursorMagnitudes.jumpFwd(stepSize);
				cursorDerivative0.jumpFwd(stepSize);
				cursorDerivative1.jumpFwd(stepSize);
			}
		}
	};

	ops().thread().chunker(chunkable, Views.flatIterable(magnitudes).size());

	// stores each Thread to execute
	final List<Callable<Void>> listCallables = new ArrayList<>();

	// compute descriptor (default 3x3, i.e. 9 channels: one channel for
	// each bin)
	final RectangleShape shape = new RectangleShape(spanOfNeighborhood, false);
	final NeighborhoodsAccessible<FloatType> neighborHood = shape.neighborhoodsRandomAccessible(angles);

	for (int i = 0; i < in.dimension(0); i++) {
		listCallables.add(new ComputeDescriptor(Views.interval(convertedIn, in), i, angles.randomAccess(),
				magnitudes.randomAccess(), (RandomAccess<FloatType>) out.randomAccess(),
				neighborHood.randomAccess()));
	}

	try {
		es.invokeAll(listCallables);
	} catch (final InterruptedException e) {
		throw new RuntimeException(e);
	}

	listCallables.clear();
}
 

@SuppressWarnings("unchecked")
private void testWithMask(final RandomAccessibleInterval<FloatType> in, final ImgLabeling<Integer, IntType> seeds) {
	// create mask which is 1 everywhere
	long[] dims = new long[in.numDimensions()];
	in.dimensions(dims);
	Img<BitType> mask = ArrayImgs.bits(dims);
	RandomAccess<BitType> raMask = mask.randomAccess();
	for (BitType b : mask) {
		b.setZero();
	}
	for (int x = 0; x < 10; x++) {
		for (int y = 0; y < 10; y++) {
			raMask.setPosition(new int[] { x, y });
			raMask.get().setOne();
		}
	}

	/*
	 * use 8-connected neighborhood
	 */
	// compute result without watersheds
	ImgLabeling<Integer, IntType> out = (ImgLabeling<Integer, IntType>) ops.run(WatershedSeeded.class, null, in,
			seeds, true, false, mask);

	assertResults(in, out, seeds, mask, false, true);

	// compute result with watersheds
	ImgLabeling<Integer, IntType> out2 = (ImgLabeling<Integer, IntType>) ops.run(WatershedSeeded.class, null, in,
			seeds, true, true, mask);

	assertResults(in, out2, seeds, mask, true, true);

	/*
	 * use 4-connected neighborhood
	 */
	// compute result without watersheds
	ImgLabeling<Integer, IntType> out3 = (ImgLabeling<Integer, IntType>) ops.run(WatershedSeeded.class, null, in,
			seeds, false, false, mask);

	assertResults(in, out3, seeds, mask, false, true);

	// compute result with watersheds
	ImgLabeling<Integer, IntType> out4 = (ImgLabeling<Integer, IntType>) ops.run(WatershedSeeded.class, null, in,
			seeds, false, true, mask);

	assertResults(in, out4, seeds, mask, true, true);
}
 

/**
 * create FFT memory, create FFT filter and run it
 */
@Override
public void compute(
	final RandomAccessibleInterval<I> input,
	final RandomAccessibleInterval<K> kernel,
	final RandomAccessibleInterval<O> output)
{

	//RandomAccessibleInterval<O> output = createOutput(input, kernel);

	final int numDimensions = input.numDimensions();

	// 1. Calculate desired extended size of the image

	final long[] paddedSize = new long[numDimensions];

	if (borderSize == null) {
		// if no border size was passed in, then extend based on kernel size
		for (int d = 0; d < numDimensions; ++d) {
			paddedSize[d] = (int) input.dimension(d) + (int) kernel.dimension(d) -
				1;
		}

	}
	else {
		// if borderSize was passed in
		for (int d = 0; d < numDimensions; ++d) {

			paddedSize[d] = Math.max(kernel.dimension(d) + 2 * borderSize[d], input
				.dimension(d) + 2 * borderSize[d]);
		}
	}

	RandomAccessibleInterval<I> paddedInput = padOp.calculate(input,
		new FinalDimensions(paddedSize));

	RandomAccessibleInterval<K> paddedKernel = padKernelOp.calculate(kernel,
		new FinalDimensions(paddedSize));

	RandomAccessibleInterval<C> fftInput = createOp.calculate(
		new FinalDimensions(paddedSize));

	RandomAccessibleInterval<C> fftKernel = createOp.calculate(
		new FinalDimensions(paddedSize));

	// TODO: in this case it is difficult to match the filter op in the
	// 'initialize' as we don't know the size yet, thus we can't create
	// memory
	// for the FFTs
	filter = createFilterComputer(paddedInput, paddedKernel, fftInput,
		fftKernel, output);

	filter.compute(paddedInput, paddedKernel, output);
}
 

@Override
public void compute(final RandomAccessibleInterval<T> input,
	final IterableInterval<DoubleType> tubeness)
{
	cancelReason = null;

	final int numDimensions = input.numDimensions();
	// Sigmas in pixel units.
	final double[] sigmas = new double[numDimensions];
	for (int d = 0; d < sigmas.length; d++) {
		final double cal = d < calibration.length ? calibration[d] : 1;
		sigmas[d] = sigma / cal;
	}

	/*
	 * Hessian.
	 */

	// Get a suitable image factory.
	final long[] gradientDims = new long[numDimensions + 1];
	final long[] hessianDims = new long[numDimensions + 1];
	for (int d = 0; d < numDimensions; d++) {
		hessianDims[d] = input.dimension(d);
		gradientDims[d] = input.dimension(d);
	}
	hessianDims[numDimensions] = numDimensions * (numDimensions + 1) / 2;
	gradientDims[numDimensions] = numDimensions;
	final Dimensions hessianDimensions = FinalDimensions.wrap(hessianDims);
	final FinalDimensions gradientDimensions = FinalDimensions.wrap(
		gradientDims);
	final ImgFactory<DoubleType> factory = ops().create().imgFactory(
		hessianDimensions);
	final Img<DoubleType> hessian = factory.create(hessianDimensions,
		new DoubleType());
	final Img<DoubleType> gradient = factory.create(gradientDimensions,
		new DoubleType());
	final Img<DoubleType> gaussian = factory.create(input, new DoubleType());

	// Handle multithreading.
	final int nThreads = Runtime.getRuntime().availableProcessors();
	final ExecutorService es = threadService.getExecutorService();

	try {
		// Hessian calculation.
		HessianMatrix.calculateMatrix(Views.extendBorder(input), gaussian,
			gradient, hessian, new OutOfBoundsBorderFactory<>(), nThreads, es,
			sigma);

		statusService.showProgress(1, 3);
		if (isCanceled()) return;

		// Hessian eigenvalues.
		final RandomAccessibleInterval<DoubleType> evs = TensorEigenValues
			.calculateEigenValuesSymmetric(hessian, TensorEigenValues
				.createAppropriateResultImg(hessian, factory, new DoubleType()),
				nThreads, es);

		statusService.showProgress(2, 3);
		if (isCanceled()) return;

		final AbstractUnaryComputerOp<Iterable<DoubleType>, DoubleType> method;
		switch (numDimensions) {
			case 2:
				method = new Tubeness2D(sigma);
				break;
			case 3:
				method = new Tubeness3D(sigma);
				break;
			default:
				System.err.println("Cannot compute tubeness for " + numDimensions +
					"D images.");
				return;
		}
		ops().transform().project(tubeness, evs, method, numDimensions);

		statusService.showProgress(3, 3);

		return;
	}
	catch (final IncompatibleTypeException | InterruptedException
			| ExecutionException e)
	{
		e.printStackTrace();
		return;
	}
}
 

/**
 * This test makes sure that a mask that is based on a lower dimension
 * image has the correct dimensionality.
 */
@Test
public void irregularRoiDimensionTest() {
	// load a 3D test image
	RandomAccessibleInterval<UnsignedByteType> img = positiveCorrelationImageCh1;
	final long width = img.dimension(0);
	final long height = img.dimension(1);
	final long slices = img.dimension(2);
	final long[] dimImg = new long[ img.numDimensions() ];
	img.dimensions(dimImg);
	// create a random noise 2D image -- set roiWidh/roiSize accordingly
	RandomAccessibleInterval<UnsignedByteType> maskSlice =
		TestImageAccessor.produceSticksNoiseImage( (int) width, (int) height, 50, 2, 10);
	RandomAccessibleInterval<BitType> mask =
			MaskFactory.createMask(dimImg, maskSlice);
	final long[] dimMask = new long[ mask.numDimensions() ];
	mask.dimensions(dimMask);
	// check the dimensions of the mask
	org.junit.Assert.assertArrayEquals(dimImg, dimMask);
	// make sure the mask actually got the same content on every slice
	final double[] offset = new double[ mask.numDimensions() ];
	Arrays.fill(offset, 0);
	double[] size = new double[ mask.numDimensions() ];
	size[0] = width;
	size[1] = height;
	size[2] = 1;
	RandomAccess<BitType> maskCursor = mask.randomAccess();
	RectangleRegionOfInterest roi = new RectangleRegionOfInterest( offset, size);
	Cursor<BitType> firstSliceCursor = roi.getIterableIntervalOverROI(mask).cursor();
	
	final long[] pos = new long[ mask.numDimensions() ];
	while (firstSliceCursor.hasNext()) {
		firstSliceCursor.fwd();
		firstSliceCursor.localize(pos);
		BitType maskValue = firstSliceCursor.get();
		// go through all slices
		for (int i=1; i<slices; ++i) {
			pos[2] = i;
			maskCursor.setPosition(pos);
			// compare the values and assume they are the same
			int cmp = maskCursor.get().compareTo(maskValue);
			assertEquals(0, cmp);
		}
	}
}
 

@Override
public void compute(RandomAccessibleInterval<I> in,
	RandomAccessibleInterval<K> kernel, RandomAccessibleInterval<O> out)
{
	// create FFT input memory if needed
	if (getFFTInput() == null) {
		setFFTInput(getCreateOp().calculate(in));
	}

	// create FFT kernel memory if needed
	if (getFFTKernel() == null) {
		setFFTKernel(getCreateOp().calculate(in));
	}

	// if a starting point for the estimate was not passed in then create
	// estimate Img and use the input as the starting point
	if (raiExtendedEstimate == null) {

		raiExtendedEstimate = createOp.calculate(in);

		copyOp.compute(in, raiExtendedEstimate);
	}

	// create image for the reblurred
	raiExtendedReblurred = createOp.calculate(in);

	// perform fft of psf
	fftKernelOp.compute(kernel, getFFTKernel());

	// -- perform iterations --

	for (int i = 0; i < getMaxIterations(); i++) {

		if (status != null) {
			status.showProgress(i, getMaxIterations());
		}

		// create reblurred by convolving kernel with estimate
		// NOTE: the FFT of the PSF of the kernel has been passed in as a
		// parameter. when the op was set up, and computed above, so we can use
		// compute
		convolverOp.compute(raiExtendedEstimate, this.raiExtendedReblurred);

		// compute correction factor
		rlCorrectionOp.compute(in, raiExtendedReblurred, raiExtendedReblurred);

		// perform update to calculate new estimate
		updateOp.compute(raiExtendedReblurred, raiExtendedEstimate);

		// apply post processing
		if (iterativePostProcessingOps != null) {
			for (UnaryInplaceOp<RandomAccessibleInterval<O>, RandomAccessibleInterval<O>> pp : iterativePostProcessingOps) {
				pp.mutate(raiExtendedEstimate);
			}
		}

		// accelerate the algorithm by taking a larger step
		if (getAccelerator() != null) {
			getAccelerator().mutate(raiExtendedEstimate);
		}
	}

	// -- copy crop padded back to original size

	final long[] start = new long[out.numDimensions()];
	final long[] end = new long[out.numDimensions()];

	for (int d = 0; d < out.numDimensions(); d++) {
		start[d] = 0;
		end[d] = start[d] + out.dimension(d) - 1;
	}

	copy2Op.compute(Views.interval(raiExtendedEstimate, new FinalInterval(start,
		end)), out);
}