cloudy: prt_lines_continuum.cpp Source File

00001 /* This file is part of Cloudy and is copyright (C)1978-2008 by Gary J. Ferland and
00002  * others.  For conditions of distribution and use see copyright notice in license.txt */
00003 /*lines_continuum put energetics, H, and He lines into line intensity stack */
00004 #include "cddefines.h"
00005 #include "taulines.h"
00006 #include "physconst.h"
00007 #include "iso.h"
00008 #include "geometry.h"
00009 #include "dense.h"
00010 #include "prt.h"
00011 #include "opacity.h"
00012 #include "coolheavy.h"
00013 #include "phycon.h"
00014 #include "rfield.h"
00015 #include "predcont.h"
00016 #include "lines_service.h"
00017 #include "radius.h"
00018 #include "continuum.h"
00019 #include "lines.h"
00020 
00021 void lines_continuum(void)
00022 {
00023 
00024         double f1, 
00025                 f2 , 
00026                 bac , 
00027                 flow;
00028         long i,nBand;
00029 
00030         DEBUG_ENTRY( "lines_continuum()" );
00031 
00032         /* code has all local emissivities zeroed out with cryptic comment about being
00033          * situation dependent.  Why?  this is option to turn back on */
00034         const bool KILL_CONT = false;
00035 
00036         i = StuffComment( "continua" );
00037         linadd( 0., (realnum)i , "####", 'i',
00038                 " start continua");
00039 
00040         /***********************************************************************
00041          * stuff in Bac ratio - continuum above the Balmer Jump 
00042          * this is trick, zeroing out saved continuum integrated so far,
00043          * and adding the current version, so that the line array gives the 
00044          * value in the final continuum 
00045          *
00046          * reflected continuum is different from others since relative to inner
00047          * radius, others for for this radius 
00048          *************************************************************************/
00049 
00051         /***************************************************************************
00052          * "Bac " , 3646,  this is residual continuum at peak of Balmer Jump
00053          * flux below - flux above                                   
00054          ***************************************************************************/
00055         /* >>chng 00 dec 02, remove opac.tmn */
00056         /* >>chng 00 dec 19, remove / radius.GeoDil */
00057         /* extrapolated continuum above head */
00058         /* >>chng 01 jul 13, from ConInterOut to ConEmitOut */
00059         f1 = (rfield.ConEmitOut[iso.ipIsoLevNIonCon[ipH_LIKE][ipHYDROGEN][2]-1] + 
00060                 rfield.ConEmitReflec[iso.ipIsoLevNIonCon[ipH_LIKE][ipHYDROGEN][2]-1]/radius.r1r0sq )/
00061           rfield.widflx[iso.ipIsoLevNIonCon[ipH_LIKE][ipHYDROGEN][2]-1]/*/radius.GeoDil /opac.tmn[iso.ipIsoLevNIonCon[ipH_LIKE][ipHYDROGEN][2]-1]*/;
00062 
00063         /* extrapolated continuum below head */
00064         /* >>chng 00 dec 19, remove / radius.GeoDil */
00065         f2 = (rfield.ConEmitOut[iso.ipIsoLevNIonCon[ipH_LIKE][ipHYDROGEN][2]-2]+ 
00066                 rfield.ConEmitReflec[iso.ipIsoLevNIonCon[ipH_LIKE][ipHYDROGEN][2]-2]/radius.r1r0sq )/
00067           rfield.widflx[iso.ipIsoLevNIonCon[ipH_LIKE][ipHYDROGEN][2]-2]/*/radius.GeoDil /opac.tmn[iso.ipIsoLevNIonCon[ipH_LIKE][ipHYDROGEN][2]-2]*/;
00068 
00069         /* convert to nuFnu units */
00070         f1 = f1*0.250*0.250*EN1RYD*radius.r1r0sq;
00071         f2 = f2*0.250*0.250*EN1RYD*radius.r1r0sq;
00072         bac = (f1 - f2);
00073 
00074         /* memory not allocated until ipass >= 0 
00075          * clear summed intrinsic and emergent intensity of following
00076          * entry - following call to linadd will enter the total and
00077          * keep entering the total but is done for each zone hence need to
00078          * keep resetting to zero*/
00079         if( LineSave.ipass > 0 )
00080         {
00081                 LineSv[LineSave.nsum].sumlin[0] = 0.;
00082                 LineSv[LineSave.nsum].sumlin[1] = 0.;
00083         }
00084 
00085         linadd(MAX2(0.,bac)/radius.dVeff,3646,"Bac ",'i',
00086                 "residual flux at head of Balmer continuum, nuFnu ");
00087         /* >>chng 03 feb 06, set to zero */
00088         /* emslin saves the per unit vol emissivity of a line, which is normally 
00089          * what goes into linadd.  We zero this unit emissivity which was set
00090          * FOR THE PREVIOUS LINE since it is so situation dependent */
00091         if( KILL_CONT && LineSave.ipass > 0 )
00092         {
00093                 LineSv[LineSave.nsum-1].emslin[0] = 0.;
00094                 LineSv[LineSave.nsum-1].emslin[1] = 0.;
00095         }
00096 
00097         /* memory not allocated until ipass >= 0 */
00098         if( LineSave.ipass > 0 )
00099         {
00100                 LineSv[LineSave.nsum].sumlin[0] = 0.;
00101                 LineSv[LineSave.nsum].sumlin[1] = 0.;
00102         }
00103 
00104         linadd(f1/radius.dVeff,3645,"nFnu",'i',
00105                 "total flux above head of Balmer continuum, nuFnu ");
00106         /* >>chng 03 feb 06, set to zero */
00107         /* emslin saves the per unit vol emissivity of a line, which is normally 
00108          * what goes into linadd.  We zero this unit emissivity which was set
00109          * FOR THE PREVIOUS LINE since it is so situation dependent */
00110         if( KILL_CONT && LineSave.ipass > 0 )
00111         {
00112                 LineSv[LineSave.nsum-1].emslin[0] = 0.;
00113                 LineSv[LineSave.nsum-1].emslin[1] = 0.;
00114         }
00115 
00116         /* memory not allocated until ipass >= 0 */
00117         if( LineSave.ipass > 0 )
00118         {
00119                 LineSv[LineSave.nsum].sumlin[0] = 0.;
00120                 LineSv[LineSave.nsum].sumlin[1] = 0.;
00121         }
00122 
00123         linadd(f2/radius.dVeff,3647,"nFnu",'i',
00124                 "total flux above head of Balmer continuum, nuFnu ");
00125         /* >>chng 03 feb 06, set to zero */
00126         /* emslin saves the per unit vol emissivity of a line, which is normally 
00127          * what goes into linadd.  We zero this unit emissivity which was set
00128          * FOR THE PREVIOUS LINE since it is so situation dependent */
00129         if( KILL_CONT && LineSave.ipass > 0 )
00130         {
00131                 LineSv[LineSave.nsum-1].emslin[0] = 0.;
00132                 LineSv[LineSave.nsum-1].emslin[1] = 0.;
00133         }
00134 
00135         /******************************************************************************
00136          * "cout" , 3646,  this is outward residual continuum at peak of Balmer Jump  *
00137          * equal to total in spherical geometry, half in opt thin open geometry       *
00138          ******************************************************************************/
00139         /* >>chng 00 dec 02, remove opac.tmn */
00140         /* >>chng 00 dec 19, remove / radius.GeoDil */
00141         f1 = rfield.ConEmitOut[iso.ipIsoLevNIonCon[ipH_LIKE][ipHYDROGEN][2]-1]/
00142                 rfield.widflx[iso.ipIsoLevNIonCon[ipH_LIKE][ipHYDROGEN][2]-1]/*/radius.GeoDil /opac.tmn[iso.ipIsoLevNIonCon[ipH_LIKE][ipHYDROGEN][2]-1]*/;
00143 
00144         /* >>chng 00 dec 19, remove / radius.GeoDil */
00145         f2 = rfield.ConEmitOut[iso.ipIsoLevNIonCon[ipH_LIKE][ipHYDROGEN][2]-2]/
00146                 rfield.widflx[iso.ipIsoLevNIonCon[ipH_LIKE][ipHYDROGEN][2]-2]/*/radius.GeoDil /opac.tmn[iso.ipIsoLevNIonCon[ipH_LIKE][ipHYDROGEN][2]-2]*/;
00147 
00148         /* net Balmer jump */
00149         bac = (f1 - f2)*0.250*0.250*EN1RYD*radius.r1r0sq;
00150 
00151         /* memory not allocated until ipass >= 0 */
00152         if( LineSave.ipass > 0 )
00153         {
00154                 LineSv[LineSave.nsum].sumlin[0] = 0.;
00155                 LineSv[LineSave.nsum].sumlin[1] = 0.;
00156         }
00157 
00158         linadd(MAX2(0.,bac)/radius.dVeff,3646,"cout",'i',
00159                 "residual flux in Balmer continuum, nuFnu ");
00160         /* >>chng 03 feb 06, set to zero */
00161         /* emslin saves the per unit vol emissivity of a line, which is normally 
00162          * what goes into linadd.  We zero this unit emissivity which was set
00163          * FOR THE PREVIOUS LINE since it is so situation dependent */
00164         if( KILL_CONT && LineSave.ipass > 0 )
00165         {
00166                 LineSv[LineSave.nsum-1].emslin[0] = 0.;
00167                 LineSv[LineSave.nsum-1].emslin[1] = 0.;
00168         }
00169 
00170         /*********************************************************************
00171          * "cref" , 3646,  this is reflected continuum at peak of Balmer Jump*
00172          * equal to zero in spherical geometry, half of total in op thin opn *
00173          *********************************************************************/
00174         /* >>chng 00 dec 02, remove opac.tmn */
00175         /* >>chng 00 dec 19, remove / radius.GeoDil */
00176         f1 = rfield.ConEmitReflec[iso.ipIsoLevNIonCon[ipH_LIKE][ipHYDROGEN][2]-1]/
00177                 rfield.widflx[iso.ipIsoLevNIonCon[ipH_LIKE][ipHYDROGEN][2]-1];
00178 
00179         f2 = rfield.ConEmitReflec[iso.ipIsoLevNIonCon[ipH_LIKE][ipHYDROGEN][2]-2]/
00180                 rfield.widflx[iso.ipIsoLevNIonCon[ipH_LIKE][ipHYDROGEN][2]-2];
00181 
00182         /* net Balmer jump */
00183         bac = (f1 - f2)*0.250*0.250*EN1RYD;
00184 
00185         /* memory not allocated until ipass >= 0 */
00186         if( LineSave.ipass > 0 )
00187         {
00188                 LineSv[LineSave.nsum].sumlin[0] = 0.;
00189                 LineSv[LineSave.nsum].sumlin[1] = 0.;
00190         }
00191 
00192         linadd(MAX2(0.,bac)/radius.dVeff,3646,"cref",'i',
00193                 "residual flux in Balmer continuum, nuFnu ");
00194         /* >>chng 03 feb 06, set to zero */
00195         /* emslin saves the per unit vol emissivity of a line, which is normally 
00196          * what goes into linadd.  We zero this unit emissivity which was set
00197          * FOR THE PREVIOUS LINE since it is so situation dependent */
00198         if( KILL_CONT && LineSave.ipass > 0 )
00199         {
00200                 LineSv[LineSave.nsum-1].emslin[0] = 0.;
00201                 LineSv[LineSave.nsum-1].emslin[1] = 0.;
00202         }
00203 
00204         /*********************************************************************
00205          * "thin" , 3646, tot optically thin continuum at peak of Balmer Jump*/
00206         if( nzone > 0 )
00207         {
00208                 /* rfield.ConEmitLocal is not defined initially, only evaluate when into model */
00209                 f1 = rfield.ConEmitLocal[nzone][iso.ipIsoLevNIonCon[ipH_LIKE][ipHYDROGEN][2]-1]/
00210                         rfield.widflx[iso.ipIsoLevNIonCon[ipH_LIKE][ipHYDROGEN][2]-1];
00211                 f2 = rfield.ConEmitLocal[nzone][iso.ipIsoLevNIonCon[ipH_LIKE][ipHYDROGEN][2]-2]/
00212                         rfield.widflx[iso.ipIsoLevNIonCon[ipH_LIKE][ipHYDROGEN][2]-2];
00213                 bac = (f1 - f2)*0.250*0.250*EN1RYD;
00214         }
00215         else
00216         {
00217                 f1 = 0.;
00218                 f2 = 0.;
00219                 bac = 0.;
00220         }
00221 
00222         linadd(MAX2(0.,bac),3646,"thin",'i',
00223                 "residual flux in Balmer continuum, nuFnu ");
00224 
00225         continuum.cn4861 /= (realnum)radius.dVeff;
00226         linadd(continuum.cn4861,4860,"Inci",'i',
00227                 "incident continuum nu*f_nu at H-beta, at illuminated face of cloud ");
00228         /* >>chng 07 jun 13, at this point nsum is the number of lines in the stack
00229          * so nsum-1 is the entry created by the above call to linadd
00230          * we want to set the emergent incident continuum to the intrinsic
00231          * incident continuum - there are reports of the incident continuum 
00232          * striking the cloud and the emergent - incident distinction does
00233          * not apply - code had left emergent intensity at zero, bug reported
00234          * by K Korista on discussion group 2007 jun 14 */
00235         if( LineSave.ipass > 0 )
00236         {
00237                 LineSv[LineSave.nsum-1].sumlin[1] = LineSv[LineSave.nsum-1].sumlin[0];
00238         }
00239 
00240         continuum.cn1216 /= (realnum)radius.dVeff;
00241         linadd(continuum.cn1216,1215,"Inci",'i',
00242                 "incident continuum nu*f_nu near Ly-alpha, at illuminated face of cloud");
00243         /* see comment concerning parallel code immediately above */
00244         if( LineSave.ipass > 0 )
00245         {
00246                 LineSv[LineSave.nsum-1].sumlin[1] = LineSv[LineSave.nsum-1].sumlin[0];
00247         }
00248 
00249         continuum.cn1216 *= (realnum)radius.dVeff;
00250         continuum.cn4861 *= (realnum)radius.dVeff;
00251 
00252         if( LineSave.ipass > 0 )
00253         {
00254                 continuum.cn4861 = 0.;
00255                 continuum.cn1216 = 0.;
00256         }
00257 
00258         flow = (iso.RadRecomb[ipH_LIKE][ipHYDROGEN][ipH2p][ipRecRad] + 
00259                 iso.RadRecomb[ipH_LIKE][ipHYDROGEN][ipH2s][ipRecRad])*
00260           iso.RadRecomb[ipH_LIKE][ipHYDROGEN][ipH2p][ipRecEsc]*
00261           dense.eden*dense.xIonDense[ipHYDROGEN][1]* 5.45e-12;
00262         linadd(flow,0,"Ba C",'i',
00263                 "integrated Balmer continuum emission");
00264 
00265         if( iso.n_HighestResolved_max[ipH_LIKE][ipHYDROGEN] >= 3 )
00266         {
00267                 flow = ( iso.RadRecomb[ipH_LIKE][ipHYDROGEN][ipH3s][ipRecRad]*
00268                         iso.RadRecomb[ipH_LIKE][ipHYDROGEN][ipH3s][ipRecEsc] +
00269                         iso.RadRecomb[ipH_LIKE][ipHYDROGEN][ipH3p][ipRecRad]*
00270                         iso.RadRecomb[ipH_LIKE][ipHYDROGEN][ipH3p][ipRecEsc] +
00271                         iso.RadRecomb[ipH_LIKE][ipHYDROGEN][ipH3d][ipRecRad]*
00272                         iso.RadRecomb[ipH_LIKE][ipHYDROGEN][ipH3d][ipRecEsc] ) * 
00273                         dense.eden*dense.xIonDense[ipHYDROGEN][1]*3.53e-12;
00274         }
00275         else
00276         {
00277                 flow = iso.RadRecomb[ipH_LIKE][ipHYDROGEN][3][ipRecRad]*
00278                         iso.RadRecomb[ipH_LIKE][ipHYDROGEN][3][ipRecEsc]*
00279                   dense.eden*dense.xIonDense[ipHYDROGEN][1]*3.53e-12;
00280         }
00281         linadd(flow,0,"PA C",'i',
00282                 "Paschen continuum emission ");
00283 
00284         /* >>chng 05 oct 16, add this output option, these are a
00285          * series of continuum bands defined in the file bands_continuum.dat
00286          * this makes it possible to enter any integrated total emission
00287          * into the emission-line stack */
00288         for( nBand=0; nBand<continuum.nContBand; ++nBand )
00289         {
00290                 double total = 0.;
00291                 /* find total emission over band - units will be erg cm-2 s-1 */
00292                 for( i=continuum.ipContBandLow[nBand]; i<continuum.ipContBandHi[nBand]; ++i )
00293                 {
00294                         double xIntenOut = 
00295                                 /* the attenuated incident continuum */
00296                                 rfield.flux[i-1] +
00297 
00298                                 /* the outward emitted continuum radiation field 
00299                                  * >>chng 06 aug 07, add geomery.covgeo here */
00300                                 (rfield.ConEmitOut[i-1] + 
00301 
00302                                 /* outward emitted lines */
00303                                 rfield.outlin[i-1])*geometry.covgeo;
00304                         /* div by opac.E2TauAbsFace[i] because ConEmitReflec has already
00305                          * been corrected for this - emergent_line will introduce a 
00306                          * further correction so this will cancel the extra factor */
00307                         double xIntenIn = 
00308                                 ((double)rfield.ConEmitReflec[i-1]/SDIV((double)opac.E2TauAbsFace[i-1]) + 
00309                                  /* outward emitted lines */
00310                                  rfield.reflin[i-1] )*geometry.covgeo;
00311 
00312                         /* the fraction of this that gets out */
00313                         total += rfield.anu[i-1] * 
00314                                 emergent_line( xIntenIn , xIntenOut , i )
00315                                 /* >>chng 06 aug 07, add opac.tmn here */
00316                                 / SDIV(opac.tmn[i-1]);
00317                 }
00318                 /* we will call lindst with an energy half way between the two
00319                  * ends of the band.  This will make an extinction correction that
00320                  * has already been applied above by multiplying by emergent_line.
00321                  * Find this factor and remove it before the call */
00322                 double corr = emergent_line( 0.5 , 0.5 , 
00323                         (continuum.ipContBandLow[nBand]+continuum.ipContBandHi[nBand])/2 );
00324                 if( corr < SMALLFLOAT )
00325                         total = 0.;
00326                 else
00327                         total /= corr;
00328 
00329                 /* convert to physical units */
00330                 total *= EN1RYD*radius.r1r0sq/radius.dVeff;
00331                 /* memory not allocated until ipass >= 0 */
00332                 if( LineSave.ipass > 0 )
00333                 {
00334                         LineSv[LineSave.nsum].sumlin[0] = 0.;
00335                         LineSv[LineSave.nsum].sumlin[1] = 0.;
00336                 }
00337 
00338                 lindst( total , 
00339                         continuum.ContBandWavelength[nBand] , 
00340                         continuum.chContBandLabels[nBand],
00341                         (continuum.ipContBandLow[nBand]+continuum.ipContBandHi[nBand])/2 ,
00342                         'i' , false,
00343                         "continuum bands defined in bands_continuum.dat");
00344         }
00345 
00346         linadd(MAX2(0.,CoolHeavy.brems_cool_net),0,"HFFc",'c',
00347                 "net free-free cooling, ALL species, free-free heating subtracted, so nearly cancels with cooling in LTE ");
00348 
00349         linadd(MAX2(0.,-CoolHeavy.brems_cool_net),0,"HFFh",'h',
00350                 "net free-free heating, nearly cancels with cooling in LTE ");
00351 
00352         linadd(CoolHeavy.brems_cool_h,0,"H FF",'i',
00353                 " H brems (free-free) cooling ");
00354 
00355         linadd(CoolHeavy.brems_heat_total,0,"FF H",'i',
00356                 "total free-free heating ");
00357 
00358         linadd(CoolHeavy.brems_cool_he,0,"HeFF",'i',
00359                 "He brems emission ");
00360 
00361         linadd(CoolHeavy.herec,0,"HeFB",'c',
00362                 "He recombination cooling ");
00363 
00364         linadd(CoolHeavy.heavfb,0,"MeFB",'c',
00365                 "heavy element recombination cooling ");
00366 
00367         linadd(CoolHeavy.brems_cool_metals,0,"MeFF",'i',
00368                 "heavy elements (metals) brems cooling, heat not subtracted ");
00369 
00370         linadd(CoolHeavy.brems_cool_h+CoolHeavy.brems_cool_he+CoolHeavy.brems_cool_metals,0,"ToFF",'i',
00371                 "total brems emission - total cooling but not minus heating ");
00372 
00373         linadd((CoolHeavy.brems_cool_h+CoolHeavy.brems_cool_he)*sexp(5.8e6/phycon.te),0,"FF X",'i',
00374                 "part of H brems, in x-ray beyond 0.5KeV ");
00375 
00376         linadd(CoolHeavy.eebrm,0,"eeff",'c',
00377                 "electron - electron brems ");
00378 
00379         /* predict emitted continuum at series of continuum points */
00380         /* variables are located in predcont.h, 
00381          * NPREDCONT number of continuum points,
00382          * next two arrays are defined in zerologic.c
00383          * EnrPredCont - energy (Ryd) where we want to predict the continuum,
00384          *               these are set in zerologic.c,
00385          * ipPreCont -   pointers to energy within continuum array, 
00386          *               these are set in ContCreatePointers
00387          *
00388          * the array of continuum energies where this will be done is
00389          * stored as EnrPredCont which lives in predcont.h
00390          * and is initialized in zerologic.c
00391          *
00392          * the entry nFnu will only be printed if the command
00393          * print diffuse continuum
00394          * is entered - 
00395          *
00396          * this code should be kept parallel with that in dopunch, where
00397          * punch continuum is produced, since two must agree */
00398         for( i=0; i < NPREDCONT; i++ )
00399         {
00400                 double SourceTransmitted , Cont_nInu;
00401                 double SourceReflected, 
00402                   DiffuseOutward,
00403                   DiffuseInward;
00404                 double renorm;
00405 
00406                 /* put wavelength in Angstroms into dummy structure, so that we can use iWavLen
00407                  * to get a proper wavelength with units, continuum energies are stored in EnrPredCont */
00408                 TauDummy.WLAng = (realnum)(RYDLAM/EnrPredCont[i]);
00409                 /*lambda = iWavLen(&TauDummy , &chUnits , &chShift );*/
00410 
00411                 /* >>chng 00 dec 02, there were three occurrences of /opac.tmn which had the
00412                  * effect of raising the summed continuum by the local opacity correction factor.
00413                  * in the case of the Lyman continuum this raised the reported value by orders
00414                  * of magnitude.  There have been commented out in the following for now. */
00415                 /* reflected total continuum (diff+incident emitted inward direction) */
00416 
00417                 /* >>chng 00 dec 08, implement the "set nFnu [SOURCE_REFLECTED] ... command, PvH */
00418                 /* >>chng 00 dec 19, remove / radius.GeoDil */
00419                 renorm = rfield.anu2[ipPredCont[i]]*EN1RYD/rfield.widflx[ipPredCont[i]]/*/radius.GeoDil*/;
00420 
00421                 /* this is the reflected diffuse continuum */
00422                 if( prt.lgDiffuseInward )
00423                 {
00424                         DiffuseInward = rfield.ConEmitReflec[ipPredCont[i]]*renorm;
00425                         /* /opac.tmn[ipPredCont[i]]*/
00426                 }
00427                 else
00428                 {
00429                         DiffuseInward = 0.;
00430                 }
00431 
00432                 /* the outward diffuse continuum */
00433                 if( prt.lgDiffuseOutward )
00434                 {
00435                         DiffuseOutward = rfield.ConEmitOut[ipPredCont[i]]*renorm*radius.r1r0sq;
00436                 }
00437                 else
00438                 {
00439                         DiffuseOutward = 0.;
00440                 }
00441 
00442                 /* reflected part of INCIDENT continuum (only incident, not diffuse, which was above) */
00443                 if( prt.lgSourceReflected )
00444                 {
00445                         SourceReflected =  rfield.ConRefIncid[ipPredCont[i]]*renorm;
00446                         /* /opac.tmn[ipPredCont[i]]*/
00447                 }
00448                 else
00449                 {
00450                         SourceReflected =  0.;
00451                 }
00452 
00453                 /* the attenuated incident continuum */
00454                 if( prt.lgSourceTransmitted )
00455                 {
00456                         SourceTransmitted = rfield.flux[ipPredCont[i]]*renorm*radius.r1r0sq;
00457                 }
00458                 else
00459                 {
00460                         SourceTransmitted = 0.;
00461                 }
00462 
00463                 /* memory has not been allocated until ipass >= 0, so must not access this element,
00464                  * this element will be used to save the following quantity */
00465                 if( LineSave.ipass > 0 )
00466                 {
00467                         LineSv[LineSave.nsum].sumlin[0] = 0.;
00468                         LineSv[LineSave.nsum].sumlin[1] = 0.;
00469                 }
00470 
00471                 linadd((DiffuseInward+SourceReflected+DiffuseOutward+SourceTransmitted)/radius.dVeff,
00472                         TauDummy.WLAng,"nFnu",'i',
00473                         "total continuum at selected energy points " );
00474 
00475                 /* emslin saves the per unit vol emissivity of a line, which is normally 
00476                  * what goes into linadd.  We zero this unit emissivity which was set
00477                  * FOR THE PREVIOUS LINE since it is so situation dependent */
00478                 if( KILL_CONT && LineSave.ipass > 0 )
00479                 {
00480                         LineSv[LineSave.nsum-1].emslin[0] = 0.;
00481                         LineSv[LineSave.nsum-1].emslin[1] = 0.;
00482                 }
00483 
00484                 /* this is the normal set to zero to trick the NEXT line into going in properly */
00485                 if( LineSave.ipass > 0 )
00486                 {
00487                         LineSv[LineSave.nsum].sumlin[0] = 0.;
00488                         LineSv[LineSave.nsum].sumlin[1] = 0.;
00489                 }
00490 
00491                 /* the nsum-1 -- emslin and nsum -- sumlin is not a bug, look above - they do
00492                  * different things to different saves */
00493                 Cont_nInu = rfield.flux[ipPredCont[i]]*renorm*radius.r1r0sq +
00494                         rfield.ConRefIncid[ipPredCont[i]]*renorm;
00495 
00496 #               if 0
00497                 /* this code can be used to create assert statements for the continuum shape */
00498                 if( !i )
00499                         fprintf(ioQQQ,"\n");
00500                 char chWL[1000];
00501                 sprt_wl( chWL , TauDummy.WLAng );
00502                 fprintf( ioQQQ,"assert line luminosity \"nInu\" %s  %.3f\n",
00503                         chWL , 
00504                         log10(SDIV(Cont_nInu/radius.dVeff)) + radius.Conv2PrtInten  );
00505 #               endif
00506 
00507                 linadd( Cont_nInu/radius.dVeff,TauDummy.WLAng,"nInu",'i',
00508                         "transmitted and reflected incident continuum at selected energy points " );
00509 
00510                 /* emslin saves the per unit volume emissivity of a line, which is normally 
00511                  * what goes into linadd.  We zero this unit emissivity since it is so situation dependent */
00512                 if( KILL_CONT && LineSave.ipass > 0 )
00513                 {
00514                         LineSv[LineSave.nsum-1].emslin[0] = 0.;
00515                         LineSv[LineSave.nsum-1].emslin[1] = 0.;
00516                 }
00517 
00518                 /* memory has not been allocated until ipass >= 0 */
00519                 if( LineSave.ipass > 0 )
00520                 {
00521                         LineSv[LineSave.nsum].sumlin[0] = 0.;
00522                         LineSv[LineSave.nsum].sumlin[1] = 0.;
00523                 }
00524 
00525                 linadd( (DiffuseInward+SourceReflected)/radius.dVeff,TauDummy.WLAng,"InwT",'i',
00526                         "total reflected continuum, total inward emission plus reflected (XXdiffuseXX) total continuum ");
00527 
00528                 if( KILL_CONT && LineSave.ipass > 0 )
00529                 {
00530                         LineSv[LineSave.nsum-1].emslin[0] = 0.;
00531                         LineSv[LineSave.nsum-1].emslin[1] = 0.;
00532                 }
00533 
00534                 /* memory has not been allocated until ipass >= 0 */
00535                 if( LineSave.ipass > 0 )
00536                 {
00537                         LineSv[LineSave.nsum].sumlin[0] = 0.;
00538                         LineSv[LineSave.nsum].sumlin[1] = 0.;
00539                 }
00540 
00541                 linadd(SourceReflected/radius.dVeff,TauDummy.WLAng,"InwC",'i',
00542                         "reflected incident continuum (only incident) ");
00543 
00544                 if( KILL_CONT && LineSave.ipass > 0 )
00545                 {
00546                         LineSv[LineSave.nsum-1].emslin[0] = 0.;
00547                         LineSv[LineSave.nsum-1].emslin[1] = 0.;
00548                 }
00549         }
00550         return;
00551 }