net.imglib2.type.numeric.RealType Java Examples

/**
 * Update a volume with the given IterableInterval.
 * This method actually populates the volume
 * @param image
 * @param name
 * @param voxelDimensions
 * @param v
 * @param <T>
 * @return a Node corresponding to the input volume
 */
public <T extends RealType<T>> Node updateVolume( IterableInterval<T> image, String name,
                                                                   float[] voxelDimensions, Volume v ) {
    List<SourceAndConverter<T>> sacs = (List<SourceAndConverter<T>>) v.getMetadata().get("sources");

    RandomAccessibleInterval<T> source = sacs.get(0).getSpimSource().getSource(0, 0);// hard coded to timepoint and mipmap 0

    Cursor<T> sCur = Views.iterable(source).cursor();
    Cursor<T> iCur = image.cursor();
    while( sCur.hasNext() ) {
        sCur.fwd();
        iCur.fwd();
        sCur.get().set(iCur.get());
    }

    v.getVolumeManager().notifyUpdate(v);
    v.getVolumeManager().requestRepaint();
    //v.getCacheControls().clear();
    //v.setDirty( true );
    v.setNeedsUpdate( true );
    //v.setNeedsUpdateWorld( true );

    return v;
}
 

public <T extends RealType<T>> RandomAccessibleInterval< T > aggregate(
		List<RandomAccessibleInterval< T >> rais,
		List<? extends ViewId> vids,
		AbstractSequenceDescription< ?, ? extends BasicViewDescription< ? >, ? > sd
		)
{
	Map<BasicViewDescription< ? >, RandomAccessibleInterval<T>> map = new HashMap<>();

	for (int i = 0; i < vids.size(); i++)
	{
		ViewId vid = vids.get( i );
		BasicViewDescription< ? > vd = sd.getViewDescriptions().get( vid );
		map.put( vd, rais.get( i ) );
	}

	for (final Action action : actions)
	{
		map = action.aggregate( map );
	}

	// return the first RAI still present
	// ideally, there should be only one left
	return map.values().iterator().next();
	
}
 

public <T extends RealType<T> & NativeType<T>, V extends AbstractVolatileRealType<T, V> & NativeType<V>>
SourceState<T, V> getRaw(
		final String name,
		final SharedQueue queue,
		final int priority) throws Exception
{
	LOG.debug("Raw data set requested. Name={}", name);
	final N5Writer           writer     = n5.get();
	final String             dataset    = this.dataset.get();
	final double[]           resolution = asPrimitiveArray(resolution());
	final double[]           offset     = asPrimitiveArray(offset());
	final N5BackendRaw<T, V> backend    = new N5BackendRaw<>(writer, dataset);
	final SourceState<T, V>  state      = new ConnectomicsRawState<>(backend, queue, priority, name, resolution, offset);
	LOG.debug("Returning raw source state {} {}", name, state);
	return state;
}
 

private static <T extends RealType<T>> RandomAccessible<RealComposite<T>>  collapseDimension(
		final RandomAccessibleInterval<T> rai,
		final int dimension,
		final long[] channels,
		final T extension
)
{
	final int lastDim = rai.numDimensions() - 1;
	final int numChannels = (int) rai.dimension(dimension);

	long[] min = Intervals.minAsLongArray(rai);
	long[] max = Intervals.maxAsLongArray(rai);

	assert LongStream.of(channels).filter(c -> c > max[dimension] && c < min[dimension]).count() == 0;

	final RandomAccessibleInterval<T> relevantRai = isFullRange(channels, numChannels)
			? rai
			: Views.stack(LongStream.of(channels).mapToObj(channel -> Views.hyperSlice(rai, dimension, channel)).collect(Collectors.toList()));

	final RandomAccessible<T> ra = Views.extendValue(lastDim == dimension
			? relevantRai
			: Views.moveAxis(relevantRai, dimension, lastDim), extension);
	return Views.collapseReal(ra, numChannels);
}
 

/**
 *
 * @param reader container
 * @param dataset dataset
 * @param transform transforms voxel data into real world coordinates
 * @param priority in fetching queue
 * @param name initialize with this name
 * @param <T> data type
 * @param <V> viewer type
 * @return {@link DataSource}
 * @throws IOException if any N5 operation throws {@link IOException}
 */
public static <T extends NativeType<T> & RealType<T>, V extends Volatile<T> & NativeType<V> & RealType<V>>
DataSource<T, V> openRawAsSource(
		final N5Reader reader,
		final String dataset,
		final AffineTransform3D transform,
		final SharedQueue queue,
		final int priority,
		final String name) throws IOException, ReflectionException {
	return openScalarAsSource(
			reader,
			dataset,
			transform,
			queue,
			priority,
			i -> i == Interpolation.NLINEAR
			     ? new NLinearInterpolatorFactory<>()
			     : new NearestNeighborInterpolatorFactory<>(),
			i -> i == Interpolation.NLINEAR
			     ? new NLinearInterpolatorFactory<>()
			     : new NearestNeighborInterpolatorFactory<>(),
			name
	                         );
}
 

private static <D extends NativeType<D> & RealType<D>, T extends RealType<T>> RealComposite<T>  createExtension(
		final D d,
		final T t,
		final Converter<D, T> converter,
		final long size,
		IntFunction<D> valueAtIndex
)
{
	LOG.debug("Creating extension with size {}", size);
	final ArrayImg<D, ?> img = new ArrayImgFactory<>(d).create(1, size);
	img.setLinkedType((D) d.getNativeTypeFactory().createLinkedType((NativeImg)img));
	final CompositeIntervalView<D, RealComposite<D>> collapsed = Views.collapseReal(img);
	RealComposite<D> extensionCopy = collapsed.randomAccess().get();
	for (int channel = 0; channel < size; ++channel)
		extensionCopy.get(channel).set(valueAtIndex.apply(channel));
	return Views.collapseReal(Converters.convert((RandomAccessibleInterval<D>)img, converter, t.createVariable())).randomAccess().get();
}
 

/**
 * Compute the partial derivative of source in a particular dimension.
 *
 * @param source
 *            source image, has to provide valid data in the interval of the
 *            gradient image plus a one pixel border in dimension.
 * @param target
 *            output image, the partial derivative of source in the
 *            specified dimension.
 * @param dimension
 *            along which dimension the partial derivatives are computed
 * @param <T> pixel type source
 * @param <S> pixel type target
 */
public static < T extends RealType< T >, S extends RealType< S > > void gradient(
		final RandomAccessible< T > source,
		final RandomAccessibleInterval< S > target,
		final int dimension )
{
	final Cursor< T > front = Views.flatIterable(
			Views.interval( source,
					Intervals.translate( target, 1, dimension ) ) ).cursor();
	final Cursor< T > back = Views.flatIterable(
			Views.interval( source,
					Intervals.translate( target, -1, dimension ) ) ).cursor();
	for( final S t : Views.flatIterable( target ) )
	{
		t.setReal( front.next().getRealDouble() - back.next().getRealDouble());
		t.mul( 0.5 );
	}
}
 

/**
 * calculate the shift between two images from the phase correlation matrix
 * @param pcm the phase correlation matrix of img1 and img2
 * @param img1 source image 1
 * @param img2 source image 2
 * @param nHighestPeaks the number of peaks in pcm to check via cross. corr.
 * @param minOverlap minimal overlap (in pixels)
 * @param subpixelAccuracy whether to do subpixel shift peak localization or not
 * @param interpolateSubpixel whether to interpolate the subpixel shift in cross. corr.
 * @param service thread pool
 * @param <R> PCM pixel type
 * @param <T> image 1 pixel type
 * @param <S> image 2 pixel type
 * @return best (highest c.c.) shift peak
 */
public static <T extends RealType<T>, S extends RealType<S>, R extends RealType<R>> PhaseCorrelationPeak2 getShift(
		RandomAccessibleInterval<R> pcm, RandomAccessibleInterval<T> img1, RandomAccessibleInterval<S> img2, int nHighestPeaks,
		long minOverlap, boolean subpixelAccuracy, boolean interpolateSubpixel, ExecutorService service)
{
	System.out.println( "PCM" );
	List<PhaseCorrelationPeak2> peaks = PhaseCorrelation2Util.getPCMMaxima(pcm, service, nHighestPeaks, subpixelAccuracy);
	//peaks = PhaseCorrelation2Util.getHighestPCMMaxima(peaks, nHighestPeaks);
	System.out.println( "expand" );
	PhaseCorrelation2Util.expandPeakListToPossibleShifts(peaks, pcm, img1, img2);
	System.out.print( "cross " );
	long t = System.currentTimeMillis();
	PhaseCorrelation2Util.calculateCrossCorrParallel(peaks, img1, img2, minOverlap, service, interpolateSubpixel);
	System.out.println( (System.currentTimeMillis() - t) );
	System.out.println( "sort" );
	Collections.sort(peaks, Collections.reverseOrder(new PhaseCorrelationPeak2.ComparatorByCrossCorrelation()));
	System.out.println( "done" );

	if (peaks.size() > 0)
		return peaks.get(0);
	else
		return null;
}
 

public static <
		D extends RealType<D> & NativeType<D>,
		T extends AbstractVolatileRealType<D, T> & NativeType<T>> N5ChannelDataSource<D, T> zeroExtended(
		final N5Meta meta,
		final AffineTransform3D transform,
		final String name,
		final SharedQueue queue,
		final int priority,
		final int channelDimension,
		final long[] channels) throws IOException, DataTypeNotSupported {

	return extended(
			meta,
			transform,
			name,
			queue,
			priority,
			channelDimension,
			channels,
			RealType::setZero,
			RealType::setZero);
}
 

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

public <T extends RealType<T>> T getMean(IterableInterval< T > img)
{
	RealSum sum = new RealSum();
	long nPix = 0;
	
	for (T t : img)
	{
		sum.add( t.getRealDouble() );
		nPix++;
	}
	
	T res = img.firstElement().createVariable();
	res.setReal( sum.getSum()/nPix );
	return res;
	
}
 

private static <R extends RealType<R> & NativeType<R>, B extends BooleanType<B>> void computeSignedDistanceTransform(
		final RandomAccessibleInterval<B> mask,
		final RandomAccessibleInterval<R> target,
		final DISTANCE_TYPE distanceType,
		final double... weights)
{
	final RandomAccessibleInterval<R> distanceOutside = target;
	final RandomAccessibleInterval<R> distanceInside = new ArrayImgFactory<>(Util.getTypeFromInterval(target)).create(target);
	DistanceTransform.binaryTransform(mask, distanceOutside, distanceType, weights);
	DistanceTransform.binaryTransform(Logical.complement(mask), distanceInside, distanceType, weights);
	LoopBuilder.setImages(distanceOutside, distanceInside, target).forEachPixel((outside, inside, result) -> {
		switch (distanceType)
		{
		case EUCLIDIAN:
			result.setReal(Math.sqrt(outside.getRealDouble()) - Math.sqrt(inside.getRealDouble()));
			break;
		case L1:
			result.setReal(outside.getRealDouble() - inside.getRealDouble());
			break;
		}
	});
}
 

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 <
		D extends RealType<D> & NativeType<D>,
		T extends AbstractVolatileRealType<D, T> & NativeType<T>> N5ChannelDataSource<D, T> valueExtended(
		final N5Meta meta,
		final AffineTransform3D transform,
		final String name,
		final SharedQueue queue,
		final int priority,
		final int channelDimension,
		final long[] channels,
		final double extension) throws IOException, DataTypeNotSupported {

	return extended(
			meta, transform,
			name,
			queue,
			priority,
			channelDimension,
			channels,
			d -> d.setReal(extension),
			t -> t.setReal(extension)
	);
}
 

public static <T extends RealType<T>, S extends RealType<S>> void copyRealImage(IterableInterval<T> source, RandomAccessibleInterval<S> dest) {
	RandomAccess<S> destRA = dest.randomAccess();
	Cursor<T> srcC = source.cursor();
	
	
	while (srcC.hasNext()){
		srcC.fwd();
		destRA.setPosition(srcC);
		destRA.get().setReal(srcC.get().getRealDouble());
	}
}
 

public static <
		D extends RealType<D> & NativeType<D>,
		T extends AbstractVolatileRealType<D, T> & NativeType<T>> N5ChannelDataSource<D, T> valueExtended(
		final N5Meta meta,
		final AffineTransform3D transform,
		final String name,
		final SharedQueue queue,
		final int priority,
		final int channelDimension,
		final long channelMin,
		final long channelMax,
		final boolean revertChannelOrder,
		final double value) throws IOException, DataTypeNotSupported {

	return extended(
			meta,
			transform,
			name,
			queue,
			priority,
			channelDimension,
			channelMin,
			channelMax,
			revertChannelOrder,
			d -> d.setReal(value),
			t -> t.setReal(value)
	);
}
 

private static <D extends NativeType<D> & RealType<D>, T extends RealType<D>> RealComposite<D>  createExtension(
		final D d,
		final long size
)
{
	return createExtension(d, d.createVariable(), new TypeIdentity<>(), size);
}
 

private static <T extends RealType<T>> Function<Interpolation, InterpolatorFactory<T, RandomAccessible<T>>>
realTypeInterpolation()
{
	return i -> i.equals(Interpolation.NLINEAR)
	            ? new NLinearInterpolatorFactory<>()
	            : new NearestNeighborInterpolatorFactory<>();
}
 

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 ComplexType<T>, S extends ComplexType<S>, R extends RealType<R>> void calculatePCM(
		RandomAccessibleInterval<T> fft1, RandomAccessibleInterval<T> fft1Copy, RandomAccessibleInterval<S> fft2, RandomAccessibleInterval<S> fft2Copy, RandomAccessibleInterval<R> pcm,
		ExecutorService service)
{
	// TODO: multithreaded & check for cursor vs randomaccess
	
	// normalize, save to copies
	PhaseCorrelation2Util.normalizeInterval(fft1, fft1Copy, service);
	PhaseCorrelation2Util.normalizeInterval(fft2, fft2Copy, service);
	// conjugate
	PhaseCorrelation2Util.complexConjInterval(fft2Copy, fft2Copy, service);
	// in-place multiplication
	PhaseCorrelation2Util.multiplyComplexIntervals(fft1Copy, fft2Copy, fft1Copy, service);
	FFT.complexToReal(fft1Copy, pcm, service);
}
 

/**
 * Adds a SourceAndConverter to the scene.
 *
 * @param sources The list of SourceAndConverter to add
 * @param name Name of the dataset
 * @param voxelDimensions Array with voxel dimensions.
 * @param <T> Type of the dataset.
 * @return THe node corresponding to the volume just added.
 */
public <T extends RealType<T>> Node addVolume(List<SourceAndConverter<T>> sources,
                                              int numTimepoints,
                                              String name,
                                              float... voxelDimensions ) {
    int setupId = 0;
    ArrayList<ConverterSetup> converterSetups = new ArrayList<>();
    for( SourceAndConverter source: sources ) {
        converterSetups.add(BigDataViewer.createConverterSetup(source, setupId++));
    }

    return addVolume(sources, converterSetups, numTimepoints, name, voxelDimensions);
}
 

public <T extends RealType<T> & NativeType<T>, V extends AbstractVolatileRealType<T, V> & NativeType<V>>
List<? extends SourceState<RealComposite<T>, VolatileWithSet<RealComposite<V>>>> getChannels(
		final String name,
		final int[] channelSelection,
		final SharedQueue queue,
		final int priority) throws Exception
{
	final N5Reader                  reader            = n5.get();
	final String                    dataset           = this.dataset.get();
	final N5Meta                    meta              = N5Meta.fromReader(reader, dataset);
	final double[]                  resolution        = asPrimitiveArray(resolution());
	final double[]                  offset            = asPrimitiveArray(offset());
	final AffineTransform3D         transform         = N5Helpers.fromResolutionAndOffset(resolution, offset);
	final long                      numChannels       = datasetAttributes.get().getDimensions()[3];

	LOG.debug("Got channel info: num channels={} channels selection={}", numChannels, channelSelection);
	final N5BackendChannel<T, V> backend = new N5BackendChannel<>(n5.get(), dataset, channelSelection, 3);
	final ConnectomicsChannelState<T, V, RealComposite<T>, RealComposite<V>, VolatileWithSet<RealComposite<V>>> state = new ConnectomicsChannelState<>(
			backend,
			queue,
			priority,
			name + "-" + Arrays.toString(channelSelection),
			resolution,
			offset);
	state.converter().setMins(i -> min().get());
	state.converter().setMaxs(i -> max().get());
	return Collections.singletonList(state);
}
 

private static <
		R extends RealType<R>,
		C extends RealComposite<R>,
		V extends Volatile<C>> void convertInverting(
		V input,
		final ARGBType output,
		final int numChannels,
		final DoubleProperty[] min,
		final double[] scaleR,
		final double[] scaleG,
		final double[] scaleB,
		final int A
) {
	double rd = 0.0;
	double gd = 0.0;
	double bd = 0.0;
	final RealComposite<? extends RealType<?>> c = input.get();
	for (int channel = 0; channel < numChannels; ++channel) {
		final double v = c.get(channel).getRealDouble() - min[channel].get();
		rd += scaleR[channel] * v;
		gd += scaleG[channel] * v;
		bd += scaleB[channel] * v;
	}
	final int r0 = (int) (rd + 0.5);
	final int g0 = (int) (gd + 0.5);
	final int b0 = (int) (bd + 0.5);
	final int r = Math.min(255, Math.max(r0, 0));
	final int g = Math.min(255, Math.max(g0, 0));
	final int b = Math.min(255, Math.max(b0, 0));
	output.set(ARGBType.rgba(r, g, b, A));
}
 

private static <D extends RealType<D> & NativeType<D>, T extends AbstractVolatileRealType<D, T> & NativeType<T>> SourceState<D, T> makeRawSourceState(
		final PainteraBaseView viewer,
		final N5Reader container,
		final String group,
		final double[] resolution,
		final double[] offset,
		final double min,
		final double max,
		final String name
) throws IOException {
	try {
		final N5BackendRaw<D, T> backend = new N5BackendRaw<>(
				N5Meta.fromReader(container, group).getWriter(),
				group);
		final ConnectomicsRawState<D, T> state =  new ConnectomicsRawState<>(
				backend,
				viewer.getQueue(),
				0,
				name,
				resolution,
				offset);
		state.converter().setMin(min);
		state.converter().setMax(max);
		return state;
	} catch (final ReflectionException e) {
		throw new IOException(e);
	}
}
 

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, ImgFactory<R> factory, R type,
		ImgFactory<C> fftFactory, C fftType, ExecutorService service) {
	
	int [] extension = new int[img1.numDimensions()];
	Arrays.fill(extension, 10);
	return calculatePCM(img1, img2, extension, factory, type, fftFactory, fftType, service);
}
 

private static <D extends RealType<D> & NativeType<D>, T extends AbstractVolatileRealType<D, T> & NativeType<T>> ChannelSourceState<D, T, ?, ?> makeChannelSourceState(
		final PainteraBaseView viewer,
		final N5Reader reader,
		final String dataset,
		final AffineTransform3D transform,
		final int channelDimension,
		final long[] channels,
		final double min,
		final double max,
		final String name
) throws IOException {
	try {
		final N5ChannelDataSource<D, T> channelSource = N5ChannelDataSource.zeroExtended(
				N5Meta.fromReader(reader, dataset),
				transform,
				name,
				viewer.getQueue(),
				0,
				channelDimension,
				channels);
		return new ChannelSourceState<>(
				channelSource,
				new ARGBCompositeColorConverter.InvertingImp0<>(channels.length, min, max),
				new ARGBCompositeAlphaAdd(),
				name);
	} catch (final ReflectionException | DataTypeNotSupported e) {
		throw new IOException(e);
	}
}
 

/**
 * Adss a SourceAndConverter to the scene.
 *
 * @param sac The SourceAndConverter to add
 * @param name Name of the dataset
 * @param voxelDimensions Array with voxel dimensions.
 * @param <T> Type of the dataset.
 * @return THe node corresponding to the volume just added.
 */
public <T extends RealType<T>> Node addVolume(SourceAndConverter<T> sac,
                                              int numTimepoints,
                                              String name,
                                              float... voxelDimensions ) {
    List<SourceAndConverter<T>> sources = new ArrayList<>();
    sources.add(sac);

    return addVolume(sources, numTimepoints, name, voxelDimensions);
}
 

public static <T extends RealType<T>, S extends RealType<S>> RandomAccessibleInterval<FloatType> dummyFuse(RandomAccessibleInterval<T> img1, RandomAccessibleInterval<S> img2, PhaseCorrelationPeak2 shiftPeak, ExecutorService service)
{
	long[] shift = new long[img1.numDimensions()];
	shiftPeak.getShift().localize(shift);
	long[] minImg1 = new long[img1.numDimensions()];
	long[] minImg2 = new long[img1.numDimensions()];
	long[] maxImg1 = new long[img1.numDimensions()];
	long[] maxImg2 = new long[img1.numDimensions()];
	long[] min = new long[img1.numDimensions()];
	long[] max = new long[img1.numDimensions()];
	
	for (int i = 0; i < img1.numDimensions(); i++){
		minImg1[i] = 0;
		maxImg1[i] = img1.dimension(i) -1;
		minImg2[i] = shiftPeak.getShift().getLongPosition(i);
		maxImg2[i] = img2.dimension(i) + minImg2[i] - 1;
		
		min[i] =  Math.min(minImg1[i], minImg2[i]);
		max[i] = Math.max(maxImg1[i], maxImg2[i]);
	}
	
	
	RandomAccessibleInterval<FloatType> res = new ArrayImgFactory<FloatType>().create(new FinalInterval(min, max), new FloatType());
	copyRealImage(Views.iterable(img1), Views.translate(res, min), service);
	copyRealImage(Views.iterable(Views.translate(img2, shift)), Views.translate(res, min), service);
	return res;	
	
}
 

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

/**
 * Add an IterableInterval as a Volume
 * @param image
 * @param <T>
 * @return a Node corresponding to the Volume
 */
public <T extends RealType<T>> Node addVolume( IterableInterval<T> image ) throws Exception {
    if( image instanceof RandomAccessibleInterval ) {
        return addVolume((RandomAccessibleInterval) image, "Volume");
    } else {
        throw new Exception("Unsupported Volume type:" + image);
    }
}