; file: findgf_chianti.pro ; init: Nov 5 2015 Rob Rutten Deil ; last: Nov 2 2019 Rob Rutten Deil ; note: enter "setssw chianti" in ed /home/rutten/bin/idl ; adapt line selection at end of this file and run per Hyper C pro findgf_chianti,wavaa,elemnr,ionstage,$ gftot,glow,gup,$ air=air,marginaa=marginaa,quiet=quiet ;+ ; NAME: ; findgf_chianti ; PURPOSE: ; determine line gf and glow, gup from Chianti database ; CALL: ; findgf_chianti,wavaa,elemnr=elemnr,ionstage=ionstage ; INPUTS: ; wavaa: wavelength in Angstrom (vacuum unless /air set) ; elemnr: element nr (H = 1) ; ionstage: ionization stage, neutral = 0 ; OPTIONAL KEYWORD INPUTS: ; air 1/0: air wavelength given instead of vacuum wavelength ; marginaa: search margin in Angstrom (default 0.1) ; quiet 1/0: less/more printout ; OUTPUTS: ; gftot = total gf value for all components together ; glow = total g for all lower levels ; gup = total g for all upper levels ; OPTIONAL OUTPUT ; print these values if quiet ne 1 ; HISTORY: ; Nov 5 2015 RR: start ;- ; answer no-parameter query if (n_params() lt 1) then begin print,' findgf_chianti,wavaa,elemnr=elemnr,ionstage=ionstage,air=air,$' print,' marginaa=marginaa,quiet=quiet)' return endif ; defaults for keywords if (not keyword_set(air)) then air=0 if (not keyword_set(marginaa)) then marginaa=1. if (not keyword_set(quiet)) then quiet=0 ; various checks if (elemnr gt 30) then begin print,' #### ABORT findgf_chianti: element not in Chianti (only 1-30)' return endif if (elemnr ne -1 and ionstage gt elemnr-1) then begin print,' #### ABORT findgf_chianti: ionstage > elemnr' return endif ; convert air to vacumm wavelength and reverse if (air) then airtovac,wavaa,vacwavaa else vacwavaa=wavaa if (not air) then vactoair,wavaa,airwavaa else airwavaa=wavaa ; define Chianti element and ion stage label arrays element=['H','He','Li','Be','B','C','N','O','F','Ne','Na',$ 'Mg','Al','Si','P','S','Cl','Ar','K','Ca','Sc','Ti',$ 'V','Cr','Mn','Fe','Co','Ni','Cu','Zn'] ;RR first 30 elements ; construct .elvlc and .wgfa file names elem=strlowcase(element[elemnr-1]) ion=elem+'_'+ntostr(ionstage+1) elvlcfile=!xuvtop+'/'+elem+'/'+ion+'/'+ion+'.elvlc' wgfafile=!xuvtop+'/'+elem+'/'+ion+'/'+ion+'.wgfa' ; print file names if (quiet ne 1) then print,' ==== species: '+ion if (quiet ne 1) then print,' ==== .elvlc file (use ew): '+elvlcfile if (quiet ne 1) then print,' ==== .wgfa file (use ew): '+wgfafile ; check Chianti file existence (not when scanning all) if (not file_test(elvlcfile)) then begin if (quiet ne 1) then print,$ ' #### ABORT findgf_chianti: no file '+elvlcfile return endif ; read elvlc file (_direct also returns mult although it doesn't say so) ; eg: ew /home/rutten/rr/ssw/packages/chianti/dbase/si/si_4/si_4.elvlc read_elvlc_direct,elvlcfile,l1,term,conf,ss,ll,spd,jj,mult,$ ecm,eryd,ecmth,erydth,ref ; read wgfa file read_wgfa2, wgfafile, lvl1, lvl2, wvl, gf, Avalue, ref ; select the line components wavsel=where(abs(wvl-vacwavaa) lt marginaa) ; get the total gf gftot=total(gf[wavsel]) ; get ranges lower and upper levels minlow=min(lvl1[wavsel]) maxlow=max(lvl1[wavsel]) minup=min(lvl2[wavsel]) maxup=max(lvl2[wavsel]) ; get total g values glow=total(mult[minlow:maxlow]) gup=total(mult[minup:maxup]) ; print results if (quiet ne 1) then begin print,' ==== number components = '+ntostr(n_elements(wavsel)) print,' ==== total gf = '+ntostr(gftot) print,' ==== glow = '+ntostr(glow)+' gup = '+ntostr(gup) endif end ; ================== main to test per IDLWAVE Hyper-c ==================== elemnr=26 ; AIA FeIX 171 ionstage=8 air=0 wavaa=171.07 margin=0.1 elemnr=26 ; AIA FeXII 193 ionstage=11 air=0 wavaa=193.51 margin=0.05 ;; elemnr=2 ; HeII 304 ;; ionstage=1 ;; air=0 ;; wavaa=304 ;; margin=0.5 ;; elemnr=1 ; HI 1216 = Lyalpha ;; ionstage=0 ;; air=0 ;; wavaa=1215.67 ;; margin=1. ;; elemnr=1 ; HI 6563 = Halpha ;; ionstage=0 ;; air=1 ;; wavaa=6563. ;; margin=1. findgf_chianti,wavaa,elemnr,ionstage,air=air,$ marginaa=margin end