implicit none integer iopt,iopt1 integer istop real*8 wmin,wmax * read theoretical emissivities from individual files or binary file "OIIlines_AB" and pick * out the relevant lines istop=0 if(istop.eq.1) stop 1 write(6,10) 10 format(/1x,'Choose one of three options:' // : '1: Generate a line list' // : '2: Extract emission coefficients for a user supplied list of' : ' lines in user'/' specified conditions' // : '3: Exit' // : 1x,'Enter an integer:') read(5,*) iopt if(iopt.ne.1.and.iopt.ne.2.and.iopt.ne.3) go to 1 if(iopt.eq.3) stop call reademiss if(iopt.eq.1) then 2 write(6,20) 20 format(//1x,'Further options:'//'1: List all lines' : ' (8889 between 367.97 and 3.90713e+07 Angstrom)'// : '2: Specify a wavelength range for the' : ' list'//1x,'Enter an integer') read(5,*) iopt1 if(iopt1.ne.1.and.iopt1.ne.2) go to 2 wmin=367.97d0 wmax= 3.907135E+07 if(iopt1.eq.2) then 3 write(6,30) 30 format(/1x,'Enter lower and upper wavelength in Angstrom') read(5,*) wmin, wmax if(wmin.gt.wmax) then write(6,31) 31 format(1x,'Lower wavelength greater than upper') go to 3 endif endif endif if(iopt.eq.1) then call opt1(wmin,wmax) elseif(iopt.eq.2) then call opt2 endif go to 1 stop end ************************************************************************* subroutine opt1(wmin,wmax) implicit none * mxl = maximum number of line emission coeffs * mxt = maximum number of temperatures * mxd = maximum number of densities * mxout = maximum number of user requested emission coeffs integer mxl,mxt,mxd,mxout parameter (mxl=10000,mxt=25,mxd=16,mxout=30) integer iopt,istop character*1 eflag(mxl),eflagc(mxl),eflagcs,ul,Case character*84 txt1(mxl),txt1c(mxl),txt1cs character*32 txt2(mxl),txt2c(mxl),txt2cs character*19 woutfile,eoutfile real*8 em(mxt,mxd,mxl),wl(mxl), : eu(mxl),el(mxl),temlg(mxt), : denslg(mxd),wlc(mxl),emc(mxl), : emcs,wmin,wmax,ws,es,emmax,wv1,ev1 integer nt,nd,il,inl,jl,nl,nw common/alltheory/wl,em,eflag,eu,el,txt1,txt2,temlg,denslg, : nt,nd,nw,Case data ul/'_'/ data nl/8889/ woutfile='OIIdata'//ul//'list'//ul//'worder' eoutfile='OIIdata'//ul//'list'//ul//'eorder' open(unit=16,status='unknown',file=woutfile) open(unit=17,status='unknown',file=eoutfile) * Load line list for canonical temperature and density into local storage * and select those in specified range inl=0 emmax=0.d0 wv1=4649.13d0 do il=1,nl if(abs(wl(il)-wv1).lt.1.d-2) then ev1=em(21,11,il) endif if(wl(il).ge.wmin.and.wl(il).le.wmax) then inl=inl+1 wlc(inl)=wl(il) emc(inl)=em(21,11,il) txt1c(inl)=txt1(il) txt2c(inl)=txt2(il) eflagc(inl)=eflag(il) if(emc(inl).gt.emmax) emmax=emc(inl) c write(6,*) inl,il,wlc(inl) endif enddo write(16,1) inl,wmin,wmax write(17,1) inl,wmin,wmax write(6,1) inl,wmin,wmax 1 format(/1x,'Found ',i5,' lines in interval ',f12.2, : ' to ',f12.2,' Angstrom'/) if(inl.eq.0) return write(6,11) eoutfile,woutfile c write(16,11) eoutfile,woutfile 11 format(1x,'File ',a19, : ' contains the line list in order of decreasing emission' : ' coefficient'// : 1x,'File ',a19,' contains the full list preceded by a list' : ' of only the strongest lines, in order of increasing' : ' wavelength') * Sort into increasing wavelength order do il=1,inl-1 ws=wlc(il) do jl=il+1,inl if(wlc(jl).lt.ws) then wlc(il)=wlc(jl) wlc(jl)=ws ws=wlc(il) emcs=emc(il) txt1cs=txt1c(il) txt2cs=txt2c(il) emc(il)=emc(jl) emc(jl)=emcs txt1c(il)=txt1c(jl) txt1c(jl)=txt1cs txt2c(il)=txt2c(jl) txt2c(jl)=txt2cs eflagcs=eflagc(il) eflagc(il)=eflagc(jl) eflagc(jl)=eflagcs endif enddo enddo write(16,4) Case 4 format(' This file is in two sections:'/ : ' Section 1 - Lines with emission coefficient at' : ' least 0.1% of strongest line in this wavelength range,'/ : 12x,' in order of increasing wavelength'/ : ' Section 2 - All lines in this wavelength range in order' : ' of increasing wavelength'// : ' All emission coefficients at a temperature of 10000K and' : ' an electron number density of 10000/cm^3 in Case 'a1/) write(16,41) 41 format(1x,'I = Line index in full list for this' : ' wavelength range (see second section below)') write(16,42) 42 format( : 1x,'X/T = Wavelength from only experimental energies (X)' : ' or at least one theoretical energy (T)'/ : 1x,' T implies wavelength not' : ' spectroscopically useful'/ : 1x,'W = Wavelength[A],'/ : 1x,'EM = Emission coefficient[cgs],'/ : 1x,'E/EMAX = EM relative to the maximum value in this' : ' wavelength range,'/ : 1x,'E/EV1 = EM rel to V1 4649.13'/ : 1x,'IPI, IPF = parent level index'/ : 1x,'NLI, NLF = valence electron n, l of initial (I)' : ' and final (F) states'/ : 1x,'2JI, 2JF = 2*J (2*total angular momentum)'/ : 1x,'PI PF = parity (0=EVEN, 1=ODD)'/) write(16,22) 22 format(' Section 1'/) write(16,21) do il=1,inl if((emc(il)/emmax).gt.1.d-3) then if(wlc(il).le.1.e5) then write(16,5) il,wlc(il),eflagc(il),emc(il),emc(il)/emmax, : emc(il)/ev1,txt1c(il) 5 format(1x,i5,f12.2,2x,a1,1x,1pe12.3,0p2f12.5,a84) else write(16,51) il,wlc(il),eflagc(il),emc(il),emc(il)/emmax, : emc(il)/ev1,txt1c(il) 51 format(1x,i5,1pe12.5,2x,a1,1x,1pe12.3,0p2f12.5,a84) endif endif enddo write(16,2) 2 format(//' Section 2'/ : ' Full list in order of increasing wavelength'//) write(16,21) 21 format(5x,'I',11x,'W',14x,'EM',6x,'E/EMAX',7x,'E/EV1',3x, : 'Initial Term',12x,'Final Term',13x, : ' IPI NLI 2JI PI',3x,' IPF NLF 2JF PF') do il=1,inl if(wlc(il).le.1.e5) then if(emc(il)/ev1.gt.1.d-5) then write(16,3) il,wlc(il),eflagc(il),emc(il),emc(il)/emmax, : emc(il)/ev1,txt1c(il) 3 format(1x,i5,f12.2,2x,a1,1x,1p3e12.3,a84) else write(16,31) il,wlc(il),eflagc(il),emc(il),emc(il)/emmax, : emc(il)/ev1,txt1c(il) 31 format(1x,i5,f12.2,2x,a1,1x,1p3e12.3,a84) endif else if(emc(il)/ev1.gt.1.d-5) then write(16,32) il,wlc(il),eflagc(il),emc(il),emc(il)/emmax, : emc(il)/ev1,txt1c(il) 32 format(1x,i5,1pe12.5,2x,a1,1x,1p3e12.3,a84) else write(16,33) il,wlc(il),eflagc(il),emc(il),emc(il)/emmax, : emc(il)/ev1,txt1c(il) 33 format(1x,i5,1pe12.5,2x,a1,1x,1p3e12.3,a84) endif endif enddo * Sort into decreasing emission coefficient order do il=1,inl-1 es=emc(il) do jl=il+1,inl if(emc(jl).gt.es) then emc(il)=emc(jl) emc(jl)=es es=emc(il) ws=wlc(il) wlc(il)=wlc(jl) wlc(jl)=ws txt1cs=txt1c(il) txt2cs=txt2c(il) txt1c(il)=txt1c(jl) txt1c(jl)=txt1cs txt2c(il)=txt2c(jl) txt2c(jl)=txt2cs eflagcs=eflagc(il) eflagc(il)=eflagc(jl) eflagc(jl)=eflagcs endif enddo enddo write(17,7) Case 7 format(1x,'Lines in order of decreasing emission coefficient'/ : 1x,'Emission coefficients at a temperature of 10000K and' : ' an electron number density of 10000/cm^3 in Case 'a1/) write(17,43) 43 format(1x,'I = Line index') write(17,42) write(17,21) c do il=1,inl c if((emc(il)/emmax).gt.1.d-3) then c if(wlc(il).le.1.e5) then c write(17,5) il,wlc(il),eflagc(il),emc(il),emc(il)/emmax, c : emc(il)/ev1,txt1c(il) c else c write(17,51) il,wlc(il),eflagc(il),emc(il),emc(il)/emmax, c : emc(il)/ev1,txt1c(il) c endif c endif c enddo c write(17,6) c 6 format(//' Full list'//' In order of decreasing emissivity'//) c write(17,21) do il=1,inl if(wlc(il).le.1.e5) then if(emc(il)/ev1.gt.1.d-5.and.emc(il)/emmax.gt.1.d-5) then write(17,5) il,wlc(il),eflagc(il),emc(il),emc(il)/emmax, : emc(il)/ev1,txt1c(il) else write(17,31) il,wlc(il),eflagc(il),emc(il),emc(il)/emmax, : emc(il)/ev1,txt1c(il) endif else if(emc(il)/ev1.gt.1.d-5.and.emc(il)/emmax.gt.1.d-5) then write(17,51) il,wlc(il),eflagc(il),emc(il),emc(il)/emmax, : emc(il)/ev1,txt1c(il) else write(17,33) il,wlc(il),eflagc(il),emc(il),emc(il)/emmax, : emc(il)/ev1,txt1c(il) endif endif enddo call flush(16) call flush(17) return end ************************************************************************* subroutine opt2 implicit none * mxl = maximum number of line emission coefficients * mxt = maximum number of temperatures * mxd = maximum number of densities * mxout = maximum number of user requested emission coefficients integer mxl,mxt,mxd,mxout parameter (mxl=10000,mxt=25,mxd=16,mxout=30) character*84 txt1(mxl) character*32 txt2(mxl) character*30 infile character*34 outfile character*22 datname character*1 eflag(mxl),ef(mxl) character*1 ul,blk,Case real*8 dens(mxd),emms(mxt,mxd,mxout),wave(mxout), : em(mxt,mxd,mxl),wl(mxl),wv(mxout), : temlg(mxt),denslg(mxd),temlog(mxt),logdens(mxd), : tin(mxt),din(mxd),eps2d(mxt,mxd,mxout),eu(mxl), : el(mxl) real*8 xn,yt,TwoDnt,d,t,dt,dd,eps(3) integer nden,idens,nline,iline,k,ic,np,nt,nd,ntdin,itdin,ip, : ntem,item,iprint,i,n,nf(mxl),kc,isline,nw, : nsline,kline,ntin,ndin,itin,idin,status,access integer iopt,istop common/obs/wave,nline common/alltheory/wl,em,eflag,eu,el,txt1,txt2,temlg,denslg, : nt,nd,nw,Case common/theory/temlog,logdens,emms,ntem,nden data ul/'_'/, blk/' '/ write(6,*) 'List of OII wavelengths is required in file', : ' OIIdata_xxxx' write(6,*) 'Enter xxxx (1-22 characters)' read(5,*) datname infile='OIIdata'//ul//datname status=access(infile,'') if(status.ne.0) then write(6,*) 'File ',infile,' does not exist' return endif k=index(infile,blk) infile=infile(1:k-1) outfile=infile(1:k-1)//ul//'out' open(unit=15,status='unknown',file=infile) open(unit=20,status='unknown',file=outfile) rewind(15) * read * number of lines for which emission coeffs required - nline read(15,*) nline write(6,*) 'number of lines = ',nline * read line wavelengths * note that lines are identified by their wavelength as * printed in the lines.. files and given to 2 d.p. do iline=1,nline read(15,*) wave(iline) write(6,10) iline,wave(iline) 10 format(1x'line number',i5,' wavelength [A] ',f12.2) enddo close(15) istop=0 if(istop.eq.1) stop ntem=nt nden=nd do item=1,ntem temlog(item)=temlg(item) enddo do idens=1,nden logdens(idens)=denslg(idens) enddo do item=1,ntem do idens=1,nden * mark all lines as not found do iline=1,nline nf(iline)=1 enddo kline=0 do iline=1,nline do isline=1,nw if(abs(wl(isline)-wave(iline)) : .lt.0.001) then kline=kline+1 if(nf(iline).eq.0) then write(6,*) 'found two lines that match ', : wave(iline) c stop endif emms(item,idens,kline)=log10(em(item,idens,isline)) ef(kline)=eflag(isline) wv(kline)=wave(iline) nf(iline)=0 endif enddo enddo * check that all lines were found do iline=1,nline if(nf(iline).gt.0) then write(6,*) 'line ',wave(iline),' not found' endif enddo enddo enddo * and print for checking * do iline=1,nline * write(6,*) wave(iline) * do item=1,ntem * write(6,*) temlog(item),(emms(item,idens,iline), * : idens=1,ndens) * enddo * enddo istop=0 if(istop.eq.1) stop 1 write(6,20) 20 format(/'Choose one of two options:' // : '1: Provide a list of temperature/density pairs in a file' : ' TDpairs_xxxx' // : '2: Provide a range of temperatures and densities defining a' : ' 1 or 2 dimensional grid in a file TDranges_xxxx'// : 'Enter an integer') read(5,*) iopt if(iopt.ne.1.and.iopt.ne.2) go to 1 write(6,*) 'Enter xxxx (1-21 characters)' read(5,*) datname write(6,21) outfile 21 format(/1x,'Output will be in file ',a34/) if(iopt.eq.1) then infile='TDpairs'//ul//datname elseif(iopt.eq.2) then infile='TDranges'//ul//datname endif k=index(infile,blk) infile=infile(1:k-1) status=access(infile,'') if(status.ne.0) then write(6,*) 'File ',infile,' does not exist' return endif 30 format(1x,'Interpolation to log10(N[/cm^3]) = ',f7.3, : ' log10(T[K]) = ',f7.3///5x,'line','wavelength[A]', : ' Interp pts ',' coefficient [cm^3/s]') open(unit=18,status='unknown',file=infile) rewind(18) if(iopt.eq.1) then write(20,24) 24 format(/1x,'I = Line index'/ : 1x,'W = Wavelength[A],'/ : 1x,'EM = Emission coefficient[cgs],'/ : 1x,'X = (air)wavelength from experimental energies'/ : 1x,'T = one or both energies from theory,' : ' wavelength inaccurate') write(20,23) 23 format(/1x,'Where necessary, the emission coefficients' : ' are interpolated from'/ : 1x,'the theory mesh using 6-point Lagrange in' : ' temperature and density'/) read(18,*) ntdin do itdin=1,ntdin read(18,*) t,d xn=log10(d) yt=log10(t) c write(6,25) t,d write(20,25) t,d,Case 25 format(/1x,'Temperature = ',1pe12.3,' Density = ', : 1pe12.3,' Case = ',a1) c write(6,26) yt,xn write(20,26) yt,xn write(20,27) 26 format(1x,'log10(T) = ',f10.3,' log10(N) = ',f7.3/) 27 format(8x,'I',12x,'W',4x,'X/T',13x,'EM') if(xn.ge.denslg(1).and.xn.le.denslg(nden).and. : yt.ge.temlg(1).and.yt.le.temlg(ntem)) then do iline=1,kline do ip=1,3 eps(ip)=10.d0**TwoDnt(iline,xn,yt,ip) enddo c write(6,40) iline,wv(iline),ef(iline),(eps(ip),ip=1,3) write(20,40) iline,wv(iline),ef(iline),eps(3) enddo else write(6,35) write(20,35) 35 format(/1x,'Temperature or density outside available', : ' theory range'/) endif enddo 40 format(5x,i4,f13.2,5x,a2,1p3e15.3) elseif(iopt.eq.2) then write(20,23) read(18,*) ntin,ndin ctest write(6,*) 'ntin,ndin',ntin,ndin if(ntin.eq.1) then read(18,*) tin(1) dt=0.d0 elseif(ntin.gt.1) then read(18,*) tin(1),dt ctest write(6,*) tin(1),dt endif if(ndin.eq.1) then read(18,*) din(1) dd=0.d0 elseif(ndin.gt.1) then read(18,*) din(1),dd endif do itin=1,ntin tin(itin)=tin(1)+dt*(itin-1) yt=tin(itin) do idin=1,ndin din(idin)=din(1)+dd*(idin-1) xn=din(idin) if(xn.ge.logdens(1).and.xn.le.logdens(nden).and. : yt.ge.temlog(1).and.yt.le.temlog(ntem)) then do iline=1,kline ip=3 eps2d(itin,idin,iline)= : 10.d0**TwoDnt(iline,xn,yt,ip) enddo else write(6,35) write(20,35) return endif enddo enddo do iline=1,kline write(20,50) iline,wv(iline),Case 50 format(//1x,'Emission coefficients [cgs] for line ',i2, : ' Wavelength ',f10.2,' Case ',a1//) write(20,60) (din(idin),idin=1,ndin) 60 format(1x,'log10(T[K])',3x,'log10(N[cgs])=', : f7.3,29(4x,f6.3)) do itin=1,ntin write(20,70) tin(itin),(eps2d(itin,idin,iline), : idin=1,ndin) 70 format(5x,f7.3,14x,1p30e10.3) enddo enddo endif call flush(20) return end ************************************************************************* SUBROUTINE LAGR(F,X,FV,XV,NP,NI,ND,IER) IMPLICIT REAL*8(A-H,O-Z) C C SUBROUTINE TO PERFORM LAGRANGIAN INTERPOLATION OF A SET OF C TABULATED FUNCTION VALUES F(X(I)) C C INPUTS C F(I) = ARRAY OF TABULATED FUNCTION VALUES C X(I) = ARRAY OF CORRESPONDING X VALUES C NI = NUMBER OF TABULATED VALUES C 2*NP = NUMBER OF POINTS TO BE USED IN INTERPOLATION C XV = VALUE OF X AT WHICH FUNCTION VALUE IS REQUIRED C C OUTPUTS C FV = INTERPOLATED VALUE C IER = ERROR CODE (.NE.0 IF XV OUTSIDE RANGE OF X) C DIMENSION F(ND),X(ND),XIQ(20) LG=2*NP IER=0 IF(XV.LT.X(1).OR.XV.GT.X(NI)) THEN IER=1 RETURN ENDIF C C LOCATE XV IN X DO 5999 I=2,NI IF(X(I).GE.XV) THEN C C SELECT SUBSET OF X FOR INTERPOLATION, TAKING CARE AT ENDS IN=NI-LG+1 ILOW=I-NP IF(ILOW.LT.1) ILOW=1 IF(ILOW.GT.IN) ILOW=IN DO 5998 ITAU=1,LG XIQ(ITAU)=X(ILOW+ITAU-1) 5998 CONTINUE C C CONSTRUCT INTERPOLATION COEFFICIENTS FV=0.D0 DO 5997 ITAU=1,LG FL=PHI(XIQ,LG,ITAU,XV)/PHI(XIQ,LG,ITAU,XIQ(ITAU)) FV=FV+FL*F(ILOW+ITAU-1) 5997 CONTINUE RETURN ENDIF 5999 CONTINUE END ******************************************************************************* DOUBLE PRECISION FUNCTION PHI(XIQ,LG,ITAU,X) IMPLICIT REAL*8(A-H,O-Z) DIMENSION XIQ(20) FK=1.D0 DO 5999 L=1,LG IF(L.NE.ITAU) THEN FK=FK*(X-XIQ(L)) ENDIF 5999 CONTINUE PHI=FK RETURN END ****************************************************************************** double precision function OneDt(il,in,yt,np) implicit none integer mxt,mxd,mxout parameter (mxout=30,mxt=25,mxd=16) real*8 temlog(mxt),logdens(mxd),emms(mxt,mxd,mxout) real*8 ft(mxt),xt(mxt),fvt,yt integer il,in,it,ntem,nden,ier,np common/theory/temlog,logdens,emms,ntem,nden * write(6,*) 'enter 1Dt, yt ',yt do it=1,ntem ft(it)=emms(it,in,il) xt(it)=temlog(it) enddo * write(6,*) 'in 1Dt ntem,ft,xt',ntem,ft,xt call lagr(ft,xt,fvt,yt,np,ntem,mxt,ier) OneDt=fvt * write(6,*) 'in 1Dt yt,fvt,1Dt',yt,fvt,OneDt return end ********************************************************************************* double precision function OneDn(il,xn,yt,np) implicit none integer mxt,mxd,mxout parameter (mxout=30,mxt=25,mxd=16) real*8 temlog(mxt),logdens(mxd),fn(mxd),an(mxd), : emms(mxt,mxd,mxout) real*8 fvn,xn,yt,OneDt integer il,np,ntem,nden,in,ier common/theory/temlog,logdens,emms,ntem,nden * write(6,*) 'enter 1Dn, il, xn, yt',il,xn,yt do in=1,nden fn(in)=OneDt(il,in,yt,np) an(in)=logdens(in) enddo * write(6,*) 'in 1Dn before call lagr, xn=',xn * write(6,*) 'fn,an',fn,an call lagr(fn,an,fvn,xn,np,nden,mxd,ier) OneDn=fvn * write(6,*) 'return' return end ******************************************************************************************** double precision function TwoDnt(il,xn,yt,np) implicit none real*8 xn,yt,OneDn integer il,np * write(6,*) 'enter TwoDnt, il,xn,yt ',il,xn,yt TwoDnt=OneDn(il,xn,yt,np) return end ******************************************************************************************* subroutine reademiss * iopt = 1 read data from binary file 'lines.bin' * iopt = 2 read data from individual lines_xxx_yyy_dr files implicit none integer mxl,mxt,mxd parameter (mxl=10000,mxt=25,mxd=16) character*1 eflag(mxl),caseA,caseB,caseC,Case,keepCase character*4 ktem(mxt),kdens(mxd) character*5 lines character*12 namedata character*25 name character*80 text character*84 txt1(mxl) character*32 txt2(mxl) real*8 em(mxt,mxd,mxl),wl(mxl), : temlg(mxt),denslg(mxd), : eu(mxl),el(mxl),temp(mxd) integer ic,it,nt,id,nd,i,icase,status,nw,access common/alltheory/wl,em,eflag,eu,el,txt1,txt2,temlg,denslg, : nt,nd,nw,Case data temlg/2.0,2.1,2.2,2.3,2.4,2.5,2.6,2.7,2.8,2.9, : 3.0,3.1,3.2,3.3,3.4,3.5,3.6,3.7,3.8,3.9, : 4.0,4.1,4.2,4.3,4.4/ data denslg/2.0,2.2,2.4,2.6,2.8,3.0,3.2,3.4,3.6,3.8,4.0, : 4.2,4.4,4.6,4.8,5.0/ data caseA/'A'/,caseB/'B'/,caseC/'C'/,keepCase/' '/ data nt/25/, nd/16/, nw/8889/ * read line list and theoretical emission coeffs from binary * file 'OIIlines_AB' and store wavelengths and emission coeffs namedata='OIIlines_ABC' status=access(namedata,'') if(status.ne.0) then write(6,*) 'Data file ',namedata,' does not exist' stop endif open(unit=12,status='unknown',file=namedata) 1 write(6,10) 10 format(1x,'Specify Case (A, B or C)') read(5,*) Case if(Case.ne.caseA.and.Case.ne.caseB.and.Case.ne.caseC) go to 1 if(Case.eq.keepCase) return write(6,*) 'open and read file ',namedata,', please wait' do ic=1,18 read(12,*) enddo do ic=1,nw read(12,1001) txt1(ic),wl(ic),eflag(ic),eu(ic),el(ic), : txt2(ic) c if(ic.eq.8889) write(6,*) txt1(ic),wl(ic) enddo 1001 format(6x,a84,f11.2,4x,a1,1x,2f11.2,a32) do ic=1,nw do icase=1,3 do i=1,3 read(12,*) enddo do it=1,nt if(icase.eq.1.and.Case.eq.caseA.or. : icase.eq.2.and.Case.eq.caseB.or. : icase.eq.3.and.Case.eq.caseC) then read(12,1002) (em(it,id,ic),id=1,nd) else read(12,1002) (temp(id),id=1,nd) endif enddo 1002 format(26x,16(e10.3)) enddo enddo write(6,40) temlg(1),temlg(nt),denslg(1),denslg(nd) c write(20,10) temlg(1),temlg(nt),denslg(1),denslg(nd) 40 format(/1x,'Theoretical emission coefficients read for '/1x, : 'Case ',a1/1x, : 'log10(Temperature[K]) = ',16x,f6.2,' to ',f6.2, : ' with step 0.1'/1x, : 'log10(Electron number density[/cm^3]) = ',f6.2, : ' to ',f6.2,' with step 0.2'/) keepCase=Case close(12) return end