OV-GSO-DC / dustem-wrapper_idl

Commit 0d19c43e3eaaf1315708450fb55dbf681f6d2714

Authored by Jean-Philippe Bernard
1 parent c24ac88b
Exists in master

improved

Showing 4 changed files with 45 additions and 25 deletions   Show diff stats
Inline Side-by-side
LabTools/IRAP/JPB/phangs_brute_force_fit_with_isrf_grid.pro
  Diff comments   View file @0d19c43
... ... @@ -47,11 +47,13 @@ restore,data_dir+use_source_name+'_isrf_min_prediction.sav',/verb
47 47 ;% RESTORE: Restored variable: SOURCE_NAME.
48 48  
49 49 ;For test
50   -class_min=13
51   -class_max=13
  50 +;class_min=13
  51 +;class_max=13
  52 +;class_min=13
  53 +;class_max=13
52 54  
53   -;class_min=10
54   -;class_max=15
  55 +class_min=10
  56 +class_max=15
55 57  
56 58 restore,data_dir+use_source_name+'_jwst_seds_muse_pixels.sav',/verb
57 59 ;% RESTORE: Restored variable: ALL_SEDS.
... ... @@ -139,7 +141,7 @@ readcol,file,Mathis_wavs,Mathis_ISRF
139 141 Mathis_ISRF=interpol(Mathis_ISRF,Mathis_wavs,isrf_wavelengths)
140 142 Mathis_1mic=interpol(Mathis_ISRF,isrf_wavelengths,1.0)
141 143  
142   -stop
  144 +;stop
143 145  
144 146 FOR isrf_class=class_min,class_max DO BEGIN
145 147 ;==== This is to use the grid for this ISRF class
... ... @@ -241,7 +243,7 @@ ENDFOR
241 243 ;FOR vid=0LL,Nvor-1 DO Ypah_map[*all_seds_indices[vid]]=Ypahs[vid]
242 244 ;FOR vid=0LL,Nvor-1 DO YVSG_map[*all_seds_indices[vid]]=Yvsgs[vid]
243 245 ;FOR vid=0LL,Nvor-1 DO chi2_map[*all_seds_indices[vid]]=chi2s[vid]
244   -stop
  246 +;stop
245 247  
246 248 ;G0 histogram
247 249 win=0L
... ... @@ -259,27 +261,27 @@ window,win,xsize=800,ysize=900 & win=win+1
259 261  
260 262 obp=[1.1,0.,1.15,1]
261 263 ;image_cont20,G0_map,Href,/square,imrange=[-2,40],axis_color_table=1,image_color_table='jpbloadct',/silent,off_bar=obp,title='G0'
262   -image_cont20,la_log10(G0_map),Href,/square,imrange=[-0.2,2],axis_color_table=1,image_color_table='jpbloadct',/silent,off_bar=obp,title='G0'
  264 +image_cont20,la_log10(G0_map),Href,/square,imrange=[-0.2,2],axis_color_table=1,image_color_table='jpbloadct',/silent,off_bar=obp,title='log10(G0)'
263 265  
264 266 window,win,xsize=800,ysize=900 & win=win+1
265   -;imrange=[-2.,6]
266   -imrange=[-3.,0.]
  267 +imrange=[-2.,6]
  268 +;imrange=[-3.,0.]
267 269 image_cont20,la_log10(Ypah_map),Href,/square,imrange=imrange,axis_color_table=1,image_color_table='jpbloadct',/silent,off_bar=obp,title='log10(Ypah)'
268 270  
269 271 window,win,xsize=800,ysize=900 & win=win+1
270   -;imrange=[-2.,2]
271   -imrange=[-3.,0.]-3. ;why ??
  272 +imrange=[-2.,2]
  273 +;imrange=[-3.,0.]-3. ;why ??
272 274 ;image_cont20,la_log10(YVSG_map),Href,/square,imrange=[-2.,6],axis_color_table=1,image_color_table='jpbloadct',/silent,off_bar=obp,title='log10(Yvsg)'
273 275 image_cont20,la_log10(YVSG_map),Href,/square,imrange=imrange,axis_color_table=1,image_color_table='jpbloadct',/silent,off_bar=obp,title='log10(Yvsg)'
274 276  
275 277 window,win,xsize=800,ysize=900 & win=win+1
276 278 imrange=[3,7]
277   -imrange=[3,7]+6.
278   -image_cont20,la_log10(chi2_map),Href,/square,imrange=imrange,image_color_table='jpbloadct',/silent,title='Chi2'
  279 +;imrange=[3,7]+6.
  280 +image_cont20,la_log10(chi2_map),Href,/square,imrange=imrange,axis_color_table=1,image_color_table='jpbloadct',/silent,title='log10(Chi2)',off_bar=obp
279 281  
280 282 window,win,xsize=800,ysize=900 & win=win+1
281   -imrange=[0,3]
282   -image_cont20,la_log10(fact_map),Href,/square,imrange=imrange,image_color_table='jpbloadct',/silent,title='fact'
  283 +imrange=[-1,2]
  284 +image_cont20,la_log10(fact_map),Href,/square,imrange=imrange,axis_color_table=1,image_color_table='jpbloadct',/silent,title='log10(fact)',off_bar=obp
283 285  
284 286 stop
285 287 jwst_4_coutours=fltarr((size(jwst_images))[1],(size(jwst_images))[2],2)
... ...
LabTools/IRAP/JPB/srun/make_phangs_grids.pro
  Diff comments   View file @0d19c43
... ... @@ -21,7 +21,8 @@ bidon=1
21 21  
22 22 ;make grid
23 23 ;make_sed_phangs_tables,/test,grid_type=2,/show_seds,isrf_class=15
24   -make_sed_phangs_tables,grid_type=2,isrf_class=16
  24 +make_sed_phangs_tables,grid_type=2,isrf_class=17
  25 +;make_sed_phangs_tables,grid_type=2,isrf_class=16
25 26 ;make_sed_phangs_tables,grid_type=2,isrf_class=15
26 27 ;make_sed_phangs_tables,grid_type=2,isrf_class=14
27 28 ;make_sed_phangs_tables,grid_type=2,isrf_class=13
... ...
src/idl/dustem_brute_force_fit.pro
  Diff comments   View file @0d19c43
... ... @@ -97,9 +97,20 @@ table_pmins=!dustem_grid.pmin_values
97 97 table_pmaxs=!dustem_grid.pmax_values
98 98  
99 99 IF keyword_set(fixed_parameters_description) THEN BEGIN
100   - dustem_marginalize_grid,fixed_parameters_description,fixed_parameters_values,seds,Nseds,Nparams,table_params,table_pmins,table_pmaxs
  100 + dustem_marginalize_grid,fixed_parameters_description,fixed_parameters_values,seds,Nseds,params,Nparams,table_params,table_pmins,table_pmaxs
101 101 ENDIF
102 102  
  103 +;help,seds,Nseds,params,Nparams,table_params,table_pmins,table_pmaxs
  104 +;SEDS STRUCT = -> <Anonymous> Array[400]
  105 +;NSEDS LONG = 400
  106 +;PARAMS STRUCT = -> <Anonymous> Array[400]
  107 +;NPARAMS LONG = 3
  108 +;TABLE_PARAMS STRING = Array[3]
  109 +;TABLE_PMINS DOUBLE = Array[3]
  110 +;TABLE_PMAXS DOUBLE = Array[3]
  111 +
  112 +;stop
  113 +
103 114 ;=== These are the full arrays needed for calculations
104 115 chi2s=dblarr(Nseds)
105 116 grid_seds=dblarr(Nseds,Nfilters)
... ... @@ -143,12 +154,13 @@ ENDIF ELSE BEGIN
143 154 chi2s[i]=total((sed.stokesI-grid_seds[i,*])^2/sed.sigmaII)
144 155 ENDFOR
145 156 ENDELSE
  157 +
146 158 ;stop
147 159  
148 160 chi2=min(chi2s)
149 161 ind=where(chi2s EQ chi2,count)
150 162  
151   -params=grid_param_values[ind[0],*] ;these are the best fit parameters
  163 +params=reform(grid_param_values[ind[0],*]) ;these are the best fit parameters
152 164 fact=facts[ind[0]] ; This is the corresponding scaling factor
153 165 Nfreedom=Nfilters-Nparams-1 ;This is the degree of freedom
154 166 rchi2=chi2/Nfreedom ;This is the reduced chi2
... ... @@ -186,8 +198,8 @@ ind=where(chi2s LE chi2+dchi2,count)
186 198 params_min=fltarr(Nparams)
187 199 params_max=fltarr(Nparams)
188 200 FOR i=0L,Nparams-1 DO BEGIN
189   - params_min[i]=min(params[ind,i])
190   - params_max[i]=max(params[ind,i])
  201 + params_min[i]=min(grid_param_values[ind[0],i])
  202 + params_max[i]=max(grid_param_values[ind[0],i])
191 203 ENDFOR
192 204 params_uncertainties=params_max-params_min
193 205  
... ...
src/idl/dustem_marginalize_grid.pro
  Diff comments   View file @0d19c43
... ... @@ -2,6 +2,7 @@ PRO dustem_marginalize_grid,fixed_parameters_description $
2 2 ,fixed_parameters_values $
3 3 ,seds $
4 4 ,Nseds $
  5 + ,params $
5 6 ,Nparams $
6 7 ,table_params $
7 8 ,table_pmins $
... ... @@ -22,6 +23,7 @@ PRO dustem_marginalize_grid,fixed_parameters_description $
22 23 ; fixed_parameters_values : fixed parameters values
23 24 ; seds : sed grid to be marginalized (also an output)
24 25 ; Nseds : Number of seds (also an output)
  26 +; params : sed grid parameter values to be marginalized (also an output)
25 27 ; Nparams : Number of parameters in sed grid (also an output)
26 28 ; table_params : description of parameters in sed grid (also an output)
27 29 ; table_pmins : min values of parameters in sed grid (also an output)
... ... @@ -70,19 +72,22 @@ FOR i=0L,N_fixed_params-1 DO BEGIN
70 72 IF count NE 0 THEN BEGIN
71 73 fixed_parameter_value=fixed_parameters_values[ind[0]]
72 74 ;compute distance between table parameter and the requested fixed value
73   - distance=abs(fixed_parameter_value-seds[*].(ind[0]))
  75 + ;distance=abs(fixed_parameter_value-seds[*].(ind[0]))
  76 + distance=abs(fixed_parameter_value-params[*].(ind[0]))
74 77 indd=where(distance EQ min(distance),ccount)
  78 + ;stop
75 79 IF ccount EQ 0 THEN BEGIN
76 80 message,'This should not happen',/continue
77 81 stop
78 82 ENDIF
79 83 seds=seds[indd]
  84 + params=params[indd]
80 85 ;Nseds=n_elements(seds) ;or ccount
81 86 Nseds=ccount
82   - table_params=table_params[indc]
83   - Nparams=nc
84   - table_pmins=table_pmins[indc]
85   - table_pmaxs=table_pmaxs[indc]
  87 + ;table_params=table_params[indc]
  88 + ;Nparams=nc
  89 + ;table_pmins=table_pmins[indc]
  90 + ;table_pmaxs=table_pmaxs[indc]
86 91 ENDIF
87 92 ENDFOR
88 93  
... ...