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 "h2.h"
00014 #include "hmi.h"
00015 #include "opacity.h"
00016 #include "radius.h"
00017 #include "atomfeii.h"
00018 #include "rt.h"
00019
00020
00021 void RT_tau_inc(void)
00022 {
00023
00024 long int i,
00025 ipHi,
00026 nelem,
00027 ipLo,
00028 ipISO;
00029
00030 double factor;
00031
00032 DEBUG_ENTRY( "RT_tau_inc()" );
00033
00034 if( trace.lgTrace )
00035 {
00036 fprintf( ioQQQ, " RT_tau_inc called.\n" );
00037 }
00038
00039
00040 RT_line_all( false , true);
00041
00042 opac.telec += (float)(radius.drad_x_fillfac*dense.eden*6.65e-25);
00043 opac.thmin += (float)(radius.drad_x_fillfac*hmi.Hmolec[ipMHm]*3.9e-17*
00044 (1. - rfield.ContBoltz[hmi.iphmin-1]/ hmi.hmidep));
00045
00046
00047 rt.dTauMase = 0;
00048 rt.mas_species = 0;
00049 rt.mas_ion = 0;
00050 rt.mas_hi = 0;
00051 rt.mas_lo = 0;
00052
00053
00054
00055 for( ipISO=ipH_LIKE; ipISO<NISO; ++ipISO )
00056 {
00057 for( nelem=ipISO; nelem < LIMELM; nelem++ )
00058 {
00059
00060
00061 int ion = nelem+1-ipISO;
00062
00063
00064
00065
00066
00067
00068 if( ion <=dense.IonHigh[nelem] && dense.xIonDense[nelem][ion] > dense.density_low_limit )
00069 {
00070
00071
00072
00073
00074
00075 factor = dense.xIonDense[nelem][ion];
00076
00077
00078
00079
00080
00081
00082
00083
00084 for( ipHi=1; ipHi < iso.numLevels_local[ipISO][nelem]; ipHi++ )
00085 {
00086 for( ipLo=0; ipLo < ipHi; ipLo++ )
00087 {
00088
00089
00090 double save = EmisLines[ipISO][nelem][ipHi][ipLo].PopOpc;
00091 EmisLines[ipISO][nelem][ipHi][ipLo].PopOpc *= factor;
00092
00093 RT_line_one_tauinc(&EmisLines[ipISO][nelem][ipHi][ipLo] ,
00094 ipISO , nelem , ipHi , ipLo );
00095
00096 EmisLines[ipISO][nelem][ipHi][ipLo].PopOpc = save;
00097 }
00098 }
00099 ipLo = 0;
00100
00101 for( ipHi=iso.quant_desig[ipISO][nelem][iso.numLevels_max[ipISO][nelem]-1].n; ipHi < iso.nLyman[ipISO]; ipHi++ )
00102
00103 {
00104
00105 iso.ExtraLymanLines[ipISO][nelem][ipHi].PopOpc =
00106
00107 EmisLines[ipISO][nelem][3][ipLo].PopOpc *factor;
00108
00109
00110 RT_line_one_tauinc(&iso.ExtraLymanLines[ipISO][nelem][ipHi] ,
00111 -1 ,ipISO , nelem , ipHi );
00112
00113
00114 iso.ExtraLymanLines[ipISO][nelem][ipHi].PopOpc =
00115 EmisLines[ipISO][nelem][3][ipLo].PopOpc;
00116 }
00117 }
00118 }
00119 }
00120
00121
00122
00123
00124
00125 for( i=1; i <= nLevel1; i++ )
00126 {
00127 RT_line_one_tauinc(&TauLines[i],
00128 -2 , -2 , -2 , i );
00129 }
00130
00131
00132 for( i=0; i < nWindLine; i++ )
00133 {
00134
00135
00136 if( TauLine2[i].IonStg < TauLine2[i].nelem+1-NISO )
00137 {
00138 RT_line_one_tauinc(&TauLine2[i] ,
00139 -3 , -3 , -3 , i );
00140 }
00141 }
00142
00143
00144 for( i=0; i < nUTA; i++ )
00145 {
00146 if( UTALines[i].Aul > 0. )
00147 {
00148
00149 UTALines[i].PopOpc = dense.xIonDense[UTALines[i].nelem-1][UTALines[i].IonStg-1];
00150 UTALines[i].PopLo = dense.xIonDense[UTALines[i].nelem-1][UTALines[i].IonStg-1];
00151 UTALines[i].PopHi = 0.;
00152 RT_line_one_tauinc(&UTALines[i], -4 , -4 , -4 , i );
00153 }
00154 }
00155
00156
00157 for( i=0; i < nHFLines; i++ )
00158 {
00159
00160 float save = dense.xIonDense[HFLines[i].nelem-1][HFLines[i].IonStg-1];
00161
00162
00163 if( save<=0. ) continue;
00164
00165
00166 dense.xIonDense[HFLines[i].nelem-1][HFLines[i].IonStg-1] *= hyperfine.HFLabundance[i];
00167
00168 RT_line_one_tauinc(&HFLines[i] , -5 , -5 , -5 , i );
00169
00170
00171 dense.xIonDense[HFLines[i].nelem-1][HFLines[i].IonStg-1] = save;
00172 }
00173 # if 0
00174 {
00175
00176
00177 enum {DEBUG_LOC=false};
00178
00179 if( DEBUG_LOC )
00180 {
00181 fprintf(ioQQQ,"isotope pops\t%.3e\t%.3e\t%.3e\t%.3e\n",
00182
00183
00184 HFLines[0].TauIn , HFLines[0].TauTot,
00185 TauLines[ipH21cm].TauIn , TauLines[ipH21cm].TauTot);
00186 }
00187 }
00188 # endif
00189
00190
00191 for( i=0; i < nCORotate; i++ )
00192 {
00193 RT_line_one_tauinc(&C12O16Rotate[i],
00194 -6 , -6 , -6 , i );
00195 RT_line_one_tauinc(&C13O16Rotate[i],
00196 -7 , -7 , -7 , i );
00197 }
00198
00199
00200 FeII_RT_TauInc();
00201
00202
00203 H2_RT_tau_inc();
00204
00205
00206 if( wind.windv == 0. )
00207 {
00208
00209
00210 for( i=0; i < NFEII; i++ )
00211 {
00212
00213 fe2ovr_la.Fe2TauLte[i] += fe2ovr_la.feopc[i]*(float)radius.drad_x_fillfac;
00214 }
00215 }
00216
00217 if( trace.lgTrace && trace.lgOptcBug )
00218 {
00219 fprintf( ioQQQ, " RT_tau_inc updated optical depths:\n" );
00220 prtmet();
00221 }
00222
00223 if( trace.lgTrace )
00224 {
00225 fprintf( ioQQQ, " RT_tau_inc returns.\n" );
00226 }
00227
00228 DEBUG_EXIT( "RT_tau_inc()" );
00229 return;
00230 }
00231