Loading Manifest.toml +55 −57 Original line number Diff line number Diff line Loading @@ -11,9 +11,9 @@ uuid = "0dad84c5-d112-42e6-8d28-ef12dabb789f" [[ArrayInterface]] deps = ["IfElse", "LinearAlgebra", "Requires", "SparseArrays", "Static"] git-tree-sha1 = "045ff5e1bc8c6fb1ecb28694abba0a0d55b5f4f5" git-tree-sha1 = "a71d224f61475b93c9e196e83c17c6ac4dedacfa" uuid = "4fba245c-0d91-5ea0-9b3e-6abc04ee57a9" version = "3.1.17" version = "3.1.18" [[Artifacts]] uuid = "56f22d72-fd6d-98f1-02f0-08ddc0907c33" Loading Loading @@ -41,9 +41,9 @@ version = "0.3.3" [[ChainRulesCore]] deps = ["Compat", "LinearAlgebra", "SparseArrays"] git-tree-sha1 = "be770c08881f7bb928dfd86d1ba83798f76cf62a" git-tree-sha1 = "f53ca8d41e4753c41cdafa6ec5f7ce914b34be54" uuid = "d360d2e6-b24c-11e9-a2a3-2a2ae2dbcce4" version = "0.10.9" version = "0.10.13" [[CodeTracking]] deps = ["InteractiveUtils", "UUIDs"] Loading @@ -59,9 +59,9 @@ version = "0.7.0" [[ColorSchemes]] deps = ["ColorTypes", "Colors", "FixedPointNumbers", "Random", "StaticArrays"] git-tree-sha1 = "c8fd01e4b736013bc61b704871d20503b33ea402" git-tree-sha1 = "ed268efe58512df8c7e224d2e170afd76dd6a417" uuid = "35d6a980-a343-548e-a6ea-1d62b119f2f4" version = "3.12.1" version = "3.13.0" [[ColorTypes]] deps = ["FixedPointNumbers", "Random"] Loading Loading @@ -119,6 +119,12 @@ git-tree-sha1 = "ee400abb2298bd13bfc3df1c412ed228061a2385" uuid = "9a962f9c-6df0-11e9-0e5d-c546b8b5ee8a" version = "1.7.0" [[DataFrames]] deps = ["Compat", "DataAPI", "Future", "InvertedIndices", "IteratorInterfaceExtensions", "LinearAlgebra", "Markdown", "Missings", "PooledArrays", "PrettyTables", "Printf", "REPL", "Reexport", "SortingAlgorithms", "Statistics", "TableTraits", "Tables", "Unicode"] git-tree-sha1 = "a19645616f37a2c2c3077a44bc0d3e73e13441d7" uuid = "a93c6f00-e57d-5684-b7b6-d8193f3e46c0" version = "1.2.1" [[DataStructures]] deps = ["Compat", "InteractiveUtils", "OrderedCollections"] git-tree-sha1 = "4437b64df1e0adccc3e5d1adbc3ac741095e4677" Loading Loading @@ -174,9 +180,9 @@ uuid = "8ba89e20-285c-5b6f-9357-94700520ee1b" [[Distributions]] deps = ["FillArrays", "LinearAlgebra", "PDMats", "Printf", "QuadGK", "Random", "SparseArrays", "SpecialFunctions", "Statistics", "StatsBase", "StatsFuns"] git-tree-sha1 = "2733323e5c02a9d7f48e7a3c4bc98d764fb704da" git-tree-sha1 = "3889f646423ce91dd1055a76317e9a1d3a23fff1" uuid = "31c24e10-a181-5473-b8eb-7969acd0382f" version = "0.25.6" version = "0.25.11" [[DocStringExtensions]] deps = ["LibGit2"] Loading @@ -185,25 +191,19 @@ uuid = "ffbed154-4ef7-542d-bbb7-c09d3a79fcae" version = "0.8.5" [[DoubleFloats]] deps = ["GenericLinearAlgebra", "GenericSchur", "LinearAlgebra", "Polynomials", "Printf", "Quadmath", "Random", "Requires", "SpecialFunctions"] git-tree-sha1 = "cfc5657a37c2881a728d76bd14ad808ca096d601" deps = ["GenericLinearAlgebra", "LinearAlgebra", "Polynomials", "Printf", "Quadmath", "Random", "Requires", "SpecialFunctions"] git-tree-sha1 = "1c962cf7e75c09a5f1fbf504df7d6a06447a1129" uuid = "497a8b3b-efae-58df-a0af-a86822472b78" version = "1.1.22" version = "1.1.23" [[Downloads]] deps = ["ArgTools", "LibCURL", "NetworkOptions"] uuid = "f43a241f-c20a-4ad4-852c-f6b1247861c6" [[ExprTools]] git-tree-sha1 = "10407a39b87f29d47ebaca8edbc75d7c302ff93e" git-tree-sha1 = "b7e3d17636b348f005f11040025ae8c6f645fe92" uuid = "e2ba6199-217a-4e67-a87a-7c52f15ade04" version = "0.1.3" [[EzXML]] deps = ["Printf", "XML2_jll"] git-tree-sha1 = "0fa3b52a04a4e210aeb1626def9c90df3ae65268" uuid = "8f5d6c58-4d21-5cfd-889c-e3ad7ee6a615" version = "1.1.0" version = "0.1.6" [[FFTW]] deps = ["AbstractFFTs", "FFTW_jll", "LinearAlgebra", "MKL_jll", "Preferences", "Reexport"] Loading Loading @@ -234,15 +234,15 @@ uuid = "7b1f6079-737a-58dc-b8bc-7a2ca5c1b5ee" [[FillArrays]] deps = ["LinearAlgebra", "Random", "SparseArrays"] git-tree-sha1 = "31939159aeb8ffad1d4d8ee44d07f8558273120a" git-tree-sha1 = "25b9cc23ba3303de0ad2eac03f840de9104c9253" uuid = "1a297f60-69ca-5386-bcde-b61e274b549b" version = "0.11.7" version = "0.12.0" [[FiniteDiff]] deps = ["ArrayInterface", "LinearAlgebra", "Requires", "SparseArrays", "StaticArrays"] git-tree-sha1 = "f6f80c8f934efd49a286bb5315360be66956dfc4" git-tree-sha1 = "8b3c09b56acaf3c0e581c66638b85c8650ee9dca" uuid = "6a86dc24-6348-571c-b903-95158fe2bd41" version = "2.8.0" version = "2.8.1" [[FixedPointNumbers]] deps = ["Statistics"] Loading @@ -262,6 +262,10 @@ git-tree-sha1 = "e2af66012e08966366a43251e1fd421522908be6" uuid = "f6369f11-7733-5829-9624-2563aa707210" version = "0.10.18" [[Future]] deps = ["Random"] uuid = "9fa8497b-333b-5362-9e8d-4d0656e87820" [[GFit]] deps = ["CMPFit", "DataStructures", "Dates", "Dierckx", "Distributions", "LsqFit", "MacroTools", "PrettyTables", "Printf", "ProgressMeter", "QuadGK", "Random", "Statistics"] path = "/home/gcalderone/.julia/dev/GFit" Loading @@ -280,12 +284,6 @@ git-tree-sha1 = "ff291c1827030ffaacaf53e3c83ed92d4d5e6fb6" uuid = "14197337-ba66-59df-a3e3-ca00e7dcff7a" version = "0.2.5" [[GenericSchur]] deps = ["LinearAlgebra", "Printf"] git-tree-sha1 = "372e48d7f3ced17fdc888a841bcce77be417ce57" uuid = "c145ed77-6b09-5dd9-b285-bf645a82121e" version = "0.5.0" [[Gnuplot]] deps = ["ColorSchemes", "ColorTypes", "Colors", "DataStructures", "REPL", "ReplMaker", "StatsBase", "StructC14N", "Test"] path = "/home/gcalderone/.julia/dev/Gnuplot" Loading Loading @@ -313,6 +311,12 @@ git-tree-sha1 = "323a38ed1952d30586d0fe03412cde9399d3618b" uuid = "d8418881-c3e1-53bb-8760-2df7ec849ed5" version = "1.5.0" [[InvertedIndices]] deps = ["Test"] git-tree-sha1 = "15732c475062348b0165684ffe28e85ea8396afc" uuid = "41ab1584-1d38-5bbf-9106-f11c6c58b48f" version = "1.0.0" [[IterTools]] git-tree-sha1 = "05110a2ab1fc5f932622ffea2a003221f4782c18" uuid = "c8e1da08-722c-5040-9ed9-7db0dc04731e" Loading Loading @@ -370,21 +374,15 @@ uuid = "29816b5a-b9ab-546f-933c-edad1886dfa8" [[Libdl]] uuid = "8f399da3-3557-5675-b5ff-fb832c97cbdb" [[Libiconv_jll]] deps = ["Artifacts", "JLLWrappers", "Libdl", "Pkg"] git-tree-sha1 = "42b62845d70a619f063a7da093d995ec8e15e778" uuid = "94ce4f54-9a6c-5748-9c1c-f9c7231a4531" version = "1.16.1+1" [[LinearAlgebra]] deps = ["Libdl"] uuid = "37e2e46d-f89d-539d-b4ee-838fcccc9c8e" [[LogExpFunctions]] deps = ["DocStringExtensions", "LinearAlgebra"] git-tree-sha1 = "1ba664552f1ef15325e68dc4c05c3ef8c2d5d885" git-tree-sha1 = "7bd5f6565d80b6bf753738d2bc40a5dfea072070" uuid = "2ab3a3ac-af41-5b50-aa03-7779005ae688" version = "0.2.4" version = "0.2.5" [[Logging]] uuid = "56ddb016-857b-54e1-b83d-db4d58db5568" Loading Loading @@ -499,9 +497,15 @@ uuid = "44cfe95a-1eb2-52ea-b672-e2afdf69b78f" [[Polynomials]] deps = ["Intervals", "LinearAlgebra", "MutableArithmetics", "RecipesBase"] git-tree-sha1 = "3685cb1e2ccbbb9973684774f956f5c6e4673bc9" git-tree-sha1 = "0bbfdcd8cda81b8144de4be8a67f5717e959a005" uuid = "f27b6e38-b328-58d1-80ce-0feddd5e7a45" version = "2.0.12" version = "2.0.14" [[PooledArrays]] deps = ["DataAPI", "Future"] git-tree-sha1 = "cde4ce9d6f33219465b55162811d8de8139c0414" uuid = "2dfb63ee-cc39-5dd5-95bd-886bf059d720" version = "1.2.1" [[Preferences]] deps = ["TOML"] Loading @@ -526,7 +530,7 @@ uuid = "92933f4c-e287-5a05-a399-4b506db050ca" version = "1.7.1" [[QSFit]] deps = ["CMPFit", "Cosmology", "DSP", "DataStructures", "Dates", "DelimitedFiles", "Dierckx", "FITSIO", "GFit", "GFitViewer", "Pkg", "Printf", "Statistics", "Unitful", "UnitfulAstro"] deps = ["CMPFit", "Cosmology", "DSP", "DataFrames", "DataStructures", "Dates", "DelimitedFiles", "Dierckx", "FITSIO", "GFit", "GFitViewer", "Gnuplot", "Pkg", "Printf", "REPL", "Statistics", "TextParse", "Unitful", "UnitfulAstro"] path = "/home/gcalderone/.julia/dev/QSFit" uuid = "fea25315-b951-4667-83c9-50e53efab241" version = "0.1.0" Loading Loading @@ -606,9 +610,9 @@ uuid = "6462fe0b-24de-5631-8697-dd941f90decc" [[SortingAlgorithms]] deps = ["DataStructures"] git-tree-sha1 = "2ec1962eba973f383239da22e75218565c390a96" git-tree-sha1 = "b3363d7460f7d098ca0912c69b082f75625d7508" uuid = "a2af1166-a08f-5f64-846c-94a0d3cef48c" version = "1.0.0" version = "1.0.1" [[SparseArrays]] deps = ["LinearAlgebra", "Random"] Loading @@ -622,15 +626,15 @@ version = "1.5.1" [[Static]] deps = ["IfElse"] git-tree-sha1 = "2740ea27b66a41f9d213561a04573da5d3823d4b" git-tree-sha1 = "62701892d172a2fa41a1f829f66d2b0db94a9a63" uuid = "aedffcd0-7271-4cad-89d0-dc628f76c6d3" version = "0.2.5" version = "0.3.0" [[StaticArrays]] deps = ["LinearAlgebra", "Random", "Statistics"] git-tree-sha1 = "745914ebcd610da69f3cb6bf76cb7bb83dcb8c9a" git-tree-sha1 = "1b9a0f17ee0adde9e538227de093467348992397" uuid = "90137ffa-7385-5640-81b9-e52037218182" version = "1.2.4" version = "1.2.7" [[Statistics]] deps = ["LinearAlgebra", "SparseArrays"] Loading Loading @@ -694,10 +698,10 @@ uuid = "e0df1984-e451-5cb5-8b61-797a481e67e3" version = "1.0.1" [[TimeZones]] deps = ["Dates", "EzXML", "LazyArtifacts", "Mocking", "Pkg", "Printf", "RecipesBase", "Serialization", "Unicode"] git-tree-sha1 = "960099aed321e05ac649c90d583d59c9309faee1" deps = ["Dates", "Future", "LazyArtifacts", "Mocking", "Pkg", "Printf", "RecipesBase", "Serialization", "Unicode"] git-tree-sha1 = "81753f400872e5074768c9a77d4c44e70d409ef0" uuid = "f269a46b-ccf7-5d73-abea-4c690281aa53" version = "1.5.5" version = "1.5.6" [[TranscodingStreams]] deps = ["Random", "Test"] Loading @@ -714,9 +718,9 @@ uuid = "4ec0a83e-493e-50e2-b9ac-8f72acf5a8f5" [[Unitful]] deps = ["ConstructionBase", "Dates", "LinearAlgebra", "Random"] git-tree-sha1 = "b3682a0559219355f1e3c8024e9f97adce2d4623" git-tree-sha1 = "a981a8ef8714cba2fd9780b22fd7a469e7aaf56d" uuid = "1986cc42-f94f-5a68-af5c-568840ba703d" version = "1.8.0" version = "1.9.0" [[UnitfulAngles]] deps = ["Dates", "Unitful"] Loading @@ -736,12 +740,6 @@ git-tree-sha1 = "28807f85197eaad3cbd2330386fac1dcb9e7e11d" uuid = "ea10d353-3f73-51f8-a26c-33c1cb351aa5" version = "0.6.2" [[XML2_jll]] deps = ["Artifacts", "JLLWrappers", "Libdl", "Libiconv_jll", "Pkg", "Zlib_jll"] git-tree-sha1 = "1acf5bdf07aa0907e0a37d3718bb88d4b687b74a" uuid = "02c8fc9c-b97f-50b9-bbe4-9be30ff0a78a" version = "2.9.12+0" [[Zlib_jll]] deps = ["Libdl"] uuid = "83775a58-1f1d-513f-b197-d71354ab007a" Loading run.jl +17 −16 Original line number Diff line number Diff line Loading @@ -11,7 +11,7 @@ using Dates using TextParse import StatsBase, JLD2 path = "output-" * readlines(`git branch --show-current`)[1] path = "output-new-" * readlines(`git branch --show-current`)[1] mkdir("$(path)") epoch_filenames = Vector{String}() for (root, dirs, files) in walkdir("AT2018zf") Loading @@ -32,8 +32,8 @@ function read_spec(epoch_id; kw...) spec = Spectrum(Val(:ASCII), file, columns=[1,2]; label="$epoch_id: $file", kw...); if isa(epoch_id, Number) spec.flux ./= all_scale[epoch_id] spec.err ./= all_scale[epoch_id] end spec.err .= 0.05 .* median(spec.flux); return spec end Loading @@ -48,13 +48,13 @@ I tried 180, 200 and 250, and the [OIII] luminosity do not change =# spec = read_spec("B03", resolution=200.) add_spec!(source, spec); res = fit(source); res = qsfit(source); viewer(res, showcomps=[:qso_cont, :galaxy, :balmer], filename="$(path)/results_boller.html") println(res.bestfit[:OIII_5007].norm.val ) # 0.0002913168348334021 println(res.bestfit[:OIII_5007].norm.val / res.bestfit[:galaxy].norm.val ) # 64.65011509608263 println(res.model[:OIII_5007].norm.val ) # 0.0002913095378291511 println(res.model[:OIII_5007].norm.val / res.model[:galaxy].norm.val ) # 61.34678931179987 Loading @@ -71,14 +71,14 @@ if !isfile("$(path)/scale_fwhm.dat") # source.options[:host_template] = "/home/gcalderone/Mbi1.30Zm1.49T00.0700_iTp0.00_baseFe_linear_FWHM_2.51" spec = read_spec(id) add_spec!(source, spec); res = fit(source); res = qsfit(source); viewer(res, showcomps=[:qso_cont, :galaxy, :balmer], filename="$(path)/results_$(id).html") # Update all_scale and store FWHM @info res.bestfit[:OIII_5007].norm.val all_scale[id] *= res.bestfit[:OIII_5007].norm.val all_fwhm[ id] = res.bestfit[:OIII_5007].fwhm.val @info res.model[:OIII_5007].norm.val all_scale[id] *= res.model[:OIII_5007].norm.val all_fwhm[ id] = res.model[:OIII_5007].fwhm.val catch end end Loading Loading @@ -150,7 +150,8 @@ for loop in 1:Nloop add_spec!(source, spec); @gp :- :zoom xr=[4750,5150] spec.λ ./ (1 + source.z) spec.flux "w l t '$(id)'" end res = multi_fit(source); res = qsfit_multi(source, ref_id=2); JLD2.jldsave(file; res) viewer(res, showcomps=[:qso_cont, :galaxy, :balmer], filename="$(path)/results_$(job)_$(loop).html") viewer(res, showcomps=[:qso_cont, :galaxy, :balmer], filename="$(path)/results_$(job)_$(loop)_rebin4.html", rebin=4) Loading @@ -166,12 +167,12 @@ for loop in 1:Nloop end end end bestfit2 = fit!(multi2, source.data, minimizer=GFit.cmpfit()) bestfit2 = fit!(res.source, multi2, res.pspecs) OIII_norm = fill(0., length(chosen_epochs)) for id in 1:length(chosen_epochs) OIII_norm[id] = bestfit2[id][:OIII_5007].norm.val all_scale[chosen_epochs[id]] *= bestfit2[id][:OIII_5007].norm.val OIII_norm[id] = multi2[id][:OIII_5007].norm.val all_scale[chosen_epochs[id]] *= multi2[id][:OIII_5007].norm.val end JLD2.jldsave(file; res, all_scale, OIII_norm) else Loading @@ -188,7 +189,7 @@ for loop in 1:Nloop file = "$(path)/results_$(job)_$(loop).dat" res, all_scale, OIII_norm = values(JLD2.load(file)) OIII_norm_evol = hcat(OIII_norm_evol, OIII_norm) println(res.bestfit.mdc.fitstat / res.bestfit.mdc.dof) println(res.fitres.fitstat / res.fitres.dof) end OIII_norm_evol = OIII_norm_evol[:,2:end] @gp "set grid" :- Loading src/CL_1ES_1927p654.jl +224 −52 Original line number Diff line number Diff line module CL_1ES_1927p654 using QSFit, DataStructures, GFit, Statistics, Dierckx import QSFit: AbstractSpectralLine, DefaultRecipe, CombinedLine, NarrowLine, BroadLine, BroadBaseLine, ComboBroadLine, default_options, qso_cont_component, line_component, line_group_name, known_spectral_lines, fit, multi_fit, spectral_coverage, reduce import QSFit: DefaultRecipe, default_options, known_spectral_lines, LineComponent, add_qso_continuum!, add_patch_functs!, guess_emission_lines!, qsfit_multi abstract type q1927p654 <: DefaultRecipe end Loading @@ -22,63 +24,233 @@ function default_options(::Type{T}) where T <: q1927p654 end function qso_cont_component(::Type{T}) where T <: q1927p654 comp = qso_cont_component(supertype(T)) function QSFit.add_qso_continuum!(source::QSO{T}, pspec::PreparedSpectrum, model::Model) where T <: q1927p654 λ = domain(model)[:] comp = QSFit.powerlaw(3000) comp.x0.val = median(λ) comp.norm.val = Spline1D(λ, pspec.data.val, k=1, bc="error")(comp.x0.val) comp.alpha.val = -1.5 comp.alpha.low = -5 comp.alpha.high = 5 return comp end model[:qso_cont] = comp push!(model[:Continuum].list, :qso_cont) evaluate!(model) end struct AsymmNarrowLine <: AbstractSpectralLine λ::Float64 end line_group_name(::Type{T}, name::Symbol, line::AsymmNarrowLine) where T <: q1927p654 = :NarrowLines function known_spectral_lines(::Type{T}) where T <: q1927p654 list = OrderedDict{Symbol, AbstractSpectralLine}() #list[:Lyb ] = CombinedLine( 1026.0 ) list[:Lya ] = CombinedLine( 1215.24 ) list[:NV_1241 ] = NarrowLine( 1240.81 ) list[:OI_1306 ] = BroadLine( 1305.53 ) list[:CII_1335 ] = BroadLine( 1335.31 ) list[:SiIV_1400 ] = BroadLine( 1399.8 ) list[:CIV_1549 ] = CombinedLine( 1549.48 ) #list[:HeII ] = BroadLine( 1640.4 ) #list[:OIII ] = BroadLine( 1665.85 ) #list[:AlIII ] = BroadLine( 1857.4 ) list[:CIII_1909 ] = BroadLine( 1908.734) list[:CII ] = BroadLine( 2326.0 ) list[:F2420 ] = BroadLine( 2420.0 ) list[:MgII_2798 ] = CombinedLine( 2799.117) #list[:NeVN ] = NarrowLine( 3346.79 ) list[:NeVI_3426 ] = NarrowLine( 3426.85 ) list[:OII_3727 ] = NarrowLine( 3729.875) list[:NeIII_3869 ] = NarrowLine( 3869.81 ) #list[:Hd ] = BroadLine( 4102.89 ) list[:Hg ] = CombinedLine( 4341.68 ) list[:bb_Hg ] = BroadBaseLine( 4341.68 ) list[:OIII_4363 ] = NarrowLine( 4363.00 ) # TODO: Check wavelength is correct list[:HeII ] = BroadLine( 4686. ) list[:Hb ] = CombinedLine( 4862.68 ) list[:bb_Hb ] = BroadBaseLine( 4862.68 ) list[:OIII_4959 ] = NarrowLine( 4960.295) list[:OIII_5007 ] = NarrowLine( 5008.240) # AsymmNarrowLine #list[:OIII_5007_bw] = NarrowLine( 5008.240) list[:HeI_5876 ] = BroadLine( 5877.30 ) list[:OI_6300 ] = NarrowLine( 6300.00 ) # TODO: Check wavelength is correct list[:OI_6364 ] = NarrowLine( 6364.00 ) # TODO: Check wavelength is correct list[:NII_6549 ] = NarrowLine( 6549.86 ) list[:Ha ] = CombinedLine( 6564.61 ) list[:bb_Ha ] = BroadBaseLine( 6564.61 ) list[:NII_6583 ] = NarrowLine( 6585.27 ) list[:SII_6716 ] = NarrowLine( 6718.29 ) list[:SII_6731 ] = NarrowLine( 6732.67 ) function known_spectral_lines(source::QSO{T}) where T <: q1927p654 list = [ # MultiCompLine( :Lyb , [BroadLine, NarrowLine]), MultiCompLine( :Lya , [BroadLine, NarrowLine]), NarrowLine( :NV_1241 ), BroadLine( custom_transition(tid=:OI_1306 , 1305.53 )), BroadLine( custom_transition(tid=:CII_1335 , 1335.31 )), BroadLine( :SiIV_1400 ), MultiCompLine( :CIV_1549 , [BroadLine, NarrowLine]), # BroadLine( :HeII_1640 ), # BroadLine( custom_transition(tid=:OIII , 1665.85 )), # BroadLine( custom_transition(tid=:AlIII , 1857.4 )), BroadLine( :CIII_1909 ), BroadLine( custom_transition(tid=:CII , 2326.0 )), BroadLine( custom_transition(tid=:F2420 , 2420.0 )), MultiCompLine( :MgII_2798 , [BroadLine, NarrowLine]), # NarrowLine( custom_transition(tid=:NeVN , 3346.79 )), NarrowLine( custom_transition(tid=:NeVI_3426 , 3426.85 )), NarrowLine( custom_transition(tid=:OII_3727 , 3729.875)), NarrowLine( custom_transition(tid=:NeIII_3869 , 3869.81 )), # BroadLine( :Hd ), MultiCompLine( :Hg , [BroadLine, NarrowLine, BroadBaseLine]), NarrowLine( :OIII_4363 ), BroadLine( :HeII_4686 ), MultiCompLine( :Hb , [BroadLine, NarrowLine, BroadBaseLine]), NarrowLine( :OIII_4959 ), NarrowLine( :OIII_5007 ), # AsymmTailLine( :OIII_5007 , :blue), BroadLine( :HeI_5876 ), NarrowLine( :OI_6300 ), NarrowLine( :OI_6364 ), NarrowLine( :NII_6549 ), MultiCompLine( :Ha , [BroadLine, NarrowLine, BroadBaseLine]), NarrowLine( :NII_6583 ), NarrowLine( :SII_6716 ), NarrowLine( :SII_6731 )] return list end function QSFit.add_patch_functs!(source::QSO{T}, pspec::PreparedSpectrum, model::Model) where T <: q1927p654 QSFit.add_patch_functs!(parent_recipe(source), pspec, model) # Force Hg and Hb lines to have the same shape as Ha model[:Ha_br].norm.low = 1.e-10 # avoid division by zero @patch! begin model[:Hg_na].norm = model[:Hg_br].norm * (model[:Ha_na].norm / model[:Ha_br].norm) model[:Hg_na].voff = model[:Ha_na].voff model[:Hg_na].fwhm = model[:Ha_na].fwhm model[:Hg_bb].norm = model[:Hg_br].norm * (model[:Ha_bb].norm / model[:Ha_br].norm) model[:Hg_bb].voff = model[:Ha_bb].voff model[:Hg_bb].fwhm = model[:Ha_bb].fwhm model[:Hg_br].voff = model[:Ha_br].voff model[:Hg_br].fwhm = model[:Ha_br].fwhm end @patch! begin model[:Hb_na].norm = model[:Hb_br].norm * (model[:Ha_na].norm / model[:Ha_br].norm) model[:Hb_na].voff = model[:Ha_na].voff model[:Hb_na].fwhm = model[:Ha_na].fwhm model[:Hb_bb].norm = model[:Hb_br].norm * (model[:Ha_bb].norm / model[:Ha_br].norm) model[:Hb_bb].voff = model[:Ha_bb].voff model[:Hb_bb].fwhm = model[:Ha_bb].fwhm model[:Hb_br].voff = model[:Ha_br].voff model[:Hb_br].fwhm = model[:Ha_br].fwhm end end function qsfit_multi(source::QSO{TRecipe}; ref_id=1) where TRecipe <: q1927p654 elapsed = time() @assert length(source.specs) > 1 @assert 1 <= ref_id <= length(source.specs) pspecs = [PreparedSpectrum(source, id=id) for id in 1:length(source.specs)] multi = MultiModel() println(logio(source), "\nFit continuum components...") for id in 1:length(pspecs) pspec = pspecs[id] model = Model(pspec.domain) model[:Continuum] = SumReducer([]) model[:main] = SumReducer([]) push!(model[:main].list, :Continuum) select_reducer!(model, :main) delete!(model.revals, :default_sum) # TODO if source.options[:instr_broadening] # TODO GFit.set_instr_response!(model[1], (l, f) -> instrumental_broadening(l, f, source.spectra[id].resolution)) # TODO end QSFit.add_qso_continuum!(source, pspec, model) QSFit.add_host_galaxy!(source, pspec, model) QSFit.add_balmer_cont!(source, pspec, model) push!(multi, model) if id != ref_id @patch! multi[id][:galaxy].norm = multi[ref_id][:galaxy].norm end end fitres = fit!(source, multi, pspecs) for id in 1:length(pspecs) model = multi[id] pspec = pspecs[id] QSFit.renorm_cont!(source, pspec, model) freeze(model, :qso_cont) haskey(model, :galaxy) && freeze(model, :galaxy) haskey(model, :balmer) && freeze(model, :balmer) evaluate!(model) end evaluate!(multi) println(logio(source), "\nFit iron templates...") for id in 1:length(pspecs) model = multi[id] pspec = pspecs[id] model[:Iron] = SumReducer([]) push!(model[:main].list, :Iron) QSFit.add_iron_uv!( source, pspec, model) QSFit.add_iron_opt!(source, pspec, model) if length(model[:Iron].list) > 0 fitres = fit!(source, model, pspec) haskey(model, :ironuv ) && freeze(model, :ironuv) haskey(model, :ironoptbr) && freeze(model, :ironoptbr) haskey(model, :ironoptna) && freeze(model, :ironoptna) end evaluate!(model) end evaluate!(multi) println(logio(source), "\nFit known emission lines...") for id in 1:length(pspecs) model = multi[id] pspec = pspecs[id] QSFit.add_emission_lines!(source, pspec, model) guess_emission_lines!(source, pspec, model) QSFit.add_patch_functs!(source, pspec, model) multi[id][:OIII_5007].norm.val = 1. multi[id][:OIII_5007].norm.fixed = true end fitres = fit!(source, multi, pspecs) for id in 1:length(pspecs) model = multi[id] pspec = pspecs[id] for lname in keys(pspec.lcs) freeze(model, lname) end end evaluate!(multi) println(logio(source), "\nFit unknown emission lines...") for id in 1:length(pspecs) model = multi[id] pspec = pspecs[id] QSFit.add_unknown_lines!(source, pspec, model) end evaluate!(multi) println(logio(source), "\nLast run with all parameters free...") for id in 1:length(pspecs) model = multi[id] pspec = pspecs[id] thaw(model, :qso_cont) haskey(model, :galaxy ) && thaw(model, :galaxy) haskey(model, :balmer ) && thaw(model, :balmer) haskey(model, :ironuv ) && thaw(model, :ironuv) haskey(model, :ironoptbr) && thaw(model, :ironoptbr) haskey(model, :ironoptna) && thaw(model, :ironoptna) for lname in keys(pspec.lcs) thaw(model, lname) end for j in 1:source.options[:n_unk] cname = Symbol(:unk, j) if model[cname].norm.val > 0 thaw(model, cname) else freeze(model, cname) end end end fitres = fit!(source, multi, pspecs) rerun = false for id in 1:length(pspecs) model = multi[id] pspec = pspecs[id] rerun = rerun || QSFit.neglect_weak_features!(source, pspec, model, fitres) end if rerun println(logio(source), "\nRe-run fit...") fitres = fit!(source, multi, pspecs) end println(logio(source)) show(logio(source), fitres) out = QSFit.QSFitMultiResults(source, pspecs, multi, fitres) elapsed = time() - elapsed println(logio(source), "\nElapsed time: $elapsed s") close_logio(source) return out end #= function line_component(::Type{T}, line::AsymmNarrowLine) where T <: q1927p654 comp = QSFit.SpecLineAsymmGauss(line.λ) Loading @@ -101,7 +273,7 @@ function line_component(::Type{T}, line::ComboBroadLine) where T <: q1927p654 end =# include("Recipe_single.jl") include("Recipe_multi.jl") # include("Recipe_single.jl") # include("Recipe_multi.jl") end Loading
Manifest.toml +55 −57 Original line number Diff line number Diff line Loading @@ -11,9 +11,9 @@ uuid = "0dad84c5-d112-42e6-8d28-ef12dabb789f" [[ArrayInterface]] deps = ["IfElse", "LinearAlgebra", "Requires", "SparseArrays", "Static"] git-tree-sha1 = "045ff5e1bc8c6fb1ecb28694abba0a0d55b5f4f5" git-tree-sha1 = "a71d224f61475b93c9e196e83c17c6ac4dedacfa" uuid = "4fba245c-0d91-5ea0-9b3e-6abc04ee57a9" version = "3.1.17" version = "3.1.18" [[Artifacts]] uuid = "56f22d72-fd6d-98f1-02f0-08ddc0907c33" Loading Loading @@ -41,9 +41,9 @@ version = "0.3.3" [[ChainRulesCore]] deps = ["Compat", "LinearAlgebra", "SparseArrays"] git-tree-sha1 = "be770c08881f7bb928dfd86d1ba83798f76cf62a" git-tree-sha1 = "f53ca8d41e4753c41cdafa6ec5f7ce914b34be54" uuid = "d360d2e6-b24c-11e9-a2a3-2a2ae2dbcce4" version = "0.10.9" version = "0.10.13" [[CodeTracking]] deps = ["InteractiveUtils", "UUIDs"] Loading @@ -59,9 +59,9 @@ version = "0.7.0" [[ColorSchemes]] deps = ["ColorTypes", "Colors", "FixedPointNumbers", "Random", "StaticArrays"] git-tree-sha1 = "c8fd01e4b736013bc61b704871d20503b33ea402" git-tree-sha1 = "ed268efe58512df8c7e224d2e170afd76dd6a417" uuid = "35d6a980-a343-548e-a6ea-1d62b119f2f4" version = "3.12.1" version = "3.13.0" [[ColorTypes]] deps = ["FixedPointNumbers", "Random"] Loading Loading @@ -119,6 +119,12 @@ git-tree-sha1 = "ee400abb2298bd13bfc3df1c412ed228061a2385" uuid = "9a962f9c-6df0-11e9-0e5d-c546b8b5ee8a" version = "1.7.0" [[DataFrames]] deps = ["Compat", "DataAPI", "Future", "InvertedIndices", "IteratorInterfaceExtensions", "LinearAlgebra", "Markdown", "Missings", "PooledArrays", "PrettyTables", "Printf", "REPL", "Reexport", "SortingAlgorithms", "Statistics", "TableTraits", "Tables", "Unicode"] git-tree-sha1 = "a19645616f37a2c2c3077a44bc0d3e73e13441d7" uuid = "a93c6f00-e57d-5684-b7b6-d8193f3e46c0" version = "1.2.1" [[DataStructures]] deps = ["Compat", "InteractiveUtils", "OrderedCollections"] git-tree-sha1 = "4437b64df1e0adccc3e5d1adbc3ac741095e4677" Loading Loading @@ -174,9 +180,9 @@ uuid = "8ba89e20-285c-5b6f-9357-94700520ee1b" [[Distributions]] deps = ["FillArrays", "LinearAlgebra", "PDMats", "Printf", "QuadGK", "Random", "SparseArrays", "SpecialFunctions", "Statistics", "StatsBase", "StatsFuns"] git-tree-sha1 = "2733323e5c02a9d7f48e7a3c4bc98d764fb704da" git-tree-sha1 = "3889f646423ce91dd1055a76317e9a1d3a23fff1" uuid = "31c24e10-a181-5473-b8eb-7969acd0382f" version = "0.25.6" version = "0.25.11" [[DocStringExtensions]] deps = ["LibGit2"] Loading @@ -185,25 +191,19 @@ uuid = "ffbed154-4ef7-542d-bbb7-c09d3a79fcae" version = "0.8.5" [[DoubleFloats]] deps = ["GenericLinearAlgebra", "GenericSchur", "LinearAlgebra", "Polynomials", "Printf", "Quadmath", "Random", "Requires", "SpecialFunctions"] git-tree-sha1 = "cfc5657a37c2881a728d76bd14ad808ca096d601" deps = ["GenericLinearAlgebra", "LinearAlgebra", "Polynomials", "Printf", "Quadmath", "Random", "Requires", "SpecialFunctions"] git-tree-sha1 = "1c962cf7e75c09a5f1fbf504df7d6a06447a1129" uuid = "497a8b3b-efae-58df-a0af-a86822472b78" version = "1.1.22" version = "1.1.23" [[Downloads]] deps = ["ArgTools", "LibCURL", "NetworkOptions"] uuid = "f43a241f-c20a-4ad4-852c-f6b1247861c6" [[ExprTools]] git-tree-sha1 = "10407a39b87f29d47ebaca8edbc75d7c302ff93e" git-tree-sha1 = "b7e3d17636b348f005f11040025ae8c6f645fe92" uuid = "e2ba6199-217a-4e67-a87a-7c52f15ade04" version = "0.1.3" [[EzXML]] deps = ["Printf", "XML2_jll"] git-tree-sha1 = "0fa3b52a04a4e210aeb1626def9c90df3ae65268" uuid = "8f5d6c58-4d21-5cfd-889c-e3ad7ee6a615" version = "1.1.0" version = "0.1.6" [[FFTW]] deps = ["AbstractFFTs", "FFTW_jll", "LinearAlgebra", "MKL_jll", "Preferences", "Reexport"] Loading Loading @@ -234,15 +234,15 @@ uuid = "7b1f6079-737a-58dc-b8bc-7a2ca5c1b5ee" [[FillArrays]] deps = ["LinearAlgebra", "Random", "SparseArrays"] git-tree-sha1 = "31939159aeb8ffad1d4d8ee44d07f8558273120a" git-tree-sha1 = "25b9cc23ba3303de0ad2eac03f840de9104c9253" uuid = "1a297f60-69ca-5386-bcde-b61e274b549b" version = "0.11.7" version = "0.12.0" [[FiniteDiff]] deps = ["ArrayInterface", "LinearAlgebra", "Requires", "SparseArrays", "StaticArrays"] git-tree-sha1 = "f6f80c8f934efd49a286bb5315360be66956dfc4" git-tree-sha1 = "8b3c09b56acaf3c0e581c66638b85c8650ee9dca" uuid = "6a86dc24-6348-571c-b903-95158fe2bd41" version = "2.8.0" version = "2.8.1" [[FixedPointNumbers]] deps = ["Statistics"] Loading @@ -262,6 +262,10 @@ git-tree-sha1 = "e2af66012e08966366a43251e1fd421522908be6" uuid = "f6369f11-7733-5829-9624-2563aa707210" version = "0.10.18" [[Future]] deps = ["Random"] uuid = "9fa8497b-333b-5362-9e8d-4d0656e87820" [[GFit]] deps = ["CMPFit", "DataStructures", "Dates", "Dierckx", "Distributions", "LsqFit", "MacroTools", "PrettyTables", "Printf", "ProgressMeter", "QuadGK", "Random", "Statistics"] path = "/home/gcalderone/.julia/dev/GFit" Loading @@ -280,12 +284,6 @@ git-tree-sha1 = "ff291c1827030ffaacaf53e3c83ed92d4d5e6fb6" uuid = "14197337-ba66-59df-a3e3-ca00e7dcff7a" version = "0.2.5" [[GenericSchur]] deps = ["LinearAlgebra", "Printf"] git-tree-sha1 = "372e48d7f3ced17fdc888a841bcce77be417ce57" uuid = "c145ed77-6b09-5dd9-b285-bf645a82121e" version = "0.5.0" [[Gnuplot]] deps = ["ColorSchemes", "ColorTypes", "Colors", "DataStructures", "REPL", "ReplMaker", "StatsBase", "StructC14N", "Test"] path = "/home/gcalderone/.julia/dev/Gnuplot" Loading Loading @@ -313,6 +311,12 @@ git-tree-sha1 = "323a38ed1952d30586d0fe03412cde9399d3618b" uuid = "d8418881-c3e1-53bb-8760-2df7ec849ed5" version = "1.5.0" [[InvertedIndices]] deps = ["Test"] git-tree-sha1 = "15732c475062348b0165684ffe28e85ea8396afc" uuid = "41ab1584-1d38-5bbf-9106-f11c6c58b48f" version = "1.0.0" [[IterTools]] git-tree-sha1 = "05110a2ab1fc5f932622ffea2a003221f4782c18" uuid = "c8e1da08-722c-5040-9ed9-7db0dc04731e" Loading Loading @@ -370,21 +374,15 @@ uuid = "29816b5a-b9ab-546f-933c-edad1886dfa8" [[Libdl]] uuid = "8f399da3-3557-5675-b5ff-fb832c97cbdb" [[Libiconv_jll]] deps = ["Artifacts", "JLLWrappers", "Libdl", "Pkg"] git-tree-sha1 = "42b62845d70a619f063a7da093d995ec8e15e778" uuid = "94ce4f54-9a6c-5748-9c1c-f9c7231a4531" version = "1.16.1+1" [[LinearAlgebra]] deps = ["Libdl"] uuid = "37e2e46d-f89d-539d-b4ee-838fcccc9c8e" [[LogExpFunctions]] deps = ["DocStringExtensions", "LinearAlgebra"] git-tree-sha1 = "1ba664552f1ef15325e68dc4c05c3ef8c2d5d885" git-tree-sha1 = "7bd5f6565d80b6bf753738d2bc40a5dfea072070" uuid = "2ab3a3ac-af41-5b50-aa03-7779005ae688" version = "0.2.4" version = "0.2.5" [[Logging]] uuid = "56ddb016-857b-54e1-b83d-db4d58db5568" Loading Loading @@ -499,9 +497,15 @@ uuid = "44cfe95a-1eb2-52ea-b672-e2afdf69b78f" [[Polynomials]] deps = ["Intervals", "LinearAlgebra", "MutableArithmetics", "RecipesBase"] git-tree-sha1 = "3685cb1e2ccbbb9973684774f956f5c6e4673bc9" git-tree-sha1 = "0bbfdcd8cda81b8144de4be8a67f5717e959a005" uuid = "f27b6e38-b328-58d1-80ce-0feddd5e7a45" version = "2.0.12" version = "2.0.14" [[PooledArrays]] deps = ["DataAPI", "Future"] git-tree-sha1 = "cde4ce9d6f33219465b55162811d8de8139c0414" uuid = "2dfb63ee-cc39-5dd5-95bd-886bf059d720" version = "1.2.1" [[Preferences]] deps = ["TOML"] Loading @@ -526,7 +530,7 @@ uuid = "92933f4c-e287-5a05-a399-4b506db050ca" version = "1.7.1" [[QSFit]] deps = ["CMPFit", "Cosmology", "DSP", "DataStructures", "Dates", "DelimitedFiles", "Dierckx", "FITSIO", "GFit", "GFitViewer", "Pkg", "Printf", "Statistics", "Unitful", "UnitfulAstro"] deps = ["CMPFit", "Cosmology", "DSP", "DataFrames", "DataStructures", "Dates", "DelimitedFiles", "Dierckx", "FITSIO", "GFit", "GFitViewer", "Gnuplot", "Pkg", "Printf", "REPL", "Statistics", "TextParse", "Unitful", "UnitfulAstro"] path = "/home/gcalderone/.julia/dev/QSFit" uuid = "fea25315-b951-4667-83c9-50e53efab241" version = "0.1.0" Loading Loading @@ -606,9 +610,9 @@ uuid = "6462fe0b-24de-5631-8697-dd941f90decc" [[SortingAlgorithms]] deps = ["DataStructures"] git-tree-sha1 = "2ec1962eba973f383239da22e75218565c390a96" git-tree-sha1 = "b3363d7460f7d098ca0912c69b082f75625d7508" uuid = "a2af1166-a08f-5f64-846c-94a0d3cef48c" version = "1.0.0" version = "1.0.1" [[SparseArrays]] deps = ["LinearAlgebra", "Random"] Loading @@ -622,15 +626,15 @@ version = "1.5.1" [[Static]] deps = ["IfElse"] git-tree-sha1 = "2740ea27b66a41f9d213561a04573da5d3823d4b" git-tree-sha1 = "62701892d172a2fa41a1f829f66d2b0db94a9a63" uuid = "aedffcd0-7271-4cad-89d0-dc628f76c6d3" version = "0.2.5" version = "0.3.0" [[StaticArrays]] deps = ["LinearAlgebra", "Random", "Statistics"] git-tree-sha1 = "745914ebcd610da69f3cb6bf76cb7bb83dcb8c9a" git-tree-sha1 = "1b9a0f17ee0adde9e538227de093467348992397" uuid = "90137ffa-7385-5640-81b9-e52037218182" version = "1.2.4" version = "1.2.7" [[Statistics]] deps = ["LinearAlgebra", "SparseArrays"] Loading Loading @@ -694,10 +698,10 @@ uuid = "e0df1984-e451-5cb5-8b61-797a481e67e3" version = "1.0.1" [[TimeZones]] deps = ["Dates", "EzXML", "LazyArtifacts", "Mocking", "Pkg", "Printf", "RecipesBase", "Serialization", "Unicode"] git-tree-sha1 = "960099aed321e05ac649c90d583d59c9309faee1" deps = ["Dates", "Future", "LazyArtifacts", "Mocking", "Pkg", "Printf", "RecipesBase", "Serialization", "Unicode"] git-tree-sha1 = "81753f400872e5074768c9a77d4c44e70d409ef0" uuid = "f269a46b-ccf7-5d73-abea-4c690281aa53" version = "1.5.5" version = "1.5.6" [[TranscodingStreams]] deps = ["Random", "Test"] Loading @@ -714,9 +718,9 @@ uuid = "4ec0a83e-493e-50e2-b9ac-8f72acf5a8f5" [[Unitful]] deps = ["ConstructionBase", "Dates", "LinearAlgebra", "Random"] git-tree-sha1 = "b3682a0559219355f1e3c8024e9f97adce2d4623" git-tree-sha1 = "a981a8ef8714cba2fd9780b22fd7a469e7aaf56d" uuid = "1986cc42-f94f-5a68-af5c-568840ba703d" version = "1.8.0" version = "1.9.0" [[UnitfulAngles]] deps = ["Dates", "Unitful"] Loading @@ -736,12 +740,6 @@ git-tree-sha1 = "28807f85197eaad3cbd2330386fac1dcb9e7e11d" uuid = "ea10d353-3f73-51f8-a26c-33c1cb351aa5" version = "0.6.2" [[XML2_jll]] deps = ["Artifacts", "JLLWrappers", "Libdl", "Libiconv_jll", "Pkg", "Zlib_jll"] git-tree-sha1 = "1acf5bdf07aa0907e0a37d3718bb88d4b687b74a" uuid = "02c8fc9c-b97f-50b9-bbe4-9be30ff0a78a" version = "2.9.12+0" [[Zlib_jll]] deps = ["Libdl"] uuid = "83775a58-1f1d-513f-b197-d71354ab007a" Loading
run.jl +17 −16 Original line number Diff line number Diff line Loading @@ -11,7 +11,7 @@ using Dates using TextParse import StatsBase, JLD2 path = "output-" * readlines(`git branch --show-current`)[1] path = "output-new-" * readlines(`git branch --show-current`)[1] mkdir("$(path)") epoch_filenames = Vector{String}() for (root, dirs, files) in walkdir("AT2018zf") Loading @@ -32,8 +32,8 @@ function read_spec(epoch_id; kw...) spec = Spectrum(Val(:ASCII), file, columns=[1,2]; label="$epoch_id: $file", kw...); if isa(epoch_id, Number) spec.flux ./= all_scale[epoch_id] spec.err ./= all_scale[epoch_id] end spec.err .= 0.05 .* median(spec.flux); return spec end Loading @@ -48,13 +48,13 @@ I tried 180, 200 and 250, and the [OIII] luminosity do not change =# spec = read_spec("B03", resolution=200.) add_spec!(source, spec); res = fit(source); res = qsfit(source); viewer(res, showcomps=[:qso_cont, :galaxy, :balmer], filename="$(path)/results_boller.html") println(res.bestfit[:OIII_5007].norm.val ) # 0.0002913168348334021 println(res.bestfit[:OIII_5007].norm.val / res.bestfit[:galaxy].norm.val ) # 64.65011509608263 println(res.model[:OIII_5007].norm.val ) # 0.0002913095378291511 println(res.model[:OIII_5007].norm.val / res.model[:galaxy].norm.val ) # 61.34678931179987 Loading @@ -71,14 +71,14 @@ if !isfile("$(path)/scale_fwhm.dat") # source.options[:host_template] = "/home/gcalderone/Mbi1.30Zm1.49T00.0700_iTp0.00_baseFe_linear_FWHM_2.51" spec = read_spec(id) add_spec!(source, spec); res = fit(source); res = qsfit(source); viewer(res, showcomps=[:qso_cont, :galaxy, :balmer], filename="$(path)/results_$(id).html") # Update all_scale and store FWHM @info res.bestfit[:OIII_5007].norm.val all_scale[id] *= res.bestfit[:OIII_5007].norm.val all_fwhm[ id] = res.bestfit[:OIII_5007].fwhm.val @info res.model[:OIII_5007].norm.val all_scale[id] *= res.model[:OIII_5007].norm.val all_fwhm[ id] = res.model[:OIII_5007].fwhm.val catch end end Loading Loading @@ -150,7 +150,8 @@ for loop in 1:Nloop add_spec!(source, spec); @gp :- :zoom xr=[4750,5150] spec.λ ./ (1 + source.z) spec.flux "w l t '$(id)'" end res = multi_fit(source); res = qsfit_multi(source, ref_id=2); JLD2.jldsave(file; res) viewer(res, showcomps=[:qso_cont, :galaxy, :balmer], filename="$(path)/results_$(job)_$(loop).html") viewer(res, showcomps=[:qso_cont, :galaxy, :balmer], filename="$(path)/results_$(job)_$(loop)_rebin4.html", rebin=4) Loading @@ -166,12 +167,12 @@ for loop in 1:Nloop end end end bestfit2 = fit!(multi2, source.data, minimizer=GFit.cmpfit()) bestfit2 = fit!(res.source, multi2, res.pspecs) OIII_norm = fill(0., length(chosen_epochs)) for id in 1:length(chosen_epochs) OIII_norm[id] = bestfit2[id][:OIII_5007].norm.val all_scale[chosen_epochs[id]] *= bestfit2[id][:OIII_5007].norm.val OIII_norm[id] = multi2[id][:OIII_5007].norm.val all_scale[chosen_epochs[id]] *= multi2[id][:OIII_5007].norm.val end JLD2.jldsave(file; res, all_scale, OIII_norm) else Loading @@ -188,7 +189,7 @@ for loop in 1:Nloop file = "$(path)/results_$(job)_$(loop).dat" res, all_scale, OIII_norm = values(JLD2.load(file)) OIII_norm_evol = hcat(OIII_norm_evol, OIII_norm) println(res.bestfit.mdc.fitstat / res.bestfit.mdc.dof) println(res.fitres.fitstat / res.fitres.dof) end OIII_norm_evol = OIII_norm_evol[:,2:end] @gp "set grid" :- Loading
src/CL_1ES_1927p654.jl +224 −52 Original line number Diff line number Diff line module CL_1ES_1927p654 using QSFit, DataStructures, GFit, Statistics, Dierckx import QSFit: AbstractSpectralLine, DefaultRecipe, CombinedLine, NarrowLine, BroadLine, BroadBaseLine, ComboBroadLine, default_options, qso_cont_component, line_component, line_group_name, known_spectral_lines, fit, multi_fit, spectral_coverage, reduce import QSFit: DefaultRecipe, default_options, known_spectral_lines, LineComponent, add_qso_continuum!, add_patch_functs!, guess_emission_lines!, qsfit_multi abstract type q1927p654 <: DefaultRecipe end Loading @@ -22,63 +24,233 @@ function default_options(::Type{T}) where T <: q1927p654 end function qso_cont_component(::Type{T}) where T <: q1927p654 comp = qso_cont_component(supertype(T)) function QSFit.add_qso_continuum!(source::QSO{T}, pspec::PreparedSpectrum, model::Model) where T <: q1927p654 λ = domain(model)[:] comp = QSFit.powerlaw(3000) comp.x0.val = median(λ) comp.norm.val = Spline1D(λ, pspec.data.val, k=1, bc="error")(comp.x0.val) comp.alpha.val = -1.5 comp.alpha.low = -5 comp.alpha.high = 5 return comp end model[:qso_cont] = comp push!(model[:Continuum].list, :qso_cont) evaluate!(model) end struct AsymmNarrowLine <: AbstractSpectralLine λ::Float64 end line_group_name(::Type{T}, name::Symbol, line::AsymmNarrowLine) where T <: q1927p654 = :NarrowLines function known_spectral_lines(::Type{T}) where T <: q1927p654 list = OrderedDict{Symbol, AbstractSpectralLine}() #list[:Lyb ] = CombinedLine( 1026.0 ) list[:Lya ] = CombinedLine( 1215.24 ) list[:NV_1241 ] = NarrowLine( 1240.81 ) list[:OI_1306 ] = BroadLine( 1305.53 ) list[:CII_1335 ] = BroadLine( 1335.31 ) list[:SiIV_1400 ] = BroadLine( 1399.8 ) list[:CIV_1549 ] = CombinedLine( 1549.48 ) #list[:HeII ] = BroadLine( 1640.4 ) #list[:OIII ] = BroadLine( 1665.85 ) #list[:AlIII ] = BroadLine( 1857.4 ) list[:CIII_1909 ] = BroadLine( 1908.734) list[:CII ] = BroadLine( 2326.0 ) list[:F2420 ] = BroadLine( 2420.0 ) list[:MgII_2798 ] = CombinedLine( 2799.117) #list[:NeVN ] = NarrowLine( 3346.79 ) list[:NeVI_3426 ] = NarrowLine( 3426.85 ) list[:OII_3727 ] = NarrowLine( 3729.875) list[:NeIII_3869 ] = NarrowLine( 3869.81 ) #list[:Hd ] = BroadLine( 4102.89 ) list[:Hg ] = CombinedLine( 4341.68 ) list[:bb_Hg ] = BroadBaseLine( 4341.68 ) list[:OIII_4363 ] = NarrowLine( 4363.00 ) # TODO: Check wavelength is correct list[:HeII ] = BroadLine( 4686. ) list[:Hb ] = CombinedLine( 4862.68 ) list[:bb_Hb ] = BroadBaseLine( 4862.68 ) list[:OIII_4959 ] = NarrowLine( 4960.295) list[:OIII_5007 ] = NarrowLine( 5008.240) # AsymmNarrowLine #list[:OIII_5007_bw] = NarrowLine( 5008.240) list[:HeI_5876 ] = BroadLine( 5877.30 ) list[:OI_6300 ] = NarrowLine( 6300.00 ) # TODO: Check wavelength is correct list[:OI_6364 ] = NarrowLine( 6364.00 ) # TODO: Check wavelength is correct list[:NII_6549 ] = NarrowLine( 6549.86 ) list[:Ha ] = CombinedLine( 6564.61 ) list[:bb_Ha ] = BroadBaseLine( 6564.61 ) list[:NII_6583 ] = NarrowLine( 6585.27 ) list[:SII_6716 ] = NarrowLine( 6718.29 ) list[:SII_6731 ] = NarrowLine( 6732.67 ) function known_spectral_lines(source::QSO{T}) where T <: q1927p654 list = [ # MultiCompLine( :Lyb , [BroadLine, NarrowLine]), MultiCompLine( :Lya , [BroadLine, NarrowLine]), NarrowLine( :NV_1241 ), BroadLine( custom_transition(tid=:OI_1306 , 1305.53 )), BroadLine( custom_transition(tid=:CII_1335 , 1335.31 )), BroadLine( :SiIV_1400 ), MultiCompLine( :CIV_1549 , [BroadLine, NarrowLine]), # BroadLine( :HeII_1640 ), # BroadLine( custom_transition(tid=:OIII , 1665.85 )), # BroadLine( custom_transition(tid=:AlIII , 1857.4 )), BroadLine( :CIII_1909 ), BroadLine( custom_transition(tid=:CII , 2326.0 )), BroadLine( custom_transition(tid=:F2420 , 2420.0 )), MultiCompLine( :MgII_2798 , [BroadLine, NarrowLine]), # NarrowLine( custom_transition(tid=:NeVN , 3346.79 )), NarrowLine( custom_transition(tid=:NeVI_3426 , 3426.85 )), NarrowLine( custom_transition(tid=:OII_3727 , 3729.875)), NarrowLine( custom_transition(tid=:NeIII_3869 , 3869.81 )), # BroadLine( :Hd ), MultiCompLine( :Hg , [BroadLine, NarrowLine, BroadBaseLine]), NarrowLine( :OIII_4363 ), BroadLine( :HeII_4686 ), MultiCompLine( :Hb , [BroadLine, NarrowLine, BroadBaseLine]), NarrowLine( :OIII_4959 ), NarrowLine( :OIII_5007 ), # AsymmTailLine( :OIII_5007 , :blue), BroadLine( :HeI_5876 ), NarrowLine( :OI_6300 ), NarrowLine( :OI_6364 ), NarrowLine( :NII_6549 ), MultiCompLine( :Ha , [BroadLine, NarrowLine, BroadBaseLine]), NarrowLine( :NII_6583 ), NarrowLine( :SII_6716 ), NarrowLine( :SII_6731 )] return list end function QSFit.add_patch_functs!(source::QSO{T}, pspec::PreparedSpectrum, model::Model) where T <: q1927p654 QSFit.add_patch_functs!(parent_recipe(source), pspec, model) # Force Hg and Hb lines to have the same shape as Ha model[:Ha_br].norm.low = 1.e-10 # avoid division by zero @patch! begin model[:Hg_na].norm = model[:Hg_br].norm * (model[:Ha_na].norm / model[:Ha_br].norm) model[:Hg_na].voff = model[:Ha_na].voff model[:Hg_na].fwhm = model[:Ha_na].fwhm model[:Hg_bb].norm = model[:Hg_br].norm * (model[:Ha_bb].norm / model[:Ha_br].norm) model[:Hg_bb].voff = model[:Ha_bb].voff model[:Hg_bb].fwhm = model[:Ha_bb].fwhm model[:Hg_br].voff = model[:Ha_br].voff model[:Hg_br].fwhm = model[:Ha_br].fwhm end @patch! begin model[:Hb_na].norm = model[:Hb_br].norm * (model[:Ha_na].norm / model[:Ha_br].norm) model[:Hb_na].voff = model[:Ha_na].voff model[:Hb_na].fwhm = model[:Ha_na].fwhm model[:Hb_bb].norm = model[:Hb_br].norm * (model[:Ha_bb].norm / model[:Ha_br].norm) model[:Hb_bb].voff = model[:Ha_bb].voff model[:Hb_bb].fwhm = model[:Ha_bb].fwhm model[:Hb_br].voff = model[:Ha_br].voff model[:Hb_br].fwhm = model[:Ha_br].fwhm end end function qsfit_multi(source::QSO{TRecipe}; ref_id=1) where TRecipe <: q1927p654 elapsed = time() @assert length(source.specs) > 1 @assert 1 <= ref_id <= length(source.specs) pspecs = [PreparedSpectrum(source, id=id) for id in 1:length(source.specs)] multi = MultiModel() println(logio(source), "\nFit continuum components...") for id in 1:length(pspecs) pspec = pspecs[id] model = Model(pspec.domain) model[:Continuum] = SumReducer([]) model[:main] = SumReducer([]) push!(model[:main].list, :Continuum) select_reducer!(model, :main) delete!(model.revals, :default_sum) # TODO if source.options[:instr_broadening] # TODO GFit.set_instr_response!(model[1], (l, f) -> instrumental_broadening(l, f, source.spectra[id].resolution)) # TODO end QSFit.add_qso_continuum!(source, pspec, model) QSFit.add_host_galaxy!(source, pspec, model) QSFit.add_balmer_cont!(source, pspec, model) push!(multi, model) if id != ref_id @patch! multi[id][:galaxy].norm = multi[ref_id][:galaxy].norm end end fitres = fit!(source, multi, pspecs) for id in 1:length(pspecs) model = multi[id] pspec = pspecs[id] QSFit.renorm_cont!(source, pspec, model) freeze(model, :qso_cont) haskey(model, :galaxy) && freeze(model, :galaxy) haskey(model, :balmer) && freeze(model, :balmer) evaluate!(model) end evaluate!(multi) println(logio(source), "\nFit iron templates...") for id in 1:length(pspecs) model = multi[id] pspec = pspecs[id] model[:Iron] = SumReducer([]) push!(model[:main].list, :Iron) QSFit.add_iron_uv!( source, pspec, model) QSFit.add_iron_opt!(source, pspec, model) if length(model[:Iron].list) > 0 fitres = fit!(source, model, pspec) haskey(model, :ironuv ) && freeze(model, :ironuv) haskey(model, :ironoptbr) && freeze(model, :ironoptbr) haskey(model, :ironoptna) && freeze(model, :ironoptna) end evaluate!(model) end evaluate!(multi) println(logio(source), "\nFit known emission lines...") for id in 1:length(pspecs) model = multi[id] pspec = pspecs[id] QSFit.add_emission_lines!(source, pspec, model) guess_emission_lines!(source, pspec, model) QSFit.add_patch_functs!(source, pspec, model) multi[id][:OIII_5007].norm.val = 1. multi[id][:OIII_5007].norm.fixed = true end fitres = fit!(source, multi, pspecs) for id in 1:length(pspecs) model = multi[id] pspec = pspecs[id] for lname in keys(pspec.lcs) freeze(model, lname) end end evaluate!(multi) println(logio(source), "\nFit unknown emission lines...") for id in 1:length(pspecs) model = multi[id] pspec = pspecs[id] QSFit.add_unknown_lines!(source, pspec, model) end evaluate!(multi) println(logio(source), "\nLast run with all parameters free...") for id in 1:length(pspecs) model = multi[id] pspec = pspecs[id] thaw(model, :qso_cont) haskey(model, :galaxy ) && thaw(model, :galaxy) haskey(model, :balmer ) && thaw(model, :balmer) haskey(model, :ironuv ) && thaw(model, :ironuv) haskey(model, :ironoptbr) && thaw(model, :ironoptbr) haskey(model, :ironoptna) && thaw(model, :ironoptna) for lname in keys(pspec.lcs) thaw(model, lname) end for j in 1:source.options[:n_unk] cname = Symbol(:unk, j) if model[cname].norm.val > 0 thaw(model, cname) else freeze(model, cname) end end end fitres = fit!(source, multi, pspecs) rerun = false for id in 1:length(pspecs) model = multi[id] pspec = pspecs[id] rerun = rerun || QSFit.neglect_weak_features!(source, pspec, model, fitres) end if rerun println(logio(source), "\nRe-run fit...") fitres = fit!(source, multi, pspecs) end println(logio(source)) show(logio(source), fitres) out = QSFit.QSFitMultiResults(source, pspecs, multi, fitres) elapsed = time() - elapsed println(logio(source), "\nElapsed time: $elapsed s") close_logio(source) return out end #= function line_component(::Type{T}, line::AsymmNarrowLine) where T <: q1927p654 comp = QSFit.SpecLineAsymmGauss(line.λ) Loading @@ -101,7 +273,7 @@ function line_component(::Type{T}, line::ComboBroadLine) where T <: q1927p654 end =# include("Recipe_single.jl") include("Recipe_multi.jl") # include("Recipe_single.jl") # include("Recipe_multi.jl") end