Skip to content

Commit b7e5cdc

Browse files
gselzergselzer
committed
Beautify Frangi vesselness op
1 parent 61ca449 commit b7e5cdc

1 file changed

Lines changed: 19 additions & 38 deletions

File tree

‎scijava-ops-image/src/main/java/org/scijava/ops/image/filter/vesselness/DefaultFrangi.java‎

Lines changed: 19 additions & 38 deletions
Original file line numberDiff line numberDiff line change
@@ -30,21 +30,18 @@
3030
package org.scijava.ops.image.filter.vesselness;
3131

3232
import java.util.ArrayList;
33-
import java.util.Collections;
3433
import java.util.Comparator;
3534

35+
import Jama.EigenvalueDecomposition;
36+
import Jama.Matrix;
3637
import net.imglib2.Cursor;
3738
import net.imglib2.RandomAccess;
3839
import net.imglib2.RandomAccessibleInterval;
3940
import net.imglib2.outofbounds.OutOfBoundsMirrorFactory;
4041
import net.imglib2.outofbounds.OutOfBoundsMirrorFactory.Boundary;
4142
import net.imglib2.type.numeric.RealType;
4243
import net.imglib2.view.Views;
43-
4444
import org.scijava.function.Computers;
45-
46-
import Jama.EigenvalueDecomposition;
47-
import Jama.Matrix;
4845
import org.scijava.ops.spi.Nullable;
4946

5047
/**
@@ -59,24 +56,13 @@
5956
public class DefaultFrangi<T extends RealType<T>, U extends RealType<U>>
6057
implements Computers.Arity3<RandomAccessibleInterval<T>, Integer, double[], RandomAccessibleInterval<U>> {
6158

62-
private double[] spacing;
63-
private int scale;
64-
6559
protected double alpha = 0.5;
6660
protected double beta = 0.5;
6761

6862
protected double minimumVesselness = Double.MIN_VALUE;
6963
protected double maximumVesselness = Double.MAX_VALUE;
7064

71-
public double getMinimumVesselness() {
72-
return minimumVesselness;
73-
}
74-
75-
public double getMaximumVesselness() {
76-
return maximumVesselness;
77-
}
78-
79-
private double getDistance(RandomAccess<T> ahead, RandomAccess<T> behind, int d) {
65+
private double getDistance(RandomAccess<T> ahead, RandomAccess<T> behind, int d, double[] spacing) {
8066
double distance = 0;
8167

8268
for (int i = 0; i < d; i++) {
@@ -97,10 +83,10 @@ private double derive(double val1, double val2, double distance) {
9783
/**
9884
* TODO
9985
*
100-
* @param input
101-
* @param scale size of vessels to search for
86+
* @param input the image containing input data
87+
* @param scale sigma for the gaussian filter, and the scale for the vesselness filter
10288
* @param spacing physical distance between data points
103-
* @param output
89+
* @param output the pre-allocated output buffer
10490
*/
10591
@Override
10692
public void compute(final RandomAccessibleInterval<T> input, final Integer scale,
@@ -111,26 +97,23 @@ public void compute(final RandomAccessibleInterval<T> input, final Integer scale
11197

11298
// set spacing if the parameter is not passed.
11399
if (spacing == null) {
114-
this.spacing = new double[input.numDimensions()];
100+
spacing = new double[input.numDimensions()];
115101
for (int i = 0; i < input.numDimensions(); i++)
116-
this.spacing[i] = 1;
117-
} else
118-
this.spacing = spacing;
119-
120-
this.scale = scale;
102+
spacing[i] = 1;
103+
}
121104

122-
frangi(input, output, this.scale);
105+
frangi(input, output, spacing, scale);
123106
}
124107

125-
private final void frangi(RandomAccessibleInterval<T> in, RandomAccessibleInterval<U> out, int step) {
108+
private void frangi(RandomAccessibleInterval<T> in,
109+
RandomAccessibleInterval<U> out, double[] spacing, int step) {
126110

127111
// create denominators used for gaussians later.
128112
double ad = 2 * alpha * alpha;
129113
double bd = 2 * beta * beta;
130114

131115
// OutOfBoundsMirrorStrategy for use when the cursor reaches the edges.
132-
OutOfBoundsMirrorFactory<T, RandomAccessibleInterval<T>> osmf = new OutOfBoundsMirrorFactory<T, RandomAccessibleInterval<T>>(
133-
Boundary.SINGLE);
116+
OutOfBoundsMirrorFactory<T, RandomAccessibleInterval<T>> osmf = new OutOfBoundsMirrorFactory<>(Boundary.SINGLE);
134117

135118
Cursor<T> cursor = Views.iterable(in).localizingCursor();
136119

@@ -164,7 +147,7 @@ private final void frangi(RandomAccessibleInterval<T> in, RandomAccessibleInterv
164147

165148
// take the derivative between the two points
166149
double derivativeA = derive(behind.get().getRealDouble(), current.get().getRealDouble(),
167-
getDistance(behind, current, in.numDimensions()));
150+
getDistance(behind, current, in.numDimensions(), spacing));
168151

169152
// move one ahead to take the other first derivative
170153
ahead.move(step, m);
@@ -173,11 +156,11 @@ private final void frangi(RandomAccessibleInterval<T> in, RandomAccessibleInterv
173156

174157
// take the derivative between the two points
175158
double derivativeB = derive(current.get().getRealDouble(), ahead.get().getRealDouble(),
176-
getDistance(current, ahead, in.numDimensions()));
159+
getDistance(current, ahead, in.numDimensions(), spacing));
177160

178161
// take the second derivative using the two first derivatives
179162
double derivative2 = derive(derivativeA, derivativeB,
180-
getDistance(behind, ahead, in.numDimensions()));
163+
getDistance(behind, ahead, in.numDimensions(), spacing));
181164

182165
hessian.set(m, n, derivative2);
183166

@@ -191,10 +174,10 @@ private final void frangi(RandomAccessibleInterval<T> in, RandomAccessibleInterv
191174
// find and sort the eigenvalues and eigenvectors of the Hessian
192175
EigenvalueDecomposition e = hessian.eig();
193176
double[] eigenvaluesArray = e.getRealEigenvalues();
194-
ArrayList<Double> eigenvaluesArrayList = new ArrayList<Double>();
177+
ArrayList<Double> eigenvaluesArrayList = new ArrayList<>();
195178
for (double d : eigenvaluesArray)
196179
eigenvaluesArrayList.add(d);
197-
Collections.sort(eigenvaluesArrayList, Comparator.comparingDouble(Math::abs));
180+
eigenvaluesArrayList.sort(Comparator.comparingDouble(Math::abs));
198181

199182
// vesselness value
200183
double v = 0;
@@ -256,9 +239,7 @@ private final void frangi(RandomAccessibleInterval<T> in, RandomAccessibleInterv
256239
v = (1 - Math.exp(an / ad)) * Math.exp(bn / bd) * (1 - Math.exp(cn / cd));
257240
}
258241

259-
} else
260-
261-
if (!Double.isNaN(v)) {
242+
} else {
262243
maximumVesselness = Math.max(v, maximumVesselness);
263244
minimumVesselness = Math.min(v, minimumVesselness);
264245
}

0 commit comments

Comments
 (0)