/* This file is part of Cloudy and is copyright (C)1978-2006 by Gary J. Ferland
 * For conditions of distribution and use see copyright notice in license.txt */
/*PrtAllTau master routine controlling printout of optical depths at end of calculation */
#include "cddefines.h"
#include "taulines.h"
#include "iso.h"
#include "opacity.h"
#include "dense.h"
#include "colden.h"
#include "elementnames.h"
#include "geometry.h"
#include "prt.h"

void PrtAllTau(void)
{
	long int i, 
		n,
	  nelem;
	float fcon, 
	  flin;

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

	/* optical depths used by code are total through model,
	 * when sphere is set, this is twice optical depth through
	 * computed structure */
	if( geometry.lgSphere )
	{
		fcon = 2.;
		if( geometry.lgStatic )
		{
			flin = 2.;
		}
		else
		{
			flin = 1.;
		}
	}
	else
	{
		fcon = 1.;
		flin = 1.;
	}

	/* print out optical depths and column densities */

	/* R(1300) is Rayleigh scattering */
	fprintf( ioQQQ, "\n Contin Optical Depths: COMP:");

	fprintf( ioQQQ,PrintEfmt("%9.2e", opac.telec));
	fprintf( ioQQQ, "    H-:");
	fprintf( ioQQQ,PrintEfmt("%9.2e",opac.thmin ));

	fprintf( ioQQQ, " R(1300):");
	fprintf( ioQQQ,PrintEfmt("%9.2e", colden.colden[ipCOL_H0]*6.71e-24));

	fprintf( ioQQQ, "  H2+:");
	fprintf( ioQQQ,PrintEfmt("%9.2e", colden.colden[ipCOL_H2p]*7e-18));

	fprintf( ioQQQ, "  Pfa:");
	if( iso.numLevels_max[ipH_LIKE][ipHYDROGEN] > 5 )
	{
		PrintE82( ioQQQ , opac.TauTotalGeo[0][EmisLines[ipH_LIKE][ipHYDROGEN][5][4].ipCont-1]/fcon);
	}
	else
	{
		PrintE82( ioQQQ , 0.);
	}
	fprintf( ioQQQ, "\n" );

	fprintf( ioQQQ, "                          Pa:");
	if( iso.numLevels_max[ipH_LIKE][ipHYDROGEN] > 4 )
	{
		fprintf( ioQQQ,PrintEfmt("%9.2e", opac.TauTotalGeo[0][EmisLines[ipH_LIKE][ipHYDROGEN][4][3].ipCont-1]/fcon));
	}
	else
	{
		PrintE82( ioQQQ , 0.);
	}

	fprintf( ioQQQ, "    Ba:");
	if( iso.numLevels_max[ipH_LIKE][ipHYDROGEN] > 3 )
	{
		fprintf( ioQQQ,PrintEfmt("%9.2e", opac.TauTotalGeo[0][EmisLines[ipH_LIKE][ipHYDROGEN][3][2].ipCont-1]/fcon));
	}
	else
	{
		PrintE82( ioQQQ , 0.);
	}

	fprintf( ioQQQ, "      Hb:");
	if( iso.numLevels_max[ipH_LIKE][ipHYDROGEN] > 4 )
	{
		fprintf( ioQQQ,PrintEfmt("%9.2e", opac.TauTotalGeo[0][EmisLines[ipH_LIKE][ipHYDROGEN][4][2].ipCont-1]/fcon));
	}
	else
	{
		PrintE82( ioQQQ , 0.);
	}

	fprintf( ioQQQ, "   La:");
	fprintf( ioQQQ,PrintEfmt("%9.2e", opac.TauTotalGeo[0][EmisLines[ipH_LIKE][ipHYDROGEN][ipH2p][ipH1s].ipCont-1]/fcon));

	fprintf( ioQQQ, "     1r:");
	PrintE93( ioQQQ , opac.TauTotalGeo[0][iso.ipIsoLevNIonCon[ipH_LIKE][ipHYDROGEN][ipH1s]-1]/fcon);

	fprintf( ioQQQ, "  1.8:");
	PrintE82( ioQQQ , opac.TauTotalGeo[0][iso.ipIsoLevNIonCon[ipHE_LIKE][ipHELIUM][0]-1]/fcon);

	fprintf( ioQQQ, " 4.:");
	PrintE93( ioQQQ , opac.TauTotalGeo[0][iso.ipIsoLevNIonCon[ipH_LIKE][1][ipH1s]-1]/fcon);
	fprintf( ioQQQ, "\n");

	if( dense.lgElmtOn[ipHELIUM] )
	{
		/* print optical depths of previous iteration with old atom */

		/* these results are for new HeI as part of he iso sequence */
		fprintf( ioQQQ, "  Line Optical Depths: 10830:");
		fprintf( ioQQQ,PrintEfmt("%9.2e", (EmisLines[ipHE_LIKE][ipHELIUM][ipHe2p3P0][ipHe2s3S].TauIn+
			EmisLines[ipHE_LIKE][ipHELIUM][ipHe2p3P1][ipHe2s3S].TauIn+
			EmisLines[ipHE_LIKE][ipHELIUM][ipHe2p3P2][ipHe2s3S].TauIn)/flin));
		fprintf( ioQQQ, "  3889:");
		fprintf( ioQQQ,PrintEfmt("%9.2e", EmisLines[ipHE_LIKE][ipHELIUM][ipHe3p3P][ipHe2s3S].TauIn/flin));
		fprintf( ioQQQ, "    5876:");
		fprintf( ioQQQ,PrintEfmt("%9.2e", (EmisLines[ipHE_LIKE][ipHELIUM][ipHe3d3D][ipHe2p3P0].TauIn+
			EmisLines[ipHE_LIKE][ipHELIUM][ipHe3d3D][ipHe2p3P1].TauIn+
			EmisLines[ipHE_LIKE][ipHELIUM][ipHe3d3D][ipHe2p3P2].TauIn)/flin));
		fprintf( ioQQQ, " 7065:");
		fprintf( ioQQQ,PrintEfmt("%9.2e", (EmisLines[ipHE_LIKE][ipHELIUM][ipHe3s3S][ipHe2p3P0].TauIn+
			EmisLines[ipHE_LIKE][ipHELIUM][ipHe3s3S][ipHe2p3P1].TauIn+
			EmisLines[ipHE_LIKE][ipHELIUM][ipHe3s3S][ipHe2p3P2].TauIn)/flin));
		fprintf( ioQQQ, "  2.06m:");
		fprintf( ioQQQ,PrintEfmt("%9.2e", EmisLines[ipHE_LIKE][ipHELIUM][ipHe2p1P][ipHe2s1S].TauIn/flin));
		fprintf( ioQQQ, "  21c:");
		fprintf( ioQQQ,PrintEfmt("%9.2e", HFLines[0].TauIn/flin));
		/* the HeII 3 cm equivalent of H 21 cm */
		fprintf( ioQQQ, "\n" );
	}

	/* print optical depths of some metal lines */
	prtmet();

	/* print H-like H, He+ optical depths */
	for( nelem=ipHYDROGEN; nelem<=ipHELIUM; ++nelem )
	{
		/* helium may be turned off */
		if( dense.lgElmtOn[nelem] )
		{
#			define NUMB_PER_LINE	8
			fprintf( ioQQQ, "\n Old, new %s%2li continuum optical depths:\n",
				elementnames.chElementSym[nelem] ,
				nelem+1);
			/* absorption continuum optical depths are energies of the h-like ion continua
			 * loop over old, then new optical depths */
			for( i=1; i>=0; --i )
			{
				/* print ground, skip t2, then do 2p */
				for( n=ipH1s; n < iso.numLevels_max[ipH_LIKE][nelem]; n++ )
				{
					if( n==ipH2s )
						continue;
					if( n%NUMB_PER_LINE ==1)
						fprintf(ioQQQ,"\n");
					/* this, combined with "continue" above, ensures that we print
					 * 1 (1s), 2(tot 2), then 3 */
					fprintf( ioQQQ , "%6ld",MAX2(1,n));
					fprintf( ioQQQ,PrintEfmt("%9.2e", opac.TauAbsGeo[i][iso.ipIsoLevNIonCon[ipH_LIKE][nelem][n]-1]/fcon));
				}
				fprintf( ioQQQ, "\n" );
			}

			/* now do h-like line optical depths */
			fprintf( ioQQQ, "\n Old, new %s%2li line optical depths:\n",
				elementnames.chElementSym[nelem] ,
				nelem+1);
			/* Lya is a special case due to 2s-2p resolution - explictly print it first */
			fprintf( ioQQQ, "%3i-%2i",2, 1 );
			fprintf( ioQQQ,PrintEfmt("%9.2e", EmisLines[ipH_LIKE][nelem][ipH2p][ipH1s].TauTot/flin ));
			/* total optical depth in 3-2s and 3-2p, is total of both so 2-1 is correct for 3-2*/
			for( n=3; n < iso.numLevels_max[ipH_LIKE][nelem]; n++ )
			{
				if( n%NUMB_PER_LINE ==1)
					fprintf(ioQQQ,"\n");
				fprintf( ioQQQ, "%3ld-%2ld",n, n-1 );
				fprintf( ioQQQ,PrintEfmt("%9.2e", EmisLines[ipH_LIKE][nelem][n][n-1].TauTot/flin ));
			}
			fprintf( ioQQQ, "\n" );

			fprintf( ioQQQ, "%3i-%2i",2, 1 );
			fprintf( ioQQQ,PrintEfmt("%9.2e", EmisLines[ipH_LIKE][nelem][ipH2p][ipH1s].TauIn/flin ));
			for( n=3; n < iso.numLevels_max[ipH_LIKE][nelem]; n++ )
			{
				if( n%NUMB_PER_LINE ==1)
					fprintf(ioQQQ,"\n");
				fprintf( ioQQQ, "%3ld-%2ld",n, n-1 );
				fprintf( ioQQQ,PrintEfmt("%9.2e", EmisLines[ipH_LIKE][nelem][n][n-1].TauIn/flin ));
			}
			fprintf( ioQQQ, "\n" );
		}
	}

	/* ================================================================================ */

	/* print helium lines if helium exists */
	if( dense.lgElmtOn[ipHELIUM] )
	{
		fprintf( ioQQQ, "\n Old He Is optical depths:" );
		for( i=1; i <= 8; i++ )
		{
			fprintf( ioQQQ, "%4ld", i );
			fprintf( ioQQQ,PrintEfmt("%9.2e", opac.TauAbsGeo[1][iso.ipIsoLevNIonCon[ipHE_LIKE][ipHELIUM][i-1]-1]/fcon) );
		}
		fprintf( ioQQQ, "\n" );

		fprintf( ioQQQ, " New He Is optical depths:" );
		for( i=1; i <= 8; i++ )
		{
			fprintf( ioQQQ, "%4ld", i );
			fprintf( ioQQQ,PrintEfmt("%9.2e", opac.TauAbsGeo[0][iso.ipIsoLevNIonCon[ipHE_LIKE][ipHELIUM][i-1]-1]/fcon ));
		}
		fprintf( ioQQQ, "\n" );

		/* ================================================================================*/


		fprintf( ioQQQ, "          Old He Is Lines:" );
		fprintf( ioQQQ, "%2d-%1d",2, 1 );
		fprintf( ioQQQ,PrintEfmt("%9.2e", EmisLines[ipHE_LIKE][ipHELIUM][ipHe2p1P][ipHe1s1S].TauTot/flin ));
		fprintf( ioQQQ, "%2d-%1d",3, 2 );
		fprintf( ioQQQ,PrintEfmt("%9.2e", EmisLines[ipHE_LIKE][ipHELIUM][ipHe3s1S][ipHe2p1P].TauTot/flin ));
		fprintf( ioQQQ, "\n" );

		fprintf( ioQQQ, "          New He Is Lines:" );
		fprintf( ioQQQ, "%2d-%1d",2, 1 );
		fprintf( ioQQQ,PrintEfmt("%9.2e", EmisLines[ipHE_LIKE][ipHELIUM][ipHe2p1P][ipHe1s1S].TauIn/flin ));
		fprintf( ioQQQ, "%2d-%1d",3, 2 );
		fprintf( ioQQQ,PrintEfmt("%9.2e", EmisLines[ipHE_LIKE][ipHELIUM][ipHe3s1S][ipHe2p1P].TauIn/flin ));
		fprintf( ioQQQ, "\n" );

		/* ================================================================================*/
	}

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

