

lle = TimeSeriesBuilder.from_lat_lle('lle',ft2_file="lle_pt.fit",
                                     lle_file="lle.fit",
                                     rsp_file="lle.rsp")
lle.set_background_interval('-100--10','150-500')

lle.set_active_time_interval('68-110')
lle.view_lightcurve(-200,500);


lle_plugin = lle.to_spectrumlike()
lle_plugin.use_effective_area_correction()

lle_plugin.display()

lle_plugin.view_count_spectrum();

"""
Explanation: Data setup
We are going to use data from Fermi-LAT, Fermi-GBM and Swift-XRT. Let's go through the process of setting up the data from each instrument. We will work from high energy to low energy.
Fermi-LAT
Once we have obtained the Fermi-LAT data, in this case, the LAT Low Energy (LLE) data, we can reduce the data into an plugin using the light curve tools provided in 3ML. LLE data is in the format of FITS event files with an associated spacecraft point history file and energy dispersion response. The TimeSeriesBuilder class has special methods for dealing with the LLE data.
End of explanation
"""


gbm_detectors = ["n4", "n7", "n8", "b0"]

for det in gbm_detectors:
    ts_cspec = TimeSeriesBuilder.from_gbm_cspec_or_ctime(
        det, cspec_or_ctime_file=f"cspec_{det}.pha", rsp_file=f"cspec_{det}.rsp2"
    )
    ts_cspec.set_background_interval("-400--10", "700-1200")

    ts_cspec.save_background(filename=f"{det}_bkg.h5", overwrite=True)

gbm_time_series = {}
gbm_plugins = {}
for det in gbm_detectors:


    ts = TimeSeriesBuilder.from_gbm_tte(
        det,
        tte_file=f"tte_{det}.fit.gz",
        rsp_file=f"cspec_{det}.rsp2",
        restore_background=f"{det}_bkg.h5",
    )

    gbm_time_series[det] = ts

    ts.view_lightcurve(-10, 200)
    ts.set_active_time_interval("68-110")
    gbm_plugins[det] = ts.to_spectrumlike()

for det, plugin in gbm_plugins.items():

    if det.startswith("b"):

        plugin.set_active_measurements("250-30000")

    else:

        plugin.set_active_measurements("10-900")
        
    if det != "n3":
        
        plugin.use_effective_area_correction()

    plugin.rebin_on_background(1)
    plugin.view_count_spectrum()

"""
Explanation: Fermi-GBM
End of explanation
"""


xrt = SwiftXRTLike(
    "xrt",
    observation="awt.pi",
    background="awtback.pi",
    arf_file="awt.arf",
    response="awt.rmf",
)
xrt.display()

xrt.remove_rebinning()
xrt.set_active_measurements('4-10')

xrt.rebin_on_background(1.)
xrt.use_effective_area_correction()
xrt.view_count_spectrum();


"""
Explanation: Swift-XRT
For Swift-XRT, we can use the normal OGIPLike plugin, but the energy resolution of the instrument is so fine that we would waste time integrating over the photon bins during forward-folding. Thus, there is a special plugin that overrides the computation of the photon integrals with a simple sum.
End of explanation
"""


all_plugins = [lle_plugin, xrt]
for _ , plugin in gbm_plugins.items():
    
    all_plugins.append(plugin)


datalist = DataList(*all_plugins)

"""
Explanation: Combining all the plugins
End of explanation
"""


sbpl = SmoothlyBrokenPowerLaw(pivot=1E3)

sbpl.alpha.prior =  Truncated_gaussian(lower_bound=-1.5, upper_bound=0, mu=-1, sigma=.5)
sbpl.beta.prior =  Truncated_gaussian(lower_bound=-3., upper_bound=-1.6, mu=-2, sigma=.5)
sbpl.break_energy.prior = Log_uniform_prior(lower_bound=1, upper_bound=1E3)
sbpl.break_scale.prior =  Log_uniform_prior(lower_bound=1E-4, upper_bound=10.)
sbpl.K.prior = Log_uniform_prior(lower_bound=1E-2, upper_bound=1E2)
sbpl.K = 1E-1
sbpl.break_energy.bounds = (0, None)

sbpl.break_scale.free=True

ps = PointSource('grb',0,0,spectral_shape=sbpl)

model = Model(ps)

bayes = BayesianAnalysis(model,datalist)
bayes.set_sampler('multinest')

for k,v in model.free_parameters.items():
    if "cons" in k:
        
        v.prior = Truncated_gaussian(lower_bound=.8, upper_bound=1.2, mu=1, sigma=.1)
        
    


#bayes.sampler.setup(dlogz=10.,frac_remain=.5)
bayes.sampler.setup(n_live_points=1000)

bayes.sample()


bayes.restore_median_fit()
#sbpl.K = 1E-1
display_spectrum_model_counts(bayes, min_rate=[-1,5,5,5,5,5]);

bayes.results.corner_plot();

plot_point_source_spectra(bayes.results, flux_unit='erg2/(cm2 s keV)',ene_max =1E5);

"""
Explanation: Fitting
Model setup
Band Function
End of explanation
"""


