dustem_read_cb19_stellar_templates.pro
5.28 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
165
166
167
168
169
170
171
172
173
174
175
FUNCTION dustem_read_cb19_stellar_templates,age_values=age_values $
,metalicity_values=metalicity_values $
,template_wav=template_wav $
,Mstars=Mstars $
,info_only=info_only $
,info_st=info_st $
,help=help
;+
; NAME:
; dustem_read_cb19_stellar_templates
; PURPOSE:
; reads EMILES stellar population templates
; CATEGORY:
; DustEMWrap, Distributed, High-Level, User Example
; CALLING SEQUENCE:
; templates=dustem_read_cb19_stellar_templates([age_values=][,metalicity_values=][,template_waves=][,Mstars=][,/info_only][,info_st=][,/help])
; INPUTS:
; age_values : array of age values
; metalicity_values : array of metalicity values
; OPTIONAL INPUT PARAMETERS:
; None
; OUTPUTS:
; templates : SSP templates (F_lambda) [Lsun Msun^-1 A^-1]
; OPTIONAL OUTPUT PARAMETERS:
; template_wav = common wavelength array for templates [Angstroem]
; Mstars = Mass of stars [Msun]
; info_st = info structure for the templates
; ACCEPTED KEY-WORDS:
; info_only : if set, read only the information structure
; help : if set, prints this help
; 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:
;
; EXAMPLES
; st=dustem_read_cb19_stellar_templates(age_values=age_values,metalicity_values=metalicity_values,template_wav=template_wav)
; MODIFICATION HISTORY:
; Written by JPB June-2023
; 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_read_cb19_stellar_templates'
templates=-1
goto,the_end
END
;=== This is what we want
age_values_wanted=[0.03, 0.05, 0.08, 0.15, 0.25, 0.40, 0.60, 1.0, 1.75, 3.0, 5.0, 8.5, 13.5]
metalicity_values_wanted=[-1.48630, -0.961400, -0.351200, +0.0600000, +0.255900, +0.397100]
NagesWanted=n_elements(age_values_wanted)
NZWanted=n_elements(metalicity_values_wanted)
;=== This is what we have
Z=[0.0001, 0.0002, 0.0005, 0.001, 0.002, 0.004, 0.006, 0.008, 0.010, 0.014, 0.017, 0.020, 0.030, 0.040, 0.060]
;la métallicité solaire pour ce set étant: Z=0.01524, et donc ces modèles vont en gros de [Z/Zo]= -2.18 à +0.59
Z0=0.01524
Z=alog10(Z/Z0)
;print,Z
; -2.18299 -1.88196 -1.48401 -1.18298 -0.881955 -0.580925 -0.404834 -0.279895 -0.182985 -0.0368569 0.0474640 0.118045 0.294136 0.419075 0.595166
NZ=n_elements(Z)
dir_templates=!phangs_data_dir+'/ISRF/SSPs/CB2019/CB19_Chabrier/Mu100/'
;dir_templates=!dustem_wrap_soft_dir+'/Data/SSPs/CB2019/CB19_Chabrier/Mu100/'
files=file_search(dir_templates+'*.fits')
Nfiles=n_elements(files)
;metalicities=fltarr(Nfiles)
;ffile=strarr(Nfiles)
;FOR i=0L,Nfiles-1 DO BEGIN
; file=files[i]
; pos=strposmulti(file,'/')
; ffile[i]=strmid(file[i],pos[n_elements(pos)-1]+1,1000)
; pos=strposmulti(ffile[i],'_')
; metalicities[i]=float(strmid(ffile[i],pos[0]+2,pos[1]-pos[0]-2))
; stop
; ;metalicities[i]=
;ENDFOR
;read one fits file for header infos
file=dir_templates+'cb2019_z001_chab_hr_xmilesi_ssp.fits'
d=mrdfits(file,1,h)
wavs=d.wavelength
tagnames=tag_names(d)
ages=tagnames[1:*]
Nages=n_elements(ages)
ages=strmid(ages,1,1000)
aages=fltarr(Nages)
FOR i=0L,Nages-1 DO BEGIN
pos=strposmulti(ages[i],'P',count)
IF count NE 0 THEN BEGIN
aages[i]=float(textoidl_str_replace(ages[i],'P','.'))/1.e9 ;Gyr
ENDIF
ENDFOR
ages=aages
Nages=n_elements(ages)
all_ssps=ptrarr(NZ,Nages)
FOR i=0L,Nfiles-1 DO BEGIN
st=mrdfits(files[i],1,h)
FOR j=0L,Nages-1 DO BEGIN
all_ssps[i,j]=ptr_new(st.(j+1))
ENDFOR
ENDFOR
;interpolate on requested age an metalicities
Ntemplates=NagesWanted*NZWanted
templates=ptrarr(Ntemplates)
FOR i=0L,Ntemplates-1 DO BEGIN
ij=index2ij([i],[NagesWanted,NZWanted])
iage=ij[0,0]
imet=ij[0,1]
age=age_values_wanted[iage]
met=metalicity_values_wanted[imet]
;==== CAUTION, this is nearest neighbour interpolation
dist_age=abs(age-ages)
dist_Z=abs(met-Z)
ind_age=where(dist_age EQ min(dist_age))
ind_Z=where(dist_Z EQ min(dist_Z))
templates[i]=all_ssps[ind_Z[0],ind_age[0]]
;print,ind_age[0],ind_Z[0]
;==== This is nearest neighbour on age and interpolation on Z
ssps_used=all_ssps[*,ind_age[0]]
xx=indgen(Nz)
toto=interpol(xx,Z,met)
u1=abs(toto-round(toto))
u2=abs(toto-fix(toto))
;print,round(toto),toto,fix(toto)
;print,1.-u1,1.-u2
vec=(*ssps_used[round(toto)])*(1.-u1) + (*ssps_used[fix(toto)])*(1.-u2)
ind=where(vec LT 0,count)
IF count NE 0 THEN BEGIN
stop
ENDIF
templates[i]=ptr_new(vec)
;stop
ENDFOR
age_values=age_values_wanted
metalicity_values=metalicity_values_wanted
;print,age_values
;print,metalicity_values
;stop
angstroem2mic=1.e6/1.e10
template_wav=wavs*angstroem2mic ;mic
;do_interpol=1L
;IF do_interpol THEN BEGIN
;stop
lambir=dustem_get_wavelengths()
ind_extrapol=where(lambir GT max(template_wav) OR lambir LT min(template_wav),count_extrapol)
FOR i=0L,Ntemplates-1 DO BEGIN
;ind=where(lambir LE max(template_wav) AND lambir GE min(template_wav),count)
vec=interpol(*templates[i],template_wav,lambir)
vec[ind_extrapol]=0. ;set extrapolated values to 0 (otherwise would be slightly negative)
templates[i]=ptr_new(vec)
ENDFOR
template_wav=lambir
;ENDIF
the_end:
RETURN,templates
END