diff --git a/LabTools/IRAP/JPB/srun/make_phangs_grids.pro b/LabTools/IRAP/JPB/srun/make_phangs_grids.pro index 8a9e2a6..31e9847 100644 --- a/LabTools/IRAP/JPB/srun/make_phangs_grids.pro +++ b/LabTools/IRAP/JPB/srun/make_phangs_grids.pro @@ -13,6 +13,9 @@ PRO make_phangs_grids,bidon ;more {$phangs_srundir}/make_phangs_grids.log | grep "iteration " ;more {$phangs_srundir}/make_phangs_grids.log | grep "Will produce" +;This won't work, because it needs the Fortran to set files to the same area ... +;defsysv,'!DUSTEM_RES','/tmp/jbernard/Dustem_nohup' + t0=systime(/sec) ;make grid diff --git a/plugin_scopes_definition.xcat b/plugin_scopes_definition.xcat index 2013490..9c9d81d 100644 --- a/plugin_scopes_definition.xcat +++ b/plugin_scopes_definition.xcat @@ -13,13 +13,13 @@ \------------------------------------------------ \This scope is for plugins adding to the output SEDs ADD SED 1 -ADD POLSED 1 +\ADD POLSED 1 \This scope is for plugins replacing the ISRF before fortran is run REPLACE ISRF 2 ADD ISRF 2 \This scope is for plugins replacing the fortran SEDs REPLACE SED 3 -REPLACE POLSED 3 +\REPLACE POLSED 3 REPLACE POLEXT 3 \This scope is for plugins which fully replace the Fortran, including by calling Fortran several times ... REPLACE FORTRAN 4 diff --git a/src/idl/dustem_activate_plugins.pro b/src/idl/dustem_activate_plugins.pro index c810bad..c16600c 100755 --- a/src/idl/dustem_activate_plugins.pro +++ b/src/idl/dustem_activate_plugins.pro @@ -78,36 +78,27 @@ ENDELSE ;stop -;Run the plugins, according to their scopes +;===== Run the plugins, according to their scopes (and the Fortran) -;normal (additive) plugins come first as they can be used by other plugins -;The ISRF plugin is next - we can always change its scope however we please to. DustEm runs here. This is with the /force_dustem_run keyword. -;The plugins that use the DustEM output come next - we can always change their 'common' scope +;The 'ADD_SED' plugins are run first. Note that their output could be used by plugins run later. +;The plugins related to 'ISRF' are run second, and the Fortrn is run at that stage, just after the corresponding plugins. +;The plugins which Replace the output of the fortran are run last. ;IF /avoid is used the specified scopes are avoided. They're considered when /avoid is missing. ;using '**' so not to write the entire scope -;=== This runs plugins with scopes including STELLAR, REPLACE and ISRF -;dustem_run_plugins, p_dim ,param_descs, param_values, param_func, ['ADD','REPLACE'], avoid=1,use_previous_fortran=use_previous_fortran -dustem_run_plugins, p_dim ,param_descs, param_values, param_func, ['ADD_SED'], avoid=1,use_previous_fortran=use_previous_fortran +;=== This runs the plugins +;dustem_run_plugins, p_dim ,param_descs, param_values, param_func, ['ADD_SED'], avoid=1,use_previous_fortran=use_previous_fortran +dustem_run_plugins, p_dim ,param_descs, param_values, param_func, 1, avoid=1,use_previous_fortran=use_previous_fortran -;=== The line below makes fawlty go Segmentation fault. st is undefined at that stage -;=== it would run fine if avoid was set to 1 -;IC: I think the segmentation fault comes from the fact that the user wants to use models that have the /pol keyword in the GRAIN.DAT file without being in polarization mode. +;dustem_run_plugins, p_dim ,param_descs, param_values, param_func, ['ADD_ISRF','REPLACE_ISRF'] ,force_dustem_run=1,st=st,use_previous_fortran=use_previous_fortran ;dustem output is available at this step. +dustem_run_plugins, p_dim ,param_descs, param_values, param_func, 2,force_dustem_run=1,st=st,use_previous_fortran=use_previous_fortran ;dustem output is available at this step. -;dustem_run_plugins, p_dim ,param_descs, param_values, param_func, ['STELLAR','ISRF'] ,force_dustem_run=1,st=st,use_previous_fortran=use_previous_fortran ;dustem output is available at this step. -dustem_run_plugins, p_dim ,param_descs, param_values, param_func, ['ADD_ISRF','REPLACE_ISRF'] ,force_dustem_run=1,st=st,use_previous_fortran=use_previous_fortran ;dustem output is available at this step. +;dustem_run_plugins, p_dim ,param_descs, param_values, param_func, ['REPLACE_SED'],use_previous_fortran=use_previous_fortran +dustem_run_plugins, p_dim ,param_descs, param_values, param_func, 3,use_previous_fortran=use_previous_fortran -;added this small condition to ensure that dustem is run -;(if no ISRF plugin is used) -; if ~isa(st) then begin -; st=dustem_run(p_dim) -; !dustem_current=ptr_new(st) -; endif -;THE ABOVE CONDITION IS UNNECESSARY SINCE DUSTEM IS CALLED ABOVE WHETHER PLUGINS ARE RUN OR NOT. - -;dustem_run_plugins, p_dim ,param_descs, param_values, param_func, ['REPLACE'],use_previous_fortran=use_previous_fortran -dustem_run_plugins, p_dim ,param_descs, param_values, param_func, ['REPLACE_SED'],use_previous_fortran=use_previous_fortran +;dustem_run_plugins, p_dim ,param_descs, param_values, param_func, ['REPLACE_SED'],use_previous_fortran=use_previous_fortran +dustem_run_plugins, p_dim ,param_descs, param_values, param_func, 4,use_previous_fortran=use_previous_fortran the_end: diff --git a/src/idl/dustem_compute_stokes.pro b/src/idl/dustem_compute_stokes.pro index ab32b3d..46c3aae 100755 --- a/src/idl/dustem_compute_stokes.pro +++ b/src/idl/dustem_compute_stokes.pro @@ -75,18 +75,25 @@ P = st.polsed.em_tot * fact ;This is Polarized intensity P Nwaves=(size(stI))[1] -frac=P/stI +;stop + +;This some times leads to division by 0 +frac=P/stI tes=where(finite(frac) eq 0) frac(tes)=0. -scopes=tag_names((*!dustem_plugin)) +scopes=(*!dustem_plugin).scope +n_plugins=n_elements(scopes) + IF scopes[0] NE 'NONE' THEN BEGIN - FOR i=0L,n_tags(*!dustem_plugin)-1 DO BEGIN - IF total(strsplit((*(*!dustem_plugin).(i).scope),'+',/extract) EQ 'REPLACE_POLSED') THEN BEGIN - Q_spec = (*(*!dustem_plugin).(i).spec)[*,1] - U_spec = (*(*!dustem_plugin).(i).spec)[*,2] - ENDIF - ENDFOR + FOR i=0L,n_plugins-1 DO BEGIN + ;IF total(strsplit((*(*!dustem_plugin).(i).scope),'+',/extract) EQ 'REPLACE_POLSED') THEN BEGIN + ;IF total(strsplit((*!dustem_plugin)[i].scope,'+',/extract) EQ 'REPLACE_POLSED') THEN BEGIN + IF total(strsplit((*!dustem_plugin)[i].scope,'+',/extract) EQ 'REPLACE_SED') THEN BEGIN + Q_spec = (*(*!dustem_plugin)[i].spec)[*,1] + U_spec = (*(*!dustem_plugin)[i].spec)[*,2] + ENDIF + ENDFOR ENDIF IF ~isa(Q_spec) && ~isa(U_spec) THEN BEGIN @@ -94,18 +101,14 @@ ENDIF ENDIF IF scopes[0] NE 'NONE' THEN BEGIN - FOR i=0L,n_tags(*!dustem_plugin)-1 DO BEGIN - - IF total(strsplit((*(*!dustem_plugin).(i).scope),'+',/extract) EQ 'ADD_POLSED') THEN BEGIN - - Q_spec+=(*(*!dustem_plugin).(i).spec)[*,1] - U_spec+=(*(*!dustem_plugin).(i).spec)[*,2] - - + FOR i=0L,n_plugins-1 DO BEGIN + ;IF total(strsplit((*(*!dustem_plugin).(i).scope),'+',/extract) EQ 'ADD_POLSED') THEN BEGIN + ;IF total(strsplit((*!dustem_plugin)[i].scope,'+',/extract) EQ 'ADD_POLSED') THEN BEGIN + IF total(strsplit((*!dustem_plugin)[i].scope,'+',/extract) EQ 'ADD_SED') THEN BEGIN + Q_spec+=(*(*!dustem_plugin)[i].spec)[*,1] + U_spec+=(*(*!dustem_plugin)[i].spec)[*,2] ENDIF - ENDFOR - ENDIF ;----------------------------------------- ind=where(finite(U_spec) NE 1 OR finite(Q_spec) NE 1,count) diff --git a/src/idl/dustem_fit_intensity_mbb_example.pro b/src/idl/dustem_fit_intensity_mbb_example.pro index 1acdaaf..e695a0f 100644 --- a/src/idl/dustem_fit_intensity_mbb_example.pro +++ b/src/idl/dustem_fit_intensity_mbb_example.pro @@ -162,7 +162,7 @@ if keyword_set(wait) then begin end ;== INITIALISE DUSTEM -dustem_init,model=use_model,polarization=use_polarization +dustem_init,model=use_model,polarization=use_polarization,/show_plots !dustem_nocatch=1 !dustem_verbose=use_verbose IF keyword_set(noobj) THEN !dustem_noobj=1 diff --git a/src/idl/dustem_init_plugins.pro b/src/idl/dustem_init_plugins.pro index b27f71e..d93e619 100644 --- a/src/idl/dustem_init_plugins.pro +++ b/src/idl/dustem_init_plugins.pro @@ -88,7 +88,13 @@ plugin_names=plugin_names[un] parameter_types=parameter_types[un] ind=where(parameter_types EQ 'PLUGIN',Nplugins) -one_plugin_st={name:'',spec:ptr_new(),scope:'',paramtag:ptr_new()} +;===== This is the dustem structure describing the plugin scopes +one_plugin_st={name:'', $ ;name of the plugin + spec:ptr_new(), $ ;output of the plugin (a spectrum) + scope:'', $ ;scope of the plugin + run_order:0L, $ ;plugin run order + paramtag:ptr_new() $ ;parameter names + } IF Nplugins EQ 0 THEN BEGIN plugin_st=replicate(one_plugin_st,1) @@ -114,6 +120,18 @@ FOR i=0L,Nplugins-1 DO BEGIN toto=execute(str) plugin_st[i].paramtag=ptr_new(paramtag) ENDFOR +;derive the run order from the scopes +file_scopes=!dustem_wrap_soft_dir+'plugin_scopes_definition.xcat' +scope_st=read_xcat(file_scopes,/silent) +possible_scopes=scope_st.objective+'_'+scope_st.target +FOR i=0L,Nplugins-1 DO BEGIN + ind=where(possible_scopes EQ plugin_st[i].scope,count) + IF count EQ 0 THEN BEGIN + message,'scope '+plugin_st[i].scope+' not recognized',/continue + stop + ENDIF + plugin_st[i].run_order=scope_st[ind[0]].run_order +ENDFOR ;creates the !dustem_plugin structure defsysv, '!dustem_plugin', ptr_new(plugin_st) diff --git a/src/idl/dustem_plugin_freefree.pro b/src/idl/dustem_plugin_freefree.pro index e1fefce..fdc8d79 100755 --- a/src/idl/dustem_plugin_freefree.pro +++ b/src/idl/dustem_plugin_freefree.pro @@ -23,8 +23,11 @@ FUNCTION dustem_plugin_freefree,key=key $ ; ; OPTIONAL INPUT PARAMETERS: ; key = input parameter number +; 1 : Tgas [K] +; 2 : free-free Amplitude [] +; 3 : polarization fraction [%] (default=0) +; 4 : polarization angle [deg] ; val = input parameter value -; ; OUTPUTS: ; freefree = free-free spectrum (on dustem wavelengths) ; @@ -77,29 +80,21 @@ IF keyword_set(key) THEN BEGIN b=where(key EQ 2,count2) c=where(key EQ 3,count3) ;default polarization fraction -newly added d=where(key EQ 4,count4) ;default polarization angle -newly added - IF count1 NE 0 then Tgas=(val(a))(0) - IF count2 NE 0 then Amplitude=(val(b))(0) - IF count3 NE 0 then smallp=val[c[0]] ;-newly added - IF count4 NE 0 then psi=val[d[0]] ;-newly added + IF count1 NE 0 then Tgas=val[a[0]] + IF count2 NE 0 then Amplitude=val[b[0]] + IF count3 NE 0 then smallp=val[c[0]] + IF count4 NE 0 then psi=val[d[0]] ENDIF -;IF !dustem_which EQ 'DESERT' THEN BEGIN -; lambir=((*!dustem_params).gemissiv.lambir) -;ENDIF ELSE BEGIN -; lambir=((*!dustem_params).lambda.lambda) -;ENDELSE -;replaced by the following to ease the life of plugins writters lambir=dustem_get_wavelengths() Nwavs=n_elements(lambir) cmic=3.e14 nu=cmic/lambir ;Hz mjy=1 ;output is in MJy/sr -lambir_ref=10000. - +lambir_ref=10000. ;in mic (this is 1cm) -output=fltarr(Nwavs,3) ;newly added +output=fltarr(Nwavs,3) -;stop ;use_method='Deschenes2008' ;use_method='WallsGabaud1998' use_method='Dickinson2003_norm' @@ -140,10 +135,7 @@ ENDCASE polar_ippsi2iqu,output[*,0],Q,U,replicate(smallp,Nwavs),replicate(psi,Nwavs) output[*,1]=Q -output[*,2]=U -;output[*,3]=replicate(smallp,Nwavs) -;output[*,4]=replicate(psi,Nwavs) - +output[*,2]=U the_end: RETURN,output diff --git a/src/idl/dustem_plugin_mbbdy.pro b/src/idl/dustem_plugin_mbbdy.pro index 096af8c..5c9c179 100755 --- a/src/idl/dustem_plugin_mbbdy.pro +++ b/src/idl/dustem_plugin_mbbdy.pro @@ -52,7 +52,8 @@ w0=100. ; Theoretical wavelength (3THz) where optical depth equals unity smallp=0. psi=0. -scope='ADD_SED+ADD_POLSED' +;scope='ADD_SED+ADD_POLSED' +scope='ADD_SED' paramtag=['Amp',textoidl('T_{MBB}')+' [K]',textoidl('\beta'),'p','Psi [deg]'] IF keyword_set(key) THEN BEGIN diff --git a/src/idl/dustem_plugin_modify_dust_pol.pro b/src/idl/dustem_plugin_modify_dust_pol.pro index e7c470e..cd7d54b 100644 --- a/src/idl/dustem_plugin_modify_dust_pol.pro +++ b/src/idl/dustem_plugin_modify_dust_pol.pro @@ -5,31 +5,33 @@ FUNCTION dustem_plugin_modify_dust_pol, key=key, val=val, scope=scope, paramtag= ; dustem_plugin_modify_dust_pol ; ; PURPOSE: -; modifies Emission Stokes Q,U values for dust polarization in st according to the given keywords and values -; +; Replaces the emission in polarization (Stokes Q,U) with a fixed polarization angle and fraction +; Q,U to be replaced are computed from the total intensity as: +; Q=I*val[0]*p*cos(2*val[1]) +; U=I*val[0]*p*sin(2*val[1]) +; where I is the total intensity in the model and p is the polarization fraction in the model ; CATEGORY: ; DustEMWrap, Plugin, Mid-level, Distributed ; ; CALLING SEQUENCE: -; a=dustem_plugin_modify_dust_pol(st,[key=][val=][,scope=][,paramtag=][,paramdefault=][,/help]) +; a=dustem_plugin_modify_dust_pol([key=][val=][,scope=][,paramtag=][,paramdefault=][,/help]) ; ; INPUTS: ; st = dustem structure ; ; OPTIONAL INPUT PARAMETERS: -; key = input parameter numbers (first = polarization fraction in %, default=1.%, second=polarization angle, default=0.) +; key = input parameter numbers +; first = multiplication factor for the polarization fraction +; second = polarization angle [deg], default=0. ; val = input parameter values -; ; OUTPUTS: -; out = array containing the stokes emission parameters associated to the dust/synchrotron component or a concatenated version for both +; out = array containing the stokes emission parameters associated to the dust/synchrotron component or a concatenated version for both ; ; OPTIONAL OUTPUT PARAMETERS: ; scope = if set, returns only the scope of the pluggin ; paramtag = if set, returns only the parameter tags -; ; ACCEPTED KEY-WORDS: ; help = if set, print this help -; ; COMMON BLOCKS: ; None ; @@ -53,20 +55,25 @@ IF keyword_set(help) THEN BEGIN ENDIF IF keyword_set(scope) THEN BEGIN + ;scope='REPLACE_POLSED' + scope='REPLACE_SED' out=0 GOTO, the_end ENDIF IF keyword_set(paramtag) THEN BEGIN + paramtag=['p',textoidl('\psi')+' [deg]'] out=0. GOTO, the_end ENDIF ;Retrieving this sytem variable (dustem output) -if ~isa(!dustem_current) then begin - out=0 +IF ~isa(!dustem_current) then begin + out=0 goto, the_end -ENDIF else st=(*!dustem_current) +ENDIF ELSE BEGIN + st=(*!dustem_current) +ENDELSE ;below are the default values for the plugin parameters smallp_fact=1. ;This is the default dust polarization factor @@ -76,12 +83,13 @@ IF keyword_set(key) THEN BEGIN ind1=where(key EQ 1,count1) ind2=where(key EQ 2,count2) IF count1 NE 0 THEN smallp_fact=val[ind1[0]] ; setting smallp from pd - this is another polarization fraction (constant) that is applied to the total dust emission - IF count2 NE 0 THEN psi=val[ind2[0]] & !dustem_psi = psi ; setting psi from pd. !dustem_psi here helps for the plotting. + IF count2 NE 0 THEN psi=val[ind2[0]] + !dustem_psi = psi ; setting psi from pd. !dustem_psi here helps for the plotting. ENDIF fact = 1.e4*(*!dustem_HCD)/(4.*!pi)/(3.e8/1.e-6/((st.polsed).wav))*1.e20/1.e7 -P=((st.polsed).em_tot)*fact ; This is the polarized emission P -I=((st.sed).em_tot)*fact ; This is the total intensity emission I +P=((st.polsed).em_tot)*fact ; This is the polarized emission P currently in the model +I=((st.sed).em_tot)*fact ; This is the total intensity emission I currently in the model Nwaves=(size(I))[1] @@ -90,7 +98,7 @@ Nwaves=(size(I))[1] frac_model=P*0. indI = where(I ne 0, countI) -if countI ne 0 then frac_model[indI] = P[indI]/I[indI] ;This is the polarization fraction in the model +if countI ne 0 then frac_model[indI] = P[indI]/I[indI] ;This is the polarization fraction currently in the model frac_used=frac_model*smallp_fact @@ -98,22 +106,15 @@ psi_used = replicate(psi,Nwaves) polar_ippsi2iqu,I,Q,U,frac_used,psi_used -out=fltarr(Nwaves,3) ; modified this. This is the only plugin that has this number of outputs. - -;degtorad = !pi/180 +out=fltarr(Nwaves,3) out[*,0]=I out[*,1]=Q out[*,2]=U -; out[*,3]=frac_used -; ;out[*,4]=0.5*atan(U,Q)/degtorad -; out[*,4]=psi_used RETURN, out the_end: -scope='REPLACE_POLSED' -paramtag=['p',textoidl('\psi')+' [deg]'] END diff --git a/src/idl/dustem_plugin_synchrotron.pro b/src/idl/dustem_plugin_synchrotron.pro index df1ff0d..bdd56fc 100755 --- a/src/idl/dustem_plugin_synchrotron.pro +++ b/src/idl/dustem_plugin_synchrotron.pro @@ -64,7 +64,8 @@ s=3 ;power spectrum of CR E distribution A=1 ;Synchrotron radiation amplitude at 10 mm psi=0. ;default polarization angle smallp=0.3;0.3 ;default polarization fraction -scope='ADD_SED+ADD_POLSED' +;scope='ADD_SED+ADD_POLSED' +scope='ADD_SED' ;paramtag=['s (plaw_index)','Amp','p','Psi (deg)'] paramtag=[textoidl('\alpha_{CR}'),'Amp','p',textoidl('\psi')+' [deg]'] paramdefault=[s,A,smallp,psi] diff --git a/src/idl/dustem_read_dustem_lv.pro b/src/idl/dustem_read_dustem_lv.pro index e294c7a..6c5c4cd 100755 --- a/src/idl/dustem_read_dustem_lv.pro +++ b/src/idl/dustem_read_dustem_lv.pro @@ -79,6 +79,7 @@ ENDWHILE ;readf,unit,str CLOSE,unit free_lun,unit + str=strcompress(str) str=strtrim(str,2) vv=str_sep(str,' ') @@ -100,12 +101,11 @@ Nlines=n_elements(wav) instruc='one_st={wav:0.,' FOR i=0L,Ngrains-1 DO BEGIN - instruc=instruc+'em_grain_'+strtrim(i+1,2)+':0.,' + instruc=instruc+'em_grain_'+strtrim(i+1,2)+':0.D0,' ENDFOR -instruc=instruc+'em_tot:0.}' +instruc=instruc+'em_tot:0.D0}' message,'Executing '+instruc,/info toto=execute(instruc) - st=replicate(one_st,Nlines) st.wav=wav @@ -113,6 +113,7 @@ FOR i=0L,Ngrains-1 DO BEGIN instruc='st.(i+1)='+'em_grain_'+strtrim(i+1,2) message,'Executing '+instruc,/info toto=execute(instruc) + ;stop ENDFOR st.em_tot=em_tot diff --git a/src/idl/dustem_run_plugins.pro b/src/idl/dustem_run_plugins.pro index aed96fa..0677694 100644 --- a/src/idl/dustem_run_plugins.pro +++ b/src/idl/dustem_run_plugins.pro @@ -1,13 +1,13 @@ -PRO dustem_run_plugins, p_dim ,$ - param_descs,$ - param_values,$ - param_func,$ - scopes,$ - st=st,$ - avoid=avoid,$ - force_dustem_run=force_dustem_run,$ - use_previous_fortran=use_previous_fortran, $ - help=help +PRO dustem_run_plugins, p_dim $ + ,param_descs $ + ,param_values $ + ,param_func $ + ,run_order $ + ,st=st $ + ,avoid=avoid $ + ,force_dustem_run=force_dustem_run $ + ,use_previous_fortran=use_previous_fortran $ + ,help=help ;+ ; NAME: @@ -28,7 +28,7 @@ PRO dustem_run_plugins, p_dim ,$ ; param_descs = parameter description vector ; param_values = current parameter values ; param-func = plugin indices array -; scopes = scopes of plugins +; run_order = run order of the plugins to be run ; ; OPTIONAL INPUT PARAMETERS: ; avoid = scopes should be avoided @@ -85,12 +85,14 @@ Nplugins=n_elements(*!dustem_plugin) must_be_run=intarr(Nplugins) FOR i=0L,Nplugins-1 DO BEGIN - FOR j=0L,n_elements(scopes)-1 DO BEGIN - pos=strposmulti((*!dustem_plugin)[i].scope,scopes[j],count) - IF count NE 0 THEN must_be_run[i]=1 - ENDFOR + IF (*!dustem_plugin)[i].run_order EQ run_order THEN BEGIN + must_be_run[i]=1 + ENDIF ENDFOR +;print,must_be_run +;stop + ;=== run plugins that need be run FOR i=0L,Nplugins-1 DO BEGIN IF must_be_run[i] THEN BEGIN -- libgit2 0.21.2