// PORT A - ADC
// PORT B - Buttons
// PORT C - LCD
// PORT D - LEDs

#include <Mega32.h>
#include <stdio.h>  
#include <stdlib.h>
#asm
    .equ __lcd_port=0x15
#endasm
#include <lcd.h>
#include <delay.h>
#define begin {
#define end }

//define subroutines
void initialize(void);
void display(void); //displays reading on LCD
void adc(void); //reads 10 bits of ADC register
void regressions(void); //calculates 2nd order regression curve

//define variables
unsigned int time1, time2, time3, time4, time5; //time counters
unsigned long total1, total2, total3, total4; //summation of ADC readings
unsigned long count; //number of readings
unsigned int temp, temp1, temp2; //ADC temporary values
unsigned int cal_count; //number of calibration readings necessary
unsigned char is_cal, state, but_state, but; //state machine
char lcd_buffer[17]; 
eeprom long int dummy; //to avoid eeprom bug
eeprom float a0,a1,a2,b0,b1,b2; //coefficients for regression curve
float ntu; //calculated value from regression curve
char u,v; //calibration array indices

//regression variables
float x[32]; //ratios
float y[32]; //ntu
float x2[32], x3[32], x4[32], xy[32], x2y[32]; //intermediate regression terms
float a,b,c,d,e,f,g,h,i,j,k,l; //regression matrix elements 
char z; //misc loop index      
float sum_x, sum_y, sum_x2, sum_x3, sum_x4, sum_xy, sum_x2y;
//more intermediate regression terms

//timer interrupt
interrupt [TIM0_COMP] void timer0_compare(void)
begin      
	time1++; // 5 ms ADC sampling period
	if (time1>5) time1=0;
	
	time2++; // switch back to default screen at "wait" state
	if (time2>3000) time2=0;
	
	time3++; // flashing LED during ADC reading
	if (time3>500) time3=0;

	time4++; // debounce buttons
	if (time4>120) time4=0;

	time5++; // button hold
end


void main(void) begin

	initialize(); 

	while (1) begin

		/* 	Button Functions:

			PINB.0	Record
			PINB.1	Calibrate
			PINB.2	Cancel
			PINB.3	+
			PINB.4	-
		*/

		switch (state) begin

		//wait state
		case 0: begin
		
			//keep current display for 3 s, then default display message
			if (time2>=3000) begin
				lcd_clear();
				lcd_gotoxy(0,0);
				lcd_putsf("Ready.");
			end

			if (PINB.0==0) begin 
				state=1; //record
				total1=0; total2=0; total3=0; total4=0; //sum of turbidity samples
				count=0; //number of samples
    			end

			if (PINB.1==0) begin 
				state=2; //calibrate
			end
	
			if (PINB.2==0) begin //display regression coefficients
  				sprintf(lcd_buffer,"%3.2f",a0);
				lcd_gotoxy(0,1);
				lcd_puts(lcd_buffer);
				delay_ms(3000);
				
  				sprintf(lcd_buffer,"%3.2f",a1);
				lcd_gotoxy(0,1);
				lcd_puts(lcd_buffer);
				delay_ms(3000);

  				sprintf(lcd_buffer,"%3.2f",a2);
				lcd_gotoxy(0,1);
				lcd_puts(lcd_buffer);
				time2=0;
			end

		break; end

		//record state
		case 1: begin

			if (time3>=500) PORTD=~PORTD; //flashing LED

			//display appropriate message on LCD
			lcd_gotoxy(0,0);
			if (is_cal==0) lcd_putsf("Recording...    ");
			else lcd_putsf("Calibrating...  ");
  			
			if (time1>=5) begin

				if (PINB.2==0) begin
					lcd_gotoxy(0,0);
					lcd_putsf("Cancelled.      ");
					time2=0;
					is_cal=0; //return to wait state
					state=0;
					break;
				end

				/*	ADC samples

					A.0 - 90  reading x
					A.1 - 180 reading x
					A.2 - 90  reading x
					A.3 - 180 reading x               
				*/
				
				adc(); //read A0  
				total1+=temp; //update sample sum

				ADMUX = 0b01100001; //prepare A.1
				ADCSR.6=1; //ready next ADC read
				while(ADCSR.6==1);

				adc(); //read A1  
				total2+=temp; //update sample sum
/* unimplemented

				ADMUX = 0b01100010; //prepare A.2
				ADCSR.6=1; //ready next ADC read
				while(ADCSR.6==1);

				adc(); //read A2  
				total3+=temp; //update sample sum

				ADMUX = 0b01100011; //prepare A.3
				ADCSR.6=1; //ready next ADC read
				while(ADCSR.6==1);

				adc(); //read A3  
				total4+=temp; //update sample sum
unimplemented */
				ADMUX = 0b01100000; //prepare A.0
				ADCSR.6=1; //ready next ADC read
				while(ADCSR.6==1);

				count++; //sample count
				
				if (count>=5000) begin

					state=is_cal; //go to wait or input state
					
					//calculate 90/180 ratios 
					x[30]=(float)total1/(float)total2; //lower amps
					x[30]=(x[30]-.5)*10;
					//normalize to avoid floating point precision bug

					//x[31]=(float)total3/(float)total4; //upper amps
					
					if (is_cal==0) begin //not calibrating

						if (1 /*x[30]<525*/ ) begin //lower range
							ntu = a0*x[30]*x[30] + a1*x[30] + a2;
						end
						else begin //upper range
							ntu = b0*x[31]*x[31] + b1*x[31] + b2;
						end

						lcd_gotoxy(0,0);
						lcd_putsf("Complete.       ");
						time2=0;																		
						display();
					end
					else begin //calibrating

						if (1 /*x[30]<525*/ ) begin //lower range
							x[u]=x[30];	
						end
						else begin //upper range
							x[v]=x[31];	
						end
					end

					PORTD=0xff; //off LED

				end	//if (count>=1500)			
				
			end //if (time1>=10)
	
		break; end
		
		//calibrate state
		case 2: begin
 
			//LCD message
			lcd_gotoxy(0,0);
			lcd_putsf("Place sample.   ");

			is_cal=3; //go to input state after recording
			
			if (cal_count>=8) begin
				state=0;
				is_cal=0;
				cal_count=0;
				regressions(); //calculate regression coefficients
				lcd_gotoxy(0,0);
				lcd_putsf("Complete.       ");
				time2=0;
			end

			if (PINB.0==0) begin //record
				state=1;
				total1=0; total2=0; total3=0; total4=0; //sum of turbidity samples
				count=0; //number of samples
				cal_count++;
			end

			if (PINB.2==0) begin //cancel
				lcd_gotoxy(0,0);
				lcd_putsf("Cancelled.      ");
				time2=0;				
				state=0;
				is_cal=0;
			end

		break; end

		//input state
		case 3: begin
		
			lcd_gotoxy(0,0);
			lcd_putsf("Adjust + place. ");

			while (~PINB.0==0) begin //while record button not pushed
			
				while (time4>=120) begin //debounce period

				switch(but_state) begin
		
					case 0: begin //’Nothing’
						if (~PINB==0) but_state=0;		
						else begin
		    					but_state=1;
		    					but=~PINB;
		    					time5=0;
						end
					break; end

					case 1: begin //’Last Push’
						if (~PINB==but) but_state=2;		
						else begin
		    					but_state=0;
		    					time5=0;
		    				end
					break; end
		
					case 2: begin //’Button Hold’
						if (time5>300) begin //held button
							but_state=0;
							if (but==0b00001000) begin
								ntu+=10;
							end
							else if (but==0b00010000) begin
								ntu-=10;
							end
							display();
							time5=0; but=~PINB;
						end
						else begin //tapped button
	    						but_state=1;
							if (but==0b00001000) begin
								ntu+=1;
							end
							else if (but==0b00010000) begin
								ntu-=1;
							end	    					
	    						display();
						end
					break; end		

				end //switch

				end //while

			end //while

			lcd_clear();
			state=2;
			
			//add user input into recorded data
			if (1 /*x[30]<525*/ ) begin //lower range
				y[u]=ntu; u++;
			end
			else begin //upper range
				y[v]=ntu; v++;	
			end			
			
		break; end

		end //switch
	
	end //while
	
end //main

void display(void)
begin

	//"xxxx NTU xx.xxxR"

  	sprintf(lcd_buffer,"%4.0f",ntu);
	lcd_gotoxy(0,1);
	lcd_puts(lcd_buffer);

	lcd_gotoxy(5,1);
	lcd_putsf("NTU");   

  	sprintf(lcd_buffer,"%2.3f",x[30]);	
	lcd_gotoxy(9,1);
	lcd_puts(lcd_buffer);

	lcd_gotoxy(15,1);
	lcd_putsf("R");   
	time2=0;
end

void adc(void)
begin
	temp1=ADCL; temp2=ADCH; //read all 10 bits of ADC register
	temp=(temp2<<2)+(temp1>>6); //concatenate ADCL and ADCH
end  

void initialize(void)
begin
	//set up the ports
	DDRB=0x00; //PORT B as button input
	PORTB=0xff; //pullup
	DDRD=0xff; //PORT D as LED output
	PORTD=0xff; //LED off
       
	//set up timer 0 for ms counter
	TIMSK=2; //turn on timer 0 cmp match ISR 
	OCR0=249; //set the compare reg to 250 time ticks
  	TCCR0=0b00001011; //prescalar to 64 and turn on clear-on-match
	time1=0; time2=0; time3=0; time4=0; time5=0; //time counters

	//set up LCD
	lcd_init(16);
	lcd_clear();
	lcd_gotoxy(0,0);
	lcd_putsf("Ready.");

	//init ADC
	//channel zero/ left adj /EXTERNAL Aref
	ADMUX = 0b01100000;   
	//enable ADC and set prescaler to 1/128*16MHz=125,000
	//and clear interupt enable
	//and start a conversion
	ADCSR = 0b11000111;

	//other variables
	state=0; count=0;
	total1=0; total2=0; total3=0; total4=0;
	temp=0; temp1=0; temp2=0;
	is_cal=0; cal_count=0;
	u=0; v=5;
 	ntu=90;   
       
    	//regression stuff                
	sum_x=0;
	sum_y=0;
	sum_x2=0;
	sum_x3=0;
	sum_x4=0;
	sum_xy=0;	
	sum_x2y=0; 
	
	#asm ("sei");
end  

//least mean square regression
void regressions(void)
begin    
	sum_x=0;
	sum_y=0;
	sum_x2=0;
	sum_x3=0;
	sum_x4=0;
	sum_xy=0;	
	sum_x2y=0;  

	for (z=0; z<u; z++) begin

		x2[z]=(x[z]*x[z]);
		x3[z]=(x2[z]*x[z]);
		x4[z]=(x3[z]*x[z]);
		xy[z]=(x[z]*y[z]);
		x2y[z]=(x2[z]*y[z]);  
	end

	for (z=0; z<u; z++) begin
	
		sum_x=sum_x+x[z];
		sum_y=sum_y+y[z];
		sum_x2=sum_x2+x2[z];
		sum_x3=sum_x3+x3[z];
		sum_x4=sum_x4+x4[z];
		sum_xy=sum_xy+xy[z];	
		sum_x2y=sum_x2y+x2y[z];  
	end
	
	a=(float)u;
	b=sum_x;
	c=sum_x2;
	d=sum_y;   
	
	b=b/a;
	c=c/a;
	d=d/a;
	a=a/a;
	
	e=sum_x;
	f=sum_x2;
	g=sum_x3;
	h=sum_xy;

  	f=f/e;
  	g=g/e;
  	h=h/e;
	e=e/e;
	
	i=sum_x2;
	j=sum_x3;
	k=sum_x4;
	l=sum_x2y;

	j=j/i;
	k=k/i;
	l=l/i;
	i=i/i;	
	
	/* 	State of matrix:
		1	b	c	d
		1	f	g	h
		1	j	k	l
	*/
  
	e=e-a;
	f=f-b;
	g=g-c;
	h=h-d;

	i=i-a;
	j=j-b;
	k=k-c;
	l=l-d;
    
	/* 	State of matrix:
		1	b	c	d
		0	f	g	h
		0	j	k	l
	*/

	g=g/f;
	h=h/f;
  	f=f/f;
	
	/* 	State of matrix:
		1	b	c	d
		0	1	g	h
		0	j	k	l
	*/
	
	c=c-g*b;
	d=d-h*b;
	b=b-f*b;
         
	k=k-g*j;
	l=l-h*j;
	j=j-f*j;
    
	/* 	State of matrix:
		1	0	c	d
		0	1	g	h	
		0	0	k	l
	*/
	
	l=l/k;
	k=1;

 	d=d-c*l;
	c=c-c*k;
    
	h=h-g*l;	
	g=g-g*k;     
	
	a0=l;
	a1=h;
	a2=d;   
	  
/* unimplemented

	sum_x=0;
	sum_y=0;
	sum_x2=0;
	sum_x3=0;
	sum_x4=0;
	sum_xy=0;	
	sum_x2y=0;  
	
	for (z=15; z<v; z++) begin

		x2[z]=(x[z]*x[z]);
		x3[z]=(x2[z]*x[z]);
		x4[z]=(x3[z]*x[z]);
		xy[z]=(x[z]*y[z]);
		x2y[z]=(x2[z]*y[z]);  
	end

	for (z=15; z<v; z++) begin
	
		sum_x=sum_x+x[z];
		sum_y=sum_y+y[z];
		sum_x2=sum_x2+x2[z];
		sum_x3=sum_x3+x3[z];
		sum_x4=sum_x4+x4[z];
		sum_xy=sum_xy+xy[z];	
		sum_x2y=sum_x2y+x2y[z];  
	end
	
	a=(float)(v-15);
	b=sum_x;
	c=sum_x2;
	d=sum_y;   
	
	b=b/a;
	c=c/a;
	d=d/a;
	a=a/a;
	
	e=sum_x;
	f=sum_x2;
	g=sum_x3;
	h=sum_xy;

  	f=f/e;
  	g=g/e;
  	h=h/e;
	e=e/e;
	
	i=sum_x2;
	j=sum_x3;
	k=sum_x4;
	l=sum_x2y;

	j=j/i;
	k=k/i;
	l=l/i;
	i=i/i;	
	
	/* 	State of matrix:
		1	b	c	d
		1	f	g	h
		1	j	k	l
	*/
    	/*
	e=e-a;
	f=f-b;
	g=g-c;
	h=h-d;

	i=i-a;
	j=j-b;
	k=k-c;
	l=l-d;
    	*/
	/* 	State of matrix:
		1	b	c	d
		0	f	g	h
		0	j	k	l
	*/
    	/*
	g=g/f;
	h=h/f;
  	f=f/f;
	
	/* 	State of matrix:
		1	b	c	d
		0	1	g	h
		0	j	k	l
	*/
	/*
	c=c-g*b;
	d=d-h*b;
	b=b-f*b;
         
	k=k-g*j;
	l=l-h*j;
	j=j-f*j;
    
	/* 	State of matrix:
		1	0	c	d
		0	1	g	h	
		0	0	k	l
	*/
	/*
	l=l/k;
	k=1;

 	d=d-c*l;
	c=c-c*k;
    
	h=h-g*l;	
	g=g-g*k; 
	
	b0=l;
	b1=h;
	b2=d;    
      */
end