diff --git a/src/idl/dustem_fit_J13.pro b/src/idl/dustem_fit_J13.pro deleted file mode 100755 index a697090..0000000 --- a/src/idl/dustem_fit_J13.pro +++ /dev/null @@ -1,266 +0,0 @@ -PRO dustem_fit_j13,help=help,noplot=noplot,mode=mode - -;+ -; NAME: -; dustem_fit -; PURPOSE: -; Runs conjointly a fit on extinction, emission and polarization by extinction -; CATEGORY: -; Dustem -; CALLING SEQUENCE: -; sed=dustem_fit([/noplot,/help]) -; INPUTS: -; None -; OPTIONAL INPUT PARAMETERS: -; noplot : do not display figures -; OUTPUTS: -; None -; OPTIONAL OUTPUT PARAMETERS: -; ACCEPTED KEY-WORDS: -; help = If set, print this help -; COMMON BLOCKS: -; None -; SIDE EFFECTS: -; None -; RESTRICTIONS: -; The dustem idl wrapper must be installed -; PROCEDURE: -; None -; EXAMPLES -; -; MODIFICATION HISTORY: -; TLS added by D. Paradis. -; Written by V. Guillet (2012) -; see evolution details on the dustem cvs maintained at CESR -; Contact J.-Ph. Bernard (Jean-Philippe.Bernard@cesr.fr) in case of problems. -;- - -IF keyword_set(mode) THEN BEGIN - use_mode=strupcase(mode) -ENDIF ELSE BEGIN - use_mode='COMPIEGNE_ETAL2010' ;Default is last dustem model -ENDELSE - -; Verbose -defsysv, '!dustem_verbose', 1 - -data_files = { $ - ; sed: !dustem_wrap_soft_dir + 'Data/SEDs/Gal_composite_spectrum.xcat', $ - sed: '/Users/paradis/Soft_Librairies/Planck/MW/HA_diffuse_COlevel.xcat' $ - ; sed: !dustem_wrap_soft_dir + 'Data/SEDs/spec_firas_chi2.xcat', $ - ; ext: !dustem_wrap_soft_dir + 'Data/EXTs/spec_sk-69228-2.xcat' $ - ; polext: !dustem_wrap_soft_dir + 'Data/POLEXTs/Serkowski_DISM.xcat', $ - ; polsed: '' $ - } - -; Local variable, to define global !fit_rchi2_weight -rchi2_weight = { $ - sed: 1. $ - ; ext: 1. $ - ; polext: 1. $ - ; polsed: 1. $ - } - -; Flags to determine what to fit -defsysv, '!dustem_data', { $ ;Data to fit - sed: ptr_new() $ - ; ext: ptr_new() $ - ; polext: ptr_new() $ -; polsed: ptr_new() $ - } - -; Define global !fit_rchi2_weight with the same structure and field order like !dustem_data -tagnames = tag_names(!dustem_data) -instr="defsysv, '!fit_rchi2_weight', {" -for i = 0, n_elements(tagnames)-1 do instr+=tagnames(i)+': rchi2_weight.'+tagnames(i)+', ' -instr+='empty: 0. } -toto=execute(instr) - -;=== Set weight for reduced chi2 ponderation between emission, extinction, polarized extinction and emission === -; Calculate total according to fileds in !dustem_data structure -total = 0 -for i = 0, n_tags(!dustem_data)-1 do total = total + !fit_rchi2_weight.(i) -; values are relative to total -for i = 0, n_tags(!dustem_data)-1 do !fit_rchi2_weight.(i) /= total - -dustem_init,mode=use_mode - -if keyword_set(noplot) then defsysv,'!dustem_show_plot', 0 - -; Check that polarisation has been run (-pol in GRAIN.DAT) if fit on polarization is asked for -if (tag_exist(!dustem_data,'polext') or tag_exist(!dustem_data,'polsed')) and !run_pol eq 0. then begin - message,"pol option must be set in GRAIN.DAT to fit polarization",/info - stop -endif - -; Check that DTLS has been run (-dtls in GRAIN.DAT) -if (tag_exist(!dustem_data,'dtls') ) and !run_tls eq 0. then begin - message,"dtls option must be set in GRAIN.DAT to use the TLS model",/info - stop -endif - -;=== Read observational data === -for i = 0, n_elements(tagnames)-1 do begin - instr='!dustem_data.('+strtrim(i,2)+') = dustem_set_data(read_xcat(data_files.'+tagnames(i)+',/silent))' - toto=execute(instr) -endfor - -;=== Parameters for fit === -pd_mass = [ $ - '(*!dustem_params).grains(0).mdust_o_mh', $ ;PAH0 mass fraction - '(*!dustem_params).grains(1).mdust_o_mh', $ ;PAH1 mass fraction - '(*!dustem_params).grains(2).mdust_o_mh' $ ;Small amCBEx - ; '(*!dustem_params).grains(3).mdust_o_mh' $ ;Large amCBEx - ] - -;tls=['(*!dustem_params).tls.a_dtls', $ -; '(*!dustem_params).tls.lc' $ -; ; '(*!dustem_params).tls.vm' $ -;] - -iv_tls=[10., 9.] - -iv_mass = [ $ - 2.23e-4, $ - 2.64e-6, $ - 3.5e-4 $ -; 5.7e-4 $ - ] ; DUST MASSES -;iv_mass = [ $ -; 9.3E-4, $ -; 0.000449, $ -; 2.2e-3 $ -;; 0.0008 $ -;; 7.8e-4 $ ;7.8e-3 -; ] ; DUST MASSES -;iv_mass = [ $ -; 1.1e-3, $ -; 1.8e-9, $ -; 6.e-3, $ -; 0.01 $ -; ] ; DUST MASSES - -pd_sdist = [ $ -;; '(*!dustem_params).grains(0).amin', $ ;PAH0 minmal size -;; '(*!dustem_params).grains(0).amax' $ ;PAH0 maximal size -; '(*!dustem_params).grains(0).alpha_o_a0', $ ;PAH0 mean size - '(*!dustem_params).grains(0).alpha_o_a0' $ ;PAH0 mean size - -;; '(*!dustem_params).grains(0).at', $ ;PAH0 sigma size -; '(*!dustem_params).grains(0).amin', $ ;PAH1 minmal size -; '(*!dustem_params).grains(0).amax' $ ;PAH1 maximal size -; '(*!dustem_params).grains(2).alpha_o_a0' $ ;PAH1 mean size -;; '(*!dustem_params).grains(1).at', $ ;PAH1 sigma size -;; '(*!dustem_params).grains(2).amin', $ ;VSG minimal size -;; '(*!dustem_params).grains(2).amax', $ ;VSG maximal size -; '(*!dustem_params).grains(3).amin', $ ;Carbon BG minimal size -;; '(*!dustem_params).grains(3).amax', $ ;Carbon BG maximal size -; '(*!dustem_params).grains(2).alpha_o_a0', $ ;Carbon BG power law slope -; '(*!dustem_params).grains(3).alpha_o_a0' $ ;Carbon BG power law slope - -; '(*!dustem_params).grains(3).at', $ ;Carbon BG threshold size (exp. decay) -; '(*!dustem_params).grains(3).ac', $ ;Carbon BG cut-off size (exp. decay) -; '(*!dustem_params).grains(3).gamma', $ ;Carbon BG exp. decay slope -; '(*!dustem_params).grains(4).amin', $ ;Silicate BG minimal size -;; '(*!dustem_params).grains(4).amax', $ ;Silicate BG maximal size -; '(*!dustem_params).grains(4).alpha_o_a0', $ ;Silicate BG power law slope -; '(*!dustem_params).grains(4).at', $ ;Silicate BG threshold size (exp. decay) -; '(*!dustem_params).grains(4).ac', $ ;Silicate BG cut-off size (exp. decay) -; '(*!dustem_params).grains(4).gamma' $ ;Silicate BG exp. decay slope - ] -iv_sdist = [ $ -;3.5000E-08, 1.2000E-07, 6.4000E-08, 1.0000E-01 $ ; PAH0 -;;3.5000E-08, 1.2000E-07, 6.4000E-08, 1.0000E-01, $ ; PAH1 -;-4.98, 3.12e-7 $ --4. $ -;1.7e-8,-3.21 $ ; VSG -;;6.0000E-07, 5.0000E-04,-2.6000E+00, 0.7000E-05, 0.7000E-05, 0.9000E+00, $ Carbon BG -;;6.0000E-07, 5.0000E-04,-3.4000E+00, 3.0000E-05, 1.0000E-05, 0.4200E+00 $ Silicate BG -;6.0000E-07, -2.6000E+00, 0.7000E-05, 0.7000E-05, 0.9000E+00, $ Carbon BG -;6.0000E-07, -3.4000E+00, 3.0000E-05, 1.0000E-05, 0.4200E+0 $ Silicate BG - ] - -pd_pol_alignment = [ $ - '(*!dustem_params).pol(4).atresh', $ ; Threshold size for grain alignment - '(*!dustem_params).pol(4).pstiff', $ ; Stiffness of transition between non-aligned and aligned grains sizes - '(*!dustem_params).pol(4).plev' $ ; Maximal alignment efficiency of grains - ] -iv_pol_alignment = [ $ - 0.1, $ - 0.5, $ - 1.0 $ - ] - -pd_G0=['(*!dustem_params).G0'] ; -iv_G0=[4.] - -pd_cont = ['dustem_create_continuum_2'] -iv_cont = [0.0039] -;stop -;pd = [ pd_pol_alignment, pd_mass ] -;iv = [ iv_pol_alignment, iv_mass ] - -; Parameters -;pd = [ pd_G0, pd_mass, pd_pol_alignment, pd_sdist, pd_cont ] -;iv = [ iv_G0, iv_mass, iv_pol_alignment, iv_sdist, iv_cont ] -;pd = [ pd_mass, pd_sdist, pd_cont ] -;iv = [ iv_mass, iv_sdist, iv_cont ] -;pd = [ pd_mass, tls,pd_cont ] -;iv = [ iv_mass, iv_tls,iv_cont ] -pd = [pd_G0, pd_mass, pd_cont,pd_sdist];,tls] -iv = [iv_G0, iv_mass, iv_cont,iv_sdist];,iv_tls] -;pd = [ pd_cont,tls,pd_G0 ] -;iv = [ iv_cont,iv_tls,iv_G0 ] - -Npar=n_elements(pd) -ulimed=replicate(0,Npar) -llimed=replicate(1,Npar) -llims=replicate(0.,Npar) -ulims=replicate(0.,Npar) - -; Upper limit on alignment efficiency -;ulimed(1+n_elements(pd_mass)+2) = 1 -;ulims(1+n_elements(pd_mass)+2) = 1. -; Upper limit on grain size -; amCBEx -;ulimed(1+n_elements(pd_mass)+n_elements(pd_pol_alignment)+9) = 1 -;ulims(1+n_elements(pd_mass)+n_elements(pd_pol_alignment)+9) = 9.9d-4 - -;Npar=n_elements(pd) -;ulimed=replicate(0,Npar) -;llimed=replicate(1,Npar) -;llims=replicate(0.,Npar) -;ulims=replicate(0.,Npar) -;stop -; ulimed=[0 , 0. ,0 ,0, 0, 0,0,1] ; aucun des parametres n'a une valeur superieure limite -;ulims=[0 , 0. ,0 ,0, 0, 0, 0,5.e-6] - llimed=[1, 1, 1., 1.,1,1.] ; tous les parametres ont une valeur inferieure limite - llims =[0.1,0, 0,0.,0.,-2];, 1.e-8];,0,0,3.] - - -;=== Fixed parameters -;fpd=['(*!dustem_params).G0'] ;Fixed parameter description -;fiv=[1000] -;dustem_init_fixed_params,fpd,fiv -fpd=['(*!dustem_params).grains(3).mdust_o_mh'];,'(*!dustem_params).grains(3).mdust_o_mh'] ;Fixed parameter description -fiv=[1.e-16];,1.e-16] -dustem_init_fixed_params,fpd,fiv - -;== SET THE FITTED PARAMETERS -dustem_init_parinfo,pd,iv,up_limited=ulimed,lo_limited=llimed,up_limits=ulims,lo_limits=llims - -;=== RUN fit -tol=1.e-13 -Nitermax=20 ;maximum number of iteration. This is the criterium which will stop the fit procedure - -loadct,13 -!y.range=[1.e-3,100.] ;This is to ajust plot range from outside the routine -t1=systime(0,/sec) -res=dustem_mpfit_data(tol=tol,Nitermax=Nitermax) ;,func_name=func_name,cf_min=cf_min) -t2=systime(0,/sec) -print,'executed in ',t2-t1,' sec' - -;st=dustem_read_all_res(!dustem_res,/silent) - -END - -- libgit2 0.21.2