dustem_make_polarization_sed_example.pro
5.49 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
137
138
139
140
141
142
143
144
145
146
147
148
149
150
151
152
153
154
155
156
157
158
159
160
161
162
163
164
165
166
167
168
169
170
171
172
173
174
175
176
177
PRO dustem_make_polarization_sed_example,model=model $
,filters=filters $
,outfile=outfile $
,help=help $
,verbose=verbose
;+
; NAME:
; dustem_make_polarization_sed_example
;
; PURPOSE:
; This is an example of how to generate an SED using DustEMWrap.
; It is meant to be an example to follow when writing your own
; programs using DustEMWrap
;
; CATEGORY:
; DustEMWrap, Distributed, High-Level, User Example
;
; CALLING SEQUENCE:
; dustem_make_polarization_sed_example[,model=][,filters=][,/help]
;
; INPUTS:
; None
;
; OPTIONAL INPUT PARAMETERS:
; filters = instrument bands for which SED is created. Default is
; ['IRAS1','IRAS2','IRAS3','IRAS4','PACS3','SPIRE1','SPIRE2','SPIRE3','HFI2','HFI3','HFI4','HFI5','HFI6','LFI1','LFI2','LFI3']
;
; OUTPUTS:
; None
;
; OPTIONAL OUTPUT PARAMETERS:
; outfile = path+filename for .xcat output of the SED. Default is './polarization_sed.xcat'
;
; ACCEPTED KEY-WORDS:
; model = if set, name of the DustEM dust model to use (default='G17_MODELD')
; help = If set print this help
; verbose = If set, subroutines will be more chatty
;
; COMMON BLOCKS:
; None
;
; SIDE EFFECTS:
; None
;
; RESTRICTIONS:
; The DustEM fortran code must be installed
; The DustEMWrap IDL code must be installed
;
; PROCEDURE:
; None
;
; EXAMPLES
; dustem_make_polarization_sed_example
; dustem_make_polarization_sed_example,model='G17_MODELA',outfile='mysed.xcat'
;
; MODIFICATION HISTORY:
; Written by JPB June-2022
; Evolution details on the DustEMWrap gitlab.
; See http://dustemwrap.irap.omp.eu/ for FAQ and help.
;-
IF keyword_set(help) THEN BEGIN
doc_library,'dustem_make_polarization_sed_example'
goto,the_end
ENDIF
dustem_define_la_common
IF keyword_set(model) THEN BEGIN
use_model=strupcase(model)
ENDIF ELSE BEGIN
use_model='G17_MODELD' ;Default is one of the Guillet et al (2017) models since they treat dust polarisation properties
ENDELSE
IF keyword_set(filters) THEN BEGIN
use_filters=strupcase(filters)
ENDIF ELSE BEGIN
use_filters=['IRAS1','IRAS2','IRAS3','IRAS4','PACS3','SPIRE1','SPIRE2','SPIRE3','HFI2','HFI3','HFI4','HFI5','HFI6','LFI1','LFI2','LFI3']
ENDELSE
IF keyword_set(outfile) THEN BEGIN
use_outfile=outfile
ENDIF ELSE BEGIN
use_outfile='./polarization_sed.xcat'
ENDELSE
exists=dustem_test_model_exists(use_model,/pol)
if exists eq 0 then $
message,'Unknown or unpolarized dust model'
use_verbose=0
use_polarization=1
if keyword_set(verbose) then use_verbose=1
dustem_init,mode=use_model,polarization=use_polarization
!dustem_verbose=use_verbose
!dustem_nocatch=1
;=== initialize filters for the SED
Nfilt=n_elements(use_filters)
sed=dustem_initialize_sed(Nfilt)
sed.filter=use_filters
sed.wave=dustem_filter2wav(use_filters)
sed.instru=dustem_filter2instru(use_filters)
;=== initialize IQU and associated errors to avoid problems when checking SED in dustem_set_data.pro
sed[*].StokesI=1.
sed.StokesQ=sed.StokesI/100.
sed.StokesU=sed.StokesI/100.
sed.SigmaII=sed.StokesI/100.
sed.SigmaQQ=sed.StokesI/100.
sed.SigmaUU=sed.StokesI/100.
sed.SigmaIQ=sed.StokesI/100.
sed.SigmaIU=sed.StokesI/100.
sed.SigmaQU=sed.StokesI/100.
;=== initialize the !dustem_data structure that will be used internally
dustem_set_data,m_fit=sed,m_show=sed
;=== Set which model parameters to use for the SED
pd = [ $
'(*!dustem_params).G0', $ ;G0
'(*!dustem_params).grains(0).mdust_o_mh',$ ;PAH0_MC10 mass fraction
'(*!dustem_params).grains(1).mdust_o_mh',$ ;amCBE_0.3333x mass fraction
'(*!dustem_params).grains(2).mdust_o_mh', $ ;aSil2001BE6pctG_0.4x mass fraction
'dustem_plugin_modify_dust_pol_2', $ ;Dust polarization angle
'dustem_plugin_synchrotron_1', $ ;Synchrotron spectra index
'dustem_plugin_synchrotron_2', $ ;Synchrotron amplitude at 10 mm
'dustem_plugin_synchrotron_3', $ ;Synchrotron polarization fraction
'dustem_plugin_synchrotron_4' $ ;Synchrotron polarization angle
]
pd_true=[1.,7.8000E-04,7.8000E-04,7.8000E-04, 48.,3., 1.e-2, 0.3, 45.]
;== SET INITIAL VALUES AND LIMITS
;== AND ACTIVATE ANY PLUGINS
dustem_init_params,use_model,pd,pd_true,pol=use_polarization
;=== compute the SED (I then QU) using the model and specified
;=== true_vals of parameters
dustem_Ised=dustem_compute_sed(pd_true,st=ssti) ;sst is not set on purpose, for dustem_compute_sed to compute the sed using fortran code
toto=dustem_compute_stokes(pd_true,st=sstqu) ;this procedure also allows for the extraction of the spectra
dustem_qsed = toto[0]
dustem_used = toto[1]
sed.stokesI=dustem_Ised
sed.stokesQ=dustem_Qsed
sed.stokesU=dustem_Used
;=== set uncertainties for I,Q,U
sed.sigmaII=abs(sed.StokesI*0.000001)
sed.sigmaQQ=abs(sed.StokesQ*0.000001)
sed.sigmaUU=abs(sed.StokesU*0.000001)
;=== set covariances to 0.
sed.sigmaIQ=0.
sed.sigmaIU=0.
sed.sigmaQU=0.
;==== fill in dependent columns of the SED.
sed=dustem_fill_sed_dependent_columns(sed)
;======== save the SED to a file
write_xcat,sed,use_outfile
message,'Wrote '+use_outfile,/continue
;======== move the file to the Data/EXAMPLE_OBSDATA/ subdirectory
;filename_final=!dustem_wrap_soft_dir+'/Data/EXAMPLE_OBSDATA/my_example_SED_###.xcat'
;str='cp '+use_outfile+' '+filename_final
;message,'Do '+str+' to make change permanent',/continue
the_end:
END