50 unsigned long StrLen =
min(strlen(chLine),ArrLen);
51 ASSERT( StrLen < LineLen );
52 unsigned long pad = (LineLen-StrLen)/2;
53 for(
unsigned long i=0; i < pad; ++i )
67 double ratio = ( q2 >
SMALLFLOAT ) ? q1/q2 : 0.;
75 vector<realnum> sclsav, scaled;
76 vector<long int> ipSortLines;
77 vector<double> xLog_line_lumin;
78 vector<short int> lgPrt;
129 for( ipEmType=0; ipEmType<nEmType; ++ipEmType )
145 char wl_str[20] =
"";
149 " >>PROBLEM Normalization line (\"%s\" %s) has small or zero intensity, its value was %.2e and its intensity was set to 1.\n"
150 " >>Please consider using another normalization line (this is set with the NORMALIZE command).\n",
152 fprintf(
ioQQQ,
" >>The relative intensities will be meaningless, and many lines may appear too faint.\n" );
158 if( ipEmType == 1 || ipEmType == 3 )
176 xLog_line_lumin[i] = ( line_lumin > 0. ) ? log10(line_lumin) : 0.;
204 for(i=0; i< 32; i++ )
206 ipNegIntensity[i] = LONG_MAX;
215 const char chIntensityType[4][100]=
216 {
" Intrinsic" ,
" Emergent" ,
"Cumulative intrinsic" ,
"Cumulative emergent" };
217 ASSERT( ipEmType==0 || ipEmType==1 || ipEmType==2 || ipEmType==3 );
228 fprintf(
ioQQQ,
" Caution: emergent intensities are not reliable on the "
229 "first iteration.\n");
249 sclsav[iprnt] = scaled[i];
258 ipNegIntensity[nNegIntenLines] = i;
259 nNegIntenLines =
MIN2(32,nNegIntenLines+1);
265 xLog_line_lumin[iprnt] = 0.;
281 xLog_line_lumin[i] = 0.;
286 sclsav[i] = scaled[i];
291 ipNegIntensity[nNegIntenLines] = i;
292 nNegIntenLines =
MIN2(32,nNegIntenLines+1);
309 for( i=0; i<iprnt; ++i )
316 sclsav[j] = sclsav[i];
319 xLog_line_lumin[j] = xLog_line_lumin[i];
320 Slines[j] = Slines[i];
329 for( i=0; i<iprnt; ++i )
376 for( i=0; i<iprnt; ++i )
395 nline = (iprnt + 2)/3;
401 nline = (iprnt + 3)/4;
411 for( i=0; i < nline; i++ )
416 for( j=i; j<iprnt; j = j + nline)
419 long ipLin = ipSortLines[j];
452 if( sclsav[ipLin] < 9999.5 )
456 else if( sclsav[ipLin] < 99999.5 )
460 else if( sclsav[ipLin] < 999999.5 )
464 else if( sclsav[ipLin] < 9999999.5 )
475 if( j+nline < iprnt )
482 if( nNegIntenLines > 0 )
484 fprintf(
ioQQQ,
" Lines with negative intensities - Linear intensities relative to normalize line.\n" );
487 for( i=0; i < nNegIntenLines; ++i )
492 scaled[ipNegIntensity[i]] );
506 ipNegIntensity[j] = i;
516 for( i=0; i < j; i++ )
533 ipNegIntensity[j] = i;
542 for( i=0; i < j; i++ )
629 int repeat = (72-len)/2;
630 for( i=0; i < repeat; ++i )
633 for( i=0; i < 72-repeat-len; ++i )
651 if( (strcmp(chCKey,
"CONT")!= 0) && !
nMatch(
"HIDE" , chCard) )
666 " *********************************> Log(U):%6.2f <*********************************\n",
672 "\n This is a beta test version of the code, and is intended for testing only.\n\n" );
680 fprintf(
ioQQQ,
" >>>>>>>>>> Warning! Warnings exist, this calculation has serious problems.\n" );
688 " >>>>>>>>>> Cautions are present.\n" );
694 " >>>>>>>>>> The calculation stopped for unintended reasons.\n" );
700 " >>>>>>>>>> Another iteration is needed. Use the ITERATE command.\n" );
706 strcpy( chGeo,
"Closed" );
710 strcpy( chGeo,
" Open" );
716 strcpy( chPlaw,
" Constant Pressure " );
721 strcpy( chPlaw,
" Constant Density " );
727 strcpy( chPlaw,
" Power Law Density " );
732 strcpy( chPlaw,
" Rapid Fluctuations " );
737 strcpy( chPlaw,
" Special Density Lw " );
742 strcpy( chPlaw,
" Hydrostatic Equlib " );
747 strcpy( chPlaw,
" Radia Driven Wind " );
752 strcpy( chPlaw,
" Dynamical Flow " );
757 strcpy( chPlaw,
" Globule " );
762 strcpy( chPlaw,
" UNRECOGNIZED CPRES " );
766 "\n Emission Line Spectrum. %20.20sModel. %6.6s geometry. Iteration %ld of %ld.\n",
776 chUnit =
"/arcsec^2";
780 sprintf( chLine,
"Flux observed at the Earth (erg/s/cm^2%s).", chUnit );
791 sprintf( chLine,
"Surface brightness (erg/s/cm^2/%s).", chUnit );
802 chUnit =
"erg/s/cm^2";
808 sprintf( chCoverage,
"with full coverage" );
810 sprintf( chCoverage,
"with a covering factor of %.1f%%",
geometry.
covgeo*100. );
813 sprintf( chLine,
"Luminosity (%s) emitted by a partial cylinder %s.", chUnit, chCoverage );
815 sprintf( chLine,
"Luminosity (%s) emitted by a shell %s.", chUnit, chCoverage );
822 chAper =
"long slit";
824 chAper =
"pencil beam";
828 sprintf( chLine,
"Observed through a %s with aperture covering factor %.1f%%.",
836 sprintf( chLine,
"Intensity (erg/s/cm^2)." );
849 static const int NWLH = 21;
851 realnum wlh[NWLH] = { 6562.85e0f, 4861.36e0f, 4340.49e0f, 4101.76e0f, 1.87511e4f, 1.28181e4f,
852 1.09381e4f, 1.00494e4f, 2.62515e4f, 2.16553e4f, 1.94456e4f, 7.45781e4f,
853 4.6525e4f, 3.73953e4f, 4.05116e4f, 7.45781e4f, 3.29609e4f, 1215.68e0f,
854 1025.73e0f, 972.543e0f, 949.749e0f };
870 for( j=0; j<NWLH; ++j )
872 if( fabs(
LineSave.
lines[i].wavelength()-wlh[j] ) <= errorwave )
885 static const int NWLHE = 20;
886 realnum wlhe[NWLHE] = {1640.e0f, 1215.23e0f, 1085.03e0f, 1025.35e0f, 4686.01e0f, 3203.30e0f,
887 2733.46e0f, 2511.36e0f, 1.01242e4f, 6560.44e0f, 5411.80e0f, 4859.57e0f,
888 1.86377e4f, 1.16270e4f, 9345.37e0f, 8237.17e0f, 303.808e0f, 256.338e0f,
889 243.046e0f, 237.350e0f};
900 for( j=0; j<NWLHE; ++j )
902 if( fabs(
LineSave.
lines[i].wavelength()-wlhe[j] ) <= errorwave )
923 fprintf(
ioQQQ,
" PrtFinal could not find H 1 4861.36 with cdLine.\n" );
943 fprintf(
ioQQQ,
" PrtFinal could not find He 1 5876 with cdLine.\n" );
957 fprintf(
ioQQQ,
" PrtFinal could not find He 2 4686 with cdLine.\n" );
971 heabun = (he4686*0.078 + he5876*0.739)/hbeta;
992 if(
cdLine(
"blnd",5007,&o5007,&absint)<=0 )
997 fprintf(
ioQQQ,
" PrtFinal could not find O 3 5007 with cdLine.\n" );
1002 if(
cdLine(
"blnd",4363,&o4363,&absint)<=0 )
1007 fprintf(
ioQQQ,
" PrtFinal could not find H 1 4363 with cdLine.\n" );
1015 if( (o4363 != 0.) && (o5007 != 0.) )
1019 bot = o5007 - o4363/o3enro;
1043 double a31 = 0.2262;
1053 ratio = ratio/1.3333/(o3cs13/o3cs12);
1055 ratio = ratio/o3enro/o3br32;
1068 if(
cdLine(
"Bac ",3646.,&bacobs,&absint)<=0 )
1070 fprintf(
ioQQQ,
" PrtFinal could not find Bac 3646 with cdLine.\n" );
1109 x = 1.e0/log10(bac);
1110 tel = -4.827969400906710 + x*(33.08493347410885 + x*(-56.08886262205931 +
1111 x*(52.44759454803217 + x*(-25.07958990094209 + x*(4.815046760060006)))));
1113 if( tel > 1. && tel < 5. )
1130 if(
cdLine(
"thin",3646.,&bacthn,&absint)<=0 )
1132 fprintf(
ioQQQ,
" PrtFinal could not find thin 3646 with cdLine.\n" );
1139 fprintf(
ioQQQ,
" PrtFinal could not find Ca B 4861.36 with cdLine.\n" );
1177 x = 1.e0/log10(bacthn);
1178 tel = -4.827969400906710 + x*(33.08493347410885 + x*(-56.08886262205931 +
1179 x*(52.44759454803217 + x*(-25.07958990094209 + x*(4.815046760060006)))));
1181 if( tel > 1. && tel < 5. )
1311 "\n IONIZE PARMET: U(1-)%8.4f U(4-):%8.4f U(sp):%6.2f "
1312 "Q(ion): %7.3f L(ion)%7.3f Q(low):%7.3f L(low)%7.3f\n",
1313 uhl, uhel, usp, qion, pion, qlow, plow );
1356 " ENERGY BUDGET: Heat:%8.3f Coolg:%8.3f Error:%5.1f%% Rec Lin:%8.3f F-F H%7.3f P(rad/tot)max ",
1467 chUnit =
"(gm/cm^2)";
1476 if(
cdTemp(
"opti",0,&THI,
"volume" ) )
1478 fprintf(
ioQQQ,
"\n>>>>>>>>>>>>>>>>> PrtFinal, impossible problem getting 21cm opt.\n");
1488 if(
cdTemp(
"21cm",0,&THI,
"radius" ) )
1490 fprintf(
ioQQQ,
"\n>>>>>>>>>>>>>>>>> PrtFinal, impossible problem getting 21cm temp.\n");
1498 if(
cdTemp(
"spin",0,&THI,
"radius" ) )
1500 fprintf(
ioQQQ,
"\n>>>>>>>>>>>>>>>>> PrtFinal, impossible problem getting 21cm spin.\n");
1583 rjeans =
exp10(rjeans);
1589 ajmass =
exp10(ajmass);
1606 ajmin =
exp10(ajmin);
1628 if( bottom > 1e-30 && top > 1e-30 )
1630 ratio = log10(top) - log10(bottom);
1631 if( ratio < 36. && ratio > -36. )
1633 ratio =
exp10(ratio);
1650 alfox = log(ratio)/log(freq_ratio);
1711 " Hatom level%3ld HatomType:%4.4s HInducImp %2c"
1712 " He tot level:%3ld He2 level: %3ld ExecTime%8.1f\n",
1725 " ConvrgError(%%) <eden>%7.3f MaxEden%7.3f <H-C>%7.2f Max(H-C)%8.2f <Press>%8.3f MaxPrs er%7.3f\n",
1734 " Continuity(%%) chng Te%6.1f elec den%6.1f n(H2)%7.1f n(CO) %7.1f\n",
1742 fprintf(
ioQQQ,
" Te Te(Ne) Te(NeNp) Te(NeHe+)Te(NeHe2+) Te(NeO+) Te(NeO2+)"
1743 " Te(H2) N(H) Ne(O2+) Ne(Np)\n" );
1744 static const char* weight[3] = {
"Radius",
"Area",
"Volume" };
1746 for(
int dim=0; dim < dmax; ++dim )
1807 fprintf(
ioQQQ,
" Be careful: grains exist. This spectrum was not corrected for reddening before analysis.\n" );
1817 double total_dust2gas = 0.;
1819 fprintf(
ioQQQ,
"\n Average Grain Properties (over radius):\n" );
1821 for(
size_t i0=0; i0 <
gv.
bin.size(); i0 += 10 )
1827 size_t i1 =
min(i0+10,
gv.
bin.size());
1830 for(
size_t nd=i0; nd < i1; nd++ )
1832 chQHeat = (char)(
gv.
bin[nd]->lgEverQHeat ?
'*' :
' ');
1838 for(
size_t nd=i0; nd < i1; nd++ )
1846 for(
size_t nd=i0; nd < i1; nd++ )
1854 for(
size_t nd=i0; nd < i1; nd++ )
1862 for(
size_t nd=i0; nd < i1; nd++ )
1870 for(
size_t nd=i0; nd < i1; nd++ )
1954 fprintf(
ioQQQ,
" Too many comments have been entered into line array; increase the value of NHOLDCOMMENTS.\n" );
1968 for( i=0; i<6; ++i )
1999 for( i=0; i <
nzone; i++ )
2019 for( i=0; i <
nzone; i++ )
2048 for( i=0; i <
nzone; i++ )
2069 for( i=0; i <
nzone; i++ )
void PrintE93(FILE *, double)
realnum WavlenErrorGet(realnum wavelength, long sig_figs)
NORETURN void TotalInsanity(void)
double widflx(size_t i) const
realnum UV_Cont_rel2_Draine_DB96_face
bool lgPrintLineCumulative
STATIC void PrintCenterLine(FILE *io, const char chLine[], size_t ArrLen, size_t LineLen)
realnum fine_opac_velocity_width
long nMatch(const char *chKey, const char *chCard)
STATIC void gett2(realnum *tsqr)
long int nTotalIoniz_start
TransitionList HFLines("HFLines",&AnonStates)
int cdTemp(const char *chLabel, long int IonStage, double *TeMean, const char *chWeight)
sys_float sexp(sys_float x)
char chVersion[INPUT_LINE_LENGTH]
bool lgHCaseBOK[2][HS_NZ]
realnum AverHeatCoolError
void PrtMeanIon(char chType, bool lgDensity, FILE *)
void cap4(char *chCAP, const char *chLab)
double anu(size_t i) const
static t_version & Inst()
t_iso_sp iso_sp[NISO][LIMELM]
bool lgPrintBlockIntrinsic
bool fp_equal(sys_float x, sys_float y, int n=3)
bool lgdBaseSourceExists[LIMELM][LIMELM+1]
void PrintE71(FILE *, double)
long ipoint(double energy_ryd)
multi_arr< double, 2 > TempEdenMean
bool lgPrintBlockEmergent
STATIC void gett2o3(realnum *tsqr)
const int INPUT_LINE_LENGTH
multi_arr< double, 4 > TempIonEdenMean
bool lgSurfaceBrightness_SR
diatomics h2("h2", 4100.,&hmi.H2_total, Yan_H2_CS)
realnum wlAirVac(double wlAir)
sys_float safe_div(sys_float x, sys_float y, sys_float res_0by0)
STATIC void PrintSpectrum()
void PrtColumns(FILE *ioMEAN)
realnum UV_Cont_rel2_Habing_TH85_face
realnum gas_phase[LIMELM]
double totlin(int chInfo)
void sprt_wl(char *chString, realnum wl)
long int cdLine(const char *chLabel, realnum wavelength, double *relint, double *absint)
double extin_mag_V_extended
multi_arr< double, 4 > TempIonMean
#define DEBUG_ENTRY(funcname)
void PrintRatio(double q1, double q2)
static vector< realnum > wavelength
multi_arr< double, 2 > TempMean
void DatabasePrintReference()
static const int NHOLDCOMMENTS
int fprintf(const Output &stream, const char *format,...)
multi_arr< double, 4 > xIonMean
sys_float SDIV(sys_float x)
char chHoldComments[NHOLDCOMMENTS][INPUT_LINE_LENGTH]
void PrintE82(FILE *, double)
bool lgSortLineWavelength
void spsort(realnum x[], long int n, long int iperm[], int kflag, int *ier)
long int StuffComment(const char *chComment)