phajy / Sesto2024

Preparation for the "From the Dolomites to the Event Horizon: Sledging Down the Black Hole Potential Well" 2024 meeting in Sesto. Will include talk, data, and calculations.
MIT License
1 stars 0 forks source link

Gradus spectral models #2

Closed fjebaker closed 2 months ago

fjebaker commented 2 months ago

Added the Gradus additive line profile model to the fit.

There's a few outstanding problems here:

Additionally the example transfer function table I've included in the data directory only has a 10x10 grid of parameters (done for brevity). Let's chat about what kind of models we want to fit, and the parameter ranges for the transfer function tables and then compute some new (denser) tables for the actual fitting.

I also ran JuliaFormatter over the script, so apologies for some unnecessary diffs!

fjebaker commented 2 months ago

The good news is for some of the XMM data (specifically index 2 and 3) we get fairly sensible fits. Both are along the lines of:

┌ CompositeModel with 2 model components:
│      a2 + a1
│ Model key and parameters:
│    a1 => LineProfile
│         K_1 ->  0.000139 ± 0.1  ∈ [ 0, Inf   ]   FREE
│         a_1 ->     0.963 ± 0.1  ∈ [ 0, 0.998 ]   FREE
│         θ_1 ->      33.6 ± 4    ∈ [ 0, Inf   ]   FREE
│       rin_1 ->         1                              FROZEN
│      rout_1 ->       100                              FROZEN
│        E₀_1 ->       6.4                              FROZEN
│    a2 => PowerLaw
│         K_2 ->   0.00796 ± 0.1  ∈ [ 0, Inf   ]   FREE
│         a_2 ->      1.83 ± 0.2  ∈ [ 0, Inf   ]   FREE
└ 

┌ FittingResult:
│   Model: CompositeModel[PowerLaw + LineProfile]
│   . u     : [0.00013859, 0.96276, 33.614, 0.0079636, 1.8282]
│   . σᵤ    : [1.2338e-05, 0.075803, 1.4552, 0.00018689, 0.015545]
│   . χ²    : 114.32 
└ 
fjebaker commented 2 months ago

AD now works! 🎉 we can remove the autodiff = :finite argument.

fjebaker commented 2 months ago

With your fix for the frozen parameter binding, we can now fit the whole thing with the Gradus models! Here's a preliminary attempt:

example

┌ MultiFittingResult:
│    Model: CompositeModel[PowerLaw + AutoCache]
│    . u     : [0.00018873, 0.95716, 34.729, 0.018017, 1.9156]
│    . σᵤ    : [9.3231e-06, 0.023654, 0.46583, 8.7455e-05, 0.0028567]
│    . χ²    : 200.40 
│    Model: CompositeModel[PowerLaw + AutoCache]
│    . u     : [0.00011664, 0.95716, 34.729, 0.011149, 1.9156]
│    . σᵤ    : [7.1866e-06, 0.023654, 0.46583, 5.9531e-05, 0.0028567]
│    . χ²    : 163.38 
│    Model: CompositeModel[PowerLaw + AutoCache]
│    . u     : [0.00016146, 0.95716, 34.729, 0.0090218, 1.9156]
│    . σᵤ    : [1.1121e-05, 0.023654, 0.46583, 6.7429e-05, 0.0028567]
│    . χ²    : 145.19 
│    Model: CompositeModel[PowerLaw + AutoCache]
│    . u     : [0.00024576, 0.95716, 34.729, 0.016804, 1.9156]
│    . σᵤ    : [3.0203e-05, 0.023654, 0.46583, 0.00017771, 0.0028567]
│    . χ²    : 348.35 
│    Model: CompositeModel[PowerLaw + AutoCache]
│    . u     : [0.00026623, 0.95716, 34.729, 0.016361, 1.9156]
│    . σᵤ    : [3.0455e-05, 0.023654, 0.46583, 0.00017768, 0.0028567]
│    . χ²    : 337.94 
│    Model: CompositeModel[PowerLaw + AutoCache]
│    . u     : [0.00020407, 0.95716, 34.729, 0.017525, 1.9156]
│    . σᵤ    : [1.2407e-05, 0.023654, 0.46583, 0.00012182, 0.0028567]
│    . χ²    : 1152.3 
│    Model: CompositeModel[PowerLaw + AutoCache]
│    . u     : [0.00020958, 0.95716, 34.729, 0.017174, 1.9156]
│    . σᵤ    : [1.2456e-05, 0.023654, 0.46583, 0.00012002, 0.0028567]
│    . χ²    : 1074.3 
│    Model: CompositeModel[PowerLaw + AutoCache]
│    . u     : [0.00011123, 0.95716, 34.729, 0.011953, 1.9156]
│    . σᵤ    : [2.1734e-05, 0.023654, 0.46583, 0.00012733, 0.0028567]
│    . χ²    : 703.13 
│    Model: CompositeModel[PowerLaw + AutoCache]
│    . u     : [9.6993e-05, 0.95716, 34.729, 0.011836, 1.9156]
│    . σᵤ    : [2.2030e-05, 0.023654, 0.46583, 0.00012830, 0.0028567]
│    . χ²    : 542.06 
└ Σχ² = 4667.1
fjebaker commented 2 months ago

Performance leaves a little to be desired:

result = @time fit(prob, LevenbergMarquadt(); verbose = true)
# looking to improve on
#  78.044531 seconds (9.23 M allocations: 1.073 GiB, 0.23% gc time, 28.93% compilation time)
#  55.780932 seconds (38.87 k allocations: 513.027 MiB, 0.24% gc time)
phajy commented 2 months ago

This looks great. I'll run through this locally and then get back to you with any comments (if any) before merging. Should we close and merge this one before looking at PR https://github.com/phajy/Sesto2024/pull/3?

fjebaker commented 2 months ago

3 is merging into this branch. But the order is unimportant, as we can change which branch to merge into. Probably makes more sense to resolve this one first, then look at #3