OV-GSO-DC / dustem-wrapper_idl

Commit ba82c5cfa325020ba327016a6f05bb4b5ef0e507

Authored by Ilyes Choubani
1 parent 7a1d28bf
Exists in master

improved with modification of polarization of spinning dust via plugin

Showing 9 changed files with 136 additions and 104 deletions   Show diff stats
Inline Side-by-side
src/idl/dustem_activate_plugins.pro
  Diff comments   View file @ba82c5c
... ... @@ -90,11 +90,7 @@ FOR i=0L,n_elements(p_min)-1 DO BEGIN
90 90 toto=execute(str)
91 91  
92 92 (*!dustem_scope).(k)=ptr_new(scope)
93   -
94   -
95   - ;print,
96   -
97   -
  93 +
98 94  
99 95 str='((*!dustem_plugin).('+strtrim(k,2)+'))=ptr_new('+ftn+'(key=index,val=value)'+')' & str=str(0)
100 96 toto=execute(str) & IF !dustem_verbose NE 0 THEN message,strupcase(strmid(ftn,7)),/info
... ...
src/idl/dustem_compute_polsed.pro
  Diff comments   View file @ba82c5c
... ... @@ -50,23 +50,76 @@ ENDIF
50 50 If keyword_set(result) THEN result=stp
51 51 fact = 1.e4*1.E20/(4.*!pi)/(3.e8/1.e-6/(stp.polsed).wav)*1.e20/1.e7
52 52 P = stp.polsed.em_tot * fact ;This is Polarized intensity P
53   -
54   -;ADDING PLUGIN TO SPECTRUM----------------
  53 +I = stp.sed.em_tot * fact ;This is Total intensity I
  54 +out=0.
  55 +Nwaves=(size(P))[1]
  56 +
  57 +frac=P/I
  58 +tes=where(finite(frac) eq 0)
  59 +frac(tes)=0.
  60 +
  61 +
  62 +;
  63 +;
  64 +; ;ADDING PLUGIN TO SPECTRUM----------------
  65 +; FOR i=0L,n_tags(*!dustem_scope)-1 DO BEGIN
  66 +; IF total(strsplit((*(*!dustem_scope).(i)),'+',/extract) EQ 'ADD_POLSED') THEN P+=(*(*!dustem_plugin).(i))[*,3] ;3 is for the polarized portion of the plugin because this is polsed
  67 +; ENDFOR
  68 +; ;-----------------------------------------
  69 +;
  70 +;
  71 +; ;MODIFYING PLUGIN----------------
  72 +; FOR i=0L,n_tags(*!dustem_scope)-1 DO BEGIN
  73 +; IF total(strsplit((*(*!dustem_scope).(i)),'+',/extract) EQ 'REPLACE_POLSED') THEN BEGIN
  74 +; P=(*(*!dustem_plugin).(i))[*,3]
  75 +; P+=(*(*!dustem_plugin).(i))[*,3] ;3 is for the polarized portion of the plugin because this is polsed
  76 +; ENDIF
  77 +; ENDFOR
  78 +; ;-----------------------------------------
  79 +
  80 +
  81 +;if Q and U have already been computed elsewhere
  82 +; FOR i=0L,n_tags(*!dustem_scope)-1 DO BEGIN
  83 +;
  84 +; IF total(strsplit((*(*!dustem_scope).(i)),'+',/extract) EQ 'MODIFY_DUST_POLAR') THEN BEGIN
  85 +;
  86 +; Q_spec=(*(*!dustem_plugin).(i))[*,1]
  87 +; U_spec=(*(*!dustem_plugin).(i))[*,2]
  88 +;
  89 +; ENDIF
  90 +; ENDFOR
  91 +
  92 +
  93 +FOR i=0L,n_tags(*!dustem_scope)-1 DO BEGIN
  94 +
  95 + IF total(strsplit((*(*!dustem_scope).(i)),'+',/extract) EQ 'REPLACE_QSED') THEN BEGIN
  96 +
  97 + Q_spec=(*(*!dustem_plugin).(i))[*,1]
  98 + ENDIF
  99 +
  100 + IF total(strsplit((*(*!dustem_scope).(i)),'+',/extract) EQ 'REPLACE_USED') THEN BEGIN
  101 + U_spec=(*(*!dustem_plugin).(i))[*,2]
  102 + ENDIF
  103 +
  104 +ENDFOR
  105 +
  106 +
  107 +IF ~isa(Q_spec) && ~isa(U_spec) THEN polar_ippsi2iqu,stI,Q_spec,U_spec,frac,replicate(0.,Nwaves)
  108 +
55 109 FOR i=0L,n_tags(*!dustem_scope)-1 DO BEGIN
56   - IF total(strsplit((*(*!dustem_scope).(i)),'+',/extract) EQ 'ADD_POLSED') THEN P+=(*(*!dustem_plugin).(i))[*,3] ;3 is for the polarized portion of the plugin because this is polsed
  110 +
  111 + IF total(strsplit((*(*!dustem_scope).(i)),'+',/extract) EQ 'ADD_QSED') THEN BEGIN
  112 + Q_spec+=(*(*!dustem_plugin).(i))[*,1]
  113 + ENDIF
  114 +
  115 + IF total(strsplit((*(*!dustem_scope).(i)),'+',/extract) EQ 'ADD_USED') THEN BEGIN
  116 + U_spec+=(*(*!dustem_plugin).(i))[*,2]
  117 + ENDIF
  118 +
57 119 ENDFOR
58   -;-----------------------------------------
59   -
60 120  
61   -;MODIFYING PLUGIN----------------
62   -FOR i=0L,n_tags(*!dustem_scope)-1 DO BEGIN
63   - IF total(strsplit((*(*!dustem_scope).(i)),'+',/extract) EQ 'REPLACE_POLSED') THEN BEGIN
64   - P=(*(*!dustem_plugin).(i))[*,3]
65   - P+=(*(*!dustem_plugin).(i))[*,3] ;3 is for the polarized portion of the plugin because this is polsed
66   - ENDIF
67   -ENDFOR
68   -;-----------------------------------------
69 121  
  122 +P=sqrt(Q_spec^2+U_spec^2)
70 123  
71 124 dustem_polsed = (*!dustem_data.sed).values * 0.
72 125  
... ...
src/idl/dustem_compute_stokes.pro
  Diff comments   View file @ba82c5c
... ... @@ -61,15 +61,15 @@ tes=where(finite(frac) eq 0)
61 61 frac(tes)=0.
62 62  
63 63 ;if Q and U have already been computed elsewhere
64   -FOR i=0L,n_tags(*!dustem_scope)-1 DO BEGIN
65   -
66   - IF total(strsplit((*(*!dustem_scope).(i)),'+',/extract) EQ 'MODIFY_DUST_POLAR_ANGLE') THEN BEGIN
67   -
68   - Q_spec=(*(*!dustem_plugin).(i))[*,0]
69   - U_spec=(*(*!dustem_plugin).(i))[*,1]
70   -
71   - ENDIF
72   -ENDFOR
  64 +; FOR i=0L,n_tags(*!dustem_scope)-1 DO BEGIN
  65 +;
  66 +; IF total(strsplit((*(*!dustem_scope).(i)),'+',/extract) EQ 'MODIFY_DUST_POLAR') THEN BEGIN
  67 +;
  68 +; Q_spec=(*(*!dustem_plugin).(i))[*,1]
  69 +; U_spec=(*(*!dustem_plugin).(i))[*,2]
  70 +;
  71 +; ENDIF
  72 +; ENDFOR
73 73  
74 74  
75 75 FOR i=0L,n_tags(*!dustem_scope)-1 DO BEGIN
... ... @@ -85,7 +85,6 @@ FOR i=0L,n_tags(*!dustem_scope)-1 DO BEGIN
85 85  
86 86 ENDFOR
87 87  
88   -
89 88  
90 89 IF ~isa(Q_spec) && ~isa(U_spec) THEN polar_ippsi2iqu,stI,Q_spec,U_spec,frac,replicate(0.,Nwaves)
91 90  
... ... @@ -122,10 +121,10 @@ ENDIF ELSE BEGIN
122 121 message,'SKIPPING color correction calculations',/info
123 122 ENDELSE
124 123  
125   -ind_polsed=where((*!dustem_data.polsed).filt_names NE 'SPECTRUM',count_polsed)
  124 +ind_polsed=where((*!dustem_data.sed).filt_names NE 'SPECTRUM',count_polsed)
126 125  
127 126 IF count_polsed NE 0 THEN BEGIN ; INCLUDE COLOR CORRECTIONS OF (Q,U) for dust & synch
128   - filter_names=((*!dustem_data.polsed).filt_names)(ind_polsed)
  127 + filter_names=((*!dustem_data.sed).filt_names)(ind_polsed)
129 128  
130 129 if isa(dustem_qsed) then begin
131 130 sqsed=dustem_cc(stp.sed.wav,Q_spec,filter_names,cc=cc)
... ... @@ -144,10 +143,10 @@ ENDIF
144 143 ;wavelengths points exist in the model (long wavelengths).
145 144  
146 145  
147   -ind_spec=where((*!dustem_data.polsed).filt_names EQ 'SPECTRUM',count_spec)
  146 +ind_spec=where((*!dustem_data.sed).filt_names EQ 'SPECTRUM',count_spec)
148 147 IF count_spec NE 0 THEN BEGIN
149   - if isa(dustem_qsed) then dustem_qsed(ind_spec)=interpol(Q_spec,stp.wav,(((*!dustem_data.polsed).wav)(ind_spec)))
150   - if isa(dustem_used) then dustem_used(ind_spec)=interpol(U_spec,stp.wav,(((*!dustem_data.polsed).wav)(ind_spec)))
  148 + if isa(dustem_qsed) then dustem_qsed(ind_spec)=interpol(Q_spec,stp.wav,(((*!dustem_data.sed).wav)(ind_spec)))
  149 + if isa(dustem_used) then dustem_used(ind_spec)=interpol(U_spec,stp.wav,(((*!dustem_data.sed).wav)(ind_spec)))
151 150 ENDIF
152 151 out_st=stp
153 152  
... ...
src/idl/dustem_fit_sed_polsed_readme.pro
  Diff comments   View file @ba82c5c
... ... @@ -74,18 +74,18 @@ sed.instru=dustem_filter2instru(filters)
74 74 sed.StokesI=sed.wave*0
75 75 sed.StokesQ=sed.wave*0
76 76 sed.StokesU=sed.wave*0
77   -sed.largeP=sed.wave*0
  77 +;sed.largeP=sed.wave*0
78 78 sed.sigmaII=sed.wave*0
79 79 sed.sigmaQQ=sed.wave*0
80 80 sed.sigmaUU=sed.wave*0
81   -sed.sigma_largep=sed.wave*0
  81 +;sed.sigma_largep=sed.wave*0
  82 +
  83 +
82 84 ;st=dustem_set_data(ext=ext,polsed=polsed,polext=polext,sed=sed,polfrac=polfrac,rchi2_weight=rchi2_weight,f_HI=f_HI,qsed=qsed,used=used,qext=qext,uext=uext)
83 85  
84 86 ;INITIALIZING THE !dustem_data tags that will be internally used
85 87 st=dustem_set_data(sed=sed,rchi2_weight=rchi2_weight,f_HI=f_HI)
86 88  
87   -
88   -
89 89 ;st=dustem_set_data() ;this should be one of the tests
90 90  
91 91  
... ... @@ -96,26 +96,18 @@ pd = [ $
96 96 '(*!dustem_params).grains(1).mdust_o_mh',$ ;PAH0 mass fraction
97 97 '(*!dustem_params).grains(2).mdust_o_mh',$ ;PAH1 mass fraction
98 98 '(*!dustem_params).grains(3).mdust_o_mh',$ ;amCBEx
99   - 'dustem_plugin_synchrotron_2', $ ;Synchrotron amplitude
100   - 'dustem_plugin_synchrotron_4', $ ;Synchrotron polarization angle
101   - 'dustem_plugin_modify_dust_polar_angle_1', $ ;Polarization angle applied to the Fortran dust model via a core plugin
102   - 'dustem_plugin_modify_spinning_sed_1', $
103   - 'dustem_plugin_modify_spinning_sed_2' $
  99 + ;'dustem_plugin_synchrotron_2', $ ;Synchrotron amplitude
  100 + ;'dustem_plugin_synchrotron_4', $ ;Synchrotron polarization angle
  101 + 'dustem_plugin_modify_dust_polarization_2', $ ;Polarization angle applied to the Fortran dust model via a core plugin
  102 + 'dustem_plugin_modify_spinning_polarization_1', $
  103 + 'dustem_plugin_modify_spinning_polarization_2' $
104 104 ]
105 105  
106   -p_truth=[1.,7.8000E-04,7.8000E-04,7.8000E-04,1.6500E-04,5e-3,20.,68.,1e-3,20]
107   -;p_truth=[1.,7.8000E-04,7.8000E-04,7.8000E-04,1.6500E-04,5e-3,20.,0.2,20]
108   -
109   -;=== set initial value of parameters you want to fit
110   -;iv=p_truth ;exact solution
111   -
  106 +;p_truth=[1.,7.8000E-04,7.8000E-04,7.8000E-04,1.6500E-04,5e-3,20.,68.,1.2,20]
  107 +p_truth=[1.,7.8000E-04,7.8000E-04,7.8000E-04,1.6500E-04,68.,4e-3,20]
  108 +iv=p_truth+[0.,8.e-4,8.e-4,8.e-4,8.e-4,-10.,2e-4,10] ;shifted from solution
112 109  
113   -;UNCOMMENT THIS BLOCK AFTER THE TESTING
114 110  
115   -iv=p_truth+[0.,8.e-4,8.e-4,8.e-4,8.e-4,8e-3,-8.,-10.,1.2,10] ;shifted from solution
116   -;iv=p_truth+[0.,8.e-4,8.e-4,8.e-4,8.e-4,8e-3,-8.,1.2,10] ;shifted from solution
117   -
118   -;TESTING THE ISRF PLUGIN
119 111  
120 112 Npar=n_elements(pd)
121 113 ulimed=replicate(0,Npar)
... ... @@ -142,8 +134,8 @@ sed.sigmaII=abs(sed.StokesI*0.000001)
142 134 ;stop
143 135  
144 136  
145   -sed.largeP=dustem_compute_polsed(p_truth,sstp)
146   -sed.sigma_largep=abs(sed.largeP*0.000001)
  137 +;sed.largeP=dustem_compute_polsed(p_truth,sstp)
  138 +;sed.sigma_largep=abs(sed.largeP*0.000001)
147 139  
148 140  
149 141 toto=dustem_compute_stokes(p_truth,sstqu,dustem_qsed=dustem_qsed,dustem_used=dustem_used) ;this procedure also allows for the extraction of the spectra
... ... @@ -154,25 +146,12 @@ sed.sigmaQQ=abs(sed.StokesQ*0.000001)
154 146 sed.StokesU=dustem_used
155 147 sed.sigmaUU=abs(sed.StokesU*0.000001)
156 148  
157   -
158   -; ;SET qext and uext data - comment if polarizarion data is to be used
159   -; qext.spec=dustem_qext
160   -; qext.error=qext.spec*0.1
161   -;
162   -;
163   -; uext.spec=dustem_uext
164   -; uext.error=uext.spec*0.1
165   -; ;
166   -
167 149 ;=== SET THE OBSERVATION STRUCTURE
168 150  
169 151 ind=where(sed.wave GE 0.1)
170 152  
171 153 st_bidon=dustem_set_data(sed=sed[ind],rchi2_weight=rchi2_weight,f_HI=f_HI) ; comment if second line is to be used
172 154  
173   -
174   -
175   -
176 155 ;=== RUN fit
177 156 tol=1.e-16
178 157 xtol=1.e-16
... ...
src/idl/dustem_mpfit_run.pro
  Diff comments   View file @ba82c5c
... ... @@ -388,7 +388,7 @@ FOR i_tag=0,count_data_tag-1 DO BEGIN
388 388 titstq='Stokes Q Polarization Intensity'
389 389 ytitstq=textoidl('\nuQ_\nu/N_H (W/m^2/sr/H)')
390 390 xtit=textoidl('\lambda (\mum)')
391   - xr=[10,8e3]
  391 + xr=[10,6e4]
392 392 normpol=Q_spec*fact
393 393 normpol1=dustem_qsed*fact1
394 394 ;normpol1=(*!dustem_data.qsed).values*fact1
... ... @@ -439,12 +439,12 @@ FOR i_tag=0,count_data_tag-1 DO BEGIN
439 439 titstu='Stokes U Polarization Intensity'
440 440 ytitstu=textoidl('\nuU_\nu/N_H (W/m^2/sr/H)')
441 441 xtit=textoidl('\lambda (\mum)')
442   - xr=[10,8e3]
  442 + xr=[10,6e4]
443 443 normpol2=U_spec*fact
444 444 normpol3=dustem_used*fact1
445 445 ;yr=[min(normpol2)-min(normpol2)/10,max(normpol2)+max(normpol2)*10]
446 446  
447   - yr=[1e-34,1e-28]
  447 + yr=[1e-39,1e-28]
448 448 indu=where(U_spec lt 0, countu)
449 449 ylog=1
450 450 if countu ne 0 then begin
... ...
src/idl/dustem_plot_polar.pro
  Diff comments   View file @ba82c5c
... ... @@ -320,8 +320,11 @@ IF ptr_valid(!dustem_data.polsed) THEN BEGIN
320 320  
321 321 ;if total(counttx) then frac_synch=p_dim(indx) else frac_synch=0.3 ;THIS LINE IS TO TAKE INTO ACCOUNT WHEN THE AMPLITUDE OF THE SYNCHROTRON IS FROZEN
322 322  
323   - cgoplot,st.polsed.wav,(*(*!dustem_plugin).synchrotron)[*,3]*(3.e8/(1.e-6*(st.polsed.wav)))/1.d40,color='Crimson',linestyle=3
324   - cgoplot,st.polsed.wav,st.polsed.em_tot*1./(4.*!pi)*1.e17/1d20+((*(*!dustem_plugin).synchrotron)[*,3])*(3.e8/(1.e-6*(st.polsed.wav)))/1.d40,color='black'
  323 + synchrotron_p=sqrt((*(*!dustem_plugin).synchrotron)[*,1]^2+(*(*!dustem_plugin).synchrotron)[*,2]^2)
  324 +
  325 +
  326 + cgoplot,st.polsed.wav,synchrotron_p*(3.e8/(1.e-6*(st.polsed.wav)))/1.d40,color='Crimson',linestyle=3
  327 + cgoplot,st.polsed.wav,st.polsed.em_tot*1./(4.*!pi)*1.e17/1d20+synchrotron_p*(3.e8/(1.e-6*(st.polsed.wav)))/1.d40,color='black'
325 328  
326 329 ENDIF
327 330  
... ...
src/idl/dustem_plugin_modify_dust_polar_angle.pro renamed to src/idl/dustem_plugin_modify_dust_polarization.pro
  Diff comments   View file @ba82c5c
1   -Function dustem_plugin_modify_dust_polar_angle, key=key, val=val, scope=scope,help=help
  1 +Function dustem_plugin_modify_dust_polarization, key=key, val=val, scope=scope,help=help
2 2  
3 3 ;+
4 4 ; NAME:
... ... @@ -38,7 +38,8 @@ IF keyword_set(help) THEN BEGIN
38 38 ENDIF
39 39  
40 40  
41   -scope='MODIFY_DUST_POLAR_ANGLE'
  41 +scope='REPLACE_QSED+REPLACE_USED'
  42 +;scope=''
42 43  
43 44  
44 45 IF ptr_valid((*!dustem_fit).CURRENT_PARAM_VALUES) THEN BEGIN
... ... @@ -51,13 +52,13 @@ ENDIF ELSE BEGIN
51 52  
52 53 ENDELSE
53 54  
54   -
55   -
56 55 psi=0.
  56 +smallp=1.
57 57 IF keyword_set(key) THEN BEGIN
58   - a=where(key EQ 1,count1)
59   - IF count1 NE 0 THEN psi=val[a[0]] ; setting psi from pd
60   -
  58 + ind1=where(key EQ 1,count1)
  59 + ind2=where(key EQ 2,count2)
  60 + IF count1 NE 0 THEN smallp=val[ind1[0]] ; setting psi from pd
  61 + IF count2 NE 0 THEN psi=val[ind2[0]] ; setting psi from pd
61 62 ENDIF
62 63  
63 64  
... ... @@ -67,17 +68,18 @@ I=st.sed.em_tot*fact
67 68  
68 69 Nwaves=(size(I))[1]
69 70  
70   -frac=P/I
  71 +frac=P/I*smallp
71 72 tes=where(finite(frac) eq 0)
72 73 frac(tes)=0.
73 74  
74 75  
75 76 polar_ippsi2iqu,I,Q,U,frac,replicate(psi,Nwaves)
76 77  
77   -out=fltarr(Nwaves,2)
  78 +out=fltarr(Nwaves,3)
78 79  
79   -out[*,0]=Q
80   -out[*,1]=U
  80 +out[*,0]=I
  81 +out[*,1]=Q
  82 +out[*,2]=U
81 83  
82 84 return, out
83 85  
... ...
src/idl/dustem_plugin_modify_spinning_sed.pro renamed to src/idl/dustem_plugin_modify_spinning_polarization.pro
  Diff comments   View file @ba82c5c
1   -FUNCTION dustem_plugin_modify_spinning_sed, key=key, val=val, scope=scope, help=help
  1 +FUNCTION dustem_plugin_modify_spinning_polarization, key=key, val=val, scope=scope, help=help
2 2  
3 3 ;+
4 4 ; NAME:
... ... @@ -32,12 +32,12 @@ FUNCTION dustem_plugin_modify_spinning_sed, key=key, val=val, scope=scope, help=
32 32 ;-
33 33  
34 34 IF keyword_set(help) THEN BEGIN
35   - doc_library,'dustem_plugin_modify_spinning_sed'
  35 + doc_library,'dustem_plugin_modify_spinning_polarization'
36 36 goto,the_end
37 37 ENDIF
38 38  
39 39  
40   -scope='REPLACE_POLSED+REPLACE_QSED+REPLACE_USED'
  40 +scope='ADD_QSED+ADD_USED'
41 41  
42 42 IF ptr_valid((*!dustem_fit).CURRENT_PARAM_VALUES) THEN BEGIN
43 43 p_dii=(*(*!dustem_fit).CURRENT_PARAM_VALUES)
... ... @@ -45,37 +45,44 @@ IF ptr_valid((*!dustem_fit).CURRENT_PARAM_VALUES) THEN BEGIN
45 45  
46 46 ENDIF ELSE BEGIN
47 47 p_dii=(*(*!dustem_fit).PARAM_INIT_VALUES)
48   - st=dustem_run(p_di)
  48 + st=dustem_run(p_dii)
49 49  
50 50 ENDELSE
51 51  
  52 +;stop
52 53  
53 54 psi=0.
54 55 smallp=0. ;supposing the spinning dust is
55 56 IF keyword_set(key) THEN BEGIN
  57 +
56 58 ind1=where(key EQ 1,count1)
57 59 ind2=where(key EQ 2,count2)
58   - IF count1 NE 0 THEN smallp=val[ind1[0]] ; setting psi from pd
  60 + IF count1 NE 0 THEN smallp=val[ind1[0]] ; setting smallp from pd
59 61 IF count2 NE 0 THEN psi=val[ind2[0]] ; setting psi from pd
  62 +
60 63 ENDIF
61 64  
62 65  
63   -fact = 1.e4*1.E20/(4.*!pi)/(3.e8/1.e-6/(st.polsed).wav)*1.e20/1.e7
  66 +fact = 1.e4*1.E20/(4.*!pi)/(3.e8/1.e-6/(st.sed).wav)*1.e20/1.e7
64 67  
65 68  
66 69 SpinningI=(st.sed.em_grain_1)+(st.sed.em_grain_2)*fact
67   -SpinningP=SpinningI*smallp
  70 +
  71 +;760
  72 +SpinningI[0:760]=0.
  73 +
68 74  
69 75 Nwaves=(size(SpinningI))[1]
70 76  
  77 +
71 78 polar_ippsi2iqu,SpinningI,SpinningQ,SpinningU,replicate(smallp,Nwaves),replicate(psi,Nwaves)
72 79  
73   -out=fltarr(Nwaves,4)
  80 +out=fltarr(Nwaves,3)
74 81  
75 82 out[*,0]=SpinningI
76 83 out[*,1]=SpinningQ
77 84 out[*,2]=SpinningU
78   -out[*,3]=SpinningP
  85 +
79 86  
80 87 return, out
81 88  
... ...
src/idl/dustem_plugin_synchrotron.pro
  Diff comments   View file @ba82c5c
... ... @@ -56,12 +56,7 @@ IF keyword_set(key) THEN BEGIN
56 56 IF count1 NE 0 then s=val[ind1[0]] ;else s=3
57 57 IF count2 NE 0 then A=val[ind2[0]] ;else A=1
58 58 IF count3 NE 0 then smallp=val[ind3[0]]
59   - IF count4 NE 0 then BEGIN
60   - psi=val[ind4[0]]
61   - print, 'psi in synchro plugin is:'
62   - print,psi
63   - ENDIF
64   -
  59 + IF count4 NE 0 then psi=val[ind4[0]]
65 60  
66 61 ENDIF
67 62  
... ... @@ -73,7 +68,7 @@ nu=cmic/lambir
73 68 lambir_ref=10000.
74 69  
75 70 ;see Deschenes et al 2008, eq 6
76   -output=fltarr(Nwavs,4)
  71 +output=fltarr(Nwavs,3)
77 72 output[*,0]=nu^(-(s+3.)/2.)
78 73 norm=interpol(output[*,0],lambir,lambir_ref)
79 74 output[*,0]=A*output[*,0]/norm
... ... @@ -83,11 +78,9 @@ polar_ippsi2iqu,output[*,0],Q,U,replicate(smallp,Nwavs),replicate(psi,Nwavs)
83 78  
84 79 output[*,1]=Q
85 80 output[*,2]=U
86   -output[*,3]=output[*,0]*smallp ;polarized part of the synchrotron emission
87 81  
88   -scope='ADD_SED+ADD_POLSED+ADD_QSED+ADD_USED'
89   -;scope=((*!dustem_scope).synchrotron)
90   -;if keyword_set(scope) then
  82 +
  83 +scope='ADD_SED+ADD_QSED+ADD_USED'
91 84  
92 85 the_end:
93 86 RETURN,output
... ...