make_phangs_ssps_isrf_prediction.pro
5.75 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
PRO make_phangs_ssps_isrf_prediction,source_name=source_name,save=save,help=help
;+
; NAME:
; make_phangs_ssps_isrf_prediction
; CALLING SEQUENCE:
; make_phangs_ssps_isrf_prediction[,source_name=][,/save][,/help]
; PURPOSE:
; predicts the minimum ISRF in voronoi bins
; INPUTS:
; None
; OPTIONAL KEYWORDS:
; save = if set, save the classes
; help = if set, print this help
; OUTPUTS:
; None
; OPTIONAL INPUT:
; source_name = source name (default='ngc0628')
; OPTIONAL OUTPUT:
; None
; PROCEDURE AND SUBROUTINE USED
; ISRF saved is the one predicted by dustem_plugin_phangs_stellar_isrf for each Muse Voronoi bin,
; scaled in amplitude to match the Muse filter data for that Voronoi bin.
; The G0 value computed is wrt the Mathis field at 1 mic.
; SIDE EFFECTS:
; produces file _isrf_min_prediction.sav
; EXAMPLE:
; make_phangs_ssps_isrf_prediction,source_name='ngc0628',/save
; MODIFICATION HISTORY:
; written by Jean-Philippe Bernard
;-
IF keyword_set(help) THEN BEGIN
doc_library,'make_phangs_ssps_isrf_prediction'
goto,the_end
ENDIF
;===== predicts the minimum ISRF in voronoi bins
;===== The prediction is scaled to the Muse filter data observations
win=0L
;window,win,xsize=900,ysize=1000 & win=win+1
;=== This is where the data is read from and the ISRFs will be stored
data_dir=!phangs_data_dir+'/ISRF/WORK/'
use_model='DBP90' ;Example with default keywords uses the DBP90 model
use_polarization=0 ; initialize Dustemwrap in no polarization mode
;== INITIALISE DUSTEM
dustem_init,model=use_model,polarization=use_polarization,show_plots=show_plots
use_source_name='ngc0628'
IF keyword_set(source_name) THEN use_source_name=source_name
CASE use_source_name OF
'ngc0628':BEGIN
object_distance=9.84 ;MPc
object_thickness=0.1 ;kpc
END
ENDCASE
dustem_define_la_common
obp=[1.1,0,1.15,1]
;dustem_init,show_plots=show_plots
;needed only for NHCO
restore,data_dir+use_source_name+'_jwst_images.sav',/verb
;% RESTORE: Restored variable: JWST_IMAGES.
;% RESTORE: Restored variable: FILTERS.
;% RESTORE: Restored variable: HREF.
;% RESTORE: Restored variable: NHCO.
;needed for stellar parameters and voronoi bin info
restore,data_dir+use_source_name+'_muse_images.sav',/verb
;% RESTORE: Restored variable: ST_TEMPLATES.
;% RESTORE: Restored variable: ST_MUSE_WEIGHTS.
;% RESTORE: Restored variable: VORONOI_ID.
;% RESTORE: Restored variable: AGE_VALUES.
;% RESTORE: Restored variable: METALICITY_VALUES.
;% RESTORE: Restored variable: BINS.
;% RESTORE: Restored variable: HREF.
;=== needed for absolute scaling of ISRF
restore,data_dir+use_source_name+'_muse_filters_data.sav',/verb
;% RESTORE: Restored variable: MUSE_IMAGES.
;% RESTORE: Restored variable: MUSE_FILTERS.
;% RESTORE: Restored variable: HREF.
restore,data_dir+use_source_name+'_muse_seds_muse_pixels.sav',/verb
;% RESTORE: Restored variable: ALL_SEDS.
;% RESTORE: Restored variable: ALL_SEDS_INDICES.
Nvor=max(voronoi_id)
use_NHmap=NHCO ;used NH map in 1.e21 (from CO)
;stop
!quiet=1
t1=systime(0,/sec)
first_vid=0L
lambir=dustem_get_wavelengths(isrf_wavelengths=isrf_wavelengths)
;show_each=10
show_each=100
Nlamb=n_elements(isrf_wavelengths)
ISRFS=dblarr(Nlamb,Nvor)
G0s=dblarr(Nvor)
;==== get Mathis field for G0 calculations
file=!dustem_soft_dir+'/data/ISRF_MATHIS.DAT'
readcol,file,Mathis_wavs,Mathis_ISRF
Mathis_ISRF=interpol(Mathis_ISRF,Mathis_wavs,isrf_wavelengths)
Mathis_1mic=interpol(Mathis_ISRF,isrf_wavelengths,1.0)
FOR vid=first_vid,Nvor-1 DO BEGIN
IF vid mod show_each EQ 0 THEN BEGIN
t2=systime(0,/sec)
delta_t_hr=(t2-t1)/60.^2 ;elapsed time [hr]
frac_perc=vid/Nvor*100
delta_t_remain=delta_t_hr/frac_perc*100.-delta_t_hr ;remaining time [hr]
message,'done '+strtrim(frac_perc)+' % in '+strtrim(delta_t_hr,2)+' hr remaining '+strtrim(delta_t_remain,2)+' hr',/continue
ENDIF
;===== get the SSP weights for the given voronoi bins
weights=phangs_binid2weights(st_muse_weights,vid,st_templates,age_values,metalicity_values,reddening=reddening)
amplitude=1. ;amplitude is set to 1, and will be corrected by G0 below
use_reddening=0. ;This sets Muse reddening to 0. Will be contsrained later
fpd=phangs_stellar_continuum_plugin_weight2params(weights,parameter_values=val,redenning=use_reddening,/force_include_reddening,amplitude=amplitude,/force_include_amplitude)
Nparams=n_elements(val)
key=intarr(Nparams)
FOR i=0L,Nparams-1 DO BEGIN
toto=dustem_parameter_description2type(fpd[i],string_name=string_name,key=one_key)
key[i]=one_key
ENDFOR
;==== compute the prediction stellar continuum for comparison to the Muse filter data.
;==== Here we use the reddening from MUSE
val_sed=val
val_sed[1]=reddening ;for the SED, we use E(B-V)
spectrum=dustem_plugin_phangs_stellar_continuum(key=key,val=val_sed)
sed=interpol(spectrum[*,0],lambir,dustem_filter2wav(muse_filters[0:2]))
muse_sed=(*all_seds[vid]).stokesI
facts=muse_sed/sed
muse_factor=la_mean(facts) ;This takes the average value over the MUSE bands
;print,G0,reddening
;stop
val_isrf=val ;for the ISRF, we use E(B-V)=0
ISRFS[*,vid]=dustem_plugin_phangs_stellar_isrf(key=key,val=val_isrf,object_distance=object_distance,object_thickness=object_thickness)
ISRFS[*,vid]=ISRFS[*,vid]*muse_factor ;=== scale the amplitude of ISRFs to the muse data
;=== compute corresponding G0
G0=interpol(ISRFS[*,vid],isrf_wavelengths,1.)/Mathis_1mic
G0s[vid]=G0
ENDFOR
IF keyword_set(save) THEN BEGIN
file_save=data_dir+use_source_name+'_isrf_min_prediction.sav'
save,ISRFS,G0s,object_distance,object_thickness,use_source_name,file=file_save,/verb
phangs_make_isrf_product,use_source_name,voronoi_id,ISRFS,G0s,object_distance,object_thickness
ENDIF
;cgplot,lindgen(Nvor),G0s,/ylog,yrange=[1.e-5,1],xtit='Voronoi ID',ytit='G0'
stop
the_end:
END