OV-GSO-DC / dustem-wrapper_idl

Commit fd152fd20bb8ba9ea15b431be9879b82e3111859

Authored by Annie Hughes
1 parent 1cf2bf59
Exists in master

cosmetic changes for consistency with tutorial text

Showing 2 changed files with 80 additions and 45 deletions   Show diff stats
Inline Side-by-side
src/idl/dustem_fit_polarization_example.pro
  Diff comments   View file @fd152fd
... ... @@ -120,12 +120,12 @@ pd = [ $
120 120 '(*!dustem_params).grains(0).mdust_o_mh',$ ;PAH0 mass fraction
121 121 '(*!dustem_params).grains(1).mdust_o_mh',$ ;amorphous carbon mass fraction
122 122 '(*!dustem_params).grains(2).mdust_o_mh', $ ;amorphous silicate mass fraction
123   - 'dustem_plugin_modify_dust_pol_2', $ ;This will set the polarization angle
  123 + 'dustem_plugin_modify_dust_pol_2', $ ;Dust polarization angle
124 124 'dustem_plugin_synchrotron_2', $ ;Synchrotron amplitude at 10 mm
125 125 'dustem_plugin_synchrotron_3', $ ;Synchrotron polarization fraction
126 126 'dustem_plugin_synchrotron_4' $ ;Synchrotron polarization angle
127 127 ]
128   -iv=[7.2e-4,7.2e-4,8.5e-4,-10.,0.015,0.25,40.]
  128 +iv=[7.2e-4,7.2e-4,8.5e-4,-10.,0.015,0.25,35.]
129 129  
130 130 Npar=n_elements(pd)
131 131 ulimed=[0,0,0,0,1,1,1]
... ... @@ -223,7 +223,7 @@ IF keyword_set(fits_save_and_restore) THEN BEGIN
223 223 dustem_read_fits_table,filename=fits_save_and_restore,dustem_spectra_st=dustem_spectra_st
224 224 ;==== plot result taken from the saved fits table
225 225 res=*(*!dustem_fit).CURRENT_PARAM_VALUES
226   - dustemwrap_plot,res,dustem_spectra_st,xr=xr,/xstyle,yr=yr,/ysty,/ylog,/xlog,title=tit+' (From Saved FITS file)'
  226 + dustemwrap_plot,res,dustem_spectra_st,xr=xr,/xstyle,yr=yr,/ysty,/ylog,/xlog,title=title+' (From saved FITS file)'
227 227 if keyword_set(wait) then begin
228 228 message,'Saved the results as FITS in the file: '+fits_save_and_restore+', and made a plot using the data in this file',/info
229 229 wait,wait
... ...
src/idl/dustem_make_polarization_sed_example.pro
  Diff comments   View file @fd152fd
1   -PRO dustem_make_polarization_sed_example,model=model,help=help
  1 +PRO dustem_make_polarization_sed_example,model=model,filters=filters,help=help,outfile=outfile
2 2  
3 3 ;+
4 4 ; NAME:
5 5 ; dustem_make_polarization_sed_example
6 6 ; PURPOSE:
7   -; This is an example of how to generate an SED using dustemwrap.
  7 +; This is an example of how to generate an SED using DustEMWrap.
8 8 ; It is meant to be an example to follow when writing your own
9   -; programs using the dustem IDL/GDL wrapper.
  9 +; programs using DustEMWrap
10 10 ; CATEGORY:
11   -; Dustem
  11 +; DustEMWrap, Distributed, High-Level, User Example
12 12 ; CALLING SEQUENCE:
13   -; dustem_make_polarization_sed_example[,model=][,/help]
  13 +; dustem_make_polarization_sed_example[,model=][,filters=][,/help]
14 14 ; INPUTS:
15 15 ; None
16 16 ; OPTIONAL INPUT PARAMETERS:
... ... @@ -18,51 +18,86 @@ PRO dustem_make_polarization_sed_example,model=model,help=help
18 18 ; OUTPUTS:
19 19 ; None
20 20 ; OPTIONAL OUTPUT PARAMETERS:
21   -; None
  21 +; outfile = path+filename for .xcat output of the SED
22 22 ; ACCEPTED KEY-WORDS:
23   -; model = if set, name of the Dustem dust model to use (default='G17_MODELC')
  23 +; model = if set, name of the DustEM dust model to use (default='G17_MODELC')
24 24 ; help = If set print this help
25 25 ; COMMON BLOCKS:
26 26 ; None
27 27 ; SIDE EFFECTS:
28 28 ; None
29 29 ; RESTRICTIONS:
30   -; The dustem idl wrapper must be installed
  30 +; The DustEM fortran code must be installed
  31 +; The DustEMWrap IDL code must be installed
31 32 ; PROCEDURE:
32 33 ; None
33 34 ; EXAMPLES
34 35 ; dustem_make_polarization_sed_example
35 36 ; dustem_make_polarization_sed_example,model='G17_MODELD'
36 37 ; MODIFICATION HISTORY:
37   -; Written by J.P. Bernard June 29th 2011
38   -; see evolution details on the dustem cvs maintained at IRAP
39   -; Contact J.-Ph. Bernard (Jean-Philippe.Bernard@cesr.fr) in case of problems.
  38 +; Written by JPB June-2022
  39 +; Evolution details on the DustEMWrap gitlab.
  40 +; See http://dustemwrap.irap.omp.eu/ for FAQ and help.
40 41 ;-
41 42  
42 43 IF keyword_set(help) THEN BEGIN
43   - doc_library,'dustem_fit_sed_polsed_readme'
  44 + doc_library,'dustem_make_polarization_sed_example'
44 45 goto,the_end
45 46 ENDIF
46 47  
47 48 dustem_define_la_common
48 49  
49   -;=== initialise dustem
50   -use_mode='G17_MODELC'
51   -dustem_init,mode=use_mode,kwords=['logn-chrg-spin','plaw-pol','plaw-pol'],/pol
  50 +IF keyword_set(model) THEN BEGIN
  51 + use_model=strupcase(model)
  52 +ENDIF ELSE BEGIN
  53 + use_model='G17_MODELD' ;Default is one of the Guillet et al (2017) models since they treat dust polarisation properties
  54 +ENDELSE
  55 +
  56 +IF keyword_set(filters) THEN BEGIN
  57 + use_filters=strupcase(filters)
  58 +ENDIF ELSE BEGIN
  59 + use_filters=['IRAS1','IRAS2','IRAS3','IRAS4','PACS3','SPIRE1','SPIRE2','SPIRE3','HFI2','HFI3','HFI4','HFI5','HFI6','LFI1','LFI2','LFI3']
  60 +ENDELSE
  61 +
  62 +IF keyword_set(outfile) THEN BEGIN
  63 + use_outfile=outfile
  64 +ENDIF ELSE BEGIN
  65 + use_outfile='./polarization_sed.xcat'
  66 +ENDELSE
  67 +
  68 +known_mdls=['MC10','DBP90','DL01','WD01_RV5P5B','DL07','J13','G17_MODELA','G17_MODELB','G17_MODELC','G17_MODELD']
  69 +pol_mdls=['G17_MODELA','G17_MODELB','G17_MODELC','G17_MODELD']
  70 +test_model = where(known_mdls eq use_model,ct)
  71 +if ct eq 0 then begin
  72 + message,'ISM dust model '+use_model+' unknown',/continue
  73 + message,'Known models are MC10,DBP90,DL01,WD01_RV5P5B,DL07,J13,G17_MODELA,G17_MODELB,G17_MODELC,G17_MODELD',/continue
  74 + stop
  75 +end
  76 +test_model = where(pol_mdls eq use_model,ct)
  77 +if ct eq 0 then begin
  78 + message,'The only models with polarisation are G17_MODELA,G17_MODELB,G17_MODELC,G17_MODELD',/continue
  79 + stop
  80 +end
  81 +
  82 +
  83 +if keyword_set(verbose) then use_verbose=1
  84 +
  85 +dustem_define_la_common
52 86  
  87 +;dustem_init,mode=use_model,grain_keywords=['logn-chrg-spin','?','plaw-pol'],polarization=1
  88 +dustem_init,mode=use_model,polarization=1
53 89 !dustem_verbose=1
54 90 !dustem_show_plot=1
55 91 !dustem_nocatch=1
56 92  
57   -;=== define filters for the SED
58   -filters=['IRAS1','IRAS2','IRAS3','IRAS4','PACS3','SPIRE1','SPIRE2','SPIRE3','HFI2','HFI3','HFI4','HFI5','HFI6','LFI1','LFI2','LFI3']
59   -Nfilt=n_elements(filters)
  93 +;=== initialize filters for the SED
  94 +Nfilt=n_elements(use_filters)
60 95 sed=dustem_initialize_sed(Nfilt)
61   -sed.filter=filters
62   -sed.wave=dustem_filter2wav(filters)
63   -sed.instru=dustem_filter2instru(filters)
  96 +sed.filter=use_filters
  97 +sed.wave=dustem_filter2wav(use_filters)
  98 +sed.instru=dustem_filter2instru(use_filters)
64 99  
65   -;=== initializing IQU and associated errors to avoid problems when checking SED in dustem_set_data.pro
  100 +;=== initialize IQU and associated errors to avoid problems when checking SED in dustem_set_data.pro
66 101 sed[*].StokesI=1.
67 102 sed.StokesQ=sed.StokesI/100.
68 103 sed.StokesU=sed.StokesI/100.
... ... @@ -73,35 +108,36 @@ sed.SigmaIQ=sed.StokesI/100.
73 108 sed.SigmaIU=sed.StokesI/100.
74 109 sed.SigmaQU=sed.StokesI/100.
75 110  
76   -;=== initializing the !dustem_data structure that will be internally used
77   -dustem_set_data,sed
  111 +;=== initialize the !dustem_data structure that will be used internally
  112 +dustem_set_data,sed,sed
78 113  
79 114 ;help,!dustem_data
80 115  
81 116 ;=== Set which model parameters to use for the SED
82 117 pd = [ $
83 118 '(*!dustem_params).G0', $ ;G0
84   - '(*!dustem_params).grains(0).mdust_o_mh',$ ;PAH0 mass fraction
85   - '(*!dustem_params).grains(1).mdust_o_mh',$ ;PAH0 mass fraction
86   - '(*!dustem_params).grains(2).mdust_o_mh', $ ;PAH1 mass fraction
87   - 'dustem_plugin_modify_dust_pol_2', $ ;This will set the polarization angle
  119 + '(*!dustem_params).grains(0).mdust_o_mh',$ ;PAH0_MC10 mass fraction
  120 + '(*!dustem_params).grains(1).mdust_o_mh',$ ;amCBE_0.3333x mass fraction
  121 + '(*!dustem_params).grains(2).mdust_o_mh', $ ;aSil2001BE6pctG_0.4x mass fraction
  122 + 'dustem_plugin_modify_dust_pol_2', $ ;Dust polarization angle
88 123 'dustem_plugin_synchrotron_1', $ ;Synchrotron spectra index
89 124 'dustem_plugin_synchrotron_2', $ ;Synchrotron amplitude at 10 mm
90 125 'dustem_plugin_synchrotron_3', $ ;Synchrotron polarization fraction
91 126 'dustem_plugin_synchrotron_4' $ ;Synchrotron polarization angle
92 127 ]
93   -p_truth=[1.,7.8000E-04,7.8000E-04,7.8000E-04,10.,3., 1.e-2, 0.3, 45.]
  128 +pd_true=[1.,7.8000E-04,7.8000E-04,7.8000E-04, 10.,3., 1.e-2, 0.3, 45.]
94 129  
95 130 Npar=n_elements(pd)
96 131 ulimed=replicate(0,Npar)
97 132 llimed=replicate(1,Npar)
98 133 llims=replicate(0,Npar)
99 134  
100   -dustem_init_parinfo,pd,p_truth,up_limited=ulimed,lo_limited=llimed,up_limits=ulims,lo_limits=llims
101   -dustem_init_plugins,pd
  135 +dustem_init_params,use_model,pd,pd_true,pol=1
  136 +;dustem_init_parinfo,pd,pd_true,up_limited=ulimed,lo_limited=llimed,up_limits=ulims,lo_limits=llims
  137 +;dustem_init_plugins,pd
102 138  
103   -dustem_Ised=dustem_compute_sed(p_truth,ssti) ;sst is not set on purpose, for dustem_compute_sed to compute the sed using fortran code
104   -toto=dustem_compute_stokes(p_truth,sstqu,dustem_Qsed,dustem_Used) ;this procedure also allows for the extraction of the spectra
  139 +dustem_Ised=dustem_compute_sed(pd_true,ssti) ;sst is not set on purpose, for dustem_compute_sed to compute the sed using fortran code
  140 +toto=dustem_compute_stokes(pd_true,sstqu,dustem_Qsed,dustem_Used) ;this procedure also allows for the extraction of the spectra
105 141  
106 142 ;print,dustem_Qsed,dustem_Used
107 143  
... ... @@ -124,16 +160,15 @@ sed.sigmaQU=0.
124 160 ;==== fill in dependent columns of the SED.
125 161 sed=dustem_fill_sed_dependent_columns(sed)
126 162  
127   -;======== save the SED
128   -filename='/tmp/fake_polarization_sed.xcat'
129   -write_xcat,sed,filename
130   -
131   -message,'Wrote '+filename,/continue
  163 +;======== save the SED to a file
  164 +write_xcat,sed,use_outfile
  165 +message,'Wrote '+use_outfile,/continue
132 166  
133   -filename_final=!dustem_wrap_soft_dir+'/Data/EXAMPLE_OBSDATA/example_sed3.xcat'
134   -str='cp '+filename+' '+filename_final
135   -message,'Do '+str+' to make change permanent',/continue
  167 +;======== move the file to the Data/EXAMPLE_OBSDATA/ subdirectory
  168 +;filename_final=!dustem_wrap_soft_dir+'/Data/EXAMPLE_OBSDATA/my_example_SED_###.xcat'
  169 +;str='cp '+use_outfile+' '+filename_final
  170 +;message,'Do '+str+' to make change permanent',/continue
136 171  
137 172 the_end:
138 173  
139   -END
140 174 \ No newline at end of file
  175 +END
... ...