import ij.*;
import ij.gui.*;
import java.awt.*;
import ij.plugin.filter.PlugInFilter;
import ij.process.*;

/** RGBfromRIS
  *
  * Gets Log Complex Color Image from Real-Imaginary-Saturation.
  *
  * This example adapted by pf from tutorial plugins named
  * ColorInverter & StackAverage found with the tutorial at
  * http://www.fh-hagenberg.at/mtd/depot/imaging/imagej
  * 
  */
public class RGBfromRIS_ implements PlugInFilter {

	protected ImagePlus imp; // defines an instance variable, protected means what?
	

	public int setup(String arg, ImagePlus imp) {
		if (arg.equals("about"))
			{showAbout(); return DONE;}
		this.imp=imp;
                return DOES_32+STACK_REQUIRED; // checks for stacked floats as input
	}

	public void run(ImageProcessor ip) {

		ImageStack stack=imp.getStack(); // sets up the image stack for access
		// get width, height and the region of interest
		int w = ip.getWidth();     
		int h = ip.getHeight();    
		int nslices = (int) stack.getSize(); // gets the number of stack slices
		int dimension = ip.getWidth()*ip.getHeight();	
		float[] pixels01 = (float[]) stack.getPixels(1); // gets the first slice pixel array
		float[] pixels02 = (float[]) stack.getPixels(2); // gets the 2nd slice pixel array
		float[] pixels03 = new float[dimension];
		if (nslices>2) {pixels03 = (float[]) stack.getPixels(3);}
	
		// create a new image with the same size and copy the pixels of the original image
		ImagePlus cRGB = NewImage.createRGBImage ("Log Complex Color", w, h, 1, NewImage.FILL_BLACK);
		ImageProcessor cc_ip = cRGB.getProcessor();
		int[] pixels = (int[]) cc_ip.getPixels(); // sets up the new HSB image pixel array

		double hue=0;
		double satL=0;
		double lightness;
		double satV=0;
		double brightness=0;
		double amplitude=0;
		double phase=0;
		// double hh,tr,tg,tb;
		int i, rint, gint, bint;
		double f, p, q, t;
		double r=0; double g=0; double b=0;
		// IJ.showMessageWithCancel("testwindow","look out");
		for (int j=0;j<dimension;j++) {
			// first determine amplitude[0..+Infinity] and phase[0..2Pi] from XYsatL
			amplitude = Math.sqrt(pixels01[j]*pixels01[j]+pixels02[j]*pixels02[j]);  // >=0
			phase = Math.atan2(pixels02[j],pixels01[j]); // between -Pi and Pi?
			// amplitude = pixels01[j];
			// phase = pixels02[j];
			// from these get hue[0..360], satL[0..1], and lightness[0..1] (HSL parameters)
			if(phase >= 0){ // already positive convert to degrees
				hue = 360.*phase/(2*Math.PI); // HSB all float between 0 and 1 except hue
			} else { // move to positive branch cut and degrees
				hue = 360.*(phase+2*Math.PI)/(2*Math.PI);
			}
			if (nslices>2) {
				satL = pixels03[j];
			} else {
				satL=1;
			}
			// if ((j>(dimension/2.))&&!IJ.showMessageWithCancel("check inputs","amplitude="+IJ.d2s(amplitude)+" satL="+IJ.d2s(satL)+" phase="+IJ.d2s(phase))) {break;}
			if(amplitude<1) { // if pixel value < 1, brightness increases with amplitude 
				if (amplitude == 0) {
					lightness = 0;
					// if ((j>(dimension/2.0))&&!IJ.showMessageWithCancel("amplitude=0","amplitude="+IJ.d2s(amplitude)+" satL="+IJ.d2s(satL)+" phase="+IJ.d2s(phase))) {break;}
				} else {
					lightness = 0.5/(1.0-Math.log(amplitude));
					// if ((j>(dimension/2.0))&&!IJ.showMessageWithCancel("amplitude<1","amplitude="+IJ.d2s(amplitude)+" satL="+IJ.d2s(satL)+" phase="+IJ.d2s(phase)+" lightness="+IJ.d2s(lightness))) {break;}
				}
			} else {    // if pixel value > 1, saturation decreases as amplitude goes up
				if (amplitude < Math.exp(100)) {
					lightness = 1 - (0.5)/(1+Math.log(amplitude));
					// if ((j>(dimension/2.0))&&!IJ.showMessageWithCancel("amplitude>1","amplitude="+IJ.d2s(amplitude)+" satL="+IJ.d2s(satL)+" phase="+IJ.d2s(phase))) {break;}
				} else {
					lightness = 1;
				}
			}
			
			// next get H[0..360]S[0..1]V[0..1] from HSL
			if (satL==0) {
				satV=0;
				brightness = lightness;
			} else if (lightness<0.5) { // cautious 1/2 in java becomes zero
				brightness = lightness*(1.+satL);
				satV = (2.*satL)/(1.+satL);
				// if ((j>(dimension/2.))&&!IJ.showMessageWithCancel("lightness<0.5","satL="+IJ.d2s(satL)+" lightness="+IJ.d2s(lightness)+" satV="+IJ.d2s(satV)+" bri="+IJ.d2s(brightness))) {break;}
				// IJ.showStatus("satV");
			} else { // lightness>0.5
				brightness = lightness + satL*(1 - lightness);
				satV = 2.*(1-lightness)*satL/brightness; // correct if lightness>1/2
				// if ((j>(dimension/2))&&!IJ.showMessageWithCancel("lightness>0.5","satL="+IJ.d2s(satL)+" lightness="+IJ.d2s(lightness)+" satV="+IJ.d2s(satV)+" bri="+IJ.d2s(brightness))) {break;}
				// satV = satL/brightness; // wrong
			}			
			// finally get R[0..1]G[0..1]B[0..1] from HSV
			if( satV == 0 ) { // achromatic (grey)
				r = brightness;
				g = brightness;
				b = brightness;
			} else {
				hue /= 60;			// sector 0 to 5
				i = (int) Math.floor( hue );
				f = hue - i;			// fractional part of h
				p = brightness * ( 1 - satV );
				q = brightness * ( 1 - satV * f );
				t = brightness * ( 1 - satV * ( 1 - f ) );
				switch( i ) {
				case 0: r = brightness; g = t; b = p; break;
				case 1: r = q; g = brightness; b = p; break;
				case 2:	r = p;	g = brightness;	b = t; break;
				case 3: r = p; g = q; b = brightness; break;
				case 4: r = t; g = p; b = brightness; break;
				default: r = brightness; g = p; b = q; break;
				}
			}

			rint = (int) (255.*r);
			gint = (int) (255.*g);
			bint = (int) (255.*b);
			// now encode these values into the new cColor image
		    	pixels[j] = ((rint & 0xff) << 16) + ((gint & 0xff) << 8) + (bint & 0xff);
		}
		// ip.resetMinAndMax();
		cRGB.show();
		cRGB.updateAndDraw();
	}

	void showAbout() {
		// called by setup if the string argument is "about"
		IJ.showMessage("RGB from Ampl-Phase-Saturation",
			"displays the log color image for a 2 or 3 slice floating point stack");
	}

}