dustem_variational_study_readme.pro
3.87 KB
1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
37
38
39
40
41
42
43
44
45
46
47
48
49
50
51
52
53
54
55
56
57
58
59
60
61
62
63
64
65
66
67
68
69
70
71
72
73
74
75
76
77
78
79
80
81
82
83
84
85
86
87
88
89
90
91
92
93
94
95
96
97
98
99
100
101
102
103
104
105
106
107
108
109
110
111
112
113
114
115
116
117
118
119
120
121
122
123
124
125
126
127
128
129
130
131
132
133
134
135
136
PRO dustem_variational_study_readme,model=model,help=help,png=png,postscript=postscript
;+
; NAME:
; dustem_variational_study_readme
; PURPOSE:
; This is an example of how to run a variational study on the dust model parameters with the dustem wrapper.
; CATEGORY:
; Dustem
; CALLING SEQUENCE:
; dustem_variational_study_readme[,model=][,/help][,/png][,/postscript]
; INPUTS:
; None
; OPTIONAL INPUT PARAMETERS:
; None
; OUTPUTS:
; None
; OPTIONAL OUTPUT PARAMETERS:
; None
; ACCEPTED KEY-WORDS:
; model = Selects one of the dust mixture used by dustem
; 'COMPIEGNE_ETAL2010' from Compiegne et al 2010 (default)
; 'DBP90' from Desert et al 1990
; 'DL01' from Draine & Li 2001
; 'DL07' from Draine & Li 2007
; postscript = if set, plot is done in DUSTEM/Docs/Figures/Last_dustem_fit.ps
; help = If set, print this help
; COMMON BLOCKS:
; None
; SIDE EFFECTS:
; None
; RESTRICTIONS:
; The dustem fortran code must be installed
; The dustem idl wrapper must be installed
; PROCEDURE:
; None
; EXAMPLES
; dustem_variational_study_readme,model='DBP90',/postscript
; dustem_variational_study_readme,model='COMPIEGNE_ETAL2010',/postscript
; MODIFICATION HISTORY:
; Written by J.P. Bernard July 15 2021
; see evolution details on the dustem gitlab
; Contact J.-Ph. Bernard (Jean-Philippe.Bernard@irap.omp.eu) in case of problems.
;-
;This Readme shows how to run dustem from the wrapper
IF keyword_set(help) THEN BEGIN
doc_library,'dustem_variational_study_readme'
goto,the_end
ENDIF
;loadct,13
IF keyword_set(model) THEN BEGIN
use_model=strupcase(model)
ENDIF ELSE BEGIN
use_model='COMPIEGNE_ETAL2010' ;Default is last dustem model
ENDELSE
;== This just initializes dustem wrapper
dustem_init,model=use_model
!dustem_verbose=1
!dustem_show_plot=1
IF keyword_set(png) THEN BEGIN
dir_png=!dustem_dat+'/Figures/'
force_mkdir,dir_png
ENDIF
IF keyword_set(postscript) THEN BEGIN
dir_ps=!dustem_dat+'/Figures/'
force_mkdir,dir_ps
ENDIF
IF not keyword_set(dir) THEN BEGIN
; CASE getenv('DUSTEM_WHICH') OF
CASE !dustem_which OF
'DESERT':dir='DESERT_POST_ISO_MORE/'
'COMPIEGNE': dir='MC_DAT/'
'VERSTRAETE': dir='d_3.5/'
'VERSTRAETE': dir='LV_DAT/'
'WEB3p8': dir_in=!dustem_soft_dir
ELSE: dir='les_DAT/'
ENDCASE
ENDIF
;== The following structure contains the inputs to the model
st_model=dustem_read_all(dir_in)
;These are the values of G0 to be explored
g0s=[1.,2.,5.,10.,100.,1000.]
n_g0s=n_elements(g0s)
win=0L
window,win,title='DUSTEMWrapper Variational study'
CASE use_model OF
'DBP90': pd = ['(*!dustem_params).G0'] ;This is the parameter description of G0 for models without gas parameters
'COMPIEGNE_ETAL2010': pd = ['(*!dustem_params).gas.G0'] ;This is the parameter description of G0 for models with gas parameters
'AJ13' : pd = ['(*!dustem_params).G0'] ;This is the parameter description of G0 for models without gas parameters
ENDCASE
FOR i=0L,n_g0s-1 DO BEGIN
iv = [g0s[i]]
dustem_init_parinfo,pd,iv
params=[g0s[i]]
st=dustem_run(params)
;do the plot
cleanplot
xtit=textoidl('\lambda (\mum)')
;=== same plot as in Compiegne et al. 2010
ytit=textoidl('\nuI_\nu^{em} (W/m^2/sr for N_H=1.e20 H/cm^2)')
tit='DUSTEMWrapper Emitted Intensity '+use_model+' G0='+strtrim(g0s[i],2)
yr=[1e-11,1.e-4]
xr=[1,5e3]
;stop
dustem_plot_nuinu_em,st,yr=yr,/ysty,xr=xr,/xsty,/xlog,/ylog,title=tit,xtit=xtit,ytit=ytit
IF keyword_set(postscript) THEN BEGIN
file_ps=dir_ps+use_model+'_variational_study_nuinuem_G0='+strtrim(g0s[i],2)+'.ps'
spawn,'rm '+file_ps
previous_device=!d.name
set_plot,'PS'
device,filename=file_ps
!p.charsize=0.8
dustem_plot_nuinu_em,st,yr=yr,/ysty,xr=xr,/xsty,/xlog,/ylog,title=tit,xtit=xtit,ytit=ytit
device,/close
set_plot,previous_device
message,'Wrote '+file_ps,/continue
ENDIF
ENDFOR
;stop
the_end:
END