/*pldata punches selected line data for all lines transferred in code */ #include #include "cddefines.h" #include "taulines.h" #include "hmi.h" #include "phycon.h" #include "ionfracs.h" #include "elementnames.h" #include "h.h" #include "tseton.h" #include "forcte.h" #include "elecden.h" #include "reccno.h" #include "chlinelbl.h" #include "printe82.h" #include "linefit.h" #include "iwavlen.h" #include "coolr.h" #include "pldata.h" void pldata(FILE * io) { char chLbl[11] , /* these are used for parts of the line label */ chLbl2[11]; long int i, iWL, j, limit; const long nskip=2; /* number of emission lines per line of output */ double tot; double CritDen; # ifdef DEBUG_FUN fputs( "<+>pldata()\n", debug_fp ); # endif /* routine punches out (on unit io) line data * for all recombination lines, and all transitions that are transferred */ /* say what is happening so we know why we stopped */ fprintf( ioQQQ, " punching line data, then stopping\n" ); /* evaluate rec coef at constant temperature if this is set, else * use 10,000K */ if( tseton.lgTSetOn ) { phycon.te = forcte.ForceTemp; } else { phycon.te = 1e4; } /* this is set of Dima's recombination lines */ LineFit(phycon.te,RecCNO.RecCoefCNO); fprintf( io, "\n Recombination lines of C, N, O\n" ); fprintf( io, " Ion WL(A) Coef Ion WL(A) Coef\n" ); for( i=0; i<471; i+=nskip) { /* nskip is set to 2 above */ limit = MIN2(471,i+nskip); fprintf( io, " " ); for( j=i; j < limit; j++ ) { fprintf( io, "%2.2s%2ld%6ld%8.3f ", elementnames.chElementSym[(long)(RecCNO.RecCoefCNO[0][j])-1], (long)(RecCNO.RecCoefCNO[0][j]-RecCNO.RecCoefCNO[1][j]+1.01), (long)(RecCNO.RecCoefCNO[2][j]+0.5), log10(MAX2(SMALLFLOAT,RecCNO.RecCoefCNO[3][j]) ) ); } fprintf( io, " \n" ); } fprintf( io, "\n\n" ); phycon.eden = 1.; phycon.hden = 1.; ElecDen.EdenHCorr = 1.; /* want very small neutral fractions so get mostly e- cs */ h.hi = 1.e-5f; hmi.htwo = 0.; h.hii = 1.; for( i=1; i <= nLevel1; i++ ) { TauLines[i].PopLo = 1.; } for( i=0; i < nWindLine; i++ ) { TauLine2[i].PopLo = 1.; } for( i=0; i < LIMELM; i++ ) { for( j=0; j < LIMELM+1; j++ ) { xIonFracs[j][i] = 1.; } } /* evaluate cooling, this forces evaluation of collision strengths */ coolr(&tot); fprintf( io, " Level 1 transferred lines\n" ); fprintf( io, "Ion label WL gl gu gf A CS n(crt)\n" ); for( i=1; i <= nLevel1; i++ ) { /* chLineLbl generates a 1 char string from the line transfer array info - * "Ne 2 128" the string is null terminated, * in following printout the first 4 char is used first, followed by * an integer, followed by the rest of the array*/ strcpy( chLbl, chLineLbl(&TauLines[i]) ); /* this is the second piece of the line label, pick up string after start */ strcpy( chLbl2, chLbl+4 ); /* following logic must be kept parallel with that in putline */ iWL = iWavLen(&TauLines[i]); fprintf( io, "%4.4s%5ld%6.6s%3ld%3ld", /* ion label, like C 1 */ chLbl , /* integer label wavelength from printout */ iWL, /* label saying something like 379m, 3456A */ chLbl2, /* lower and upper stat weights */ (long)(TauLines[i].gLo), (long)(TauLines[i].gHi) ); /* oscillator strength */ fprintf( io,PrintEfmt("%9.2e", TauLines[i].gf )); /* einstein A for transition */ fprintf( io,PrintEfmt("%9.2e", TauLines[i].Aul) ); /* check because some CS are huge, others small, but avoid exp notation */ if( TauLines[i].cs > 1. ) { /* collision strength in f foramt */ fprintf( io, "%7.3f",TauLines[i].cs ); } else { fprintf( io, "%7.4f",TauLines[i].cs ); } /* now print critical density, * the 100 is the sqrt of Te */ if( TauLines[i].cs> 0. ) { CritDen = TauLines[i].Aul * TauLines[i].gHi*100. / (TauLines[i].cs*8.629e-6); CritDen = log10(CritDen); } else { CritDen = 0.; } fprintf( io, "%7.3f\n",CritDen ); } fprintf( io, "\n\n\n" ); fprintf( io, " end level 1, start level 2\n" ); fprintf( io, "Ion label WL gl gu gf A CS n(crt)\n" ); for( i=0; i < nWindLine; i++ ) { iWL = iWavLen(&TauLine2[i]); /* iWavLen = WindLine( ipLnWlAng ,i) * following logic must be kept parallel with that in putline * if(iWavLen.gt.1 000 000 ) then * iWavLen = iWavLen / 10000 * else if(iWavLen.gt.10 000 ) then * iWavLen = iWavLen / 1000 * else if( iWavLen.lt.10 ) then * iWavLen = WindLine( ipLnWlAng ,i) * 10. * endif */ strcpy( chLbl, chLineLbl(&TauLine2[i]) ); /* this is the second piece of the line label, pick up string after start */ strcpy( chLbl2, chLbl+4 ); fprintf( io, "%4.4s%5ld%6.6s%3ld%3ld", /* ion label, like C 1 */ chLbl , /* integer label wavelength from printout */ iWL, /* label saying something like 379m, 3456A */ chLbl2, /* lower and upper stat weights */ (long)(TauLine2[i].gLo), (long)(TauLine2[i].gHi) ); /* oscillator strength */ fprintf( io,PrintEfmt("%9.2e", TauLine2[i].gf)); /* einstein A for transition */ fprintf( io,PrintEfmt("%9.2e", TauLine2[i].Aul)); /* use different formats depending on size of collision strength */ if( TauLine2[i].cs > 1. ) { fprintf( io, "%7.3f", TauLine2[i].cs ); } else { fprintf( io, "%7.4f", TauLine2[i].cs ); } /* now print critical density, * the 100 is the sqrt of Te */ if( TauLine2[i].cs> 0. ) { CritDen = TauLine2[i].Aul * TauLine2[i].gHi*100. / (TauLine2[i].cs*8.629e-6); CritDen = log10(CritDen); } else { CritDen = 0.; } fprintf( io, "%7.3f\n",CritDen ); } /* stop when done, we have done some serious damage to the code */ puts( "[Normal top in pldata]" ); fprintf( ioQQQ, " Cloudy ends:\n" ); exit(1); }