00001
00002
00003
00004 #include "cddefines.h"
00005 #include "taulines.h"
00006 #include "iso.h"
00007 #include "rfield.h"
00008 #include "trace.h"
00009 #include "dense.h"
00010 #include "hyperfine.h"
00011 #include "wind.h"
00012 #include "prt.h"
00013 #include "conv.h"
00014 #include "h2.h"
00015 #include "mole.h"
00016 #include "hmi.h"
00017 #include "opacity.h"
00018 #include "cooling.h"
00019 #include "thermal.h"
00020 #include "radius.h"
00021 #include "atomfeii.h"
00022 #include "rt.h"
00023 #include "doppvel.h"
00024 #include "mole.h"
00025
00026
00027 void RT_tau_inc(void)
00028 {
00029
00030 long int i,
00031 ipHi,
00032 ipLo;
00033
00034 DEBUG_ENTRY( "RT_tau_inc()" );
00035
00036 if( trace.lgTrace )
00037 {
00038 fprintf( ioQQQ, " RT_tau_inc called.\n" );
00039 }
00040
00041
00042 ASSERT( !conv.lgFirstSweepThisZone );
00043 conv.lgLastSweepThisZone = true;
00044 RT_line_all( );
00045
00046
00047 CoolEvaluate( &thermal.ctot );
00048
00049 if( nzone <=1 )
00050 {
00051 opac.telec = (realnum)(radius.drad_x_fillfac*dense.eden*6.65e-25);
00052 opac.thmin = (realnum)(radius.drad_x_fillfac*findspecieslocal("H-")->den*3.9e-17*
00053 (1. - rfield.ContBoltz[hmi.iphmin-1]/ hmi.hmidep));
00054 }
00055 else
00056 {
00057 opac.telec += (realnum)(radius.drad_x_fillfac*dense.eden*6.65e-25);
00058 opac.thmin += (realnum)(radius.drad_x_fillfac*findspecieslocal("H-")->den*3.9e-17*
00059 (1. - rfield.ContBoltz[hmi.iphmin-1]/ hmi.hmidep));
00060 }
00061
00062
00063 rt.dTauMase = 0;
00064 rt.mas_species = 0;
00065 rt.mas_ion = 0;
00066 rt.mas_hi = 0;
00067 rt.mas_lo = 0;
00068
00069
00070 for( long ipISO=ipH_LIKE; ipISO<NISO; ++ipISO )
00071 {
00072 for( long nelem=ipISO; nelem < LIMELM; nelem++ )
00073 {
00074
00075
00076 int ion = nelem+1-ipISO;
00077
00078 if( ion <=dense.IonHigh[nelem] && dense.xIonDense[nelem][ion] > dense.density_low_limit )
00079 {
00080
00081 for( ipLo=0; ipLo < iso_sp[ipISO][nelem].numLevels_local; ipLo++ )
00082 {
00083 if( iso_ctrl.lgDielRecom[ipISO] )
00084 RT_line_one_tauinc(SatelliteLines[ipISO][nelem][ipSatelliteLines[ipISO][nelem][ipLo]], ipISO, nelem, -1, ipLo,
00085 GetDopplerWidth(dense.AtomicWeight[nelem]) );
00086 }
00087
00088 for( ipHi=1; ipHi < iso_sp[ipISO][nelem].numLevels_local; ipHi++ )
00089 {
00090 for( ipLo=0; ipLo < ipHi; ipLo++ )
00091 {
00092 if( iso_sp[ipISO][nelem].trans(ipHi,ipLo).ipCont() <= 0 )
00093 continue;
00094
00095
00096 RT_line_one_tauinc(iso_sp[ipISO][nelem].trans(ipHi,ipLo), ipISO, nelem, ipHi, ipLo,
00097 GetDopplerWidth(dense.AtomicWeight[nelem]) );
00098 }
00099 }
00100 ipLo = 0;
00101
00102 for( ipHi=iso_sp[ipISO][nelem].st[iso_sp[ipISO][nelem].numLevels_local-1].n()+1; ipHi < iso_ctrl.nLyman[ipISO]; ipHi++ )
00103 {
00104 TransitionList::iterator tr = ExtraLymanLines[ipISO][nelem].begin()+ipExtraLymanLines[ipISO][nelem][ipHi];
00105 (*tr).Emis().PopOpc() = iso_sp[ipISO][nelem].st[0].Pop();
00106
00107
00108 RT_line_one_tauinc(*tr, -1 ,ipISO, nelem, ipHi,
00109 GetDopplerWidth(dense.AtomicWeight[nelem]) );
00110 }
00111 }
00112 }
00113 }
00114
00115
00116
00117
00118 for( i=1; i <= nLevel1; i++ )
00119 {
00120 RT_line_one_tauinc(TauLines[i], -2, -2, -2, i, GetDopplerWidth(dense.AtomicWeight[(*TauLines[i].Hi()).nelem()-1]) );
00121 }
00122
00123
00124 for( i=0; i < nWindLine; i++ )
00125 {
00126
00127
00128 if( (*TauLine2[i].Hi()).IonStg() < (*TauLine2[i].Hi()).nelem()+1-NISO )
00129 {
00130 RT_line_one_tauinc(TauLine2[i], -3, -3, -3, i, GetDopplerWidth(dense.AtomicWeight[(*TauLine2[i].Hi()).nelem()-1]) );
00131 }
00132 }
00133
00134
00135 for( i=0; i < nUTA; i++ )
00136 {
00137
00138 UTALines[i].Emis().PopOpc() = dense.xIonDense[(*UTALines[i].Hi()).nelem()-1][(*UTALines[i].Hi()).IonStg()-1];
00139 (*UTALines[i].Lo()).Pop() = dense.xIonDense[(*UTALines[i].Hi()).nelem()-1][(*UTALines[i].Hi()).IonStg()-1];
00140 (*UTALines[i].Hi()).Pop() = 0.;
00141 RT_line_one_tauinc(UTALines[i], -4 , -4 , -4 , i, GetDopplerWidth(dense.AtomicWeight[(*UTALines[i].Hi()).nelem()-1]) );
00142 }
00143
00144
00145 for( i=0; i < nHFLines; i++ )
00146 {
00147
00148 realnum save = dense.xIonDense[(*HFLines[i].Hi()).nelem()-1][(*HFLines[i].Hi()).IonStg()-1];
00149
00150
00151 if( save<=0. ) continue;
00152
00153
00154 dense.xIonDense[(*HFLines[i].Hi()).nelem()-1][(*HFLines[i].Hi()).IonStg()-1] *= hyperfine.HFLabundance[i];
00155
00156 RT_line_one_tauinc(HFLines[i] , -5 , -5 , -5 , i, GetDopplerWidth(dense.AtomicWeight[(*HFLines[i].Hi()).nelem()-1]) );
00157
00158
00159 dense.xIonDense[(*HFLines[i].Hi()).nelem()-1][(*HFLines[i].Hi()).IonStg()-1] = save;
00160 }
00161
00162
00163 FeII_RT_TauInc();
00164
00165
00166 for( diatom_iter diatom = diatoms.begin(); diatom != diatoms.end(); ++diatom )
00167 (*diatom)->H2_RT_tau_inc();
00168
00169
00170 for( long ipSpecies=0; ipSpecies<nSpecies; ipSpecies++ )
00171 {
00172 if( dBaseSpecies[ipSpecies].lgActive )
00173 {
00174 realnum DopplerWidth = GetDopplerWidth( dBaseSpecies[ipSpecies].fmolweight );
00175 for (TransitionList::iterator tr=dBaseTrans[ipSpecies].begin();
00176 tr != dBaseTrans[ipSpecies].end(); ++tr)
00177 {
00178 int ipHi = (*tr).ipHi();
00179 if (ipHi >= dBaseSpecies[ipSpecies].numLevels_local || (*tr).ipCont() <= 0)
00180 continue;
00181 int ipLo = (*tr).ipLo();
00182
00183 RT_line_one_tauinc( *tr, -10, ipSpecies, ipHi, ipLo, DopplerWidth );
00184 }
00185 }
00186 }
00187
00188
00189 if( wind.lgStatic() )
00190 {
00191
00192 t_fe2ovr_la::Inst().tau_inc();
00193 }
00194
00195 if( trace.lgTrace && trace.lgOptcBug )
00196 {
00197 fprintf( ioQQQ, " RT_tau_inc updated optical depths:\n" );
00198 prtmet();
00199 }
00200
00201 if( trace.lgTrace )
00202 fprintf( ioQQQ, " RT_tau_inc returns.\n" );
00203
00204 return;
00205 }