OV-GSO-DC / dustem-wrapper_idl

Commit 0cf486f6ed8b410ef06649a0485cbd8587250372

Authored by Annie Hughes
2 parents 7ecea017 efd9ac73
Exists in master

Merge branch 'master' of https://gitlab.irap.omp.eu/OV-GSO-DC/dustem-wrapper_idl

Showing 5 changed files with 62 additions and 18 deletions   Show diff stats
Inline Side-by-side
src/idl/dustem_fit_intensity_example.pro
  Diff comments   View file @0cf486f
... ... @@ -370,7 +370,9 @@ end
370 370 ;I think modifications will have to be applied to dustemwrap_plot...
371 371 ;idea: create a system variable that signals if the fits file has been successfully saved and use it in the dustemwrap_plot(obj/noobj) routine(s)
372 372 ;to add the '(From saved FITS file)' string to the title
373   -;
  373 +
  374 +;stop
  375 +
374 376 IF keyword_set(fits_save_and_restore) THEN BEGIN
375 377 message,'Writing out results structure: '+fits_save_and_restore,/info
376 378 dustem_write_fits_table,filename=fits_save_and_restore,help=help
... ...
src/idl/dustem_init.pro
  Diff comments   View file @0cf486f
... ... @@ -349,7 +349,7 @@ defsysv,'!dustem_show_plot',1 ;0=show no plot
349 349 ;== This is for the default mode
350 350 dustem_inputs={GRAIN:'GRAIN.DAT',ISRF:'ISRF.DAT',ALIGN:'ALIGN.DAT'}
351 351 !dustem_inputs=ptr_new(dustem_inputs)
352   -
  352 +use_model='DEFAULT'
353 353  
354 354 defsysv, '!dustem_HCD', ptr_new(1.00E20) ;I do not see the reason behind using a pointer in this case.
355 355  
... ... @@ -357,10 +357,6 @@ defsysv, '!dustem_HCD', ptr_new(1.00E20) ;I do not see the reason behind using a
357 357 ;===This is the Desert option
358 358  
359 359 IF keyword_set(model) THEN BEGIN
360   -
361   -defsysv, '!dustem_model', model
362   -
363   -
364 360 CASE model OF
365 361 'MC10':BEGIN
366 362 (*!dustem_inputs).grain='GRAIN_MC10.DAT'
... ... @@ -414,8 +410,11 @@ defsysv, '!dustem_model', model
414 410 stop
415 411 END
416 412 ENDCASE
  413 + use_model=model
417 414 ENDIF
418 415  
  416 +defsysv, '!dustem_model', use_model
  417 +
419 418 ;Initialize !dustem_params structure
420 419 IF not keyword_set(dir) THEN BEGIN
421 420 ; CASE getenv('DUSTEM_WHICH') OF
... ...
src/idl/dustem_initialize_internal_sed.pro 0 → 100644
  Diff comments   View file @0cf486f
... ... @@ -0,0 +1,9 @@
  1 +FUNCTION dustem_initialize_internal_sed,Nsed
  2 +
  3 +;initializes an SED in the !dustem_data internal format
  4 +
  5 +sed={instru_names:replicate('',Nsed),filt_names:replicate('',Nsed),wav:replicate(la_undef(5),Nsed),values:replicate(la_undef(5),Nsed),sigma:replicate(la_undef(5),Nsed)}
  6 +
  7 +RETURN,sed
  8 +
  9 +END
0 10 \ No newline at end of file
... ...
src/idl/dustem_read_fits_table.pro
  Diff comments   View file @0cf486f
... ... @@ -60,7 +60,20 @@ str_predicted_SED=mrdfits(file,unit_predicted_sed,header_predicted_SED)
60 60 ;==== add the predicted total spectrum in extension 3
61 61 ;str_predicted_spectrum_tot=mrdfits(file,3,header_predicted_spectrum_tot)
62 62  
  63 +used_model=strtrim(sxpar(header_predicted_SED,'MODEL'),2)
  64 +IF used_model NE 'DEFAULT' THEN BEGIN
  65 + dustem_init,model=used_model
  66 +ENDIF ELSE BEGIN
  67 + dustem_init,model=0
  68 +ENDELSE
  69 +!dustem_model=used_model
  70 +
  71 +;stop
  72 +;==== This is to create (*!dustem_data).sed in case it does not exist already
63 73 Nsed=n_elements(str_input_SED)
  74 +sed=dustem_initialize_internal_sed(Nsed)
  75 +(*!dustem_data).sed=ptr_new(sed)
  76 +
64 77 (*(*!dustem_data).sed).instru_names=dustem_filter2instru(str_input_SED.filter)
65 78 (*(*!dustem_data).sed).filt_names=str_input_SED.filter
66 79 (*(*!dustem_data).sed).wav=str_input_SED.wavelength
... ... @@ -77,29 +90,34 @@ best_fit_chi2=sxpar(header_predicted_SED,'CHI2')
77 90 best_fit_rchi2=sxpar(header_predicted_SED,'RCHI2')
78 91  
79 92 ;===== Fill in header with parameter values
80   -parameters_desc=*(*!dustem_fit).PARAM_DESCS
81   -parameters_init_values=*(*!dustem_fit).PARAM_INIT_VALUES
82   -best_parameters=*(*!dustem_fit).CURRENT_PARAM_VALUES
83   -best_parameters_uncertainties=*(*!dustem_fit).CURRENT_PARAM_ERRORS
84   -best_fit_chi2=(*!dustem_fit).CHI2
85   -best_fit_rchi2=(*!dustem_fit).RCHI2
  93 +;parameters_desc=*(*!dustem_fit).PARAM_DESCS
  94 +;parameters_init_values=*(*!dustem_fit).PARAM_INIT_VALUES
  95 +;best_parameters=*(*!dustem_fit).CURRENT_PARAM_VALUES
  96 +;best_parameters_uncertainties=*(*!dustem_fit).CURRENT_PARAM_ERRORS
  97 +;best_fit_chi2=(*!dustem_fit).CHI2
  98 +;best_fit_rchi2=(*!dustem_fit).RCHI2
86 99  
87 100 Nparams=sxpar(header_predicted_SED,'NPARAMS')
88 101 parameters_desc=strarr(Nparams)
89 102 best_parameters=dblarr(Nparams)
90 103 best_parameters_init_values=dblarr(Nparams)
91 104 best_parameters_uncertainties=dblarr(Nparams)
  105 +parameters_func_values=dblarr(Nparams)
92 106 FOR i=0L,Nparams-1 DO BEGIN
93 107 parameters_desc[i]=sxpar(header_predicted_SED,'PARN'+strtrim(i+1,2))
94 108 best_parameters[i]=sxpar(header_predicted_SED,'PARV'+strtrim(i+1,2))
95 109 best_parameters_uncertainties[i]=sxpar(header_predicted_SED,'PARU'+strtrim(i+1,2))
96 110 best_parameters_init_values[i]=sxpar(header_predicted_SED,'PARI'+strtrim(i+1,2))
  111 + parameters_func_values[i]=sxpar(header_predicted_SED,'PARI'+strtrim(i+1,2))
97 112 ENDFOR
98 113  
99   -*(*!dustem_fit).PARAM_DESCS=parameters_desc
100   -*(*!dustem_fit).CURRENT_PARAM_VALUES=best_parameters
101   -*(*!dustem_fit).CURRENT_PARAM_ERRORS=best_parameters_uncertainties
102   -*(*!dustem_fit).PARAM_INIT_VALUES=best_parameters_init_values
  114 +(*!dustem_fit).PARAM_DESCS=ptr_new(parameters_desc)
  115 +(*!dustem_fit).PARAM_INIT_VALUES=ptr_new(best_parameters_init_values)
  116 +(*!dustem_fit).CURRENT_PARAM_VALUES=ptr_new(best_parameters)
  117 +(*!dustem_fit).CURRENT_PARAM_ERRORS=ptr_new(best_parameters_uncertainties)
  118 +(*!dustem_fit).CHI2=best_fit_chi2
  119 +(*!dustem_fit).RCHI2=best_fit_rchi2
  120 +(*!dustem_fit).PARAM_FUNC=ptr_new(parameters_func_values)
103 121  
104 122 ;=== Below is only to get the same output form as dustem_compute_sed
105 123 ;Note: If the same model is used, the outputs should be the same as what is in the fits file ....
... ...
src/idl/dustem_write_fits_table.pro
  Diff comments   View file @0cf486f
... ... @@ -47,6 +47,8 @@ IF not keyword_set(function_name) THEN BEGIN
47 47 ENDIF ELSE BEGIN
48 48 ENDELSE
49 49  
  50 +used_model=!dustem_model
  51 +
50 52 Ngrains=n_tags(dustem_spectra_st.sed)-2
51 53  
52 54 ;===== define structures containing results
... ... @@ -89,6 +91,7 @@ best_parameters=*(*!dustem_fit).CURRENT_PARAM_VALUES
89 91 best_parameters_uncertainties=*(*!dustem_fit).CURRENT_PARAM_ERRORS
90 92 best_fit_chi2=(*!dustem_fit).CHI2
91 93 best_fit_rchi2=(*!dustem_fit).RCHI2
  94 +parameters_func_values=*(*!dustem_fit).PARAM_FUNC
92 95  
93 96 Nparams=n_elements(parameters_desc)
94 97  
... ... @@ -169,6 +172,12 @@ FOR i=0L,Nparams-1 DO BEGIN
169 172 sxaddpar,header_predicted_SED,'PARN'+strtrim(i+1,2),parameters_desc[i],' DustemWrap parameter name'
170 173 ENDFOR
171 174 sxaddpar,header_predicted_SED,'COMMENT','===================================================',after='PARN'+strtrim(i,2)
  175 +
  176 +;====== write dustem wrap parameter best fit values
  177 +sxaddpar,header_predicted_SED,'COMMENT','======= Dustemwrap dust model used ======',after='COMMENT'
  178 +sxaddpar,header_predicted_SED,'MODEL',used_model,'Dust model used'
  179 +sxaddpar,header_predicted_SED,'COMMENT','===================================================',after='MODEL'
  180 +
172 181 ;====== write dustem wrap parameter best fit values
173 182 sxaddpar,header_predicted_SED,'COMMENT','======= Dustemwrap parameter best fit values ======',after='COMMENT'
174 183 FOR i=0L,Nparams-1 DO BEGIN
... ... @@ -181,14 +190,21 @@ sxaddpar,header_predicted_SED,'COMMENT','======= Dustemwrap best fit parameter u
181 190 FOR i=0L,Nparams-1 DO BEGIN
182 191 sxaddpar,header_predicted_SED,'PARU'+strtrim(i+1,2),best_parameters_uncertainties[i],'Param. '+parameters_desc[i]
183 192 ENDFOR
184   -sxaddpar,header_predicted_SED,'COMMENT','===================================================',after='PARIU'+strtrim(i,2)
  193 +sxaddpar,header_predicted_SED,'COMMENT','===================================================',after='PARU'+strtrim(i,2)
185 194  
186 195 ;====== write dustem wrap parameter initial fit values
187 196 sxaddpar,header_predicted_SED,'COMMENT','======= Dustemwrap parameter initial fit values ======',after='COMMENT'
188 197 FOR i=0L,Nparams-1 DO BEGIN
189 198 sxaddpar,header_predicted_SED,'PARI'+strtrim(i+1,2),parameters_init_values[i],'Param. '+parameters_desc[i]
190 199 ENDFOR
191   -sxaddpar,header_predicted_SED,'COMMENT','===================================================',after='PARIV'+strtrim(i,2)
  200 +sxaddpar,header_predicted_SED,'COMMENT','===================================================',after='PARI'+strtrim(i,2)
  201 +
  202 +;====== write dustem wrap parameter func values (not sure exactly what this is)
  203 +sxaddpar,header_predicted_SED,'COMMENT','======= Dustemwrap parameter func values ======',after='COMMENT'
  204 +FOR i=0L,Nparams-1 DO BEGIN
  205 + sxaddpar,header_predicted_SED,'PARF'+strtrim(i+1,2),parameters_func_values[i],'Param. '+parameters_desc[i]
  206 +ENDFOR
  207 +sxaddpar,header_predicted_SED,'COMMENT','===================================================',after='PARF'+strtrim(i,2)
192 208  
193 209 ;===== complement extension header for predicted spectrum total
194 210 sxaddpar,header_predicted_spectrum_tot,'TITLE','OUTPUT BEST FIT TOTAL SPECTRUM BY DUSTEMWRAP'
... ...