//-----------------------------------------------------------------------------
	//		applet calculates Figure3 of Basin inversion paper
	//
	//		coded by Mike Sandiford, December 7-13, 1998
	//
	//-----------------------------------------------------------------------------

	import java.awt.*;
	import java.applet.*;
	import Jlab.*;					//my own package of matlab like procedures
 	import java.util.*;

	public class inversion extends java.applet.Applet{

		Plot plot1 ;
		Plot plot2 ;

     	String  s[]  = {"H(s)","H(c)","k(s)","k(c)","h(r)", "q(c)","q(m)","rho(s)","rho(c)" };	//  labels of text fields in order that they will appear								
		double enteredValues[] =  { 1.3, 2.5, 3., 3.,12., 30, 30., 2350., 2750.};				//  values of text fields in order that they will appear								
		double minEnteredValues[] =  { 0, 1, 1, 1,.5, 0, 0, 2000., 2550.};						//  values of text fields in order that they will appear								
		double maxEnteredValues[] =  {8, 12, 5, 6, 34, 120, 80, 3000., 3100.};					//  values of text fields in order that they will appear								
		double factors[] =  { 1e-6,1e-6,1,1,1e3,1e-3,1e-3,1,1};									//  conversion factors that relate entered values to real values 
	 	double  Hs, Hc,ks, kc, hr, qc, qm, rc, rs;
		PanelObject panelObject;
	 	Jlab values = new Jlab(enteredValues).multiply(factors);
	    Color colorArray[] ={Color.black,Color.blue, Color.green, Color.red};

 		double delz= 100;
		int zlength = (int) (40*1000/delz)  +1;
		Jlab z = new Jlab(zlength).multiply(delz);
	
	
		Jlab kmz  = new Jlab(zlength).negative().divide(10);
	  	Jlab 	temp0  , temp4 , temp8 , temp12 ;
	
		//-----------------------------------------------------------------------------
		//				init method
		//-----------------------------------------------------------------------------

		public void init() {

 			setBackground(Color.lightGray);
			panelObject = new PanelObject(this, s,  enteredValues,3);
			
			}

		//-----------------------------------------------------------------------------
		//				event handler
		//-----------------------------------------------------------------------------

		public boolean action(Event e, Object arg) {

	 	  	panelObject.update(enteredValues,minEnteredValues ,maxEnteredValues) ;
	 	
	 		// **** define panel dependencies here
	 		if (e.target == panelObject.texts[1])  {
				enteredValues[4] = enteredValues[5] / enteredValues[1] ;
				panelObject.texts[4].setText(String.valueOf((double) (int)(enteredValues[4] *1e5 ) /1e5  ));
				}	
			else if (e.target == panelObject.texts[4])  {
				enteredValues[1] = enteredValues[5] / enteredValues[4] ;
				panelObject.texts[1].setText(String.valueOf((double) (int)(enteredValues[1] *1e3 ) /1e3  ));
		 		}
			else if (e.target == panelObject.texts[5])  {
				enteredValues[1] = enteredValues[5] / enteredValues[4] ;
				panelObject.texts[1].setText(String.valueOf((double) (int)(enteredValues[1] *1e3 ) /1e3 ));
				}
	 		 if (enteredValues[4]> 34) {
				//	 		enteredValues[4] = 34;
				//	 		texts[4].setText(String.valueOf((double) (int)(enteredValues[4] *1e5 ) /1e5  ));	 		
				enteredValues[1] = enteredValues[5] / enteredValues[4] ;
				panelObject.texts[1].setText(String.valueOf((double) (int)(enteredValues[1] *1e3 ) /1e3 ));
		 		}
		 	//  *** end panel dependencies  
	  		repaint();	 	
		 	return super.action(e, arg);
			}

		//-----------------------------------------------------------------------------
		//				calculation routines
		//-----------------------------------------------------------------------------
	
		private double[]  calculateTemps( Jlab z,  double zs,  double delz  ) {
					 
	 		double zci = 35.*1000;	 	
			double rp = (rc-rs)/(3350-rc);		
		    double temp[] = new double[z.length];
	 	    double 	b = zci/(zci - (1 + rp)*zs);
		    int ns = (int) (zs/delz)+1;
	        int  nh1 = (int)  (hr/delz)+1;
	        int nh = (int) (nh1/b)+ns;
	  		int zc = (int) (((zci/delz)+1)/b) + ns;	 
			for(int i=0; i<=ns; i++) {			
				temp[i] = -(Hs*z.index[i]*z.index[i])/(2*ks) - (z.index[i]*(-(Hc*hr) - qm*b - Hs*zs*b))/ (ks*b);
				}
			for(int i=ns+1; i<= nh; i++) {
				temp[i]  = -(Hc*z.index[i]*z.index[i])/(2*kc) - (z.index[i]*(-(Hc*hr) - qm*b - Hc*zs*b))/ (kc*b) - (-2*Hc*hr*kc*zs + 2*Hc*hr*ks*zs - 2*kc*qm*zs*b + 2*ks*qm*zs*b - Hs*kc*zs*zs*b + Hc*ks*zs*zs*b)/(2*kc*ks*b);
				}
			for(int i=nh+1; i<z.length; i++) {			
				temp[i] = ((2*ks*qm*(z.index[i] - zs) + kc*zs*(2*qm + Hs*zs))*b*b + Hc*hr*(hr*ks + 2*kc*zs*b))/(2*kc*ks*b*b);
				}
			return temp;		
			}

		private int  calculateMoho(   double zs,   double delz  ) {

	   		values =  new Jlab(enteredValues).multiply(factors);
	 	 	rs 	=   values.index[7];
	 		rc 	=   values.index[8];
			double zci = 35.*1000;	 
			double rp = (rc-rs)/(3350-rc);
	 	    double 	b = zci/(zci - (1 + rp)*zs);
		    int ns = (int) (zs/delz)+1;  
		    int zc = (int) (Math.round(zci/delz/b) + ns);    
	        return  zc-1;
			}
	
		//-----------------------------------------------------------------------------
		//				plot methods does  2-D plot ->  uses Jlab
		//-----------------------------------------------------------------------------
	 	
		public void makePlot1 () { 
			
			Jlab tempa = temp0.join(temp4,temp8, temp12);
			int moho0 	= calculateMoho(     0,    delz ) ;
			int moho4 	= calculateMoho(  4000,    delz ) ;
			int moho8 	= calculateMoho(  8000,    delz ) ;
			int moho12 	= calculateMoho( 12000,    delz ) ;
		    double yint[]  = {kmz.index[moho0],kmz.index[moho4],kmz.index[moho8], kmz.index[moho12]} ;
			double xint[]  = {temp0.index[moho0],temp4.index[moho4],temp8.index[moho8],temp12.index[moho12]};

			plot1 = new Plot(this); 
			plot1.subplot(2,1,1 );
		  	plot1.options.put("Color", colorArray);
   	 	  	plot1.plotYX(tempa, kmz  );
	     	plot1.point( xint,yint,2. )	;
	     	plot1.text( "T(zc) = " + (int)temp0.index[moho0], 	 .45 ,	.9,   	colorArray[0]);
	      	plot1.text( "T(zc) = " + (int) temp4.index[moho4], 	 .45 ,	.83,   	colorArray[1]);
 	     	plot1.text( "T(zc) = " + (int) temp8.index[moho8], 	 .05 ,	.13,   	colorArray[2]);
	   		plot1.text( "T(zc) = " + (int) temp12.index[moho12],  .05 ,	.05,  	colorArray[3]);
	     	plot1.xlabel("T   (C)" );
	     	plot1.ylabel("Z" );
			}

		public void makePlot2 () { 
			Jlab tempb = temp0.subtract(temp0).join(temp4.subtract(temp0),temp8.subtract(temp0),temp12.subtract(temp0));
			plot2 = new Plot(this); 
			plot2.subplot(2,1,2);		
	 	    plot2.options.put("Color", colorArray);
			plot2.plotYX(tempb, kmz ) ;
   	   	    plot2.xlabel("del T (C) "   );	
   	   	    plot1.title("geotherms");     
   	   	    plot2.title("delta T");    
   	   		}
		 	  	
		//-----------------------------------------------------------------------------
		//				paint method does  2-D plot ->  uses Jlab
		//-----------------------------------------------------------------------------
		 	
		public void paint(Graphics g){
		        
	   		values 	=  	new Jlab(enteredValues).multiply(factors);
	 		Hs 		= 	values.index[0];
	 		Hc 		= 	values.index[1];
	 		ks 		= 	values.index[2]; 
	 		kc 		= 	values.index[3];
	 		hr 		=   values.index[4];
	 		qc 		=   values.index[5];
	 		qm 		= 	values.index[6];
	 	 	rs 		=   values.index[7];
	 		rc 		=   values.index[8];
	
			temp0 	= 	new Jlab(calculateTemps( z,  0,    	 delz ));
			temp4 	= 	new Jlab(calculateTemps( z,  4000,   delz ));
			temp8 	= 	new Jlab(calculateTemps( z,  8000,   delz ));
			temp12 	=  	new Jlab(calculateTemps( z,  12000,  delz ));
			
			makePlot1();
			makePlot2();
 	  	    plot1.render(g,this); 	
 		  	plot2.render(g, this);
 	   	}
 	 	   	
 	  	    
		}