<SYSTEM_TASK:>
Given the following text description, write Python code to implement the functionality described below step by step
<END_TASK>
<USER_TASK:>
Description:
 Example joint fit between GBM and Swift BAT
This demonstrates 3ML's ability to to joint fit between two instruments with different likelihoods. Here, we have GBM data which obey a Poisson-Gaussian profile likelihoog (<a href=http
Step1:  Plugin setup
We have data from the same time interval from Swift BAT and a GBM NAI and BGO detector. We simply read in the files and select the energy range that we want to perform the fit over.
Step2:  Model setup
We setup up or spectral and likelihood model and combine the data. 3ML will automatically assign the proper likelihood to each data set.
Step3:  Spectral fitting
Now we simply fit the data
Step4:  It seems that the effective area between GBM and BAT do not agree!
Step5:  Let's add an effective area correction between the detectors to investigate the problem
Step6:  Now we have a much better fit to all data sets
Step7:  Examining the differences
We stored the Analysis Results of each analysis. We can save them to disk, but also, we can use them to see the differences between the resulting spectral fit.
Step8:  Let's plot the fits in model space and see how different the resulting models are.
<ASSISTANT_TASK:>
Python Code:
%matplotlib inline
%matplotlib notebook
from threeML import *
import os
Explanation: Example joint fit between GBM and Swift BAT
This demonstrates 3ML's ability to to joint fit between two instruments with different likelihoods. Here, we have GBM data which obey a Poisson-Gaussian profile likelihoog (<a href=http://heasarc.gsfc.nasa.gov/docs/xanadu/xspec/manual/node293.html> PGSTAT</a> in XSPEC lingo) and Swift BAT which data which are the result of a "fit" via a coded mask and hence obey a Gaussian ( $\chi^2$ ) likelihood.
3ML automatically probes the data and selects the proper likelihood for each data set. Details for the basics of plugin setup can be found in the basic_tests example. Here we concentrate on joint fitting.
End of explanation
gbm_dir = "gbm"
bat_dir = "bat"
bat = OGIPLike('BAT',
               observation=os.path.join(bat_dir,'gbm_bat_joint_BAT.pha'),
               response=os.path.join(bat_dir,'gbm_bat_joint_BAT.rsp'))
bat.set_active_measurements('15-150')
bat.view_count_spectrum()
nai6 = OGIPLike('n6',
                os.path.join(gbm_dir,'gbm_bat_joint_NAI_06.pha'),
                os.path.join(gbm_dir,'gbm_bat_joint_NAI_06.bak'),
                os.path.join(gbm_dir,'gbm_bat_joint_NAI_06.rsp'),
                spectrum_number=1)
nai6.set_active_measurements('8-900')
nai6.view_count_spectrum()
bgo0 = OGIPLike('b0',
                os.path.join(gbm_dir,'gbm_bat_joint_BGO_00.pha'),
                os.path.join(gbm_dir,'gbm_bat_joint_BGO_00.bak'),
                os.path.join(gbm_dir,'gbm_bat_joint_BGO_00.rsp'),
                spectrum_number=1)
bgo0.set_active_measurements('250-10000')
bgo0.view_count_spectrum()
Explanation: Plugin setup
We have data from the same time interval from Swift BAT and a GBM NAI and BGO detector. We simply read in the files and select the energy range that we want to perform the fit over.
End of explanation
band = Band()
model = Model(PointSource('joint_fit',0,0,spectral_shape=band))
data_list = DataList(bat,nai6,bgo0)
jl = JointLikelihood(model, data_list)
Explanation: Model setup
We setup up or spectral and likelihood model and combine the data. 3ML will automatically assign the proper likelihood to each data set.
End of explanation
_=jl.fit()
no_eac_results = jl.results
Explanation: Spectral fitting
Now we simply fit the data
End of explanation
_=display_spectrum_model_counts(jl,step=False)
Explanation: It seems that the effective area between GBM and BAT do not agree!
End of explanation
# turn on the effective area correction and set it's bounds
nai6.use_effective_area_correction(.2,1.8)
bgo0.use_effective_area_correction(.2,1.8)
# refit the data
_=jl.fit()
with_eac_res = jl.results
Explanation: Let's add an effective area correction between the detectors to investigate the problem:
End of explanation
_=display_spectrum_model_counts(jl,step=False)
Explanation: Now we have a much better fit to all data sets
End of explanation
no_eac_results.display()
with_eac_res.display()
Explanation: Examining the differences
We stored the Analysis Results of each analysis. We can save them to disk, but also, we can use them to see the differences between the resulting spectral fit.
End of explanation
plot_point_source_spectra(no_eac_results,with_eac_res,flux_unit='erg2/(keV s cm2)',equal_tailed=True)
Explanation: Let's plot the fits in model space and see how different the resulting models are.
End of explanation
<END_TASK>