26 static const int NISM = 23;
35 static const double tnuHM96[
NHM96]={-8,-1.722735683,-0.351545683,-0.222905683,-0.133385683,
37 -0.127655683,-0.004575683,0.297544317,0.476753,0.476756,0.588704317,
38 0.661374317,1.500814317,2.245164317};
41 static const double fnuHM96[
NHM96]={-32.53342863,-19.9789,-20.4204,-20.4443,-20.5756,-20.7546,
42 -21.2796,-21.6256,-21.8404,-21.4823,-22.2102,-22.9263,-23.32,-24.2865};
65 const double RESETFACTOR = 0.98;
74 power = (fluxlog[1] - fluxlog[0] ) / log10( tnu[1]/tnu[0] );
77 fluxlog[0] = fluxlog[1] + power * log10( tnu[0] /tnu[1] );
83 power = log10( fluxlog[1]/fluxlog[0]) / log10( tnu[1]/tnu[0] );
86 fluxlog[0] = log10(fluxlog[1]) + power * log10( tnu[0] /tnu[1] );
88 fluxlog[0] =
exp10(fluxlog[0]);
116 bool lgNoContinuum =
false,
130 fprintf(
ioQQQ,
" %ld is too many spectra entered. Increase LIMSPC\n Sorry.\n",
140 lgQuoteFound =
false;
144 bool lgKeyword =
false;
149 chFile =
"akn120.sed";
157 chFile =
"CrabDavidson.sed";
161 chFile =
"CrabHester.sed";
174 chFile =
"Rubin.sed";
180 chFile =
"Trapezium.sed";
203 bool lgHM05 =
false, lgHM12 =
false;
209 for( i=0; i <
NAGN; i++ )
230 fprintf(
ioQQQ,
" There must be a number for the break.\n Sorry.\n" );
237 fprintf(
ioQQQ,
" The break must be greater than 0.2 Ryd.\n Sorry.\n" );
244 ConBreak = 0.0912/ConBreak;
250 ConBreak =
exp10(ConBreak);
253 else if( ConBreak == 0. )
255 fprintf(
ioQQQ,
" An energy of 0 is not allowed.\n Sorry.\n" );
261 fprintf(
ioQQQ,
" The energy of the break cannot be greater than%10.2e Ryd.\n Sorry.\n",
268 fprintf(
ioQQQ,
" The energy of the break cannot be less than%10.2e Ryd.\n Sorry.\n",
295 for( j=0; j <
NHM96; j++ )
307 scale = log10(scale);
312 fprintf(
ioQQQ,
" %ld is too many continua entered. Increase LIMSPC\n Sorry.\n",
322 fprintf(
ioQQQ,
" This command has come between a previous ordered pair of continuum shape and luminosity commands.\n Reorder the commands to complete each continuum specification before starting another.\n" );
353 bool lgQuasar = p.
nMatch(
"QUAS");
354 if( lgHM12 && lgQuasar )
356 fprintf(
ioQQQ,
" The QUASAR option is not supported on the TABLE HM12 command.\n" );
368 scale = log10(scale);
370 UNUSED double zlow, zhigh;
371 int version = lgHM05 ? 2005 : 2012;
385 fprintf(
ioQQQ,
" %ld is too many continua entered. Increase LIMSPC\n Sorry.\n", p.
m_nqh );
394 fprintf(
ioQQQ,
" This command has come between a previous ordered pair of"
395 " continuum shape and luminosity commands.\n Reorder the commands"
396 " to complete each continuum specification before starting another.\n" );
414 else if( p.
nMatch(
" ISM") )
421 for( i=6; i <
NISM; i++ )
435 if( scale > 0. && !p.
nMatch(
" LOG"))
436 scale = log10(scale);
441 fprintf(
ioQQQ,
" %4ld is too many continua entered. Increase LIMSPC\n Sorry.\n",
451 fprintf(
ioQQQ,
" This command has come between a previous ordered pair of continuum shape and luminosity commands.\n Reorder the commands to complete each continuum specification before starting another.\n" );
491 else if( p.
nMatch(
"DRAI") )
507 if( scale > 0. && !p.
nMatch(
" LOG") )
508 scale = log10(scale);
513 fprintf(
ioQQQ,
" %4ld is too many continua entered. Increase LIMSPC\n Sorry.\n",
523 fprintf(
ioQQQ,
" This command has come between a previous ordered pair of continuum shape and luminosity commands.\n Reorder the commands to complete each continuum specification before starting another.\n" );
565 else if( p.
nMatch(
"LINES") )
572 lgNoContinuum =
true;
576 fprintf(
ioQQQ,
" sorry, only one table line per input stream\n");
582 if( lgQuoteFound && chFile.length() > 0 )
592 fprintf(
ioQQQ,
"\n DISASTER PROBLEM ParseTable could not find "
594 fprintf(
ioQQQ,
" Please check the spelling of the file name and that it "
595 "is in either the local or data directory.\n\n");
605 else if( p.
nMatch(
"POWE") )
632 else if( brakmm == 0. )
641 else if( brakmm < 0. )
645 brakmm =
exp10(brakmm);
650 brakmm = RYDLAM / (1e4*brakmm);
657 " Check the order of parameters on this table power law command - the low-energy break of %f Ryd seems high to me.\n",
683 else if( brakxr == 0. )
691 brakxr =
exp10(brakxr);
702 if( brakmm >= brakxr )
704 fprintf(
ioQQQ,
" HELP!! The lower energy for the power law, %f, is greater than the upper energy, %f. This is not possible.\n Sorry.\n",
727 else if( p.
nMatch(
"READ") )
734 fprintf(
ioQQQ,
" Name of file must appear on TABLE READ.\n");
749 if( scale > 0. && !p.
nMatch(
" LOG"))
750 scale = log10(scale);
754 fprintf(
ioQQQ,
" %ld is too many continua entered. Increase LIMSPC\n Sorry.\n",
763 fprintf(
ioQQQ,
" This command has come between a previous ordered pair of continuum shape and luminosity commands.\n Reorder the commands to complete each continuum specification before starting another.\n" );
804 fprintf(
ioQQQ,
" The TABLE TLUSTY command is no longer supported.\n" );
805 fprintf(
ioQQQ,
" Please use TABLE STAR TLUSTY instead. See Hazy for details.\n" );
809 else if( p.
nMatch(
"SED") || lgKeyword )
816 fprintf(
ioQQQ,
"PROBLEM in TABLE SED: No quotes were found.\n The SED table file must be designated in quotes.\n" );
825 strcpy( chPath,
"SED" );
827 strcat( chPath, chFile.c_str() );
833 long int fileLength = 0;
834 const char *chEnergyUnits=NULL;
835 bool lgUnitsSet =
false, lgNuFnu =
false, lgExtrapolate =
false;
837 while(
read_whole_line( chLine , (
int)
sizeof(chLine) , ioDATA ) != NULL )
847 if( chLine[0]==
'\n' || chLine[0]==
'\0' )
849 fprintf(
ioQQQ,
"PROBLEM in TABLE SED: Encountered unexpected empty line.\n");
855 strcpy( chLineCAPS , chLine );
859 if(
nMatch(
"UNIT", chLineCAPS ) != 0)
866 if(
nMatch(
"NUFN" , chLineCAPS ) != 0)
871 if(
nMatch(
"EXTRAP" , chLineCAPS ) != 0)
872 lgExtrapolate =
true;
877 double *tnuInput = NULL;
878 double *tslopInput = NULL;
879 tnuInput = (
double *)
MALLOC( (
size_t)fileLength*
sizeof(double ) );
880 tslopInput = (
double *)
MALLOC( (
size_t)fileLength*
sizeof(double ) );
884 long int entryNum = 0;
886 while(
read_whole_line( chLine , (
int)
sizeof(chLine) , ioDATA ) != NULL )
895 tnuInput[entryNum] =
FFmtRead(chLine,&i,
sizeof(chLine),&lgEOL);
898 fprintf(
ioQQQ,
"PROBLEM in TABLE SED: frequency and flux pairs must be entered. The first number was not found"
899 " on this line:\n%s\n", chLine);
902 tslopInput[entryNum] =
FFmtRead(chLine,&i,
sizeof(chLine),&lgEOL);
908 fprintf(
ioQQQ,
"PROBLEM in TABLE SED: frequency and flux pairs must be entered. Two numbers were not found"
909 "on this line:\n%s\n", chLine);
913 if( tslopInput[entryNum] <= 0. )
915 fprintf(
ioQQQ,
"PROBLEM flux in TABLE SED file entry %li is not positive\nphoton energy %e where the value is %e.\n",
916 entryNum , tnuInput[entryNum] , tslopInput[entryNum] );
921 if( tnuInput[entryNum] <= 0. )
923 fprintf(
ioQQQ,
"PROBLEM photon energy in TABLE SED file entry %li is not positive, photon energy is %e.\n",
924 entryNum , tnuInput[entryNum] );
929 if( entryNum>1 && (tnuInput[entryNum-2]-tnuInput[entryNum-1])*
930 (tnuInput[entryNum-1]-tnuInput[entryNum])<=0 )
932 fprintf(
ioQQQ,
"PROBLEM photon energies in TABLE SED file must be in increasing order\n");
933 fprintf(
ioQQQ,
"The following photon energy sequence is not monotonically changing:\n");
934 for(
long ij=entryNum-2; ij<=entryNum; ++ij )
939 tslopInput[entryNum] = log10(tslopInput[entryNum]);
946 for(
long i=0; i < entryNum; i++ )
948 unitChange.
set(tnuInput[i], chEnergyUnits );
949 tnuInput[i] = unitChange.
Ryd();
951 ASSERT( tnuInput[i]>tnuInput[i-1]);
958 for(
long i=0; i < entryNum; i++ )
959 tslopInput[i] -= log10(tnuInput[i]);
966 for(
long i=0; i < entryNum; i++ )
980 else if( p.
nMatch(
"STAR") )
982 char chMetalicity[5] =
"", chODFNew[8], chVaryFlag[7] =
"";
983 bool lgHCa =
false, lgHNi =
false;
985 double Tlow = -1., Thigh = -1.;
1001 sprintf(chVaryFlag,
"%1ld-DIM",ndim);
1007 static const char table[][2][10] = {
1039 strncpy( chMetalicity,
"p00",
sizeof(chMetalicity) );
1040 for (i=0; i < (long)(
sizeof(table)/
sizeof(*table)); ++i)
1044 strncpy( chVaryFlag, table[i][0],
sizeof(chVaryFlag) );
1045 strncpy( chMetalicity, table[i][1],
sizeof(chMetalicity) );
1053 bool lgHalo = p.
nMatch(
"HALO");
1054 bool lgSolar = ( !lgHalo && strcmp( chMetalicity,
"p00" ) == 0 );
1059 bool lgList = p.
nMatch(
"LIST");
1068 for( nval=0; nval <
MDIM; nval++ )
1075 if( nval == 0 && !lgList )
1077 fprintf(
ioQQQ,
" The stellar temperature MUST be entered.\n" );
1085 val[0] =
exp10(val[0]);
1088 fprintf(
ioQQQ,
" DISASTER the log of the parameter was specified, "
1089 "but the numerical value is too large.\n Sorry.\n\n");
1096 nstar =
GridInterpolate(val,&nval,&ndim,chFile.c_str(),lgList,&Tlow,&Thigh);
1099 else if( p.
nMatch(
"ATLA") )
1106 strncpy( chODFNew,
"_odfnew",
sizeof(chODFNew) );
1108 strncpy( chODFNew,
"",
sizeof(chODFNew) );
1111 nstar =
AtlasInterpolate(val,&nval,&ndim,chMetalicity,chODFNew,lgList,&Tlow,&Thigh);
1114 else if( p.
nMatch(
"COST") )
1131 if( val[1] < 1. || val[1] > 7. )
1133 fprintf(
ioQQQ,
" The second number must be the id sequence number, 1 to 7.\n" );
1138 else if( p.
nMatch(
"ZAMS") )
1143 fprintf(
ioQQQ,
" A second number, the age of the star, must be present.\n" );
1147 else if( p.
nMatch(
" AGE") )
1152 fprintf(
ioQQQ,
" A second number, the ZAMS mass of the star, must be present.\n" );
1174 if( ! ( lgSolar || lgHalo ) )
1176 fprintf(
ioQQQ,
" Please choose SOLAR or HALO abundances.\n" );
1183 else if( p.
nMatch(
"KURU") )
1188 else if( p.
nMatch(
"MIHA") )
1193 else if( p.
nMatch(
"RAUC") )
1195 if( ! ( lgSolar || lgHalo ) )
1197 fprintf(
ioQQQ,
" Please choose SOLAR or HALO abundances.\n" );
1209 else if( p.
nMatch(
"HYDR") )
1213 else if( p.
nMatch(
"HELI") )
1217 else if( p.
nMatch(
"H+HE") )
1225 else if( p.
nMatch(
"CO W") )
1236 else if( p.
nMatch(
"TLUS") )
1244 else if( p.
nMatch(
"BSTA") )
1250 else if( p.
nMatch(
"OSTA") )
1258 fprintf(
ioQQQ,
" There must be a third key on TABLE STAR TLUSTY;" );
1259 fprintf(
ioQQQ,
" the options are OBSTAR, BSTAR, OSTAR.\n" );
1264 else if( p.
nMatch(
"WERN") )
1271 else if( p.
nMatch(
"WMBA") )
1280 fprintf(
ioQQQ,
" There must be a second key on TABLE STAR;" );
1281 fprintf(
ioQQQ,
" the options are ATLAS, KURUCZ, MIHALAS, RAUCH, WERNER, and WMBASIC.\n" );
1301 for( i=1; i < nval; i++ )
1321 else if( p.
nMatch(
"COST") )
1349 else if( p.
nMatch(
"KURU") )
1354 "TABLE STAR KURUCZ %f LOG" );
1357 else if( p.
nMatch(
"MIHA") )
1362 "TABLE STAR MIHALAS %f LOG" );
1365 else if( p.
nMatch(
"RAUC") )
1374 else if( p.
nMatch(
"HELI") )
1376 else if( p.
nMatch(
"H+HE") )
1380 else if( p.
nMatch(
"CO W") )
1385 if( ( lgHCa || lgHNi ) && ndim == 2 )
1411 else if( p.
nMatch(
"BSTA") )
1413 else if( p.
nMatch(
"OSTA") )
1424 else if( p.
nMatch(
"WERN") )
1429 "TABLE STAR WERNER %f LOG %f" );
1432 else if( p.
nMatch(
"WMBA") )
1437 "TABLE STAR WMBASIC %f LOG %f %f" );
1450 for( i=1; i < nval; i++ )
1464 fprintf(
ioQQQ,
" There MUST be a keyword on this line. The keys are:"
1465 "AGN, DRAINE, HM96, HM05, HM12, _ISM, LINES, POWERlaw, "
1466 "READ, SED, and STAR.\n Sorry.\n" );
1514 fprintf(
ioQQQ,
" Table for this continuum; TNU, TFAC, TSLOP\n" );
1587 ++i;
tsldrn[i] = -17.7575;
1588 ++i;
tsldrn[i] = -17.7268;
1589 ++i;
tsldrn[i] = -17.7036;
1590 ++i;
tsldrn[i] = -17.6953;
1591 ++i;
tsldrn[i] = -17.7182;
1592 ++i;
tsldrn[i] = -17.7524;
1593 ++i;
tsldrn[i] = -17.8154;
1594 ++i;
tsldrn[i] = -17.9176;
1595 ++i;
tsldrn[i] = -18.1675;
1596 ++i;
tsldrn[i] = -18.1690;
1597 ++i;
tsldrn[i] = -18.2477;
1598 ++i;
tsldrn[i] = -18.4075;
1599 ++i;
tsldrn[i] = -18.5624;
1600 ++i;
tsldrn[i] = -18.7722;
1694 if( nLINE_TABLE == 0 )
1701 for(
long n=0; n < nLINE_TABLE; ++n )
1706 fprintf(
ioQQQ,
"lines_table in parse_table.cpp did not find line %s ",chLabel[n].c_str());
1714 fprintf(
ioQQQ ,
" BOTCHED MONITORS!!! Botched Monitors!!! lines_table could not find a total of %li lines\n\n", miss );
1744 fprintf(
ioQQQ,
" the input continuum file has been truncated.\n" );
1753 char *first = strrchr(chLine,
'/');
1756 unit = string(first+1);
1764 fprintf(
ioQQQ,
" the input continuum file has been truncated.\n" );
1772 sscanf( chLine,
"%ld", &version );
1776 " the input continuum file has the wrong version number, I expected %li and found %li.\n",
1783 fprintf(
ioQQQ,
" the input continuum file has been truncated.\n" );
1789 double mesh_lo, mesh_hi;
1798 strncpy( md5sum, chLine,
NMD5 );
1802 fprintf(
ioQQQ,
" the input continuum file has been truncated.\n" );
1810 sscanf( chLine,
"%x %x", &u.i[0], &u.i[1] );
1812 sscanf( chLine,
"%x %x", &u.i[1], &u.i[0] );
1817 fprintf(
ioQQQ,
" the input continuum file has been truncated.\n" );
1825 sscanf( chLine,
"%x %x", &u.i[0], &u.i[1] );
1827 sscanf( chLine,
"%x %x", &u.i[1], &u.i[0] );
1832 fprintf(
ioQQQ,
" the input continuum file has been truncated.\n" );
1840 fprintf(
ioQQQ,
" the input continuum file has an incompatible energy grid.\n" );
1848 sscanf( chLine,
"%x %x", &u.i[0], &u.i[1] );
1850 sscanf( chLine,
"%x %x", &u.i[1], &u.i[0] );
1855 fprintf(
ioQQQ,
" the input continuum file has been truncated.\n" );
1862 sscanf( chLine,
"%ld", &nflux );
1866 fprintf(
ioQQQ,
" the input continuum file has been truncated.\n" );
1873 fprintf(
ioQQQ,
" the input continuum file has been truncated.\n" );
1883 sscanf( chLine,
"%lf%lf ", &help[0], &help[1] );
1887 fprintf(
ioQQQ,
"\n\nPROBLEM in TABLE READ continuum frequencies: point %li should "
1888 "have %e and has %e\n", (
unsigned long)i,
rfield.
anu(i), help[0] );
1915 fprintf(
ioQQQ,
" the input continuum file has the wrong number of points: %ld\n",
1924 fprintf(
ioQQQ,
" please recreate this file using the SAVE TRANSMITTED CONTINUUM command.\n" );
bool nMatch(const char *chKey) const
void prt_wl(FILE *ioOUT, realnum wl)
FILE * open_data(const char *fname, const char *mode, access_scheme scheme)
long RauchInterpolateHydr(double val[], long *nval, long *ndim, bool lgList, double *Tlow, double *Thigh)
static double tsldrn[NDRAINE]
const int FILENAME_PATH_LENGTH_2
NORETURN void TotalInsanity(void)
STATIC void ReadTable(const char *fnam)
long findline(const char *chLabel, realnum wavelength)
string mesh_md5sum() const
long CoStarInterpolate(double val[], long *nval, long *ndim, IntMode imode, bool lgHalo, bool lgList, double *val0_lo, double *val0_hi)
static const long VERSION_TRNCON
long TlustyInterpolate(double val[], long *nval, long *ndim, tl_grid tlg, const char *chMetalicity, bool lgList, double *Tlow, double *Thigh)
void ParseTable(Parser &p)
long nMatch(const char *chKey, const char *chCard)
long HaardtMadauInterpolate(double val, int version, bool lgQuasar, double *zlow, double *zhigh)
static const double tnuHM96[NHM96]
int GetQuote(string &chLabel)
realnum varang[LIMPAR][2]
static double tnuism[NISM]
bool fp_equal_tol(sys_float x, sys_float y, sys_float tol)
long RauchInterpolateHCa(double val[], long *nval, long *ndim, bool lgHalo, bool lgList, double *Tlow, double *Thigh)
bool nMatchErase(const char *chKey)
vector< Energy > tNu[LIMSPC]
vector< realnum > tFluxLog[LIMSPC]
static double tnudrn[NDRAINE]
long int cdGetLineList(const char chFile[], vector< string > &chLabels, vector< realnum > &wl)
char chVarFmt[LIMPAR][FILENAME_PATH_LENGTH_2]
long WMBASICInterpolate(double val[], long *nval, long *ndim, bool lgList, double *Tlow, double *Thigh)
realnum vparm[LIMEXT][LIMPAR]
vector< realnum > tslop[LIMSPC]
double anu(size_t i) const
static string chLINE_LIST
STATIC void resetBltin(double *tnu, double *fluxlog, bool lgLog)
size_t ipointC(double anu) const
long RauchInterpolateHpHe(double val[], long *nval, long *ndim, bool lgList, double *Tlow, double *Thigh)
long MihalasInterpolate(double val[], long *nval, long *ndim, bool lgList, double *Tlow, double *Thigh)
Illuminate::IlluminationType Illumination[LIMSPC]
const int INPUT_LINE_LENGTH
NORETURN void NoNumb(const char *chDesc) const
const char * strchr_s(const char *s, int c)
static const double fnuHM96[NHM96]
static double tnuagn[NAGN]
STATIC void ZeroContin(void)
long RauchInterpolateCOWD(double val[], long *nval, long *ndim, bool lgList, double *Tlow, double *Thigh)
const char * StandardEnergyUnit(const char *chCard)
long WernerInterpolate(double val[], long *nval, long *ndim, bool lgList, double *Tlow, double *Thigh)
#define DEBUG_ENTRY(funcname)
static const unsigned int NMD5
long RauchInterpolateHelium(double val[], long *nval, long *ndim, bool lgList, double *Tlow, double *Thigh)
long GridInterpolate(double val[], long *nval, long *ndim, const char *FileName, bool lgList, double *Tlow, double *Thigh)
int fprintf(const Output &stream, const char *format,...)
char * read_whole_line(char *chLine, int nChar, FILE *ioIN)
static double tslagn[NAGN]
long Kurucz79Interpolate(double val[], long *nval, long *ndim, bool lgList, double *Tlow, double *Thigh)
long RauchInterpolatePG1159(double val[], long *nval, long *ndim, bool lgList, double *Tlow, double *Thigh)
long RauchInterpolateHNi(double val[], long *nval, long *ndim, bool lgHalo, bool lgList, double *Tlow, double *Thigh)
long AtlasInterpolate(double val[], long *nval, long *ndim, const char *chMetalicity, const char *chODFNew, bool lgList, double *Tlow, double *Thigh)
static double fnuism[NISM]
double FFmtRead(const char *chCard, long int *ipnt, long int last, bool *lgEOL)