Blame view

src/idl/dustem_mpfit_run.pro 13.1 KB
40597848   Ilyes Choubani   Plotting of numbe...
1
FUNCTION dustem_mpfit_run,x,p_min,niter=niter_completed,_EXTRA=_extra
427f1205   Jean-Michel Glorian   version 4.2 merged
2
3
4

;+
; NAME:
36e8b879   Jean-Michel Glorian   update from deborah
5
;    dustem_mpfit_run_vg
427f1205   Jean-Michel Glorian   version 4.2 merged
6
7
8
9
10
; PURPOSE:
;    function used to fit dustem model with dustem_mpfit_sed
; CATEGORY:
;    Dustem
; CALLING SEQUENCE:
36e8b879   Jean-Michel Glorian   update from deborah
11
;    res=dustem_mpfit_run_vg(x,params[,_EXTRA=extra][,cc_sed=cc_sed][,/help])
427f1205   Jean-Michel Glorian   version 4.2 merged
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
; INPUTS:
;    x         = necessary for interface but unused
;    params    = Dustem Model parameters normalized to initial values
; OPTIONAL INPUT PARAMETERS:
;    cc_sed    = SED color correction factors
;    _extra    = extra parameters for the plot routine
; OUTPUTS:
;    None
; OPTIONAL OUTPUT PARAMETERS:
;    None
; ACCEPTED KEY-WORDS:
;    help      = If set, print this help
; COMMON BLOCKS:
;    None
; SIDE EFFECTS:
;    SED and model are plotted if !dustem_show_plot
; RESTRICTIONS:
;    The dustem idl wrapper must be installed
; PROCEDURE:
;    None
; EXAMPLES
36e8b879   Jean-Michel Glorian   update from deborah
33
;    res=dustem_mpfit_run_vg(x,params[,_EXTRA=extra][,cc_sed=cc_sed][,/help])
4750086c   Ilyes Choubani   nouvelle philosph...
34
;  ODIFICATION HISTORY:
427f1205   Jean-Michel Glorian   version 4.2 merged
35
36
37
38
39
40
41
42
;    Adapted to extinction and polarization by Vincent Guillet (IAS)
;    Written by J.-Ph. Bernard
;    see evolution details on the dustem cvs maintained at CESR
;    Contact J.-Ph. Bernard (Jean-Philippe.Bernard@cesr.fr) in case of problems.
;-

;note that x is not used.

7cc6c3de   Jean-Philippe Bernard   added and removed...
43
44
;stop

427f1205   Jean-Michel Glorian   version 4.2 merged
45
IF keyword_set(help) THEN BEGIN
9be94157   Jean-Philippe Bernard   improved
46
  doc_library,'dustem_mpfit_run'
427f1205   Jean-Michel Glorian   version 4.2 merged
47
48
49
  sed=0.
  goto,the_end
ENDIF
427f1205   Jean-Michel Glorian   version 4.2 merged
50
51

p_dim=p_min*(*(*!dustem_fit).param_init_values)
427f1205   Jean-Michel Glorian   version 4.2 merged
52

427f1205   Jean-Michel Glorian   version 4.2 merged
53
IF !run_ionfrac gt 0. THEN toto = dustem_create_ionfrac()
452c334e   Ilyes Choubani   Implementation Of...
54

427f1205   Jean-Michel Glorian   version 4.2 merged
55
;RUN THE MODEL AND GET THE SPECTRUM
45da29c6   Ilyes Choubani   updating mpfit_run
56

45da29c6   Ilyes Choubani   updating mpfit_run
57
58
dustem_activate_plugins,p_dim/(*(*!dustem_fit).param_init_values),st ;st is an output...

36e8b879   Jean-Michel Glorian   update from deborah
59
st_model=dustem_read_all(!dustem_soft_dir)
427f1205   Jean-Michel Glorian   version 4.2 merged
60

45da29c6   Ilyes Choubani   updating mpfit_run
61

427f1205   Jean-Michel Glorian   version 4.2 merged
62
63
chi2_sed      = 0  
chi2_ext      = 0  
452c334e   Ilyes Choubani   Implementation Of...
64
65
66
67
chi2_qsed      =0
chi2_used      =0 
chi2_qext      =0
chi2_uext      =0 
18e4331f   Ilyes Choubani   general update (f...
68

427f1205   Jean-Michel Glorian   version 4.2 merged
69
70
rchi2_sed     = 0
rchi2_ext     = 0
452c334e   Ilyes Choubani   Implementation Of...
71
72
73
74
rchi2_qsed     =0
rchi2_used     =0
rchi2_qext     =0
rchi2_uext     =0
18e4331f   Ilyes Choubani   general update (f...
75

427f1205   Jean-Michel Glorian   version 4.2 merged
76
77
wrchi2_sed    = 0
wrchi2_ext    = 0
452c334e   Ilyes Choubani   Implementation Of...
78
79
80
81
wrchi2_qsed    =0
wrich2_used    =0
wrchi2_qext    =0
wrich2_uext    =0
18e4331f   Ilyes Choubani   general update (f...
82

427f1205   Jean-Michel Glorian   version 4.2 merged
83
84
n_sed         = 0
n_ext         = 0
452c334e   Ilyes Choubani   Implementation Of...
85
86
87
88
n_qsed         =0
n_used         =0
n_qext         =0
n_uext         =0
36e8b879   Jean-Michel Glorian   update from deborah
89
90

win = 0
427f1205   Jean-Michel Glorian   version 4.2 merged
91

0068116a   Ilyes Choubani   General update + ...
92
tagnames=tag_names((*!dustem_data))
427f1205   Jean-Michel Glorian   version 4.2 merged
93
dustem_out = [0]
36e8b879   Jean-Michel Glorian   update from deborah
94
95
ind_data_tag=dustem_data_tag(count=count_data_tag)

40597848   Ilyes Choubani   Plotting of numbe...
96

b5ccb706   Jean-Philippe Bernard   improved to fit p...
97
FOR i_tag=0,count_data_tag-1 DO BEGIN
427f1205   Jean-Michel Glorian   version 4.2 merged
98

36e8b879   Jean-Michel Glorian   update from deborah
99
    CASE tagnames(ind_data_tag(i_tag)) OF
266ae799   Ilyes Choubani   General update
100
        
b74c4452   Ilyes Choubani   Treating extincti...
101
		'SED'	: BEGIN 
5f04fa07   Ilyes Choubani   general update
102
	        dustem_sed = dustem_compute_sed(p_dim,st,SED_spec)
0068116a   Ilyes Choubani   General update + ...
103
104
	        chi2_sed   = total(((*(*!dustem_data).sed).values-dustem_sed)^2/(*(*!dustem_data).sed).sigma^2)
			n_sed = n_elements((*(*!dustem_data).sed).values)
b5ccb706   Jean-Philippe Bernard   improved to fit p...
105
106
107
	       	rchi2_sed  = chi2_sed / n_sed
	       	wrchi2_sed  = rchi2_sed * !fit_rchi2_weight.sed
	       	message,"chi2 SED = "+strtrim(chi2_sed,2)+"  rchi2 SED = "+strtrim(rchi2_sed,2),/continue;," weighted rchi2 SED=",wrchi2_sed
0068116a   Ilyes Choubani   General update + ...
108
109
110
111
112
113
114
115
			;print,'wl ',(*(*!dustem_data).sed).wav
			;print,"Delta SED per wl: ",(*(*!dustem_data).sed).values-dustem_sed
			;print,"error per wl: ",(*(*!dustem_data).sed).sigma
			;print,"chi2 per wl: ",((*(*!dustem_data).sed).values-dustem_sed)^2/(*(*!dustem_data).sed).sigma^2
	        alpha = total(((*(*!dustem_data).sed).values*dustem_sed)/(*(*!dustem_data).sed).sigma^2) $
	              		/ total(((*(*!dustem_data).sed).values)^2/(*(*!dustem_data).sed).sigma^2)
	        beta = total(((*(*!dustem_data).sed).values*dustem_sed)/(*(*!dustem_data).sed).sigma^2) $
	             		/ total(dustem_sed^2/(*(*!dustem_data).sed).sigma^2)
427f1205   Jean-Michel Glorian   version 4.2 merged
116
		
b5ccb706   Jean-Philippe Bernard   improved to fit p...
117
	        print,"Best SED chi2 obtained if f_HI = ",alpha, " or if we multiply model by ",beta
b74c4452   Ilyes Choubani   Treating extincti...
118
119
120
121
122
             
            ;These two lines were put here so that the chi2 is updated everytime we plot
            ;Moving it to the bottom so that all chi2 are taken into account  
;  	        (*!dustem_fit).chi2=chi2_sed
;  	        (*!dustem_fit).rchi2=rchi2_sed    
fcb6eade   Ilyes Choubani   general update - ...
123
	        
5a2643a4   Ilyes Choubani   cleaning/improvin...
124
	        dustem_out = [ dustem_out, dustem_sed ]
266ae799   Ilyes Choubani   General update
125
            
36e8b879   Jean-Michel Glorian   update from deborah
126
127
		END

b74c4452   Ilyes Choubani   Treating extincti...
128
129
		'EXT' 	: BEGIN     		
            dustem_ext = dustem_compute_ext(p_dim,st,EXT_spec)   
0068116a   Ilyes Choubani   General update + ...
130
131
	        chi2_ext   = total(((*(*!dustem_data).ext).values-dustem_ext)^2/(*(*!dustem_data).ext).sigma^2)
			n_ext = n_elements((*(*!dustem_data).ext).values)
b5ccb706   Jean-Philippe Bernard   improved to fit p...
132
			rchi2_ext  = chi2_ext / n_ext
b74c4452   Ilyes Choubani   Treating extincti...
133
134
	        wrchi2_ext  = rchi2_ext * !fit_rchi2_weight.ext
	        print,"chi2 EXT = ",chi2_ext,"  rchi2 EXT = ",rchi2_ext;," weighted rchi2 EXT=",wrchi2_ext
0068116a   Ilyes Choubani   General update + ...
135
        	;print,"chi2 per wl: ",((*(*!dustem_data).ext).values-dustem_ext)^2/(*(*!dustem_data).ext).sigma^2
b5ccb706   Jean-Philippe Bernard   improved to fit p...
136
137
138
139
140

		 	dustem_out = [ dustem_out, dustem_ext ]

		 	END

b74c4452   Ilyes Choubani   Treating extincti...
141
		'POLEXT': BEGIN 
b5ccb706   Jean-Philippe Bernard   improved to fit p...
142

759a527d   Ilyes Choubani   general update
143
			
b74c4452   Ilyes Choubani   Treating extincti...
144
			;Haven't investigated !run_circ=1 yet... 
759a527d   Ilyes Choubani   general update
145
			;---------------the procedure within this block has not been updated for a vert long time
b5ccb706   Jean-Philippe Bernard   improved to fit p...
146
147
148
149
150
			IF !dustem_show_plot NE 0 and !run_circ then begin
	                   win+=1
	                 
				dustem_plot_circ,st,win=win,xr=[0.1,5],aligned=st_model.align.grains.aligned,yr=[-0.03,0.03]
			endif
759a527d   Ilyes Choubani   general update
151
152
153
			;--------------------------------------------------------------------------------------
			
			
b5ccb706   Jean-Philippe Bernard   improved to fit p...
154
			IF !run_lin then begin
759a527d   Ilyes Choubani   general update
155
                    ;compute_polext needs updating. It's been months. 
67bd858a   Ilyes Choubani   Replication of th...
156
                    dustem_polext = dustem_compute_polext(p_dim,st,POLEXT_spec,SPEXT_spec,dustem_fpolext,dustem_ext=dustem_ext)
0068116a   Ilyes Choubani   General update + ...
157
	        		;chi2_polext   = total(((*(*!dustem_data).polext).values-dustem_polext)^2/(*(*!dustem_data).polext).sigma^2)
36e8b879   Jean-Michel Glorian   update from deborah
158

b5ccb706   Jean-Philippe Bernard   improved to fit p...
159
			ENDIF
b5ccb706   Jean-Philippe Bernard   improved to fit p...
160
			END
427f1205   Jean-Michel Glorian   version 4.2 merged
161

b74c4452   Ilyes Choubani   Treating extincti...
162
        'QSED': BEGIN 
45da29c6   Ilyes Choubani   updating mpfit_run
163
        		    
fcb6eade   Ilyes Choubani   general update - ...
164
		        toto = dustem_compute_stokes(p_dim,st,dustem_qsed,dustem_used,Q_spec,U_spec,PSI_spec,dustem_psi_em) 
28a7a1a7   Ilyes Choubani   Remaining fixes. ...
165
		      
45da29c6   Ilyes Choubani   updating mpfit_run
166
167
		        ;taking care of the color correction for spass ; SHOULD NOT LEAVE THIS ON
		        ;dustem_qsed(-1)=dustem_qsed(-2)
0068116a   Ilyes Choubani   General update + ...
168
; 		        testspass1 = ((*(*!dustem_data).qsed).filt_names)(-1) eq 'SPASS1' 
45da29c6   Ilyes Choubani   updating mpfit_run
169
170
;  			    if testspass1 then dustem_qsed(-1)=dustem_qsed(-2)
;  			    
0068116a   Ilyes Choubani   General update + ...
171
172
			    chi2_qsed = total(((*(*!dustem_data).qsed).values-dustem_qsed)^2/(*(*!dustem_data).qsed).sigma^2) 
			    n_qsed = n_elements((*(*!dustem_data).qsed).values)
45da29c6   Ilyes Choubani   updating mpfit_run
173
174
			    rchi2_qsed  = chi2_qsed / n_qsed
	        	wrchi2_qsed  = rchi2_qsed * !fit_rchi2_weight.qsed 
5a2643a4   Ilyes Choubani   cleaning/improvin...
175
	            dustem_out = [ dustem_out, dustem_qsed ] ; classic command
45da29c6   Ilyes Choubani   updating mpfit_run
176
	            message,"chi2 QSED = "+strtrim(chi2_qsed,2)+"  rchi2 QSED = "+strtrim(rchi2_qsed,2),/continue;," weighted rchi2 POLSED=",wrchi2_polsed
0068116a   Ilyes Choubani   General update + ...
177
	        	print,"chi2 per wl: ",((*(*!dustem_data).qsed).values-dustem_qsed)^2/(*(*!dustem_data).qsed).sigma^2
45da29c6   Ilyes Choubani   updating mpfit_run
178
179
		        
		    END   
1355825c   Ilyes Choubani   General update
180

b74c4452   Ilyes Choubani   Treating extincti...
181
182
        'USED': BEGIN 
    		     
0068116a   Ilyes Choubani   General update + ...
183
184
 			    chi2_used  = total(((*(*!dustem_data).used).values-dustem_used)^2/(*(*!dustem_data).used).sigma^2) 
 			    n_used = n_elements((*(*!dustem_data).used).values)
45da29c6   Ilyes Choubani   updating mpfit_run
185
186
 			    rchi2_used  = chi2_used / n_used
  	        	wrchi2_used  = rchi2_used * !fit_rchi2_weight.used
5a2643a4   Ilyes Choubani   cleaning/improvin...
187
 	            dustem_out = [ dustem_out, dustem_used] ; classic command
45da29c6   Ilyes Choubani   updating mpfit_run
188
 	            message,"chi2 USED = "+strtrim(chi2_used,2)+"  rchi2 USED = "+strtrim(rchi2_used,2),/continue;," weighted rchi2 POLSED=",wrchi2_polsed
0068116a   Ilyes Choubani   General update + ...
189
  	        	print,"chi2 per wl: ",((*(*!dustem_data).used).values-dustem_used)^2/(*(*!dustem_data).used).sigma^2
45da29c6   Ilyes Choubani   updating mpfit_run
190
191
192
                
  		END
  		
b74c4452   Ilyes Choubani   Treating extincti...
193
        'POLSED': BEGIN
45da29c6   Ilyes Choubani   updating mpfit_run
194
195
196
 			
             	  IF !run_pol and !run_lin THEN BEGIN
  
67bd858a   Ilyes Choubani   Replication of th...
197
            	       	dustem_polsed = dustem_compute_polsed(p_dim,st,P_spec,SP_spec,dustem_polfrac,dustem_sed=dustem_sed)
bffa7018   Ilyes Choubani   update
198
            	      	
45da29c6   Ilyes Choubani   updating mpfit_run
199
200
         		  ENDIF
		END
67bd858a   Ilyes Choubani   Replication of th...
201
        	
452c334e   Ilyes Choubani   Implementation Of...
202
 		
b74c4452   Ilyes Choubani   Treating extincti...
203
204
205
  		'QEXT': BEGIN ;UPDATED  
                
                toto = dustem_compute_stokext(p_dim,st,dustem_qext,dustem_uext,QEXT_spec,UEXT_spec,PSIEXT_spec,dustem_psi_ext) 
0068116a   Ilyes Choubani   General update + ...
206
207
                chi2_qext   =total(((*(*!dustem_data).qext).values-dustem_qext)^2/(*(*!dustem_data).qext).sigma^2) 
   			    n_qext = n_elements((*(*!dustem_data).qext).values)
452c334e   Ilyes Choubani   Implementation Of...
208
209
210
211
   			    rchi2_qext  = chi2_qext / n_qext
  	        	wrchi2_qext  = rchi2_qext * !fit_rchi2_weight.qext
   	            dustem_out = [ dustem_out, dustem_qext] ; classic command
   	            message,"chi2 QEXT = "+strtrim(chi2_qext,2)+"  rchi2 QEXT = "+strtrim(rchi2_qext,2),/continue;," weighted rchi2 POLSED=",wrchi2_polsed
0068116a   Ilyes Choubani   General update + ...
212
  	        	print,"chi2 per wl: ",((*(*!dustem_data).qext).values-dustem_qext)^2/(*(*!dustem_data).qext).sigma^2 
452c334e   Ilyes Choubani   Implementation Of...
213
214
215
     	           
  		END
;  		
b74c4452   Ilyes Choubani   Treating extincti...
216
217
  		'UEXT': BEGIN ;UPDATED
  		
0068116a   Ilyes Choubani   General update + ...
218
219
 			    chi2_uext   =total(((*(*!dustem_data).uext).values-dustem_uext)^2/(*(*!dustem_data).uext).sigma^2) 
   			    n_uext = n_elements((*(*!dustem_data).uext).values)
452c334e   Ilyes Choubani   Implementation Of...
220
221
222
223
   			    rchi2_uext  = chi2_uext / n_uext
  	        	wrchi2_uext  = rchi2_uext * !fit_rchi2_weight.uext
   	            dustem_out = [ dustem_out, dustem_uext] ; classic command
   	            message,"chi2 UEXT = "+strtrim(chi2_uext,2)+"  rchi2 UEXT = "+strtrim(rchi2_uext,2),/continue;," weighted rchi2 POLSED=",wrchi2_polsed
0068116a   Ilyes Choubani   General update + ...
224
  	        	print,"chi2 per wl: ",((*(*!dustem_data).uext).values-dustem_uext)^2/(*(*!dustem_data).uext).sigma^2
b74c4452   Ilyes Choubani   Treating extincti...
225
            
452c334e   Ilyes Choubani   Implementation Of...
226
227
 	           
  		END
67bd858a   Ilyes Choubani   Replication of th...
228
      	
b74c4452   Ilyes Choubani   Treating extincti...
229
230
	    ELSE : ;Do nothing
	    
427f1205   Jean-Michel Glorian   version 4.2 merged
231
232
    ENDCASE

b5ccb706   Jean-Philippe Bernard   improved to fit p...
233
ENDFOR
b74c4452   Ilyes Choubani   Treating extincti...
234

9eb53891   Jean-Philippe Bernard   modified call to ...
235
;=== computing chi^2
36e8b879   Jean-Michel Glorian   update from deborah
236
237
238
239
240
241
242
243
244
245
246
247
248
249
250
DOF = n_elements(p_min)
total_chi2 = 0
total_rchi2 = 0
total_n = 0
if !fit_rchi2_weight.ext ne 0 then begin
	total_n += n_ext
	total_chi2 +=  chi2_ext
	total_rchi2 += rchi2_ext
endif
if !fit_rchi2_weight.sed ne 0 then begin
	total_n += n_sed
	total_chi2 += chi2_sed
	total_rchi2 += rchi2_sed
endif
if !run_pol then begin
1355825c   Ilyes Choubani   General update
251
252
253
254
255
256
257
258
259
260
261
262
263
264
265
266
267
268
269
270
271

    if !fit_rchi2_weight.qsed ne 0 then begin
    	total_n += n_qsed
    	total_chi2 += chi2_qsed
    	total_rchi2 += rchi2_qsed
    endif
    if !fit_rchi2_weight.used ne 0 then begin
    	total_n += n_used
    	total_chi2 += chi2_used
    	total_rchi2 += rchi2_used
    endif
    if !fit_rchi2_weight.qext ne 0 then begin
    	total_n += n_qext
    	total_chi2 += chi2_qext
    	total_rchi2 += rchi2_qext
    endif
    if !fit_rchi2_weight.uext ne 0 then begin
    	total_n += n_uext
    	total_chi2 += chi2_uext
    	total_rchi2 += rchi2_uext
    endif
452c334e   Ilyes Choubani   Implementation Of...
272
273
274

endif

36e8b879   Jean-Michel Glorian   update from deborah
275
276
277
278

total_wchi2 = chi2_ext * !fit_rchi2_weight.ext $
            + chi2_sed * !fit_rchi2_weight.sed 
if !run_pol then begin
18e4331f   Ilyes Choubani   general update (f...
279
	total_wchi2 	+=  $;chi2_polext  * !fit_rchi2_weight.polext $
3ff9cac3   Ilyes Choubani   general update
280
281
;              		 + chi2_polsed  * !fit_rchi2_weight.polsed $
;              		 + chi2_polfrac * !fit_rchi2_weight.polfrac $
452c334e   Ilyes Choubani   Implementation Of...
282
283
284
285
             		 + chi2_qsed * !fit_rchi2_weight.qsed $
             		 + chi2_used * !fit_rchi2_weight.used $
             		 + chi2_qext * !fit_rchi2_weight.qext $
             		 + chi2_uext * !fit_rchi2_weight.uext  
36e8b879   Jean-Michel Glorian   update from deborah
286
287
endif

b74c4452   Ilyes Choubani   Treating extincti...
288
289
290
291
292
;Using this expression for now for the weighting. 
(*!dustem_fit).chi2 = total_wchi2
(*!dustem_fit).rchi2 = total_wchi2/(total_n-DOF)


9ccf7615   Jean-Philippe Bernard   modified to fit u...
293
message,'======================================================================',/continue
36e8b879   Jean-Michel Glorian   update from deborah
294
295
296
print,"Total weighted chi2 driving mpfitfun = ", total_wchi2
print,"Nb of data points = ", total_n
print,"DOF = ", DOF
18e4331f   Ilyes Choubani   general update (f...
297
298
299
print,"Reduced chi2 EXT/SED/POLEXT/POLSED/POLFRAC",rchi2_ext,rchi2_sed;,rchi2_polext;,rchi2_polsed,rchi2_polfrac
print,"Weighted reduced chi2 EXT/SED/POLEXT/POLSED/POLFRAC",wrchi2_ext,wrchi2_sed;,wrchi2_polext;,wrchi2_polsed,wrchi2_polfrac
if !run_pol then 	print,"Weights",!fit_rchi2_weight.ext,!fit_rchi2_weight.sed $;,!fit_rchi2_weight.polext $;,!fit_rchi2_weight.polsed,!fit_rchi2_weight.polfrac $ 
36e8b879   Jean-Michel Glorian   update from deborah
300
301
302
else 			print,"Weights",!fit_rchi2_weight.ext,!fit_rchi2_weight.sed
print,"Mean weighted reduced chi2 = ", total_wchi2/(total_n-DOF)
print,"Unweighted reduced Total chi2 = ", total_chi2/(total_n-DOF)
b74c4452   Ilyes Choubani   Treating extincti...
303
print,"Unweighted mean reduced chi2 = ", total_rchi2/4*total_n/(total_n-DOF);total_rchi2/4*total_n/(total_n-DOF) original ; I do not understand this line. Maybe redo some computation.
9ccf7615   Jean-Philippe Bernard   modified to fit u...
304
message,'======================================================================',/continue
386ad2de   Jean-Philippe Bernard   cant remeber what...
305
306
;the following is needed under fawlty (not IDL) which message,/reset does not work like in IDL
!Error_State.msg=''
9eb53891   Jean-Philippe Bernard   modified call to ...
307
308

;=== plotting the fit
74043007   Jean-Philippe Bernard   included tet on d...
309
310
311
312
313
314
if !dustem_iter.act NE !dustem_iter.prv or !dustem_iter.prv EQ 1 then begin
    IF !dustem_noobj EQ 0 THEN BEGIN
      dustemwrap_plot,p_dim,st,dustem_sed,SED_spec,dustem_qsed,Q_spec,dustem_used,U_spec,dustem_polsed,P_spec,dustem_polfrac,SP_spec, dustem_psi_em, PSI_spec,dustem_ext,EXT_spec,dustem_qext,QEXT_spec,dustem_uext,UEXT_spec,dustem_polext,POLEXT_spec,dustem_fpolext,SPEXT_spec,dustem_psi_ext,PSIEXT_spec,_extra=_extra
    ENDIF ELSE BEGIN
      dustemwrap_plot_noobj,p_dim,st,dustem_sed,SED_spec,dustem_qsed,Q_spec,dustem_used,U_spec,dustem_polsed,P_spec,dustem_polfrac,SP_spec, dustem_psi_em, PSI_spec,dustem_ext,EXT_spec,dustem_qext,QEXT_spec,dustem_uext,UEXT_spec,dustem_polext,POLEXT_spec,dustem_fpolext,SPEXT_spec,dustem_psi_ext,PSIEXT_spec,_extra=_extra
    ENDELSE
9eb53891   Jean-Philippe Bernard   modified call to ...
315
316
317
318
319
320
    if !dustem_iter.act EQ 1 then begin
        !dustem_iter.prv=0
        !dustem_iter.act=0 
    endif else !dustem_iter.prv=!dustem_iter.act
endif

d1b8cc75   Jean-Philippe Bernard   add then commente...
321
;stop
36e8b879   Jean-Michel Glorian   update from deborah
322
323
324
325
326
327
328
329
330

;print,"Total reduced chi2 = ", (wrchi2_sed+wrchi2_ext+wrchi2_polext+wrchi2_polsed+wrchi2_polfrac)/(n_sed+n_ext+n_polext+n_polsed+n_polfrac-DOF)
; Reduced chi2 is like if there were the same nb of data points for each curve
;print,"Mean reduced chi2 = ", (wrchi2_sed+wrchi2_ext+wrchi2_polext+wrchi2_polsed+wrchi2_polfrac) * total_n/(total_n-DOF)

; Plot size distribution
win+=1
xr=[0.3,5e3]
yr=[0.1,100]
9ccf7615   Jean-Philippe Bernard   modified to fit u...
331
;IF !dustem_show_plot NE 0 THEN dustem_plot_sdist,xr=xr,yr=yr,ps=filename,win=win;,/donotclose
427f1205   Jean-Michel Glorian   version 4.2 merged
332
333
334
335

; Suppress the artificial [0] at the beginning 
dustem_out = dustem_out[1:*]

9ccf7615   Jean-Philippe Bernard   modified to fit u...
336
337
338
339
340
341
342
;IF defined(extra) THEN BEGIN
;  IF keyword_set(extra.plot_it) THEN BEGIN
;;    message,'params='+strtrim(p_min,2),/info,/continue
;;    stop
;    print,'dustem_mpfit_run_vg: Current parameters values:',p_dim
;  ENDIF
;ENDIF
427f1205   Jean-Michel Glorian   version 4.2 merged
343
344
345
346
347
348

the_end:

RETURN,dustem_out

END