This post shows how to create and apply non-linear transformations with ImgLib2, specifically using DisplacementFieldTransforms and PositionFieldTransforms.

The last example of this post will show how to use a PositionFieldTransform to create a 2D image from a 1D signal by a transformation of coordinates.

A 2D image created after transformation	of this 1D function

The code folded below sets up dependencies, imports classes, and makes it convenient to show image output as shown in this post.

Code

%mavenRepo scijava.public https://maven.scijava.org/content/groups/public

%maven net.imagej:ij:1.54f
%maven org.scijava:scijava-common:2.97.0
%maven net.imglib2:imglib2:6.2.0
%maven net.imglib2:imglib2-ij:2.0.1
%maven net.imglib2:imglib2-realtransform:4.0.1
%maven org.knowm.xchart:xchart:3.5.4
    
import java.util.function.*;
import java.util.stream.*;
import java.awt.image.BufferedImage;

import io.github.spencerpark.jupyter.kernel.display.common.*;
import io.github.spencerpark.jupyter.kernel.display.mime.*;

import org.knowm.xchart.*;
import org.knowm.xchart.style.Styler.LegendPosition;
import org.knowm.xchart.style.markers.SeriesMarkers;

import ij.*;
import net.imglib2.*;
import net.imglib2.view.*;
import net.imglib2.util.*;
import net.imglib2.realtransform.*;
import net.imglib2.position.*;
import net.imglib2.type.numeric.*;
import net.imglib2.type.numeric.real.*;
import net.imglib2.type.numeric.integer.*;
import net.imglib2.interpolation.randomaccess.*;
import net.imglib2.img.display.imagej.*;
import net.imglib2.img.imageplus.ImagePlusImgs;
import net.imglib2.img.array.ArrayImgs;

public static <T extends NumericType<T>> BufferedImage toBuffferedImage(RandomAccessibleInterval rai) {

    return ImageJFunctions.wrap((RandomAccessibleInterval<T>)rai, "").getBufferedImage();
}

getKernelInstance().getRenderer().createRegistration(RandomAccessibleInterval.class)
        .preferring(MIMEType.IMAGE_PNG)
        .supporting(MIMEType.IMAGE_JPEG, MIMEType.IMAGE_GIF)
        .register((rai, context) -> Image.renderImage(toBuffferedImage(rai),context));

var printWriter = new PrintWriter(System.out,true);
char rarrow = '\u2192';

public static void printTformPts(RealPoint x, RealPoint y) {
    
    printWriter.println(String.format("%s %s %s", x, rarrow, y));
};

Displacement fields

A displacement field is the most common way to represent and store non-linear transformations. Imagine an image where a vector is stored at every position \(\mathbf{x}\). That vector describes how much and in what direction the position \(\mathbf{x}\) should be translated, or displaced.

\[ \begin{align} \mathbf{x} &: \text{a position } \\ \mathbf{v}( \mathbf{x} ) &: \text{the vector at position } \mathbf{x}\\ \mathbf{x} + \mathbf{v}( \mathbf{x} ) &: \text{the output of the transformation} \end{align} \]

For example, let’s consider an example vector field with the vector \([50, -25]\) at every position. This will be equivalent to a simple global translation and is not the recommended way to represent a translation, but will be instructive thanks to its simplicity.

A FunctionRealRandomAccessible is one way to make a image containing a constant value (ConstantUtils is another way). This image can be passed directly to a DisplacementFieldTransform.

var displacementVector = new double[]{50, -25};        // the constant vector

// an image that takes the value of displacementVector everywhere
var constantVector = new FunctionRealRandomAccessible<>(2,
        (x, v) -> { v.setPosition(displacementVector); }, // set the vector 
        () -> { return DoubleType.createVector(2); });     // make a 2d vector

// a constant displacement field
var dfieldConstant = new DisplacementFieldTransform(constantVector);

Applying the transformation to any point translates that point by the same amount \([50, -25]\) , as expected:

// two ways of initializing a point (0.0, 0.0)
var x = new RealPoint(0.0, 0.0); // give the coordinates as arguments
var y = new RealPoint(2);        // give the number of dimensions as an argument

// transform x, store the result in y
dfieldConstant.apply(x, y);
 // [0, 0] is transformed to [50, -25]
printTformPts(x, y);

// [-50, 25] is transformed to [0, 0]
x = new RealPoint(-50, 25);
dfieldConstant.apply(x, y);
printTformPts(x, y);

// [pi, 100k] is transformed to [pi + 50, 100k - 25]
x = new RealPoint(Math.PI, 100000);
dfieldConstant.apply(x, y);
printTformPts(x, y);

(0.0,0.0) → (50.0,-25.0)
(-50.0,25.0) → (0.0,0.0)
(3.141592653589793,100000.0) → (53.1415926535898,99975.0)

Let’s make a slightly more interesting displacement field with four vectors - one at each corner. Our image will be:

----------------
|[1, 1]  [2, 2]|
|[3, 3]  [4, 4]|
----------------

i.e. the vector at position (0,1) is \([3,3]\).

Our data needs to be 3D, where the first dimension holds the two vector components, and the other two hold the spatial dimensions. We can use an ArrayImg to store these data.

var dfieldData = ArrayImgs.doubles(
        new double[]{1, 1, 2, 2, 3, 3, 4, 4}, // the data
        2, 2, 2);  // the dimensions

// the displacement field
var dfieldCorners = new DisplacementFieldTransform(dfieldData);

For example, we expect the transformation to take the point \([0,1]\) to \([3,4] = [0,1] + [3,3]\), and that is indeed what we see:

// [0,0] is transformed to [1,1] = [0,0] + [1,1]
var x = new RealPoint(0.0, 0.0);
var y = new RealPoint(2);
dfieldCorners.apply(x, y);
printTformPts(x, y);

// [1,0] is transformed to [3,2] = [1,0] + [2,2]
var x = new RealPoint(1.0, 0.0);
dfieldCorners.apply(x, y);
printTformPts(x, y);

// [0,1] is transformed to [3,4] = [0,1] + [3,3]
var x = new RealPoint(0.0, 1.0);
dfieldCorners.apply(x, y);
printTformPts(x, y);

// [1,1] is transformed to [5,5] = [1,1] + [4,4]
var x = new RealPoint(1.0, 1.0);
dfieldCorners.apply(x, y);
printTformPts(x, y);

(0.0,0.0) → (1.0,1.0)
(1.0,0.0) → (3.0,2.0)
(0.0,1.0) → (3.0,4.0)
(1.0,1.0) → (5.0,5.0)

What happens if we try to apply the transformation “in between” the discrete values of the array, or out-of-bounds of the array?

// [0.5, 0.5] is transformed to [3.0, 3.0]
var x = new RealPoint(0.5, 0.5);
var y = new RealPoint(2);
dfieldCorners.apply(x, y);
printTformPts(x, y);

// [-100, -100] is transformed to [-99.0,-99.0]
var x = new RealPoint(-100, -100);
dfieldCorners.apply(x, y);
printTformPts(x, y);

(0.5,0.5) → (3.0,3.0)
(-100.0,-100.0) → (-99.0,-99.0)

By default, linear interpolation determines the displacement within the array but off of grid values. The nearest value on the border of the array determines the displacement outside the array’s bounds. If, instead, we want nearest-neighbor interpolation and all out-of-bounds displacements to be zero, it can be acheived by:

// interpret our 3d image as a 2d image of vectors
var vectorsFrom3D = Views.collapseReal(Views.moveAxis(dfieldData, 0, dfieldData.numDimensions() - 1));

// make the displacement field with custom interpolation and boundary extension
var dfieldCornersCustom = new DisplacementFieldTransform(
        Views.interpolate(
                Views.extendZero(vectorsFrom3D),            // zeros out-of-bounds
                new NearestNeighborInterpolatorFactory())); // nearest-neighbor interpolation

// [0.4, 0.4] is transformed to [1.4, 1.4]
var x = new RealPoint(0.4, 0.4);
var y = new RealPoint(2);
dfieldCornersCustom.apply(x, y);
printTformPts(x, y);

// [-100, -100] is transformed to [-100, -100]
var x = new RealPoint(-100, -100);
dfieldCornersCustom.apply(x, y);
printTformPts(x, y);

(0.4,0.4) → (1.4,1.4)
(-100.0,-100.0) → (-100.0,-100.0)

Another common requirement is to specify the spacing and offset of the displacement field’s discrete grid. For example, suppose we want the displacements above to be applied to the field of view \([-10, 10] \times [-20, 20]\) rather than the array coordinates \([0, 1] \times [0, 1]\).

(-10,-20)  ------------------------------------------------  (10,-20)
           |[1,1]                                    [2,2]|
           |                                              |                                          
           |                                              |
           |                                              |
           |            (interpolated values)             |
           |                                              |
           |                                              |
           |                                              |
           |[3,3]                                    [4,4]|
(-10, 20)  ------------------------------------------------ (10,20)

where the values outside the box in parenthesis are the coordinates of the corners.

This is possible by passing the desired spacing and offset directly to the DisplacementFieldTransform.

/*
 * the spacing and offset below map the array origin [0,0] to [-10,-20]
 * and the element at array index [1,1] to [10,20]
 */
var spacing = new double[]{20, 40};
var offset = new double[]{-10, -20};
var dfieldCornersSpacingOffset = new DisplacementFieldTransform(dfieldData, spacing, offset);

// [-10,-20] is displaced by the vector [1,1], so goes to [-9,-19]
var x = new RealPoint(-10, -20);
var y = new RealPoint(2);
dfieldCornersSpacingOffset.apply(x, y);
printTformPts(x, y);

// similarly, [10,20] is displaced by the vector [4,4], so goes to [14,24]
var x = new RealPoint(10, 20);
dfieldCornersSpacingOffset.apply(x, y);
printTformPts(x, y);

(-10.0,-20.0) → (-9.0,-19.0)
(10.0,20.0) → (14.0,24.0)

Transforming images

Displacement fields can also be applied to images. First, let’s open an image:

var img = ImageJFunctions.wrap(IJ.openImage("https://mirror.imagej.net/ij/images/boats.gif"));
img

Now let’s transform it using our constant valued displacement field. Remember this field has the vector \([50, -25]\) everwhere.

First, we’ll define two functions for transforming images:

/*
 * Transform an image of arbitrary size (RandomAccessible) and output an image of size
 * given by resultInterval.
 */ 
public RandomAccessibleInterval transformImageInterval(
        RandomAccessible img,
        RealTransform transform,
        Interval resultInterval) {
    
    // use n-linear interpolation
    var interpImg = Views.interpolate(img, new NLinearInterpolatorFactory());

    // transform the image
    var transformedImg = new RealTransformRandomAccessible(interpImg, transform);

    // rasterize and set bounding box
    return Views.interval(Views.raster(transformedImg), resultInterval);
}

public RandomAccessibleInterval transformImage(RandomAccessibleInterval img, RealTransform transform) {
    
    // default out-of-bounds extension and output interval
    return transformImageInterval(Views.extendZero(img), transform, img);
}

transformImage(img, dfieldConstant);

Notice that the image is shifted down (+y direction) and to the left (-x direction) : the opposite direction of the displacement vector. This is because the “inverse” transformation (from output coordinates to input coordinates) is needed to transform images. Learn why here.

Now let’s use what we learned above to make a displacement field whose corners are the corners of the image. We do this by setting the spacing of the displacement field’s coordinate grid as we learned above. This means displacement vector stored at at the array coordinate [1,1] to a new position: [imageWidth, imageHeight].

// some vector field
var dfieldData = ArrayImgs.doubles( 
        new double[]{-45, -50, 35, -35, -25, 25, 70, 75},  // the vector data
        2, 2, 2); // the dimensions
    
// use the image width and height as the vector fields spacing
var dfield = new DisplacementFieldTransform(
        dfieldData, 
        img.dimension(0),
        img.dimension(1));

// transform the image
transformImage(img, dfield);

Finally, let’s make a field of random displacements at one-tenth the resolution of the image and apply it.

// create image of random vectors
var randDfieldData = ArrayImgs.doubles(2, (img.dimension(0) / 20), (img.dimension(1) / 20));
randDfieldData.forEach(x -> {x.set(15 * (Math.random() - 0.5));});
    
// spacing of 20
var dfield = new DisplacementFieldTransform(randDfieldData, 20, 20);

// transform the image
transformImage(img, dfield);

Position fields

A position, or coordinate field is similar to a displacement field in that it is also represented by a field of vectors, but those vectors represent positions directly, rather than displacements of the current position.

\[ \begin{align} \mathbf{x} &: \text{a position } \\ \mathbf{v}( \mathbf{x} ) &: \text{the vector at position } \mathbf{x}\\ \mathbf{v}( \mathbf{x} ) &: \text{the output of the transformation} \end{align} \]

If we use a constant vector field to make a PositionFieldTransform the output will always be the same:

var coordinateVector = new double[]{1, 2};
var constantVector = new FunctionRealRandomAccessible<>(
        2,
        (x, v) -> {v.setPosition( coordinateVector );},  // set the vector 
        () -> {return DoubleType.createVector(2);});     // make a 2d vector

var pfieldConstant = new PositionFieldTransform(constantVector);

// [0, 0] is transformed to [1, 2]
var x = new RealPoint(0, 0);
var y = new RealPoint(2);
pfieldConstant.apply(x, y);
printTformPts(x, y);

// [9999, 9999] is also transformed to [1, 2]
var x = new RealPoint(9999, 9999);
var y = new RealPoint(2);
pfieldConstant.apply(x, y);
printTformPts(x, y);

(0.0,0.0) → (1.0,2.0)
(9999.0,9999.0) → (1.0,2.0)

The Localizables class has some convenience methods for producing images of coordinates. Creating a position field with this image gives the identity. FunctionRandomAccessibles can also be used to generate coordinate images.

Challenge: Try using Converters to extract the x and y coordinates from the coordinates variable.

// make an image of coordinates. The argument "2" specifies the number of dimensions.
var coordinates = Localizables.realRandomAccessible(2);

// make an image whose intensities are the x coordinate
var xCoordinates = Views.interval(
        new FunctionRandomAccessible<>(
                2,
                (x, v) -> {v.set(x.getIntPosition(0));}, 
                UnsignedByteType::new),
        new FinalInterval(255, 255));

// make an image whose intensities are the x coordinate
var yCoordinates = Views.interval(
        new FunctionRandomAccessible<>(
                2,
                (x, v) -> {v.set(x.getIntPosition(1));}, 
                UnsignedByteType::new),
        new FinalInterval(255, 255));

display("x coordinates");
display(xCoordinates );
display(" ");

display("y coordinates");
yCoordinates

x coordinates

y coordinates

// the identity transformation
var identityPfield = new PositionFieldTransform(Localizables.realRandomAccessible(2));

// [0, 0] is transformed to [0, 0]
var x = new RealPoint(0, 0);
var y = new RealPoint(2);
identityPfield.apply(x, y);
printTformPts(x, y);

// [9999, 9999] is transformed to [9999, 9999]
var x = new RealPoint(9999, 9999);
identityPfield.apply(x, y);
printTformPts(x, y);

(0.0,0.0) → (0.0,0.0)
(9999.0,9999.0) → (9999.0,9999.0)

For the above examples, notice that the identityPfield does indeed behaves as the identity transformation.

Let’s make another position field that stretches and compresses the image in a non-linear way.

/*
 * The details of this vector field are not so important,
 * but notice that it is some non-linear function.
 */ 
var pFieldVectorField = new FunctionRealRandomAccessible<>(
        2,
        (p, v) -> {
            var cx = img.dimension(0) / 2.0;
            var cy = img.dimension(1) / 2.0;
            var ex = Math.exp(-0.020 * (p.getDoublePosition(0) - cx));
            var ey = Math.exp(-0.020 * (p.getDoublePosition(1) - cy));
            v.setPosition(10 + 700 / (1 + ex), 0);
            v.setPosition(500 / ( 1 + ey), 1);
        },
        () -> {return DoubleType.createVector(2);});
    
// spacing of 10
var spacing = new double[]{1, 1};
var offset = new double[]{img.dimension(0) / 2.0, img.dimension(1) / 2.0};
var pfield = new PositionFieldTransform(pFieldVectorField);

// transform the image
transformImage(img, pfield);

Position fields can have different input and output dimensions. The length of the vector (first) dimension of the position field indicates its output dimension. The domain of the field indicates the input dimension. So generally, an N+1 dimensional image has N-dimensional inputs since one dimension is the vector dimension. Let’s consider an example.

var randPfieldData = ArrayImgs.doubles(2, 20, 30, 40, 50); // make a 5d array

var pfieldRand = new PositionFieldTransform(randPfieldData);
System.out.println("input dimensions:  " + pfieldRand.numSourceDimensions());
System.out.println("output dimensions: " + pfieldRand.numTargetDimensions());

input dimensions:  4
output dimensions: 2

Make a 1D Function and plot it:

// a 1d signal
var sin1d = new FunctionRandomAccessible<>(
        1,
        (p, v) -> {
            double t = p.getDoublePosition(0);
            if (t > 250) v.set(64);
            else v.setReal(150 + ((t/5) * Math.sin(t/10)));
        },
        UnsignedByteType::new);

Code

// plot the 1d function
var xpts = new double[256];
var ypts = new double[256];

var access = sin1d.randomAccess();
IntStream.range(0, 256).forEach(
        i -> {
            xpts[i] = i;
            access.fwd(0);
            ypts[i] = access.get().getRealDouble();
        });

var chart = new XYChartBuilder().width(800).height(600).build();
chart.getStyler().setChartTitleVisible(false);
chart.getStyler().setLegendVisible(false);
chart.getStyler().setChartTitleBoxBackgroundColor(java.awt.Color.white);
chart.getStyler().setChartBackgroundColor(java.awt.Color.white);

var series = chart.addSeries("function", xpts, ypts);
series.setMarker(SeriesMarkers.NONE);
BitmapEncoder.getBufferedImage(chart)

/* 
 * Create the position field (1D -> 2D)
 *
 * The details of this vector field are not so important,
 * but notice:
 *    1) it is some non-linear function.
 *    2) its input p is a 2D; notice "p.getDoublePosition(i)"
 *    3) its output v is 1D (scalar)
 */ 
var pFieldTransition = new PositionFieldTransform( 
    new FunctionRealRandomAccessible<>(
            2, 
            (p, v) -> {
                double x = p.getDoublePosition(0) - 300;
                double y = p.getDoublePosition(1) - 300;
                double r = Math.sqrt(x*x + y*y);
                double tht = Math.atan2(x, y);
                v.setPosition(1.2 * r +  32 * tht, 0); 
            },
            () -> {return DoubleType.createVector(1);})
);

// transform the 1D signal, make a 2D image of size 512 x 512
transformImageInterval(sin1d, pFieldTransition, new FinalInterval(512, 512));