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

public void initialize(RandomAccessibleInterval<T> yk_iterated) {
	if (yk_prediction == null) {

		long[] temp = new long[yk_iterated.numDimensions()];
		yk_iterated.dimensions(temp);

		FinalDimensions dims = new FinalDimensions(temp);

		yk_prediction = create.calculate(dims);
		xkm1_previous = create.calculate(dims);
		yk_prediction = create.calculate(dims);
		gk = create.calculate(dims);
		hk_vector = create.calculate(dims);

	}

}
 

/**
 * Inverts an image.
 *
 * @param <T> The images data type.
 * @param image The image to convert.
 * @return The inverted image.
 */
public static <T extends RealType<T> & NativeType<T>> RandomAccessibleInterval<T> invertImage(
		RandomAccessibleInterval<T> image) {
	Cursor<T> imgCursor = Views.iterable(image).localizingCursor();
	// invert the image
	long[] dim = new long[ image.numDimensions() ];
	image.dimensions(dim);
	ArrayImgFactory<T> imgFactory = new ArrayImgFactory<T>();
	RandomAccessibleInterval<T> invImg = imgFactory.create(
			dim, image.randomAccess().get().createVariable() ); // "Inverted " + image.getName());
	RandomAccess<T> invCursor = invImg.randomAccess();

	while (imgCursor.hasNext()) {
		imgCursor.fwd();
		invCursor.setPosition(imgCursor);
		invCursor.get().setReal( imgCursor.get().getMaxValue() - imgCursor.get().getRealDouble() );
	}

	return invImg;
}
 

/**
 * Tests a PredicateCursor by checking if all visited values are "true".
 * @throws MissingPreconditionException
 */
@Test
public void predicateCursorTest() throws MissingPreconditionException {
	// load a 3D test image
	RandomAccessibleInterval<UnsignedByteType> img = positiveCorrelationImageCh1;
	long[] roiOffset = createRoiOffset(img);
	long[] roiSize = createRoiSize(img);
	long[] dim = new long[ img.numDimensions() ];
	img.dimensions(dim);
	RandomAccessibleInterval<BitType> mask = MaskFactory.createMask(dim,
			roiOffset, roiSize);

	// create cursor to walk an image with respect to a mask
	final Predicate<BitType> predicate = new MaskPredicate();
	Cursor<BitType> roiCursor
		= new PredicateCursor<BitType>(
				Views.iterable(mask).localizingCursor(), predicate);

	// test if all visited voxels are "true"
	while (roiCursor.hasNext()) {
		roiCursor.fwd();
		assertTrue(roiCursor.get().get());
	}
}
 

/**
 * Converts an arbitrary image to a black/white version of it.
 * All image data lower or equal the splitValue will get black,
 * the rest will turn white.
 */
public static <T extends RealType<T> & NativeType<T>> RandomAccessibleInterval<T> makeBinaryImage(
		RandomAccessibleInterval<T> image, T splitValue) {
	Cursor<T> imgCursor = Views.iterable(image).localizingCursor();
	// make a new image of the same type, but binary
	long[] dim = new long[ image.numDimensions() ];
	image.dimensions(dim);
	ArrayImgFactory<T> imgFactory = new ArrayImgFactory<T>();
	RandomAccessibleInterval<T> binImg = imgFactory.create( dim,
			image.randomAccess().get().createVariable() ); // "Binary image of " + image.getName());
	RandomAccess<T> invCursor = binImg.randomAccess();

	while (imgCursor.hasNext()) {
		imgCursor.fwd();
		invCursor.setPosition(imgCursor);
		T currentValue = invCursor.get();
		if (currentValue.compareTo(splitValue) > 0)
			currentValue.setReal(  currentValue.getMaxValue() );
		else
			currentValue.setZero();
	}

	return binImg;
}
 

/**
 * Compute the inverse Hessian matrix from the the steepest descent images.
 * @param descent descent image
 * @param <T> pixel type
 * @return Hessian
 */
public static <T extends RealType< T >> double[][] computeInverseHessian(
		final RandomAccessibleInterval< T > descent)
{
	final int n = descent.numDimensions() - 1;
	final int numParameters = (int) descent.dimension( n );
	final long[] dim = new long[n + 1];
	descent.dimensions( dim );
	dim[n] = 1;
	final LocalizingIntervalIterator pos = new LocalizingIntervalIterator( dim );
	final RandomAccess< T > r = descent.randomAccess();
	final double[] deriv = new double[numParameters];
	final double[][] H = new double[numParameters][numParameters];
	while ( pos.hasNext() )
	{
		pos.fwd();
		r.setPosition( pos );
		for ( int p = 0; p < numParameters; ++p )
		{
			deriv[p] = r.get().getRealDouble();
			r.fwd( n );
		}
		for ( int i = 0; i < numParameters; ++i )
			for ( int j = 0; j < numParameters; ++j )
				H[i][j] += deriv[i] * deriv[j];
	}
	return new Matrix( H ).inverse().getArray();
}
 

/**
 * Creates a ROI offset array with a distance of 1/4 to the origin
 * in each dimension.
 */
protected <T extends RealType<T>> long[] createRoiOffset(RandomAccessibleInterval<T> img) {
	final long[] offset = new long[ img.numDimensions() ];
	img.dimensions(offset);
	for (int i=0; i<offset.length; i++) {
		offset[i] = Math.max(1, img.dimension(i) / 4);
	}
	return offset;
}
 

@Test
public void regularRoiPixelCountTest() {
	// load a 3D test image
	RandomAccessibleInterval<UnsignedByteType> img = positiveCorrelationImageCh1;
	final long[] roiOffset = createRoiOffset(img);
	final long[] roiSize = createRoiSize(img);
	final long width = img.dimension(0);
	final long height = img.dimension(1);

	RandomAccessibleInterval<UnsignedByteType> maskImg
		= TestImageAccessor.createRectengularMaskImage(width, height, roiOffset, roiSize);
	final long[] dim = new long[ img.numDimensions() ];
	img.dimensions(dim);
	RandomAccessibleInterval<BitType> mask =
			MaskFactory.createMask(dim, maskImg);

	TwinCursor<UnsignedByteType> cursor = new TwinCursor<UnsignedByteType>(
			img.randomAccess(),
			img.randomAccess(),
			Views.iterable(mask).localizingCursor());
	// calculate volume of mask bounding box
	long roiVolume = roiSize[0] * roiSize[1] * img.dimension(2);
	long count = 0;
	while (cursor.hasNext()) {
		cursor.fwd();
		count++;
	}

	assertEquals(roiVolume, count);
}
 

@Override
public void compute(final RandomAccessibleInterval<T> op,
	final RandomAccessibleInterval<T> r)
{
	final IterableInterval<T> iterOp = Views.flatIterable(op);
	final IterableInterval<T> iterR = Views.flatIterable(r);

	long[] dim = new long[r.numDimensions()];
	r.dimensions(dim);
	Cursor<T> rc = iterR.cursor();
	Cursor<T> opc = iterOp.localizingCursor();
	// Fill with non background marker
	while (rc.hasNext()) {
		rc.next().setOne();
	}

	rc.reset();
	boolean border;
	// Flood fill from every background border voxel
	while (rc.hasNext()) {
		rc.next();
		opc.next();
		if (rc.get().get() && !opc.get().get()) {
			border = false;
			for (int i = 0; i < r.numDimensions(); i++) {
				if (rc.getLongPosition(i) == 0 || rc.getLongPosition(i) == dim[i] -
					1)
				{
					border = true;
					break;
				}
			}
			if (border) {
				floodFillComp.compute(op, rc, r);
			}
		}
	}
}
 

/**
 * Creates a new {@link DataContainer} for a specific set of image and
 * channel combination. It will give access to the image according to
 * the region of interest (ROI) passed. Default thresholds, min, max and
 * mean will be set according to the ROI as well.
 *
 * @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
 * @param offset The offset of the ROI in each dimension
 * @param size The size of the ROI in each dimension
 */
public DataContainer(RandomAccessibleInterval<T> src1,
		RandomAccessibleInterval<T> src2, int ch1, int ch2,
		String name1, String name2,
		final long[] offset, final long size[])
		throws MissingPreconditionException {
	sourceImage1 = src1;
	sourceImage2 = src2;
	sourceImage1Name = name1;
	sourceImage2Name = name2;
	
	final int numDims = src1.numDimensions();
	final long[] dim = new long[numDims];
	src1.dimensions(dim);
	long[] roiOffset = new long[numDims];
	long[] roiSize = new long[numDims];

	adjustRoiOffset(offset, roiOffset, dim);
	adjustRoiSize(size, roiSize, dim, roiOffset);

	// create a mask that is valid everywhere
	mask = MaskFactory.createMask(dim, roiOffset, roiSize);
	maskBBOffset = roiOffset.clone();
	maskBBSize = roiSize.clone();
	// this constructor only supports regular masks
	maskType = MaskType.Regular;

	this.ch1 = ch1;
	this.ch2 = ch2;

	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();
}
 

/**
 * Counts the number of foreground boxes in a grid over the given interval
 *
 * @param sizes Sizes of the interval's dimensions in pixels
 * @param translation Translation of the box grid in each dimension
 * @param boxSize Size of a box (n * n * ... n)
 */
private long countForegroundBoxes(final RandomAccessibleInterval<B> input,
		final long boxSize, final long[] sizes, final long[] translation) {
	final long[] dimensions = new long[input.numDimensions()];
	input.dimensions(dimensions);
	long boxCount = 1;
	for (int i = 0; i < sizes.length; i++) {
		boxCount *= Math.ceil((1.0 * sizes[i] - translation[i]) / boxSize);
	}
	int nThreads = countNThreads(boxCount, boxSize);
	final long boxesPerThread = Math.max(1, boxCount / nThreads);
	final long remainder = boxCount % nThreads;

	final ExecutorService pool = Executors.newFixedThreadPool(nThreads);
	final List<Future<Long>> futures = new ArrayList<>(nThreads);
	for (int i = 0; i < nThreads; i++) {
		final long first = i * boxesPerThread;
		final long boxes = i < nThreads - 1 ? boxesPerThread :
				boxesPerThread + remainder;
		final Callable<Long> task = createCalculationTask(dimensions,
				boxSize, translation, first, input, boxes);
		futures.add(pool.submit(task));
	}

	long foregroundBoxes = 0;
	for (Future<Long> future : futures) {
		try {
			foregroundBoxes += future.get();
		} catch (InterruptedException | ExecutionException e) {
			e.printStackTrace();
		}
	}
	return foregroundBoxes;
}
 

protected static <T extends RealType<T>> LabelRegion<String> createLabelRegion(
	final RandomAccessibleInterval<T> interval, final float min, final float max, long... dims)
{
	if (dims == null || dims.length == 0) {
		dims = new long[interval.numDimensions()];
		interval.dimensions(dims);
	}
	final ImgLabeling<String, IntType> labeling = 
		new ImgLabeling<>(ArrayImgs.ints(dims));

	final RandomAccess<LabelingType<String>> ra = labeling.randomAccess();
	final RandomAccessibleIntervalCursor<T> c = new RandomAccessibleIntervalCursor<>(interval);		
	final long[] pos = new long[labeling.numDimensions()];
	while (c.hasNext()) {
		final T item = c.next();
		final float value = item.getRealFloat();
		if (value >= min && value <= max) {
			c.localize(pos);
			ra.setPosition(pos);
			ra.get().add("1");
		}
	}
	final LabelRegions<String> labelRegions = new LabelRegions<>(labeling);

	return labelRegions.getLabelRegion("1");

}
 

/**
 * This test creates first an "always true" mask and count the data
 * values. There should be as many as the number of vocels in total.
 * After that an "always false" mask is created. The predicate cursor
 * there should not return any values.
 */
@Test
public void simpleMaskCreationTest() {
	final RandomAccessibleInterval<UnsignedByteType> img = positiveCorrelationImageCh1;
	// first, create an always true mask
	final long[] dim = new long[ img.numDimensions() ];
	img.dimensions(dim);
	RandomAccessibleInterval<BitType> mask = MaskFactory.createMask(dim, true);
	final Predicate<BitType> predicate = new MaskPredicate();
	Cursor<BitType> cursor
		= new PredicateCursor<BitType>(
				Views.iterable(mask).localizingCursor(), predicate);
	// iterate over mask and count values
	long count = 0;
	while (cursor.hasNext()) {
		cursor.fwd();
		count++;
		assertTrue(cursor.get().get());
	}
	assertEquals(ImageStatistics.getNumPixels(mask), count);

	// second, create an always false mask
	mask = MaskFactory.createMask(dim, false);
	cursor = new PredicateCursor<BitType>(
			Views.iterable(mask).localizingCursor(), predicate);
	// iterate over mask and count values
	count = 0;
	while (cursor.hasNext()) {
		cursor.fwd();
		count++;
	}
	assertEquals(0, count);
}
 

/**
 * 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();
}
 

/**
 * Create a new mask based on a threshold condition for two images.
 */
public static<T extends RealType< T >> RandomAccessibleInterval<BitType> createMask(
		RandomAccessibleInterval<T> ch1, RandomAccessibleInterval<T> ch2,
		T threshold1, T threshold2, ThresholdMode tMode, CombinationMode cMode) {
	
	final long[] dims = new long[ ch1.numDimensions() ];
	ch1.dimensions(dims);
	RandomAccessibleInterval<BitType> mask = createMask(dims);
	Cursor<T> cursor1 = Views.iterable(ch1).cursor();
	Cursor<T> cursor2 = Views.iterable(ch2).cursor();
	Cursor<BitType> maskCursor = Views.iterable(mask).cursor();
	
	while (cursor1.hasNext() && cursor2.hasNext() && maskCursor.hasNext()) {
		cursor1.fwd();
		cursor2.fwd();
		maskCursor.fwd();
		
		boolean ch1Valid, ch2Valid;
		
		T data1 = cursor1.get();
		T data2 = cursor2.get();
		
		// get relation to threshold
		if (tMode == ThresholdMode.Above) {
			ch1Valid = data1.compareTo(threshold1) > 0;
			ch2Valid = data2.compareTo(threshold2) > 0;
		} else if (tMode == ThresholdMode.Below) {
			ch1Valid = data1.compareTo(threshold1) < 0;
			ch2Valid = data2.compareTo(threshold2) < 0;
		} else {
			throw new UnsupportedOperationException();
		}
		
		BitType maskData = maskCursor.get();
		
		// combine the results into mask
		if (cMode == CombinationMode.AND) {
			maskData.set( ch1Valid && ch2Valid );
		} else if (cMode == CombinationMode.OR) {
			maskData.set( ch1Valid || ch2Valid );
		} else if (cMode == CombinationMode.NONE) {
			maskData.set( !(ch1Valid || ch2Valid) );
		} else {
			throw new UnsupportedOperationException();
		}
	}
	
	return mask;
}
 

public static < S extends RealType< S > > Img< S > computeMaxProjection( final RandomAccessibleInterval< S > avgPSF, final ImgFactory< S > factory, int minDim )
{
	final long[] dimensions = new long[ avgPSF.numDimensions() ];
	avgPSF.dimensions( dimensions );
	
	if ( minDim < 0 )
	{
		long minSize = dimensions[ 0 ];
		minDim = 0;
		
		for ( int d = 0; d < dimensions.length; ++d )
		{
			if ( avgPSF.dimension( d ) < minSize )
			{
				minSize = avgPSF.dimension( d );
				minDim = d;
			}
		}
	}
	
	final long[] projDim = new long[ dimensions.length - 1 ];
	
	int dim = 0;
	long sizeProjection = 0;
	
	// the new dimensions
	for ( int d = 0; d < dimensions.length; ++d )
		if ( d != minDim )
			projDim[ dim++ ] = dimensions[ d ];
		else
			sizeProjection = dimensions[ d ];
	
	final Img< S > proj = factory.create( projDim, Views.iterable( avgPSF ).firstElement() );
	
	final RandomAccess< S > psfIterator = avgPSF.randomAccess();
	final Cursor< S > projIterator = proj.localizingCursor();
	
	final int[] tmp = new int[ avgPSF.numDimensions() ];
	
	while ( projIterator.hasNext() )
	{
		projIterator.fwd();

		dim = 0;
		for ( int d = 0; d < dimensions.length; ++d )
			if ( d != minDim )
				tmp[ d ] = projIterator.getIntPosition( dim++ );

		tmp[ minDim ] = -1;
		
		double maxValue = -Double.MAX_VALUE;
		
		psfIterator.setPosition( tmp );
		for ( int i = 0; i < sizeProjection; ++i )
		{
			psfIterator.fwd( minDim );
			final double value = psfIterator.get().getRealDouble();
			
			if ( value > maxValue )
				maxValue = value;
		}
		
		projIterator.get().setReal( maxValue );
	}
	
	return proj;
}
 

/**
 * This test test if a regular mask is created by the MaskFactory
 * correctly. First, the dimensions are checked, they must be the
 * same as the original images ones. Then it is checked if all
 * values in the mask image have the value they should have. For
 * a regular ROI this is easy to tell as one can calculate it out
 * of the position.
 */
@Test
public void regularMaskCreationTest() throws MissingPreconditionException {
	// load a 3D test image
	RandomAccessibleInterval<UnsignedByteType> img = positiveCorrelationImageCh1;
	final long[] roiOffset = createRoiOffset(img);
	final long[] roiSize = createRoiSize(img);
	final long[] dim = new long[ img.numDimensions() ];
	img.dimensions(dim);
	RandomAccessibleInterval<BitType> mask = MaskFactory.createMask(dim,
			roiOffset, roiSize);

	// is the number of dimensions the same as in the image?
	final long[] dimMask = new long[ mask.numDimensions() ];
	mask.dimensions(dimMask);
	assertTrue( Arrays.equals(dim, dimMask) );

	// go through the mask and check if all valid points are in the ROI
	long[] pos = new long[ img.numDimensions() ];
	final Cursor<BitType> cursor = Views.iterable(mask).localizingCursor();
	while ( cursor.hasNext() ) {
		cursor.fwd();
		cursor.localize(pos);
		// get values in mask image
		boolean onInMask = cursor.get().get();
		// calculate value that the current point *should* have
		boolean onInROI = true;
		for(int i=0; i<pos.length; ++i)
			onInROI &= pos[i] >= roiOffset[i] && pos[i] < (roiOffset[i] + roiSize[i]);
		// both values must match
		assertTrue(onInMask == onInROI);
	}

	/* go once more trough the image wrt. the mask to build a
	 * bounding box
	 */
	// create cursor to walk an image with respect to a mask
	final Predicate<BitType> predicate = new MaskPredicate();
	Cursor<BitType> roiCursor
		= new PredicateCursor<BitType>(
				Views.iterable(mask).localizingCursor(), predicate);
	long[] min = new long[ mask.numDimensions() ];
	long[] max = new long[ mask.numDimensions() ];
	Arrays.fill(min, Integer.MAX_VALUE);
	Arrays.fill(max, Integer.MIN_VALUE);
	while (roiCursor.hasNext()) {
		roiCursor.fwd();
		roiCursor.localize(pos);
		for (int i=0; i<pos.length; i++) {
			if (pos[i] < min[i])
				min[i] = pos[i];
			if (pos[i] > max[i])
				max[i] = pos[i];
		}
	}
	// the bounding box min should equal the ROI offset
	assertTrue(Arrays.equals(min, roiOffset));
	// create theoretical bounding box max and check it
	long[] roiMax = roiOffset.clone();
	for (int i=0; i<roiMax.length; i++)
		roiMax[i] += roiSize[i] - 1;
	assertTrue(Arrays.equals(max, roiMax));
}
 

/**
 * Tests if a masked walk over an image refers to the correct data
 * by copying the data to a separate image and then compare it with
 * the original image data. The position data in the original image
 * is calculated based on the ROI offset and the relative position
 * in the copied ROI image.
 * @throws MissingPreconditionException
 */
@Test
public void maskContentTest() throws MissingPreconditionException {
	// load a 3D test image
	final RandomAccessibleInterval<UnsignedByteType> img = positiveCorrelationImageCh1;
	final long[] roiOffset = createRoiOffset(img);
	final long[] roiSize = createRoiSize(img);
	final long[] dim = new long[ img.numDimensions() ];
	img.dimensions(dim);
	final RandomAccessibleInterval<BitType> mask =
			MaskFactory.createMask(dim, roiOffset, roiSize);

	// create cursor to walk an image with respect to a mask
	TwinCursor<UnsignedByteType> cursor = new TwinCursor<UnsignedByteType>(
			img.randomAccess(), img.randomAccess(),
			Views.iterable(mask).localizingCursor());

	// create an image for the "clipped ROI"
	ImgFactory<UnsignedByteType> maskFactory =
			new ArrayImgFactory<UnsignedByteType>();
	RandomAccessibleInterval<UnsignedByteType> clippedRoiImage =
			maskFactory.create( roiSize, new UnsignedByteType() ); //  "Clipped ROI" );
	RandomAccess<UnsignedByteType> outputCursor =
			clippedRoiImage.randomAccess();

	// copy ROI data to new image
	long[] pos = new long[ clippedRoiImage.numDimensions() ];
	while (cursor.hasNext()) {
		cursor.fwd();
		cursor.localize(pos);
		// shift position by offset
		for (int i=0; i<pos.length; i++) {
			pos[i] = pos[i] - roiOffset[i];
		}
		outputCursor.setPosition(pos);
		outputCursor.get().set( cursor.getFirst() );
	}

	/* go through the clipped ROI and compare the date to offset values
	 * of the original data.
	 */
	Cursor<UnsignedByteType> roiCopyCursor =
			Views.iterable(clippedRoiImage).localizingCursor();
	RandomAccess<UnsignedByteType> imgCursor =
			img.randomAccess();
	// create variable for summing up and set it to zero
	double sum = 0;
	pos = new long [ clippedRoiImage.numDimensions() ];
	while (roiCopyCursor.hasNext()) {
		roiCopyCursor.fwd();
		roiCopyCursor.localize(pos);
		// shift position by offset
		for (int i=0; i<pos.length; i++) {
			pos[i] = pos[i] + roiOffset[i];
		}
		// set position in original image
		imgCursor.setPosition(pos);
		// get ROI and original image data
		double roiData = roiCopyCursor.get().getRealDouble();
		double imgData = imgCursor.get().getRealDouble();
		// sum up the difference
		double diff = roiData - imgData;
		sum += diff * diff;
	}

	// check if sum is zero
	assertTrue("The sum of squared differences was " + sum + ".", Math.abs(sum) < 0.00001);
}
 

/**
 * 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(final RandomAccessibleInterval<I> input, final RandomAccessibleInterval<O> output) {

	final long[] size = new long[input.numDimensions()];
	input.dimensions(size);

	final RandomAccess<O> outputRA = output.randomAccess();
	final Cursor<I> inputCursor = Views.iterable(input).localizingCursor();
	final long[] currentPos = new long[input.numDimensions()];
	final long[] neighborhoodPos = new long[input.numDimensions()];
	final long[] neighborhoodMin = new long[input.numDimensions()];
	final long[] neighborhoodMax = new long[input.numDimensions()];
	Neighborhood<I> neighborhood;
	Cursor<I> neighborhoodCursor;
	final RectangleNeighborhoodFactory<I> fac = RectangleNeighborhood.factory();
	while (inputCursor.hasNext()) {
		inputCursor.fwd();
		inputCursor.localize(currentPos);
		double distance;
		inputCursor.localize(neighborhoodMin);
		inputCursor.localize(neighborhoodMax);
		neighborhoodMin[0] = Math.max(0, neighborhoodMin[0] - radius);
		neighborhoodMin[1] = Math.max(0, neighborhoodMin[1] - radius);
		neighborhoodMax[0] = Math.min(input.max(0), neighborhoodMax[0] + radius);
		neighborhoodMax[1] = Math.min(input.max(1), neighborhoodMax[1] + radius);
		final Interval interval = new FinalInterval(neighborhoodMin, neighborhoodMax);
		neighborhood = fac.create(currentPos, neighborhoodMin, neighborhoodMax, interval, input.randomAccess());
		neighborhoodCursor = neighborhood.localizingCursor();
		double weight, v = 0.0;
		double w = 0.0;
		do {
			neighborhoodCursor.fwd();
			neighborhoodCursor.localize(neighborhoodPos);
			distance = getDistance(currentPos, neighborhoodPos);
			weight = gauss(distance, sigmaS);// spatial kernel

			distance = Math.abs(inputCursor.get().getRealDouble() - neighborhoodCursor.get().getRealDouble());// intensity
																												// difference
			weight *= gauss(distance, sigmaR);// range kernel, then exponent addition

			v += weight * neighborhoodCursor.get().getRealDouble();
			w += weight;
		} while (neighborhoodCursor.hasNext());
		outputRA.setPosition(currentPos);
		outputRA.get().setReal(v / w);
	}

}
 

@SuppressWarnings("unchecked")
@Override
public void compute(final RandomAccessibleInterval<I> input,
	final O output)
{

	// List to store all directions occuring within the image on borders
	ArrayList<DoubleType> dirList = new ArrayList<>();

	// Dimension of input region
	long[] dims = new long[input.numDimensions()];
	input.dimensions(dims);

	// create image for derivations in x and y direction
	Img<I> derX = imgCreator.calculate(input);
	Img<I> derY = imgCreator.calculate(input);

	// calculate derivations in x and y direction
	PartialDerivative.gradientCentralDifference2(Views.extendMirrorSingle(
		input), derX, 0);
	PartialDerivative.gradientCentralDifference2(Views.extendMirrorSingle(
		input), derY, 1);

	// calculate theta at each position: theta = atan(dX/dY) + pi/2
	Cursor<I> cX = derX.cursor();
	Cursor<I> cY = derY.cursor();

	// for each position calculate magnitude and direction
	while (cX.hasNext()) {
		cX.next();
		cY.next();

		double dx = cX.get().getRealDouble();
		double dy = cY.get().getRealDouble();

		double dir = 0.0;
		double mag = 0.0;

		mag = Math.sqrt(dx * dx + dy * dy);

		if (dx != 0 && mag > 0.0) {
			dir = Math.atan(dy / dx) + Math.PI / 2;
			dirList.add(new DoubleType(dir));
		}
	}

	// No directions: output is zero
	if (dirList.isEmpty()) {
		output.setReal(0.0);
	}
	// Otherwise compute histogram over all occuring directions
	// and calculate inverse second moment on it as output
	else {
		Histogram1d<Integer> hist = histOp.calculate(dirList);
		double std = stdOp.calculate(hist).getRealDouble();
		output.setReal(1 / std);
	}
}