00001
00002
00003
00004
00005
00006
00007
00008
00009
00010
00011
00012
00013
00014
00015
00016
00017 #include "cddefines.h"
00018 #include "lines_service.h"
00019 #include "dense.h"
00020 #include "geometry.h"
00021 #include "hydrogenic.h"
00022 #include "ipoint.h"
00023 #include "iso.h"
00024 #include "lines.h"
00025 #include "opacity.h"
00026 #include "physconst.h"
00027 #include "radius.h"
00028 #include "rfield.h"
00029 #include "rt.h"
00030 #include "taulines.h"
00031
00032 STATIC void normalizeProfile( double *profile, long numPoints, realnum *eArray, double E_centroid );
00033
00034
00035 double eina(double gf,
00036 double enercm,
00037 double gup)
00038 {
00039 double eina_v;
00040
00041 DEBUG_ENTRY( "eina()" );
00042
00043
00044
00045
00046
00047 eina_v = (gf/gup)*TRANS_PROB_CONST*POW2(enercm);
00048 return( eina_v );
00049 }
00050
00051
00052 double GetGF(double trans_prob,
00053 double enercm,
00054 double gup)
00055 {
00056 double GetGF_v;
00057
00058 DEBUG_ENTRY( "GetGF()" );
00059
00060 ASSERT( enercm > 0. );
00061 ASSERT( trans_prob > 0. );
00062 ASSERT( gup > 0.);
00063
00064
00065
00066
00067
00068 GetGF_v = trans_prob*gup/TRANS_PROB_CONST/POW2(enercm);
00069 return( GetGF_v );
00070 }
00071
00072
00073 double abscf(double gf,
00074 double enercm,
00075 double gl)
00076 {
00077 double abscf_v;
00078
00079 DEBUG_ENTRY( "abscf()" );
00080
00081 ASSERT(gl > 0. && enercm > 0. && gf > 0. );
00082
00083
00084
00085
00086 abscf_v = 1.4974e-6*(gf/gl)*(1e4/enercm);
00087 return( abscf_v );
00088 }
00089
00090
00091
00092 double RefIndex(double EnergyWN )
00093 {
00094 double RefIndex_v,
00095 WaveMic,
00096 xl,
00097 xn;
00098
00099 DEBUG_ENTRY( "RefIndex()" );
00100
00101 ASSERT( EnergyWN > 0. );
00102
00103
00104 WaveMic = 1.e4/EnergyWN;
00105
00106
00107 if( WaveMic > 0.2 )
00108 {
00109
00110
00111 xl = 1.0/WaveMic/WaveMic;
00112
00113
00114
00115 xn = 255.4/(41. - xl);
00116 xn += 29498.1/(146. - xl);
00117 xn += 64.328;
00118 RefIndex_v = xn/1.e6 + 1.;
00119 }
00120 else
00121 {
00122 RefIndex_v = 1.;
00123 }
00124 ASSERT( RefIndex_v >= 1. );
00125 return( RefIndex_v );
00126 }
00127
00128
00129
00130
00131
00132
00133 realnum WavlenErrorGet( realnum wavelength )
00134 {
00135 double a;
00136 realnum errorwave;
00137
00138 DEBUG_ENTRY( "WavlenErrorGet()" );
00139
00140 ASSERT( LineSave.sig_figs <= 6 );
00141
00142 if( wavelength > 0. )
00143 {
00144
00145 a = log10( wavelength+FLT_EPSILON);
00146 a = floor(a);
00147 }
00148 else
00149 {
00150
00151
00152 a = 0.;
00153 }
00154
00155 errorwave = 5.f * (realnum)pow( 10., a - (double)LineSave.sig_figs );
00156 return errorwave;
00157 }
00158
00159
00160 STATIC void lincom(
00161 double xInten,
00162 realnum wavelength,
00163 const char *chLab,
00164
00165 long int ipnt,
00166 char chInfo,
00167
00168
00169 const char *chComment,
00170
00171 bool lgAdd)
00172 {
00173 DEBUG_ENTRY( "lincom()" );
00174
00175
00176
00177
00178
00179
00180
00181 if( LineSave.ipass > 0 )
00182 {
00183
00184
00185 if (lgAdd || xInten > 0.)
00186 {
00187
00188
00189
00190
00191 LineSv[LineSave.nsum].SumLine[0] += xInten*radius.dVeffAper;
00192
00193
00194 LineSv[LineSave.nsum].emslin[0] = xInten;
00195 }
00196
00197 if (lgAdd)
00198 {
00199 if (wavelength > 0 )
00200 {
00201
00202
00203
00204
00205
00206 LineSv[LineSave.nsum].emslin[1] = LineSv[LineSave.nsum].emslin[0];
00207 LineSv[LineSave.nsum].SumLine[1] = LineSv[LineSave.nsum].SumLine[0];
00208 }
00209 }
00210 else
00211 {
00212 if ( xInten > 0. && ipnt <= rfield.nflux )
00213 {
00214
00215
00216 const double saveemis = emergent_line(
00217 xInten*rt.fracin , xInten*(1.-rt.fracin) , ipnt );
00218 LineSv[LineSave.nsum].emslin[1] = saveemis;
00219 LineSv[LineSave.nsum].SumLine[1] += saveemis*radius.dVeffAper;
00220 }
00221 }
00222 }
00223
00224 else if( LineSave.ipass == 0 )
00225 {
00226
00227
00228 ASSERT( (chInfo == 'c') || (chInfo == 'h') || (chInfo == 'i') || (chInfo == 'r' ) );
00229
00230 LineSv[LineSave.nsum].chSumTyp = (char)chInfo;
00231 LineSv[LineSave.nsum].emslin[0] = 0.;
00232 LineSv[LineSave.nsum].emslin[1] = 0.;
00233 LineSv[LineSave.nsum].chComment = chComment;
00234
00235
00236 ASSERT( strlen( chLab )<5 );
00237 strcpy( LineSv[LineSave.nsum].chALab, chLab );
00238
00239 if (lgAdd)
00240 {
00241 LineSv[LineSave.nsum].wavelength = wavelength;
00242 }
00243 else
00244 {
00245
00246
00247 LineSv[LineSave.nsum].wavelength = fabs(wavelength);
00248 LineSv[LineSave.nsum].SumLine[0] = 0.;
00249 LineSv[LineSave.nsum].SumLine[1] = 0.;
00250
00251
00252
00253
00254
00255 ASSERT( ipnt > 0 );
00256 # ifndef NDEBUG
00257 double error = MAX2(0.1*rfield.AnuOrg[ipnt-1] , rfield.widflx[ipnt-1] );
00258 ASSERT( wavelength<=0 ||
00259 fabs( rfield.AnuOrg[ipnt-1] - RYDLAM / wavelength) < error );
00260 # endif
00261 }
00262 }
00263
00264
00265 ++LineSave.nsum;
00266
00267
00268
00269 }
00270
00271
00272 void linadd(
00273 double xInten,
00274 realnum wavelength,
00275 const char *chLab,
00276 char chInfo,
00277
00278 const char *chComment )
00279 {
00280 DEBUG_ENTRY( "linadd()" );
00281
00282
00283 const long int ipnt = LONG_MAX;
00284
00285 lincom( xInten, wavelength, chLab, ipnt, chInfo, chComment, true );
00286 }
00287
00288
00289
00290
00291 double emergent_line(
00292
00293 double emissivity_in ,
00294
00295 double emissivity_out ,
00296
00297 long int ipCont )
00298 {
00299
00300 double emergent_in , emergent_out;
00301 long int i = ipCont-1;
00302
00303 DEBUG_ENTRY( "emergent_line()" );
00304
00305 ASSERT( i >= 0 && i < rfield.nupper-1 );
00306
00307
00308
00309
00310
00311
00312
00313
00314
00315
00316 if( iteration == 1 )
00317 {
00318
00319
00320 emergent_in = emissivity_in*opac.E2TauAbsFace[i];
00321 emergent_out = emissivity_out;
00322 }
00323 else
00324 {
00325 if( geometry.lgSphere )
00326 {
00327
00328
00329
00330 emergent_in = emissivity_in * opac.E2TauAbsFace[i] *opac.E2TauAbsTotal[i];
00331
00332
00333 emergent_out = emissivity_out * opac.E2TauAbsOut[i];
00334 }
00335 else
00336 {
00337
00338
00339
00340 double reflected = emissivity_out * opac.albedo[i] * (1.-opac.E2TauAbsOut[i]);
00341
00342 emergent_in = (emissivity_in + reflected) * opac.E2TauAbsFace[i];
00343
00344 emergent_out = (emissivity_out - reflected) * opac.E2TauAbsOut[i];
00345 }
00346 }
00347
00348 return( emergent_in + emergent_out );
00349 }
00350
00351
00352 void outline_base(double dampXvel, double damp, bool lgTransStackLine, long int ip, double phots, realnum inwd,
00353 double nonScatteredFraction)
00354 {
00355 DEBUG_ENTRY( "outline_base()" );
00356
00357 #define DO_PROFILE false
00358
00359 if( !DO_PROFILE )
00360 outline_base_bin(lgTransStackLine, ip, phots, inwd, nonScatteredFraction);
00361 else
00362 {
00363 ASSERT( damp > 0. );
00364 double LineWidth = dampXvel/damp;
00365 LineWidth = MIN2( 0.1 * SPEEDLIGHT, LineWidth );
00366 double sigma = (LineWidth/SPEEDLIGHT);
00367 long ip3SigmaRed = ipoint( MAX2( rfield.emm, rfield.anu[ip] - 3.*sigma*rfield.anu[ip] ) );
00368 long ip3SigmaBlue = ipoint( MIN2( rfield.egamry, rfield.anu[ip] + 3.*sigma*rfield.anu[ip] ) );
00369 ASSERT( ip3SigmaBlue >= ip3SigmaRed );
00370 long numBins = ip3SigmaBlue - ip3SigmaRed + 1;
00371
00372 if( numBins < 3 )
00373 outline_base_bin(lgTransStackLine, ip, phots, inwd, nonScatteredFraction);
00374 else
00375 {
00376 valarray<double> profile(numBins);
00377
00378 for( long ipBin=ip3SigmaRed; ipBin<=ip3SigmaBlue; ipBin++ )
00379 {
00380 double x = (rfield.anu[ip] - rfield.anu[ipBin])/rfield.anu[ip]/sigma;
00381 profile[ipBin-ip3SigmaRed] = vfun( damp, x );
00382 }
00383
00384 normalizeProfile( &profile[0], numBins, rfield.anu+ip3SigmaRed, rfield.anu[ip] );
00385
00386 for( long ipBin=ip3SigmaRed; ipBin<=ip3SigmaBlue; ipBin++ )
00387 outline_base_bin(lgTransStackLine, ipBin, phots*profile[ipBin-ip3SigmaRed], inwd, nonScatteredFraction);
00388 }
00389 }
00390 }
00391
00392
00393 void outline_base_bin(bool lgTransStackLine, long int ip, double phots, realnum inwd,
00394 double nonScatteredFraction)
00395 {
00396 DEBUG_ENTRY( "outline_base_bin()" );
00397
00398 if (lgTransStackLine)
00399 {
00400 rfield.DiffuseLineEmission[ip] +=
00401 (realnum)phots;
00402
00403
00404 rfield.reflin[0][ip] +=
00405 (realnum)(inwd*phots*radius.BeamInIn);
00406
00407
00408 rfield.outlin[0][ip] +=
00409 (realnum)(inwd*phots*radius.BeamInOut*opac.tmn[ip]*nonScatteredFraction);
00410
00411
00412 rfield.outlin[0][ip] +=
00413 (realnum)((1.-inwd)*phots*radius.BeamOutOut*opac.tmn[ip]*nonScatteredFraction);
00414 }
00415 else
00416 {
00417 rfield.reflin[0][ip] +=
00418 (realnum)(phots*radius.dVolReflec);
00419
00420 rfield.outlin[0][ip] +=
00421 (realnum)(phots*radius.dVolOutwrd*opac.ExpZone[ip]);
00422 }
00423 }
00424
00425 STATIC void normalizeProfile( double *profile, long numPoints, realnum *eArray, double E_centroid )
00426 {
00427 double psum = 0.;
00428 double pEsum = 0.;
00429 ASSERT( numPoints >= 1 );
00430 if( numPoints==1 )
00431 {
00432 profile[0] = 1.;
00433 return;
00434 }
00435
00436 for( long i=0; i<numPoints; i++ )
00437 {
00438 pEsum += eArray[i] * profile[i];
00439 psum += profile[i];
00440 }
00441
00442 for( long i=0; i<numPoints; i++ )
00443 {
00444 ASSERT( pEsum > SMALLFLOAT );
00445 profile[i] *= E_centroid/pEsum;
00446 }
00447
00448 for( long i=0; i<numPoints; i++ )
00449 {
00450 pEsum += eArray[i] * profile[i];
00451 psum += profile[i];
00452 }
00453
00454 return;
00455 }
00456
00457
00458 void lindst(
00459
00460 double xInten,
00461
00462 realnum wavelength,
00463
00464 const char *chLab,
00465
00466 long int ipnt,
00467
00468 char chInfo,
00469
00470 bool lgOutToo,
00471
00472 const char *chComment )
00473 {
00474 DEBUG_ENTRY( "lindst()" );
00475
00476
00477 ASSERT( !lgOutToo || chInfo!='i' );
00478
00479 lincom(xInten, wavelength, chLab, ipnt, chInfo, chComment, false );
00480
00481 if( LineSave.ipass > 0 )
00482 {
00483
00484
00485 if (lgOutToo && xInten > 0.)
00486 {
00487
00488
00489
00490
00491
00492
00493 const bool lgTransStackLine = false;
00494 const long int ip = ipnt - 1;
00495 const double phots = xInten/(rfield.anu[ipnt-1]*EN1RYD);
00496 const realnum inwd = (realnum)(1.0-(1.+geometry.covrt)/2.);
00497 const double nonScatteredFraction = 1.;
00498
00499 outline_base_bin(lgTransStackLine, ip, phots, inwd, nonScatteredFraction);
00500 }
00501 }
00502 }
00503
00504
00505 void lindst(
00506 double dampXvel,
00507 double damp,
00508
00509 double xInten,
00510
00511 realnum wavelength,
00512
00513 const char *chLab,
00514
00515 long int ipnt,
00516
00517 char chInfo,
00518
00519 bool lgOutToo,
00520
00521 const char *chComment )
00522 {
00523 DEBUG_ENTRY( "lindst()" );
00524
00525
00526 ASSERT( !lgOutToo || chInfo!='i' );
00527
00528 lincom(xInten, wavelength, chLab, ipnt, chInfo, chComment, false );
00529
00530 if( LineSave.ipass > 0 )
00531 {
00532
00533
00534 if (lgOutToo && xInten > 0.)
00535 {
00536
00537
00538
00539
00540
00541
00542 const bool lgTransStackLine = false;
00543 const long int ip = ipnt - 1;
00544 const double phots = xInten/(rfield.anu[ipnt-1]*EN1RYD);
00545 const realnum inwd = (realnum)(1.0-(1.+geometry.covrt)/2.);
00546 const double nonScatteredFraction = 1.;
00547
00548 outline_base(dampXvel, damp, lgTransStackLine, ip, phots, inwd, nonScatteredFraction);
00549 }
00550 }
00551 }
00552
00553
00554 void lindst(
00555 transition *t,
00556
00557 const char *chLab,
00558
00559 char chInfo,
00560
00561 bool lgOutToo,
00562
00563 const char *chComment )
00564 {
00565 DEBUG_ENTRY( "lindst()" );
00566
00567 lindst( t->Emis->dampXvel, t->Emis->damp, t->Emis->xIntensity, t->WLAng, chLab, t->ipCont, chInfo,
00568 lgOutToo, chComment );
00569
00570 }
00571
00572
00573 void PntForLine(
00574
00575 double wavelength,
00576
00577 const char *chLabel,
00578
00579
00580 long int *ipnt)
00581 {
00582
00583
00584
00585
00586
00587 const int MAXFORLIN = 1000;
00588 static long int ipForLin[MAXFORLIN]={0};
00589
00590
00591 static long int nForLin;
00592
00593 DEBUG_ENTRY( "PntForLine()" );
00594
00595
00596 ASSERT( wavelength >= 0. );
00597
00598 if( wavelength == 0. )
00599 {
00600
00601 nForLin = 0;
00602 }
00603 else
00604 {
00605
00606 if( LineSave.ipass > 0 )
00607 {
00608
00609 *ipnt = ipForLin[nForLin];
00610 }
00611 else if( LineSave.ipass == 0 )
00612 {
00613
00614 if( nForLin >= MAXFORLIN )
00615 {
00616 fprintf( ioQQQ, "PROBLEM %5ld lines is too many for PntForLine.\n",
00617 nForLin );
00618 fprintf( ioQQQ, " Increase the value of maxForLine everywhere in the code.\n" );
00619 cdEXIT(EXIT_FAILURE);
00620 }
00621
00622
00623 const double EnergyRyd = RYDLAM/wavelength;
00624 ipForLin[nForLin] = ipLineEnergy(EnergyRyd,chLabel , 0);
00625 *ipnt = ipForLin[nForLin];
00626 }
00627 else
00628 {
00629
00630 *ipnt = 0;
00631 }
00632 ++nForLin;
00633 }
00634 return;
00635 }
00636
00637
00638 void EmLineJunk( emission *t )
00639 {
00640
00641 DEBUG_ENTRY( "EmLineJunk()" );
00642
00643
00644 t->TauCon = -FLT_MAX;
00645
00646
00647 t->TauIn = -FLT_MAX;
00648 t->TauTot = -FLT_MAX;
00649
00650
00651 t->iRedisFun = INT_MIN;
00652
00653
00654 t->ipFine = -10000;
00655
00656
00657 t->FracInwd = -FLT_MAX;
00658
00659
00660 t->pump = -FLT_MAX;
00661
00662
00663 t->xIntensity = -FLT_MAX;
00664
00665
00666 t->phots = -FLT_MAX;
00667
00668
00669 t->gf = -FLT_MAX;
00670
00671
00672 t->Pesc = -FLT_MAX;
00673 t->Pdest = -FLT_MAX;
00674 t->Pelec_esc = -FLT_MAX;
00675
00676
00677 t->dampXvel = -FLT_MAX;
00678 t->damp = -FLT_MAX;
00679
00680
00681 t->ColOvTot = -FLT_MAX;
00682
00683
00684 t->AutoIonizFrac = -FLT_MAX;
00685
00686
00687 t->opacity = -FLT_MAX;
00688
00689 t->PopOpc = -FLT_MAX;
00690
00691
00692 t->Aul = -FLT_MAX;
00693
00694
00695 t->ots = -FLT_MAX;
00696 return;
00697 }
00698
00699
00700 void CollisionJunk( collision * t )
00701 {
00702
00703 DEBUG_ENTRY( "CollisionJunk()" );
00704
00706 t->ColUL = -FLT_MAX;
00707
00708
00709 t->cool = -FLT_MAX;
00710 t->heat = -FLT_MAX;
00711
00712
00713 t->col_str = -FLT_MAX;
00714
00715 for( long i=0; i<ipNCOLLIDER; i++ )
00716 t->col_stri[i] = -FLT_MAX;
00717
00718 return;
00719 }
00720
00721
00722 void StateJunk( quantumState * t )
00723 {
00724
00725 DEBUG_ENTRY( "StateJunk()" );
00726
00727 t->chLabel[0] = '\0';
00728
00730 t->g = -FLT_MAX;
00731
00733 t->Pop = -FLT_MAX;
00734
00736 t->IonStg = -10000;
00737
00739 t->nelem = -10000;
00740
00741
00742 t->sp = NULL;
00743
00744 return;
00745 }
00746
00747
00748 void EmLineZero( emission * t )
00749 {
00750
00751 DEBUG_ENTRY( "EmLineZero()" );
00752
00753
00754
00755 t->TauCon = opac.taumin;
00756
00757
00758
00759 t->TauIn = opac.taumin;
00760
00761
00762 t->TauTot = 1e20f;
00763
00764 t->TauTrack.clear();
00765
00766
00767
00768 t->FracInwd = 1.;
00769
00770
00771 t->pump = 0.;
00772
00773
00774 t->xIntensity = 0.;
00775
00776
00777 t->phots = 0.;
00778
00779
00780
00781 t->Pesc = 1.;
00782 t->Pdest = 0.;
00783 t->Pelec_esc = 0.;
00784
00785
00786 t->ColOvTot = 1.;
00787
00788
00789 t->PopOpc = 0.;
00790
00791
00792 t->ots = 0.;
00793 return;
00794 }
00795
00796
00797 void CollisionZero( collision * t )
00798 {
00799
00800 DEBUG_ENTRY( "CollisionZero()" );
00801
00802 t->ColUL = 0.;
00803
00804 t->cool = 0.;
00805 t->heat = 0.;
00806 return;
00807 }
00808
00809
00810 void StateZero( quantumState * t )
00811 {
00812
00813 DEBUG_ENTRY( "StateZero()" );
00814
00816 t->Pop = 0.;
00817 return;
00818 }
00819
00820
00821 double ConvRate2CS( realnum gHi , realnum rate )
00822 {
00823
00824 double cs;
00825
00826 DEBUG_ENTRY( "ConvRate2CS()" );
00827
00828
00829
00830
00831 cs = rate * gHi / dense.cdsqte;
00832
00833
00834
00835 ASSERT( cs >= 0. );
00836 return cs;
00837 }
00838
00839
00840 double ConvCrossSect2CollStr( double CrsSectCM2, double gLo, double E_ProjectileRyd, double reduced_mass_grams )
00841 {
00842 double CollisionStrength;
00843
00844 DEBUG_ENTRY( "ConvCrossSect2CollStr()" );
00845
00846 ASSERT( CrsSectCM2 >= 0. );
00847 ASSERT( gLo >= 0. );
00848 ASSERT( E_ProjectileRyd >= 0. );
00849 ASSERT( reduced_mass_grams >= 0. );
00850
00851 CollisionStrength = CrsSectCM2 * gLo * E_ProjectileRyd / (PI*BOHR_RADIUS_CM*BOHR_RADIUS_CM);
00852
00853
00854 #if 0
00855 CollisionStrength *= reduced_mass_grams / ELECTRON_MASS;
00856 #endif
00857
00858 ASSERT( CollisionStrength >= 0. );
00859 return CollisionStrength;
00860 }
00861
00862
00863 double totlin(
00864
00865
00866
00867
00868 int chInfo)
00869 {
00870 long int i;
00871 double totlin_v;
00872
00873 DEBUG_ENTRY( "totlin()" );
00874
00875
00876
00877
00878
00879
00880 if( (chInfo != 'i' && chInfo != 'r') && chInfo != 'c' )
00881 {
00882 fprintf( ioQQQ, " TOTLIN does not understand chInfo=%c\n",
00883 chInfo );
00884 cdEXIT(EXIT_FAILURE);
00885 }
00886
00887
00888
00889 totlin_v = 0.;
00890 for( i=0; i < LineSave.nsum; i++ )
00891 {
00892 if( LineSv[i].chSumTyp == chInfo )
00893 {
00894 totlin_v += LineSv[i].SumLine[0];
00895 }
00896 }
00897 return( totlin_v );
00898 }
00899
00900
00901
00902 void FndLineHt(long int *level,
00903
00904 long int *ipStrong,
00905 double *Strong)
00906 {
00907 long int i;
00908
00909 DEBUG_ENTRY( "FndLineHt()" );
00910
00911 *Strong = 0.;
00912 *level = 0;
00913
00914
00915 for( i=1; i <= nLevel1; i++ )
00916 {
00917
00918 if( TauLines[i].Coll.heat > *Strong )
00919 {
00920 *ipStrong = i;
00921 *level = 1;
00922 *Strong = TauLines[i].Coll.heat;
00923 }
00924 }
00925
00926
00927 for( i=0; i < nWindLine; i++ )
00928 {
00929 if( TauLine2[i].Hi->IonStg < TauLine2[i].Hi->nelem+1-NISO )
00930 {
00931
00932 if( TauLine2[i].Coll.heat > *Strong )
00933 {
00934 *ipStrong = i;
00935 *level = 2;
00936 *Strong = TauLine2[i].Coll.heat;
00937 }
00938 }
00939 }
00940
00941
00942 for( i=0; i < nHFLines; i++ )
00943 {
00944
00945 if( HFLines[i].Coll.heat > *Strong )
00946 {
00947 *ipStrong = i;
00948 *level = 3;
00949 *Strong = HFLines[i].Coll.heat;
00950 }
00951 }
00952
00953
00954 for( i=0; i <linesAdded2; i++)
00955 {
00956
00957 if( dBaseLines[i].tran->Coll.heat > *Strong )
00958 {
00959 *ipStrong = i;
00960 *level = 4;
00961 *Strong = dBaseLines[i].tran->Coll.heat;
00962 }
00963 }
00964
00965 fixit();
00966
00967
00968 return;
00969 }
00970
00971 quantumState *AddState2Stack( void )
00972 {
00973 DEBUG_ENTRY( "AddState2Stack()" );
00974
00975 ASSERT( !lgStatesAdded );
00976
00977 currentState = new quantumState;
00978
00979 StateJunk( currentState );
00980
00981 if( statesAdded == 0 )
00982 {
00983 GenericStates = currentState;
00984 GenericStates->next = NULL;
00985 lastState = GenericStates;
00986 }
00987 else
00988 {
00989 StateZero( currentState );
00990 lastState->next = currentState;
00991 lastState = lastState->next;
00992 }
00993
00994 statesAdded++;
00995
00996 return currentState;
00997 }
00998
00999 emission *AddLine2Stack( bool lgRadiativeTrans )
01000 {
01001 DEBUG_ENTRY( "AddLine2Stack()" );
01002
01003 if( !lgRadiativeTrans )
01004 {
01005 return &DummyEmis;
01006 }
01007 else
01008 {
01009 ASSERT( lgLinesAdded == false );
01010
01011 currentLine = new emission;
01012
01013 EmLineJunk( currentLine );
01014
01015 if( linesAdded == 0 )
01016 {
01017 GenericLines = currentLine;
01018 GenericLines->next = NULL;
01019 lastLine = GenericLines;
01020 }
01021 else
01022 {
01023
01024
01025
01026
01027 EmLineZero( currentLine );
01028
01029 lastLine->next = currentLine;
01030 lastLine = lastLine->next;
01031 }
01032
01033 linesAdded++;
01034 return currentLine;
01035 }
01036 }