OV-GSO-DC / dustem-wrapper_idl

Commit cbdb57bc9e5248841fabc33e14a8d135c90357d0

Authored by Jean-Philippe Bernard
1 parent 46d50ffd
Exists in master

changed

Showing 9 changed files with 913 additions and 295 deletions   Show diff stats
Inline Side-by-side
LabTools/IRAP/JPB/calz_guess_ebv.pro
  Diff comments   View file @cbdb57b
1   -FUNCTION calz_guess_ebv,wave,flux_ratio,R_v=R_v
  1 +FUNCTION calz_guess_ebv,wave,flux_ratio,R_v=R_v,tau=tau
2 2  
3 3 ;+
4 4 ; NAME:
... ... @@ -56,6 +56,9 @@ IF (wave GE 912) AND (wave LT 6300) THEN $
56 56 ;ebv=alog10(funred/flux)/0.4/klam
57 57 ebv=alog10(flux_ratio)/0.4/klam
58 58  
  59 +;1.085736 is alog10(exp(1.d0))/0.4
  60 +tau=ebv*klam/1.085736
  61 +
59 62 RETURN,ebv
60 63  
61 64 END
62 65 \ No newline at end of file
... ...
LabTools/IRAP/JPB/phangs_brute_force_fit_with_isrf_grid.pro
  Diff comments   View file @cbdb57b
... ... @@ -315,6 +315,7 @@ restore,file,/verb
315 315 ;% RESTORE: Restored variable: ALL_SEDS.
316 316 astrosat_seds=all_seds
317 317 astrosat_filters=(*astrosat_seds[0]).filter
  318 +N_astrosat_filters=n_elements(astrosat_filters)
318 319 file=data_dir+use_source_name+'_herschel_seds_muse_pixels'+reso_str+'.sav'
319 320 st_info=file_info(file)
320 321 IF st_info.exists NE 1 THEN BEGIN
... ... @@ -409,6 +410,9 @@ ENDFOR
409 410 F0s=fltarr(Nvor)+la_undef()
410 411 G0s=fltarr(Nvor)+la_undef()
411 412 EBVs=fltarr(Nvor)+la_undef()
  413 +EBVs_uv=fltarr(Nvor,N_astrosat_filters)+la_undef()
  414 +tausV=fltarr(Nvor)+la_undef()
  415 +taus_uv=fltarr(Nvor,N_astrosat_filters)+la_undef()
412 416 observed_seds=ptrarr(Nvor) ;All observed SEDs
413 417 best_fit_seds=ptrarr(Nvor) ;All best fit SEDs
414 418 Ygrain1=fltarr(Nvor)+la_undef()
... ... @@ -515,33 +519,58 @@ NH_min=la_abs(la_min(NHmap[ind]))
515 519 ;stop
516 520  
517 521 ;=== stuff for stellar subtraction
518   -;starlight_filter_names=['F200W'] ;These are the filters that will be used to normalize stellar light (ie assumed to be dominated by stellar emission)
519   -starlight_filter_names=['F300M']
  522 +starlight_filter_names=['F200W'] ;These are the filters that will be used to normalize stellar light (ie assumed to be dominated by stellar emission)
  523 +;starlight_filter_names=['F300M']
520 524 N_starlight_filters=n_elements(starlight_filter_names)
521 525 use_starlight_filters=dustem_filter_names2filters(starlight_filter_names)
522 526 wav_star=dustem_filter2wav(use_starlight_filters)
523 527  
524 528 ;starlight_ebv_names=['F148W'] ;Astrosat1, used to compute extinction
525   -starlight_ebv_names=['sdss_g'] ;MUSE filter shortest wavelength, used to compute extinction
  529 +starlight_ebv_names=['sdss_g'] ;MUSE filter shortest wavelength, used to compute VIS extinction
526 530 N_ebv_filters=n_elements(starlight_ebv_names)
527 531 use_ebv_filters=dustem_filter_names2filters(starlight_ebv_names)
528 532 wav_ebv_star=dustem_filter2wav(use_ebv_filters)
529 533  
  534 +;starlight_ebv_uv_names=['F148W','F148Wa','F154W','F169M','F172M','N219M','N242W','N245M','N263M','N279N']
  535 +;starlight_ebv_uv_names=['F148W'] ;Astrosat filter used to compute UV extinction
  536 +;use_ebv_uv_filters=dustem_filter_names2filters(starlight_ebv_uv_names)
  537 +;N_ebv_uv_filters=n_elements(starlight_ebv_uv_names)
  538 +;=== This is to use all astrosat filters for extinction analysis
  539 +use_ebv_uv_filters=astrosat_filters
  540 +N_ebv_uv_filters=N_astrosat_filters
  541 +wav_ebv_uv_star=dustem_filter2wav(use_ebv_uv_filters)
  542 +EBV_uv=fltarr(N_ebv_uv_filters)
  543 +tau_uv=fltarr(N_ebv_uv_filters)
  544 +
  545 +;stop
  546 +
530 547 indstar=[-1L]
531   -indebv=[-1L]
532 548 FOR i=0L,N_starlight_filters-1 DO BEGIN
533 549 ind=where(all_filters EQ use_starlight_filters[i],count)
534 550 IF count NE 0 THEN BEGIN
535 551 indstar=[indstar,ind[0]]
536 552 ENDIF
  553 +ENDFOR
  554 +indstar=indstar[1:*]
  555 +
  556 +indebv=[-1L]
  557 +FOR i=0L,N_ebv_filters-1 DO BEGIN
537 558 ind=where(all_filters EQ use_ebv_filters[i],count)
538 559 IF count NE 0 THEN BEGIN
539 560 indebv=[indebv,ind[0]]
540 561 ENDIF
541 562 ENDFOR
542   -indstar=indstar[1:*]
543 563 indebv=indebv[1:*]
544 564  
  565 +indebv_uv=[-1L]
  566 +FOR i=0L,N_ebv_uv_filters-1 DO BEGIN
  567 + ind=where(all_filters EQ use_ebv_uv_filters[i],count)
  568 + IF count NE 0 THEN BEGIN
  569 + indebv_uv=[indebv_uv,ind[0]]
  570 + ENDIF
  571 +ENDFOR
  572 +indebv_uv=indebv_uv[1:*]
  573 +
545 574 one_st={isrf_class:0L, $
546 575 Nvor:0L, $
547 576 min_G0:0.,med_G0:0.,max_G0:0., $
... ... @@ -590,6 +619,10 @@ IF use_yvsg_option EQ 'MAP' THEN BEGIN
590 619 ENDIF
591 620 ;stop
592 621  
  622 +wwav_star=all_wavs[indebv[0]] ;mic
  623 +;wwav_uv=all_wavs[indebv_uv[0]] ;mic
  624 +wwav_uv=all_wavs[indebv_uv] ;mic
  625 +
593 626 FOR isrf_class=class_min,class_max DO BEGIN ;isrf_class runs from 1
594 627 ;print,isrf_class
595 628 ;stop
... ... @@ -691,11 +724,22 @@ FOR isrf_class=class_min,class_max DO BEGIN ;isrf_class runs from 1
691 724 sed_stars=interpol(ISRFs[*,vid],isrf_wavelengths,all_wavs)/fact_mjy2ergs/NH_value
692 725 sed_stars=la_mul(sed_stars,norm_stars) ;This is the stellar Full SED
693 726  
694   - wwav=all_wavs[indebv[0]]*1.e4 ;AA
  727 + ;wwav=all_wavs[indebv[0]]*1.e4 ;AA
695 728 ;flux_ratio=(sed_norm.stokesI/sed_stars)[indebv[0]] ;flux ratio between unredened and redenned (observed) flux
696 729 ;flux_ratio=sed_norm[indebv[0]].stokesI/sed_stars[indebv[0]] ;flux ratio between unredened and redenned (observed) flux
697 730 flux_ratio=sed_stars[indebv[0]]/sed_norm[indebv[0]].stokesI ;flux ratio between unredened and redenned (observed) flux
698   - EBV=calz_guess_ebv(wwav,flux_ratio,R_v=R_v)
  731 + EBV=calz_guess_ebv(wwav_star*1.e4,flux_ratio,R_v=R_v,tau=tau_calz)
  732 + tauV=tau_calz
  733 + ;stop
  734 + ;wwav_uv=all_wavs[indebv_uv[0]]*1.e4
  735 + ;flux_ratio_uv=sed_stars[indebv_uv[0]]/sed_norm[indebv_uv[0]].stokesI ;flux ratio between unredened and redenned (observed) flux
  736 + ;;stop
  737 + ;EBV_uv=calz_guess_ebv(wwav_uv*1.e4,flux_ratio_uv,R_v=R_v)
  738 + FOR k=0L,N_ebv_uv_filters-1 DO BEGIN
  739 + flux_ratio_uv=sed_stars[indebv_uv[k]]/sed_norm[indebv_uv[k]].stokesI ;flux ratio between unredened and redenned (observed) flux
  740 + EBV_uv[k]=calz_guess_ebv(wwav_uv[k]*1.e4,flux_ratio_uv,R_v=R_v,tau=tau_calz)
  741 + tau_uv[k]=tau_calz
  742 + ENDFOR
699 743 ;stop
700 744 sed_stars_reddened=sed_stars
701 745 ;We want to reden, so, use -1.*EBV
... ... @@ -793,6 +837,9 @@ FOR isrf_class=class_min,class_max DO BEGIN ;isrf_class runs from 1
793 837 G0s_predicted[vid]=la_mul(F0s_predicted[vid],this_G0prime)
794 838 IF keyword_set(subtract_star_light) THEN BEGIN
795 839 EBVs[vid]=EBV
  840 + EBVs_uv[vid,*]=EBV_uv
  841 + tausV[vid]=tauV
  842 + taus_uv[vid,*]=tau_UV
796 843 ENDIF
797 844 observed_seds[vid]=ptr_new(sed_norm.stokesI)
798 845 best_fit_seds[vid]=ptr_new(total_predicted_sed.stokesI)
... ... @@ -886,7 +933,7 @@ stats_st=st
886 933 IF keyword_set(save) THEN BEGIN
887 934 dustem_grid_saved=!dustem_grid
888 935 file_save=data_dir+use_source_name+use_model+'_logn_JWST_G0_YPAH_YVSG_on-voronoi_with-ISRFclasses'+norm_str+mathis_str+reso_str+fit_G0_str+set_yvsg_str+'.sav'
889   - save,all_filters,use_fit_filters,dustem_grid_saved,stellar_normalizations,stats_st,hit_st,hit_st_perc,EBVs,F0s,F0s_predicted,G0s,G0s_predicted, $
  936 + save,all_filters,use_fit_filters,dustem_grid_saved,stellar_normalizations,stats_st,hit_st,hit_st_perc,EBVs,EBVs_uv,tausV,taus_UV,F0s,F0s_predicted,G0s,G0s_predicted, $
890 937 Ygrain1,Ygrain2,Ygrain3,facts,chi2s,rchi2s,dF0s,dYgrain1,dYgrain2,dYgrain3, $
891 938 F0_hit,Ygrain1_hit,Ygrain2_hit,Ygrain3_hit,NH_values,best_fit_seds,observed_seds,file=file_save
892 939 message,'Wrote '+file_save,/continue
... ... @@ -895,6 +942,9 @@ ENDIF
895 942 ;stop
896 943 ;=== make maps of parameters
897 944 EBV_map=fltarr(sxpar(href,'NAXIS1'),sxpar(href,'NAXIS2'))+la_undef()
  945 +EBV_uv_map=fltarr(sxpar(href,'NAXIS1'),sxpar(href,'NAXIS2'),N_ebv_uv_filters)+la_undef()
  946 +tauV_map=fltarr(sxpar(href,'NAXIS1'),sxpar(href,'NAXIS2'))+la_undef()
  947 +tau_uv_map=fltarr(sxpar(href,'NAXIS1'),sxpar(href,'NAXIS2'),N_ebv_uv_filters)+la_undef()
898 948 chi2_map=fltarr(sxpar(href,'NAXIS1'),sxpar(href,'NAXIS2'))+la_undef()
899 949 rchi2_map=fltarr(sxpar(href,'NAXIS1'),sxpar(href,'NAXIS2'))+la_undef()
900 950 F0_map=fltarr(sxpar(href,'NAXIS1'),sxpar(href,'NAXIS2'))+la_undef()
... ... @@ -937,6 +987,11 @@ FOR vid=0LL,Nvor-1 DO BEGIN
937 987 ENDFOR
938 988 ENDIF
939 989 EBV_map[*all_seds_indices[vid]]=EBVs[vid]
  990 + tauV_map[*all_seds_indices[vid]]=tausV[vid]
  991 + FOR k=0L,N_ebv_uv_filters-1 DO BEGIN
  992 + EBV_uv_map[*all_seds_indices[vid],k]=EBVs_uv[vid,k]
  993 + tau_uv_map[*all_seds_indices[vid],k]=taus_uv[vid,k]
  994 + ENDFOR
940 995 F0_map[*all_seds_indices[vid]]=F0s[vid]
941 996 G0_map[*all_seds_indices[vid]]=G0s[vid]
942 997 F0_predicted_map[*all_seds_indices[vid]]=F0s_predicted[vid]
... ... @@ -971,12 +1026,16 @@ IF keyword_set(save) THEN BEGIN
971 1026 ; F0_hit,Ygrain1_hit,Ygrain2_hit,Ygrain3_hit,Yvsg_hit,NH_values,best_fit_seds,observed_seds,file=file_save
972 1027 ; message,'Wrote '+file_save,/continue
973 1028 file_save=data_dir+use_source_name+use_model+'_logn_JWST_G0_YPAH_YVSG_maps_with-ISRFclasses'+norm_str+mathis_str+reso_str+fit_G0_str+set_yvsg_str+'.sav'
974   - save,NHCO,NHI,NHmap,norm_stars_map,EBV_map,F0_map,G0_map,F0_predicted_map,G0_predicted_map,Yg1_map,Yg2_map,Yg3_map,fact_map,chi2_map,rchi2_map,dF0_map,dYg1_map,dYg2_map,dYg3_map,F0_hit_map,Yg1_hit_map,Yg2_hit_map,Yg3_hit_map,NH_values_map,residual_maps,residual_maps_perc,file=file_save
  1029 + save,NHCO,NHI,NHmap,norm_stars_map,EBV_map,EBV_uv_map,tauV_map,tau_uv_map,N_ebv_uv_filters,wwav_star,wwav_uv,F0_map,G0_map,F0_predicted_map,G0_predicted_map,Yg1_map,Yg2_map,Yg3_map,fact_map,chi2_map,rchi2_map,dF0_map,dYg1_map,dYg2_map,dYg3_map,F0_hit_map,Yg1_hit_map,Yg2_hit_map,Yg3_hit_map,NH_values_map,residual_maps,residual_maps_perc,file=file_save
975 1030 message,'Wrote '+file_save,/continue
976 1031 ENDIF
977 1032  
978 1033 from_restore:
979 1034  
  1035 +message,'DONE. The /from_restore option is obsolete. Use phangs_make_figures_from_restore.pro to show maps and figures ...',/continue
  1036 +stop
  1037 +goto,the_end
  1038 +
980 1039 ;Below could be replaced by a call to
981 1040 ; phangs_make_figures_from_restore,source_name=source_name $
982 1041 ; ,model=model $
... ... @@ -1590,7 +1649,7 @@ IF keyword_set(do_postscripts) THEN BEGIN
1590 1649 postscript=ps_file,orientation='portrait',/nologo
1591 1650 message,'Wrote '+ps_file,/continue
1592 1651 ENDIF
1593   -stop
  1652 +;stop
1594 1653  
1595 1654 window,win,xsize=800,ysize=900 & win=win+1
1596 1655 imrange=[-1.,+1.0]
... ...
LabTools/IRAP/JPB/phangs_extract_seds.pro
  Diff comments   View file @cbdb57b
... ... @@ -2,17 +2,18 @@ PRO phangs_extract_seds,source_name=source_name $
2 2 ,resolution_filter=resolution_filter $
3 3 ,resolution_value=resolution_value $
4 4 ,reset_seds=reset_seds $
  5 + ,do_all_seds=do_all_seds $
5 6 ,help=help
6 7  
7 8 ;+
8 9 ; NAME:
9 10 ; phangs_extract_seds
10 11 ; PURPOSE:
11   -; extract SEDs in Muse Voronoi bins
  12 +; extract SEDs in Muse Voronoi bins in the JWST,MUSE,ASTROSAT,HERSCHEL data for a given phangs object
12 13 ; CATEGORY:
13 14 ; PHANGS ISRF
14 15 ; CALLING SEQUENCE:
15   -; phangs_extract_seds[,source_name=][,/save][,/reset_seds]
  16 +; phangs_extract_seds[,source_name=][,/save][,/reset_seds][,/do_all_seds]
16 17 ; INPUTS:
17 18 ; None
18 19 ; OPTIONAL INPUT PARAMETERS:
... ... @@ -25,6 +26,7 @@ PRO phangs_extract_seds,source_name=source_name $
25 26 ; help = If set, print this help
26 27 ; nostop = if set, does not stop
27 28 ; reset_seds = if set, removes all SED saved files from disk and recomputes all.
  29 +; do_all_seds= if set, will also produce the seds aggregated with all data at once (JWST,MUSE,ASTROSAT,HERSCHEL)
28 30 ; COMMON BLOCKS:
29 31 ; None
30 32 ; SIDE EFFECTS:
... ... @@ -34,7 +36,7 @@ PRO phangs_extract_seds,source_name=source_name $
34 36 ; _muse_seds_muse_pixels.sav
35 37 ; _astrosat_seds_muse_pixels.sav
36 38 ; _herschel_seds_muse_pixels.sav
37   -; _all_seds_muse_pixels.sav
  39 +; _all_seds_muse_pixels.sav (only if do_all_seds=1)
38 40 ; RESTRICTIONS:
39 41 ; None
40 42 ; PROCEDURE:
... ... @@ -66,6 +68,7 @@ ENDIF
66 68  
67 69 data_dir=!phangs_data_dir+'/ISRF/WORK/'
68 70  
  71 +;==== If needed, remove existing SED files
69 72 IF keyword_set(reset_seds) THEN BEGIN
70 73 files_2b_removed=data_dir+use_source_name+'*seds_muse_pixels*'+reso_str+'.sav'
71 74 str='rm '+files_2b_removed
... ... @@ -81,6 +84,7 @@ ENDIF
81 84  
82 85 extract_seds:
83 86  
  87 +;==== read the MUSE data
84 88 file=data_dir+use_source_name+'_muse_data'+reso_str+'.sav'
85 89 st_info=file_info(file)
86 90 IF st_info.exists NE 1 THEN BEGIN
... ... @@ -95,6 +99,8 @@ restore,file,/verb
95 99 ;% RESTORE: Restored variable: METALICITY_VALUES.
96 100 ;% RESTORE: Restored variable: BINS.
97 101 ;% RESTORE: Restored variable: HREF.
  102 +
  103 +;==== read the JWST data
98 104 file=data_dir+use_source_name+'_jwst_images'+reso_str+'.sav'
99 105 st_info=file_info(file)
100 106 IF st_info.exists NE 1 THEN BEGIN
... ... @@ -106,6 +112,8 @@ restore,file,/verb
106 112 ;% RESTORE: Restored variable: FILTERS.
107 113 ;% RESTORE: Restored variable: HREF.
108 114 ;% RESTORE: Restored variable: NHCO.
  115 +
  116 +;==== read the ASTROSAT data
109 117 file=data_dir+use_source_name+'_astrosat_images'+reso_str+'.sav'
110 118 st_info=file_info(file)
111 119 IF st_info.exists NE 1 THEN BEGIN
... ... @@ -116,6 +124,8 @@ restore,file,/verb ;changed to astrosat_images
116 124 ;% RESTORE: Restored variable: ASTROSAT_IMAGES.
117 125 ;% RESTORE: Restored variable: ASTROSAT_FILTERS.
118 126 ;% RESTORE: Restored variable: HREF.
  127 +
  128 +;==== read the MUSE filter data
119 129 file=data_dir+use_source_name+'_muse_filters_images'+reso_str+'.sav'
120 130 st_info=file_info(file)
121 131 IF st_info.exists NE 1 THEN BEGIN
... ... @@ -126,6 +136,8 @@ restore,file,/verb ;changed to muse_images
126 136 ;% RESTORE: Restored variable: MUSE_IMAGES.
127 137 ;% RESTORE: Restored variable: MUSE_FILTERS.
128 138 ;% RESTORE: Restored variable: HREF.
  139 +
  140 +;==== read the HERSCHEL filter data
129 141 file=data_dir+use_source_name+'_herschel_images'+reso_str+'.sav'
130 142 st_info=file_info(file)
131 143 IF st_info.exists NE 1 THEN BEGIN
... ... @@ -148,6 +160,7 @@ file_save_indices=data_dir+use_source_name+'_seds_indices'+reso_str+'.sav'
148 160 file_save=data_dir+use_source_name+'_jwst_seds_muse_pixels'+reso_str+'.sav'
149 161 res=file_info(file_save)
150 162  
  163 +;==== If the JWST SED file does not already exist, do it. This will also save all_seds_indices needed for other extractions
151 164 IF not res.exists THEN BEGIN
152 165 ;stop
153 166 all_seds=extract_all_muse_phangs_seds(use_source_name,voronoi_id,jwst_images,filters,indices=all_seds_indices,counts=counts)
... ... @@ -160,12 +173,11 @@ ENDIF ELSE BEGIN
160 173 stop
161 174 ENDELSE
162 175  
163   -;stop
164   -
  176 +;restore the all_seds_indices variable for further SED extractions
165 177 restore,file_save_indices,/verb ;just to get the all_seds_indices variable
166 178 ;% RESTORE: Restored variable: ALL_SEDS_INDICES
167 179  
168   -;==== extract Muse SEDs, using previously computed all_seds_indices
  180 +;==== If the MUSE SED file does not already exist, do it. Extract MUSE SEDs, using previously computed all_seds_indices
169 181 file_save=data_dir+use_source_name+'_muse_seds_muse_pixels'+reso_str+'.sav'
170 182 res=file_info(file_save)
171 183 IF not res.exists THEN BEGIN
... ... @@ -177,6 +189,7 @@ ENDIF ELSE BEGIN
177 189 stop
178 190 ENDELSE
179 191  
  192 +;==== If the ASTROSAT SED file does not already exist, do it. Extract ASTROSAT SEDs, using previously computed all_seds_indices
180 193 file_save=data_dir+use_source_name+'_astrosat_seds_muse_pixels'+reso_str+'.sav'
181 194 res=file_info(file_save)
182 195 IF not res.exists THEN BEGIN
... ... @@ -188,6 +201,7 @@ ENDIF ELSE BEGIN
188 201 stop
189 202 ENDELSE
190 203  
  204 +;==== If the HERSCHEL SED file does not already exist, do it. Extract HERSCHEL SEDs, using previously computed all_seds_indices
191 205 file_save=data_dir+use_source_name+'_herschel_seds_muse_pixels'+reso_str+'.sav'
192 206 res=file_info(file_save)
193 207 IF not res.exists THEN BEGIN
... ... @@ -201,29 +215,30 @@ ENDELSE
201 215  
202 216 ;==== The following is not really useful, can be done through SED aggregation later
203 217 ;==== However, not equivalent to the above, because of undefined pixels in some of the images, and the way this is handled by sed extractor
204   -;stop
205   -
206   -file_save=data_dir+use_source_name+'_all_seds_muse_pixels'+reso_str+'.sav'
207   -res=file_info(file_save)
208   -IF not res.exists THEN BEGIN
209   - use_filters=[astrosat_filters,muse_filters,filters,herschel_filters]
210   - use_Nfilters=n_elements(use_filters)
211   - use_images=fltarr([sxpar(href,'NAXIS1'),sxpar(href,'NAXIS2'),2,use_Nfilters])
212   - from=0L & to=from+n_elements(filters)-1
213   - use_images[*,*,*,from:to]=jwst_images
214   - from=to+1 & to=from+n_elements(muse_filters)-1
215   - use_images[*,*,*,from:to]=muse_images
216   - from=to+1 & to=from+n_elements(astrosat_filters)-1
217   - use_images[*,*,*,from:to]=astrosat_images
218   - from=to+1 & to=from+n_elements(herschel_filters)-1
219   - use_images[*,*,*,from:to]=herschel_images
220   - all_seds=extract_all_muse_phangs_seds(use_source_name,voronoi_id,use_images,use_filters,indices=all_seds_indices,counts=counts)
221   - save,all_seds,file=file_save
222   - message,'Saved '+file_save,/info
223   -ENDIF ELSE BEGIN
224   - message,'File '+file_save+' already exist',/continue
225   - stop
226   -ENDELSE
  218 +;==== so ....
  219 +IF keyword_set(do_all_seds) THEN BEGIN
  220 + file_save=data_dir+use_source_name+'_all_seds_muse_pixels'+reso_str+'.sav'
  221 + res=file_info(file_save)
  222 + IF not res.exists THEN BEGIN
  223 + use_filters=[astrosat_filters,muse_filters,filters,herschel_filters]
  224 + use_Nfilters=n_elements(use_filters)
  225 + use_images=fltarr([sxpar(href,'NAXIS1'),sxpar(href,'NAXIS2'),2,use_Nfilters])
  226 + from=0L & to=from+n_elements(filters)-1
  227 + use_images[*,*,*,from:to]=jwst_images
  228 + from=to+1 & to=from+n_elements(muse_filters)-1
  229 + use_images[*,*,*,from:to]=muse_images
  230 + from=to+1 & to=from+n_elements(astrosat_filters)-1
  231 + use_images[*,*,*,from:to]=astrosat_images
  232 + from=to+1 & to=from+n_elements(herschel_filters)-1
  233 + use_images[*,*,*,from:to]=herschel_images
  234 + all_seds=extract_all_muse_phangs_seds(use_source_name,voronoi_id,use_images,use_filters,indices=all_seds_indices,counts=counts)
  235 + save,all_seds,file=file_save
  236 + message,'Saved '+file_save,/info
  237 + ENDIF ELSE BEGIN
  238 + message,'File '+file_save+' already exist',/continue
  239 + stop
  240 + ENDELSE
  241 +ENDIF
227 242  
228 243 the_end:
229 244  
... ...
LabTools/IRAP/JPB/phangs_isrf_pipeline.pro
  Diff comments   View file @cbdb57b
... ... @@ -151,6 +151,8 @@ phangs_smooth_muse_isrf,'ngc0628',resolution_filter='PACS3',/save,/show,/nostop
151 151  
152 152 make_phangs_isrf_classes,bidon,source_name='ngc0628',/save,resolution_filter='PACS3'
153 153  
  154 +phangs_make_intensity_mask,source_name='ngc0628',resolution_filter='PACS3',/show,/save
  155 +
154 156 ;===test new tables
155 157 dustem_init
156 158  
... ... @@ -160,35 +162,41 @@ class_min_max=[15,20]
160 162 phangs_brute_force_fit_with_isrf_grid,model='DL07',source_name='ngc0628',/include_herschel,/normalize,/fit_G0,resolution_filter='PACS3' $
161 163 ,/subtract_star_light,/show,/nostop,table_name_root=table_name_root,class_min_max=class_min_max,/save
162 164  
163   -;ISRF classes
  165 +;fit map using MUSE derived ISRF
164 166 table_name_root='_logn_MuseISRF_JWST_adaptedG0_Ypah1added_Yvsgadded_4Phangs_noionis_F0upext' ;Uses extended F0 range ....
165 167 class_min_max=0
166 168 show_sed=0
167 169 phangs_brute_force_fit_with_isrf_grid,model='DL07',source_name='ngc0628',/include_herschel,/normalize,/fit_G0,/save,resolution_filter='PACS3' $
168 170 ,/subtract_star_light,show_sed=show_sed,/nostop,table_name_root=table_name_root,class_min_max=class_min_max
169 171  
170   -;idem from restore and do postscripts
171   -table_name_root='_logn_MuseISRF_JWST_adaptedG0_Ypah1added_Yvsgadded_4Phangs_noionis_F0upext' ;Uses extended F0 range ....
172   -class_min_max=0
  172 +;doing postscripts
  173 +table_name_root='_logn_MuseISRF_JWST_adaptedG0_Ypah1added_Yvsgadded_4Phangs_noionis_F0upext'
173 174 postscripts_dir=!phangs_data_dir+'/ISRF/WORK/PS/ngc0628_PACS3/'
174   -phangs_brute_force_fit_with_isrf_grid,model='DL07',source_name='ngc0628',/include_herschel,/normalize,/fit_G0,/save,resolution_filter='PACS3' $
175   - ,/subtract_star_light,/nostop,table_name_root=table_name_root,class_min_max=class_min_max $
176   - ,/from_restore,/do_postscripts,postscripts_dir=postscripts_dir
  175 +phangs_make_figures_from_restore,model='DL07',source_name='ngc0628',/normalize,fit_G0=1,resolution_filter='PACS3' $
  176 + ,/nostop,table_name_root=table_name_root $
  177 + ,/do_postscripts,postscripts_dir=postscripts_dir
  178 +
  179 +
  180 +
  181 +
  182 +
  183 +
  184 +
177 185  
178 186 ;Mathis ISRF
179   -table_name_root='_logn_MuseISRF_JWST_adaptedG0_Ypah1added_Yvsgadded_4Phangs_noionis' ;rem: They also have YPAH0 added ....
  187 +table_name_root='_logn_MuseISRF_JWST_adaptedG0_Ypah1added_Yvsgadded_4Phangs_noionis_F0upext' ;rem: They also have YPAH0 added ....
180 188 class_min_max=0
181 189 show_sed=0
182 190 phangs_brute_force_fit_with_isrf_grid,model='DL07',source_name='ngc0628',/include_herschel,/normalize,/fit_G0,/save,resolution_filter='PACS3' $
183 191 ,/subtract_star_light,show_sed=show_sed,/nostop,table_name_root=table_name_root,class_min_max=class_min_max,/force_mathis
184 192  
185 193 ;idem from restore and do postscripts
186   -table_name_root='_logn_MuseISRF_JWST_adaptedG0_Ypah1added_Yvsgadded_4Phangs_noionis' ;rem: They also have YPAH0 added ....
  194 +table_name_root='_logn_MuseISRF_JWST_adaptedG0_Ypah1added_Yvsgadded_4Phangs_noionis_F0upext' ;rem: They also have YPAH0 added ....
187 195 class_min_max=0
188 196 postscripts_dir=!phangs_data_dir+'/ISRF/WORK/PS/ngc0628_PACS3_Mathis/'
189   -phangs_brute_force_fit_with_isrf_grid,model='DL07',source_name='ngc0628',/include_herschel,/normalize,/fit_G0,/save,resolution_filter='PACS3' $
190   - ,/subtract_star_light,/show,/nostop,table_name_root=table_name_root,class_min_max=class_min_max $
191   - ,/from_restore,/do_postscripts,/force_mathis,postscripts_dir=postscripts_dir
  197 +phangs_make_figures_from_restore,model='DL07',source_name='ngc0628',/normalize,fit_G0=1,resolution_filter='PACS3' $
  198 + ,/nostop,table_name_root=table_name_root $
  199 + ,/do_postscripts,postscripts_dir=postscripts_dir,/force_mathis
192 200  
193 201 ;This is option A1 (using default value):
194 202 ;This is when not trying to fit G0, ie, assuming same average relation between predicted F0 and F0
... ... @@ -595,45 +603,16 @@ phangs_brute_force_fit_with_isrf_grid,model='DL07',source_name='ngc0628',/includ
595 603 ,/subtract_star_light,show_sed=show_sed,/nostop,table_name_root=table_name_root,class_min_max=class_min_max $
596 604 ,yvsg_option='MAP',yvsg_fits_map=Yg3_map,yvsg_fits_header=href,marginalize_method=marginalize_method,/silent
597 605  
598   -;same as above, doing postscripts
599   -table_name_root='_logn_MuseISRF_JWST_adaptedG0_Ypah1added_Yvsgadded_4Phangs_noionis_F0upext'
600   -postscripts_dir=!phangs_data_dir+'/ISRF/WORK/PS/ngc0628_1arcsec_B1_NEAREST'
601   -phangs_brute_force_fit_with_isrf_grid,model='DL07',source_name='ngc0628',/normalize,fit_G0=1,set_yvsg=1,resolution_value=1./60.^2 $
602   - ,/nostop,table_name_root=table_name_root $
603   - ,/from_restore,/do_postscripts,postscripts_dir=postscripts_dir
604   -
605   -;same as above, doing postscripts, using hangs_make_figures_from_restore
  606 +;doing postscripts, using phangs_make_figures_from_restore
606 607 table_name_root='_logn_MuseISRF_JWST_adaptedG0_Ypah1added_Yvsgadded_4Phangs_noionis_F0upext'
607 608 postscripts_dir=!phangs_data_dir+'/ISRF/WORK/PS/ngc0628_1arcsec_B1_NEAREST'
608 609 phangs_make_figures_from_restore,model='DL07',source_name='ngc0628',/normalize,fit_G0=1,set_yvsg=1,resolution_value=1./60.^2 $
609 610 ,/nostop,table_name_root=table_name_root $
610   - ,/do_postscripts,postscripts_dir=postscripts_dir
  611 + ,/do_postscripts,postscripts_dir=postscripts_dir,/no_MC
611 612  
612   -;======== test
613   -;ngc0628_1arcsec_B1_NEAREST (1arcsec_B1_NEAREST)
614   -table_name_root='_logn_MuseISRF_JWST_adaptedG0_Ypah1added_Yvsgadded_4Phangs_noionis_F0upext' ;rem: They also have YPAH0 added ....
615   -;table_name_root='_logn_MuseISRF_JWST_adaptedG0_Ypah1added_Yvsgadded_4Phangs_noionis' ;rem: They also have YPAH0 added ....
616   -;class_min_max=[19,20]
617   -data_dir=!phangs_data_dir+'/ISRF/WORK/'
618   -source_name='ngc0628'
619   -model='DL07'
620   -norm_str='_normalized'
621   -mathis_str=''
622   -reso_str='_PACS3'
623   -fit_G0_str=''
624   -set_yvsg_str=''
625   -file_save=data_dir+source_name+model+'_JWST_G0_YPAH_YVSG_maps_with-ISRFclasses'+norm_str+mathis_str+reso_str+fit_G0_str+set_yvsg_str+'.sav'
626   -restore,file_save,/verb
627   -restore,data_dir+source_name+'_ref_header'+reso_str+'.sav',/verb
628   -nostop=1
629   -show_sed=0
630   -marginalize_method='NEAREST'
631   -phangs_brute_force_fit_with_isrf_grid,model='DL07',source_name='ngc0628',/include_herschel,/normalize,fit_G0=1,set_yvsg=1,/save,resolution_value=1./60.^2 $
632   - ,/subtract_star_light,show_sed=show_sed,/nostop,table_name_root=table_name_root,class_min_max=class_min_max $
633   - ,yvsg_option='MAP',yvsg_fits_map=Yg3_map,yvsg_fits_header=href,marginalize_method=marginalize_method,/silent
634 613  
635   -;ngc0628_1arcsec_B1M_NEAREST (1arcsec_B1M_NEAREST)
636   -table_name_root='_logn_MuseISRF_JWST_adaptedG0_Ypah1added_Yvsgadded_4Phangs_noionis'
  614 +;ngc0628_1arcsec_B1M_NEAREST (1arcsec_B1M_NEAREST) Mathis
  615 +table_name_root='_logn_MuseISRF_JWST_adaptedG0_Ypah1added_Yvsgadded_4Phangs_noionis_F0upext'
637 616 class_min_max=0
638 617 data_dir=!phangs_data_dir+'/ISRF/WORK/'
639 618 source_name='ngc0628'
... ... @@ -656,16 +635,16 @@ phangs_brute_force_fit_with_isrf_grid,model='DL07',source_name='ngc0628',/includ
656 635 ,yvsg_option='MAP',yvsg_fits_map=Yg3_map,yvsg_fits_header=href,marginalize_method=marginalize_method,force_mathis=force_mathis $
657 636 ,silent=silent
658 637  
659   -;same as above, doing postscripts
660   -table_name_root='_logn_MuseISRF_JWST_adaptedG0_Ypah1added_Yvsgadded_4Phangs_noionis'
  638 +;doing postscripts
  639 +table_name_root='_logn_MuseISRF_JWST_adaptedG0_Ypah1added_Yvsgadded_4Phangs_noionis_F0upext'
661 640 postscripts_dir=!phangs_data_dir+'/ISRF/WORK/PS/ngc0628_1arcsec_B1M_NEAREST'
662   -force_mathis=1
663   -phangs_brute_force_fit_with_isrf_grid,model='DL07',source_name='ngc0628',/normalize,fit_G0=1,set_yvsg=1,resolution_value=1./60.^2 $
  641 +phangs_make_figures_from_restore,model='DL07',source_name='ngc0628',/normalize,fit_G0=1,set_yvsg=1,resolution_value=1./60.^2,/force_mathis $
664 642 ,/nostop,table_name_root=table_name_root $
665   - ,/from_restore,/do_postscripts,postscripts_dir=postscripts_dir,force_mathis=force_mathis
  643 + ,/do_postscripts,postscripts_dir=postscripts_dir
666 644  
667 645 ;J'en suis la
668 646  
  647 +
669 648 ;ngc0628_1arcsec_B1M_NEAREST (1arcsec_B1M_NEAREST)
670 649 table_name_root='_MuseISRF_JWST_adaptedG0_Ypah1added_Ypah2added_adaptedYvsgadded_4Phangs_noionis'
671 650 class_min_max=0
... ... @@ -1050,17 +1029,20 @@ phangs_smooth_muse_isrf,'ngc4321',resolution_filter='PACS3',/save,/show,/nostop
1050 1029  
1051 1030 make_phangs_isrf_classes,bidon,source_name='ngc4321',/save,resolution_filter='PACS3',/do_postscripts
1052 1031  
  1032 +phangs_make_intensity_mask,source_name='ngc4321',resolution_filter='PACS3',/show,/save
  1033 +
1053 1034 ;brute-force fit with ISRF classes
1054   -table_name_root='_logn_MuseISRF_JWST_adaptedG0_Ypah1added_Yvsgadded_4Phangs_noionis' ;rem: They also have YPAH0 added ....
  1035 +table_name_root='_logn_MuseISRF_JWST_adaptedG0_Ypah1added_Yvsgadded_4Phangs_noionis_F0upext' ;rem: They also have YPAH0 added ....
1055 1036 class_min_max=0
1056 1037 show_sed=0
1057 1038 phangs_brute_force_fit_with_isrf_grid,model='DL07',source_name='ngc4321',/include_herschel,/normalize,/fit_G0,/save,resolution_filter='PACS3' $
1058 1039 ,/subtract_star_light,show_sed=show_sed,/nostop,table_name_root=table_name_root,class_min_max=class_min_max
1059 1040 dustem_init
1060 1041 postscripts_dir=!phangs_data_dir+'/ISRF/WORK/PS/ngc4321_PACS3/'
1061   -phangs_brute_force_fit_with_isrf_grid,model='DL07',source_name='ngc4321',/include_herschel,/normalize,/fit_G0,resolution_filter='PACS3' $
1062   - ,/subtract_star_light,show_sed=show_sed,/nostop,class_min_max=class_min_max $
1063   - ,/from_restore,/do_postscripts,postscripts_dir=postscripts_dir
  1042 +phangs_make_figures_from_restore,model='DL07',source_name='ngc4321',/normalize,fit_G0=1,resolution_filter='PACS3' $
  1043 + ,/nostop,table_name_root=table_name_root $
  1044 + ,/do_postscripts,postscripts_dir=postscripts_dir,/no_HII,/no_MC
  1045 +
1064 1046 ;J'en SUIS LA
1065 1047  
1066 1048 ;============================================================================================================
... ... @@ -1087,6 +1069,8 @@ phangs_smooth_muse_isrf,'ngc4321',resolution_value=1./60.^2,/save,/show,/nostop
1087 1069  
1088 1070 make_phangs_isrf_classes,bidon,source_name='ngc4321',/save,resolution_value=1./60.^2,/do_postscripts
1089 1071  
  1072 +phangs_make_intensity_mask,source_name='ngc4321',resolution_value=1./60.^2,/show,/save
  1073 +
1090 1074 ;ngc4321_1arcsec_B1_NEAREST (1arcsec_B1_NEAREST)
1091 1075 table_name_root='_logn_MuseISRF_JWST_adaptedG0_Ypah1added_Yvsgadded_4Phangs_noionis' ;rem: They also have YPAH0 added ....
1092 1076 class_min_max=0
... ... @@ -1108,12 +1092,13 @@ phangs_brute_force_fit_with_isrf_grid,model='DL07',source_name='ngc4321',/includ
1108 1092 ,/subtract_star_light,show_sed=show_sed,/nostop,table_name_root=table_name_root,class_min_max=class_min_max $
1109 1093 ,yvsg_option='MAP',yvsg_fits_map=Yg3_map,yvsg_fits_header=href,marginalize_method=marginalize_method,/silent
1110 1094  
1111   -
1112 1095 table_name_root='_logn_MuseISRF_JWST_adaptedG0_Ypah1added_Yvsgadded_4Phangs_noionis'
1113 1096 postscripts_dir=!phangs_data_dir+'/ISRF/WORK/PS/ngc4321_1arcsec_B1_NEAREST'
1114   -phangs_brute_force_fit_with_isrf_grid,model='DL07',source_name='ngc4321',/normalize,fit_G0=1,set_yvsg=1,resolution_value=1./60.^2 $
  1097 +phangs_make_figures_from_restore,model='DL07',source_name='ngc4321',/normalize,fit_G0=1,set_yvsg=1,resolution_value=1./60.^2 $
1115 1098 ,/nostop,table_name_root=table_name_root $
1116   - ,/from_restore,/do_postscripts,postscripts_dir=postscripts_dir
  1099 + ,/do_postscripts,postscripts_dir=postscripts_dir,/no_HII,/no_MC
  1100 +
  1101 +;J'en SUIS LA
1117 1102  
1118 1103 ;========
1119 1104 ;== NGC1566
... ... @@ -1251,30 +1236,215 @@ phangs_brute_force_fit_with_isrf_grid,model='DL07',source_name='ngc1566',/normal
1251 1236 ;========
1252 1237 ;== NGC3627
1253 1238 ;========
1254   -;J'en SUIS LA
1255   -;There is no data in /data/projects/phangs/phangs_drive/Archive/PHANGS_JWST/DR1p1/ngc3627 on alma1
1256 1239 phangs_make_jwst_images,source_name='ngc3627',/save,/show,/nostop
1257 1240 ;produces _muse_filters_data.sav,_muse_images.sav, _muse_NH.sav
1258   -phangs_make_muse_filters_images,source_name='ngc3351',/save,/show,/nostop
  1241 +phangs_make_muse_filters_images,source_name='ngc3627',/save,/show,/nostop
1259 1242 ;produces _astrosat_data.sav
1260   -phangs_make_astrosat_images,source_name='ngc3351',/save,/show,/nostop
  1243 +phangs_make_astrosat_images,source_name='ngc3627',/save,/show,/nostop
1261 1244 ;produces _herschel_images.sav
1262   -phangs_make_herschel_images,source_name='ngc3351',/save,/show,/nostop
1263   -phangs_make_hi_images,source_name='ngc3351',/save,/show,/nostop
1264   -phangs_make_co_images,source_name='ngc3351',/save,/show,/nostop
  1245 +phangs_make_herschel_images,source_name='ngc3627',/save,/show,/nostop
  1246 +phangs_make_hi_images,source_name='ngc3627',/save,/show,/nostop
  1247 +phangs_make_co_images,source_name='ngc3627',/save,/show,/nostop
1265 1248  
1266   -phangs_make_muse_filters_data,source_name='ngc3351',/save,/show,/nostop
  1249 +phangs_make_muse_filters_data,source_name='ngc3627',/save,/show,/nostop
1267 1250  
1268   -phangs_extract_seds,source_name='ngc3351'
  1251 +;phangs_extract_seds,source_name='ngc3627'
1269 1252  
1270 1253 ;weirdly here, the muse factors are not close to 1, when they are for NGC0628 ... why ??
1271   -make_phangs_ssps_isrf_prediction,source_name='ngc3351',/save
  1254 +make_phangs_ssps_isrf_prediction,source_name='ngc3627',/save
1272 1255  
1273   -make_phangs_isrf_classes,source_name='ngc3351',/save,/do_postscripts
  1256 +make_phangs_isrf_classes,source_name='ngc3627',/save,/do_postscripts
  1257 +;J'en SUIS LA
1274 1258  
1275   -table_name_root='_logn_MuseISRF_JWST_adaptedG0_Ypah1added_Yvsgadded_4Phangs_noionis'
  1259 +table_name_root='_logn_MuseISRF_JWST_adaptedG0_Ypah1added_Yvsgadded_4Phangs_noionis_F0upext'
1276 1260 show_sed=0
1277 1261 phangs_brute_force_fit_with_isrf_grid,model='DL07',table_name_root=table_name_root,source_name='ngc3351',/include_herschel,/normalize,/fit_G0,show_sed=show_sed,/subtract_star_light,/nostop
1278 1262  
  1263 +;=========================
  1264 +;== NGC3627 @ PACS3 resolution
  1265 +;========================
  1266 +phangs_make_jwst_images,source_name='ngc3627',/save,/show,/nostop,resolution_filter='PACS3'
  1267 +phangs_make_hi_images,source_name='ngc3627',/save,/show,/nostop,resolution_filter='PACS3'
  1268 +phangs_make_co_images,source_name='ngc3627',/save,/show,/nostop,resolution_filter='PACS3'
  1269 +phangs_make_muse_filters_images,source_name='ngc3627',/save,/show,/nostop,resolution_filter='PACS3'
  1270 +phangs_make_astrosat_images,source_name='ngc3627',/save,/show,/nostop,resolution_filter='PACS3'
  1271 +phangs_make_herschel_images,source_name='ngc3627',/save,/show,/nostop,resolution_filter='PACS3'
  1272 +phangs_make_muse_filters_data,source_name='ngc3627',/save,/show,/nostop,resolution_filter='PACS3'
  1273 +
  1274 +phangs_extract_seds,source_name='ngc3627',resolution_filter='PACS3'
  1275 +
  1276 +phangs_smooth_muse_isrf,'ngc3627',resolution_filter='PACS3',/save,/show,/nostop
  1277 +
  1278 +make_phangs_isrf_classes,bidon,source_name='ngc3627',/save,resolution_filter='PACS3',/do_postscripts
  1279 +
  1280 +phangs_make_intensity_mask,source_name='ngc3627',resolution_filter='PACS3',/show,/save
  1281 +
  1282 +;brute-force fit with ISRF classes
  1283 +table_name_root='_logn_MuseISRF_JWST_adaptedG0_Ypah1added_Yvsgadded_4Phangs_noionis_F0upext' ;rem: They also have YPAH0 added ....
  1284 +class_min_max=0
  1285 +show_sed=0
  1286 +phangs_brute_force_fit_with_isrf_grid,model='DL07',source_name='ngc3627',/include_herschel,/normalize,/fit_G0,/save,resolution_filter='PACS3' $
  1287 + ,/subtract_star_light,show_sed=show_sed,/nostop,table_name_root=table_name_root,class_min_max=class_min_max
  1288 +dustem_init
  1289 +postscripts_dir=!phangs_data_dir+'/ISRF/WORK/PS/ngc3627_PACS3/'
  1290 +phangs_make_figures_from_restore,model='DL07',source_name='ngc3627',/normalize,fit_G0=1,resolution_filter='PACS3' $
  1291 + ,/nostop,table_name_root=table_name_root $
  1292 + ,/do_postscripts,postscripts_dir=postscripts_dir,/no_HII,/no_MC
  1293 +;J'en SUIS LA
  1294 +
  1295 +;========
  1296 +;== NGC3627 @ 1"
  1297 +;========
  1298 +phangs_make_jwst_images,source_name='ngc3627',/save,/show,/nostop,resolution_value=1./60.^2
  1299 +phangs_make_hi_images,source_name='ngc3627',/save,/show,/nostop,resolution_value=1./60.^2
  1300 +phangs_make_co_images,source_name='ngc3627',/save,/show,/nostop,resolution_value=1./60.^2
  1301 +phangs_make_muse_filters_images,source_name='ngc3627',/save,/show,/nostop,resolution_value=1./60.^2
  1302 +phangs_make_astrosat_images,source_name='ngc3627',/save,/show,/nostop,resolution_value=1./60.^2
  1303 +phangs_make_herschel_images,source_name='ngc3627',/save,/show,/nostop,resolution_value=1./60.^2
  1304 +phangs_make_muse_filters_data,source_name='ngc3627',/save,/show,/nostop,resolution_value=1./60.^2
  1305 +
  1306 +phangs_extract_seds,source_name='ngc3627',resolution_value=1./60.^2
  1307 +
  1308 +phangs_smooth_muse_isrf,'ngc3627',resolution_value=1./60.^2,/save,/show,/nostop
  1309 +
  1310 +make_phangs_isrf_classes,bidon,source_name='ngc3627',/save,resolution_value=1./60.^2,/do_postscripts
  1311 +
  1312 +phangs_make_intensity_mask,source_name='ngc3627',resolution_value=1./60.^2,/show,/save
  1313 +
  1314 +;ngc3627_1arcsec_B1_NEAREST (1arcsec_B1_NEAREST)
  1315 +table_name_root='_logn_MuseISRF_JWST_adaptedG0_Ypah1added_Yvsgadded_4Phangs_noionis_F0upext' ;rem: They also have YPAH0 added ....
  1316 +class_min_max=0
  1317 +data_dir=!phangs_data_dir+'/ISRF/WORK/'
  1318 +source_name='ngc3627'
  1319 +model='DL07'
  1320 +norm_str='_normalized'
  1321 +mathis_str=''
  1322 +reso_str='_PACS3'
  1323 +fit_G0_str=''
  1324 +set_yvsg_str=''
  1325 +file_save=data_dir+source_name+model+'_logn_JWST_G0_YPAH_YVSG_maps_with-ISRFclasses'+norm_str+mathis_str+reso_str+fit_G0_str+set_yvsg_str+'.sav'
  1326 +restore,file_save,/verb
  1327 +restore,data_dir+source_name+'_ref_header'+reso_str+'.sav',/verb
  1328 +nostop=1
  1329 +show_sed=0
  1330 +marginalize_method='NEAREST'
  1331 +phangs_brute_force_fit_with_isrf_grid,model='DL07',source_name='ngc3627',/include_herschel,/normalize,fit_G0=1,set_yvsg=1,/save,resolution_value=1./60.^2 $
  1332 + ,/subtract_star_light,show_sed=show_sed,/nostop,table_name_root=table_name_root,class_min_max=class_min_max $
  1333 + ,yvsg_option='MAP',yvsg_fits_map=Yg3_map,yvsg_fits_header=href,marginalize_method=marginalize_method,/silent
  1334 +
  1335 +table_name_root='_logn_MuseISRF_JWST_adaptedG0_Ypah1added_Yvsgadded_4Phangs_noionis_F0upext'
  1336 +postscripts_dir=!phangs_data_dir+'/ISRF/WORK/PS/ngc3627_1arcsec_B1_NEAREST'
  1337 +phangs_make_figures_from_restore,model='DL07',source_name='ngc3627',/normalize,fit_G0=1,set_yvsg=1,resolution_value=1./60.^2 $
  1338 + ,/nostop,table_name_root=table_name_root $
  1339 + ,/do_postscripts,postscripts_dir=postscripts_dir,/no_HII,/no_MC
  1340 +;J'en SUIS LA
  1341 +
  1342 +;========
  1343 +;== NGC4254
  1344 +;========
  1345 +phangs_make_jwst_images,source_name='ngc4254',/save,/show,/nostop
  1346 +phangs_make_muse_filters_images,source_name='ngc4254',/save,/show,/nostop
  1347 +;produces _astrosat_data.sav
  1348 +phangs_make_astrosat_images,source_name='ngc4254',/save,/show,/nostop
  1349 +;produces _herschel_images.sav
  1350 +phangs_make_herschel_images,source_name='ngc4254',/save,/show,/nostop
  1351 +phangs_make_hi_images,source_name='ngc4254',/save,/show,/nostop
  1352 +phangs_make_co_images,source_name='ngc4254',/save,/show,/nostop
  1353 +
  1354 +phangs_make_muse_filters_data,source_name='ngc4254',/save,/show,/nostop
  1355 +
  1356 +phangs_extract_seds,source_name='ngc4254'
  1357 +
  1358 +;weirdly here, the muse factors are not close to 1, when they are for NGC0628 ... why ??
  1359 +make_phangs_ssps_isrf_prediction,source_name='ngc4254',/save
  1360 +
  1361 +make_phangs_isrf_classes,source_name='ngc4254',/save,/do_postscripts
  1362 +;J'en SUIS LA
  1363 +
  1364 +table_name_root='_logn_MuseISRF_JWST_adaptedG0_Ypah1added_Yvsgadded_4Phangs_noionis_F0upext'
  1365 +show_sed=0
  1366 +phangs_make_figures_from_restore,model='DL07',source_name='ngc4254',/normalize,fit_G0=1,set_yvsg=1,resolution_value=1./60.^2 $
  1367 + ,/nostop,table_name_root=table_name_root $
  1368 + ,/do_postscripts,postscripts_dir=postscripts_dir,/no_HII,/no_MC
  1369 +
  1370 +;=========================
  1371 +;== NGC4254 @ PACS3 resolution
  1372 +;========================
  1373 +phangs_make_jwst_images,source_name='ngc4254',/save,/show,/nostop,resolution_filter='PACS3'
  1374 +phangs_make_hi_images,source_name='ngc4254',/save,/show,/nostop,resolution_filter='PACS3'
  1375 +phangs_make_co_images,source_name='ngc4254',/save,/show,/nostop,resolution_filter='PACS3'
  1376 +phangs_make_muse_filters_images,source_name='ngc4254',/save,/show,/nostop,resolution_filter='PACS3'
  1377 +phangs_make_astrosat_images,source_name='ngc4254',/save,/show,/nostop,resolution_filter='PACS3'
  1378 +phangs_make_herschel_images,source_name='ngc4254',/save,/show,/nostop,resolution_filter='PACS3'
  1379 +phangs_make_muse_filters_data,source_name='ngc4254',/save,/show,/nostop,resolution_filter='PACS3'
  1380 +
  1381 +phangs_extract_seds,source_name='ngc4254',resolution_filter='PACS3'
  1382 +
  1383 +phangs_smooth_muse_isrf,'ngc4254',resolution_filter='PACS3',/save,/show,/nostop
  1384 +
  1385 +make_phangs_isrf_classes,bidon,source_name='ngc4254',/save,resolution_filter='PACS3',/do_postscripts
  1386 +
  1387 +phangs_make_intensity_mask,source_name='ngc4254',resolution_filter='PACS3',/show,/save
  1388 +
  1389 +;brute-force fit with ISRF classes
  1390 +table_name_root='_logn_MuseISRF_JWST_adaptedG0_Ypah1added_Yvsgadded_4Phangs_noionis_F0upext' ;rem: They also have YPAH0 added ....
  1391 +class_min_max=0
  1392 +show_sed=0
  1393 +phangs_brute_force_fit_with_isrf_grid,model='DL07',source_name='ngc4254',/include_herschel,/normalize,/fit_G0,/save,resolution_filter='PACS3' $
  1394 + ,/subtract_star_light,show_sed=show_sed,/nostop,table_name_root=table_name_root,class_min_max=class_min_max
  1395 +dustem_init
  1396 +postscripts_dir=!phangs_data_dir+'/ISRF/WORK/PS/ngc4254_PACS3/'
  1397 +phangs_make_figures_from_restore,model='DL07',source_name='ngc4254',/normalize,fit_G0=1,resolution_filter='PACS3' $
  1398 + ,/nostop,table_name_root=table_name_root $
  1399 + ,/do_postscripts,postscripts_dir=postscripts_dir,/no_HII,/no_MC
  1400 +;J'en SUIS LA
  1401 +
  1402 +;========
  1403 +;== NGC4254 @ 1"
  1404 +;========
  1405 +
  1406 +phangs_make_jwst_images,source_name='ngc4254',/save,/show,/nostop,resolution_value=1./60.^2
  1407 +phangs_make_hi_images,source_name='ngc4254',/save,/show,/nostop,resolution_value=1./60.^2
  1408 +phangs_make_co_images,source_name='ngc4254',/save,/show,/nostop,resolution_value=1./60.^2
  1409 +phangs_make_muse_filters_images,source_name='ngc4254',/save,/show,/nostop,resolution_value=1./60.^2
  1410 +phangs_make_astrosat_images,source_name='ngc4254',/save,/show,/nostop,resolution_value=1./60.^2
  1411 +phangs_make_herschel_images,source_name='ngc4254',/save,/show,/nostop,resolution_value=1./60.^2
  1412 +phangs_make_muse_filters_data,source_name='ngc4254',/save,/show,/nostop,resolution_value=1./60.^2
  1413 +
  1414 +phangs_extract_seds,source_name='ngc4254',resolution_value=1./60.^2
  1415 +
  1416 +phangs_smooth_muse_isrf,'ngc4254',resolution_value=1./60.^2,/save,/show,/nostop
  1417 +
  1418 +make_phangs_isrf_classes,bidon,source_name='ngc4254',/save,resolution_value=1./60.^2,/do_postscripts
  1419 +
  1420 +phangs_make_intensity_mask,source_name='ngc4254',resolution_value=1./60.^2,/show,/save
  1421 +
  1422 +;ngc4254_1arcsec_B1_NEAREST (1arcsec_B1_NEAREST)
  1423 +table_name_root='_logn_MuseISRF_JWST_adaptedG0_Ypah1added_Yvsgadded_4Phangs_noionis_F0upext' ;rem: They also have YPAH0 added ....
  1424 +class_min_max=0
  1425 +data_dir=!phangs_data_dir+'/ISRF/WORK/'
  1426 +source_name='ngc4254'
  1427 +model='DL07'
  1428 +norm_str='_normalized'
  1429 +mathis_str=''
  1430 +reso_str='_PACS3'
  1431 +fit_G0_str=''
  1432 +set_yvsg_str=''
  1433 +file_save=data_dir+source_name+model+'_logn_JWST_G0_YPAH_YVSG_maps_with-ISRFclasses'+norm_str+mathis_str+reso_str+fit_G0_str+set_yvsg_str+'.sav'
  1434 +restore,file_save,/verb
  1435 +restore,data_dir+source_name+'_ref_header'+reso_str+'.sav',/verb
  1436 +nostop=1
  1437 +show_sed=0
  1438 +marginalize_method='NEAREST'
  1439 +phangs_brute_force_fit_with_isrf_grid,model='DL07',source_name='ngc4254',/include_herschel,/normalize,fit_G0=1,set_yvsg=1,/save,resolution_value=1./60.^2 $
  1440 + ,/subtract_star_light,show_sed=show_sed,/nostop,table_name_root=table_name_root,class_min_max=class_min_max $
  1441 + ,yvsg_option='MAP',yvsg_fits_map=Yg3_map,yvsg_fits_header=href,marginalize_method=marginalize_method,/silent
  1442 +
  1443 +table_name_root='_logn_MuseISRF_JWST_adaptedG0_Ypah1added_Yvsgadded_4Phangs_noionis_F0upext'
  1444 +postscripts_dir=!phangs_data_dir+'/ISRF/WORK/PS/ngc4254_1arcsec_B1_NEAREST'
  1445 +phangs_make_figures_from_restore,model='DL07',source_name='ngc4254',/normalize,fit_G0=1,set_yvsg=1,resolution_value=1./60.^2 $
  1446 + ,/nostop,table_name_root=table_name_root $
  1447 + ,/do_postscripts,postscripts_dir=postscripts_dir,/no_HII,/no_MC
  1448 +;J'en SUIS LA
1279 1449  
1280 1450 END
... ...
LabTools/IRAP/JPB/phangs_make_figures_from_restore.pro
  Diff comments   View file @cbdb57b
... ... @@ -12,6 +12,8 @@ PRO phangs_make_figures_from_restore,source_name=source_name $
12 12 ,marginalize_method=marginalize_method $
13 13 ,do_postscripts=do_postscripts $
14 14 ,postscripts_dir=postscripts_dir $
  15 + ,no_HII=no_HII $
  16 + ,no_MC=no_MC $
15 17 ,silent=silent
16 18  
17 19 obp=[1.1,0.,1.15,1]
... ... @@ -99,20 +101,37 @@ IF st_info.exists NE 1 THEN BEGIN
99 101 stop
100 102 ENDIF
101 103 restore,file,/verb
  104 +;stop
  105 +; % RESTORE: Restored variable: ALL_FILTERS.
  106 +; % RESTORE: Restored variable: USE_FIT_FILTERS.
  107 +; % RESTORE: Restored variable: DUSTEM_GRID_SAVED.
  108 +; % RESTORE: Restored variable: STELLAR_NORMALIZATIONS.
102 109 ; % RESTORE: Restored variable: STATS_ST.
  110 +; % RESTORE: Restored variable: HIT_ST.
  111 +; % RESTORE: Restored variable: HIT_ST_PERC.
  112 +; % RESTORE: Restored variable: EBVS.
  113 +; % RESTORE: Restored variable: EBVS_UV.
103 114 ; % RESTORE: Restored variable: F0S.
104 115 ; % RESTORE: Restored variable: F0S_PREDICTED.
105   -; % RESTORE: Restored variable: YPAHS.
106   -; % RESTORE: Restored variable: YVSGS.
  116 +; % RESTORE: Restored variable: G0S.
  117 +; % RESTORE: Restored variable: G0S_PREDICTED.
  118 +; % RESTORE: Restored variable: YGRAIN1.
  119 +; % RESTORE: Restored variable: YGRAIN2.
  120 +; % RESTORE: Restored variable: YGRAIN3.
107 121 ; % RESTORE: Restored variable: FACTS.
108 122 ; % RESTORE: Restored variable: CHI2S.
109 123 ; % RESTORE: Restored variable: RCHI2S.
110 124 ; % RESTORE: Restored variable: DF0S.
111   -; % RESTORE: Restored variable: DYPAHS.
112   -; % RESTORE: Restored variable: DYVSGS.
  125 +; % RESTORE: Restored variable: DYGRAIN1.
  126 +; % RESTORE: Restored variable: DYGRAIN2.
  127 +; % RESTORE: Restored variable: DYGRAIN3.
113 128 ; % RESTORE: Restored variable: F0_HIT.
114   -; % RESTORE: Restored variable: YPAH_HIT.
115   -; % RESTORE: Restored variable: YVSG_HIT.
  129 +; % RESTORE: Restored variable: YGRAIN1_HIT.
  130 +; % RESTORE: Restored variable: YGRAIN2_HIT.
  131 +; % RESTORE: Restored variable: YGRAIN3_HIT.
  132 +; % RESTORE: Restored variable: NH_VALUES.
  133 +; % RESTORE: Restored variable: BEST_FIT_SEDS.
  134 +; % RESTORE: Restored variable: OBSERVED_SEDS.
116 135 file=data_dir+use_source_name+use_model+'_logn_JWST_G0_YPAH_YVSG_maps_with-ISRFclasses'+norm_str+mathis_str+reso_str+fit_G0_str+set_yvsg_str+'.sav'
117 136 st_info=file_info(file)
118 137 IF st_info.exists NE 1 THEN BEGIN
... ... @@ -120,11 +139,18 @@ IF st_info.exists NE 1 THEN BEGIN
120 139 stop
121 140 ENDIF
122 141 restore,file,/verb
  142 +;stop
123 143 ; % RESTORE: Restored variable: NHCO.
124 144 ; % RESTORE: Restored variable: NHI.
125 145 ; % RESTORE: Restored variable: NHMAP.
126 146 ; % RESTORE: Restored variable: NORM_STARS_MAP.
127 147 ; % RESTORE: Restored variable: EBV_MAP.
  148 +; % RESTORE: Restored variable: EBV_UV_MAP.
  149 +; % RESTORE: Restored variable: TAUV_MAP.
  150 +; % RESTORE: Restored variable: TAU_UV_MAP.
  151 +; % RESTORE: Restored variable: N_EBV_UV_FILTERS.
  152 +; % RESTORE: Restored variable: WWAV_STAR.
  153 +; % RESTORE: Restored variable: WWAV_UV.
128 154 ; % RESTORE: Restored variable: F0_MAP.
129 155 ; % RESTORE: Restored variable: G0_MAP.
130 156 ; % RESTORE: Restored variable: F0_PREDICTED_MAP.
... ... @@ -153,7 +179,8 @@ IF st_info.exists NE 1 THEN BEGIN
153 179 stop
154 180 ENDIF
155 181 restore,file,/verb
156   -;% RESTORE: Restored variable: VOR_CLASS.
  182 +; % RESTORE: Restored variable: VOR_CLASS.
  183 +; % RESTORE: Restored variable: G0PRIME.
157 184 file=data_dir+use_source_name+'_ref_header'+reso_str+'.sav'
158 185 st_info=file_info(file)
159 186 IF st_info.exists NE 1 THEN BEGIN
... ... @@ -166,6 +193,11 @@ restore,file,/verb
166 193 IF defined(dustem_grid_saved) THEN defsysv,'!dustem_grid',dustem_grid_saved
167 194 ;=== get the stellar mask
168 195 file=data_dir+use_source_name+'_star_mask'+reso_str+'.sav'
  196 +st_info=file_info(file)
  197 +IF st_info.exists NE 1 THEN BEGIN
  198 + message,'Could not find '+file,/continue
  199 + stop
  200 +ENDIF
169 201 restore,file,/verb
170 202 ;% RESTORE: Restored variable: STAR_MASK.
171 203 ;% RESTORE: Restored variable: HSTARMASK.
... ... @@ -175,6 +207,11 @@ restore,file,/verb
175 207  
176 208 ;==== restore mask of low brightness regions
177 209 file=data_dir+use_source_name+'_brightness_mask'+reso_str+'.sav'
  210 +st_info=file_info(file)
  211 +IF st_info.exists NE 1 THEN BEGIN
  212 + message,'Could not find '+file,/continue
  213 + stop
  214 +ENDIF
178 215 restore,file,/verb
179 216 ;% RESTORE: Restored variable: BRIGHTNESS_MASK.
180 217 ;stop
... ... @@ -196,6 +233,9 @@ IF count_masked NE 0 THEN BEGIN
196 233 NHmap[ind_masked]=la_undef()
197 234 norm_stars_map[ind_masked]=la_undef()
198 235 EBV_map[ind_masked]=la_undef()
  236 + EBV_uv_map[ind_masked,*]=la_undef()
  237 + tauV_map[ind_masked]=la_undef()
  238 + tau_uv_map[ind_masked,*]=la_undef()
199 239 F0_map[ind_masked]=la_undef()
200 240 G0_map[ind_masked]=la_undef()
201 241 F0_predicted_map[ind_masked]=la_undef()
... ... @@ -217,6 +257,7 @@ IF count_masked NE 0 THEN BEGIN
217 257 NH_values_map[ind_masked]=la_undef()
218 258 ENDIF
219 259  
  260 +Xions=la_div(Ygrain2,la_add(Ygrain1,Ygrain2))
220 261 Xion_map=la_div(Yg2_map,la_add(Yg1_map,Yg2_map))
221 262  
222 263 cleanplot
... ... @@ -332,174 +373,266 @@ cont_st.(0).thick=2.
332 373 cont_st.(1).color=100
333 374 cont_st.(1).thick=1.0
334 375  
  376 +;===== get molecular clouds mask
  377 +IF not keyword_set(no_MC) THEN BEGIN
  378 + dir=!phangs_data_dir+'/phangs_drive/Archive/Products/cloud_catalogs/team_release/by_galaxy/'+use_source_name+'/'
  379 + file=dir+use_source_name+'_12m+7m+tp_co21_native_binmask2D.fits.bz2'
  380 + mcs=mrdfits(file,0,hmcs)
  381 + image_cont20,mcs,hmcs,/square,/silent
  382 +ENDIF
  383 +
335 384 ;===== get HII region mask
336   -obj_name=strupcase(source_name)
337   -IF obj_name EQ 'NGC0628' THEN obj_name='NGC628'
338   -dir=!phangs_data_dir+'/phangs_drive/Archive/Products/Nebulae_catalogs/Nebulae_catalogue_v1/spatial_masks/'
339   -file=dir+obj_name+'_HIIreg_mask.fits'
340   -hii=readfits(file,hhii)
341   -ind=where(finite(hii) NE 1,count)
342   -IF count NE 0 THEN hii[ind]=la_undef()
343   -hhii=cd2astro_header(hhii)
344   -sxaddpar,hhii,'EQUINOX',2000.
345   -image_cont20,hii,hhii,/square,/silent
346   -hii=project2(hhii,hii,href,/silent,/round)
347   -window,win,xsize=wxs,ysize=wys & win=win+1
348   -image_cont20,hii,href,/square,/silent,image_color_table='jpbloadct,/white',axis_color_table=1,off_bar=obp
349   -;=== do a simple HII region mask (1 if an HII, 0 otherwise) and a contour structure
350   -is_hii=hii
351   -ind=where(hii NE -1 and hii NE la_undef(),count)
352   -is_hii[ind]=1
353   -is_hii[ind_bad]=la_undef()
354   -;stop
355   -hii_m=fltarr(Nx,Ny,2)
356   -hii_m[*,*,0]=is_hii
357   -hii_m[*,*,1]=is_hii
358   -levs_HII=[0.5,1]
359   -hii_st=create_contour_st(hii_m,lev1=levs_HII,lev2=levs_HII)
360   -hii_st.(0).color=0
361   -hii_st.(0).thick=2.
362   -hii_st.(1).color=100
363   -hii_st.(1).thick=1.0
364   -;hii[*]=la_undef()
365   -image_cont20,hii,href,/square,/silent,levels=hii_st,image_color_table='jpbloadct,/no',axis_color_table=1,off_bar=obp
366   -;image_cont20,is_hii,href,/square,/silent,levels=hii_st,image_color_table='jpbloadct',axis_color_table=1
367   -;stop
368   -;=== do a mask per HII region and a surounding mask
369   -N_hii=long(max(hii))
370   -hii_indices=ptrarr(N_hii)
371   -hii_count=lonarr(N_hii)
372   -hii_numbers=lonarr(N_hii)
373   -FOR i=0L,N_hii-1 DO BEGIN
374   - ind_hii=where(hii EQ i,count)
  385 +IF not keyword_set(no_HII) THEN BEGIN
  386 + obj_name=strupcase(use_source_name)
  387 + CASE obj_name OF
  388 + 'NGC0628': obj_name='NGC628'
  389 + ENDCASE
  390 + dir=!phangs_data_dir+'/phangs_drive/Archive/Products/Nebulae_catalogs/Nebulae_catalogue_v1/spatial_masks/'
  391 + file=dir+obj_name+'_HIIreg_mask.fits'
  392 + hii=readfits(file,hhii)
  393 + ind=where(finite(hii) NE 1,count)
  394 + IF count NE 0 THEN hii[ind]=la_undef()
  395 + hhii=cd2astro_header(hhii)
  396 + sxaddpar,hhii,'EQUINOX',2000.
  397 + image_cont20,hii,hhii,/square,/silent
  398 + hii=project2(hhii,hii,href,/silent,/round)
  399 + window,win,xsize=wxs,ysize=wys & win=win+1
  400 + image_cont20,hii,href,/square,/silent,image_color_table='jpbloadct,/white',axis_color_table=1,off_bar=obp
  401 + ;=== do a simple HII region mask (1 if an HII, 0 otherwise) and a contour structure
  402 + is_hii=hii
  403 + ind=where(hii NE -1 and hii NE la_undef(),count)
  404 + is_hii[ind]=1
  405 + is_hii[ind_bad]=la_undef()
  406 + ;stop
  407 + hii_m=fltarr(Nx,Ny,2)
  408 + hii_m[*,*,0]=is_hii
  409 + hii_m[*,*,1]=is_hii
  410 + levs_HII=[0.5,1]
  411 + hii_st=create_contour_st(hii_m,lev1=levs_HII,lev2=levs_HII)
  412 + hii_st.(0).color=0
  413 + hii_st.(0).thick=2.
  414 + hii_st.(1).color=100
  415 + hii_st.(1).thick=1.0
  416 + ;hii[*]=la_undef()
  417 + image_cont20,hii,href,/square,/silent,levels=hii_st,image_color_table='jpbloadct,/no',axis_color_table=1,off_bar=obp
  418 + ;image_cont20,is_hii,href,/square,/silent,levels=hii_st,image_color_table='jpbloadct',axis_color_table=1
  419 + ;stop
  420 + ;=== do a mask per HII region and a surounding mask
  421 + N_hii=long(max(hii))
  422 + hii_indices=ptrarr(N_hii)
  423 + hii_count=lonarr(N_hii)
  424 + hii_numbers=lonarr(N_hii)
  425 + FOR i=0L,N_hii-1 DO BEGIN
  426 + ind_hii=where(hii EQ i,count)
  427 + IF count NE 0 THEN BEGIN
  428 + hii_indices[i]=ptr_new(ind_hii)
  429 + ENDIF
  430 + hii_count[i]=count
  431 + hii_numbers[i]=i+1
  432 + ENDFOR
  433 + ;=== remove unfound HII regions
  434 + ind=where(hii_count NE 0,count)
375 435 IF count NE 0 THEN BEGIN
376   - hii_indices[i]=ptr_new(ind_hii)
  436 + hii_count=hii_count[ind]
  437 + hii_indices=hii_indices[ind]
  438 + hii_numbers=hii_numbers[ind] ;Original HII region number
  439 + N_hii=n_elements(hii_indices) ;new number of HII regions
377 440 ENDIF
378   - hii_count[i]=count
379   - hii_numbers[i]=i+1
380   -ENDFOR
381   -;=== remove unfound HII regions
382   -ind=where(hii_count NE 0,count)
383   -IF count NE 0 THEN BEGIN
384   - hii_count=hii_count[ind]
385   - hii_indices=hii_indices[ind]
386   - hii_numbers=hii_numbers[ind] ;Original HII region number
387   - N_hii=n_elements(hii_indices) ;new number of HII regions
388   -ENDIF
389   -;=== order HII regions by size
390   -order=sort(hii_count)
391   -hii_count=hii_count[order]
392   -hii_indices=hii_indices[order]
393   -hii_numbers=hii_numbers[order]
394   -;=== do a surounding mask per HII region
395   -around_hii_indices=ptrarr(N_hii)
396   -Npix=5 ;number of pixels by which to extend
397   -count_env=lonarr(N_hii)
398   -maskk=intarr(Nx,Ny)
399   -ind=where(is_hii EQ 1)
400   -maskk[ind]=1
401   -FOR i=0L,N_hii-1 DO BEGIN
402   - indd=enlarge(*hii_indices[i],Npix,Nx,Ny,/remove_seed,remove_these=maskk,count=count) ;This enlarges but does remove coliding hii regions.
403   - IF count NE 0 THEN around_hii_indices[i]=ptr_new(indd)
404   - count_env[i]=count
405   -ENDFOR
406   -stop
  441 + ;=== order HII regions by size
  442 + order=sort(hii_count)
  443 + hii_count=hii_count[order]
  444 + hii_indices=hii_indices[order]
  445 + hii_numbers=hii_numbers[order]
  446 + ;=== do a surounding mask per HII region
  447 + around_hii_indices=ptrarr(N_hii)
  448 + Npix=5 ;number of pixels by which to extend
  449 + count_env=lonarr(N_hii)
  450 + maskk=intarr(Nx,Ny)
  451 + ind=where(is_hii EQ 1)
  452 + maskk[ind]=1
  453 + FOR i=0L,N_hii-1 DO BEGIN
  454 + indd=enlarge(*hii_indices[i],Npix,Nx,Ny,/remove_seed,remove_these=maskk,count=count) ;This enlarges but does remove coliding hii regions.
  455 + IF count NE 0 THEN around_hii_indices[i]=ptr_new(indd)
  456 + count_env[i]=count
  457 + ENDFOR
  458 + ;stop
  459 +
  460 + ;=== compute global statistics of fitted parameters in various regions
  461 + one_stat={region:'',count:0L,med_G0:la_undef(4),med_YG1:la_undef(4),med_YG2:la_undef(4),med_YG3:la_undef(4),med_Xion:la_undef(4),med_rchi2:la_undef(4)}
  462 + regions=['ALL','HII','not-HII']
  463 + Nregions=n_elements(regions)
  464 + regional_stats=replicate(one_stat,Nregions)
  465 + regional_stats.region=regions
  466 + FOR i=0L,Nregions-1 DO BEGIN
  467 + CASE regions[i] OF
  468 + 'ALL': BEGIN
  469 + ind=where(G0_map NE la_undef(),count)
  470 + END
  471 + 'HII': BEGIN
  472 + ind=where(G0_map NE la_undef() AND is_hii NE 0 AND is_hii NE la_undef(),count)
  473 + END
  474 + 'not-HII': BEGIN
  475 + ind=where(G0_map NE la_undef() AND is_hii EQ 1 AND is_hii NE la_undef(),count)
  476 + END
  477 + ENDCASE
  478 + regional_stats[i].med_G0=la_median(G0_map[ind])
  479 + regional_stats[i].count=count
  480 + regional_stats[i].med_YG1=la_median(Yg1_map[ind])
  481 + regional_stats[i].med_YG2=la_median(Yg2_map[ind])
  482 + regional_stats[i].med_YG3=la_median(Yg3_map[ind])
  483 + regional_stats[i].med_Xion=la_median(Xion_map[ind])
  484 + regional_stats[i].med_rchi2=la_median(rchi2_map[ind])
  485 + ENDFOR
407 486  
408   -;=== compute statistics of fitted parameters in various regions
409   -one_stat={region:'',count:0L,med_G0:la_undef(4),med_YG1:la_undef(4),med_YG2:la_undef(4),med_YG3:la_undef(4),med_Xion:la_undef(4),med_rchi2:la_undef(4)}
410   -regions=['ALL','HII','not-HII']
411   -Nregions=n_elements(regions)
412   -regional_stats=replicate(one_stat,Nregions)
413   -regional_stats.region=regions
414   -FOR i=0L,Nregions-1 DO BEGIN
415   - CASE regions[i] OF
416   - 'ALL': BEGIN
417   - ind=where(G0_map NE la_undef(),count)
418   - END
419   - 'HII': BEGIN
420   - ind=where(G0_map NE la_undef() AND is_hii NE 0 AND is_hii NE la_undef(),count)
421   - END
422   - 'not-HII': BEGIN
423   - ind=where(G0_map NE la_undef() AND is_hii EQ 1 AND is_hii NE la_undef(),count)
424   - END
425   - ENDCASE
426   - regional_stats[i].med_G0=la_median(G0_map[ind])
427   - regional_stats[i].count=count
428   - regional_stats[i].med_YG1=la_median(Yg1_map[ind])
429   - regional_stats[i].med_YG2=la_median(Yg2_map[ind])
430   - regional_stats[i].med_YG3=la_median(Yg3_map[ind])
431   - regional_stats[i].med_Xion=la_median(Xion_map[ind])
432   - regional_stats[i].med_rchi2=la_median(rchi2_map[ind])
433   -ENDFOR
  487 + write_xcat,regional_stats,'/tmp/bidon',/nowrite
  488 +
  489 + ;=== compute local statistics of fitted parameters locally around HII regions
  490 + one_stat={region:'',med_NH:la_undef(4),std_NH:la_undef(4),med_G0:la_undef(4),std_G0:la_undef(4),med_YG1:la_undef(4),med_YG2:la_undef(4),med_YG3:la_undef(4),med_Xion:la_undef(4),med_rchi2:la_undef(4)}
  491 + regions=['HII']
  492 + Nregions=n_elements(regions)
  493 + hii_stats=replicate(one_stat,Nregions)
  494 + hii_stats.region=regions
  495 + med_NH_in=fltarr(N_hii) & med_NH_out=fltarr(N_hii)
  496 + med_G0_in=fltarr(N_hii) & med_G0_out=fltarr(N_hii)
  497 + med_YG1_in=fltarr(N_hii) & med_YG1_out=fltarr(N_hii)
  498 + med_YG2_in=fltarr(N_hii) & med_YG2_out=fltarr(N_hii)
  499 + med_YG3_in=fltarr(N_hii) & med_YG3_out=fltarr(N_hii)
  500 + med_Xion_in=fltarr(N_hii) & med_Xion_out=fltarr(N_hii)
  501 + med_rchi2_in=fltarr(N_hii) & med_rchi2_out=fltarr(N_hii)
  502 + FOR i=0L,N_hii-1 DO BEGIN
  503 + med_NH_in[i]=la_median(NH_values_map[*hii_indices[i]])
  504 + med_NH_out[i]=la_median(NH_values_map[*around_hii_indices[i]])
  505 + med_G0_in[i]=la_median(G0_map[*hii_indices[i]])
  506 + med_G0_out[i]=la_median(G0_map[*around_hii_indices[i]])
  507 + med_YG1_in[i]=la_median(YG1_map[*hii_indices[i]])
  508 + med_YG1_out[i]=la_median(YG1_map[*around_hii_indices[i]])
  509 + med_YG2_in[i]=la_median(YG2_map[*hii_indices[i]])
  510 + med_YG2_out[i]=la_median(YG2_map[*around_hii_indices[i]])
  511 + med_YG3_in[i]=la_median(YG3_map[*hii_indices[i]])
  512 + med_YG3_out[i]=la_median(YG3_map[*around_hii_indices[i]])
  513 + med_Xion_in[i]=la_median(Xion_map[*hii_indices[i]])
  514 + med_Xion_out[i]=la_median(Xion_map[*around_hii_indices[i]])
  515 + med_rchi2_in[i]=la_median(rchi2_map[*hii_indices[i]])
  516 + med_rchi2_out[i]=la_median(rchi2_map[*around_hii_indices[i]])
  517 + ENDFOR
  518 + ;=== compute average HII statistics per HII region size
  519 + ind=where(med_NH_in NE la_undef() and med_NH_out NE la_undef(),count)
  520 + hii_count=hii_count[ind]
  521 + med_NH_in=med_NH_in[ind] & med_NH_out=med_NH_out[ind]
  522 + med_G0_in=med_G0_in[ind] & med_G0_out=med_G0_out[ind]
  523 + med_YG1_in=med_YG1_in[ind] & med_YG1_out=med_YG1_out[ind]
  524 + med_YG2_in=med_YG2_in[ind] & med_YG2_out=med_YG2_out[ind]
  525 + med_YG3_in=med_YG3_in[ind] & med_YG3_out=med_YG3_out[ind]
  526 + med_Xion_in=med_Xion_in[ind] & med_Xion_out=med_Xion_out[ind]
  527 + med_rchi2_in=med_rchi2_in[ind] & med_rchi2_out=med_rchi2_out[ind]
  528 +
  529 + order=sort(hii_count)
  530 + toto=hii_count[order]
  531 + un=uniq(toto)
  532 + un_sizes=toto[un]
  533 + Nun_sizes=n_elements(un_sizes)
  534 + one_stat={size:0.,NHII:0L, $
  535 + med_d_NH:la_undef(4),std_d_NH:la_undef(4), $
  536 + med_d_G0:la_undef(4),std_d_G0:la_undef(4), $
  537 + med_d_YG1:la_undef(4),std_d_YG1:la_undef(4), $
  538 + med_d_YG2:la_undef(4),std_d_YG2:la_undef(4), $
  539 + med_d_YG3:la_undef(4),std_d_YG3:la_undef(4), $
  540 + med_d_Xion:la_undef(4),std_d_Xion:la_undef(4), $
  541 + med_d_rchi2:la_undef(4),std_d_rchi2:la_undef(4)}
  542 + size_hii_stat=replicate(one_stat,Nun_sizes)
  543 + FOR i=0L,Nun_sizes-1 DO BEGIN
  544 + ind=where(hii_count EQ un_sizes[i],count)
  545 + size_hii_stat[i].size=la_median(hii_count[ind])
  546 + size_hii_stat[i].NHII=count
  547 + vv=(med_NH_in[ind]-med_NH_out[ind])/med_NH_out[ind]*100.
  548 + size_hii_stat[i].med_d_NH=la_median(vv) & size_hii_stat[i].std_d_NH=la_sigma(vv)
  549 + vv=(med_G0_in[ind]-med_G0_out[ind])/med_G0_out[ind]*100.
  550 + size_hii_stat[i].med_d_G0=la_median(vv) & size_hii_stat[i].std_d_G0=la_sigma(vv)
  551 + vv=(med_YG1_in[ind]-med_YG1_out[ind])/med_YG1_out[ind]*100.
  552 + size_hii_stat[i].med_d_YG1=la_median(vv) & size_hii_stat[i].std_d_YG1=la_sigma(vv)
  553 + vv=(med_YG2_in[ind]-med_YG2_out[ind])/med_YG2_out[ind]*100.
  554 + size_hii_stat[i].med_d_YG2=la_median(vv) & size_hii_stat[i].std_d_YG2=la_sigma(vv)
  555 + vv=(med_YG3_in[ind]-med_YG3_out[ind])/med_YG3_out[ind]*100.
  556 + size_hii_stat[i].med_d_YG3=la_median(vv) & size_hii_stat[i].std_d_YG3=la_sigma(vv)
  557 + vv=(med_Xion_in[ind]-med_Xion_out[ind])/med_Xion_out[ind]*100.
  558 + size_hii_stat[i].med_d_Xion=la_median(vv) & size_hii_stat[i].std_d_Xion=la_sigma(vv)
  559 + vv=(med_rchi2_in[ind]-med_rchi2_out[ind])/med_rchi2_out[ind]*100.
  560 + size_hii_stat[i].med_d_rchi2=la_median(vv) & size_hii_stat[i].std_d_rchi2=la_sigma(vv)
  561 + ENDFOR
434 562  
435   -write_xcat,regional_stats,'/tmp/bidon',/nowrite
436   -
437   -;=== compute statistics of fitted parameters locally around HII regions
438   -one_stat={region:'',med_NH:la_undef(4),std_NH:la_undef(4),med_G0:la_undef(4),std_G0:la_undef(4),med_YG1:la_undef(4),med_YG2:la_undef(4),med_YG3:la_undef(4),med_Xion:la_undef(4),med_rchi2:la_undef(4)}
439   -regions=['HII']
440   -Nregions=n_elements(regions)
441   -hii_stats=replicate(one_stat,Nregions)
442   -hii_stats.region=regions
443   -med_NH_in=fltarr(N_hii) & med_NH_out=fltarr(N_hii)
444   -med_G0_in=fltarr(N_hii) & med_G0_out=fltarr(N_hii)
445   -med_YG1_in=fltarr(N_hii) & med_YG1_out=fltarr(N_hii)
446   -med_YG2_in=fltarr(N_hii) & med_YG2_out=fltarr(N_hii)
447   -med_YG3_in=fltarr(N_hii) & med_YG3_out=fltarr(N_hii)
448   -med_Xion_in=fltarr(N_hii) & med_Xion_out=fltarr(N_hii)
449   -med_rchi2_in=fltarr(N_hii) & med_rchi2_out=fltarr(N_hii)
450   -FOR i=0L,N_hii-1 DO BEGIN
451   - med_NH_in[i]=la_median(NH_values_map[*hii_indices[i]])
452   - med_NH_out[i]=la_median(NH_values_map[*around_hii_indices[i]])
453   - med_G0_in[i]=la_median(G0_map[*hii_indices[i]])
454   - med_G0_out[i]=la_median(G0_map[*around_hii_indices[i]])
455   - med_YG1_in[i]=la_median(YG1_map[*hii_indices[i]])
456   - med_YG1_out[i]=la_median(YG1_map[*around_hii_indices[i]])
457   - med_YG2_in[i]=la_median(YG2_map[*hii_indices[i]])
458   - med_YG2_out[i]=la_median(YG2_map[*around_hii_indices[i]])
459   - med_YG3_in[i]=la_median(YG3_map[*hii_indices[i]])
460   - med_YG3_out[i]=la_median(YG3_map[*around_hii_indices[i]])
461   - med_Xion_in[i]=la_median(Xion_map[*hii_indices[i]])
462   - med_Xion_out[i]=la_median(Xion_map[*around_hii_indices[i]])
463   - med_rchi2_in[i]=la_median(rchi2_map[*hii_indices[i]])
464   - med_rchi2_out[i]=la_median(rchi2_map[*around_hii_indices[i]])
465   -ENDFOR
466   -window,win,xsize=wxs,ysize=wys & win=win+1
467   -!p.multi=[0,3,3]
468   -ind=where(med_NH_in NE la_undef() and med_NH_out NE la_undef(),count)
469   -vv=(med_NH_in[ind]-med_NH_out[ind])/med_NH_out[ind]*100.
470   -xrange=[1,max(hii_count)]
471   -cgplot,hii_count[ind],vv,tit='NH',xrange=xrange,/xlog & cgoplot,10^!x.crange,[0,0.],linestyle=2,color='red'
472   -hii_stats.med_NH=la_median(vv) & hii_stats.std_NH=hii_stats.med_NH/la_sigma(vv)
473   -ind=where(med_G0_in NE la_undef() and med_G0_out NE la_undef(),count)
474   -vv=(med_G0_in[ind]-med_G0_out[ind])/med_G0_out[ind]*100.
475   -cgplot,hii_count[ind],vv,tit='G0',xrange=xrange,/xlog & cgoplot,10^!x.crange,[0,0.],linestyle=2,color='red'
476   -hii_stats.med_G0=la_median(vv) & hii_stats.std_G0=hii_stats.med_G0/la_sigma(vv)
477   -ind=where(med_YG1_in NE la_undef() and med_YG1_out NE la_undef(),count)
478   -vv=(med_YG1_in[ind]-med_YG1_out[ind])/med_YG1_out[ind]*100.
479   -cgplot,hii_count[ind],vv,tit='YG1',xrange=xrange,/xlog & cgoplot,10^!x.crange,[0,0.],linestyle=2,color='red'
480   -hii_stats.med_YG1=la_median(vv)
481   -ind=where(med_YG2_in NE la_undef() and med_YG2_out NE la_undef(),count)
482   -vv=(med_YG2_in[ind]-med_YG2_out[ind])/med_YG2_out[ind]*100.
483   -cgplot,hii_count[ind],vv,tit='YG2',xrange=xrange,/xlog & cgoplot,10^!x.crange,[0,0.],linestyle=2,color='red'
484   -hii_stats.med_YG2=la_median(vv)
485   -ind=where(med_YG3_in NE la_undef() and med_YG3_out NE la_undef(),count)
486   -vv=(med_YG3_in[ind]-med_YG3_out[ind])/med_YG3_out[ind]*100.
487   -cgplot,hii_count[ind],vv,tit='YG3',xrange=xrange,/xlog & cgoplot,10^!x.crange,[0,0.],linestyle=2,color='red'
488   -hii_stats.med_YG3=la_median(vv)
489   -ind=where(med_Xion_in NE la_undef() and med_Xion_out NE la_undef(),count)
490   -vv=(med_Xion_in[ind]-med_Xion_out[ind])/med_Xion_out[ind]*100.
491   -cgplot,hii_count[ind],vv,tit='Xion',xrange=xrange,/xlog & cgoplot,10^!x.crange,[0,0.],linestyle=2,color='red'
492   -hii_stats.med_Xion=la_median(vv)
493   -ind=where(med_rchi2_in NE la_undef() and med_rchi2_out NE la_undef(),count)
494   -vv=(med_rchi2_in[ind]-med_rchi2_out[ind])/med_rchi2_out[ind]*100.
495   -cgplot,hii_count[ind],vv,tit='rchi2',xrange=xrange,/xlog,/ylog & cgoplot,10^!x.crange,[1,1.],linestyle=2,color='red'
496   -hii_stats.med_rchi2=la_median(vv)
497   -
498   -write_xcat,hii_stats,'/tmp/bidon',/nowrite
  563 + write_xcat,size_hii_stat,/nowrite
  564 + ;stop
  565 +
  566 + window,win,xsize=wxs*2,ysize=wys & win=win+1
  567 + !p.multi=[0,3,2]
  568 + xrange=[1,max(hii_count)]
  569 +
  570 + ind=where(med_NH_in NE la_undef() and med_NH_out NE la_undef(),count)
  571 + vv=(med_NH_in[ind]-med_NH_out[ind])/med_NH_out[ind]*100.
  572 + cgplot,hii_count,vv,tit='NH',xrange=xrange,/xlog
  573 + cgoplot,10^!x.crange,[0,0.],linestyle=2,color='red'
  574 + cgoplot,size_hii_stat.size,size_hii_stat.med_d_NH,color='blue',psym='filled circle'
  575 + hii_stats.med_NH=la_median(size_hii_stat.med_d_NH)
  576 + cgoplot,10^!x.crange,replicate(hii_stats.med_NH,2),linestyle=2,color='blue'
  577 + hii_stats.std_NH=hii_stats.med_NH/la_sigma(vv)
  578 +
  579 + ind=where(med_G0_in NE la_undef() and med_G0_out NE la_undef(),count)
  580 + vv=(med_G0_in[ind]-med_G0_out[ind])/med_G0_out[ind]*100.
  581 + cgplot,hii_count[ind],vv,tit='G0',xrange=xrange,/xlog
  582 + cgoplot,10^!x.crange,[0,0.],linestyle=2,color='red'
  583 + cgoplot,size_hii_stat.size,size_hii_stat.med_d_G0,color='blue',psym='filled circle'
  584 + hii_stats.med_G0=la_median(size_hii_stat.med_d_G0)
  585 + cgoplot,10^!x.crange,replicate(la_median(size_hii_stat.med_d_G0),2),linestyle=2,color='blue'
  586 + hii_stats.std_G0=hii_stats.med_G0/la_sigma(vv)
  587 +
  588 + ind=where(med_YG1_in NE la_undef() and med_YG1_out NE la_undef(),count)
  589 + vv=(med_YG1_in[ind]-med_YG1_out[ind])/med_YG1_out[ind]*100.
  590 + cgplot,hii_count[ind],vv,tit='YG1',xrange=xrange,/xlog
  591 + cgoplot,10^!x.crange,[0,0.],linestyle=2,color='red'
  592 + cgoplot,size_hii_stat.size,size_hii_stat.med_d_YG1,color='blue',psym='filled circle'
  593 + hii_stats.med_YG1=la_median(size_hii_stat.med_d_YG1)
  594 + cgoplot,10^!x.crange,replicate(hii_stats.med_YG1,2),linestyle=2,color='blue'
  595 +
  596 + ind=where(med_YG2_in NE la_undef() and med_YG2_out NE la_undef(),count)
  597 + vv=(med_YG2_in[ind]-med_YG2_out[ind])/med_YG2_out[ind]*100.
  598 + cgplot,hii_count[ind],vv,tit='YG2',xrange=xrange,/xlog
  599 + cgoplot,10^!x.crange,[0,0.],linestyle=2,color='red'
  600 + cgoplot,size_hii_stat.size,size_hii_stat.med_d_YG2,color='blue',psym='filled circle'
  601 + hii_stats.med_YG2=la_median(size_hii_stat.med_d_YG2)
  602 + cgoplot,10^!x.crange,replicate(hii_stats.med_YG2,2),linestyle=2,color='blue'
  603 +
  604 + goto,jump_this_one
  605 + ind=where(med_YG3_in NE la_undef() and med_YG3_out NE la_undef(),count)
  606 + vv=(med_YG3_in[ind]-med_YG3_out[ind])/med_YG3_out[ind]*100.
  607 + cgplot,hii_count[ind],vv,tit='YG3',xrange=xrange,/xlog
  608 + cgoplot,10^!x.crange,[0,0.],linestyle=2,color='red'
  609 + cgoplot,size_hii_stat.size,size_hii_stat.med_d_YG3,color='blue',psym='filled circle'
  610 + hii_stats.med_YG3=la_median(size_hii_stat.med_d_YG3)
  611 + cgoplot,10^!x.crange,replicate(hii_stats.med_YG3,2),linestyle=2,color='blue'
  612 + jump_this_one:
  613 +
  614 + ind=where(med_Xion_in NE la_undef() and med_Xion_out NE la_undef(),count)
  615 + vv=(med_Xion_in[ind]-med_Xion_out[ind])/med_Xion_out[ind]*100.
  616 + cgplot,hii_count[ind],vv,tit='Xion',xrange=xrange,/xlog
  617 + cgoplot,10^!x.crange,[0,0.],linestyle=2,color='red'
  618 + cgoplot,size_hii_stat.size,size_hii_stat.med_d_Xion,color='blue',psym='filled circle'
  619 + hii_stats.med_Xion=la_median(size_hii_stat.med_d_Xion)
  620 + cgoplot,10^!x.crange,replicate(hii_stats.med_Xion,2),linestyle=2,color='blue'
  621 +
  622 + ind=where(med_rchi2_in NE la_undef() and med_rchi2_out NE la_undef(),count)
  623 + vv=(med_rchi2_in[ind]-med_rchi2_out[ind])/med_rchi2_out[ind]*100.
  624 + cgplot,hii_count[ind],vv,tit='rchi2',xrange=xrange,/xlog,yrange=[-100,500],/ysty;,/ylog
  625 + cgoplot,10^!x.crange,[1,1.],linestyle=2,color='red'
  626 + cgoplot,size_hii_stat.size,size_hii_stat.med_d_rchi2,color='blue',psym='filled circle'
  627 + hii_stats.med_rchi2=la_median(size_hii_stat.med_d_rchi2)
  628 + cgoplot,10^!x.crange,replicate(hii_stats.med_rchi2,2),linestyle=2,color='blue'
  629 +
  630 + write_xcat,hii_stats,'/tmp/bidon',/nowrite
499 631  
500   -cleanplot
  632 + cleanplot
  633 +ENDIF
501 634  
502   -stop
  635 +;stop
503 636  
504 637 ;===== plot nhits vs isrf_class
505 638 Ncl=n_elements(stats_st)
... ... @@ -594,9 +727,9 @@ cgoplot,replicate(alog10(parval_min),2),10^(!y.crange),linestyle=2,color='blue'
594 727 cgoplot,replicate(alog10(parval_max),2),10^(!y.crange),linestyle=2,color='red'
595 728  
596 729 ind=where(Ygrain1 NE la_undef())
597   -default_value=dustem_get_param_values_default('(*!dustem_params).grains[0].MDUST_O_MH')
598   -parval_min=((*!dustem_grid).PMIN_VALUES)[1]
599   -parval_max=((*!dustem_grid).PMAX_VALUES)[1]
  730 +default_value=dustem_get_param_values_default('(*!dustem_params).grains[0].MDUST_O_MH')/2. ;/2. is because in grain.dat, this value is for the sum of PAH+ and PAH0
  731 +parval_min=((*!dustem_grid).PMIN_VALUES)[1]/2. ;/2. is because in grain.dat, this value is for the sum of PAH+ and PAH0
  732 +parval_max=((*!dustem_grid).PMAX_VALUES)[1]/2. ;/2. is because in grain.dat, this value is for the sum of PAH+ and PAH0
600 733 ;xrange=alog10([default_value/range_fact,default_value*range_fact])
601 734 xrange=enlarge_range(alog10([parval_min,parval_max]),0.1)
602 735 yrange=[1,max(res)]
... ... @@ -608,9 +741,9 @@ cgoplot,replicate(alog10(parval_max),2),10^(!y.crange),linestyle=2,color='red'
608 741  
609 742 ind=where(Ygrain2 NE la_undef(),count)
610 743 IF count NE 0 THEN BEGIN
611   - default_value=dustem_get_param_values_default('(*!dustem_params).grains[1].MDUST_O_MH')
612   - parval_min=((*!dustem_grid).PMIN_VALUES)[2]
613   - parval_max=((*!dustem_grid).PMAX_VALUES)[2]
  744 + default_value=dustem_get_param_values_default('(*!dustem_params).grains[1].MDUST_O_MH')/2. ;/2. is because in grain.dat, this value is for the sum of PAH+ and PAH0
  745 + parval_min=((*!dustem_grid).PMIN_VALUES)[2]/2. ;/2. is because in grain.dat, this value is for the sum of PAH+ and PAH0
  746 + parval_max=((*!dustem_grid).PMAX_VALUES)[2]/2. ;/2. is because in grain.dat, this value is for the sum of PAH+ and PAH0
614 747 ;xrange=alog10([default_value/range_fact,default_value*range_fact])
615 748 xrange=enlarge_range(alog10([parval_min,parval_max]),0.1)
616 749 yrange=[1,max(res)]
... ... @@ -670,9 +803,9 @@ IF keyword_set(do_postscripts) THEN BEGIN
670 803 device, filename=ps_file
671 804  
672 805 ind=where(Ygrain1 NE la_undef())
673   - default_value=dustem_get_param_values_default('(*!dustem_params).grains[0].MDUST_O_MH')
674   - parval_min=((*!dustem_grid).PMIN_VALUES)[1]
675   - parval_max=((*!dustem_grid).PMAX_VALUES)[1]
  806 + default_value=dustem_get_param_values_default('(*!dustem_params).grains[0].MDUST_O_MH')/2. ;/2. is because in grain.dat, this value is for the sum of PAH+ and PAH0
  807 + parval_min=((*!dustem_grid).PMIN_VALUES)[1]/2. ;This is because in grain.dat, this value is for the sum of PAH+ and PAH0
  808 + parval_max=((*!dustem_grid).PMAX_VALUES)[1]/2. ;This is because in grain.dat, this value is for the sum of PAH+ and PAH0
676 809 ;xrange=alog10([default_value/range_fact,default_value*range_fact])
677 810 xrange=enlarge_range(alog10([parval_min,parval_max]),0.1)
678 811 yrange=[1,max(res)]
... ... @@ -690,9 +823,9 @@ IF keyword_set(do_postscripts) THEN BEGIN
690 823 device, filename=ps_file
691 824  
692 825 ind=where(Ygrain2 NE la_undef())
693   - default_value=dustem_get_param_values_default('(*!dustem_params).grains[1].MDUST_O_MH')
694   - parval_min=((*!dustem_grid).PMIN_VALUES)[2]
695   - parval_max=((*!dustem_grid).PMAX_VALUES)[2]
  826 + default_value=dustem_get_param_values_default('(*!dustem_params).grains[1].MDUST_O_MH')/2. ;/2. is because in grain.dat, this value is for the sum of PAH+ and PAH0
  827 + parval_min=((*!dustem_grid).PMIN_VALUES)[2]/2. ;This is because in grain.dat, this value is for the sum of PAH+ and PAH0
  828 + parval_max=((*!dustem_grid).PMAX_VALUES)[2]/2. ;This is because in grain.dat, this value is for the sum of PAH+ and PAH0
696 829 ;xrange=alog10([default_value/range_fact,default_value*range_fact])
697 830 xrange=enlarge_range(alog10([parval_min,parval_max]),0.1)
698 831 yrange=[1,max(res)]
... ... @@ -730,6 +863,169 @@ IF keyword_set(do_postscripts) THEN BEGIN
730 863 ENDIF
731 864 ;stop
732 865  
  866 +;==== plot tau_uv/tau_vis vs Xion
  867 +dustem_init,model='DL07'
  868 +st_model=dustem_read_all(!dustem_soft_dir)
  869 +dustem_write_all,st_model,!dustem_dat
  870 +st=dustem_run()
  871 +;compute standard extinction curves
  872 +tau_pah0=reform((st.ext.ABS_GRAIN)[0,*]+(st.ext.SCA_GRAIN)[0,*])
  873 +tau_pahp=reform((st.ext.ABS_GRAIN)[1,*]+(st.ext.SCA_GRAIN)[1,*])
  874 +tau_vsg=reform((st.ext.ABS_GRAIN)[2,*]+(st.ext.SCA_GRAIN)[2,*])
  875 +tau_bg=reform((st.ext.ABS_GRAIN)[3,*]+(st.ext.SCA_GRAIN)[3,*])+reform((st.ext.ABS_GRAIN)[4,*]+(st.ext.SCA_GRAIN)[4,*])
  876 +tau_tot=tau_bg+tau_vsg+tau_pahp+tau_pah0
  877 +;compute default model abundances
  878 +Ypah0_s=(dustem_get_param_values_default('(*!dustem_params).grains[0].MDUST_O_MH'))[0]/2. ;/2. is because in grain.dat, this value is for the sum of PAH+ and PAH0
  879 +Ypahp_s=(dustem_get_param_values_default('(*!dustem_params).grains[1].MDUST_O_MH'))[0]/2. ;/2. is because in grain.dat, this value is for the sum of PAH+ and PAH0
  880 +Yvsg_s=(dustem_get_param_values_default('(*!dustem_params).grains[2].MDUST_O_MH'))[0]
  881 +Ybg_s=(dustem_get_param_values_default('(*!dustem_params).grains[3].MDUST_O_MH'))[0]+ $
  882 + (dustem_get_param_values_default('(*!dustem_params).grains[4].MDUST_O_MH'))[0]
  883 +;=== This is with Desert value
  884 +;Ypah0_s=4.3e-4 ;devided by 2 because assuming Xion=0.5 in Desert et al.
  885 +;YBG_s=6.4e-3
  886 +
  887 +;== compute model predicted tau values on voronoi pixels, using fitted emission
  888 +Nvor=n_elements(Ygrain1)
  889 +wwavs=dustem_get_wavelengths()
  890 +Nwwavs=n_elements(wwavs)
  891 +tau_uv_model=fltarr(Nvor,N_ebv_uv_filters)
  892 +;stop
  893 +FOR k=0L,N_ebv_uv_filters-1 DO BEGIN
  894 + tau_uv_model[*,k]=(interpol(tau_pah0,wwavs,wwav_uv[k])/Ypah0_s*Ygrain1+ $
  895 + interpol(tau_pahp,wwavs,wwav_uv[k])/Ypahp_s*Ygrain2+ $
  896 + interpol(tau_vsg,wwavs,wwav_uv[k])/Yvsg_s*Ygrain3+ $
  897 + interpol(tau_bg,wwavs,wwav_uv[k])/Ybg_s*Ybg_s)*NH_values
  898 +ENDFOR
  899 +tau_vis_model=(interpol(tau_pah0,wwavs,wwav_star)/Ypah0_s*Ygrain1+ $
  900 + interpol(tau_pahp,wwavs,wwav_star)/Ypahp_s*Ygrain2+ $
  901 + interpol(tau_vsg,wwavs,wwav_star)/Yvsg_s*Ygrain3+ $
  902 + interpol(tau_bg,wwavs,wwav_star)/Ybg_s*Ybg_s)*NH_values
  903 +
  904 +ind=where(Xions NE la_undef() and finite(tausV),count)
  905 +xion_min=0.5
  906 +ind_low_Xion=where(Xions NE la_undef() and finite(tausV) and Xions LE Xion_min,count)
  907 +ind_high_Xion=where(Xions NE la_undef() and finite(tausV) and Xions GT Xion_min,count)
  908 +Xion_low=la_median(Xions[ind_low_Xion])
  909 +Xion_high=la_median(Xions[ind_high_Xion])
  910 +
  911 +window,win & win=win+1
  912 +k=0L
  913 +yrange=[1.5,10]
  914 +ytit='IR fit predicted tau_uv / tau_V []'
  915 +xtit='Voronoi bin []'
  916 +cgplot,Xions[ind],tau_uv_model[ind,0]/tau_vis_model[ind],/ysty,psym='filled circle',/nodata,yrange=yrange,xtit=xtit,ytit=ytit
  917 +color_uv_filters=round(range_gen(N_ebv_uv_filters,[200,50]))
  918 +avg_predicted_tauuv_o_tauv=fltarr(N_ebv_uv_filters)
  919 +FOR k=0L,N_ebv_uv_filters-1 DO BEGIN
  920 + cgoplot,Xions[ind],tau_uv_model[ind,k]/tau_vis_model[ind],color=color_uv_filters[k],psym='filled circle'
  921 + avg_predicted_tauuv_o_tauv[k]=la_mean(tau_uv_model[ind,k]/tau_vis_model[ind])
  922 +ENDFOR
  923 +avg_predicted_tauv=la_mean(tau_vis_model[ind])
  924 +
  925 +;compute the average observed stellar extinction for 1.e21 H/cm2
  926 +star_abs_tau_uv=fltarr(N_ebv_uv_filters)
  927 +FOR k=0L,N_ebv_uv_filters-1 DO BEGIN
  928 + star_abs_tau_uv[k]=la_mean(taus_uv[ind,k]/NH_values[ind])*10. ;10, because NH_values is in 1e20 H/cm2
  929 +ENDFOR
  930 +star_abs_tau_v=la_mean(tausV[ind]/NH_values[ind])*10. ;10, because NH_values is in 1e20 H/cm2
  931 +
  932 +star_abs_tau_uv_lowXion=fltarr(N_ebv_uv_filters)
  933 +FOR k=0L,N_ebv_uv_filters-1 DO BEGIN
  934 + star_abs_tau_uv_lowXion[k]=la_mean(taus_uv[ind_low_Xion,k]/NH_values[ind_low_Xion])*10. ;10, because NH_values is in 1e20 H/cm2
  935 +ENDFOR
  936 +star_abs_tau_v_lowXion=la_mean(tausV[ind_low_Xion]/NH_values[ind_low_Xion])*10. ;10, because NH_values is in 1e20 H/cm2
  937 +star_abs_tau_uv_highXion=fltarr(N_ebv_uv_filters)
  938 +FOR k=0L,N_ebv_uv_filters-1 DO BEGIN
  939 + star_abs_tau_uv_highXion[k]=la_mean(taus_uv[ind_high_Xion,k]/NH_values[ind_high_Xion])*10. ;10, because NH_values is in 1e20 H/cm2
  940 +ENDFOR
  941 +star_abs_tau_v_highXion=la_mean(tausV[ind_high_Xion]/NH_values[ind_high_Xion])*10. ;10, because NH_values is in 1e20 H/cm2
  942 +
  943 +;== plot the standard model extinction curve and predicted (dust) and observed (stars) extinction SEDs
  944 +tit='solid=model prediction, open=stellar extinction'
  945 +ytit='tau (for 1.e21 H/cm2)'
  946 +xtit='1/wav [mic-1]'
  947 +window,win & win=win+1
  948 +xrange=[0,10.]
  949 +cgplot,1/wwavs,tau_tot,xrange=xrange,/xsty,tit=tit,xtit=xtit,ytit=ytit
  950 +tau_tot_v=interpol(tau_tot,wwavs,wwav_star)
  951 +cgoplot,1/wwavs,tau_bg,linestyle=0,color='red'
  952 +cgoplot,1/wwavs,tau_vsg,linestyle=0,color='orange'
  953 +cgoplot,1/wwavs,tau_pahp,linestyle=0,color='violet'
  954 +cgoplot,1/wwavs,tau_pah0,linestyle=0,color='blue'
  955 +FOR k=0L,N_ebv_uv_filters-1 DO BEGIN
  956 + cgoplot,replicate(1./wwav_uv[k],2),!y.crange,linestyle=1,color=color_uv_filters[k]
  957 + cgoplot,1./wwav_uv[k],avg_predicted_tauuv_o_tauv[k]*tau_tot_v,psym='filled circle',color=color_uv_filters[k]
  958 +ENDFOR
  959 +cgoplot,replicate(1./wwav_star,2),!y.crange,linestyle=1
  960 +cgoplot,1./wwav_star,tau_tot_v,psym='filled circle',color='black'
  961 +FOR k=0L,N_ebv_uv_filters-1 DO BEGIN
  962 + cgoplot,1./wwav_uv[k],star_abs_tau_uv[k],psym='open circle',color=color_uv_filters[k]
  963 +ENDFOR
  964 +cgoplot,1./wwav_star,star_abs_tau_v,psym='open circle',color='black'
  965 +fact=tau_tot_v/star_abs_tau_v
  966 +FOR k=0L,N_ebv_uv_filters-1 DO BEGIN
  967 + cgoplot,1./wwav_uv[k],star_abs_tau_uv[k]*fact,psym='open square',color=color_uv_filters[k]
  968 +ENDFOR
  969 +cgoplot,1./wwav_star,star_abs_tau_v*fact,psym='open square',color='black'
  970 +print,fact
  971 +FOR k=0L,N_ebv_uv_filters-1 DO BEGIN
  972 + cgoplot,1./wwav_uv[k],star_abs_tau_uv_lowXion[k],psym='Open Down Triangle',color=color_uv_filters[k]
  973 +ENDFOR
  974 +cgoplot,1./wwav_star,star_abs_tau_v_lowXion,psym='Open Down Triangle',color='black'
  975 +FOR k=0L,N_ebv_uv_filters-1 DO BEGIN
  976 + cgoplot,1./wwav_uv[k],star_abs_tau_uv_highXion[k],psym='Open Up Triangle',color=color_uv_filters[k]
  977 +ENDFOR
  978 +cgoplot,1./wwav_star,star_abs_tau_v_highXion,psym='Open Up Triangle',color='black'
  979 +
  980 +
  981 +stop
  982 +
  983 +;=== print UV to VIS tau ratios for various dust components
  984 +print,interpol(tau_pah0,wwavs,wwav_uv)/interpol(tau_pah0,wwavs,wwav_star)
  985 +print,interpol(tau_pahp,wwavs,wwav_uv)/interpol(tau_pahp,wwavs,wwav_star)
  986 +print,interpol(tau_vsg,wwavs,wwav_uv)/interpol(tau_vsg,wwavs,wwav_star)
  987 +print,interpol(tau_bg,wwavs,wwav_uv)/interpol(tau_bg,wwavs,wwav_star)
  988 +
  989 +tau_bar_uv0=interpol(tau_pah0,wwavs,wwav_uv)/Ypah0_s[0]
  990 +tau_bar_v0=interpol(tau_bg,wwavs,wwav_star)/Ybg_s[0]
  991 +tau_bar_uv_bg=interpol(tau_bg,wwavs,wwav_uv)/Ybg_s[0]
  992 +;=== This is with Desert value
  993 +;tau_bar_uv0=0.2/Ypah0_s
  994 +;tau_bar_v0=0.42/YBG_s
  995 +
  996 +print,tau_bar_uv0,tau_bar_v0,tau_bar_uv0/tau_bar_v0
  997 +print,wwav_uv,wwav_star
  998 +print,1./wwav_uv,1./wwav_star
  999 +
  1000 +stop
  1001 +
  1002 +window,win & win=win+1
  1003 +xtit='Xion'
  1004 +ytit='tau_uv / tau_v'
  1005 +xrange=[0,1] & yrange=[1.e-2,1.e2]
  1006 +tau_ratio=la_div(taus_uv,tausV)
  1007 +tau_rap=Xions/(1.-Xions)*(tau_ratio*YBG_s/Ygrain1*(tau_bar_v0/tau_bar_uv0-tau_bar_uv_bg/tau_bar_uv0)-1.)
  1008 +ind=where(Xions NE la_undef() AND tau_ratio NE la_undef())
  1009 +cgplot,Xions[ind],tau_ratio[ind], psym='Filled circle',color='black',xtit=xtit,ytit=ytit,xrange=xrange,/xsty,/ylog,yrange=yrange,/ysty
  1010 +cgoplot,Xions[ind],tau_rap[ind],psym='open circle',color='red'
  1011 +ind=where(tau_rap NE la_undef() and finite(tau_rap) EQ 1)
  1012 +mean_tau_rap=la_mean(tau_rap[ind])
  1013 +cgoplot,!x.crange,replicate(mean_tau_rap,2),psym=2,color='red',linestyle=2
  1014 +print,la_mean(tau_rap[ind])
  1015 +
  1016 +stop
  1017 +
  1018 +; exp_tau_ratio=tau_ratio
  1019 +; ind=where(tau_ratio NE la_undef())
  1020 +; exp_tau_ratio=exp(-1.*ebvs_uv[ind])/exp(-1.*ebvs[ind])
  1021 +; xtit='Xion'
  1022 +; ytit='exp(-tau_uv) / exp(-tau_v)'
  1023 +; xrange=[0,1]
  1024 +; ind=where(Xions NE la_undef() AND exp_tau_ratio NE la_undef())
  1025 +; cgplot,Xions[ind],exp_tau_ratio[ind], psym='Filled circle',color='black',xtit=xtit,ytit=ytit,xrange=xrange,/xsty
  1026 +
  1027 +
  1028 +
733 1029 ;==== show maps
734 1030 obp=[1.1,0.,1.15,1]
735 1031 cleanplot
... ... @@ -741,6 +1037,9 @@ CASE use_source_name OF
741 1037 'ngc4321':imrange=[-1,2.]
742 1038 'ngc3351':imrange=[-1,1.5]
743 1039 'ngc1566':imrange=[-1,1.5]
  1040 + 'ngc4254':imrange=[-1,1.5]
  1041 + 'ngc4254':imrange=[-1,1.5]
  1042 + 'ngc3627':imrange=[-1,1.5]
744 1043 ENDCASE
745 1044 title='NH map'
746 1045 image_cont20,la_log10(NHmap),Href,/square,imrange=imrange,axis_color_table=1,image_color_table=image_color_table,/silent,off_bar=obp,title=title,levels=cont_st
... ... @@ -761,6 +1060,8 @@ CASE use_source_name OF
761 1060 'ngc4321':imrange=1.+[-0.5,+0.5]
762 1061 'ngc3351':imrange=2.+[-0.5,+0.5]
763 1062 'ngc1566':imrange=1.+[-0.5,+0.5]
  1063 + 'ngc4254':imrange=1.+[-0.5,+0.5]
  1064 + 'ngc3627':imrange=1.+[-0.5,+0.5]
764 1065 ENDCASE
765 1066 ;imrange=1.+[-0.5,+0.5]
766 1067 title='stellar normalization'
... ... @@ -772,12 +1073,15 @@ IF keyword_set(do_postscripts) THEN BEGIN
772 1073 message,'Wrote '+ps_file,/continue
773 1074 ENDIF
774 1075  
  1076 +;==== EBV stars
775 1077 window,win,xsize=wxs,ysize=wys & win=win+1
776 1078 CASE use_source_name OF
777 1079 'ngc0628':imrange=[0.1,0.25]
778 1080 'ngc4321':imrange=[0.1,0.25]
779 1081 'ngc3351':imrange=[0.1,0.6]
780 1082 'ngc1566':imrange=[0.1,0.25]
  1083 + 'ngc4254':imrange=[0.1,0.25]
  1084 + 'ngc3627':imrange=[0.1,0.25]
781 1085 ENDCASE
782 1086 title='E(B-V)'
783 1087 image_cont20,EBV_map,Href,/square,imrange=imrange,axis_color_table=1,image_color_table=image_color_table,/silent,off_bar=obp,title=title,levels=cont_st
... ... @@ -789,6 +1093,60 @@ IF keyword_set(do_postscripts) THEN BEGIN
789 1093 ENDIF
790 1094 ;stop
791 1095  
  1096 +;==== tauV stars
  1097 +window,win,xsize=wxs,ysize=wys & win=win+1
  1098 +CASE use_source_name OF
  1099 + 'ngc0628':imrange=[0.1,0.25]*4.75
  1100 + 'ngc4321':imrange=[0.1,0.25]*4.75
  1101 + 'ngc3351':imrange=[0.1,0.6]*4.75
  1102 + 'ngc1566':imrange=[0.1,0.25]*4.75
  1103 + 'ngc4254':imrange=[0.1,0.25]*4.75
  1104 + 'ngc3627':imrange=[0.1,0.25]*4.75
  1105 +ENDCASE
  1106 +title='tauV stars'
  1107 +image_cont20,tauV_map,Href,/square,imrange=imrange,axis_color_table=1,image_color_table=image_color_table,/silent,off_bar=obp,title=title,levels=cont_st
  1108 +IF keyword_set(do_postscripts) THEN BEGIN
  1109 + ps_file=ps_dir+use_source_name+'_tauVstars_map'+reso_str+fit_G0_str+mathis_str+'.ps'
  1110 + image_cont20,tauV_map,Href,/square,imrange=imrange,axis_color_table=1,image_color_table=image_color_table,/silent,off_bar=obp,title=title,levels=cont_st, $
  1111 + postscript=ps_file,orientation='portrait',/nologo
  1112 + message,'Wrote '+ps_file,/continue
  1113 +ENDIF
  1114 +
  1115 +window,win,xsize=wxs,ysize=wys & win=win+1
  1116 +title='E(B-V) UV'
  1117 +image_cont20,EBV_uv_map,Href,/square,imrange=imrange,axis_color_table=1,image_color_table=image_color_table,/silent,off_bar=obp,title=title,levels=cont_st
  1118 +IF keyword_set(do_postscripts) THEN BEGIN
  1119 + ps_file=ps_dir+use_source_name+'_EBV_uv_map'+reso_str+fit_G0_str+mathis_str+'.ps'
  1120 + image_cont20,EBV_uv_map,Href,/square,imrange=imrange,axis_color_table=1,image_color_table=image_color_table,/silent,off_bar=obp,title=title,levels=cont_st, $
  1121 + postscript=ps_file,orientation='portrait',/nologo
  1122 + message,'Wrote '+ps_file,/continue
  1123 +ENDIF
  1124 +
  1125 +window,win,xsize=wxs,ysize=wys & win=win+1
  1126 +title='tau UV map'
  1127 +image_cont20,tau_uv_map,Href,/square,imrange=imrange,axis_color_table=1,image_color_table=image_color_table,/silent,off_bar=obp,title=title,levels=cont_st
  1128 +IF keyword_set(do_postscripts) THEN BEGIN
  1129 + ps_file=ps_dir+use_source_name+'_tau_uv_map'+reso_str+fit_G0_str+mathis_str+'.ps'
  1130 + image_cont20,tau_uv_map,Href,/square,imrange=imrange,axis_color_table=1,image_color_table=image_color_table,/silent,off_bar=obp,title=title,levels=cont_st, $
  1131 + postscript=ps_file,orientation='portrait',/nologo
  1132 + message,'Wrote '+ps_file,/continue
  1133 +ENDIF
  1134 +
  1135 +stop
  1136 +
  1137 +window,win,xsize=wxs,ysize=wys & win=win+1
  1138 +title='tauV/NH'
  1139 +imrange=[0,1]
  1140 +tauv_o_NH_map=la_div(tauv_map,NHmap)
  1141 +image_cont20,tauv_o_NH_map,Href,/square,imrange=imrange,axis_color_table=1,image_color_table=image_color_table,/silent,off_bar=obp,title=title,levels=cont_st
  1142 +IF keyword_set(do_postscripts) THEN BEGIN
  1143 + ps_file=ps_dir+use_source_name+'_tauV_o_NH_map'+reso_str+fit_G0_str+mathis_str+'.ps'
  1144 + image_cont20,tauv_o_NH_map,Href,/square,imrange=imrange,axis_color_table=1,image_color_table=image_color_table,/silent,off_bar=obp,title=title,levels=cont_st, $
  1145 + postscript=ps_file,orientation='portrait',/nologo
  1146 + message,'Wrote '+ps_file,/continue
  1147 +ENDIF
  1148 +
  1149 +
792 1150 window,win,xsize=wxs,ysize=wys & win=win+1
793 1151 imrange=[-1.5,1.2]
794 1152 title='log10(F0_predicted)'
... ... @@ -855,6 +1213,8 @@ CASE use_source_name OF
855 1213 'ngc4321':imrange=[-3.8,-2.9]
856 1214 'ngc3351':imrange=[-3.8,-2.9]-0.3
857 1215 'ngc1566':imrange=[-3.8,-2.9]
  1216 + 'ngc4254':imrange=[-3.8,-2.9]
  1217 + 'ngc3627':imrange=[-3.8,-2.9]
858 1218 ENDCASE
859 1219 title='log10(YG1)'
860 1220 image_cont20,la_log10(Yg1_map),Href,/square,imrange=imrange,axis_color_table=1,image_color_table=image_color_table,/silent,off_bar=obp,title=title,levels=cont_st
... ... @@ -872,7 +1232,7 @@ IF keyword_set(do_postscripts) THEN BEGIN
872 1232 postscript=ps_file,orientation='portrait',/nologo
873 1233 message,'Wrote '+ps_file,/continue
874 1234 ENDIF
875   -stop
  1235 +;stop
876 1236  
877 1237 window,win,xsize=wxs,ysize=wys & win=win+1
878 1238 CASE use_source_name OF
... ... @@ -880,6 +1240,8 @@ CASE use_source_name OF
880 1240 'ngc4321':imrange=[-3.8,-2.9]
881 1241 'ngc3351':imrange=[-3.8,-2.9]-0.3
882 1242 'ngc1566':imrange=[-3.8,-2.9]
  1243 + 'ngc4254':imrange=[-3.8,-2.9]
  1244 + 'ngc3627':imrange=[-3.8,-2.9]
883 1245 ENDCASE
884 1246 title='log10(YG2)'
885 1247 image_cont20,la_log10(Yg2_map),Href,/square,imrange=imrange,axis_color_table=1,image_color_table=image_color_table,/silent,off_bar=obp,title=title,levels=cont_st
... ... @@ -986,6 +1348,8 @@ CASE use_source_name OF
986 1348 'ngc4321':imrange=[0.5-0.2,0.5+0.2]
987 1349 'ngc3351':imrange=[0.5-0.05,0.5+0.05]
988 1350 'ngc1566':imrange=[0.5-0.2,0.5+0.2]
  1351 + 'ngc4254':imrange=[0.5-0.2,0.5+0.2]
  1352 + 'ngc3627':imrange=[0.5-0.2,0.5+0.2]
989 1353 ENDCASE
990 1354 tit='Xion'
991 1355 image_cont20,Xion_map,Href,/square,imrange=imrange,axis_color_table=1,image_color_table=image_color_table,/silent,title=tit,off_bar=obp,levels=cont_st
... ... @@ -1002,7 +1366,7 @@ IF keyword_set(do_postscripts) THEN BEGIN
1002 1366 postscript=ps_file,orientation='portrait',/nologo
1003 1367 message,'Wrote '+ps_file,/continue
1004 1368 ENDIF
1005   -stop
  1369 +;stop
1006 1370  
1007 1371 window,win,xsize=wxs,ysize=wys & win=win+1
1008 1372 CASE use_source_name OF
... ... @@ -1010,6 +1374,8 @@ CASE use_source_name OF
1010 1374 'ngc4321':imrange=[-1.0,1.5]
1011 1375 'ngc3351':imrange=[-0.1,1.5]
1012 1376 'ngc1566':imrange=[-0.1,1.5]
  1377 + 'ngc4254':imrange=[-0.1,1.5]
  1378 + 'ngc3627':imrange=[-0.1,1.5]
1013 1379 ENDCASE
1014 1380 rap_F0=la_div(F0_map,F0_predicted_map)
1015 1381 tit='log(G0/G0_predicted)'
... ... @@ -1022,6 +1388,9 @@ IF keyword_set(do_postscripts) THEN BEGIN
1022 1388 ENDIF
1023 1389 ;stop
1024 1390  
  1391 +
  1392 +win=0L ;avoid "no more windows"
  1393 +
1025 1394 ;==== Show residual maps
1026 1395 window,win,xsize=wxs,ysize=wys & win=win+1
1027 1396 residual_range=[-100,100.]/2.
... ...
LabTools/IRAP/JPB/phangs_make_intensity_mask.pro
  Diff comments   View file @cbdb57b
... ... @@ -77,6 +77,6 @@ IF keyword_set(save) THEN BEGIN
77 77 message,'Wrote '+file_save,/continue
78 78 ENDIF
79 79  
80   -stop
  80 +;stop
81 81  
82 82 END
83 83 \ No newline at end of file
... ...
LabTools/IRAP/JPB/phangs_objectname2reso_str.pro
  Diff comments   View file @cbdb57b
... ... @@ -7,6 +7,8 @@ CASE use_object OF
7 7 'NGC3351': reso_str='1.05asec'
8 8 'NGC4321': reso_str='1.16asec'
9 9 'NGC1566': reso_str='0.80asec'
  10 + 'NGC4254': reso_str='0.89asec'
  11 + 'NGC3627': reso_str='1.05asec'
10 12 ELSE: BEGIN
11 13 reso_str='1asec'
12 14 message,'Using default value of '+reso_str+' for object '+use_object,/continue
... ...
LabTools/IRAP/RM/phangs_fits_seds_sizedist.pro
  Diff comments   View file @cbdb57b
... ... @@ -36,9 +36,9 @@ Nfilters=n_elements(filters)
36 36  
37 37 array_seds = fltarr(Nseds,Nfilters)
38 38  
39   -;fixed_param_desc= ['dustem_plugin_phangs_class_isrf_2'];,'(*!dustem_params).grains['+vsg_id_str+'].MDUST_O_MH']
  39 +fixed_param_desc= ['dustem_plugin_phangs_class_isrf_2'];,'(*!dustem_params).grains['+vsg_id_str+'].MDUST_O_MH']
40 40 ;;fixed_param_val= [1.];,0.0016600000]
41   -;fixed_param_val= [0.81113082];,0.0016600000]
  41 +fixed_param_val= [0.81113082];,0.0016600000]
42 42  
43 43 FOR j=0L,Nseds-1 DO BEGIN
44 44  
... ...
src/idl/dustem_marginalize_grid.pro
  Diff comments   View file @cbdb57b
... ... @@ -87,7 +87,7 @@ CASE use_method OF
87 87 'NEAREST': BEGIN
88 88 FOR i=0L,N_fixed_params-1 DO BEGIN
89 89 ind=where(table_params EQ fixed_parameters_description[i],count,complement=indc,ncomplement=nc)
90   - stop
  90 + ;stop
91 91 IF count NE 0 THEN BEGIN
92 92 ;fixed_parameter_value=fixed_parameters_values[ind[0]]
93 93 fixed_parameter_value=fixed_parameters_values[i]
... ...