Java Code Examples for net.imglib2.img.array.ArrayImgs#floats()

@Override
public RandomAccessibleInterval< FloatType > getFloatImage(
		int timepointId, boolean normalize, ImgLoaderHint... hints )
{
	final Random rnd = new Random( setupId );

	final Img< FloatType > img = ArrayImgs.floats( 512, 512, 86 );

	final float scale;
	
	if ( normalize )
		scale = 1;
	else
		scale = 20000;
	
	for ( final FloatType t : img )
		t.set( rnd.nextFloat() * scale);

	return img;
}
 

public static void main( String[] args )
{
	new ImageJ();
	
	Img< FloatType > img = ArrayImgs.floats( 500, 500 );
	BlendingRealRandomAccess blend = new BlendingRealRandomAccess(
			img,
			new float[]{ 100, 0 },
			new float[]{ 12, 150 } );
	
	Cursor< FloatType > c = img.localizingCursor();
	
	while ( c.hasNext() )
	{
		c.fwd();
		blend.setPosition( c );
		c.get().setReal( blend.get().getRealFloat() );
	}
	
	ImageJFunctions.show( img );
}
 

public static void main( String[] args )
{
	new ImageJ();
	
	Img< FloatType > img = ArrayImgs.floats( 500, 500 );
	BlendingRealRandomAccess blend = new BlendingRealRandomAccess(
			img,
			new float[]{ 100, 0 },
			new float[]{ 12, 150 } );
	
	Cursor< FloatType > c = img.localizingCursor();
	
	while ( c.hasNext() )
	{
		c.fwd();
		blend.setPosition( c );
		c.get().setReal( blend.get().getRealFloat() );
	}
	
	ImageJFunctions.show( img );
}
 

public static void main( String[] args )
{
	new ImageJ();
	
	Img< FloatType > img = ArrayImgs.floats( 500, 500 );
	BlendingRealRandomAccess blend = new BlendingRealRandomAccess(
			img,
			new float[]{ 100, 0 },
			new float[]{ 12, 150 } );
	
	Cursor< FloatType > c = img.localizingCursor();
	
	while ( c.hasNext() )
	{
		c.fwd();
		blend.setPosition( c );
		c.get().setReal( blend.get().getRealFloat() );
	}
	
	ImageJFunctions.show( img );
}
 

public static Img< FloatType > copyChannel( final ImagePlus imp, final int channel )
{
	final int w, h, d;

	Img< FloatType > img = ArrayImgs.floats( w = imp.getWidth(), h = imp.getHeight(), d = imp.getNSlices() );
	
	final Cursor< FloatType > c = img.cursor();

	for ( int z = 0; z < d; ++z )
	{
		final int[] pixels = (int[])imp.getStack().getProcessor( z + 1 ).getPixels();
		
		for ( int i = 0; i < w*h; ++i )
		{
			if ( channel == 0 )
				c.next().set( ( pixels[ i ] & 0xff0000) >> 16 );
			else if ( channel == 1 )
				c.next().set( ( pixels[ i ] & 0xff00 ) >> 8 );
			else
				c.next().set( pixels[ i ] & 0xff );
		}
	}
	
	return img;
}
 

public static void main(String[] args)
{
	final ImagePlus imp = IJ.openImage( "/Users/david/Desktop/stable HelaK-GFP-H2A.Z20000.tif" );
	new ImageJ();

	RandomAccessibleInterval< ? extends RealType > img = ImageJFunctions.wrapReal( imp );
	ArrayImg< FloatType, FloatArray > f = ArrayImgs.floats( 1024, 1024 );
	ArrayImg< FloatType, FloatArray > w = ArrayImgs.floats( 1024, 1024 );
	RandomAccessibleInterval< FloatType > interp = (RandomAccessibleInterval< FloatType >) getLinearInterpolation( img, new FloatType(), new float[] {0.5f,0.5f}, Executors.newSingleThreadExecutor() ).getA();
	
	RandomAccessibleInterval< FloatType > weight = new ArrayImgFactory( new FloatType() ).create( interp );
	applyWeights( interp, weight, new float[] {0.5f,0.5f}, new float[] {0,0}, new float[] {20,20}, false, Executors.newSingleThreadExecutor() );
	addTranslated( Views.iterable( interp ), f, new int[] {500, 700}, Executors.newSingleThreadExecutor() );
	addTranslated( Views.iterable( interp ), f, new int[] {400, 500}, Executors.newSingleThreadExecutor() );
	addTranslated( Views.iterable( weight ), w, new int[] {500, 700}, Executors.newSingleThreadExecutor() );
	addTranslated( Views.iterable( weight ), w, new int[] {400, 500}, Executors.newSingleThreadExecutor() );

	normalizeWeights( f, w, Executors.newSingleThreadExecutor() );
	
	ImageJFunctions.show( f );
}
 

private static Img<FloatType> generateFloatImg(final float[] values) {

		final float[] array =
			new float[(int) Intervals.numElements(new FinalInterval(dims))];

		if (array.length != values.length) {
			throw new RuntimeException("Number of values doesn't match dimmensions");
		}

		for (int i = 0; i < array.length; i++) {
			array[i] = values[i];
		}

		return ArrayImgs.floats(array, dims);
	}
 

/**
 * Function is called once with original images. Thereafter, each call is with
 * a shuffled version of the first image.
 * 
 * @param expectedPValue
 * @param expectedColocValue
 * @param result
 */
private void assertColoc(double expectedPValue, double expectedColocValue, double[] expectedColocValuesArray, double... result) {
	Img<FloatType> ch1 = ArrayImgs.floats(1); // NB: Images will be ignored.
	Img<FloatType> ch2 = ch1;

	// Mock the underlying op.
	final int[] count = { 0 };
	BinaryFunctionOp<Iterable<FloatType>, Iterable<FloatType>, Double> op = //
		op((input1, input2) -> {
			Double r;
			synchronized(this) {
			r = result[count[0]++];
			}
			return r;
		});
	
	PValueResult output = new PValueResult();
	output = ops.coloc().pValue(output, ch1, ch2, op, result.length - 1);
	Double actualPValue = output.getPValue();
	Double actualColocValue = output.getColocValue();
	double[] actualColocValuesArray = output.getColocValuesArray();
	assertEquals(expectedPValue, actualPValue, 0.0);
	assertEquals(expectedColocValue, actualColocValue, 0.0);
	Arrays.sort(actualColocValuesArray);
	for (int i = 0; i < expectedColocValuesArray.length; i++) {
		assertEquals(expectedColocValuesArray[i], actualColocValuesArray[i], 0.0);
	}
}
 

@Test
public void test() {
	long[] dims = { 15, 30 };
	// create input image
	Img<FloatType> input = ArrayImgs.floats(dims);
	MersenneTwisterFast random = new MersenneTwisterFast(SEED);
	for (FloatType b : input) {
		b.setReal(random.nextDouble());
	}

	// create 3 seeds
	Img<BitType> bits = ArrayImgs.bits(dims);
	RandomAccess<BitType> ra = bits.randomAccess();
	ra.setPosition(new int[] { 0, 0 });
	ra.get().set(true);
	ra.setPosition(new int[] { 4, 6 });
	ra.get().set(true);
	ra.setPosition(new int[] { 10, 20 });
	ra.get().set(true);

	// compute labeled seeds
	final ImgLabeling<Integer, IntType> labeledSeeds = ops.labeling().cca(bits, StructuringElement.EIGHT_CONNECTED);

	testWithoutMask(input, labeledSeeds);

	testWithMask(input, labeledSeeds);
}
 

final private static Img< FloatType > createGaussianKernel( final double[] sigmas )
{
	final int numDimensions = sigmas.length;

	final long[] imageSize = new long[ numDimensions ];
	final double[][] kernel = new double[ numDimensions ][];

	for ( int d = 0; d < numDimensions; ++d )
	{
		kernel[ d ] = Util.createGaussianKernel1DDouble( sigmas[ d ], true );
		imageSize[ d ] = kernel[ d ].length;
	}

	final Img< FloatType > kernelImg = ArrayImgs.floats( imageSize );

	final Cursor< FloatType > cursor = kernelImg.localizingCursor();
	final int[] position = new int[ numDimensions ];

	while ( cursor.hasNext() )
	{
		cursor.fwd();
		cursor.localize( position );

		double value = 1;

		for ( int d = 0; d < numDimensions; ++d )
			value *= kernel[ d ][ position[ d ] ];

		cursor.get().set( ( float ) value );
	}

	return kernelImg;
}
 

@Test
public void regressionTest() throws Exception {

	// load in input image and expected output image.
	Img<FloatType> inputImg = (Img<FloatType>) ops.run(
		net.imagej.ops.image.equation.DefaultEquation.class,
		"Math.tan(0.3*p[0]) + Math.tan(0.1*p[1])");
	Img<FloatType> expectedOutput = ((Img<FloatType>) openFloatImg(
		"Result.tif"));

	// create ouput image
	long[] dims = new long[inputImg.numDimensions()];
	inputImg.dimensions(dims);
	Img<FloatType> actualOutput = ArrayImgs.floats(dims);

	// scale over which the filter operates (sensitivity)
	int scale = 1;

	// physical spacing between data points (1,1 since I got it from the
	// computer)
	double[] spacing = { 1, 1 };

	// run the op
	ops.run(net.imagej.ops.filter.vesselness.DefaultFrangi.class, actualOutput,
		inputImg, spacing, scale);

	// compare the output image data to that stored in the file.
	Cursor<FloatType> cursor = Views.iterable(actualOutput).localizingCursor();
	RandomAccess<FloatType> actualRA = actualOutput.randomAccess();
	RandomAccess<FloatType> expectedRA = expectedOutput.randomAccess();

	while (cursor.hasNext()) {
		cursor.fwd();
		actualRA.setPosition(cursor);
		expectedRA.setPosition(cursor);
		assertEquals(expectedRA.get().get(), actualRA.get().get(), 0);
	}
}
 

@Override
public RandomAccessibleInterval< FloatType > getFloatImage( final ViewId view, final boolean normalize )
{
	try
	{
		final MultipageTiffReader r = new MultipageTiffReader( mmFile );

		final ArrayImg< FloatType, ? > img = ArrayImgs.floats( r.width(), r.height(), r.depth() );
		final BasicViewDescription< ? > vd = sequenceDescription.getViewDescriptions().get( view );

		populateImage( img, vd, r );

		if ( normalize )
			normalize( img );

		updateMetaDataCache( view, r.width(), r.height(), r.depth(), r.calX(), r.calY(), r.calZ() );

		r.close();

		return img;
	}
	catch ( Exception e )
	{
		IOFunctions.printlnSafe( "Failed to load viewsetup=" + view.getViewSetupId() + " timepoint=" + view.getTimePointId() + ": " + e );
		e.printStackTrace();
		return null;
	}
}
 

@Override
public RandomAccessibleInterval< FloatType > getFloatImage( final ViewId view, final boolean normalize )
{
	final BasicViewDescription< ? > vd = sd.getViewDescriptions().get( view );
	final Dimensions d = vd.getViewSetup().getSize();
	final VoxelDimensions dv = vd.getViewSetup().getVoxelSize();

	final ArrayImg< FloatType, ? > img = ArrayImgs.floats( d.dimension( 0 ), d.dimension( 1 ), d.dimension(  2 ) );

	final String ampOrPhaseDir;

	if ( vd.getViewSetup().getAttribute( Channel.class ).getId() == ampChannelId )
		ampOrPhaseDir = amplitudeDir;
	else if ( vd.getViewSetup().getAttribute( Channel.class ).getId() ==  phaseChannelId )
		ampOrPhaseDir = phaseDir;
	else
		throw new RuntimeException( "viewSetupId=" + view.getViewSetupId() + " is not Amplitude nor phase." );

	populateImage( img, directory, stackDir, ampOrPhaseDir, zPlanes, timepoints.get( view.getTimePointId() ), extension );

	if ( normalize )
		normalize( img );

	updateMetaDataCache( view, (int)d.dimension( 0 ), (int)d.dimension( 1 ), (int)d.dimension( 2 ), dv.dimension( 0 ), dv.dimension( 1 ), dv.dimension( 2 ) );

	return img;
}
 

public static void main( String[] args )
{
	// define the blocksize so that it is one single block
	final RandomAccessibleInterval< FloatType > block = ArrayImgs.floats( 384, 384 );
	final long[] blockSize = new long[ block.numDimensions() ];
	block.dimensions( blockSize );

	final RandomAccessibleInterval< FloatType > image = ArrayImgs.floats( 1024, 1024 );
	final long[] imgSize = new long[ image.numDimensions() ];
	image.dimensions( imgSize );

	// whatever the kernel size is (extra size/2 in general)
	final long[] kernelSize = new long[]{ 16, 32 };

	final BlockGeneratorFixedSizePrecise blockGenerator = new BlockGeneratorFixedSizePrecise( blockSize );
	final Block[] blocks = blockGenerator.divideIntoBlocks( imgSize, kernelSize );

	int i = 0;

	for ( final Block b : blocks )
	{
		// copy data from the image to the block (including extra space for outofbounds/real image data depending on kernel size)
		b.copyBlock( Views.extendMirrorDouble( image ), block );

		// do something with the block (e.g. also multithreaded, cluster, ...)
		for ( final FloatType f : Views.iterable( block ) )
			f.set( i );

		++i;

		// write the block back (use a temporary image if multithreaded or in general not all are copied first)
		b.pasteBlock( image, block );
	}

	ImageJFunctions.show( image );
}
 

@Override
public boolean process()
{
	// do not operate at the edge, 80% of the memory is a good idea I think
	final long memAvail = Math.round( cudaDevice.getFreeDeviceMemory() * ( percentGPUMem / 100.0 ) );
	final long imgBytes = numPixels() * 4 * 2; // float, two images on the card at once

	final long[] numBlocksDim = net.imglib2.util.Util.int2long( computeNumBlocksDim( memAvail, imgBytes, percentGPUMem, img.numDimensions(), "CUDA-Device " + cudaDevice.getDeviceId() ) );
	final BlockGenerator< Block > generator;

	if ( accurate )
		generator = new BlockGeneratorVariableSizePrecise( numBlocksDim );
	else
		generator = new BlockGeneratorVariableSizeSimple( numBlocksDim );

	final Block[] blocks = generator.divideIntoBlocks( getImgSize( img ), getKernelSize( sigma ) );

	if ( !accurate && blocks.length == 1 && ArrayImg.class.isInstance( img ) )
	{
		IOFunctions.println( "Conovlving image as one single block." );
		long time = System.currentTimeMillis();

		// copy the only directly into the result
		blocks[ 0 ].copyBlock( img, result );
		long copy = System.currentTimeMillis();
		IOFunctions.println( "Copying data took " + ( copy - time ) + "ms" );

		// convolve
		final float[] resultF = ((FloatArray)((ArrayImg< net.imglib2.type.numeric.real.FloatType, ? > )result).update( null ) ).getCurrentStorageArray();
		cudaconvolve.gauss( resultF, getImgSizeInt( result ), sigma, OutOfBounds.EXTEND_BORDER_PIXELS, 0 );
		IOFunctions.println( "Convolution took " + ( System.currentTimeMillis() - copy ) + "ms using device=" + cudaDevice.getDeviceName() + " (id=" + cudaDevice.getDeviceId() + ")" );

		// no copy back required
	}
	else
	{
		final RandomAccessible< net.imglib2.type.numeric.real.FloatType > input;
		
		if ( accurate )
			input = Views.extendMirrorSingle( img );
		else
			input = img;
		
		for( final Block block : blocks )
		{
			//long time = System.currentTimeMillis();
			final ArrayImg< net.imglib2.type.numeric.real.FloatType, FloatArray > imgBlock = ArrayImgs.floats( block.getBlockSize() );

			// copy the block
			block.copyBlock( input, imgBlock );
			//long copy = System.currentTimeMillis();
			//IOFunctions.println( "Copying block took " + ( copy - time ) + "ms" );

			// convolve
			final float[] imgBlockF = ((FloatArray)((ArrayImg< net.imglib2.type.numeric.real.FloatType, ? > )imgBlock).update( null ) ).getCurrentStorageArray();
			cudaconvolve.gauss( imgBlockF, getImgSizeInt( imgBlock ), sigma, OutOfBounds.EXTEND_BORDER_PIXELS, 0 );
			//long convolve = System.currentTimeMillis();
			//IOFunctions.println( "Convolution took " + ( convolve - copy ) + "ms using device=" + cudaDevice.getDeviceName() + " (id=" + cudaDevice.getDeviceId() + ")" );

			// no copy back required
			block.pasteBlock( result, imgBlock );
			//IOFunctions.println( "Pasting block took " + ( System.currentTimeMillis() - convolve ) + "ms" );
		}
	}

	return true;
}
 

@Test
public void testRandomPeaksMT()
{
	Img< FloatType > img = ArrayImgs.floats( 50, 50 );
	Random rnd = new Random( seed );
	
	for( FloatType t : img )
		t.set( rnd.nextFloat() );
	
	ArrayList< Pair< Localizable, Double > > correct = new ArrayList< Pair<Localizable,Double> >(); 
	
	RandomAccess<FloatType> ra = img.randomAccess();
	
	for ( int i = 0; i < 10; ++i ){
		for (int d = 0; d< img.numDimensions(); d++){
			ra.setPosition((int) (rnd.nextDouble() * 50), d);
		}
		ra.get().set((float) (rnd.nextDouble() + 1.0));
		
		correct.add(new ValuePair<Localizable, Double>(new Point(ra), new Double(ra.get().get())));
	}
	
	// sort the peaks in descending order
	Collections.sort(correct, new Comparator<Pair< Localizable, Double >>() {
		@Override
		public int compare(Pair<Localizable, Double> o1, Pair<Localizable, Double> o2) {
			return Double.compare(o2.getB(), o1.getB());
		}
	});
	
	int nMax = 5;
	ArrayList< Pair< Localizable, Double > > found = FourNeighborhoodExtrema.findMaxMT(Views.extendPeriodic(img), img, nMax, Executors.newFixedThreadPool(Runtime.getRuntime().availableProcessors()));
	assertEquals(nMax, found.size());
	
	
	long[] posCorrect = new long[img.numDimensions()];
	long[] posFound = new long[img.numDimensions()];
	
	for (int i = 0; i<found.size(); i++){			
		assertEquals(correct.get(i).getB(), found.get(i).getB());
		
		correct.get(i).getA().localize(posCorrect);
		found.get(i).getA().localize(posFound);			
		assertArrayEquals(posCorrect, posFound);
	}
}
 

@SuppressWarnings("unchecked")
@Override
protected boolean mapIntensities(final Patch p, final ImagePlus imp) {
	
	final ImagePlus coefficients = new Opener().openImage(
		getUNUIdFolder() +
		"trakem2.its/" +
		createIdPath(Long.toString(p.getId()), "it", ".tif"));

	if (coefficients == null)
		return false;
	
	final ImageProcessor ip = imp.getProcessor();
	
	@SuppressWarnings({"rawtypes"})
	final LinearIntensityMap<FloatType> map =
			new LinearIntensityMap<FloatType>(
					(FloatImagePlus)ImagePlusImgs.from(coefficients));

	@SuppressWarnings("rawtypes")
	Img img;

	final long[] dims = new long[]{imp.getWidth(), imp.getHeight()};
	switch (p.getType()) {
	case ImagePlus.GRAY8:
	case ImagePlus.COLOR_256:		// this only works for continuous color tables
		img = ArrayImgs.unsignedBytes((byte[])ip.getPixels(), dims);
		break;
	case ImagePlus.GRAY16:
		img = ArrayImgs.unsignedShorts((short[])ip.getPixels(), dims);
		break;
	case ImagePlus.COLOR_RGB:
		img = ArrayImgs.argbs((int[])ip.getPixels(), dims);
		break;
	case ImagePlus.GRAY32:
		img = ArrayImgs.floats((float[])ip.getPixels(), dims);
		break;
	default:
		img = null;
	}
	
	if (img == null)
		return false;

	map.run(img);
	
	return true;
}
 

public static void main(String[] args) {
	
	new ImageJ();
	Img<FloatType> img1 = ImgLib2Util.openAs32Bit(new File("src/main/resources/img1.tif"));
	long[] dims = new long[img1.numDimensions()];

	BlendedExtendedMirroredRandomAccesible2<FloatType> ext = new BlendedExtendedMirroredRandomAccesible2<FloatType>(img1, new int[]{100, 100, 10});
	
	ext.getExtInterval().dimensions(dims);		
	Img<FloatType> img2 = ArrayImgs.floats(dims);

	// TODO: For efficiency reasons we should now also iterate the BlendedExtendedMirroredRandomAccesible and not the image
	long start = System.currentTimeMillis();
	
	Cursor<FloatType> c = img2.cursor();
	for (FloatType e : Views.iterable(Views.interval(ext, ext.getExtInterval())))
	{
		c.fwd();
		c.get().set(e);
	}
	
	long end = System.currentTimeMillis();		
	System.out.println(end-start);
	
	ImageJFunctions.show(img2);	
	
	//RandomAccessibleInterval<FloatType> img3 = Views.interval(ext, ext.getExtInterval());		
	//ImageJFunctions.show(img3);
	

	

}
 

@Test
public void testPC() {
	
	// TODO: very large shifts (nearly no overlap) lead to incorrect shift determination (as expected)
	// maybe we can optimize behaviour in this situation
	Img< FloatType > img = ArrayImgs.floats( 200, 200 );
	Random rnd = new Random( seed );
	
	for( FloatType t : img )
		t.set( rnd.nextFloat());
	
	long shiftX = 28;
	long shiftY = 0;
	
	FinalInterval interval1 = new FinalInterval(new long[] {50, 50});
	FinalInterval interval2 = Intervals.translate(interval1, shiftX, 0);
	interval2 = Intervals.translate(interval2, shiftY, 1);

	
	int [] extension = new int[img.numDimensions()];
	Arrays.fill(extension, 10);
	
	RandomAccessibleInterval<FloatType> pcm = PhaseCorrelation2.calculatePCM(Views.interval(img, interval1), Views.zeroMin(Views.interval(img, interval2)), extension, new ArrayImgFactory<FloatType>(), 
			new FloatType(), new ArrayImgFactory<ComplexFloatType>(), new ComplexFloatType(), Executors.newFixedThreadPool(Runtime.getRuntime().availableProcessors()));
	
	PhaseCorrelationPeak2 shiftPeak = PhaseCorrelation2.getShift(pcm, Views.interval(img, interval1), Views.zeroMin(Views.interval(img, interval2)), 20, 0, false, Executors.newFixedThreadPool(Runtime.getRuntime().availableProcessors()));
	
	long[] expected = new long[]{shiftX, shiftY};
	long[] found = new long[img.numDimensions()];
	
	
	
	shiftPeak.getShift().localize(found);
	
	assertArrayEquals(expected, found);
	
}
 

@Test
public void testPCRealShift() {
	
	// TODO: very large shifts (nearly no overlap) lead to incorrect shift determination (as expected)
	// maybe we can optimize behaviour in this situation
	Img< FloatType > img = ArrayImgs.floats( 200, 200 );
	Random rnd = new Random( seed );
	
	for( FloatType t : img )
		t.set( rnd.nextFloat());
	
	double shiftX = -20.9;
	double shiftY = 1.9;
	
	// to test < 0.5 px off
	final double eps = 0.5;
	
	FinalInterval interval2 = new FinalInterval(new long[] {50, 50});
	
	

	AffineRandomAccessible< FloatType, AffineGet > imgTr = RealViews.affine( Views.interpolate( Views.extendZero( img ), new NLinearInterpolatorFactory<>() ), new Translation2D( shiftX, shiftY ));
	IntervalView< FloatType > img2 = Views.interval( Views.raster( imgTr ), interval2);
	
	int [] extension = new int[img.numDimensions()];
	Arrays.fill(extension, 10);
	
	RandomAccessibleInterval<FloatType> pcm = PhaseCorrelation2.calculatePCM(Views.zeroMin(img2), Views.zeroMin(Views.interval(img, interval2)), extension, new ArrayImgFactory<FloatType>(), 
			new FloatType(), new ArrayImgFactory<ComplexFloatType>(), new ComplexFloatType(), Executors.newFixedThreadPool(Runtime.getRuntime().availableProcessors()));
	
	PhaseCorrelationPeak2 shiftPeak = PhaseCorrelation2.getShift(pcm, Views.zeroMin(img2), Views.zeroMin(Views.interval(img, interval2)), 20, 0, true, Executors.newFixedThreadPool(Runtime.getRuntime().availableProcessors()));
	
	
	double[] expected = new double[]{shiftX, shiftY};
	double[] found = new double[img.numDimensions()];
	
	
	
	
	shiftPeak.getSubpixelShift().localize(found);
	
	System.out.println( Util.printCoordinates( found ) );
	
	
	for (int d = 0; d < expected.length; d++)
		assertTrue( Math.abs( expected[d] - found[d] ) < eps );
	
}