content: /*
   * Copyright (C) 2010-2018 Mathias Unberath
   * CONRAD is developed as an Open Source project under the GNU General Public License (GPL).
  */
  package edu.stanford.rsl.conrad.angio.preprocessing.segmentation.hessian.tools;
  
  import ij.IJ;
  import ij.ImagePlus;
  import ij.ImageStack;
  import ij.measure.Calibration;
  import ij.process.FloatProcessor;
  
  public abstract class HessianEvalueProcessor implements GaussianGenerationCallback {
  
  	public void proportionDone(double d) {
  		IJ.showProgress(d);
  	}
  
  	public abstract float measureFromEvalues2D( float [] evalues );
  	public abstract float measureFromEvalues3D( float [] evalues );
  
  	protected boolean normalize = false;
  	protected double sigma = 1.0;
  	protected boolean useCalibration = false;
  
  	public void setSigma( double newSigma ) {
  		sigma = newSigma;
  	}
  
  	public ImagePlus generateImage(ImagePlus original) {
  
  		Calibration calibration=original.getCalibration();
  
  		float sepX = 1, sepY = 1, sepZ = 1;
  		if( useCalibration && (calibration!=null) ) {
  			sepX = (float)calibration.pixelWidth;
  			sepY = (float)calibration.pixelHeight;
  			sepZ = (float)calibration.pixelDepth;
  		}
  
  		double minimumSeparation = Math.min(sepX,
  						    Math.min(sepY,sepX));
  
  		ComputeCurvatures c = new ComputeCurvatures(original, sigma, this, useCalibration);
  		IJ.showStatus("Convolving with Gaussian \u03C3="+sigma+" (min. pixel separation: "+minimumSeparation+")...");
  		c.run();
  
  		int width = original.getWidth();
  		int height = original.getHeight();
  		int depth = original.getStackSize();
  
  		ImageStack stack = new ImageStack(width, height);
  
  		float[] evalues = new float[3];
  
  		IJ.showStatus("Calculating Hessian eigenvalues at each point...");
  
  		float minResult = Float.MAX_VALUE;
  		float maxResult = Float.MIN_VALUE;
  
  		if( depth == 1 ) {
  
  			float[] slice = new float[width * height];
  
  			for (int y = 1; y < height - 1; ++y) {
  				for (int x = 1; x < width - 1; ++x) {
  
  					boolean real = c.hessianEigenvaluesAtPoint2D(x, y,
  										     true, // order absolute
  										     evalues,
  										     normalize,
  										     false,
  										     sepX,
  										     sepY);
  					int index = y * width + x;
  					float value = 0;
  					if( real )
  						value = measureFromEvalues2D(evalues);
  					slice[index] = value;
  					if( value < minResult )
  						minResult = value;
  					if( value > maxResult )
  						maxResult = value;
  
  				}
  				IJ.showProgress(1 / (double) height);
  			}
  
  			FloatProcessor fp = new FloatProcessor(width, height);
  			fp.setPixels(slice);
  			stack.addSlice(null, fp);
  
  		} else {
  
  			for (int z = 0; z < depth; ++z) {
  
  				float[] slice = new float[width * height];
  
  				if ((z >= 1) && (z < depth - 1)) {
  					for (int y = 1; y < height - 1; ++y) {
  						for (int x = 1; x < width - 1; ++x) {
  
  							boolean real = c.hessianEigenvaluesAtPoint3D(x, y, z,
  												     true, // order absolute
  												     evalues,
  												     normalize,
  												     false,
  												     sepX,
  												     sepY,
  												     sepZ);
  							int index = y * width + x;
  							float value = 0;
  							if( real )
  								value = measureFromEvalues3D(evalues);
  							slice[index] = value;
  							if( value < minResult )
  								minResult = value;
  							if( value > maxResult )
  								maxResult = value;
  						}
  					}
  				}
  
  				FloatProcessor fp = new FloatProcessor(width, height);
  				fp.setPixels(slice);
  				stack.addSlice(null, fp);
  				IJ.showProgress(z / (double) depth);
  			}
  		}
  
  		IJ.showProgress(1.0);
  
  		ImagePlus result=new ImagePlus("processed " + original.getTitle(), stack);
  		result.setCalibration(calibration);
  
  
  		result.getProcessor().setMinAndMax(minResult,maxResult);
  		result.updateAndDraw();
  
  		return result;
  	}
  }