dustem_set_data.pro
12.2 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
178
179
180
181
182
183
184
185
186
187
188
189
190
191
192
193
194
195
196
197
198
199
200
201
202
203
204
205
206
207
208
209
210
211
212
213
214
215
216
217
218
219
220
221
222
223
224
225
226
227
228
229
230
231
232
233
234
235
236
237
238
239
240
241
242
243
244
245
246
247
248
249
250
251
252
253
254
255
256
257
258
259
260
261
262
263
264
265
266
267
268
269
270
271
272
273
274
275
276
277
278
279
280
281
282
283
284
285
286
287
288
289
290
291
292
293
294
295
296
297
298
299
300
301
302
303
304
305
306
307
308
309
310
311
312
313
314
315
316
317
318
319
320
321
322
323
324
325
326
327
328
329
330
331
332
333
334
335
336
337
338
339
340
341
342
343
344
345
346
347
PRO dustem_set_data,m_fit=m_fit $
,m_show=m_show $
,x_fit=x_fit $
,x_show=x_show $
,rchi2_weight=rchi2_weight $
,f_HI=f_HI $
,help=help
;+
; NAME:
; dustem_set_data
;
; PURPOSE:
; Set the *!dustem data (and *!dustem_show) datasets in emission and/or extinction from input data structures
; *!dustem_data is important for the fitting of data and its general handling by DustEMWrap.
; CATEGORY:
; DustEMWrap, Mid-Level, Distributed, Initialization
;
; CALLING SEQUENCE:
; dustem_set_data,[,m_fit=,x_fit=][,m_show=,x_show=],[,rchi2_weight=][,f_HI=][,/help]
;
; INPUTS:
; None
;
; OPTIONAL INPUT PARAMETERS:
; m_fit : DustEMWrap structure contaning emission data to be fitted
; m_show : DustEMWrap structure containing emission data to be displayed
; x_fit : DusteEMWrap structure containing extinction data to be fitted
; x_show : DustEMWrap structure containing extinctiond data to be displayed
; rchi2_weight : DustEMWrap structure containing the wieghts for chi2 calculation
; f_HI : multiplicative factor for the SED values (extinction and emission)
;
; OUTPUTS:
; None/ALL OF THE ABOVE?
;
; OPTIONAL OUTPUT PARAMETERS:
; None
;
; ACCEPTED KEY-WORDS:
; help = If set print this help
;
; COMMON BLOCKS:
; None
;
; SIDE EFFECTS:
; initializes !dustem_data and !dustem_show
; wavelengths are forced to be the central filter wavelength for photometric data
;
; RESTRICTIONS:
; The DustEMWrap IDL code must be installed
;
; PROCEDURE:
; None
;
; EXAMPLES
; dustem_init
; dir=!dustem_wrap_soft_dir+'/Data/EXAMPLE_OBSDATA/'
; file=dir+'example_SED_1.xcat'
; spec=read_xcat(file,/silent)
; dustem_set_data,m_fit=spec,m_show=spec
; nb: in this example *!dustem_data and *!dustem_show end up identical
;
; MODIFICATION HISTORY:
; Written by J.-Ph. Bernard 2007
; V. Guillet (2012) : dustem_set_data is turned into a function which can be used
; indifferently for sed, ext, polext
; I. Choubani (2022): dustem_set_data is turned into a procedure for emission/extinction with/without polarization.
; 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_set_data'
goto,the_end
ENDIF
;Cleaninig fitting/showing structure(s)
for i=0L,n_tags(*!dustem_data)-1 do begin
(*!dustem_data).(i) = ptr_new()
(*!dustem_show).(i) = ptr_new()
endfor
;Outputs have the format required by !dustem_data and the subsequent fitting process.
dustem_check_data, m_fit,x_fit,sed=sed,ext=ext,polext=polext, $
polsed=polsed,polfrac=polfrac,psi_em=psi_em, $
qsed=qsed,used=used,qext=qext,uext=uext, $
psi_ext=psi_ext,fpolext=fpolext
;FITTING STRUCTURE INITIALIZATION
IF isa(sed) THEN (*!dustem_data).sed = ptr_new(sed)
IF isa(ext) THEN (*!dustem_data).ext = ptr_new(ext)
;####################
;RELIC. For Annie to judge if this is necessary
; If f_HI is specified, multiply the data by f_HI
If keyword_set(f_HI) and ptr_valid((*!dustem_data).sed) then begin
if f_HI gt 0 then (*(*!dustem_data).sed).values *= f_HI else f_HI = 1.
endif else f_HI = 1.
if keyword_set(f_HI) and ptr_valid((*!dustem_data).ext) then begin
if f_HI gt 0 then (*(*!dustem_data).ext).values *= f_HI else f_HI = 1.
endif else f_HI = 1.
defsysv,'!dustem_f_HI',f_HI
;####################
;Problem here is that the user can activate !run_pol without wanting to be in polarization mode.
;reinitialze !run_pol here? IC - January 2023
; March 2023 -- Annie commented following lines since they prevent using a model to
; generate SEDs.
;; If !run_pol THEN BEGIN
;; If tag_exist(*!dustem_show,'QSED') THEN BEGIN
;; IF ~isa((*!dustem_data).qsed) THEN !run_pol=0
;; ENDIF ELSE !run_pol=0
;; ENDIF
if !run_pol then begin
if keyword_set(m_fit) then begin
teststks = isa(qsed) and isa(used)
IF ~teststks then message, 'Stokes parameters are not set correctly. Please check the input structure. Aborting...'
message, 'Polarization component(s) to fit: QSED and USED' ,/continue
IF isa(qsed) THEN BEGIN
(*!dustem_data).qsed = ptr_new(qsed)
;RELIC
if keyword_set(f_HI) then if f_HI gt 0 then $
(*(*!dustem_data).qsed).values *= f_HI
ENDIF
if isa(used) then BEGIN
(*!dustem_data).used = ptr_new(used)
;RELIC
if keyword_set(f_HI) then if f_HI gt 0 then $
(*(*!dustem_data).used).values *= f_HI
ENDIF
if isa(polsed) then BEGIN
(*!dustem_data).polsed = ptr_new(polsed)
;RELIC
if keyword_set(f_HI) then if f_HI gt 0 then $
(*(*!dustem_data).polsed).values *= f_HI
ENDIF
if isa(polfrac) then BEGIN
(*!dustem_data).polfrac = ptr_new(polfrac)
;RELIC
if keyword_set(f_HI) then if f_HI gt 0 then $
(*(*!dustem_data).polfrac).values *= f_HI
ENDIF
if isa(psi_em) then BEGIN
(*!dustem_data).psi_em = ptr_new(psi_em)
;RELIC
if keyword_set(f_HI) then if f_HI gt 0 then $
(*(*!dustem_data).psi_em).values *= f_HI
ENDIF
endif
if keyword_set(x_fit) then begin
testexstks = isa(qext) and isa(uext)
IF ~testexstks then message, 'Extinction stokes parameters/Polarization data is not set correctly. Please check the input structure. Aborting...'
;EXTINCTION -
message, 'Extinction polarization component(s) to fit: QEXT and UEXT',/continue
if isa(polext) then BEGIN
(*!dustem_data).polext = ptr_new(polext)
if keyword_set(f_HI) then if f_HI gt 0 then $
(*(*!dustem_data).polext).values *= f_HI
ENDIF
IF isa(qext) THEN BEGIN
(*!dustem_data).qext = ptr_new(qext)
if keyword_set(f_HI) then if f_HI gt 0 then $
(*(*!dustem_data).qext).values *= f_HI
ENDIF
if isa(uext) then BEGIN
(*!dustem_data).uext = ptr_new(uext)
if keyword_set(f_HI) then if f_HI gt 0 then $
(*(*!dustem_data).uext).values *= f_HI
ENDIF
if isa(psi_ext) then BEGIN
(*!dustem_data).psi_ext = ptr_new(psi_ext)
if keyword_set(f_HI) then if f_HI gt 0 then $
(*(*!dustem_data).psi_ext).values *= f_HI
ENDIF
;added for completeness
if isa(fpolext) then BEGIN
(*!dustem_data).fpolext = ptr_new(fpolext)
if keyword_set(f_HI) then if f_HI gt 0 then $
(*(*!dustem_data).fpolext).values *= f_HI
ENDIF
endif
endif
;SHOWING STRUCTURE / IT IS NEEDED TO BE ABLE TO LOCATE THE HIDDEN DATA POINTS WHEN COMPARING WITH !dustem_data.
;when creating the HUGE system variable JP wants to create. We'll be able to omit it from DusteEMWrap for good.
if keyword_set(m_fit) and not keyword_set(m_show) then m_show=m_fit
if keyword_set(x_fit) and not keyword_set(x_show) then x_show = x_fit
dustem_check_data, m_show,x_show,sed=sed,ext=ext,polext=polext, $
polsed=polsed,polfrac=polfrac,psi_em=psi_em, $
qsed=qsed,used=used,qext=qext,uext=uext, $
psi_ext=psi_ext,fpolext=fpolext
IF isa(sed) THEN (*!dustem_show).sed = ptr_new(sed)
IF isa(ext) THEN (*!dustem_show).ext = ptr_new(ext)
;####################
;RELIC
; If f_HI is specified, multiply the data by f_HI
If keyword_set(f_HI) and ptr_valid((*!dustem_show).sed) then begin
if f_HI gt 0 then (*(*!dustem_show).sed).values *= f_HI else f_HI = 1.
endif else f_HI = 1.
if keyword_set(f_HI) and ptr_valid((*!dustem_show).ext) then begin
if f_HI gt 0 then (*(*!dustem_show).ext).values *= f_HI else f_HI = 1.
endif else f_HI = 1.
defsysv,'!dustem_f_HI',f_HI
;####################
if !run_pol then begin
if isa(m_show) then begin ; because sometimes this keyword isn't set.
teststks = isa(qsed) and isa(used)
IF (~teststks and ~isa(polsed)) or (~teststks and ~isa(polfrac)) then message, 'Stokes parameters/Polarization showing data is not set correctly. Please check the input structure. Aborting...'
If isa(polsed) then BEGIN ;Apply changes to !dustem_data here.
(*!dustem_show).polsed = ptr_new(polsed)
if keyword_set(f_HI) then if f_HI gt 0 then $
(*(*!dustem_show).polsed).values *= f_HI
ENDIF
If isa(polfrac) then BEGIN
(*!dustem_show).polfrac = ptr_new(polfrac)
if keyword_set(f_HI) then if f_HI gt 0 then $
(*(*!dustem_show).polfrac).values *= f_HI
ENDIF
IF isa(qsed) THEN BEGIN
(*!dustem_show).qsed = ptr_new(qsed)
if keyword_set(f_HI) then if f_HI gt 0 then $
(*(*!dustem_show).qsed).values *= f_HI
ENDIF
if isa(used) then BEGIN
(*!dustem_show).used = ptr_new(used)
if keyword_set(f_HI) then if f_HI gt 0 then $
(*(*!dustem_show).used).values *= f_HI
ENDIF
if isa(psi_em) then BEGIN
(*!dustem_show).psi_em = ptr_new(psi_em)
if keyword_set(f_HI) then if f_HI gt 0 then $
(*(*!dustem_show).psi_em).values *= f_HI
ENDIF
endif
if isa(x_show) then begin
teststks = isa(qext) and isa(uext)
IF (~teststks) then message, 'Extinction stokes parameters/Polarization data is not set correctly. Please check the input structure. Aborting...'
;EXTINCTION
if isa(polext) then BEGIN
(*!dustem_show).polext = ptr_new(polext)
if keyword_set(f_HI) then if f_HI gt 0 then $
(*(*!dustem_show).polext).values *= f_HI
ENDIF
IF isa(qext) THEN BEGIN
(*!dustem_show).qext = ptr_new(qext)
if keyword_set(f_HI) then if f_HI gt 0 then $
(*(*!dustem_show).qext).values *= f_HI
ENDIF
if isa(uext) then BEGIN
(*!dustem_show).uext = ptr_new(uext)
if keyword_set(f_HI) then if f_HI gt 0 then $
(*(*!dustem_show).uext).values *= f_HI
ENDIF
if isa(psi_ext) then BEGIN
(*!dustem_show).psi_ext = ptr_new(psi_ext)
if keyword_set(f_HI) then if f_HI gt 0 then $
(*(*!dustem_show).psi_ext).values *= f_HI
ENDIF
;Added for completeness
if isa(fpolext) then BEGIN
(*!dustem_show).fpolext = ptr_new(fpolext)
if keyword_set(f_HI) then if f_HI gt 0 then $
(*(*!dustem_show).fpolext).values *= f_HI
ENDIF
endif
endif
IF keyword_set(rchi2_weight) THEN BEGIN
(*!fit_rchi2_weight).sed = rchi2_weight.sed
(*!fit_rchi2_weight).ext = rchi2_weight.ext
IF !run_pol THEN BEGIN
IF isa(qsed) THEN (*!fit_rchi2_weight).qsed = rchi2_weight.qsed
IF isa(used) THEN (*!fit_rchi2_weight).used = rchi2_weight.used
IF isa(qext) THEN (*!fit_rchi2_weight).qext = rchi2_weight.qext
IF isa(uext) THEN (*!fit_rchi2_weight).uext = rchi2_weight.uext
ENDIF
ENDIF ELSE BEGIN
(*!fit_rchi2_weight).sed=1.
(*!fit_rchi2_weight).ext=1.
If !run_pol THEN BEGIN
IF isa(qsed) THEN (*!fit_rchi2_weight).qsed=1.
IF isa(used) THEN (*!fit_rchi2_weight).used=1.
IF isa(qext) THEN (*!fit_rchi2_weight).qext=1.
IF isa(uext) THEN (*!fit_rchi2_weight).uext=1.
ENDIF
ENDELSE
; ;Setting the data and wavelength tags in '!dustem_fit'
; ;ask about it because !dustem_fit's tag was initially empty...s
; ((*!dustem_fit).data)=ptr_new(!dustem_data)
IF not ptr_valid((*!dustem_data).sed) AND not ptr_valid((*!dustem_data).ext) THEN BEGIN
message,'dustemwrap system variable was not set properly.',/continue
message,'Check the existence of at least some non undefined values in StokesI and sigmaII.',/continue
stop
ENDIF
the_end:
END