/* 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 */
/* rfield.h */

/* wavelength of V filter in Angstroms */
#define WL_V_FILT	5500.

/* wavelength of B filter in Angstroms */
#define WL_B_FILT	4400.

/* zero out rfield arrays between certain limits, code in zero.c */
void rfield_opac_zero( long lo , long ihi );

/* set true when malloced, init to false */
extern int lgRfieldMalloced;

EXTERN struct t_rfield {

	/* ================================================================================= */
	/* the following define the continuum energy scale and its limits */

	/* nflux is number of continuum points needed to get to high energy
	 * end of this continuum.  this is initally set to nupper in ContCreatePointers,
	 * and then trimed down in ContSetIntensity.  It is finally reset in ConvInitSolution,
	 * to make sure continuum includes all possible line and continuum emission
	 * produced by the code */
	long int nflux;

	/* number of frequency cells needed to get full energy resolution of code */
	long int nupper;

	/* the energy of the lower limit of the continuum */
	float emm;

	/* the energy of the upper limit of the continuum */
	float egamry;

	/*faintest high energy flux to consider, set with set flxfnt command */
	float FluxFaint;

	/* energy in Ryd of center of cell */
	float *anu/*[NC_ELL]*/;

	/* original set of energy units, always the same no matter what the
	 * shape of the continuum - actual energy units change with shape */
	float *AnuOrg/*[NC_ELL]*/;

	/* width of cells in ryd */
	float *widflx/*[NC_ELL]*/;

	/* used to keep track of number of lines per freq interval */
	long int *line_count;

	/* brems occupation number - only used for two photon */
	float *OccNumbBremsCont;

	/* outward emitted continuum */
	float *OccNumbContEmitOut;

	/* these are log, sqrt, square, and cube of anu array */
	float *anulog/*[NC_ELL]*/, 
	  *anusqr/*[NC_ELL]*/, 
	  *anu2/*[NC_ELL]*/, 
	  *anu3/*[NC_ELL]*/;

	/* ================================================================================= */
	/* the following are the arrays containing the local radiation field */

	/*flux is photons per cell N.B. width of cells vary with energy, given by widflx */
	float *flux/*[NC_ELL]*/;

	/* the accumulated flux, sum from this energy to infinity */
	float *flux_accum/*[NC_ELL]*/;

	/* option to not do line transfer, set false with no line transfer command */
	int lgDoLineTrans;

	/* says whether to constantly reevaluate opacities, normally true,
	 * set false with no opacity reevaluate command */
	int lgOpacityReevaluate;

	/* this flag says that CMB has been set */
	int lgCMB_set;

	/* says whether to constantly reevaluate ionization, normally true,
	 * set false with no ionization reevaluate command */
	int lgIonizReevaluate;

	/* convoc is the conversion factor from rfield to OccNumbIncidCont */
	float *convoc/*[NC_ELL]*/;

	/* OccNumbIncidCont is the continuum occupation number, attenuated incident ONLY */
	float *OccNumbIncidCont/*[NC_ELL]*/;

	/* OccNumbDiffCont is the continuum occupation number, for local diffuse continuum */
	float *OccNumbDiffCont/*[NC_ELL]*/;

	/* array of Boltzmann factors for the continuum energy grid and current temperature */
	/* >>>chng 00 apr 04 from float to double so that corr for stim emis H- not go negative */
	double *ContBoltz/*[NC_ELL]*/;

	/* ConEmitLocal is local diffuse continuum, per unit vol,
	 * for current zone, evaluated in RT_diffuse */
	float *ConEmitLocal/*[NC_ELL]*/;

	/* reflected diffuse emission continuum */
	float *ConEmitReflec/*[NC_ELL]*/;

	/* outward diffuse emission continuum (not the interactive one),
	 * this is incremented in radinc because of interplay between absorption
	 * and emission - get the outward brems right */
	float *ConEmitOut/*[NC_ELL]*/;

	/* this is set in RT_diffuse and carries interactive continua */
	float *ConInterOut/*[NC_ELL]*/;

	/* ConRefIncid is reflected portion of incident continuum */
	float *ConRefIncid/*[NC_ELL]*/;

	/* these are energy-by-energy sums of various arrays, used to save time in
	 * evaluating rate integrals */
	float *SummedCon/*[NC_ELL]*/;
	float *SummedDif/*[NC_ELL]*/;
	float *SummedOcc/*[NC_ELL]*/;
	float *SummedDifSave/*[NC_ELL]*/;

	/* this will contail array of locally destroyed continuum photons,
	 * zeroed and evaluated in RT_diffuse, currently only two-photon */
	float *ConOTS_local_photons/*[NC_ELL]*/,
		/* the local photoionization rate corresponding to above photons */
		*ConOTS_local_OTS_rate/*[NC_ELL]*/;

	/* saves total two photon continuum for debugging, set in RT_diffuse */
	float *TotDiff2Pht;

	/* otsline and ots con - local ots fields for line and continua
	 * outlin outward line fields */
	float 
		/* the local ots line rates */
		*otslin/*[NC_ELL]*/, 
		/* the local ots continuum rates */
		*otscon/*[NC_ELL]*/, 

		/* these two are the new rates, that will be entered into the above once
		 * RTOTSUpdata has been called */
		 *otslinNew/*[NC_ELL]*/,
		 *otsconNew/*[NC_ELL]*/,

		*outlin/*[NC_ELL]*/, 
		*outlin_noplot/*[NC_ELL]*/, 
		**otssav/*[NC_ELL][2]*/;

	/* reflected line */
	float *reflin/*[NC_ELL]*/;

	/* save incident continuum for later iterations */
	float *FluxSave/*[NC_ELL]*/;

	/* this is zero or one depending whether pumping by diffuse fields
	 * is turned on (1) or turned off (0) */
	float DiffPumpOn;

	/* string identifing the first line that occured at this energy */
	char **chLineLabel/*[NC_ELL][5]*/;

	/* string identifing the first continuum edge that occured at this energy */
	char **chContLabel/*[NC_ELL][5]*/;

	/* free free gaunt factor for all charges */
	/* First index is EFFECTIVE CHARGE!	*/
	float **gff/*[LIMELM][NC_ELL]*/;

	/* flag which, if set to false, causes gauntff.dat, if it exists, to be read in,
	 * or if true, causes the file gauntff.dat to be created. rfield.gff[][] is
	 * filled by interpolation on the values in this file	*/
	int lgCompileGauntFF;

	/* method for transferring diffuse continuum
	 * either 'ots' or 'oux' where x is n for which */
	char chDffTrns[4];

	/* another flag, true if outward only - 
	 * used to multipy the ConInterOut continuum for creating the interactive
	 * continuum - when false, no outward only, this is not added */
	int lgOutOnly;

	/* ipEnergyBremsThin is index for lowest energy thin to ff abs and plasma frequency
	 * EnergyBremsThin is energy there, ryd */
	long int ipEnergyBremsThin;
	float EnergyBremsThin;

	/* index of highest cell with positive Boltzmann factor */
	long int ipMaxBolt;

	/* turn off continuum pumping, set with 'no induced processes' command */
	int lgInducProcess;

	/* saves for the upward and downward compton shifts */
	double *comup/*[NC_ELL]*/, 
	  *comdn/*[NC_ELL]*/;

	/* array indices for centers of B and V filters */
	long int ipB_filter , ipV_filter;

	/* these are the lower and upper bounds for the G0 radiation field
	 * used by Tielens & Hollenbach in their PDR work */
	long int ipG0_TH85_lo , ipG0_TH85_hi;

	/* these are the lower and upper bounds for the G0 radiation field
	 * used by Tielens & Hollenbach in their PDR work */
	long int ipG0_DB96_lo , ipG0_DB96_hi;

	/* these are the lower and upper bounds for the special G0 radiation field */
	long int ipG0_spec_lo , ipG0_spec_hi;

	/* this is the wavelength where Bertoldi & Draine estimate the Habing field */
	long int ip1000A;

	/* extinction in magnitudes at B and V filters for a point source,
	 * this does not discount forward scattering */
	double extin_mag_B_point , extin_mag_V_point;

	/* extinction in magnitudes at B and V filters for a resolved source,
	 * this does discount forward scattering, so is appropriate for a resolved source */
	double extin_mag_B_extended , extin_mag_V_extended;

	/* these are total opacities at these wavelengths, used to stop at exact Av */
	double opac_mag_B_point, opac_mag_V_point, opac_mag_B_extended , opac_mag_V_extended;

	/* coefficients for fitting Tarter expressions for compton 
	 * heating and cooling, over full energy range of code */
	float *csigh/*[NC_ELL]*/, 
	  *csigc/*[NC_ELL]*/;

	double comtot, 
	  cmheat, 
	  cmcool, 
	  comoff, 
	  cinrat;

	/* set true if compton cooling underflows */
	int lgComUndr;

	/* parameters to do with incident continuum */
	/* limit to number of spectra that can be entered */
#	define	LIMSPC	100

	/* the limit to the size of the continuum array */
#	define	NCELL	130000

	double totpow[LIMSPC], 
	  slope[LIMSPC], 
	  cutoff[LIMSPC][3], 
	  spfac[LIMSPC];

	long int nspec, 
	  ipspec;

	/* these are used for interpolate and table commands, 
	 * this is the linear energy Rydbergs of continuum point */
	/* >>chng 01 oct 16, made into explicit arrays again so that
	 * table commands will work before continuum defined.*/
	/* >>chng 05 jun 30, change from declared arrays to malloced pointers */
	float **tNuRyd /*[LIMSPC][NCELL]*/ , 
	  **tslop /*[LIMSPC][NCELL]*/, 
	  /* this is the log of f_nu of continuum point */
	  **tFluxLog /*[LIMSPC][NCELL]*/;

	/* this flag is set true if we malloced out the three previous arrays -
	 * these can be returned once the continuum is generated */
	int *lgContMalloc;

	/* range for particular continuum source */
	double range[LIMSPC][2];

	/* save array for table read command, to save continuum previously punched */
	float *ConTabRead/*[NC_ELL]*/;

	/* chSpNorm says how spectrum was normlized - are or lum
	 * chrspec says whether per unit anea or tot 4 pi  */
	char chSpNorm[LIMSPC][5], 
	  chRSpec[LIMSPC][5], 
	  chSpType[LIMSPC][6];

	/* pointer to file for table read command, only 1 possible */
	FILE *ioTableRead[LIMSPC];

	/* these are total numbers of photons over various energy ranges */
	float qhtot, 
	  qhe, 
	  qheii, 
	  qbal, 
	  qrad, 
	  qtot;

	/* hydrogen ionization parameter */
	float uh;

	/* helium ion ionization parameter */
	float uheii;

	/* lgUSphON flag set when we hit stromgren radius in spherical geometry	*/
	int lgUSphON;
	/* the stromgren radius var set to get u spherical */
	float rstrom;

	/* flag set if incident radiation field less than
	 * 10x the Habing ISM field */
	int lgHabing;

	/* >>chng 03 feb 12, add fine array to keep track of line overlap */
	/* this is the velocity width of a call in the fine opacity scale,
	 * 1e5 cm/s == 1 km/s
#	define FINE_OPAC_VELOCITY_WIDTH	(1e5) */
	/* >>chng 05 feb 27, make this a variable depending on which is the heaviest
	 * element, and the lowest temperature */
	/* heaviest element to be considered - default is iron */
	long int fine_opac_nelem;

	/* number of resolution elements over width of this element, default is 4 */
	long int fine_opac_nresolv;

	/* number of cm/s of each cell in fine mesh */
	float fine_opac_velocity_width;

	/* these are transmission coeffficients for fine opacity degraded to coarse continuum */
	float *trans_coef;
	/* and array indices for lower and upper bounds of each coarse continuum mappped onto the
	 * fine continuum 0 (false) if fine continuum does not extend to these cells */
	long int *ipnt_coarse_2_fine;

	/* low and high bounds of fine array - only must include resonance lines */
	float fine_ener_lo, fine_ener_hi;
	/* the number of cells actually used in this continuum */
	long nfine;
	/* the number of cells malloced - will be several million */
	long nfine_malloc;
	/* the resolution of the array - initally 1 km/s */
	double fine_resol;
	/* the opacity array itself */
	float *fine_opac;
	/* total optical depth array */
	float *fine_opt_depth;
	/* energies at center of each bin */
	float *fine_anu;

	/* option to turn off fine opacities with no fine opacity command */
	int lgOpacityFine;
	/* says that fine optical depths will be punched, so save them */
	int lgPunchOpacityFine;

	/* flag saying to update fine opacities - not done all the time */
	int lgFine_opac_update;

	/* flag saying whether to include H1 Lya ots - usually TRUE, set to FALSE - 0 -
	 * with no lya ots command */
	int lgLyaOTS;

	/* flag saying whether to include HeII Lya and rec cont ots - usually TRUE, set to FALSE - 0 -
	 * with no heii ots command */
	int lgHeIIOTS;

	/* flag saying whether to kill outward only lines */
	int lgKillOutLine;

	/* flag saying whether to kill outward only continuum */
	int lgKillOutCont;

	/* flag saying whether to kill ots lines */
	int lgKillOTSLine;

	/* flag to prevent more than one table read command */
	int lgTableRead;

	/* set true if plasma freq enters energy array */
	int lgPlasNu;

	/* plasma frequency for current locatino */
	float plsfrq, 

		/* store highest energy plasma frequency encountered*/
		plsfrqmax ;

	/* pointer to current plasma frequency */
	long int ipPlasma, 

		/* pointer to largest plasma freq encountered */
	  ipPlasmax;

	/* these are series of flags that say whether different parts of the
	 * continuum where entered ok */
	int lgMMok, 
	  lgHPhtOK, 
	  lgXRayOK, 
	  lgGamrOK;

	/* lowest energy counted as gamma rays, 100 keV	*/
	float EnerGammaRay;
	long int ipEnerGammaRay;

	/* lgHionRad set to .true. if no hydrogen ionizing radiation */
	int lgHionRad;

	/* these store photon occupation numbers at various energies in the continuum */
	float occmax, 
	  occmnu, 
	  tbrmax, 
	  tbrmnu, 
	  tbr4nu, 
	  occ1nu;

	/* flag saying that occupation number at 1 ryd > 1 */
	int lgOcc1Hi;

}	rfield;

