impreoved for extinction + emission fits

Jean-Philippe Bernard
1 parent bafa720f
Showing 3 changed files with 30 additions and 12 deletions Show diff stats
src/idl/dustem_fit_sed_readme.pro
src/idl/dustem_grains_colors.pro
src/idl/dustem_plot_ext.pro
@@ -92,6 +92,8 @@ dustem_init,model=use_model
 ;dir=getenv('DUSTEM_SOFT_DIR')+'/src/dustem3.8_web/SEDs/'
 dir=!dustem_wrap_soft_dir+'/Data/SEDs/'
 file=dir+'Gal_composite_spectrum.xcat'
+;file=dir+'Mathis1990_DISM.xcat'
+
 ;file=dir+'Gal_composite_spectrum_wAKARI.xcat'
 ;stop
 spec=read_xcat(file,/silent)
@@ -12,7 +12,7 @@ ENDIF ELSE BEGIN
 ENDELSE
  
  
-ls = [ [ 0,4,3,5,1 ], [ 0,4,3,5,1 ],[ 0,4,3,5,1 ] ]   ;JPB: unsure what this is ....
+ls = [ [ 0,4,3,5,1 ], [ 0,4,3,5,1 ],[ 0,4,3,5,1 ] ]   ;JPB: This is the line style for curves
  
 IF keyword_set(ps) THEN BEGIN
 	; Postscript
@@ -98,8 +98,11 @@ endif else if keyword_set(albedo) then begin
 	yr=[0,1]
 	xtit=textoidl('1/\lambda (\mum^{-1})')
 	ytit=textoidl('albedo')
-	plot,xrange,/nodata,xtit=xtit,ytit=ytit,tit=tit,xlog=xlog,ylog=ylog,yr=yr,ys=ys,xr=xr,xsty=xsty
-	oplot,xrange,sca_tot/(sca_tot+abs_tot)
+	stop
+	;plot,xrange,/nodata,xtit=xtit,ytit=ytit,tit=tit,xlog=xlog,ylog=ylog,yr=yr,ys=ys,xr=xr,xsty=xsty
+	;oplot,xrange,sca_tot/(sca_tot+abs_tot)
+	cgplot,xrange,/nodata,xtit=xtit,ytit=ytit,tit=tit,xlog=xlog,ylog=ylog,yr=yr,ys=ys,xr=xr,xsty=xsty
+	cgoplot,xrange,sca_tot/(sca_tot+abs_tot)
  
  
 endif else begin
@@ -114,8 +117,12 @@ endif else begin
 	endelse
 	endif
  
-	Ngrains=n_tags(st.sed)-2
-	colors=dustem_grains_colors(Ngrains,ls,ps=ps)
+	;Ngrains=n_tags(st.sed)-2
+	;colors=dustem_grains_colors(Ngrains,ls,ps=ps)
+	Ngrains=(*!dustem_params).Ngrains
+	colors=['dark',dustem_grains_colors(Ngrains,ls,/cgplot,ps=ps)]
+	use_lines=replicate(0,Ngrains)
+
 ;	stop
 	ytit=textoidl('\sigma_{ext}/N_H (cm^2/H)') 
 	if not keyword_set(ps) then tit='Dustem extinction'
@@ -147,21 +154,27 @@ endif else begin
         	xrange=st.ext.wav
 		xlog = 1
 	endelse
-	if multi eq 0 then plot,xrange,/nodata,xtit=xtit,ytit=ytit,tit=tit,xlog=xlog,ylog=ylog,yr=yr,ys=ys,xr=xr,xsty=xsty
-	if multi eq 1 then plot,xrange,/nodata,xtit='',ytit=ytit,tit=tit,xlog=xlog,ylog=ylog,yr=yr,/ys,xr=xr,xsty=xsty,position=[0.17,0.35,0.95,0.95],xtickformat='(A1)',yticks=3,yminor=10
-	if multi eq 2 then plot,xrange,/nodata,xtit=xtit,xlog=xlog,ylog=ylog,ys=ys,xsty=xsty,ytit='Normalized',tit='',xr=xr,yr=[0,2],position=[0.17,0.14,0.95,0.35],/noerase,yticks=2,yminor=5,xticklen=0.1
+	;if multi eq 0 then plot,xrange,/nodata,xtit=xtit,ytit=ytit,tit=tit,xlog=xlog,ylog=ylog,yr=yr,ys=ys,xr=xr,xsty=xsty
+	;if multi eq 1 then plot,xrange,/nodata,xtit='',ytit=ytit,tit=tit,xlog=xlog,ylog=ylog,yr=yr,/ys,xr=xr,xsty=xsty,position=[0.17,0.35,0.95,0.95],xtickformat='(A1)',yticks=3,yminor=10
+	;if multi eq 2 then plot,xrange,/nodata,xtit=xtit,xlog=xlog,ylog=ylog,ys=ys,xsty=xsty,ytit='Normalized',tit='',xr=xr,yr=[0,2],position=[0.17,0.14,0.95,0.35],/noerase,yticks=2,yminor=5,xticklen=0.1
+	if multi eq 0 then cgplot,xrange,/nodata,xtit=xtit,ytit=ytit,tit=tit,xlog=xlog,ylog=ylog,yr=yr,ys=ys,xr=xr,xsty=xsty
+	if multi eq 1 then cgplot,xrange,/nodata,xtit='',ytit=ytit,tit=tit,xlog=xlog,ylog=ylog,yr=yr,/ys,xr=xr,xsty=xsty,position=[0.17,0.35,0.95,0.95],xtickformat='(A1)',yticks=3,yminor=10
+	if multi eq 2 then cgplot,xrange,/nodata,xtit=xtit,xlog=xlog,ylog=ylog,ys=ys,xsty=xsty,ytit='Normalized',tit='',xr=xr,yr=[0,2],position=[0.17,0.14,0.95,0.35],/noerase,yticks=2,yminor=5,xticklen=0.1
 	;VG
  
 	;VG
 	; Plot model
 	FOR i=0L,Ngrains-1 DO BEGIN
 		if i le 1 and keyword_set(mergePAH) then begin
-  			if i eq 0 then oplot,xrange,(st.ext.abs_grain(0)+st.ext.sca_grain(0)+st.ext.abs_grain(1)+st.ext.sca_grain(1))/norm,color=colors(1),psym=psym,line=ls(1) 
+  			;if i eq 0 then oplot,xrange,(st.ext.abs_grain(0)+st.ext.sca_grain(0)+st.ext.abs_grain(1)+st.ext.sca_grain(1))/norm,color=colors(1),psym=psym,line=ls(1) 
+  			if i eq 0 then cgoplot,xrange,(st.ext.abs_grain(0)+st.ext.sca_grain(0)+st.ext.abs_grain(1)+st.ext.sca_grain(1))/norm,color=colors(1),psym=psym,line=ls(1) 
 		endif else begin
-  			oplot,xrange,(st.ext.abs_grain(i)+st.ext.sca_grain(i))/norm,color=colors(i+1),psym=psym,line=ls(i+1)
+  			;oplot,xrange,(st.ext.abs_grain(i)+st.ext.sca_grain(i))/norm,color=colors(i+1),psym=psym,line=ls(i+1)
+  			cgoplot,xrange,(st.ext.abs_grain(i)+st.ext.sca_grain(i))/norm,color=colors(i+1),psym=psym,line=ls(i+1)
 		endelse
 	ENDFOR
-	oplot,xrange,(abs_tot+sca_tot)/norm,psym=psym
+	;oplot,xrange,(abs_tot+sca_tot)/norm,psym=psym
+	cgoplot,xrange,(abs_tot+sca_tot)/norm,psym=psym
  
 	if !dustem_show_plot EQ 2 or multi eq 2 then begin
 		norm = interpol(st.ext.ext_tot,st.ext.wav,(*!dustem_data.ext).wav)
@@ -173,8 +186,11 @@ endif else begin
 	if data_exists then begin
    		if ptr_valid(!dustem_data.ext) then begin
       			if keyword_set(UV) then xdata = 1./(*!dustem_data.ext).wav else xdata = (*!dustem_data.ext).wav
-      			if keyword_set(noerrbars) then oplot,xdata,(*!dustem_data.ext).values/norm,psym=4 $
+      			;if keyword_set(noerrbars) then oplot,xdata,(*!dustem_data.ext).values/norm,psym=4 $
+      			if keyword_set(noerrbars) then cgoplot,xdata,(*!dustem_data.ext).values/norm,psym=4 $
+        		;else oploterror,xdata,(*!dustem_data.ext).values/norm,xdata*0.,(*!dustem_data.ext).sigma/norm,psym=4
         		else oploterror,xdata,(*!dustem_data.ext).values/norm,xdata*0.,(*!dustem_data.ext).sigma/norm,psym=4
+        		;oploterror,xdata,(*!dustem_data.ext).values/norm,xdata*0.,(*!dustem_data.ext).sigma/norm,psym=4
    		endif
 	endif