/*
	This C language subroutine computes the Variance
	of its floating point input.

	The Run Time Library for the C Language.

	Gordon A. Sterling	(617) 461 - 3076
	DSP Development Tools Engineering

	Created on 05/19/90 
	Updated on 5/94 by AS

	#include <stats.h>			<Header>
	float var(float in[], int length);	<Prototype>

*/

#include "lib_glob.h"
#include "sts_glob.h"

.SEGMENT/CODE	Code_Space_Name;
.FILE 		RTL_FILENAME;

.GLOBAL	    	_var;
.GLOBAL		___varD;

_var:
___varD:	put(R7);
		R7=R7-R7, put(R11);		/* Clear summation register	*/
	    
		R2=PASS R8, dm_ptr=R4;		/* Point to input array		*/
		F12=FLOAT R2, dm_lnt=dm_0;	/* Convert to float		*/

		R2=R2-1, F0=DM(dm_ptr,dm_1);
		IF LE JUMP error_return;	/*Must have > 1 point*/
		LCNTR=R2, DO sumd UNTIL LCE;
sumd:		    F7=F7+F0, F0=DM(dm_ptr,dm_1);
		CALL (PC, ___float_divide) (DB);
		F7=F7+F0;
		F11=2.0;

		F8=FLOAT R8, dm_ptr=R4;		/* Point to input array		*/
		F12=PASS F7, F0=DM(dm_ptr,1);
		R7=R7-R7;			/* Clear summation again	*/
		LCNTR=R2, DO varianced UNTIL LCE;
		    F0=F0-F12;		    	/* Compute X(n) - mean 		*/
		    F0=F0*F0;		    	/* Compute (X(n) - mean) **2	*/
varianced:	    F7=F7+F0, F0=DM(dm_ptr,1);	/* Compute sum += (X(n)-mean)**2*/
		F0=F0-F12;		    	/* Compute X(n) - mean 		*/
	        CALL (PC, ___float_divide) (DB);
		F0=F0*F0;
		F7=F7+F0, F12=F8;		/* Compute above sum/N		*/

		F0=PASS F7, FETCH_RETURN

restore_state:	get(R11,1);
		get(R7,2);
		RETURN (DB);
		RESTORE_STACK
		RESTORE_FRAME

___varP:	put(R7);
		R7=R7-R7, put(R11);		/* Clear summation register	*/
	    
		R2=PASS R8, pm_ptr=R4;		/* Point to input array		*/
		F12=FLOAT R2, pm_lnt=pm_0;	/* Convert to float		*/

		R2=R2-1, F0=PM(pm_ptr,pm_1);
		IF LE JUMP error_return;	/*Must have > 1 point*/
		LCNTR=R2, DO sump UNTIL LCE;
sump:		    F7=F7+F0, F0=PM(pm_ptr,pm_1);
		CALL (PC, ___float_divide) (DB);
		F7=F7+F0;
		F11=2.0;

		F8=FLOAT R8, pm_ptr=R4;		/* Point to input array		*/
		F12=PASS F7, F0=PM(pm_ptr,1);
		R7=R7-R7;			/* Clear summation again	*/
		LCNTR=R2, DO variancep UNTIL LCE;
		    F0=F0-F12;		    	/* Compute X(n) - mean 		*/
		    F0=F0*F0;		    	/* Compute (X(n) - mean) **2	*/
variancep:	    F7=F7+F0, F0=PM(pm_ptr,1);	/* Compute sum += (X(n)-mean)**2*/
		F0=F0-F12;		    	/* Compute X(n) - mean 		*/
	        CALL (PC, ___float_divide) (DB);
		F0=F0*F0;
		F7=F7+F0, F12=F8;		/* Compute above sum/N		*/

		JUMP (PC, restore_state) (DB);
		FETCH_RETURN
		F0=PASS F7;

error_return:	JUMP (PC, restore_state) (DB);
		R0=0;
		FETCH_RETURN

.ENDSEG;