dustem_set_data.pro
14.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
348
349
350
351
352
353
354
355
356
357
358
359
360
361
362
363
364
365
366
367
368
369
370
371
372
373
PRO dustem_set_data,m_fit,m_show,x_fit,x_show,rchi2_weight=rchi2_weight,f_HI=f_HI,help=help
;+
; NAME:
; dustem_set_data
;
; PURPOSE:
; Initializes the dustem data structure (*!dustem_data) used in the fitting
; Initializes the !dustem_show structure which is either a copy or a more detailed version of !dustem_data
;
; 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
; x_fit : dustemwrap structure containing extinction data to be fitted
; m_show : dustemwrap structure containing emission data to be displayed
; 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
;
; 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,spec,x_fit,spec,x_show)
; #in this example !dustem_data and !dustem_show are 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_check_data is introduced and the data structure format of the latest release is respected.
; dustem_set_data is turned into a procedure.
; 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
;#FITTING STRUCTURE
;To me this is unnecessary... as dustem_check_data is always needed as it shapes the output structures (that would otherwise not exist)
;in order to prepare them for the mpfit run.
;the *!dustem_data tags are always cleaned so data will never exist if this keyword is used.
;commenting the use of nocheck_data for now
;IF not keyword_set(nocheck_data) THEN BEGIN
dustem_check_data,m_fit,x_fit,sed,ext,polext,polsed,polfrac,psi_em,qsed,used,qext,uext,psi_ext,fpolext
;ENDIF
;outputs whatever the input structure has, into organized seperate structures formatted for the fitting process.
;FITTING STRUCTURE INITIALIZATION
IF isa(sed) THEN (*!dustem_data).sed = ptr_new(sed)
IF isa(ext) THEN (*!dustem_data).ext = ptr_new(ext)
;####################
;OLD VARIABLE. KEPT IT.
; 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
;####################
if !run_pol then begin
if keyword_set(m_fit) then begin
teststks = isa(qsed) and isa(used)
;stop
IF ~teststks then message, 'Stokes parameters are not set correctly. Please check the input structure. Aborting...'
;☆EMISSION
message, 'Polarization component(s) to fit: QSED and USED' ,/continue
;☆Mandatory presence in polar mode
IF isa(qsed) THEN BEGIN
(*!dustem_data).qsed = ptr_new(qsed)
if keyword_set(f_HI) then if f_HI gt 0 then $
(*(*!dustem_data).qsed).values *= f_HI
ENDIF
;☆Mandatory presence in polar mode
if isa(used) then BEGIN
(*!dustem_data).used = ptr_new(used)
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)
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)
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)
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
;Two tests that determine what the user will eventually fit
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 - allowed for the default fitting of both because I do not know if the same formalism should be applied to the extinction data because it also implies deviding extinction data into two components.
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
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
;Setting of the structure that will be displayed
;IF not keyword_set(nocheck_data) THEN BEGIN
dustem_check_data,m_show,x_show,sed,ext,polext,polsed,polfrac,psi_em,qsed,used,qext,uext,psi_ext,fpolext
;ENDIF
IF isa(sed) THEN (*!dustem_show).sed = ptr_new(sed)
IF isa(ext) THEN (*!dustem_show).ext = ptr_new(ext)
;####################
;OLD VARIABLE. KEPT IT.
; 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 teststks and isa(polsed) and ~(isa(polfrac)) and ~(isa(psi_em)) then nowshowing = 'QSED, USED and POLSED'
if teststks and isa(polsed) and ~(isa(polfrac)) and isa(psi_em) then nowshowing = 'QSED, USED, POLSED and PSI_EM'
if teststks and ~(isa(polsed)) and isa(polfrac) and ~(isa(psi_em)) then nowshowing = 'QSED, USED and POLFRAC'
if teststks and ~(isa(polsed)) and isa(polfrac) and isa(psi_em) then nowshowing = 'QSED, USED, POLFRAC and PSI_EM'
if teststks and ~(isa(polsed)) and ~(isa(polfrac)) and ~(isa(psi_em)) then nowshowing = 'QSED and USED'
if teststks and isa(polsed) and isa(polfrac) and ~(isa(psi_em)) then nowshowing = 'QSED, USED, POLSED and POLFRAC'
if teststks and isa(polsed) and isa(polfrac) and isa(psi_em) then nowshowing = 'QSED, USED, POLFRAC, POLSED and PSI_EM'
if teststks and ~(isa(polsed)) and ~(isa(polfrac)) and isa(psi_em) then nowshowing = 'QSED, USED and PSI_EM'
message, 'Polarization component(s) to show: '+nowshowing ,/continue
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...'
if teststks and isa(polext) and ~(isa(polextfrac)) and ~(isa(psi_ext)) then nowshowing = 'QEXT, UEXT and POLEXT'
if teststks and isa(polext) and ~(isa(polextfrac)) and isa(psi_ext) then nowshowing = 'QEXT, UEXT, POLEXT and PSI_EXT'
if teststks and ~(isa(polext)) and isa(polextfrac) and ~(isa(psi_ext)) then nowshowing = 'QEXT, UEXT and POLEXTFRAC'
if teststks and ~(isa(polext)) and isa(polextfrac) and isa(psi_ext) then nowshowing = 'QEXT, UEXT, POLEXTFRAC and PSI_EXT'
if teststks and ~(isa(polext)) and ~(isa(polextfrac)) and ~(isa(psi_ext)) then nowshowing = 'QEXT and UEXT'
if teststks and isa(polext) and isa(polextfrac) and ~(isa(psi_ext)) then nowshowing = 'QEXT, UEXT, POLEXT and POLEXTFRAC'
if teststks and isa(polext) and isa(polextfrac) and isa(psi_ext) then nowshowing = 'QEXT, UEXT, POLEXTFRAC, POLEXT and PSI_EXT'
if teststks and ~(isa(polext)) and ~(isa(polextfrac)) and isa(psi_ext) then nowshowing = 'QEXT, UEXT and PSI_EXT'
message, 'Polarization component(s) to show: '+nowshowing ,/continue
;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:
;
; return,0.
END