/*PunCool punch coolants */
#include <float.h>
#include <string.h>
#include <limits.h>
#include "cddefines.h"
#include "sign.h"
#include "coolants.h"
#include "tcool.h"
#include "heat.h"
#include "converge.h"
#include "zonecnt.h"
#include "wkhtcl.h"
#include "puncool.h"
/* this is limit to number of coolants to save */
#define	IPSAVE	100
/* this is limit to number of coolants to print out */
#define  IPRINT 100

void PunCool(FILE * io)
{
	char chL2[IPSAVE][5], 
	  chLSav[IPSAVE][5];

	long int i, 
	  imax, 
	  ip, 
	  is, 
	  isav[IPSAVE], 
	  isav2[IPSAVE], 
	  j;

	double csav[IPSAVE], 
	  csav2[IPSAVE], 
	  cset, 
	  sgn2[IPSAVE], 
	  sgnsav[IPSAVE], 
	  xmax;

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

	for( i=0; i<IPSAVE; ++i )
	{
		csav[i] = -DBL_MAX;
		strcpy( chLSav[i] ,  "crap" );
		isav[i] = INT_MIN;
		sgnsav[i] = -DBL_MAX;
	}

	/* cset will be weakest cooling to consider
	 * WeakHeatCool set with 'set weakheatcool' command
	 * default is 0.05 */
	cset = tcool.ctot*WkHtCl.WeakHeatCool;

	/* first find all strong lines, both + and - sign */
	ip = 0;
	for( i=0; i < Coolants.ncltot; i++ )
	{
		/* only remember the most important heating/ cooling agents */
		if( (Coolants.cooling[i] > cset || Coolants.heatnt[i] > cset) )
		{
			strcpy( chLSav[ip], chCoolants.chClntLab[i] );
			isav[ip] = Coolants.collam[i];
			csav[ip] = MAX2(Coolants.cooling[i],Coolants.heatnt[i])/
			  tcool.ctot;

			/* save sign to remember if heating or cooling line */
			if( Coolants.heatnt[i] == 0. )
			{
				sgnsav[ip] = 1.;
			}
			else
			{
				sgnsav[ip] = -1.;
			}

			/* this will now be one too many, which is ok since following
			 * loop is < ip */
			++ip;

			/* IPSAVE is dimension of save arrays, do not exceed this */
			if( ip >= IPSAVE ) break;
		}
	}

	/* will be pointer to largest value */
	imax = 0;

	/* order strongest to weakest */
	for( i=0; i < ip; i++ )
	{
		xmax = 0.;
		for( j=0; j < ip; j++ )
		{
			if( csav[j] > xmax )
			{
				xmax = csav[j];
				imax = j;
			}
		}
		strcpy( chL2[i], chLSav[imax] );
		isav2[i] = isav[imax];
		csav2[i] = csav[imax];
		sgn2[i] = sgnsav[imax];
		csav[imax] = -1.;
	}

	/* warn if tcovergence failure occurred */
	if( !conv.lgConvTemp )
	{
		fprintf( io, " >>>>  Temperature notconverged.\n" );
	}
	else if( !conv.lgConvEden )
	{
		fprintf( io, " >>>>  Electron density not convergenced.\n" );
	}
	else if( !conv.lgConvIoniz )
	{
		fprintf( io, " >>>>  Ionization not converged.\n" );
	}
	else if( !conv.lgConvPres )
	{
		fprintf( io, " >>>>  Pressure not converged.\n" );
	}

	/* convert very long wavelengths to micros */
	for( i=0; i < ip; i++ )
	{
		if( isav2[i] > 100000 )
			isav2[i] /= 10000;
	}

	/* start the printout with zone number, total heating and total cooling */
	fprintf( io, "%4ld%10.2e%10.2e", ZoneCnt.nzone, heat.htot, 
	  tcool.ctot );

	fprintf( io, " " );

	/* print only up to IPRINT, which is defined above */
	ip = MIN2( ip , IPRINT );

	/* now print the coolants 
	 * keep sign of coolant, for strong negative cooling 
	 * order is ion, wavelength, fraction of total */
	for( is=0; is < ip; is++ )
	{
		fprintf( io, "%4.4s%6ld%6.3f ", chL2[is], isav2[is], 
		  sign(csav2[is],sgn2[is]) );
	}
	fprintf( io, " \n" );

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