Blame view

src/idl/dustem_fit_intensity_example.pro 16.4 KB
a1839067   Jean-Philippe Bernard   First commit
1
2
3
4
5
PRO dustem_fit_intensity_example,model=model $
                                ,sed_file=sed_file $
                                ,Nitermax=Nitermax $
                                ,postscript=postscript $
                                ,fits_save_and_restore=fits_save_and_restore $
a6bd0451   Annie Hughes   made prettier for...
6
7
8
                                 ,help=help $
                                 ,wait=wait $
                                 ,verbose=verbose
a1839067   Jean-Philippe Bernard   First commit
9

a6bd0451   Annie Hughes   made prettier for...
10
11
;This routine is an example of how to fit an observational SED (Stokes
;I only) with DustEM and DustEMWrap. The objective is to illustrate how to use DustEMWrap
a1839067   Jean-Philippe Bernard   First commit
12
13
14
15
16
17
18
19
; (and not to do science -- the fit obtained by running this example is
; likely to be poor!)
;  
;For this example, the code uses the SED in the file example_SED_1.xcat,
;which is distributed in the Data/EXAMPLE_OBSDATA/ directory 
;
;The example SED has Stokes I photometric data points from SPITZER
;IRAC and MIPS and IRAS. Examples illustrating running DustEMWrap to
a6bd0451   Annie Hughes   made prettier for...
20
21
;fit spectral data, polarisation data and extinction data
;are provided in other _example routines in the src/idl/ directory.
a1839067   Jean-Philippe Bernard   First commit
22
23
24
25
26
  
;+
; NAME:
;    dustem_fit_intensity_example  
; PURPOSE:
a6bd0451   Annie Hughes   made prettier for...
27
28
29
30
31
32
33
34
35
36
37
38
39
;
;  This routine is an example of how to fit an observational SED
;  (Stokes I only) using DustEM and DustEMWrap. The objective is to
;  illustrate how to use DustEMWrap (and not to do science -- the fit
;  obtained by running this example is likely to be poor!)
;  
;  For this example, the code uses the SED in the file example_SED_1.xcat,
;  which is distributed in the Data/EXAMPLE_OBSDATA/ directory 
;
;  The example SED has Stokes I photometric data points from SPITZER
;  IRAC and MIPS and IRAS. Examples illustrating running DustEMWrap to
;  fit spectral data, polarisation data and extinction data
;  are provided in other _example routines in the src/idl/ directory.
a1839067   Jean-Philippe Bernard   First commit
40
41
42
43
44
45
46
47
48
49
50
; CATEGORY:
;    DustEMWrap, Distributed, High-Level, User Example
; CALLING SEQUENCE:
;    dustem_fit_intensity_example[,model=][sed_file=][,postscript=][,Nitermax=][,fits_save_and_restore=][,/help]
; INPUTS:
;    None
; OPTIONAL INPUT PARAMETERS:
;    None
; OUTPUTS:
;    None
; OPTIONAL OUTPUT PARAMETERS:
a6bd0451   Annie Hughes   made prettier for...
51
;    Plots, results structure in binary FITS table format
a1839067   Jean-Philippe Bernard   First commit
52
53
54
55
56
57
58
59
60
61
62
63
64
65
66
67
68
69
; ACCEPTED KEY-WORDS:
;    model = specifies the interstellar dust mixture used by DustEM
;           'MC10' model from Compiegne et al 2010 (default)
;           'DBP90' model from Desert et al 1990
;           'DL01' model from Draine & Li 2001
;           'WD01_RV5p5B' model from Weingartner & Draine 2002 with Rv=5.5
;           'DL07' model from Draine & Li 2007
;           'J13' model from Jones et al 2013, as updated in
;                 Koehler et al 2014
;           'G17_ModelA' model A from Guillet et al (2018). Includes
;                 polarisation. See Tables 2 and 3 of that paper for details.
;           'G17_ModelB' model B from Guillet et al (2018)
;           'G17_ModelC' model C from Guillet et al (2018)
;           'G17_ModelD' model A from Guillet et al (2018)
;    sed_file = string naming the path to text file in .xcat format that
;          describes the observational SED. If not set, the file
;          'Data/SEDs/sample_SED.xcat' is used.  
;    postscript = if set, final plot is saved as postscript in the
a6bd0451   Annie Hughes   made prettier for...
70
;                 current working directory
a1839067   Jean-Philippe Bernard   First commit
71
;    Nitermax = maximum number of fit iterations. Default is 5.
a6bd0451   Annie Hughes   made prettier for...
72
73
74
;    fits_save_and_restore = if set, save the fit results in a binary
;               FITS file. The code then restore this file and plots
;               the results using the saved results information.
a1839067   Jean-Philippe Bernard   First commit
75
;    help      = if set, print this help
a6bd0451   Annie Hughes   made prettier for...
76
77
78
;    wait      = if set, wait this many seconds between each step of
;                the code (for illustration purposes)
;    verbose      = if set, subroutines will run in verbose mode
a1839067   Jean-Philippe Bernard   First commit
79
80
81
82
83
84
85
86
87
88
89
90
; COMMON BLOCKS:
;    None
; SIDE EFFECTS:
;    None
; RESTRICTIONS:
;    The DustEM fortran code must be installed
;    The DustEMWrap IDL code must be installed
; PROCEDURES AND SUBROUTINES USED:
;    *** COMMENT AH --> is this really NONE? ****
; EXAMPLES
;    dustem_fit_intensity_example
;    dustem_fit_intensity_example,Nitermax=1,fits_save_and_restore='/tmp/mysavefile.fits'
16adf9ea   Jean-Philippe Bernard   improved
91
;    dustem_fit_intensity_example,model='DBP90'
a1839067   Jean-Philippe Bernard   First commit
92
93
94
95
96
97
98
99
100
101
102
103
104
105
; MODIFICATION HISTORY:
;    Written by JPB Apr-2011
;    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_fit_intensity_example'
  goto,the_end
END

IF keyword_set(model) THEN BEGIN
  use_model=strupcase(model)
ENDIF ELSE BEGIN
a6bd0451   Annie Hughes   made prettier for...
106
  use_model='MC10'    ;Default is the Compeigne et al (2011) model
a1839067   Jean-Philippe Bernard   First commit
107
108
ENDELSE

a6bd0451   Annie Hughes   made prettier for...
109
110
111
112
use_polarization=0   ; initialize Dustemwrap in no polarization mode 
use_window=2          ; default graphics window number to use for plotting the results
use_verbose=0
if keyword_set(verbose) then use_verbose=1
a1839067   Jean-Philippe Bernard   First commit
113

a6bd0451   Annie Hughes   made prettier for...
114
115
116
117
;=== Set the (model-dependent) parameters  that you want to fit (pd)
;    their initial values (iv)
;    and whether they are bounded (ulimed,llimed,llims,ulims)
;    fixed parameters (fpd) and their values (fiv) are also set here
a1839067   Jean-Philippe Bernard   First commit
118
119
120
121
122
123
124
125
;=== Refer to the DustEM and DustEMWrap userguides for an explanation
;    of the different grain types

CASE use_model OF
      'DBP90':BEGIN
         pd = [ $
             '(*!dustem_params).G0', $      ;G0
             '(*!dustem_params).grains(0).mdust_o_mh',$  ;PAH0 mass fraction
a6bd0451   Annie Hughes   made prettier for...
126
127
             '(*!dustem_params).grains(1).mdust_o_mh',$  ;VSG mass fraction
             '(*!dustem_params).grains(2).mdust_o_mh', $ ;BG mass fraction
a1839067   Jean-Philippe Bernard   First commit
128
129
130
131
132
133
             'dustem_plugin_continuum_2']                ;Intensity of NIR continuum
        iv =   [1.0, 4.3e-4, 4.7e-4,6.4e-3,0.001]
        Npar=n_elements(pd)
        ulimed=replicate(0,Npar)
        llimed=replicate(1,Npar)
        llims=replicate(0.,Npar)
a6bd0451   Annie Hughes   made prettier for...
134
        fpd=[] & fiv=[]
a1839067   Jean-Philippe Bernard   First commit
135
136
137
138
     END
      'DL01':BEGIN
         pd = [ $
             '(*!dustem_params).G0', $      ;G0
a6bd0451   Annie Hughes   made prettier for...
139
             '(*!dustem_params).grains(0).mdust_o_mh',$     ;PAH0 mass fraction
a1839067   Jean-Philippe Bernard   First commit
140
141
142
143
144
145
146
147
148
149
             '(*!dustem_params).grains(1).mdust_o_mh',$     ;PAH1 mass fraction
             '(*!dustem_params).grains(2).mdust_o_mh', $    ;Gra
             '(*!dustem_params).grains(3).mdust_o_mh', $    ;Gra
             '(*!dustem_params).grains(4).mdust_o_mh', $    ;aSil
             'dustem_plugin_continuum_2']                ;Intensity of NIR continuum
        iv =   [1.0,5.4e-4, 5.4e-4,1.8e-4,2.33e-3,8.27e-3,0.001]
        Npar=n_elements(pd)
        ulimed=replicate(0,Npar)
        llimed=replicate(1,Npar)
        llims=replicate(0.,Npar)
a6bd0451   Annie Hughes   made prettier for...
150
        fpd=[] & fiv=[]
a1839067   Jean-Philippe Bernard   First commit
151
152
      END
      'WD01_RV5p5B':BEGIN
a6bd0451   Annie Hughes   made prettier for...
153
154
155
156
157
158
159
160
161
162
163
164
165
166
         pd = [ $
             '(*!dustem_params).G0', $      ;G0
             '(*!dustem_params).grains(0).mdust_o_mh',$     ;PAH0 mass fraction
             '(*!dustem_params).grains(1).mdust_o_mh',$     ;PAH1 mass fraction
             '(*!dustem_params).grains(2).mdust_o_mh', $    ;Gra
             '(*!dustem_params).grains(3).mdust_o_mh', $    ;Gra
             '(*!dustem_params).grains(4).mdust_o_mh', $    ;aSil
             'dustem_plugin_continuum_2']                ;Intensity of NIR continuum
        iv =   [1.0,5.4e-4, 5.4e-4,1.8e-4,1.64e-3,6.99e-3,0.001]
        Npar=n_elements(pd)
        ulimed=replicate(0,Npar)
        llimed=replicate(1,Npar)
        llims=replicate(0.,Npar)
        fpd=[] & fiv=[]
a1839067   Jean-Philippe Bernard   First commit
167
168
169
170
       END
      'DL07':BEGIN
          pd = [ $
              '(*!dustem_params).G0', $      ;G0
a6bd0451   Annie Hughes   made prettier for...
171
              '(*!dustem_params).grains(0).mdust_o_mh',$     ;PAH0 mass fraction
a1839067   Jean-Philippe Bernard   First commit
172
173
174
175
176
              '(*!dustem_params).grains(1).mdust_o_mh',$     ;PAH1 mass fraction
              '(*!dustem_params).grains(2).mdust_o_mh', $    ;Gra
              '(*!dustem_params).grains(3).mdust_o_mh', $    ;Gra
              '(*!dustem_params).grains(4).mdust_o_mh', $    ;aSil
              'dustem_plugin_continuum_2']                ;Intensity of NIR continuum
a6bd0451   Annie Hughes   made prettier for...
177
        iv =   [1.0,4.97e-4, 4.97e-4,1.66e-4,2.21e-3,7.64e-3,0.001]
a1839067   Jean-Philippe Bernard   First commit
178
179
180
181
        Npar=n_elements(pd)
        ulimed=replicate(0,Npar)
        llimed=replicate(1,Npar)
        llims=replicate(0.,Npar)
a6bd0451   Annie Hughes   made prettier for...
182
        fpd=[] & fiv=[]
a1839067   Jean-Philippe Bernard   First commit
183
184
      END
      'MC10':BEGIN
a6bd0451   Annie Hughes   made prettier for...
185
186
; parameters of the model to fit
         pd = [ $
a1839067   Jean-Philippe Bernard   First commit
187
188
189
190
191
192
              '(*!dustem_params).grains(0).mdust_o_mh',$     ;PAH0 mass fraction
              '(*!dustem_params).grains(1).mdust_o_mh',$     ;PAH1 mass fraction
              '(*!dustem_params).grains(2).mdust_o_mh', $    ;amCBEx
              '(*!dustem_params).grains(3).mdust_o_mh', $    ;amCBEx
              '(*!dustem_params).grains(4).mdust_o_mh' $     ;aSilx
              ]
a6bd0451   Annie Hughes   made prettier for...
193
194
; initial values of the model parameters that are to be varied during the fit
         iv =   [ 7.8e-4, 7.8e-4, 1.65e-4, 1.45e-3, 6.7e-3]
a1839067   Jean-Philippe Bernard   First commit
195
        Npar=n_elements(pd)
a6bd0451   Annie Hughes   made prettier for...
196
; flags for whether the parameters are upper/lower bounded
a1839067   Jean-Philippe Bernard   First commit
197
198
        ulimed=replicate(0,Npar)
        llimed=replicate(1,Npar)
a6bd0451   Annie Hughes   made prettier for...
199
; lower bound value
a1839067   Jean-Philippe Bernard   First commit
200
        llims=replicate(0.,Npar)
a6bd0451   Annie Hughes   made prettier for...
201
; parameters of the model that are fixed
a1839067   Jean-Philippe Bernard   First commit
202
        fpd=[ $
a6bd0451   Annie Hughes   made prettier for...
203
          '(*!dustem_params).gas.G0' ,  $   ; ISRF intensity (multiplicative factor x G0)
a1839067   Jean-Philippe Bernard   First commit
204
205
          'dustem_plugin_continuum_2' $     ;intensity of NIR continuum
            ]
a6bd0451   Annie Hughes   made prettier for...
206
207
; values of the fixed model parameters
        fiv=[2.4, 3.e-3]
a1839067   Jean-Philippe Bernard   First commit
208
209
210
211
      END
     'J13':BEGIN
        pd = [ $
             '(*!dustem_params).G0', $      ;G0
a6bd0451   Annie Hughes   made prettier for...
212
213
214
215
             '(*!dustem_params).grains(0).mdust_o_mh',$     ;CM20 -- power law size distribution
             '(*!dustem_params).grains(1).mdust_o_mh',$     ;CM20 -- logN size distribution
             '(*!dustem_params).grains(2).mdust_o_mh', $    ;aPyM5
             '(*!dustem_params).grains(3).mdust_o_mh', $    ;aOlM5
a1839067   Jean-Philippe Bernard   First commit
216
             'dustem_plugin_continuum_2']                ;Intensity of NIR continuum
a6bd0451   Annie Hughes   made prettier for...
217
        iv =   [1.0, 1.7e-3, 6.3e-4,2.55e-3,2.55e-3,0.001]
a1839067   Jean-Philippe Bernard   First commit
218
219
220
221
        Npar=n_elements(pd)
        ulimed=replicate(0,Npar)
        llimed=replicate(1,Npar)
        llims=replicate(0.,Npar)
a6bd0451   Annie Hughes   made prettier for...
222
223
        fpd=[] & fiv=[]
     END
a1839067   Jean-Philippe Bernard   First commit
224
      'G17_MODELA':BEGIN
a1839067   Jean-Philippe Bernard   First commit
225
226
227
         pd = [ $
             '(*!dustem_params).G0', $      ;G0
             '(*!dustem_params).grains(0).mdust_o_mh',$  ;PAH0 mass fraction
a6bd0451   Annie Hughes   made prettier for...
228
229
             '(*!dustem_params).grains(1).mdust_o_mh',$     ;amCB_0.3333x
             '(*!dustem_params).grains(2).mdust_o_mh', $    ;aSil2001_0.3333_p20B
a1839067   Jean-Philippe Bernard   First commit
230
             'dustem_plugin_continuum_2']                ;Intensity of NIR continuum
a6bd0451   Annie Hughes   made prettier for...
231
        iv =   [1.0, 7.094483e-4, 1.319896e-3, 5.588265e-3, 0.001]
a1839067   Jean-Philippe Bernard   First commit
232
233
234
235
        Npar=n_elements(pd)
        ulimed=replicate(0,Npar)
        llimed=replicate(1,Npar)
        llims=replicate(0.,Npar)
a6bd0451   Annie Hughes   made prettier for...
236
237
        use_polarization=1
        fpd=[] & fiv=[]
a1839067   Jean-Philippe Bernard   First commit
238
239
240
      END
      'G17_MODELB':BEGIN
         pd = [ $
a6bd0451   Annie Hughes   made prettier for...
241
            '(*!dustem_params).G0', $      ;G0
a1839067   Jean-Philippe Bernard   First commit
242
             '(*!dustem_params).grains(0).mdust_o_mh',$  ;PAH0 mass fraction
a6bd0451   Annie Hughes   made prettier for...
243
244
             '(*!dustem_params).grains(1).mdust_o_mh',$     ;amCB_0.3333x
             '(*!dustem_params).grains(2).mdust_o_mh', $    ;aSil_0.3333_p20B
a1839067   Jean-Philippe Bernard   First commit
245
             'dustem_plugin_continuum_2']                ;Intensity of NIR continuum
a6bd0451   Annie Hughes   made prettier for...
246
        iv =   [1.0, 6.480655e-4, 1.327672e-3, 5.267745e-3, 0.001]
a1839067   Jean-Philippe Bernard   First commit
247
248
249
250
251
        Npar=n_elements(pd)
        ulimed=replicate(0,Npar)
        llimed=replicate(1,Npar)
        llims=replicate(0.,Npar)
        use_polarization=1
a6bd0451   Annie Hughes   made prettier for...
252
        fpd=[] & fiv=[]
a1839067   Jean-Philippe Bernard   First commit
253
254
     END
      'G17_MODELC':BEGIN
a1839067   Jean-Philippe Bernard   First commit
255
256
257
         pd = [ $
             '(*!dustem_params).G0', $      ;G0
             '(*!dustem_params).grains(0).mdust_o_mh',$  ;PAH0 mass fraction
a6bd0451   Annie Hughes   made prettier for...
258
259
             '(*!dustem_params).grains(1).mdust_o_mh',$     ;amCBE_0.3333x
             '(*!dustem_params).grains(2).mdust_o_mh', $    ;aSil2001_0.3333_p20B
a1839067   Jean-Philippe Bernard   First commit
260
             'dustem_plugin_continuum_2']                ;Intensity of NIR continuum
a6bd0451   Annie Hughes   made prettier for...
261
        iv =   [1.0, 6.979060e-4, 1.555919e-3, 4.366849e-3, 0.001]
a1839067   Jean-Philippe Bernard   First commit
262
263
264
265
266
        Npar=n_elements(pd)
        ulimed=replicate(0,Npar)
        llimed=replicate(1,Npar)
        llims=replicate(0.,Npar)
        use_polarization=1
a6bd0451   Annie Hughes   made prettier for...
267
        fpd=[] & fiv=[]
a1839067   Jean-Philippe Bernard   First commit
268
269
270
      END
      'G17_MODELD':BEGIN
         pd = [ $
a6bd0451   Annie Hughes   made prettier for...
271
            '(*!dustem_params).G0', $      ;G0
a1839067   Jean-Philippe Bernard   First commit
272
             '(*!dustem_params).grains(0).mdust_o_mh',$  ;PAH0 mass fraction
a6bd0451   Annie Hughes   made prettier for...
273
274
             '(*!dustem_params).grains(1).mdust_o_mh',$     ;amCBE_0.3333x
             '(*!dustem_params).grains(2).mdust_o_mh', $    ;aSil2001BE6pctG_0.4x
a1839067   Jean-Philippe Bernard   First commit
275
             'dustem_plugin_continuum_2']                ;Intensity of NIR continuum
a6bd0451   Annie Hughes   made prettier for...
276
        iv =   [1.0, 7.235580e-4, 1.222725e-3, 5.324418e-3, 0.001]
a1839067   Jean-Philippe Bernard   First commit
277
278
279
280
281
        Npar=n_elements(pd)
        ulimed=replicate(0,Npar)
        llimed=replicate(1,Npar)
        llims=replicate(0.,Npar)
        use_polarization=1
a6bd0451   Annie Hughes   made prettier for...
282
        fpd=[] & fiv=[]
a1839067   Jean-Philippe Bernard   First commit
283
284
285
     END
      'ELSE':BEGIN
         message,'model '+model+' unknown',/continue
a6bd0451   Annie Hughes   made prettier for...
286
         message,'Known models are MC10,DBP90,DL01,WD01_RV5p5B,DL07,J13,G17_MODELA,G17_MODELB,G17_MODELC,G17_MODELD',/continue
a1839067   Jean-Philippe Bernard   First commit
287
288
289
         stop
      END
ENDCASE
a6bd0451   Annie Hughes   made prettier for...
290
291
292
293
294
Nfix=n_elements(fpd)
if n_elements(fiv) ne Nfix then begin
   message,'Number of fixed parameters (fpd) does not equal number of initial values of fixed parameters (fiv)',/info
   stop
end
a1839067   Jean-Philippe Bernard   First commit
295

a6bd0451   Annie Hughes   made prettier for...
296
297
298
299
if keyword_set(wait) then begin
   message,'Finished setting dust model: '+model,/info
   wait,wait
end
a1839067   Jean-Philippe Bernard   First commit
300
301
302


;== INITIALISE DUSTEM
a6bd0451   Annie Hughes   made prettier for...
303
dustem_init,model=use_model,polarization=use_polarization
a1839067   Jean-Philippe Bernard   First commit
304
!dustem_nocatch=1
a6bd0451   Annie Hughes   made prettier for...
305
!dustem_verbose=use_verbose
a1839067   Jean-Philippe Bernard   First commit
306
307
308
!dustem_show_plot=1


a6bd0451   Annie Hughes   made prettier for...
309
;=== READ EXAMPLE SED DATA
a1839067   Jean-Philippe Bernard   First commit
310
311
312
313
314
dir=!dustem_wrap_soft_dir+'/Data/EXAMPLE_OBSDATA/'
file=dir+'example_SED_1.xcat'
IF keyword_set(sed_file) THEN file=sed_file
sed=read_xcat(file,/silent)

a6bd0451   Annie Hughes   made prettier for...
315
316
317
318
319
320
321
if keyword_set(wait) then begin
   message,'Finished reading SED data: '+file,/info
   wait,wait
end

;;=== ADJUST THE UNCERTAINTIES (FOR ILLUSTRATION)
ind=where(sed.sigmaII LT (0.2*sed.StokesI)^2,count)
a1839067   Jean-Philippe Bernard   First commit
322
IF count NE 0 THEN sed[ind].sigmaII=(0.2*sed[ind].StokesI)^2
a1839067   Jean-Philippe Bernard   First commit
323
324

;== SET THE OBSERVATIONAL STRUCTURE
a6bd0451   Annie Hughes   made prettier for...
325
326
327
328
;== sed is passed twice -- the first occurrence is the SED that you
;== wish to fit, the second occurrence is the
;== SED that you wish to visualise / plot. These are often, but not
;== always!, the same SED.
80783331   Jean-Philippe Bernard   changed
329
dustem_set_data,sed,sed
a1839067   Jean-Philippe Bernard   First commit
330
331
332
333
334
335
336

;== SET INITIAL VALUES AND LIMITS OF THE PARAMETERS THAT WILL BE
;== ADJUSTED DURING THE FIT
dustem_init_parinfo,pd,iv,up_limited=ulimed,lo_limited=llimed,up_limits=ulims,lo_limits=llims

;== INITIALIZE ANY PLUGINS
dustem_init_plugins,pd,fpd 
a6bd0451   Annie Hughes   made prettier for...
337

a1839067   Jean-Philippe Bernard   First commit
338
;== INITIALIZE ANY FIXED PARAMETERS FOR PLUGINS
a6bd0451   Annie Hughes   made prettier for...
339
340
341
342
343
344
if Nfix gt 0 then dustem_init_fixed_params,fpd,fiv

if keyword_set(wait) then begin
   message,'Finished initializing DustEmWrap, including plugins and fixed parameters',/info
   wait,wait
end
a1839067   Jean-Philippe Bernard   First commit
345

a6bd0451   Annie Hughes   made prettier for...
346
;=== INFORMATION TO RUN THE FIT
a1839067   Jean-Philippe Bernard   First commit
347
348
349
350
tol=1.e-16     ;fit tolerence
use_Nitermax=5 ;maximum number of iterations. 
IF keyword_set(Nitermax) THEN use_Nitermax=Nitermax

a6bd0451   Annie Hughes   made prettier for...
351
352
353
354
355
356
;=== INFORMATION TO MAKE THE PLOT
yr=[1.00e-4,1.00E2] ; y-axis limits
xr=[1.00E0,6.00e4] ; x-axis limits
tit='SED -- Example 1' ; plot title
ytit=textoidl('I_\nu (MJy/sr) for N_H=10^{20} H/cm^2') ; y-axis title
xtit=textoidl('\lambda (\mum)') ; x-axis title
a1839067   Jean-Philippe Bernard   First commit
357

a6bd0451   Annie Hughes   made prettier for...
358
;===  RUN THE FIT
a1839067   Jean-Philippe Bernard   First commit
359
360
361
362
363
364
365
t1=systime(0,/sec)
res=dustem_mpfit_data(tol=tol,Nitermax=use_Nitermax,gtol=gtol $
                      ,/xlog,/ylog,xr=xr,yr=yr,xtit=xtit,ytit=ytit,title=tit $
                      ,legend_xpos=legend_xpos,legend_ypos=legend_ypos $
                      ,errors=errors,chi2=chi2,rchi2=rchi2)
t2=systime(0,/sec)

a6bd0451   Annie Hughes   made prettier for...
366
367
368
369
370
371
if keyword_set(wait) then begin
   message,'Finished running DustEmWrap, using Niters: '+strtrim(string(use_Nitermax),2),/info
   message,'Time taken [sec]: '+sigfig(t2-t1,2,/sci),/info
   wait,wait
end
;=== MAKE THE FINAL PLOT
a1839067   Jean-Philippe Bernard   First commit
372
373
374
375
376
377
IF keyword_set(postscript) THEN BEGIN
    dir_ps='./'
    set_plot,'PS'
    ps_file=dir_ps+postscript
    device,filename=ps_file,/color
ENDIF
a6bd0451   Annie Hughes   made prettier for...
378
dustemwrap_plot,*(*!dustem_fit).CURRENT_PARAM_VALUES,dummy,xr=xr,/xstyle,yr=yr,/ysty,/ylog,/xlog,title=tit+' (Final fit)'
a1839067   Jean-Philippe Bernard   First commit
379
380
381
382
383
384
IF keyword_set(postscript) THEN BEGIN
  set_plot,'X'
  device,/close
  message,'Wrote '+ps_file,/info
ENDIF

a6bd0451   Annie Hughes   made prettier for...
385
386
387
388
389
if keyword_set(wait) then begin
   message,'Made the plot of the final results',/info
   wait,wait
end

a1839067   Jean-Philippe Bernard   First commit
390
IF keyword_set(fits_save_and_restore) THEN BEGIN
a6bd0451   Annie Hughes   made prettier for...
391
392
393
394
   message,'Writing out results structure: '+fits_save_and_restore,/info
   dustem_write_fits_table,filename=fits_save_and_restore,help=help
;=== at this point, you could erase all dustem system variables, or exit idl... all the
; information needed has been saved in the FITS table
a1839067   Jean-Philippe Bernard   First commit
395
396
397
  dustem_read_fits_table,filename=fits_save_and_restore,dustem_spectra_st=dustem_spectra_st
  ;==== plot result taken from the saved fits table
  res=*(*!dustem_fit).CURRENT_PARAM_VALUES
a6bd0451   Annie Hughes   made prettier for...
398
399
400
401
402
  dustemwrap_plot,res,dustem_spectra_st,xr=xr,/xstyle,yr=yr,/ysty,/ylog,/xlog,title=tit+' (From Saved FITS file)'
  if keyword_set(wait) then begin
     message,'Saved the results as FITS in the file: '+fits_save_and_restore+', and made a plot using the data in this file',/info
     wait,wait
  end
a1839067   Jean-Philippe Bernard   First commit
403
404
ENDIF

a6bd0451   Annie Hughes   made prettier for...
405
message,'Finished dustem_fit_intensity_example',/info
a1839067   Jean-Philippe Bernard   First commit
406
407
408
409

the_end:

END