Package gov.nih.mipav.model.algorithms
Class AlgorithmAntigradient2
java.lang.Object
java.lang.Thread
gov.nih.mipav.model.algorithms.AlgorithmBase
gov.nih.mipav.model.algorithms.AlgorithmAntigradient2
- All Implemented Interfaces:
ActionListener,WindowListener,Runnable,EventListener
This is a port of the files anitgradient2.m and antigradient2.c created by Gunnar Farneback in the
Spatial domain toolbox at http://www.imt.liu.se/mi/Tools. Email permission to do so was granted by
Gunnar Farneback:
On 03/15/2013 08:24 PM, Gandler, William (NIH/CIT) [E] wrote:
> Dear Gunnar Farneback:
>
> I am a Java programmer for the MIPAV software project at the National
> Institutes of Health. I would like to port your inverse gradient C and
> MATLAB routines in the Spatial Domain Toolbox to Java for inclusion in
> the MIPAV imaging analysis package. However, while MIPAV is open
> source, it does not distribute its code under the GNU General public
> license and has no plans to do so. Our open source license is:
>
> http://mipav.cit.nih.gov/license/license.html
>
> Would it be possible to obtain permission to distribute a port of your
> inverse gradient code in MIPAV without the GNU General public license
> as long as we give full credit to you?
Yes, go ahead.
/Gunnar
On 04/29/2013 11:32 PM, Gandler, William (NIH/CIT) [E] wrote:
> Gunnar,
>
> I have ported Antigradient2.c to MIPAV and found its performance to be very impressive.
Would there be any reason to port Antigradient.c as well? Does Antigradient.c have any capability
not possessed by Antigradient2.c?
>
>
> Sincerely,
>
>
> William Gandler
The difference between antigradient.c and antigradient2.c is that the former doesn't allow
a mask and therefore only can handle rectangular domains. Not having a mask is a big limitation
but it does simplify the algorithm drastically, making the code simpler and faster. On the other hand,
for rectangular domains it's entirely possible that there are even faster transform-based methods.
I mostly implemented antigradient.c to learn the multigrid approach.
Computing inverse gradients on irregular domains was central to my research at the time.
I'm happy to hear that it performs well.
/Gunnar
function f = antigradient2(g, mask, mu, n)
% ANTIGRADIENT2 Reconstruction from gradients in 2D and 3D.
%
% This differs from antigradient in that it can handle arbitrarily
% shaped domains.
%
% f = antigradient2(g) computes the function f which has a gradient
% as close as possible to g. The vector field g must either be an
% MxNx2 array or an MxNxPx3 array.
%
% f = antigradient2(g, mask) does the same thing but restricted to
% points where mask is nonzero. For meaningful results the domain must
% be simply connected (otherwise run it multiple times for each
% connected component). Certain shapes may give poor convergence; in
% general almost convex shapes can be expected to have the best
% properties. The mask should have size MxN or MxNxP. It can also be
% set to an empty matrix or omitted, in which case all points are
% included in the domain.
%
% f = antigradient2(g, mask, mu) sets the unrecoverable
% mean of f to mu. When omitted, the mean of f is set to zero.
%
% f = antigradient2(g, mask, mu, n) uses n multigrid iterations. The
% default value of 2 iterations is fast but only moderately accurate.
% The algorithm converges quickly, however, so usually 10-15
% iterations are sufficient to reach full convergence. This may depend
% on the shape of the domain, however.
%
% This function is based on transforming the inverse gradient problem to a
% PDE - a Poisson equation with Neumann boundary conditions. This equation
% is solved by an implementation of the full multigrid algorithm.
%
% The implementation is most efficient when all sides are of approximately
% the same size. The worst case is when two sides are large and the
% third side much smaller. The reason is that only uniform
% downsampling is implemented and when the smallest side becomes
% smaller than 4, a direct solution is applied. It makes no big
% difference whether the sides are odd or even, although sides which
% are of the form 2^k+1 give somewhat faster convergence.
%
% Author: Gunnar Farneback
% Medical Informatics
% Linkoping University, Sweden
% gunnar@imt.liu.se
On 04/23/2013 12:08 AM, Gandler, William (NIH/CIT) [E] wrote:
> Dear Gunnar:
>
> In Antigradient2.c I note the 2 differences in poisson_multigrid2D and poisson_multigrid3D:
>
> 1.) poisson_multigrid2D has:
> // Initialize solution.
> memcpy(f_out, f, M * N * sizeof(*f_out)); but no corresponding line
> is found in poisson_multigrid3D.
Looks like an artefact from some earlier development phase of the code or a not completely removed experiment.
In all cases when poisson_multigrid2D is called, f_out and f point to the same memory, so the memcpy call doesn't do anything.
> 2.) poisson_multigrid2D has if (1) at the end but poisson_multigrid3D has if (0) at the end.
That's definitely leftovers from frustrated attempts to work around a bug that long plagued the code.
I never checked if it still had any effect after finally tracking the bug down. It's not supposed to make a difference.
/Gunnar
Reference:
"Efficient Computation of the Inverse Gradient on Irregular Domains" by Gunnar Farneback,
Joakim Rydell, Tino Ebbers, Mats Andersson, and Hans Knutsson.
-
Nested Class Summary
Nested classes/interfaces inherited from class java.lang.Thread
Thread.Builder, Thread.State, Thread.UncaughtExceptionHandler -
Field Summary
FieldsModifier and TypeFieldDescriptionprivate dataStructprivate booleanprivate doubleprivate double[]private BitSetstatic final intprivate intprivate static final intprivate ModelImage(package private) booleanprivate intprivate intprivate intFields inherited from class gov.nih.mipav.model.algorithms.AlgorithmBase
destFlag, destImage, image25D, maxProgressValue, minProgressValue, multiThreadingEnabled, nthreads, progress, progressModulus, progressStep, runningInSeparateThread, separable, srcImage, threadStoppedFields inherited from class java.lang.Thread
MAX_PRIORITY, MIN_PRIORITY, NORM_PRIORITY -
Constructor Summary
ConstructorsConstructorDescriptionAlgorithmAntigradient2(ModelImage destImage, ModelImage srcImg, boolean entireImage, double fMean, int number_of_iterations) AlgorithmAntigradient2(ModelImage destImage, ModelImage srcImg, boolean entireImage, double fMean, int number_of_iterations, boolean testMode) -
Method Summary
Modifier and TypeMethodDescriptionprivate voidprivate voidprivate voidcompute_residual2D(double[] r, double[] A, double[] d, double[] f, int M, int N) private voidcompute_residual3D(double[] r, double[] A, double[] d, double[] f, int M, int N, int P) private voiddownsample2D(double[] rhs, int M, int N, double[] rhs_coarse, int Mhalf, int Nhalf, double[] weight, double[] coarse_weight) private voiddownsample3D(double[] rhs, int M, int N, int P, double[] rhs_coarse, int Mhalf, int Nhalf, int Phalf, double[] weight, double[] coarse_weight) private voidgalerkin2D(int level, int M, int N, int Mhalf, int Nhalf, double[] weight, double[] coarse_weight) private voidgalerkin3D(int level, int M, int N, int P, int Mhalf, int Nhalf, int Phalf, double[] weight, double[] coarse_weight) private voidgauss_seidel2D(double[] f, double[] A, double[] d, int M, int N) private voidgauss_seidel3D(double[] f, double[] A, double[] d, int M, int N, int P) private voidpoisson_full_multigrid2D(int level, double[] rhs, double[] weight, int number_of_iterations, int M, int N, double[] f_out) private voidpoisson_full_multigrid3D(int level, double[] rhs, double[] weight, int number_of_iterations, int M, int N, int P, double[] f_out) private voidpoisson_multigrid2D(double[] f, int level, double[] rhs, double[] weight, int n1, int n2, int nm, double[] f_out, int M, int N, int[] directly_solved) private voidpoisson_multigrid3D(int level, double[] rhs, double[] weight, int n1, int n2, int nm, double[] f_out, int M, int N, int P, int[] directly_solved) voidvoidvoidActually runs the algorithm.private voidsolve_directly2D(double[] lhs, double[] rhs, double[] f_out, int M, int N) private voidsolve_directly3D(double[] lhs, double[] rhs, double[] f_out, int M, int N, int P) private voidupsample2D(int M, int N, double[] v, int Mhalf, int Nhalf, double[] f_out, double[] weight, double[] coarse_weight) private voidupsample3D(int M, int N, int P, double[] v, int Mhalf, int Nhalf, int Phalf, double[] f_out, double[] weight, double[] coarse_weight) Methods inherited from class gov.nih.mipav.model.algorithms.AlgorithmBase
actionPerformed, addListener, addProgressChangeListener, calculateImageSize, calculatePrincipleAxis, computeElapsedTime, computeElapsedTime, convertIntoFloat, delinkProgressToAlgorithm, delinkProgressToAlgorithmMulti, displayError, errorCleanUp, finalize, fireProgressStateChanged, fireProgressStateChanged, fireProgressStateChanged, fireProgressStateChanged, fireProgressStateChanged, generateProgressValues, getDestImage, getElapsedTime, getMask, getMaxProgressValue, getMinProgressValue, getNumberOfThreads, getProgress, getProgressChangeListener, getProgressChangeListeners, getProgressModulus, getProgressStep, getProgressValues, getSrcImage, isCompleted, isImage25D, isMultiThreadingEnabled, isRunningInSeparateThread, isThreadStopped, linkProgressToAlgorithm, linkProgressToAlgorithm, makeProgress, notifyListeners, removeListener, removeProgressChangeListener, run, setCompleted, setImage25D, setMask, setMaxProgressValue, setMinProgressValue, setMultiThreadingEnabled, setNumberOfThreads, setProgress, setProgressModulus, setProgressStep, setProgressValues, setProgressValues, setRunningInSeparateThread, setSrcImage, setStartTime, setThreadStopped, startMethod, windowActivated, windowClosed, windowClosing, windowDeactivated, windowDeiconified, windowIconified, windowOpenedMethods inherited from class java.lang.Thread
activeCount, checkAccess, clone, countStackFrames, currentThread, dumpStack, enumerate, getAllStackTraces, getContextClassLoader, getDefaultUncaughtExceptionHandler, getId, getName, getPriority, getStackTrace, getState, getThreadGroup, getUncaughtExceptionHandler, holdsLock, interrupt, interrupted, isAlive, isDaemon, isInterrupted, isVirtual, join, join, join, join, ofPlatform, ofVirtual, onSpinWait, resume, setContextClassLoader, setDaemon, setDefaultUncaughtExceptionHandler, setName, setPriority, setUncaughtExceptionHandler, sleep, sleep, sleep, start, startVirtualThread, stop, suspend, threadId, toString, yield
-
Field Details
-
RECURSION_SIZE_LIMIT
private static final int RECURSION_SIZE_LIMIT- See Also:
-
MAX_LEVELS
public static final int MAX_LEVELS- See Also:
-
entireImage
private boolean entireImage -
mask
-
fMean
private double fMean -
number_of_iterations
private int number_of_iterations -
xDim
private int xDim -
yDim
private int yDim -
zDim
private int zDim -
g
private double[] g -
data
-
testImage
-
testMode
boolean testMode
-
-
Constructor Details
-
AlgorithmAntigradient2
public AlgorithmAntigradient2(ModelImage destImage, ModelImage srcImg, boolean entireImage, double fMean, int number_of_iterations) -
AlgorithmAntigradient2
public AlgorithmAntigradient2(ModelImage destImage, ModelImage srcImg, boolean entireImage, double fMean, int number_of_iterations, boolean testMode) -
AlgorithmAntigradient2
public AlgorithmAntigradient2()
-
-
Method Details
-
run2DSelfTest
public void run2DSelfTest() -
run3DSelfTest
public void run3DSelfTest() -
runAlgorithm
public void runAlgorithm()Description copied from class:AlgorithmBaseActually runs the algorithm. Implemented by inheriting algorithms.- Specified by:
runAlgorithmin classAlgorithmBase
-
antigradient2D
private void antigradient2D() -
poisson_full_multigrid2D
private void poisson_full_multigrid2D(int level, double[] rhs, double[] weight, int number_of_iterations, int M, int N, double[] f_out) -
poisson_multigrid2D
private void poisson_multigrid2D(double[] f, int level, double[] rhs, double[] weight, int n1, int n2, int nm, double[] f_out, int M, int N, int[] directly_solved) -
solve_directly2D
private void solve_directly2D(double[] lhs, double[] rhs, double[] f_out, int M, int N) -
gauss_seidel2D
private void gauss_seidel2D(double[] f, double[] A, double[] d, int M, int N) -
compute_residual2D
private void compute_residual2D(double[] r, double[] A, double[] d, double[] f, int M, int N) -
downsample2D
private void downsample2D(double[] rhs, int M, int N, double[] rhs_coarse, int Mhalf, int Nhalf, double[] weight, double[] coarse_weight) -
galerkin2D
private void galerkin2D(int level, int M, int N, int Mhalf, int Nhalf, double[] weight, double[] coarse_weight) -
upsample2D
private void upsample2D(int M, int N, double[] v, int Mhalf, int Nhalf, double[] f_out, double[] weight, double[] coarse_weight) -
antigradient3D
private void antigradient3D() -
poisson_full_multigrid3D
private void poisson_full_multigrid3D(int level, double[] rhs, double[] weight, int number_of_iterations, int M, int N, int P, double[] f_out) -
poisson_multigrid3D
private void poisson_multigrid3D(int level, double[] rhs, double[] weight, int n1, int n2, int nm, double[] f_out, int M, int N, int P, int[] directly_solved) -
gauss_seidel3D
private void gauss_seidel3D(double[] f, double[] A, double[] d, int M, int N, int P) -
compute_residual3D
private void compute_residual3D(double[] r, double[] A, double[] d, double[] f, int M, int N, int P) -
solve_directly3D
private void solve_directly3D(double[] lhs, double[] rhs, double[] f_out, int M, int N, int P) -
downsample3D
private void downsample3D(double[] rhs, int M, int N, int P, double[] rhs_coarse, int Mhalf, int Nhalf, int Phalf, double[] weight, double[] coarse_weight) -
galerkin3D
private void galerkin3D(int level, int M, int N, int P, int Mhalf, int Nhalf, int Phalf, double[] weight, double[] coarse_weight) -
upsample3D
private void upsample3D(int M, int N, int P, double[] v, int Mhalf, int Nhalf, int Phalf, double[] f_out, double[] weight, double[] coarse_weight)
-