pro for_pop_photo_calc,iz,ion,outdir=outdir
;+
; Project     : ISSI - Cavities
;
; Name        :  for_pop_photo_calc
;
; Purpose:
; 		This program calculates the level populations of an ion as a function
; 		of density, temperature, and distance from the limb, including photoexcitation.
; 		Results are saved in a file.
;
; 		It reaches 3 solar radii, with density between Log Ne=5 and 14 and temperatures
; 		between 4.0 and 7.5.
;
; Inputs:
;	iz - atomic number
;	ion
;
;  Keyword Outputs:
;		Outdir - output directory
;
; Calls: emiss_calc_enrico, findfile
;
; History: Based on a program by Enrico Landi.
;		12-Nov-2013 Adapted by T. Kucera
;		15-Dec-2015 Replaced () by [] for arrays. Also added keyword OUTDIR T. Kucera
;		20-Jun-2016 Bug fix
;-

default,outdir,''
if not dir_exist(outdir) then begin
		mk_dir,outdir,err=err
		if err ne '' then begin
			message,/cont,'Returning without performing calculation'
			return
		endif
endif


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']


filesave=element(iz-1)+'_'+strcompress(string(ion,'(i5)'),/remove_all)+'_pop_photo'
								;radial distance from surface
rphot=[1.001,1.01,1.02,1.05,1.10,1.25,1.5,2.0,2.5,3.0]
nphot=n_elements(rphot)
ndens=10
dens=5.0+findgen(ndens)
ntemp=71.
temp=4.0+0.05*findgen(ntemp)
radtemp=5770.0
ioneq=!xuvtop+'/ioneq/chianti.ioneq'
					;writes out 'populations.sav'
emiss=emiss_calc_enrico(iz,ion,dens=dens,temp=temp,ioneq_file=ioneq,/quiet)
	
restore,'populations.save'
n_levels=n_elements(pop[0,0,*])

pops=dblarr(ntemp,ndens,nphot,n_levels)

for i=0,nphot-1 do begin
   emiss=emiss_calc_enrico(iz,ion,dens=dens,rphot=rphot(i),temp=temp,ioneq_file=ioneq,/quiet)
   restore,'populations.save'
   pops[*,*,i,*]=pop[*,*,*]

endfor

ff = findfile(concat_dir(!xuvtop,'VERSION'))
chianti_version = ' '
IF  ff(0)  NE ''  THEN BEGIN 
   openr, 1, ff(0)
   readf, 1, chianti_version
   close, 1
ENDIF ELSE BEGIN 
   print,  'Please update your CHIANTI database version, which is older than 4.0'
   chianti_version = ' '
END 


outstruct={element:element(iz-1),pops:pops,temp:temp,dens:dens,radtemp:radtemp,radius:rphot,ioneq:ioneq,cversion:chianti_version}
;should add ioneq
savegen,file=concat_dir(outdir,filesave),outstruct
;print,'Fatto lo ione: ',iz,ion,format='(a16,2i3)'


return
end