00001
00002
00003
00004
00005
00006 #include "cddefines.h"
00007 #include "taulines.h"
00008 #include "atomfeii.h"
00009 #include "dense.h"
00010 #include "conv.h"
00011 #include "atoms.h"
00012 #include "rfield.h"
00013 #include "iso.h"
00014 #include "h2.h"
00015 #include "opacity.h"
00016 #include "trace.h"
00017 #include "lines_service.h"
00018 #include "atmdat.h"
00019 #include "hydrogenic.h"
00020 #include "rt.h"
00021
00022 void RT_line_all(
00023
00024
00025 bool lgDoEsc ,
00026
00027 bool lgUpdateFineOpac )
00028 {
00029 long int i,
00030 ion,
00031 ipISO,
00032 nelem;
00033 long ipHi , ipLo;
00034 double factor,
00035 coloi;
00036 bool lgTOnSave;
00037 double SaveLyaPesc[NISO][LIMELM] ,
00038 SaveLyaPdest[NISO][LIMELM];
00039
00040 bool lgRT_update = lgUpdateFineOpac || lgDoEsc || conv.lgIonStageTrimed;
00041
00042 DEBUG_ENTRY( "RT_line_all()" );
00043
00044 if( trace.lgTrace )
00045 fprintf( ioQQQ, " RT_line_all called\n" );
00046
00047
00048 if( !rfield.lgDoLineTrans )
00049 {
00050 if( conv.nPres2Ioniz )
00051 {
00052 DEBUG_EXIT( "RT_line_all()" );
00053 return;
00054 }
00055 }
00056
00057
00058
00059 rfield.lgFine_opac_update = lgUpdateFineOpac;
00060 if( rfield.lgFine_opac_update )
00061 {
00062
00063 memset(rfield.fine_opac_zone , 0 , (unsigned long)rfield.nfine*sizeof(float) );
00064 }
00065
00066
00067
00068 lgTOnSave = opac.lgTauOutOn;
00069 opac.lgTauOutOn = true;
00070
00071 # if 0
00072
00073
00074
00075 {
00076 #include "doppvel.h"
00077 double dTau , aa;
00078 ipISO = 0; nelem = 0;ipLo = 0;
00079 ipHi = 23;
00080 dTau = EmisLines[ipISO][nelem][ipHi][ipLo].PopOpc *
00081 EmisLines[ipISO][nelem][ipHi][ipLo].opacity /
00082 DoppVel.doppler[nelem] + opac.opacity_abs[EmisLines[ipISO][nelem][ipHi][ipLo].ipCont-1];
00083 dTau *= radius.drad_x_fillfac_mean;
00084 aa = log(1. + dTau ) / SDIV(dTau);
00085 fprintf(ioQQQ,"DEBUG dTau\t%.2f\t%.5e\t%.5e\t%.5e\n",fnzone,dTau,
00086 radius.drad_x_fillfac_mean,
00087 aa);
00088 }
00089 # endif
00090
00091
00092
00093 if( lgDoEsc )
00094 {
00095 for( nelem=0; nelem < LIMELM; nelem++ )
00096 {
00097
00098 if( dense.IonHigh[nelem] == (nelem + 1) )
00099 {
00100
00101 HLineTransOpacSet( nelem );
00102 }
00103 }
00104 }
00105
00106
00107 RT_stark();
00108
00109
00110
00111 for( ipISO=ipH_LIKE; ipISO<NISO; ++ipISO )
00112 {
00113
00114 for( nelem=ipISO; nelem<LIMELM; ++nelem )
00115 {
00116
00117
00118 ion = nelem+1-ipISO;
00119
00120 if( ion <=dense.IonHigh[nelem] )
00121 {
00122
00123
00124
00125
00126
00127 if( ipISO<=nelem )
00128 {
00129 SaveLyaPesc[ipISO][nelem] = EmisLines[ipISO][nelem][iso.nLyaLevel[ipISO]][0].Pesc;
00130 SaveLyaPdest[ipISO][nelem] = EmisLines[ipISO][nelem][iso.nLyaLevel[ipISO]][0].Pdest;
00131 }
00132
00133
00134 if( dense.xIonDense[nelem][ion] > 1e-30 )
00135 {
00136 factor = dense.xIonDense[nelem][ion];
00137 }
00138 else
00139 {
00140
00141
00142 factor = 1.;
00143 }
00144
00145
00146 if( ipISO==ipH_LIKE )
00147 opac.lgTauOutOn = true;
00148
00149
00150 for( ipLo=0; ipLo < (iso.numLevels_local[ipISO][nelem]-1); ++ipLo )
00151 {
00152
00153 for( ipHi=ipLo+1; ipHi < iso.numLevels_local[ipISO][nelem]; ++ipHi )
00154 {
00155
00156
00157 if( EmisLines[ipISO][nelem][ipHi][ipLo].ipCont < 0 )
00158 continue;
00159
00160
00161 EmisLines[ipISO][nelem][ipHi][ipLo].PopOpc *= factor;
00162
00163
00164
00165
00166
00167
00168
00169
00170
00171 RT_line_one(&EmisLines[ipISO][nelem][ipHi][ipLo] , lgDoEsc , lgUpdateFineOpac,true);
00172 {
00173
00174
00175 enum {DEBUG_LOC=false};
00176
00177 if( DEBUG_LOC && nelem==1&& ipLo==0 )
00178 {
00179 fprintf(ioQQQ,"DEBUG pdest\t%3li\t%.2f\t%.3e\n",
00180 ipHi ,
00181 fnzone,
00182 EmisLines[ipISO][nelem][ipHi][ipLo].Pdest);
00183 }
00184 }
00185
00186
00187 EmisLines[ipISO][nelem][ipHi][ipLo].PopOpc /= factor;
00188 }
00189
00190
00191 opac.lgTauOutOn = lgTOnSave;
00192 }
00193
00194
00195
00196 if( ipISO > 0 && lgDoEsc )
00197 {
00198
00199 EmisLines[ipISO][nelem][iso.nLyaLevel[ipISO]][0].Pesc *=
00200 opac.ExpmTau[EmisLines[ipISO][nelem][iso.nLyaLevel[ipISO]][0].ipCont-1];
00201 }
00202
00203
00204 atmdat_2phot_rate(ipISO , nelem);
00205
00206
00207
00208 if( opac.lgCaseB_no_pdest )
00209 {
00210 ipLo = 0;
00211
00212 for( ipHi=ipLo+1; ipHi < iso.numLevels_max[ipISO][nelem]; ++ipHi )
00213 {
00214
00215 EmisLines[ipISO][nelem][ipHi][ipLo].Pdest = SMALLFLOAT;
00216 }
00217 }
00218
00219 ipLo = 0;
00220
00221
00222
00223 if( dense.xIonDense[nelem][ion] > 1e-30 && lgUpdateFineOpac )
00224 {
00225
00226 factor = EmisLines[ipISO][nelem][3][ipLo].PopOpc *dense.xIonDense[nelem][ion];
00227 for( ipHi=iso.quant_desig[ipISO][nelem][iso.numLevels_max[ipISO][nelem]-1].n; ipHi < iso.nLyman[ipISO]; ipHi++ )
00228 {
00229
00230 iso.ExtraLymanLines[ipISO][nelem][ipHi].PopOpc = factor;
00231 iso.ExtraLymanLines[ipISO][nelem][ipHi].PopLo =
00232 EmisLines[ipISO][nelem][3][ipLo].PopLo;
00233
00234
00235 RT_line_one(&iso.ExtraLymanLines[ipISO][nelem][ipHi] , lgDoEsc , lgUpdateFineOpac,true);
00236
00237
00238 iso.ExtraLymanLines[ipISO][nelem][ipHi].PopOpc =
00239 EmisLines[ipISO][nelem][3][ipLo].PopOpc;
00240 }
00241 }
00242
00243 }
00244 }
00245 }
00246
00247
00248
00249
00250 if( lgTauGood( &EmisLines[ipH_LIKE][ipHYDROGEN][iso.nLyaLevel[ipH_LIKE]][0] ) )
00251 {
00252 if( lgDoEsc )
00253 {
00254
00255
00256
00257 for( ipLo=ipH1s; ipLo < (iso.numLevels_max[ipH_LIKE][ipHYDROGEN] - 1); ipLo++ )
00258 {
00259
00260 for( ipHi=MAX2((long)3,ipLo+1); ipHi < iso.numLevels_max[ipH_LIKE][ipHYDROGEN]; ipHi++ )
00261 {
00262
00263
00264 if( EmisLines[ipH_LIKE][ipHYDROGEN][ipHi][ipLo].Pesc<1. )
00265 {
00266 EmisLines[ipH_LIKE][ipHYDROGEN][ipHi][ipLo].Pesc = MIN2(1.f,
00267 EmisLines[ipH_LIKE][ipHYDROGEN][ipHi][ipLo].Pesc+
00268 (float)hydro.pestrk[ipLo][ipHi]);
00269 }
00270 }
00271 }
00272 }
00273
00274
00275 if( EmisLines[ipH_LIKE][ipHYDROGEN][ipH2p][ipH1s].Pesc < 1. )
00276 {
00277 EmisLines[ipH_LIKE][ipHYDROGEN][ipH2p][ipH1s].Pesc = MIN2(1.f,
00278 EmisLines[ipH_LIKE][ipHYDROGEN][ipH2p][ipH1s].Pesc+
00279 (float)hydro.pestrk[ipH1s][ipH2p]);
00280 }
00281
00282
00283 atom_oi_calc(&coloi);
00284
00285 EmisLines[ipH_LIKE][ipHYDROGEN][3][ipH1s].Pesc = (float)(atoms.pmph31/
00286 EmisLines[ipH_LIKE][ipHYDROGEN][3][ipH1s].Aul);
00287
00288 if( trace.lgTrace && trace.lgIsoTraceFull[ipH_LIKE] )
00289 {
00290 fprintf( ioQQQ, " RT_line_all calls P8446 who found pmph31=%10.2e\n",
00291 atoms.pmph31 );
00292 }
00293
00294
00295
00296
00297
00298
00299
00300
00301 atoms_fe2ovr();
00302
00303
00304
00305
00306
00307
00308 EmisLines[ipH_LIKE][ipHYDROGEN][ipH2p][ipH1s].Pdest += hydro.dstfe2lya;
00309 }
00310
00311
00312 if( nzone > 1 )
00313 {
00314 for( ipISO=ipH_LIKE; ipISO<NISO; ++ipISO )
00315 {
00316
00317 for( nelem=ipISO; nelem<LIMELM; ++nelem )
00318 {
00319 ion = nelem+1-ipISO;
00320 if( ion <=dense.IonHigh[nelem] && ipISO<=nelem )
00321 {
00322
00323
00324 if( EmisLines[ipISO][nelem][iso.nLyaLevel[ipISO]][0].iRedisFun==ipLY_A &&
00325 lgTauGood( &EmisLines[ipISO][nelem][iso.nLyaLevel[ipISO]][0] ) )
00326 {
00327
00328
00329
00330
00331 EmisLines[ipISO][nelem][iso.nLyaLevel[ipISO]][0].Pdest =
00332 ((float)SaveLyaPdest[ipISO][nelem] +
00333 EmisLines[ipISO][nelem][iso.nLyaLevel[ipISO]][0].Pdest) / 2.f;
00334 EmisLines[ipISO][nelem][iso.nLyaLevel[ipISO]][0].Pesc =
00335 ((float)SaveLyaPesc[ipISO][nelem] +
00336 EmisLines[ipISO][nelem][iso.nLyaLevel[ipISO]][0].Pesc) / 2.f;
00337
00338 }
00339 }
00340 }
00341 }
00342 }
00343
00344
00345
00346
00347
00348
00349
00350 if( !hydro.lgLymanPumping )
00351 {
00352 ipISO = ipH_LIKE;
00353 nelem = ipHYDROGEN;
00354 ipLo = 0;
00355 for( ipHi=ipLo+1; ipHi < iso.numLevels_max[ipISO][nelem]; ++ipHi )
00356 {
00357
00358
00359 EmisLines[ipISO][nelem][ipHi][ipLo].pump = 0.;
00360 }
00361 }
00362
00363 {
00364
00365
00366 enum {DEBUG_LOC=false};
00367
00368 if( DEBUG_LOC && nzone>433 )
00369 {
00370
00371 DumpLine(&EmisLines[ipH_LIKE][ipHYDROGEN][2][0] );
00372 # if 0
00373 fprintf(ioQQQ,"debugggg\t%.3f\t%i\t%.2e\t%.2e\t%.2e\t%.2e\t%.2e\t%.2e\t%.2e\t%.2e\t%.2e\n",
00374 fnzone,
00375 lgDoEsc ,
00376 EmisLines[ipH_LIKE][ipHYDROGEN][2][0].Pesc ,
00377 EmisLines[ipH_LIKE][ipHYDROGEN][2][0].Pdest,
00378 EmisLines[ipH_LIKE][ipHYDROGEN][2][0].pump ,
00379 EmisLines[ipH_LIKE][ipHYDROGEN][2][0].TauCon ,
00380 EmisLines[ipH_LIKE][ipHYDROGEN][2][0].TauIn ,
00381 EmisLines[ipH_LIKE][ipHYDROGEN][2][0].TauTot,
00382 EmisLines[ipH_LIKE][ipHYDROGEN][3][2].Pesc ,
00383 EmisLines[ipH_LIKE][ipHYDROGEN][3][2].Pdest ,
00384 EmisLines[ipH_LIKE][ipHYDROGEN][3][2].pump );
00385 # endif
00386 }
00387 }
00388
00389
00390 for( i=1; i <= nLevel1; i++ )
00391 {
00392 RT_line_one(&TauLines[i] , lgDoEsc , lgUpdateFineOpac,true);
00393 }
00394
00395 for( i=0; i < nCORotate; i++ )
00396 {
00397 RT_line_one(&C12O16Rotate[i] , lgDoEsc , lgUpdateFineOpac,true);
00398 RT_line_one(&C13O16Rotate[i] , lgDoEsc , lgUpdateFineOpac,true);
00399 }
00400 for( i=0; i < nHFLines; i++ )
00401 {
00402 RT_line_one(&HFLines[i] , lgDoEsc , lgUpdateFineOpac,true);
00403 }
00404
00405
00406 if( lgRT_update )
00407 {
00408 for( i=0; i < nUTA; i++ )
00409 {
00410
00411
00412 if( UTALines[i].Aul > 0. )
00413 {
00414
00415 UTALines[i].PopOpc = dense.xIonDense[UTALines[i].nelem-1][UTALines[i].IonStg-1];
00416 UTALines[i].PopLo = dense.xIonDense[UTALines[i].nelem-1][UTALines[i].IonStg-1];
00417 UTALines[i].PopHi = 0.;
00418 RT_line_one(&UTALines[i] , lgDoEsc , lgUpdateFineOpac,true);
00419 }
00420 }
00421 }
00422
00423
00424
00425
00426
00427
00428
00429
00430
00431 for( i=0; i < nWindLine; i++ )
00432 {
00433
00434
00435 if( TauLine2[i].IonStg < TauLine2[i].nelem+1-NISO )
00436 {
00437 RT_line_one(&TauLine2[i] , lgDoEsc , lgUpdateFineOpac,true);
00438 }
00439 }
00440
00441
00442 H2_RTMake( lgDoEsc , lgUpdateFineOpac);
00443
00444
00445 FeII_RT_Make( lgDoEsc , lgUpdateFineOpac);
00446
00447 DEBUG_EXIT( "RT_line_all()" );
00448 return;
00449 }
00450