;+
; Project     : ISSI - Cavities
;
; Name		 : for_spec_calc
;
; Purpose    :  precalculate spectra used to calculate the responses for the imaging instruments. 
;                     calls for_1spec_calc (also in this file) for a series of different densities.
; Inputs Keywords:
;				 Abund: abundance file(s)
;				 IonEq - Ionization Equilibrium file
;				 OutDir - directory for output
; Outputs: 
;				units - units of spectra
;
; History:
; 	16-Nov-2011 T. Kucera
; 	24-Mar-2013 Output file name now specifies abundance and ion equilibrium 
;				files and chianti version number. TAK
; 	7-Apr-2013 Default abundance now sun_coronal_feldman_1992_ext.abund (different name for sun_coronal_ext.abund)
; 	18-Dec-2015 Added OutDir keyword. 
;	                      Allow the abundance file to have a specified directory other than the Chianti dbase abundance directory
;						  Increased efficiency by allowing multiple abundance files.
;						  TAK
;  29-Feb-2016. Bug fixes TAK
;  16-Jun-2016  Changing the wavelength range to 1-901 A, Fixed bugs TAK
;-


pro for_1spec_calc,Emissivity,lambda_s,density,LogT,ioneq=ioneq0,abund=abund0,units,outdir=outdir,outfile=outfile0
;
;  Purpose: Calculate a spectrum for a range of temperatures
;
;Inputs
;	      density
;Input Keywords
;           Abund: abundance file(s)
;          IonEq - Ionization Equilibrium file
;          OutDir - directory for output
;          OutFile - output filename
;Outputs
;          Emissivity
;          lambda_s
;          LogT
;          density (if not an inout)
;          units - units of emissivity

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 

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

default,density,1e10
default,ioneq0,'chianti.ioneq'
default,abund0,'sun_coronal_1992_feldman_ext.abund'
break_file,ioneq0, ioneq_disk_log, ioneq_dir, ioneq_nm, ioneq_ext
if n_elements(outfile0) gt 0 then if n_elements(outfile0) ne NAbund then $
		message,'Number of different outfile names does not match number of abundance files to be used'
NAbund=n_elements(abund0)
NTemp=201
NLam=901
LogT=findgen(NTemp)/(NTemp-1)*4+4.

lambda_s=findgen(NLam)+1
WMin=min(lambda_s)
WMax=max(lambda_s)


if ioneq_dir eq '' then ioneq_dir=concat_dir(!xuvtop,'ioneq') else ioneq_dir=ioneq_disk_log+ioneq_dir
ioneq=concat_dir(ioneq_dir,ioneq_nm+ioneq_ext)

PRINT,'NAbund: ',NAbund,' ', abund0
PRINT,'ioneq: ',ioneq
PRINT,'Wavelength Range: ',Wmin,' - ',WMax,' A'
abund_nms=strarr(NAbund)

Emissivity =fltarr(n_elements(lambda_s),n_elements(LogT),Nabund)
								
for i=0,NTemp-1 do begin
	if (LogT[i] mod 1) eq 0 then print,'LogT: ',LogT[i],systime()
	ch_synthetic,WMin,WMax,logt=LogT[i],density=density,output=syn_struct,$
		ioneq=ioneq
	for j=0,nabund-1 do begin
			break_file,abund0[j], abund_disk_log, abund_dir0, abund_nm, abund_ext
			if abund_dir0 eq '' then abund_dir=concat_dir(!xuvtop,'abundance') else abund_dir=abund_disk_log+abund_dir0
			abund=concat_dir(abund_dir,abund_nm+abund_ext)
			abund_nms[j]=abund_nm
			print,'Abund: ',abund
		    make_chianti_spec,syn_struct, lambda_s,spec_struct, /CONTINUUM,$
					BIN_SIZE=1.,abund_name=abund,/photons
		    Emissivity[*,i,j]=spec_struct.spectrum
	endfor  ;end j (abundance loop)
endfor      ;end i (temperature loop)
units=spec_struct.units[1]

;save,file=outfile,Emissivity,lambda_s,density,LogT,ioneq,abund,units,chianti_version
for j=0,nabund-1 do begin
	case 1 of
		n_elements(outfile0) eq 1: outfile=outfile0
		n_elements(outfile0) gt 1: outfile=outfile0[j]
		else: outfile='spectrum_'+trim(fix(WMin))+'t'+trim(fix(WMax))+'A_ne'+trim(fix(alog10(density)))+'_'+$
							str_replace(ioneq_nm,'_')+'_'+str_replace(abund_nms[j],'_')+'_cv'+chianti_version
	endcase
		spec={Emissivity:reform(Emissivity[*,*,j]),lambda:lambda_s,density:density,LogT:LogT,$
			ioneq:ioneq,abund:abund,units:units,chianti_version:chianti_version}
	savegen,file=concat_dir(outdir,outfile),spec
endfor
end

pro for_spec_calc,ioneq=ioneq,abund=abund,outdir=outdir,dens=dens

	default,dens,10.^[7,8,9,10,11.]
	for i=0,n_elements(dens)-1 do for_1spec_calc,Emissivity,lambda_s,dens[i],ioneq=ioneq,abund=abund,outdir=outdir
end