Java Code Examples for net.imglib2.img.ImgFactory#create()

/**
 * 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>> RandomAccessibleInterval<T> gaussianSmooth(
		RandomAccessibleInterval<T> img, double[] sigma) {
	Interval interval = Views.iterable(img);

	ImgFactory<T> outputFactory = new ArrayImgFactory<T>(Util.getTypeFromInterval(img));
	final long[] dim = new long[ img.numDimensions() ];
	img.dimensions(dim);
	RandomAccessibleInterval<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<FloatType>(new FloatType());
	RandomAccessible<T> input = Views.extendMirrorSingle(img);
	Gauss.inFloat(sigma, input, interval, output, origin, tempFactory);

	return output;
}
 

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

public static <T extends RealType<T>, S extends RealType<S>, R extends RealType<R>, C extends ComplexType<C>> RandomAccessibleInterval<R> calculatePCM(
		RandomAccessibleInterval<T> img1, RandomAccessibleInterval<S> img2, int[] extension,
		ImgFactory<R> factory, R type, ImgFactory<C> fftFactory, C fftType, ExecutorService service){

	
	// TODO: Extension absolute per dimension in pixels, i.e. int[] extension
	// TODO: not bigger than the image dimension because the second mirroring is identical to the image
	
	Dimensions extSize = PhaseCorrelation2Util.getExtendedSize(img1, img2, extension);
	long[] paddedDimensions = new long[extSize.numDimensions()];
	long[] fftSize = new long[extSize.numDimensions()];
	FFTMethods.dimensionsRealToComplexFast(extSize, paddedDimensions, fftSize);
	
	RandomAccessibleInterval<C> fft1 = fftFactory.create(fftSize, fftType);
	RandomAccessibleInterval<C> fft2 = fftFactory.create(fftSize, fftType);
	
	FFT.realToComplex(Views.interval(PhaseCorrelation2Util.extendImageByFactor(img1, extension), 
			FFTMethods.paddingIntervalCentered(img1, new FinalInterval(paddedDimensions))), fft1, service);
	FFT.realToComplex(Views.interval(PhaseCorrelation2Util.extendImageByFactor(img2, extension), 
			FFTMethods.paddingIntervalCentered(img2, new FinalInterval(paddedDimensions))), fft2, service);
	
	RandomAccessibleInterval<R> pcm = calculatePCMInPlace(fft1, fft2, factory, type, service);
	return pcm;
	
}
 

public static < T extends RealType< T > > RandomAccessibleInterval< T > simple2x( final RandomAccessibleInterval<T> input, final ImgFactory< T > imgFactory, final boolean[] downsampleInDim )
{
	RandomAccessibleInterval< T > src = input;

	for ( int d = 0; d < input.numDimensions(); ++d )
		if ( downsampleInDim[ d ] )
		{
			final long dim[] = new long[ input.numDimensions() ];

			for ( int e = 0; e < input.numDimensions(); ++e )
			{
				if ( e == d )
					dim[ e ] = src.dimension( e ) / 2;
				else
					dim[ e ] = src.dimension( e );
			}

			final Img< T > img = imgFactory.create( dim, Views.iterable( input ).firstElement() );
			simple2x( src, img, d );
			src = img;
		}

	return src;
}
 

/**
 * This method creates a noise image that has a specified mean. Every pixel
 * has a value uniformly distributed around mean with the maximum spread
 * specified.
 *
 * @return IllegalArgumentException if specified means and spreads are not
 *         valid
 */
public static <T extends RealType<T> & NativeType<T>>
	Img<T> produceMeanBasedNoiseImage(T type, int width,
		int height, double mean, double spread, double[] smoothingSigma,
		long seed) throws IllegalArgumentException
{
	if (mean < spread || (mean + spread) > type.getMaxValue()) {
		throw new IllegalArgumentException(
			"Mean must be larger than spread, and mean plus spread must be smaller than max of the type");
	}
	// create the new image
	ImgFactory<T> imgFactory = new ArrayImgFactory<>(type);
	Img<T> noiseImage = imgFactory.create(width, height);

	Random r = new Random(seed);
	for (T value : Views.iterable(noiseImage)) {
		value.setReal(mean + ((r.nextDouble() - 0.5) * spread));
	}

	// TODO: call Ops filter.gauss instead
	return gaussianSmooth(noiseImage, smoothingSigma);
}
 

final protected static Thread getCUDAThread2(
		final AtomicInteger ai, final ImgFactory< FloatType > blockFactory, final Block[] blocks, final int[] blockSize,
		final Img< FloatType > image, final Img< FloatType > result, final int deviceId, final Img< FloatType > kernel2 )
{
	final Thread cudaThread2 = new Thread( new Runnable()
	{
		public void run()
		{
			final Img< FloatType > block = blockFactory.create( Util.int2long( blockSize ), new FloatType() );

			int i;

			while ( ( i = ai.getAndIncrement() ) < blocks.length )
				convolve2BlockCUDA( blocks[ i ], deviceId, image, result, block, kernel2 );
		}
	});
	
	return cudaThread2;
}
 

final protected static Thread getCUDAThread1(
		final AtomicInteger ai, final ImgFactory< FloatType > blockFactory, final Block[] blocks, final int[] blockSize,
		final Img< FloatType > image, final Img< FloatType > result, final int deviceId, final Img< FloatType > kernel1 )
{
	final Thread cudaThread1 = new Thread( new Runnable()
	{
		public void run()
		{
			final Img< FloatType > block = blockFactory.create( Util.int2long( blockSize ), new FloatType() );

			int i;

			while ( ( i = ai.getAndIncrement() ) < blocks.length )
				convolve1BlockCUDA( blocks[ i ], deviceId, image, result, block, kernel1, i );
		}
	});
	
	return cudaThread1;
}
 

/**
 * This method duplicates the given images, but respects ROIs if present.
 * Meaning, a sub-picture will be created when source images are
 * ROI/MaskImages.
 * 
 * @throws MissingPreconditionException
 */
protected RandomAccessibleInterval<T> createMaskImage(
	final RandomAccessibleInterval<T> image,
	final RandomAccessibleInterval<BitType> mask, final long[] offset,
	final long[] size) throws MissingPreconditionException
{
	final long[] pos = new long[image.numDimensions()];
	// sanity check
	if (pos.length != offset.length || pos.length != size.length) {
		throw new MissingPreconditionException(
			"Mask offset and size must be of same dimensionality like image.");
	}
	// use twin cursor for only one image
	final TwinCursor<T> cursor = new TwinCursor<>(image.randomAccess(), //
		image.randomAccess(), Views.iterable(mask).localizingCursor());
	// prepare output image
	final ImgFactory<T> maskFactory = new ArrayImgFactory<>();
	// Img<T> maskImage = maskFactory.create( size, name );
	final RandomAccessibleInterval<T> maskImage = maskFactory.create(size, //
		image.randomAccess().get().createVariable());
	final RandomAccess<T> maskCursor = maskImage.randomAccess();
	// go through the visible data and copy it to the output
	while (cursor.hasNext()) {
		cursor.fwd();
		cursor.localize(pos);
		// shift coordinates by offset
		for (int i = 0; i < pos.length; ++i) {
			pos[i] = pos[i] - offset[i];
		}
		// write out to correct position
		maskCursor.setPosition(pos);
		maskCursor.get().set(cursor.getFirst());
	}

	return maskImage;
}
 

/**
  * To randomize blockwise we enumerate the blocks, shuffle that list and
  * write the data to their new position based on the shuffled list.
  */
 protected Img<T> generateRandomImageStack(Img<T> img, int[] blockDimensions) {
int numberOfDimensions = Math.min(img.numDimensions(), blockDimensions.length);
int numberOfBlocks = 0;
long[] numberOfBlocksPerDimension = new long[numberOfDimensions];

for (int i = 0 ; i<numberOfDimensions; i++){
	if (img.dimension(i) % blockDimensions[i] != 0){
		System.out.println("sorry, for now image dims must be divisable by block size");
		return null;
	}
	numberOfBlocksPerDimension[i] = img.dimension(i) / blockDimensions[i];
	numberOfBlocks *= numberOfBlocksPerDimension[i];
}
List<Integer> allTheBlocks = new ArrayList<Integer>(numberOfBlocks);
for (int i = 0; i<numberOfBlocks; i++){
	allTheBlocks.add(new Integer(i));
}
Collections.shuffle(allTheBlocks, new Random());
Cursor<T> cursor = img.cursor();

// create factories for new image stack
//ContainerFactory containerFactory = new ImagePlusContainerFactory();
ImgFactory<T> imgFactory = new ArrayImgFactory<T>();
//new ImageFactory<T>(cursor.getType(), containerFactory);

// create a new stack for the random images
final long[] dim = new long[ img.numDimensions() ];
img.dimensions(dim);
Img<T> randomStack = imgFactory.create(dim, img.firstElement().createVariable());

// iterate over image data
while (cursor.hasNext()) {
	cursor.fwd();
	T type = cursor.get();
	// type.getRealDouble();
}

return randomStack;
 }
 

protected Img< FloatType > approximateEntropy(
		final RandomAccessibleInterval< FloatType > input,
		final ImgFactory< ComplexFloatType > imgFactory,
		final double[] sigma1,
		final double[] sigma2 )

{
	// the result
	ImgFactory<FloatType> f;
	try { f = imgFactory.imgFactory( new FloatType() ); } catch (IncompatibleTypeException e) { f = new ArrayImgFactory< FloatType >(); }
	
	final Img< FloatType > conv = f.create( input, new FloatType() );
	
	// compute I*sigma1
	FFTConvolution< FloatType > fftConv = new FFTConvolution<FloatType>( input, createGaussianKernel( sigma1 ), conv, imgFactory );
	fftConv.convolve();
	
	// compute ( I - I*sigma1 )^2
	final Cursor< FloatType > c = conv.cursor();
	final RandomAccess< FloatType > r = input.randomAccess();
	
	while ( c.hasNext() )
	{
		c.fwd();
		r.setPosition( c );
		
		final float diff = c.get().get() - r.get().get();
		c.get().set( diff * diff );
	}
	
	// compute ( ( I - I*sigma1 )^2 ) * sigma2
	fftConv = new FFTConvolution<FloatType>( conv, createGaussianKernel( sigma2 ), imgFactory );
	fftConv.convolve();

	// normalize to [0...1]
	FusionHelper.normalizeImage( conv );

	return conv;
}
 

/**
 * Generates a Perlin noise image. It is based on Ken Perlin's
 * reference implementation (ImprovedNoise class) and a small
 * bit of Kas Thomas' sample code (http://asserttrue.blogspot.com/).
 */
public static <T extends RealType<T> & NativeType<T>> RandomAccessibleInterval<T> producePerlinNoiseImage(T type, int width,
		int height, double z, double scale) {
	// create the new image
	ImgFactory<T> imgFactory = new ArrayImgFactory<T>();
	RandomAccessibleInterval<T> noiseImage = imgFactory.create( new int[] {width, height}, type);
	Cursor<T> noiseCursor = Views.iterable(noiseImage).localizingCursor();

	double xOffset = Math.random() * (width*width);
	double yOffset = Math.random() * (height*height);

	while (noiseCursor.hasNext()) {
		noiseCursor.fwd();
		double x = (noiseCursor.getDoublePosition(0) + xOffset) * scale;
		double y = (noiseCursor.getDoublePosition(1) + yOffset) * scale;

		float t = (float)ImprovedNoise.noise( x, y, z);

		// ImprovedNoise.noise returns a float in the range [-1..1],
		// whereas we want a float in the range [0..1], so:
                       t = (1 + t) * 0.5f;

                       noiseCursor.get().setReal(t);
	}

	//return gaussianSmooth(noiseImage, imgFactory, smoothingSigma);
	return noiseImage;
}
 

public static void main( String[] args )
{
	new ImageJ();
	
	// test blending
	ImgFactory< FloatType > f = new ArrayImgFactory< FloatType >();
	Img< FloatType > img = f.create( new int[] { 400, 400 }, new FloatType() ); 
	
	Cursor< FloatType > c = img.localizingCursor();
	final int numDimensions = img.numDimensions();
	final double[] tmp = new double[ numDimensions ];
	
	// for blending
	final long[] dimensions = new long[ numDimensions ];
	img.dimensions( dimensions );
	final float percentScaling = 0.2f;
	final double[] border = new double[ numDimensions ];
				
	while ( c.hasNext() )
	{
		c.fwd();
		
		for ( int d = 0; d < numDimensions; ++d )
			tmp[ d ] = c.getFloatPosition( d );
		
		c.get().set( (float)BlendingPixelFusion.computeWeight( tmp, dimensions, border, percentScaling ) );
	}
	
	ImageJFunctions.show( img );
	Log.debug( "done" );
}
 

public static <T extends ComplexType<T> & NativeType<T>, S extends ComplexType<S> & NativeType <S>, R extends RealType<R>> RandomAccessibleInterval<R> calculatePCM(
		RandomAccessibleInterval<T> fft1, RandomAccessibleInterval<S> fft2, ImgFactory<R> factory, R type, ExecutorService service){
	
	long[] paddedDimensions = new long[fft1.numDimensions()];
	long[] realSize = new long[fft1.numDimensions()];
	
	FFTMethods.dimensionsComplexToRealFast(fft1, paddedDimensions, realSize);
	RandomAccessibleInterval<R> res = factory.create(realSize, type);
	
	final T typeT = Views.iterable(fft1).firstElement().createVariable();
	final S typeS = Views.iterable(fft2).firstElement().createVariable();
	RandomAccessibleInterval< T > fft1Copy;
	RandomAccessibleInterval< S > fft2Copy;

	try
	{
		fft1Copy = factory.imgFactory( typeT ).create(fft1, typeT );
		fft2Copy = factory.imgFactory( typeS ).create(fft2, typeS );
	}
	catch ( IncompatibleTypeException e )
	{
		throw new RuntimeException( "Cannot instantiate Img for type " + typeS.getClass().getSimpleName() + " or " + typeT.getClass().getSimpleName() );
	}
	
	
	calculatePCM(fft1, fft1Copy, fft2, fft2Copy, res, service);
	
	return res;
}
 

protected static Img< FloatType > loadInitialImage(
		final String fileName,
		final boolean checkNumbers,
		final float minValue,
		final Dimensions dimensions,
		final ImgFactory< FloatType > imageFactory )
{
	IOFunctions.println( "Loading image '" + fileName + "' as start for iteration." );

	final ImagePlus impPSI = LegacyStackImgLoaderIJ.open( new File( fileName ) );

	if ( impPSI == null )
	{
		IOFunctions.println( "Could not load image '" + fileName + "'." );
		return null;
	}

	final long[] dimPsi = impPSI.getStack().getSize() == 1 ? 
			new long[]{ impPSI.getWidth(), impPSI.getHeight() } : new long[]{ impPSI.getWidth(), impPSI.getHeight(), impPSI.getStack().getSize() };
	final Img< FloatType > psi = imageFactory.create( dimPsi, new FloatType() );
	LegacyStackImgLoaderIJ.imagePlus2ImgLib2Img( impPSI, psi, false );

	if ( psi == null )
	{
		IOFunctions.println( "Could not load image '" + fileName + "'." );
		return null;
	}
	else
	{
		boolean dimensionsMatch = true;

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

		for ( int d = 0; d < psi.numDimensions(); ++d )
		{
			if ( psi.dimension( d ) != dimensions.dimension( d ) )
				dimensionsMatch = false;

			dim[ d ] = dimensions.dimension( d );
		}

		if ( !dimensionsMatch )
		{
			IOFunctions.println(
					"Dimensions of '" + fileName + "' do not match: " +
					Util.printCoordinates( dimPsi ) + " != " + Util.printCoordinates( dim ) );
			return null;
		}

		if ( checkNumbers )
		{
			IOFunctions.println(
					"Checking values of '" + fileName + "' you can disable this check by setting " +
					"spim.process.fusion.deconvolution.MVDeconvolution.checkNumbers = false;" );

			boolean smaller = false;
			boolean hasZerosOrNeg = false;
			
			for ( final FloatType v : psi )
			{
				if ( v.get() < minValue )
					smaller = true;

				if ( v.get() <= 0 )
				{
					hasZerosOrNeg = true;
					v.set( minValue );
				}
			}

			if ( smaller )
				IOFunctions.println(
						"Some values '" + fileName + "' are smaller than the minimal value of " +
						minValue + ", this can lead to instabilities." );

			if ( hasZerosOrNeg )
				IOFunctions.println(
						"Some values '" + fileName + "' were smaller or equal to zero," +
						"they have been replaced with the min value of " + minValue );
		}
	}

	return psi;
}
 

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

public WeightNormalizer( final List< RandomAccessibleInterval< FloatType > > weights, final ImgFactory< FloatType > factory )
{
	this.weights = weights;
	this.sumWeights = factory.create( weights.get( 0 ), new FloatType() );
}
 

@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;
	}
}
 

public static < T extends RealType< T > & NativeType< T > > ImagePlus createReRegisteredSeries( final T targetType, final ImagePlus imp, final ArrayList<InvertibleBoundable> models, final int dimensionality )
{
	final int numImages = imp.getNFrames();

	// the size of the new image
	final int[] size = new int[ dimensionality ];
	// the offset relative to the output image which starts with its local coordinates (0,0,0)
	final double[] offset = new double[ dimensionality ];

	final int[][] imgSizes = new int[ numImages ][ dimensionality ];
	
	for ( int i = 0; i < numImages; ++i )
	{
		imgSizes[ i ][ 0 ] = imp.getWidth();
		imgSizes[ i ][ 1 ] = imp.getHeight();
		if ( dimensionality == 3 )
			imgSizes[ i ][ 2 ] = imp.getNSlices();
	}
	
	// estimate the boundaries of the output image and the offset for fusion (negative coordinates after transform have to be shifted to 0,0,0)
	Fusion.estimateBounds( offset, size, imgSizes, models, dimensionality );
			
	// for output
	final ImgFactory< T > f = new ImagePlusImgFactory< T >();
	// the composite
	final ImageStack stack = new ImageStack( size[ 0 ], size[ 1 ] );

	for ( int t = 1; t <= numImages; ++t )
	{
		for ( int c = 1; c <= imp.getNChannels(); ++c )
		{
			final Img<T> out = f.create( size, targetType );
			final Img< FloatType > in = ImageJFunctions.convertFloat( Hyperstack_rearranger.getImageChunk( imp, c, t ) );

			fuseChannel( out, Views.interpolate( Views.extendZero( in ), new NLinearInterpolatorFactory< FloatType >() ), offset, models.get( t - 1 ) );

			try
			{
				final ImagePlus outImp = ((ImagePlusImg<?,?>)out).getImagePlus();
				for ( int z = 1; z <= out.dimension( 2 ); ++z )
					stack.addSlice( imp.getTitle(), outImp.getStack().getProcessor( z ) );
			} 
			catch (ImgLibException e) 
			{
				Log.error( "Output image has no ImageJ type: " + e );
			}				
		}
	}
	
	//convertXYZCT ...
	ImagePlus result = new ImagePlus( "registered " + imp.getTitle(), stack );
	
	// numchannels, z-slices, timepoints (but right now the order is still XYZCT)
	if ( dimensionality == 3 )
	{
		result.setDimensions( size[ 2 ], imp.getNChannels(), imp.getNFrames() );
		result = OverlayFusion.switchZCinXYCZT( result );
		return CompositeImageFixer.makeComposite( result, CompositeImage.COMPOSITE );
	}
	//Log.info( "ch: " + imp.getNChannels() );
	//Log.info( "slices: " + imp.getNSlices() );
	//Log.info( "frames: " + imp.getNFrames() );
	result.setDimensions( imp.getNChannels(), 1, imp.getNFrames() );
	
	if ( imp.getNChannels() > 1 )
		return CompositeImageFixer.makeComposite( result, CompositeImage.COMPOSITE );
	return result;
}
 

/**
 * Create a new mask image without any specific content, but with
 * a defined size.
 */
public static RandomAccessibleInterval<BitType> createMask(long[] dim) {
	ImgFactory< BitType > imgFactory = new ArrayImgFactory< BitType >();
	return imgFactory.create(dim, new BitType());
}
 

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