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

public class AlgorithmAntigradient2 extends AlgorithmBase
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.
  • 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

      private BitSet 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

      private dataStruct data
    • testImage

      private ModelImage 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: AlgorithmBase
      Actually runs the algorithm. Implemented by inheriting algorithms.
      Specified by:
      runAlgorithm in class AlgorithmBase
    • 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)