imglib / imglib2

A generic next-generation Java library for image processing
http://imglib2.net/
Other
293 stars 93 forks source link

feat: add RadialKDTreeInterpolatorFactory #360

Closed bogovicj closed 1 month ago

bogovicj commented 4 months ago

This PR adds an interpolator for KDTrees that takes a custom function of the (squared) radius to produce a RealRandomAccessible. This is similar to but more general than the existing InverseDistanceWeightingInterpolatorFactory which uses a particular class of functions for the radius. Another difference - InverseDistanceWeightingInterpolatorFactory uses RealType, the classes provided by this PR can use any NumericType.

The demo below produces this output: Screenshot from 2024-03-08 19-33-26

// generate points
int N = 12; // number of points
final double bigCircleRadius = 10;
final double angleDelta = 2 * Math.PI / N;
final List<RealPoint> pts = DoubleStream.iterate(0.0, x -> x + angleDelta).limit(N).mapToObj(tht -> {
    return new RealPoint(
            bigCircleRadius * Math.cos(tht),
            bigCircleRadius * Math.sin(tht));
}).collect(Collectors.toList());    

// generate values
final FloatType one = new FloatType(1f);
final ArrayList<FloatType> vals = new ArrayList<>();
IntStream.range(0, pts.size()).forEach( x -> { vals.add(one); });

// build a tree
final KDTree<FloatType> tree = new KDTree<>(vals, pts);

// a function of the squared radius for each rendered point
final double value = 2.0;
final double maxRadius = 2.0;
final double maxSqrRadius = maxRadius*maxRadius;
final DoubleUnaryOperator sqrRadiusFunction = r2 -> value * (1 + Math.cos(2 * Math.PI * r2 / maxRadius));

// create an image from the points and the radial function
final RealRandomAccessible<FloatType> image = Views.interpolate(
        tree, 
        new RadialKDTreeInterpolatorFactory<>(squaredRadiusFunction, maxRadius, new FloatType()));

// show the image
BdvStackSource<FloatType> bdv = BdvFunctions.show( 
        image, 
        Intervals.createMinSize(-20,-20, 40, 40),
        "rendered points",
        BdvOptions.options().is2D());

bdv.setDisplayRange(0, 2);
bogovicj commented 4 months ago

Another point: existing code does a nearest neighbor search but this PR uses a radius search.

InverseDistanceWeightingInterpolatorFactory implements InterpolatorFactory< T, KNearestNeighborSearch<T>> emphasis on the KNearestNeighborSearch.

In this PR: RadialKDTreeInterpolatorFactory implements InterpolatorFactory<T, KDTree<T>> but internally uses a RadiusNeighborSearchOnKDTree

we could make RadialKDTreeInterpolatorFactory implement InterpolatorFactory<T, RadiusNeighborSearchOnKDTree<T>>, if that's preferable.

StephanPreibisch commented 4 months ago

Hi, what is the math that defines what happens when the points overlap? I think there are some similarities with what I have written in STIM ... would be great to catch up on that (https://github.com/PreibischLab/STIM/blob/master/src/main/java/render/Render.java) ...

bogovicj commented 4 months ago

@StephanPreibisch ,

Yea, the code in this PR is similar to what you have in STIM, in the same ways that it's similar to what is here in imglib2. It's also different from what's in STIM in the ways its different from the imglib2 code. Some of your code there could do well in imglib2 core in my view. Summary of the differences

  1. general function of the squared radius
  2. uses radius search instead of neighbor (or combined radius-neighbor search
  3. Can work with NumericTypes (current impls use RealTypes

The math is:

$$ \sum_{p \in P} f (r^2) v(p) $$

where $P$ is all the points within some radius. and $v(p)$ is whatever value the point $p$ has (needed to build the KDTree). and $f$ is a developer-provided function of the squared radius.