tblgen.c File Reference

Prepapre look-up tables from monochromatic absorption spectra. More...

#include "jurassic.h"

Go to the source code of this file.

Macros
#define	MAXNF   20000
#define	MAXNPTS   10000000

Functions
int	main (int argc, char *argv[])

Detailed Description

Prepapre look-up tables from monochromatic absorption spectra.

Definition in file tblgen.c.

Macro Definition Documentation

MAXNF

#define MAXNF   20000

Definition at line 32 of file tblgen.c.

MAXNPTS

#define MAXNPTS   10000000

Definition at line 35 of file tblgen.c.

Function Documentation

main()

int main	(	int	argc,
		char *	argv[]
	)

Definition at line 41 of file tblgen.c.

                {
 
  FILE *in;
 
  static char *line = NULL;
 
  static double dnu, abs[MAXNPTS], epsold, f, filt[MAXNF],
    nu, nu0, nu1, nuf[MAXNF], press, temp, u;
 
  static int i, idx, nf, npts;
 
  size_t line_buf_size = 0;
 
  /* Read command line arguments... */
  if (argc != 5)
    ERRMSG("Give parameters: <press> <temp> <spec> <filter>");
  sscanf(argv[1], "%lg", &press);
  sscanf(argv[2], "%lg", &temp);
 
  /* Compute column density [molec/cm^2] (1 km path length, 1 ppmv)... */
  u = 1e-6 * press * 100 / (1.380658e-23 * temp) * 1000 / 1e4;
 
  /* Read filter function... */
  if (!(in = fopen(argv[4], "r")))
    ERRMSG("Cannot open filter file!");
  while (getline(&line, &line_buf_size, in) != -1)
    if (sscanf(line, "%lg %lg", &nuf[nf], &filt[nf]) == 2)
      if (++nf >= MAXNF)
        ERRMSG("Too many points in filter function");
  fclose(in);
 
  /* Read spectrum... */
  if (!(in = fopen(argv[3], "r")))
    ERRMSG("Cannot open spectrum!");
  for (i = 0; i < 4; ++i)
    if (getline(&line, &line_buf_size, in) == -1)
      ERRMSG("Error while reading spectrum!");
  sscanf(line, "%d %lg %lg %lg", &npts, &nu0, &dnu, &nu1);
  if (npts > MAXNPTS)
    ERRMSG("Too many points in optical depth spectrum!");
  i = 0;
  while (getline(&line, &line_buf_size, in) != -1) {
    char *tok = strtok(line, " \t\n");
    while (tok != NULL) {
      if (i >= MAXNPTS)
        ERRMSG("Too many data points in spectrum file!");
      sscanf(tok, "%lg", &abs[i]);
      abs[i] /= u;
      ++i;
      tok = strtok(NULL, " \t\n");
    }
  }
  fclose(in);
 
  /* Set grid spacing... */
  dnu = (nu1 - nu0) / ((double) npts - 1.0);
  const double r0 = (nuf[0] - nu0) / (nu1 - nu0) * (double) npts;
  const int i0 = (int) r0;
 
  /* Loop over column densities... */
  for (u = 1.0; u <= 1e30; u *= 1.122) {
 
    /* Integrate... */
    double epssum = 0, fsum = 0;
    for (i = i0; i < npts; i++) {
      nu = nu0 + dnu * (double) i;
      if (nu < nuf[0])
        continue;
      else if (nu > nuf[nf - 1])
        break;
      else {
        if (nu < nuf[idx] || nu > nuf[idx + 1])
          idx = locate_irr(nuf, nf, nu);
        f = LIN(nuf[idx], filt[idx], nuf[idx + 1], filt[idx + 1], nu);
        fsum += f;
        epssum += f * exp(-abs[i] * u);
      }
    }
    epssum = 1 - epssum / fsum;
 
    /* Write output... */
    if (epssum >= 1e-6 && epssum <= 0.999999 && epssum > epsold)
      printf("%g %g %g %g\n", press, temp, u, epssum);
    epsold = epssum;
 
    /* Check for termination... */
    if (epssum > 0.999999)
      return EXIT_SUCCESS;
  }
 
  return EXIT_SUCCESS;
}

Here is the call graph for this function: