52 fprintf(
ioQQQ,
" FFmtRead ParseDrive entered. Enter number.\n" );
58 fprintf(
ioQQQ,
" ParseDrive.dat error getting magic number\n" );
62 fac =
FFmtRead(chInput,&i,
sizeof(chInput),&lgEOL);
65 fprintf(
ioQQQ,
" FFmtRead hit the EOL with no value, return=%10.2e\n",
75 fprintf(
ioQQQ,
" FFmtRead returned with value%11.4e\n",
78 fprintf(
ioQQQ,
" Enter 0 to stop, or another value.\n" );
80 fprintf(
ioQQQ,
" FFmtRead ParseDrive exits.\n" );
83 else if( p.
nMatch(
"CASE") )
89 else if( p.
nMatch(
"CDLI") )
92 trace.lgDrv_cdLine =
true;
95 else if( p.
nMatch(
" E1 ") )
100 for( i=0; i<50; ++i )
102 fprintf(
ioQQQ,
"tau\t%.3e\t exp-tau\t%.5e\t e1 tau\t%.5e \t e2 "
103 "\t%.5e \te2n %.5e \t e3\t%.5e \t e4\t%.5e \n",
106 tau = pow( 10. , ((
double)i/4. - 9.) );
111 else if( p.
nMatch(
"ESCA") )
117 else if( p.
nMatch(
"HYAS") )
123 else if( p.
nMatch(
"GAUN") )
129 else if( p.
nMatch(
"POIN") )
132 fprintf(
ioQQQ,
" Define entire model first, later I will ask for energies.\n" );
133 trace.lgPtrace =
true;
136 else if( p.
nMatch(
"PUMP") )
140 fprintf(
ioQQQ,
" Continuum pump ParseDrive entered - Enter log tau\n" );
145 fprintf(
ioQQQ,
" ParseDrive.dat error getting magic number\n" );
150 fac =
FFmtRead(chInput,&i,
sizeof(chInput),&lgEOL);
154 fprintf(
ioQQQ,
" Tau =%11.4e\n", fac );
155 (*TauDummy).Emis().TauIn() = (
realnum)fac;
159 fprintf(
ioQQQ,
" ContPump =%11.4e\n", fac );
160 fprintf(
ioQQQ,
" Enter null to stop, or another value.\n" );
162 fprintf(
ioQQQ,
" ContPump ParseDrive exits.\n" );
165 else if( p.
nMatch(
"STAR") )
169 for( i=0; i < 40; i++ )
171 zed = ((double)i+1.)/4. + 0.01;
172 sprintf( chInput,
"starburst, zed=%10.4f", zed );
175 fprintf(
ioQQQ,
"%8.1e", zed );
178 fprintf(
ioQQQ,
"\n" );
182 else if( p.
nMatch(
"VOIGT") )
185 FILE *ioVOIGT = fopen(
"voigt.dat" ,
"w");
186 fprintf(ioVOIGT,
"x \\ a");
187 const realnum DampLogMin = -4., DampLogMax = 4.01;
188 for(
realnum damplog=DampLogMin; damplog<DampLogMax; ++damplog)
189 fprintf(ioVOIGT,
"\tlog a=%.3e",pow(10.,damplog));
190 fprintf(ioVOIGT ,
"\n");
192 for(
realnum x=-2.; x<5.;x+=0.05)
195 fprintf(ioVOIGT,
"%.3e",xlin);
196 for(
realnum damplog=DampLogMin; damplog<DampLogMax; ++damplog)
203 fprintf(ioVOIGT ,
"\t%.3e",yval[0]);
205 fprintf(ioVOIGT ,
"\n");
214 " Unrecognized key; keys are CASE, CDLIne, E1 , ESCApe, FFMTread, GAUNt, "
215 "HYAS, POINt, PUMP, STAR, and VOIGt. Sorry.\n" );
233 fprintf(
ioQQQ,
" Enter the log of the one-sided optical depth; line with no number to stop.\n" );
244 tau =
FFmtRead(chCard,&i,
sizeof(chCard),&lgEOL);
251 fprintf(
ioQQQ,
"tau was %e\n", tau );
279 fprintf(
ioQQQ,
" I will get needed H data files. This will take a second.\n");
285 enum {DEBUG_LOC=
false};
289 double xLyman , alpha;
294 for( ipHi=3; ipHi<25; ++ipHi )
296 double photons = (1./
POW2(ipHi-1.)-1./
POW2((
double)ipHi) ) /(1.-1./ipHi/ipHi );
298 alpha =
atmdat_HS_caseB(ipHi-1,ipHi, nelem,Temperature , Density ,
'A' );
299 fprintf(
ioQQQ,
"%li\t%.3e\t%.3e\n", ipHi, xLyman/alpha*photons, photons );
310 fprintf(
ioQQQ,
" Enter atomic number of species, either 1(H) or 2(He).\n" );
313 fprintf(
ioQQQ,
" error getting species \n" );
318 nelem = (
long int)
FFmtRead(chCard,&i,
sizeof(chCard),&lgEOL);
319 if( lgEOL || nelem< 1 || nelem > 2 )
321 fprintf(
ioQQQ,
" must be either 1 or 2!\n" );
329 fprintf(
ioQQQ,
" In the following temperatures <10 are log, >=10 linear.\n");
330 fprintf(
ioQQQ,
" The density is always a log.\n");
331 fprintf(
ioQQQ,
" The order of the quantum numbers do not matter.\n");
332 fprintf(
ioQQQ,
" The smallest must not be smaller than 2,\n");
333 fprintf(
ioQQQ,
" and the largest must not be larger than 25.\n");
334 fprintf(
ioQQQ,
" Units of emissivity are erg cm^3 s^-1\n\n");
335 fprintf(
ioQQQ,
" The limits of the HS tables are 2 <= n <= 25.\n");
341 fprintf(
ioQQQ,
" Enter 4 numbers, temperature, density, 2 quantum numbers, null line stop.\n" );
344 fprintf(
ioQQQ,
" Thanks for interpolating on the Hummer & Storey data set!\n" );
349 Temperature =
FFmtRead(chCard,&i,
sizeof(chCard),&lgEOL);
352 fprintf(
ioQQQ,
" error getting temperature!\n" );
357 if( Temperature < 10. )
359 Temperature = pow(10., Temperature );
361 fprintf(
ioQQQ,
" Temperature is %g\n", Temperature );
363 Density =
FFmtRead(chCard,&i,
sizeof(chCard),&lgEOL);
366 fprintf(
ioQQQ,
" error getting density!\n" );
369 Density = pow(10., Density );
370 fprintf(
ioQQQ,
" Density is %g\n", Density );
373 n1 = (long)
FFmtRead(chCard,&i,
sizeof(chCard),&lgEOL);
376 fprintf(
ioQQQ,
" error getting quantum number!\n" );
380 n2 = (long)
FFmtRead(chCard,&i,
sizeof(chCard),&lgEOL);
383 fprintf(
ioQQQ,
" error getting quantum number!\n" );
387 if(
MAX2( n1 , n2 ) > 25 )
389 fprintf(
ioQQQ,
" The limits of the HS tables are 2 <= n <= 25. Sorry.\n");
394 " 4pJ(%ld,%ld)/n_e n_p=%11.3e\n",
402 double tempTable[33] = {
403 11870.,12490.,12820.,
404 11060.,17740.,12560.,
405 16390.,16700.,11360.,
406 10240.,20740.,12030.,
407 14450.,19510.,12550.,
408 16470.,16560.,12220.,
409 15820.,12960.,10190.,
410 12960.,14060.,12560.,
411 11030.,10770.,13360.,
412 10780.,11410.,11690.,
413 12500.,13190.,21120. };
414 double edenTable[33] = {
415 10.,270.,80.,10.,70.,
416 110.,200.,10.,40.,90.,
417 340.,80.,60.,340.,30.,
418 120.,10.,50.,450.,30.,
419 180.,20.,170.,60.,20.,
420 40.,30.,20.,100.,130.,
424 for( j=0; j<33; j++ )
426 double halpha, hbeta, hgamma;
432 fprintf(
ioQQQ,
"%e\t%e\t%e\t%e\n",
442 fprintf(
ioQQQ,
" Thanks for interpolating on the Hummer & Storey data set!\n" );
452 long int i, nHi, lHi, nLo, lLo;
463 fprintf(
ioQQQ,
" Enter four quantum numbers (n, l, n', l'), null line to stop.\n" );
466 fprintf(
ioQQQ,
" error getting drvhyas \n" );
471 nHi = (
long int)
FFmtRead(chCard,&i,
sizeof(chCard),&lgEOL);
475 lHi = (
long int)
FFmtRead(chCard,&i,
sizeof(chCard),&lgEOL);
478 fprintf(
ioQQQ,
" must be four numbers!\n" );
484 fprintf(
ioQQQ,
" l must be less than n!\n" );
488 nLo = (
long int)
FFmtRead(chCard,&i,
sizeof(chCard),&lgEOL);
491 fprintf(
ioQQQ,
" must be four numbers!\n" );
495 lLo = (
long int)
FFmtRead(chCard,&i,
sizeof(chCard),&lgEOL);
498 fprintf(
ioQQQ,
" must be four numbers!\n" );
504 fprintf(
ioQQQ,
" l must be less than n!\n" );
521 fprintf(
ioQQQ,
" A(%3ld,%3ld->%3ld,%3ld)=%11.3e\n",
526 fprintf(
ioQQQ,
" Driver exits, enter next line.\n" );
542 double loggamma2, logu;
551 fprintf(
ioQQQ,
" Enter 0 to input temp, energy, and net charge, or 1 for u and gamma**2.\n" );
554 fprintf(
ioQQQ,
" dgaunt error getting magic number\n" );
558 inputflag = (int)
FFmtRead(chCard,&i,
sizeof(chCard),&lgEOL);
562 fprintf(
ioQQQ,
" Enter the temperature (log if <=10), energy (Ryd), and net charge. Null line to stop.\n" );
570 fprintf(
ioQQQ,
" dgaunt error getting magic number\n" );
578 fprintf(
ioQQQ,
" Gaunt driver exits, enter next line.\n" );
589 TeNew = pow(10.,
phycon.alogte);
594 if( lgEOL || enerlin[0] == 0. )
596 fprintf(
ioQQQ,
" Sorry, but there should be two more numbers, energy and charge.\n" );
599 z = (double)
FFmtRead(chCard,&i,
sizeof(chCard),&lgEOL);
600 if( lgEOL || z == 0. )
602 fprintf(
ioQQQ,
" Sorry, but there should be a third number, charge.\n" );
608 fprintf(
ioQQQ,
" Using my routine, Gff= \t" );
609 fprintf(
ioQQQ,
"%11.3e\n", mygaunt );
618 fprintf(
ioQQQ,
" Enter log u and log gamma2. Null line to stop.\n" );
624 fprintf(
ioQQQ,
" dgaunt error getting magic number\n" );
628 logu =
FFmtRead(chCard,&i,
sizeof(chCard),&lgEOL);
632 fprintf(
ioQQQ,
" Gaunt driver exits, enter next line.\n" );
636 loggamma2 =
FFmtRead(chCard,&i,
sizeof(chCard),&lgEOL);
639 fprintf(
ioQQQ,
" Sorry, but there should be another numbers, log gamma2.\n" );
647 fprintf(
ioQQQ,
" Using my routine, Gaunt factor is" );
648 fprintf(
ioQQQ,
"%11.3e\n", mygaunt );
void abund_starburst(Parser &p)
double atmdat_HS_caseB(long int iHi, long int iLo, long int iZ, double TempIn, double DenIn, char chCase)
sys_float sexp(sys_float x)
const int INPUT_LINE_LENGTH
double FFmtRead(const char *chCard, long int *ipnt, long int last, bool *lgEOL)
char * read_whole_line(char *chLine, int nChar, FILE *ioIN)
#define DEBUG_ENTRY(funcname)
bool nMatch(const char *chKey) const
void setline(const char *const card)
double cont_gaunt_calc(double temp, double z, double photon)
double DrvContPump(const TransitionProxy &t, realnum DopplerWidth)
double H_Einstein_A(long int n, long int l, long int np, long int lp, long int iz)
STATIC void DrvCaseBHS(void)
STATIC void DrvEscP(void)
STATIC void DrvHyas(void)
void ParseDrive(Parser &p)
UNUSED const double TE1RYD
double esc_PRD_1side(double tau, double a)
double esca0k2(double taume)
double esc_CRDwing_1side(double tau, double a)
TransitionProxy::iterator TauDummy
void TempChange(double TempNew, bool lgForceUpdate)
double expn(int n, double x)
void VoigtH(realnum a, const realnum v[], realnum y[], int n)