/*fill define the continuum energy grid over a specified range */
#include "cddefines.h"
#include "trace.h"
#include "rfield.h"
#include "fillnu.h"
#include "fill.h"

void fill(double fenlo, 
  double fenhi, 
  double resolv, 
  long int *n0, 
  long int *ipnt)
{
	long int i, 
	  nbin;
	float widtot;
	double aaa , bbb;

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

	/* nrange is number of ranges  */
	*ipnt += 1;
	if( *ipnt > 1 && fabs(1.-fenlo/fillnu.filbnd[*ipnt-1]) > 1e-4 )
	{
		fprintf( ioQQQ, " FILL improper bounds.\n" );
		fprintf( ioQQQ, " ipnt=%3ld fenlo=%11.4e filbnd(ipnt)=%11.4e\n", 
		  *ipnt, fenlo, fillnu.filbnd[*ipnt-1] );
		puts( "[Stop in fill]" );
		exit(1);
	}
	fillnu.nrange = MAX2(fillnu.nrange,*ipnt);

	if( *ipnt > NFILL )
	{
		fprintf( ioQQQ, " The limit to the number of FILL ranges is%4ld\n", 
		  NFILL );
		puts( "[Stop in fill]" );
		exit(1);
	}

	fillnu.ifill0[*ipnt-1] = *n0 - 1;
	fillnu.filbnd[*ipnt-1] = (float)fenlo;
	fillnu.filbnd[*ipnt] = (float)fenhi;
	nbin = (long int)(log(10.)*log10(fenhi/fenlo)/resolv + 1);
	fillnu.fildel[*ipnt-1] = (float)(log10(fenhi/fenlo)/nbin);

	if( fillnu.fildel[*ipnt-1] < 0.01 )
	{
		fillnu.filres[*ipnt-1] = (float)(log(10.)*fillnu.fildel[*ipnt-1]);
	}
	else
	{
		fillnu.filres[*ipnt-1] = (float)((pow(10.,2.*fillnu.fildel[*ipnt-1]) - 
		  1.)/2./pow(10.,fillnu.fildel[*ipnt-1]));
	}
	if( *n0 + nbin > NCELL )
	{
		fprintf( ioQQQ, " Fill would need%4ld cells to get to an energy of%10.3e\n", 
		  *n0 + nbin, fenhi );
		fprintf( ioQQQ, " Increase NCELL in all parameter statements.  The current value is%6ld\n", 
		  NCELL );
		puts( "[Stop in fill]" );
		exit(1);
	}

	widtot = 0.;
	for( i=0; i < nbin; i++ )
	{
		bbb = fillnu.fildel[*ipnt-1]*((float)(i) + 0.5) ;
		aaa = pow( 10. , bbb );
		rfield.anu[i+fillnu.ifill0[*ipnt-1]] = (float)(fenlo*aaa);
		rfield.widflx[i+fillnu.ifill0[*ipnt-1]] = rfield.anu[i+fillnu.ifill0[*ipnt-1]]*
		  fillnu.filres[*ipnt-1];
		widtot += rfield.widflx[i+fillnu.ifill0[*ipnt-1]];
	}

	*n0 += nbin;
	if( trace.lgTrace && (trace.lgConBug || trace.lgPtrace) )
	{
		fprintf( ioQQQ, " FILL range%2ld from%10.3e to%10.3eR in%4ld cell; ener res=%10.3e WIDTOT=%10.3e\n", 
		  *ipnt, rfield.anu[fillnu.ifill0[*ipnt-1]] - rfield.widflx[fillnu.ifill0[*ipnt-1]]/
		  2., rfield.anu[fillnu.ifill0[*ipnt-1]+nbin-1] + rfield.widflx[fillnu.ifill0[*ipnt-1]+
		  nbin-1]/2., nbin, fillnu.filres[*ipnt-1], widtot );
		fprintf( ioQQQ, " The requested range was%10.3e%10.3e The requested resolution was%10.3e\n", 
		  fenlo, fenhi, resolv );
	}


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