/* 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 */
/*
 * cddrive.h
 * this file contains definitions for the suite of routines that
 * allow the code to be driven as a subroutine.
 * These routines set up model parameters, 
 * control the execution of Cloudy, and obtain results once complete
 * these are the only "public" routines, and only these should
 * be accessed when controlling Cloudy
 */

/* 
 * DRIVING CLOUDY FROM A FORTRAN PROGRAM:
 * This should not be too hard - the recommended approach is to use
 * the cfortran.h file described at http://www-zeus.desy.de/~burow/cfortran/ 
 */

/*
 * A note on return error conditions:
 *
 * Many functions return an int to indicate success or failure.
 * I try to follow standard Unix/C convention.
 * A return value of zero usually indicates that the routine was successful,
 * and a non-zero value (not always 1) indicates failure.  This is conventional
 * in both C and Unix.  So the way to call Cloudy and test for success is
 *
 * if( cdDdrive() )
 * {
 *	   printf(" Cloudy failed.\n");
 * }
 *
 * Although I try to follow this, there ARE exceptions.
 */

/*
 * cdInit 
 * This routine must be called before any of the others - 
 * it initialialies many variables, and must be called before any
 * of the other routines.  In a large grid of calculations it must be repeatedly called
 * before the start of the new calculation and after all results have
 * been obtained from the previous model */
void cdInit(void);

/*
 * cdTalk 
 * tells the code whether or not to produce any of its normal ouput, 
 * If the argument is TRUE (1) it produces output, produces no
 * output if it is FALSE (0)
 * the default, if this routine is not called, is to produce the normal output */
void cdTalk(int);

/*
 * cdOutp 
 * This tells the code where to send output.  The argument must be
 * a previously oppened file handle.  All output will go to this file.
 * If this routine is not called then all output will go to 
 * stdout, the standard c output */
void cdOutp(FILE* );

/* cdInp - redirect input stream to open file */
void cdInp(FILE* );

/*
 * cdDepth_depth 
 * returns depth structure of previous model */
void cdDepth_depth( double cdDepth[] );

/*
 * cdnZone 
 * returns number of zones */
long int cdnZone(void );

/*
 * cdB21cm
 * returns B as measured by 21 cm 
 * assumes tangled field weighted by n(H0)/T */
double cdB21cm( void );

/*
 * cdRead 
 * This sends commands to the code.  The normal set of commands
 * described in Part I of Hazy must be entered into a null-terminated
 * string.  These strings are then fed to Cloudy with this command.  The function
 * returns the number of commands that can still be entered before the command
 * stack is full.  The code will stop if you try to continue giving it commands
 * after the command has returned zero. This return value is the opposite of the
 * standard - a non-zero return is normal */
int cdRead(char* );

/*
 * cdNoExec 
 * This provides nn option to have the code prepare the initial conditions for a model,
 * but not actually try to compute the model.  I use this when setting up a large
 * grid so that I can quickly run through the full grid as a check that the commands
 * are entered properly and the parameters I am going to vary do so properly.
 * The command is then commented out when the grid is properly set up. */
void cdNoExec(void);

/*
 * cdDrive 
 * This command actually computes a model.
 * It returns 0 if the calculation was successful, and 1 if an error
 * condition was encountered */
int cdDrive(void);


/* The next two routines confirm that the previous calculation was ok 
 * or produce a list of error conditions */

/*
 * cdErrors
 * After the calculation is completed, a summary of all error messages can be
 * be generated by calling this routine.  The argument is the output file */
void cdErrors(FILE* );

/*
 * cdNwcns
 * This command returns the number of warnings, cautions, notes, surprises, 
 * assorted types of failures found the last computed model */
void cdNwcns(
  /* abort status, if non-zero then big problems happened */
  int *lgAbort ,
  /* the number of warnings */
  long int *NumberWarnings, 
  /* the number of cautions */
  long int *NumberCautions, 
  /* the number of notes */
  long int *NumberNotes, 
  /* the number of surprises */
  long int *NumberSurprises, 
  /* the number of temperature convergence failures */
  long int *NumberTempFailures, 
  /* the number of pressure convergence failures */
  long int *NumberPresFailures,
  /* the number of ionzation convergence failures */
  long int *NumberIonFailures, 
  /* the number of electron density convergence failures */
  long int *NumberNeFailures );

/* This prints the reason why the model stopped, and the model geometry, on 
 * the io file pointed to by the file handle */
void cdReasonGeo(FILE*);

/*
 * These routines produce lists of warnings, cautions, notes, surprises
 * These are the routines that are called by cdErrors.  Normally
 * cdErrors and not these routines would be called by the user. */
void cdWarnings(FILE*);
void cdCautions(FILE*);
void cdSurprises(FILE*);
void cdNotes(FILE*);

/***********************************************************
 *
 * The next routines examine the preditions of the previous model
 *
 ***********************************************************/

/* 
 * cdLine
 * This routine finds the predicted intensity of any line in the standard output.
 *
 * arguments:
 * 1st parameter is the 4-char label + null terminated label, as it appears in the output.
 * 2nd parameter is the float wavelength in Angstroms, not how it appears in printout.
 *     The first four digits must agree with the number in the printout, but the units must be Angstroms.
 * 3rd parameter is the predicted intensity relative to the normalization line.
 * 4th par is the log of the predicted luminosity or intensity of the line (ergs).
 * 5th is pointer to relative intensity, a double that is returned
 * 6th is pointer to log of absolute intensity
 *
 * return value:
 * The routine returns an index (>0) of the arrray element within stack if it finds the line, 
 * It returns the negative of the total number of lines if it could not find the line.
 * (this is a debugging aid)
 * note that this returns a long int since there are LOTS of lines
 * this also IS NOT the standard C convention for success or failure */
long int cdLine(
	const char *chLabel, 
	float wavelength, 
	double *relint, 
	double *absint);

/*************************************************************************
 *
 *cdLine_ip get the predicted line intensity, using index for line in stack 
 *
 ************************************************************************/

void cdLine_ip(long int ipLine, 
	  /* linear intensity relative to normalization line*/
	  double *relint, 
	  /* log of luminosity or intensity of line */
	  double *absint );

/*************************************************************************
 *
 *cdDLine get the predicted emergent line intensity, also index for line in stack 
 *
 ************************************************************************/

 /* returns array index for line in array stack if we found the line, 
 * or FALSE==0 if we did not find the line */
long int cdDLine(char *chLabel, 
		/* wavelength of line as printed by code*/
	  float wavelength, 
	  /* linear intensity relative to normalization line*/
	  double *relint, 
	  /* log of luminosity or intensity of line */
	  double *absint );

/*************************************************************************
 *
 *cdLineListPunch create a file with a list of all emission lines punched, 
 *and their index within the emission line stack
 *
 ************************************************************************/

/* returns total number of lines in the list */
long int cdLineListPunch( 
	/* a file handle pointing to a file that is read for writing -
	 * the calling routine must close it */
	 FILE* io );

/*
 * cdColm 
 * This obtains the column density of a species in the previously computed model.
 * The first parameter is a 4 character + NULL terminated string which gives
 * the first 4 char of the element name as spelled by Cloudy, either upper or lower case.
 * The second parameter is the stage of ionization, 1 for atom, 2 for first ion, etc; 0 is special.
 *
 * examples: 
 * column density of atomic carbon
 * cdColm( "carb" , 1 , &col1 );
 *
 * doubly ionized helium
 * cdColm( "heli" , 3 , &col3 );
 *
 * molecular hydrogen
 * cdColm("H2  " , 0 , &col2 );
 *
 * If the ion stage is zero then the routine will check the first label
 * for the values "H2  ", "OH  ", "CO  " and "CII* *", 
 * and will return the H2, OH, CO or CII* column density in this case
 *
 * The column density is returned as the third parameter in all cases
 *
 * The function returns 0 if it found the species, 1 if it failed */
int cdColm(const char*, long, double* );

/*cdH2_colden return column density in H2, returns -1 if cannot find state,
 * header is in cddrive, source in h2.c */
double cdH2_colden( long iVib , long iRot );

/*cdCO_colden return column density in H2, negative -1 if cannot find state,
 * header is cddrive */
double cdCO_colden( long isotope , long iRot );

/*
 * cdEmis
 * This routine finds the local emissivity for any line.
 * The first argument is the 4 character (null terminated) label as it appears in the 
 * standard output. 
 * The second argument is float wavelength as it appears in the standard output.
 * The emissivity (erg /cm^3 /s) is returned as the last parameter.
 * cdEms returns the index of the line in the stack if it was found, 
 * the negative of the total number of lines in the stack if it could not find the line */
long int cdEmis(
	/* return value will be index of line within stack */
	/* 4 char null terminated string label */
	char *chLabel,
	/* line wavelength */
	float wavelength, 
	/* the vol emissivity of this line in last computed zone */
	double *emiss );


/* cdEms_ip obtain the local emissivity for a line with known index*/

void cdEmis_ip(
	/* index of the line in the stack */
	long int ipLine, 
	/* the vol emissivity of this line in last computed zone */
	double *emiss );

/*
 * cdCooling_last
 * The returns the total cooling (erg cm^-3 s^-1) for the last computed zone */
double cdCooling_last(void);

/*
 * cdHeating_last
 * returns the total heating (erg cm^-3 s^-1) for the last computed zone */
double cdHeating_last(void);

/*
 * cdPressure_last
 * This returns the pressure and its constituents for the last computed zone. */
void cdPressure_last(
	/* total pressure, all forms*/
	double *TotalPressure,			
	/* gas pressure */
	double *GasPressure,				
	/* radiation pressure */
	double *RadiationPressure);	

/*
 * cdPressure_depth
 * This returns the pressure and its constituents for the last iteration. 
 * space was allocated in the calling routine for the vectors - 
 * before calling this, cdnZone should have been called to get the number of
 * zones, then space allocated for the arrays */
void cdPressure_depth(
	/* total pressure, all forms*/
	double TotalPressure[],			
	/* gas pressure */
	double GasPressure[],				
	/* radiation pressure */
	double RadiationPressure[]);

/*
 * cdTemp_last
 * returns the temperature of the last zone on last iteration */
double cdTemp_last(void);

/*
 * cdIonFrac
 * This returns the ionization fraction for any element included in the calculation. 
 * The first parameter is 4 char null terminad string giving the first 4 letters of
 * element name as spelled by Cloudy.  
 * The second parameter is an integer giving the ionization stage, 
 * 1 for atom, 2 for first ion, etc.
 * The third parameter returns the predicted ionization fraction of that ion stage.
 * The last paramger is an 8 character + null string that says either "volume" or "radius",
 * to specify whether the average should be weighted by volume or radius.
 * The return value is 0 if the routine could find the species and
 * non-zero if it failed to find the element */
int cdIonFrac(
	/* four char string, null terminzed, giving the element name */
	const char *chLabel, 
	/* IonStage is ionization stage, 1 for atom, up to N+1 where N is atomic number */
	long int IonStage, 
	/* will be fractional ionization */
	double *fracin, 
	/* how to weight the average, must be "VOLUME" or "RADIUS" */
	const char *chWeight ,
	/* if true then weighting also has electron density, if false then only volume or radius */
	int lgDensity );

/*
 * cdVersion
 * The argument is a string with at least 8 characters that will receive a null terminated
 * string with the version number of the code. */
void cdVersion(char chString[] ) ;

/*
 * cdDate
 * The argument is a string with at least 8 char that will receive a null terminated
 * string with the date of the current version of the code. */
void cdDate(char chString[] ) ;

/* The followingh pairs of routines can keep track of the execution time for one model -
 * cdSetExecTime called first (in cdInit, not by the user) to initialize timer.
 * When cdExecTime is called it will return the ellapsed time in seconds
 * since cdInit called cdSetExecTime*/

/* normally called by cdInit, this routine sets initial variables for times */
void cdSetExecTime(void);

/* cdExecTime returns the cpu time (sec) that has elapsed since cdInit called cdSetExecTime.
 * The source is in setexectime */
double cdExecTime(void);

/* cdGetLineList will read in a list of emission line labels and wavelengths
 * from a file.  I use it for generating loc grids.  
 * Two files (cdGetLineList and cdGetLineList) are included in the main data 
 * distribution and have list of strong emission lines for high and low density gas.
 * other files can be created by the user.
 *
 * The first argument is the name of the file to read.
 * It it is void ("") then the routine will open LineList_BLR.dat 
 *
 * The next two arguments are pointers that will become arrays giving the
 * list of labels and wavelengths.  The routine will allocate the
 * needed space, but the pointer is defined in the calling routine.  
 * in the calling routine the two variable should be declared like this:
 * char **chLabels;
 * float *wavelength;
 * They would appear as follows in the call to the routine:
 * chGetLineList("", &chLabels , &wavelength );
 *
 * cdGetLineList returns the number of lines it found in the file if it was successful,
 * and -1 if it could not open the file.
 * */
long int cdGetLineList(
	/* chFile is optional filename, if void then use LineList_BLR.dat,
	 * if not void then use file specified */
	 char chFile[] ,
	 	/* 2d array of null term strings giving line labels char chLabels[nLines][10] */
	char ***chLabels ,
	/* a 1-d array of line wavelengths */
	float **wl );

/* 
 * cdTimescales returns the longest thermal, recombination, and H2 foramtion 
 * timescales that occurred in the previous model */
void cdTimescales(
	/* the thermal cooling timescale */
	double *TTherm , 
	/* the hydrogen recombination timescale */
	double *THRecom , 
	/* the H2 formation timescale */
	double *TH2 );

/* print line wavelengths in Angstroms in the standard format - just a wrapper */
void cdPrtWL( FILE *io , float wl );

/* ******************************************************************
 *
 * next part deals with FeII bands.  There are two types, the tabulated
 * band that are defined in bands_Fe2.dat, and the psuedo-continuum bins
 * that are generatedby the code in FeIIContCreate.
 * nFeIIConBins is number of continuum bins in FeII_Cont 
 * nFeIIBands is number of bands in bands_Fe2.dat, and are saved in FeII_Bands
 * the bands are created by hand and the entries in bands_Fe2.dat are
 * meant to be created by a person  */

/* the declarations for the next four are in FeIILevelPops.c */
/* this is the number of bands read in from bands_Fe2.dat */
extern long int nFeIIBands;

/* number of bands in continuum array */
extern long int nFeIIConBins;

/* band wavelength, lower and upper bounds, in vacuum Angstroms */
/* FeII.bands[n][3], where n is the number of bands in fe2bands.dat
 * these bands are defined in fe2bands.dat and read in at startup
 * of calculation */
extern float **FeII_Bands; 

/* continuum wavelengths, lower and upper bounds, in vacuum Angstroms
 * third is integrated intensity */
/* FeII_Cont[n][3], where n is the number of cells needed
 * these bands are defined in cdGetFeIIBands */
extern float **FeII_Cont; 


/*cdMPI sets flag telling exit handler to call MPI_Finalize, This routine
 * must be called after cdInit when MPI is used */
void cdMPI(void); /* set flag so that exit handler will clean up MPI */

/* 
 * this routine returns the spectrum needed for Keith Arnaud's XSPEC
 * X-Ray analysis code.  It should be called after cdDrive has successfully
 * computed a model.  The calling routine must ensure that the  vectors
 * have enough space to store the resulting spectrum, 
 * given the bounds and energy resolution */
void cdSPEC( 
	/* option - the type of spectrum to be returned
	* 1			the incident continuum 4\pi nuJ_nu, , erg cm-2 s-1
	*
	* 2			the attenuated incident continuum, same units
	* 3			the reflected continuum, same units
	*
	* 4			diffuse continuous emission outward direction
	* 5			diffuse continuous emission, reflected
	*
	* 6			collisional+recombination lines, outward
	* 7			collisional+recombination lines, reflected
	*
	*			all lines and continuum emitted by the cloud assume full coverage of 
	*			continuum source */
	int Option ,

	/* the energy of the lower edge of each cell 
	 * (in Ryd to be consistent with all the rest of Cloudy) */
	double EnergyLow[] , 

	/* the number of cells + 1*/
	long int nEnergy ,

	/* the returned spectrum, same size is two energy spectra (see option), returns nEnergy -1 pts */
	double ReturnedSpectrum[] );

void cdSPEC2( 
	/* option - the type of spectrum to be returned
	* 1			the incident continuum 4\pi nuJ_nu, , erg cm-2 s-1
	*
	* 2			the attenuated incident continuum, same units
	* 3			the reflected continuum, same units
	*
	* 4			diffuse continuous emission outward direction
	* 5			diffuse continuous emission, reflected
	*
	* 6			collisional+recombination lines, outward
	* 7			collisional+recombination lines, reflected
	*
	*			all lines and continuum emitted by the cloud assume full coverage of 
	*			continuum source */
	int Option ,

	/* the number of cells + 1*/
	long int nEnergy ,

	/* the returned spectrum, same size is two energy spectra (see option), returns nEnergy -1 pts */
	float ReturnedSpectrum[] );

/* cdTemp 
 * This routine finds the mean electron temperature for any ionization stage 
 * It returns 0 if it could find the species, 1 if it could not find the species.
 * The first argument is a null terminated 4 char string that gives the element
 * name as spelled by Cloudy.  
 * The second argument is ion stage, 1 for atom, 2 for first ion, etc
 * This third argument will be returned as result,
 * Last parameter is either "VOLUME" or "RADIUS" to give weighting 
 *
 * if the ion stage is zero then the element label will have a special meaning.
 * The string "21CM" is will return the 21 cm temperature.
 * The string "H2  " will return the temperature weighted wrt the H2 density */
int cdTemp(
	/* four char string, null terminzed, giving the element name */
	const char *chLabel, 
	/* IonStage is ionization stage, 1 for atom, up to N+1 where N is atomic number */
	long int IonStage, 
	/* will be temperature */
	double *TeMean, 
	/* how to weight the average, must be "VOLUME" or "RADIUS" */
	const char *chWeight );

/* cdPrintCommands( FILE *)
 * This routine prints all input commands into file whose hanlde is the argument */
void cdPrintCommands( FILE * );

/* wrapper to close all punch files */
void cdClosePunchFiles( void );

 /*cdH2_Line returns 1 if we found the line, 
  * or FALSE==0 if we did not find the line because ohoto-para transition
  * or upper level has lower energy than lower level 
  * NB - this is in mole_h2_io.c */
long int cdH2_Line(
	  /* indices for the upper level */
	  long int iElecHi, 
	  long int iVibHi ,
	  long int iRotHi ,
	  /* indices for lower level */
	  long int iElecLo, 
	  long int iVibLo ,
	  long int iRotLo ,
	  /* linear intensity relative to normalization line*/
	  double *relint, 
	  /* log of luminosity or intensity of line */
	  double *absint );

/* none of the following are generally needed */

/* this is the value that will be set true when cdInit is called.  
 * Other routines will check that this is TRUE when they are called, 
 * to verify that cdInit was called first.  Definition is in cdDefines.c */
extern int lgcdInitCalled;

/* some macros that will be of usein writing a new main.  In Cloudy itself
 * these may or many not already be defined, hence the if */
#ifndef TRUE
#define TRUE	1
#endif
/* */

#ifndef FALSE
#define FALSE	0
#endif
/* */

#ifndef MAX2
/* MAX2 takes two arguments, returns the larger of the two */
#define MAX2(a,b) (((a)>(b)) ? (a) : (b))
#endif
/* */

#ifndef MIN2
/* MIN2 takes two arguments, returns the smaller of the two */
#define MIN2(a,b) (((a)<(b)) ? (a) : (b))
#endif
/* */

