/*aver compute average of various quantities over the computed geometry */
#include <string.h>
#include "cddefines.h"
#include "radius.h"
#include "sphere.h"
#include "printe82.h"
#include "showme.h"
#include "aver.h"

void aver(
			/* 4 char + null string giving job, prin, zero, zone, doit*/
			char *chWhat, 
			/* the quantity to average, like the temperature */
			double quan, 
			/* what we weight against, like the o++ abundance */
			double weight, 
			/* 10 char + null string giving title for this average */
			char *chLabl)
{
	/* NAVER is the limit to the number than can be averaged */
#define	NAVER	20
	static char chLabavr[NAVER][11];
	long int i, 
	  ioff, 
	  j;
	static long int n;
	double raver[NAVER]={0.}, 
	  vaver[NAVER]={0.};
	static double aversv[4][NAVER];

#	ifdef DEBUG_FUN
	fputs( "<+>aver()\n", debug_fp );
#	endif


	if( strcmp(chWhat,"zero") == 0 )
	{
		/* chWhat='zero' - zero out the save array */
		for( i=0; i < NAVER; i++ )
		{
			for( j=0; j < 4; j++ )
			{
				aversv[j][i] = 0.;
			}
		}
		n = 0;
	}
	else if( strcmp(chWhat,"zone") == 0 )
	{
		n = 0;
	}
	else if( strcmp(chWhat,"doit") == 0 )
	{
		/* chWhat='doit' - enter values to average */
		if( n >= NAVER )
		{
			fprintf( ioQQQ, " Too many values entered into AVER, increase NAVER\n" );
			puts( "[Stop in aver]" );
			exit(1);
		}
		aversv[0][n] += weight*quan*radius.dReff;
		aversv[1][n] += weight*radius.dReff;
		aversv[2][n] += weight*quan*radius.dVeff;
		aversv[3][n] += weight*radius.dVeff;

		/* this is the label for this particular average, like " TeNe "*/
		strcpy( chLabavr[n], chLabl );
		n += 1;
	}
	else if( strcmp(chWhat,"prin") == 0 )
	{
		/* set things up to get answers */
		ioff = n*10/2 + 1;

		/* space out to center the label on the numbers */
		fprintf( ioQQQ,"\n");
		for( i=0; i < ioff; i++ )
		{
			/*chTitAvr[i] = ' ';*/
			fprintf( ioQQQ, " " );
		}

		/* now print header with cr */
		fprintf( ioQQQ, "Averaged Quantities\n " );

		fprintf( ioQQQ, "        " );
		for( i=0; i < n; i++ )
		{
			fprintf( ioQQQ, "%10.10s", chLabavr[i] );
		}
		fprintf( ioQQQ, " \n" );

		for( i=0; i < n; i++ )
		{
			if( aversv[1][i] > 0. )
			{
				raver[i] = aversv[0][i]/aversv[1][i];
			}
			else
			{
				raver[i] = 0.;
			}
			if( aversv[3][i] > 0. )
			{
				vaver[i] = aversv[2][i]/aversv[3][i];
			}
			else
			{
				vaver[i] = 0.;
			}
		}

		fprintf( ioQQQ, " Radius:" );
		for( i=0; i < n; i++ )
		{
			/*fprintf( ioQQQ, "%11.3e", raver[i] );*/
			fprintf( ioQQQ, " ");
			fprintf(ioQQQ,PrintEfmt("%9.2e", raver[i] ) );
		}
		fprintf( ioQQQ, "\n" );

		/* only print vol aver is lgSphere is set */
		if( sphere.lgSphere )
		{
			fprintf( ioQQQ, " Volume:" );
			for( i=0; i < n; i++ )
			{
				/*fprintf( ioQQQ, "%11.3e", vaver[i] );*/
				fprintf( ioQQQ, " ");
				fprintf(ioQQQ,PrintEfmt("%9.2e", vaver[i] ) );
			}
			fprintf( ioQQQ, "\n" );
		}
	}
	else
	{
		fprintf( ioQQQ, " AVER does not understand chWhat=%4.4s\n", 
		  chWhat );
		ShowMe();
		puts( "[Stop in aver]" );
		exit(1);
	}

#	ifdef DEBUG_FUN
	fputs( " <->aver()\n", debug_fp );
#	endif
	return;
}