/*tav computes average of old and new optical depths for new scale at end of iter,
 * called by update, also pop371 */
#include <string.h>
#include "cddefines.h"
#include "taulines.h"
#include "itercnt.h"
#include "taumin.h"
#include "sphere.h"
#include "double.h"
#include "prtmaser.h"
#include "chlinelbl.h"
#include "tav.h"

void tav(EmLine * t, 
  double WeightNew)
{
	char chSave[11];

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

	/* option to print masing lines, set with print maser */
	if( PrtMaser.lgPrtMaser && 
		( t->TauTot < -0.01 || t->TauIn <  -0.01) )
	{
		strcpy( chSave, chLineLbl(t) );
		fprintf( ioQQQ, " Masing line:%10.10s t(in, out)=%10.2e%10.2e\n", 
		  chSave, t->TauIn, t->TauTot );
	}

	if( IterCnt.iter == 1 && (!sphere.lgStatic) )
	{
		/* first shot at model, set T(ipLnTauTot) to T(1)
		 * DoubleTau normally 1, set to 2 by DoubleTau command in order
		 * to simulate two-sided photoionization */
		t->TauTot = t->TauIn*double_.DoubleTau;
		t->TauIn = (float)MIN2(tauminCom.taumin,t->TauTot/2.);
	}
	else if( sphere.lgSphere && sphere.lgStatic )
	{
		/* static sphere, both sides interact */
		if( IterCnt.iter == 1 )
		{
			/* on first iteration TauIn was 0 for first zone */
			t->TauTot = (float)(2.*t->TauIn);
		}
		else
		{
			t->TauTot = (float)(t->TauIn*WeightNew + t->TauTot*
			  (1. - WeightNew));
		}
		t->TauIn = (float)(t->TauTot/2.);
	}
	else
	{
		/* take some sort of mean of old and new limiting optical depths */
		if( t->TauIn > 0. && t->TauTot > 0. )
		{
			t->TauTot = (float)pow(10.,log10(t->TauTot)*(1. - WeightNew) + 
				log10(t->TauIn*double_.DoubleTau)*WeightNew);
			t->TauIn = (float)MIN2(tauminCom.taumin,t->TauTot/ 2.);
		}
		else
		{
			t->TauTot = (float)((t->TauIn*double_.DoubleTau + 
			  t->TauTot)/2.);
			t->TauIn = (float)MIN2(tauminCom.taumin,t->TauTot/2.);
		}
	}

	/* this is escape prob */
	t->Pesc = (float)(0.5*(1. + 1./MAX2(1.,t->TauTot)));

	/* this is fraction inward */
	t->FracInwd = (float)MIN2(1.,1.5-t->Pesc);

	/* this is destruction probability
	 * >>chng 96 sep 4, was not here, needed to add since now taking
	 * mean of old and new dest prob */
	t->Pdest = 0.;

	/* line opacity */
	t->dTau = 0.;

	/* optical depth to the continuum source */
	t->TauCon = tauminCom.taumin;

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