Blame view

src/idl/dustem_compute_polext.pro 5.12 KB
07372e26   Ilyes Choubani   Corrected huge pl...
1
FUNCTION dustem_compute_polext ,p_dim,sti,POLEXT_spec,SPEXT_spec,dustem_fpolext,dustem_ext=dustem_ext,out_st=out_st,_extra=extra
e7938fa3   Ilyes Choubani   Corrected02: Impl...
2
3
4
5
6
7
8

;CREATED AS A CONSEQUENCE OF THE INCLUSION OF THE PLUGINS


;written to homogonize DustEM

IF not keyword_set(sti) THEN BEGIN
f2caf3f1   Ilyes Choubani   Debugging some er...
9
10
  dustem_activate_plugins,p_dim/(*(*!dustem_fit).param_init_values),sti
  
e7938fa3   Ilyes Choubani   Corrected02: Impl...
11
12
ENDIF

e7938fa3   Ilyes Choubani   Corrected02: Impl...
13

b74c4452   Ilyes Choubani   Treating extincti...
14
15
POLEXT_spec = sti.polext.ext_tot * (*!dustem_HCD)/1.0e21  ; 

4b745987   Ilyes Choubani   minor update
16
;Hard-coded test to set to potential zero negative values in the polarized extinction arrays:
f2caf3f1   Ilyes Choubani   Debugging some er...
17
18
19
20
21
;so far this has been seen in absorption arrays  
;This test needs to include all grain species that polarize

IF !run_pol THEN BEGIN
    
f2caf3f1   Ilyes Choubani   Debugging some er...
22
23
24
25
26
27
    ind_ngtv_plxt = where (POLEXT_spec LT 0, ct_ngtv_plxt)
   
    IF ct_ngtv_plxt NE 0 THEN (POLEXT_spec)[ind_ngtv_plxt] = 0.
        
        
ENDIF 
b74c4452   Ilyes Choubani   Treating extincti...
28

afde94a3   Ilyes Choubani   general update
29
EXT_spec = sti.ext.ext_tot * (*!dustem_HCD)/1.0e21   ;This is Total intensity I
67bd858a   Ilyes Choubani   Replication of th...
30
31
32
33
34
35
36

out=0.
Nwaves=(size(POLEXT_spec))[1]

frac=POLEXT_spec/EXT_spec
tes=where(finite(frac) eq 0)
frac(tes)=0.
b74c4452   Ilyes Choubani   Treating extincti...
37

afde94a3   Ilyes Choubani   general update
38

b74c4452   Ilyes Choubani   Treating extincti...
39
40
scopes=tag_names((*!dustem_plugin))
IF scopes[0] NE 'NONE' THEN BEGIN
afde94a3   Ilyes Choubani   general update
41
42
43
44
45
46
47
48
;IF ptr_valid(!dustem_plugin) THEN BEGIN
  for i=0L,n_tags(*!dustem_plugin)-1 do begin
    if total(strsplit((*(*!dustem_plugin).(i).scope),'+',/extract) eq 'ADD_EXT') then EXT_spec+=(*(*!dustem_plugin).(i).spec)[*,0]
  endfor
ENDIF


IF scopes[0] NE 'NONE' THEN BEGIN
67bd858a   Ilyes Choubani   Replication of th...
49
50
51
52
53
54
55

    FOR i=0L,n_tags(*!dustem_plugin)-1 DO BEGIN   
        
        IF total(strsplit((*(*!dustem_plugin).(i).scope),'+',/extract) EQ 'REPLACE_POLEXT') THEN BEGIN
            
            QEXT_spec=(*(*!dustem_plugin).(i).spec)[*,1]
            UEXT_spec=(*(*!dustem_plugin).(i).spec)[*,2]
67bd858a   Ilyes Choubani   Replication of th...
56
57
58
            
        ENDIF 
        
887b421b   Ilyes Choubani   correcting some c...
59
60
    ENDFOR   
ENDIF 
67bd858a   Ilyes Choubani   Replication of th...
61
62
63
64
65
    
    IF ~isa(QEXT_spec) && ~isa(UEXT_spec) THEN BEGIN
        polar_ippsi2iqu,EXT_spec,QEXT_spec,UEXT_spec,frac,replicate(!dustem_psi_ext,Nwaves)
        SPEXT_spec = frac 
    ENDIF
887b421b   Ilyes Choubani   correcting some c...
66
67
68

IF scopes[0] NE 'NONE' THEN BEGIN

67bd858a   Ilyes Choubani   Replication of th...
69
70
71
72
73
74
    FOR i=0L,n_tags(*!dustem_plugin)-1 DO BEGIN
        
        IF total(strsplit((*(*!dustem_plugin).(i).scope),'+',/extract) EQ 'ADD_POLEXT') THEN BEGIN
        
            QEXT_spec+=(*(*!dustem_plugin).(i).spec)[*,1]
            UEXT_spec+=(*(*!dustem_plugin).(i).spec)[*,2]
67bd858a   Ilyes Choubani   Replication of th...
75
76
77
78
            
        ENDIF
      
    ENDFOR
67bd858a   Ilyes Choubani   Replication of th...
79

887b421b   Ilyes Choubani   correcting some c...
80
81
82
ENDIF    

POLEXT_spec=sqrt(QEXT_spec^2+UEXT_spec^2)
afde94a3   Ilyes Choubani   general update
83
SPEXT_spec = POLEXT_SPEC/EXT_spec
887b421b   Ilyes Choubani   correcting some c...
84

e7938fa3   Ilyes Choubani   Corrected02: Impl...
85

0068116a   Ilyes Choubani   General update + ...
86
dustem_polext = (*(*!dustem_data).qext).values * 0.
67bd858a   Ilyes Choubani   Replication of th...
87
88
89
90
91
92
93
94
95
96
97
98
99

    
if not isarray(POLEXT_spec) THEN stop ;I don't understand this test. The only thing it can indicate is a problem with the dust parameters or the fortran executable. But in my opinion these tests would be  


;NO COLOR CORRECTION FOR EXTINCTION. FOR NOW ...

; IF !dustem_do_cc NE 0 AND !dustem_never_do_cc EQ 0 THEN BEGIN
;   message,'DOING color correction calculations',/info
; ENDIF ELSE BEGIN
;   message,'SKIPPING color correction calculations',/info
; ENDELSE
; 
0068116a   Ilyes Choubani   General update + ...
100
; ind_polsed=where((*(*!dustem_data).qsed).filt_names NE 'SPECTRUM',count_polsed)  
67bd858a   Ilyes Choubani   Replication of th...
101
102
;   
; IF count_polsed NE 0 THEN BEGIN
0068116a   Ilyes Choubani   General update + ...
103
;   filter_names=((*(*!dustem_data).qsed).filt_names)(ind_polsed)
67bd858a   Ilyes Choubani   Replication of th...
104
105
106
107
108
109
110
111
112
113
;   spolsed=dustem_cc(st.polsed.wav,P_spec,filter_names,cc=cc)
;   dustem_polsed[ind_polsed]=spolsed
; ENDIF


;For spectrum data points, interpolate in log-log.
;Linear interpolation leads to wrong values, in particular where few
;wavelengths points exist in the model (long wavelengths).
;I believe this to be automatic because the provided dustem spectra are already sampled on a log-log grid.

0068116a   Ilyes Choubani   General update + ...
114
115
ind_spec=where((*(*!dustem_data).qext).filt_names EQ 'SPECTRUM',count_spec) 
IF count_spec NE 0 THEN dustem_polext(ind_spec)=interpol(POLEXT_spec,sti.polext.wav,(((*(*!dustem_data).qext).wav)(ind_spec)))
67bd858a   Ilyes Choubani   Replication of th...
116

afde94a3   Ilyes Choubani   general update
117
out_st=sti
67bd858a   Ilyes Choubani   Replication of th...
118
119
120
121
122

;GENERATING THE INTERPOLATES FOR POLFRAC (dustem_polfrac)
if !run_lin then begin
      		
    
afde94a3   Ilyes Choubani   general update
123
    if not keyword_set(dustem_ext) then dustem_ext = dustem_compute_ext(p_dim,sti,EXT_spec) 
67bd858a   Ilyes Choubani   Replication of th...
124
125
126
127
128
129
130
131
132
133
134
135
    
    If n_elements(dustem_ext) ne n_elements(dustem_polext) then begin                            
        if n_elements(dustem_ext) gt n_elements(dustem_polext) then begin 
      
            dustem_polext_x = dustem_polext
            dustem_ext_x = dustem_polext ;meaning dustem_sed needs to be modified
      
            nwaves = n_elements(dustem_polext)
      
      
              for i=0L,nwaves-1 do begin
      
0068116a   Ilyes Choubani   General update + ...
136
                  j=where((*(*!dustem_data).ext).wav EQ (*(*!dustem_data).fpolext).wav(i),testwav)
67bd858a   Ilyes Choubani   Replication of th...
137
138
139
140
141
142
143
144
145
146
147
148
                  if testwav ne 0 then dustem_ext_x(i) = dustem_ext(j(0))                                    
      
              endfor
      
        endif ELSE begin
        
            dustem_polext_x = dustem_ext
            dustem_ext_x = dustem_ext ;meaning dustem_polsed needs to be modified
        
            nwaves = n_elements(dustem_ext)
        
            for i=0L,nwaves-1 do begin
0068116a   Ilyes Choubani   General update + ...
149
                j=where((*(*!dustem_data).ext).wav EQ (*(*!dustem_data).fpolext).wav(i),testwav)
67bd858a   Ilyes Choubani   Replication of th...
150
151
152
153
154
155
156
157
158
159
160
161
162
163
164
165
                if testwav ne 0 then dustem_polext_x(i) = dustem_polext(j(0))                                    
            endfor
          
        
        ENDELSE
    
    endif ELSE BEGIN
        
        dustem_polext_x = dustem_polext 
        dustem_ext_x = dustem_ext
    
    ENDELSE
    
    dustem_fpolext = dustem_polext_x/dustem_ext_x 
	       
ENDIF
e7938fa3   Ilyes Choubani   Corrected02: Impl...
166

e7938fa3   Ilyes Choubani   Corrected02: Impl...
167
168
169
return, dustem_polext


759a527d   Ilyes Choubani   general update
170

e7938fa3   Ilyes Choubani   Corrected02: Impl...
171
END