/* This file is part of Cloudy and is copyright (C) 1978-2004 by Gary J. Ferland.
 * For conditions of distribution and use, see copyright notice in license.txt */
/*lines_molecules put energetics, H, and He lines into line intensity stack */
#include "cddefines.h"
#include "taulines.h"
#include "physconst.h"
#include "coolheavy.h"
#include "thermal.h"
#include "dense.h"
#include "hmi.h"
#include "phycon.h"
#include "h2.h"
#include "co.h"
#include "lines_service.h"
#include "radius.h"
#include "lines.h"

void lines_molecules(void)
{
	long int i;

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

	/* molecules */
	i = StuffComment( "molecules" );
	linadd( 0., (float)i , "####", 'i');

	/* >>refer	H2	rot	Lepp, S., & Shull, J.M., 1983, ApJ, 270, 578-582
	 * roughly two microns */
	/* cooling due H2 rotation lines from simple model */
	linadd(CoolHeavy.h2line,20000.,"H2 l",'c');
	/* remember largest fraction of H2 cooling for possible comment */
	hmi.h2line_cool_frac = (float)MAX2( CoolHeavy.h2line/thermal.ctot , hmi.h2line_cool_frac );

	/* HD rotation cooling */
	linadd(CoolHeavy.HD,0,"HDro",'c');

	/* molecular hydrogen heating */
	hmi.h2dtot += hmi.HeatH2Dish_used*(float)radius.dVeff;
	hmi.h2dfrc = (float)MAX2(hmi.HeatH2Dish_used/thermal.htot,hmi.h2dfrc);

	/* largest fraction of heating due to photo dissoc of H2+ */
	hmi.h2pmax = MAX2(hmi.h2pmax,(float)(thermal.heating[0][16]/thermal.htot));

	/* heating by H2 dissociation by Lyman continuum */
	linadd(hmi.HeatH2Dish_used,0,"H2dH",'h');

	/*remember largest fraction of heating due to H2 vib deexcitation */
	hmi.HeatH2DexcMax = MAX2((float)(hmi.HeatH2Dexc_used/thermal.htot),hmi.HeatH2DexcMax);

	/*remember largest fraction of cooling due to H2 coolingn */
	hmi.CoolH2DexcMax = MAX2((float)(-hmi.HeatH2Dexc_used/thermal.htot),hmi.CoolH2DexcMax);

	/* heating by coll deexcit of vib-excited H2 */
	linadd( MAX2(0.,hmi.HeatH2Dexc_used),0,"H2vH",'h');

	/* cooling by coll deexcit of vib-excited H2 */
	linadd( MAX2(0.,-hmi.HeatH2Dexc_used) ,0,"H2vC",'c');

	/* line emission by vib-excited H2 */
	if( h2.lgH2ON )
	{
		/* when large molecule is turned on do not print this simple estimate */
		/* line emission by vib-excited H2 */
		linadd( 0. ,0,"H2 v",'i');
	}
	else
	{
		/* H2 vib-excited lines from Tielens & Hollenbach 1985 */
		linadd( hmi.Hmolec[ipMH2s]*2e-7*4.17e-12,0,"H2 v",'i');
	}

	/* add in explicit lines from the large H2 molecule */
	H2_LinesAdd();

	/* neg H ion free-bound emission, H + e -> H- + hnu */
	linadd(hmi.hmicol,0,"H-FB",'c');

	/* neg H ion free-free emission */
	linadd(CoolHeavy.brems_cool_hminus,0,"H-FF",'i');

	/* H-alpha produce by H- mutual neutralization */
	linadd(hmi.HalphaHmin*3.032e-12,6563,"H-CT",'i');

	/* remember total heating */
	hmi.hmitot += hmi.hmihet*radius.dVeff;

	/* H- heating */
	linadd(MAX2(0.,hmi.hmihet),0,"H- H",'h');

	/* induced H- cooling */
	linadd(MAX2(0.,-hmi.hmihet),0,"H-Hc",'c');

	/* H+ + H => H2+ + photon continuum cooling */
	linadd(CoolHeavy.H2PlsCool,0,"H2+ ",'c');

	/* H2+ photo dissoc heating */
	linadd(hmi.h2plus_heat,0,"H2+p",'h');

	/* HeH+ formation cooling */
	linadd(MAX2(3.27e-12+phycon.te*BOLTZMANN,0.)*dense.xIonDense[ipHYDROGEN][1]*dense.xIonDense[ipHELIUM][0]*1e-20+
	  (1.76e-11+phycon.te*BOLTZMANN)*dense.xIonDense[ipHYDROGEN][0]*dense.xIonDense[ipHELIUM][1]*1e-16,0,"HEH+",'i');

	/* carbon monoxide heating */
	co.codtot += co.CODissHeat*(float)radius.dVeff;
	co.codfrc = (float)MAX2(co.CODissHeat/thermal.htot,co.codfrc);

	/* carbon monoxide co photodissociation */
	linadd(co.CODissHeat,0,"COdh",'h');

	/* total c12o16 cooling */
	linadd(CoolHeavy.C12O16Rot,0,"CO12",'c');

	/* total c13o16 cooling */
	linadd(CoolHeavy.C13O16Rot,0,"CO13",'c');

	/* remember most amount of CO cooling */
	co.COCoolBigFrac = MAX2( co.COCoolBigFrac, (CoolHeavy.C12O16Rot+CoolHeavy.C13O16Rot)/thermal.ctot);

	/* add in CO carbon monoxide lines */
	for( i=0; i < nCORotate; i++ )
	{
		PutLine(&C12O16Rotate[i]);
	}
	for( i=0; i < nCORotate; i++ )
	{
		PutLine(&C13O16Rotate[i]);
	}

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

