;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;; ; Document name: atomrd.pro ; Created by: Phil &, High Altitude Observatory/NCAR, Boulder CO, November 23, 1994 ; ; Last Modified: Mon May 1 10:53:24 1995 by judge (Phil &) on pika ;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;; ; PRO atomrd,file,coldata,head,plot=plot,info=info,make=make, describe=describe, $ help=help ;+ ; PROJECT: ; HAOS-DIAPER ; ; NAME: ; atomrd ; ; PURPOSE: ; Reads atomic data in standard MULTI format from file ; file, including collisional data and header information ; ; EXPLANATION: ; ; CALLING SEQUENCE: ; atomrd,file,coldata,head ; ; INPUTS: ; ; OPTIONAL INPUTS: ; None. ; ; OUTPUTS: ; Data are stored in common block catom (accessed via @catom) ; ; OPTIONAL OUTPUTS: ; coldata: structure containing collisional parameters ; head: structure containing header information ; ; KEYWORD PARAMETERS: ; /help. Will call doc_library and list header, and return ; plot=plot make a term diagram ; info=info lists header info ; make=make make the file if it does not exist with ground + next ion ; stage only ; /describe will describe variables in the file ; CALLS: ; None. ; ; COMMON BLOCKS: ; None. ; ; RESTRICTIONS: ; None. ; ; SIDE EFFECTS: ; None. ; ; CATEGORY: ; ; PREVIOUS HISTORY: ; Written November 23, 1994, by Phil &, High Altitude Observatory/NCAR, Boulder CO ; ; MODIFICATION HISTORY: ; M. Carlsson circa 1988: original coding ; P. Judge circa 1990: double precision variables used ; P. Judge 21-Jan-1993: collisional file is read if coldata in arg list ; P. Judge 29-Jul-1993: Bug corrected (TRAD=TRAD2 was commented out earlier) ; P. Judge 14-Jun-1994: common block replaced by @catom ; P. Judge 05-Jul-1994: header input added, order of jrad, irad input ; made unimportant, SDP prolog added ; P. Judge 04-Aug-1994: Bug corrected (jrad, irad for photoionization continua ; were set to input values minus 2, caused by 05-Jul-1994 fix) ; P. Judge 19-Aug-1994: Search done on SDP_ATOM directories if file ; : does not exist in current directory, also ; : if no '.' exists the input filename is assumed to ; : consist of the extension (e.g., atomrd,'o1') ; P. Judge November 23, 1994: /make option added ; P. Judge November 30, 1994: Bug fixed so that gencolrd is called ; when n_params(0) ge 2 (previously ; nparams(0) eq 2) ; P. Judge May 1, 1995: Added describe option ; VERSION: ; Version 1, November 23, 1994 ;- ; @catom IF(N_ELEMENTS(help) GT 0) THEN BEGIN doc_library,'atomrd' RETURN ENDIF ; ;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;; ; May 1, 1995 Modifications begin IF(N_ELEMENTS(describe) GT 0) THEN BEGIN PRINT,'ATOMIC PARAMETERS (SEE CARLSSON 1986 FPR MORE DETAILS)' PRINT,'NK : NUMBER OF LEVELS, CONTINUUM LEVELS INCLUDED' PRINT,'NRAD : NUMBER OF RADIATIVE TRANSITIONS' PRINT,'NLINE : NUMBER OF BOUND-BOUND TRANSITIONS' PRINT,'NCONT : NUMBER OF BOUND-FREE TRANSITIONS I DETAIL' PRINT,'LABEL : e.g., ''N IV 1S2 2S2 1SE 0'' ' PRINT,'G : STATISTICAL WEIGHT OF LEVEL' PRINT,'EV : ENERGY IN EV (INPUT IN CM-1)' PRINT,'F : FOR LINES : *ABSORPTION* OSCILLATOR STRENGTH' PRINT,' FOR CONTINUA: CROSS-SECTION AT LIMIT' PRINT,' THE WAVELENGTH DEPENDENCE OF THE CROSS-SECTION IS' PRINT,' ASSUMED TO BE A=A0*(NY0/NY)**3 IF NOT GIVEN EXPLICITLY' PRINT,'A : EINSTEIN A COEFFICIENT' PRINT,'B : B' PRINT,'' PRINT,'IWIDE : =.TRUE. FOR FREQUENCY DEPENDENT TRANSITIONS' PRINT,' (CONTINUA AND WIDE LINES)' PRINT,'KTRANS(K) TRANSITION NUMBER K IS FREQUENCY DEPENDENT' PRINT,' TRANSITION NUMBER KTRANS' PRINT,'NRFIX : NUMBER OF TRANSITIONS WITH FIXED RATES' PRINT,'FOR FIXED RATES VARIABLES SEE ROUTINE FIXRAD' PRINT,'IRAD(K): LOWER LEVEL FOR RADIATIVE TRANSITION K' PRINT,'JRAD(K): UPPER "' PRINT,'ALAMB : VACUUM WAVELENGTH IN ANGSTROM' PRINT,' IN PRINTOUT ROUTINES IT IS PRINTED AS EITHER VACUUM OR AIR' PRINT,' DEPENDING ON THE VALUE (.LT.2000 VACUUM, .GT.2000 AIR)' PRINT,'HNY4P : H*NY/4PI, NY IN UNITS OF A TYPICAL DOPPLER WIDTH' PRINT,'' ENDIF ; ; May 1, 1995 Modifications end ;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;; if (n_params(0) eq 0) then begin print,'atomrd,file [,coldata,head]' return endif ; ;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;; ; 19-Aug-1994 modifications begin ; top: fil=file fil1=fil j=findfile(fil,count=count) if(count ne 0) then goto, ok k=strpos(fil,'.') if(k eq -1) then fil='atom.'+fil j=findfile(fil,count=count) fil1=fil if(count eq 0) then begin atomfnd,fil,full,status = status if(status eq 0) then begin IF(N_ELEMENTS(make) NE 0) THEN begin k=STRPOS(fil,'.') l=STRLEN(fil) ext=STRMID(fil,k+1,l) l=STRLEN(ext) el_ion,ext,el,ion strt_atom,el,FIX(ion),status = status IF(status GT 0) then GOTO,top ELSE return ENDIF ELSE BEGIN PRINT,'atomrd: error- file does not exist in current or SDP_ATOM ' + $ 'directory' RETURN ENDELSE endif fil=full print,'atomrd: file not in current directory, reading file ', fil endif else begin fil=j(0) endelse ok: ; ; 19-Aug-1994 modifications end ;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;; ; ; header? ; if(n_params(0) eq 3) then rd_head,fil,head ,/dir ; ; strip comment lines from input file ; text=' ' openr,lu1,/get_lun,fil,err=err if(err ne 0) then begin print,'atomrd: error opening file',fil return endif openw,lu2,/get_lun,'dums.dat' while not eof(lu1) do begin readf,lu1,text if (strmid(text,0,1) ne '*') then begin k0=strpos(text,"'") if (k0 ne -1) then begin k1=strpos(text,"'",k0+1) k2=strlen(text) printf,lu2,strmid(text,0,k0) ; text up to label printf,lu2,strmid(text,k0+1,k1-k0-1) ; label printf,lu2,strmid(text,k1+1,k2-k1) ; text after label endif else begin printf,lu2,text endelse endif endwhile free_lun,lu1 close,lu2 openr,lu2,'dums.dat' ;read,'Give QNORM',qnorm qnorm=5.0 one=double(1.) EE=one*1.602189E-12 HH=one*6.626176E-27 CC=one*2.99792458E10 EM=one*9.109534E-28 UU=one*1.6605655E-24 BK=one*1.380662E-16 PI=one*3.14159265359 HCE=HH*CC/EE*1.E8 HC2=2.*HH*CC *1.E24 HCK=HH*CC/BK*1.E8 EK=EE/BK HNY4P=HH*CC/QNORM/4./PI*1.E-5 atomid='' readf,lu2, ATOMID atomid=strtrim(atomid,2) abnd = 0.0 awgt = 0.0 READF,LU2, ABND,AWGT AWGT=AWGT*UU nk=0 nline=0L ; longword ncont=0 nrfix=0 READF,LU2, NK,NLINE,NCONT,NRFIX NRAD=NLINE+NCONT if(float(nk)*float(nrad>1) gt 1.e4) then begin print,'This is a big atom- please be patient...' endif ev=dblarr(nk) g=fltarr(nk) label=strarr(nk) ion=intarr(nk) krad=intarr(nk,nk) lab='' for I=0,NK-1 do begin READF,LU2, ev2,g2 readf,lu2, lab readf,lu2, ion2 G(I)=g2 LABEL(I)=lab ION(I)=ion2 EV(I)=EV2*CC *HH/EE for J=0,NK-1 do begin KRAD(I,J)=0. endfor endfor ;C ;C BOUND-BOUND TRANSITIONS IN DETAIL ;C CALCULATE LAMBDA, A AND B ;C IF QMAX OR Q0.LT.0 FREQUENCY POINTS IN DOPPLER UNITS ARE READ ;C KT=-1 mq=100 if (nrad gt 0) then begin f=fltarr(nrad) nq=fltarr(nrad) qmax=fltarr(nrad) q0=fltarr(nrad) iwide=intarr(nrad) ga=fltarr(nrad) gw=fltarr(nrad) gq=fltarr(nrad) irad=intarr(nrad) jrad=intarr(nrad) alamb=dblarr(nrad) a=fltarr(nrad) b=fltarr(nk,nk) ktrans=intarr(nrad) q=fltarr(mq,nrad) endif if (nline gt 0) then begin for KR=0L,NLINE-1 do begin readf,lu2, j,i,f2,nq2,qmax2,q02,io2,ga2,gw2,gq2 up = j > i dn = j < i i=dn-1 j=up-1 f(kr)=f2 nq(kr)=nq2 qmax(kr)=qmax2 q0(kr)=q02 iwide(kr)=io2 ga(kr)=ga2 gw(kr)=gw2 gq(kr)=gq2 IF(QMAX(KR) LT 0.0) OR (Q0(KR) LT 0.0) THEN begin dum=fltarr(nq2) READF,LU2, dum q(0,kr)=dum ENDIF ; ; PGJ modification 27 Jan 1992 ; ; KRAD(I,J)=KR ; KRAD(J,I)=KR KRAD(I,J)=KR+1 KRAD(J,I)=KR+1 ; PGJ end modification 27 Jan 1992 ; IRAD(KR)=I+1 JRAD(KR)=J+1 ALAMB(KR)=HCE/(EV(J)-EV(I)) A(KR)=F(KR)*6.671E15*G(I)/G(J)/ALAMB(KR)/ALAMB(KR) B(J,I)=ALAMB(KR)^3/HC2*A(KR) B(I,J)=G(J)/G(I)*B(J,I) endfor endif ;C ;C BOUND-FREE TRANSITIONS IN DETAIL ;C IF QMAX.LT.0.0 FREQUENCY POINTS IN ANGSTROM (STARTING WITH ;C THRESHOLD AND DECREASING) AND CROSSECTIONS IN CM2 ARE READ. ;C UNIT CONVERSION IN ROUTINE FREQC ;C mq=100 if (nrad gt nline) then begin alfac=fltarr(mq+1,nrad-nline) frq=alfac for KR=NLINE,NRAD-1 do begin READF,LU2, J,I,f2,nq2,qmax2 up = j > i dn = j < i i=dn-1 j=up-1 ; ; PGJ 04-AUG-94 bug affecting photoionization: begin ; continua ; i=i-1 ; j=j-1 ; ; PGJ 04-AUG-94 bug affecting photoionization: end ; kt=kt+1 F(KR)=f2 alfac(0,kr-nline)=f2 NQ(KR)=nq2 QMAX(KR)=qmax2 IF(QMAX(KR) LT 0.0) THEN begin tab=fltarr(2,nq2) READF,LU2, tab tab=transpose(tab) q(0,kr)=hh*cc/ee/(ev(j)-ev(i))*1.e8 ; wavelength at edge q(1,kr)=tab(*,0) ; wavelength table frq(0,kt)=cc*1.e8/q(0:nq(kr),kr) ; frequencies alfac(1,kt)=tab(*,1) ENDIF else begin q(0,kr)=hh*cc/ee/(ev(j)-ev(i))*1.e8 ; wavelength at edge frq(0,kt)=cc*1.e8/q(0,kr) ; frequency at edge frq(nq(kr),kt)=cc*1.e8/qmax(kr) ; frequency at shortest wavelength frq(0,kt)=findgen(nq(kr)+1)/nq(kr)*(frq(nq(kr),kt)-frq(0,kt))+frq(0,kt) q(0,kr)=cc*1.e8/frq(0:nq(kr),kt) alfac(0,kt)=alfac(0,kt)*(frq(0,kt)/frq(0:nq(kr),kt))^3 endelse KTRANS(KR)=KR-NLINE+KT IRAD(KR)=I+1 JRAD(KR)=J+1 KRAD(I,J)=KR KRAD(J,I)=KR GA(KR)=0. GW(KR)=0. GQ(KR)=0. ALAMB(KR)=HCE/(EV(J)-EV(I)) endfor mq=max(nq)+1 alfac=alfac(0:mq,*) frq=frq(0:mq,*) endif ;C ;C INPUT FIXED TRANSITIONS PARAMETERS ;C if (nrfix gt 0) then begin jfx=intarr(nrfix) ifx=intarr(nrfix) ipho=intarr(nrfix) a0=fltarr(nrfix) trad=fltarr(nrfix) itrad=intarr(nrfix) for KF=0,NRFIX-1 do begin readf,lu2, jfx2,ifx2,ipho2,a02,trad2,itrad2 up = jfx2 > ifx2 dn = jfx2 < ifx2 ifx2=dn-1 jfx2=up-1 JFX(KF)=jfx2 IFX(KF)=ifx2 IPHO(KF)=ipho2 A0(KF)=a02 TRAD(KF)=trad2 ITRAD(KF)=itrad2 IF(ITRAD(KF) EQ 4) THEN begin IF(IPHO(KF) EQ 0) THEN begin readf,lu2,j,i,ff2,nqf2,qmaxf2,q0f2,io2,gaf2,gwf2,gqf2 IF(qmaxf2 LT 0.0) OR (q0f2 LT 0.0) THEN begin tab=fltarr(nqf2) READF,LU2, tab ENDIF endif ELSE begin READF,LU2, J,I,FF2,NQF2,QMAXF2 IF(QMAXF2 LT 0.0) THEN begin tab=fltarr(2*nqf2) READF,LU2, tab ENDIF ENDelse endif endfor endif free_lun,lu2 if(n_params(0) ge 2) then gencolrd,coldata,filen=fil if(n_elements(q) gt 0) then q=q(0:max(nq)-1,*) ; ; plotting? ; if(n_elements(plot) ne 0) then termdiag if(n_elements(info) ne 0) then atom_info,fil1 RETURN END ;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;; ; End of 'atomrd.pro'. ;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;