/*plankf evaluate Planck function for any cell at current electron temperature */
#include "cddefines.h"
#include "physconst.h"
#include "rfield.h"
#include "plankf.h"

double plankf(long int ip)
{
	double plankf_v;

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


	/* evaluate Planck function
	 * argument is pointer to cell energy in ANU
	 * return photon flux for cell IP */
	if( rfield.ContBoltz[ip] <= 0. )
	{
		plankf_v = 1e-35;
	}
	else
	{
		plankf_v = 6.991e-21*POW2(FR1RYD*rfield.anu[ip])/
			(1./rfield.ContBoltz[ip] - 1.)*FR1RYD*4.;
	}


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