Python DataPHA.set_arf方法代码示例

# 需要导入模块: from sherpa.astro.data import DataPHA [as 别名]
# 或者: from sherpa.astro.data.DataPHA import set_arf [as 别名]
def test_rsp_normf_error(analysis):
    """Check that an error is raised on set_analysis

    """

    exposure = 200.1

    # rdata is only used to define the grids
    rdata = create_non_delta_rmf()
    specresp = create_non_delta_specresp()
    adata = create_arf(rdata.energ_lo,
                       rdata.energ_hi,
                       specresp,
                       exposure=exposure)

    nchans = rdata.e_min.size
    channels = np.arange(1, nchans + 1, dtype=np.int16)
    counts = np.ones(nchans, dtype=np.int16)
    pha = DataPHA('test-pha', channel=channels, counts=counts,
                  exposure=exposure)
    pha.set_arf(adata)

    with pytest.raises(DataErr) as exc:
        pha.set_analysis(analysis)

    emsg = "response incomplete for dataset test-pha, " + \
           "check the instrument model"
    assert str(exc.value) == emsg

# 需要导入模块: from sherpa.astro.data import DataPHA [as 别名]
# 或者: from sherpa.astro.data.DataPHA import set_arf [as 别名]
def test_arfmodelpha_call(ignore):
    """What happens calling an arf with a pha?

    The ignore value indicates what channel to ignore (0 means
    nothing is ignored). The aim is to check edge effects,
    and as there are only a few channels, it was decided to
    test all channels.
    """

    # Note: the exposure is set in the PHA and ARF, but should not be
    #       used when evaluating the model; it's value has been
    #       set to a value that the test will fail it it is.
    #
    exposure = 200.1
    estep = 0.01
    egrid = np.arange(0.01, 0.06, estep)
    svals = [1.1, 1.2, 1.3, 1.4]
    specresp = np.asarray(svals)
    adata = create_arf(egrid[:-1], egrid[1:], specresp,
                       exposure=exposure)

    constant = 2.3
    mdl = Const1D('flat')
    mdl.c0 = constant

    channels = np.arange(1, 5, dtype=np.int16)
    counts = np.asarray([10, 5, 12, 7], dtype=np.int16)
    pha = DataPHA('test-pha', channel=channels, counts=counts,
                  exposure=exposure)
    pha.set_arf(adata)

    # force energy units (only needed if ignore is set)
    pha.set_analysis('energy')

    if ignore is not None:
        de = estep * 0.9
        e0 = egrid[ignore]
        pha.notice(lo=e0, hi=e0 + de, ignore=True)

        # The assert are intended to help people reading this
        # code rather than being a useful check that the code
        # is working.
        mask = [True, True, True, True]
        mask[ignore] = False
        assert (pha.mask == mask).all()

    wrapped = ARFModelPHA(adata, pha, mdl)

    # The model is evaluated on the ARF grid, not whatever
    # is sent in. It is also integrated across the bins,
    # which is why there is a multiplication by the
    # grid width (for this constant model).
    #
    # Note that the filter doesn't change the grid.
    #
    de = egrid[1:] - egrid[:-1]
    expected = constant * np.asarray(svals) * de
    out = wrapped([4, 5])
    assert_allclose(out, expected)

# 需要导入模块: from sherpa.astro.data import DataPHA [as 别名]
# 或者: from sherpa.astro.data.DataPHA import set_arf [as 别名]
def test_rsp1d_matrix_pha_zero_energy_bin():
    """What happens when the first bin starts at 0, with replacement.

    Unlike test_rsp1d_delta_pha_zero_energy_bin this directly
    calls Response1D to create the model.
    """

    ethresh = 1.0e-5

    rdata = create_non_delta_rmf()

    # hack the first bin to have 0 energy
    rdata.energ_lo[0] = 0.0

    # PHA and ARF have different exposure ties
    exposure_arf = 0.1
    exposure_pha = 2.4

    specresp = create_non_delta_specresp()

    with warnings.catch_warnings(record=True) as ws:
        warnings.simplefilter("always")
        adata = create_arf(rdata.energ_lo,
                           rdata.energ_hi,
                           specresp,
                           exposure=exposure_arf,
                           ethresh=ethresh)

    validate_zero_replacement(ws, 'ARF', 'user-arf', ethresh)

    nchans = rdata.e_min.size
    channels = np.arange(1, nchans + 1, dtype=np.int16)
    counts = np.ones(nchans, dtype=np.int16)
    pha = DataPHA('test-pha', channel=channels, counts=counts,
                  exposure=exposure_pha)
    pha.set_rmf(rdata)
    pha.set_arf(adata)

    pha.set_analysis('energy')

    mdl = MyPowLaw1D()

    rsp = Response1D(pha)
    wrapped = rsp(mdl)

    # Evaluate the statistic / model. The value was calculated using
    # commit a65fb94004664eab219cc09652172ffe1dad80a6 on a linux
    # system (Ubuntu 17.04).
    #
    f = Fit(pha, wrapped)
    ans = f.calc_stat()
    assert ans == pytest.approx(37971.8716151947)

# 需要导入模块: from sherpa.astro.data import DataPHA [as 别名]
# 或者: from sherpa.astro.data.DataPHA import set_arf [as 别名]
def test_rspmodelpha_matrix_call_xspec():
    """Check XSPEC constant is invariant to wavelength/energy setting.

    As XSPEC models internally convert from Angstrom to keV,
    do a simple check here.
    """

    exposure = 200.1
    rdata = create_non_delta_rmf()
    specresp = create_non_delta_specresp()
    adata = create_arf(rdata.energ_lo,
                       rdata.energ_hi,
                       specresp,
                       exposure=exposure)

    constant = 2.3
    mdl = XSconstant('flat')
    mdl.factor = constant

    nchans = rdata.e_min.size
    channels = np.arange(1, nchans + 1, dtype=np.int16)
    counts = np.ones(nchans, dtype=np.int16)
    pha = DataPHA('test-pha', channel=channels, counts=counts,
                  exposure=exposure)

    # The set_arf call isn't necessary, but leave in
    pha.set_arf(adata)
    pha.set_rmf(rdata)

    # The XSPEC models are evaluated on an energy grid, even when
    # the analysis setting is wavelength. Also, unlike the Sherpa
    # Constant model, the XSPEC XSconstant model is defined
    # over the integrated bin, so no correction is needed for the
    # bin width.
    #
    modvals = constant * specresp
    matrix = get_non_delta_matrix()
    expected = np.matmul(modvals, matrix)

    wrapped = RSPModelPHA(adata, rdata, pha, mdl)

    pha.set_analysis('wave')
    out_wl = wrapped([4, 5])
    assert_allclose(out_wl, expected)

    pha.set_analysis('energy')
    out_en = wrapped([4, 5])
    assert_allclose(out_en, expected)

# 需要导入模块: from sherpa.astro.data import DataPHA [as 别名]
# 或者: from sherpa.astro.data.DataPHA import set_arf [as 别名]
def test_rsp1d_delta_pha_zero_energy_bin():
    "What happens when the first bin starts at 0, with replacement"

    ethresh = 2.0e-7

    # PHA and ARF have different exposure ties
    exposure1 = 0.1
    exposure2 = 2.4
    egrid = np.asarray([0.0, 0.1, 0.2, 0.4, 0.5, 0.7, 0.8])
    elo = egrid[:-1]
    ehi = egrid[1:]
    specresp = np.asarray([10.2, 9.8, 10.0, 12.0, 8.0, 10.0])

    with warnings.catch_warnings(record=True) as ws:
        warnings.simplefilter("always")
        adata = create_arf(elo, ehi, specresp, exposure=exposure1,
                           ethresh=ethresh)

    validate_zero_replacement(ws, 'ARF', 'user-arf', ethresh)

    with warnings.catch_warnings(record=True) as ws:
        warnings.simplefilter("always")
        rdata = create_delta_rmf(elo, ehi, ethresh=ethresh)

    validate_zero_replacement(ws, 'RMF', 'delta-rmf', ethresh)

    channels = np.arange(1, 7, dtype=np.int16)
    counts = np.ones(6, dtype=np.int16)
    pha = DataPHA('test-pha', channel=channels, counts=counts,
                  exposure=exposure2)
    pha.set_rmf(rdata)
    pha.set_arf(adata)

    pha.set_analysis('energy')

    mdl = MyPowLaw1D()
    tmdl = PowLaw1D()

    wrapped = RSPModelPHA(adata, rdata, pha, mdl)

    out = wrapped([0.1, 0.2])

    elo[0] = ethresh
    expected = specresp * tmdl(elo, ehi)

    assert_allclose(out, expected)
    assert not np.isnan(out[0])

# 需要导入模块: from sherpa.astro.data import DataPHA [as 别名]
# 或者: from sherpa.astro.data.DataPHA import set_arf [as 别名]
def test_arf1d_pha_zero_energy_bin():
    "What happens when the first bin starts at 0, with replacement"

    ethresh = 1.0e-10

    # Note: the two exposures are different to check which is
    #       used (the answer is neither, which seems surprising)
    #
    exposure1 = 0.1
    egrid = np.asarray([0.0, 0.1, 0.2, 0.4, 0.5, 0.7, 0.8])
    elo = egrid[:-1]
    ehi = egrid[1:]
    specresp = np.asarray([10.2, 9.8, 10.0, 12.0, 8.0, 10.0])

    with warnings.catch_warnings(record=True) as ws:
        warnings.simplefilter("always")
        adata = create_arf(elo, ehi, specresp, exposure=exposure1,
                           ethresh=ethresh)

    validate_zero_replacement(ws, 'ARF', 'user-arf', ethresh)

    arf = ARF1D(adata)

    exposure2 = 2.4
    channels = np.arange(1, 7, dtype=np.int16)
    counts = np.ones(6, dtype=np.int16)
    pha = DataPHA('test-pha', channel=channels, counts=counts,
                  exposure=exposure2)
    pha.set_arf(adata)

    pha.set_analysis('energy')

    mdl = MyPowLaw1D()
    tmdl = PowLaw1D()

    wrapped = ARFModelPHA(arf, pha, mdl)

    out = wrapped([0.1, 0.2])
    elo[0] = ethresh
    expected = specresp * tmdl(elo, ehi)

    assert_allclose(out, expected)
    assert not np.isnan(out[0])

# 需要导入模块: from sherpa.astro.data import DataPHA [as 别名]
# 或者: from sherpa.astro.data.DataPHA import set_arf [as 别名]
def test_rsp_normf_call(arfexp, phaexp):
    """Check out Response1D with no RMF.

    analysis is the analysis setting
    arfexp determines whether the arf has an exposure time
    phaexp determines whether the PHA has an exposure time

    This only uses the channel setting
    """

    # Chose different exposure times for ARF and PHA to see which
    # gets picked up.
    #
    if arfexp:
        arf_exposure = 200.1
    else:
        arf_exposure = None

    if phaexp:
        pha_exposure = 220.9
    else:
        pha_exposure = None

    if phaexp:
        exposure = pha_exposure
        mdl_label = '({} * flat)'.format(exposure)
    elif arfexp:
        exposure = arf_exposure
        mdl_label = '({} * flat)'.format(exposure)
    else:
        exposure = 1.0
        mdl_label = 'flat'

    # rdata is only used to define the grids
    rdata = create_non_delta_rmf()
    specresp = create_non_delta_specresp()
    adata = create_arf(rdata.energ_lo,
                       rdata.energ_hi,
                       specresp,
                       exposure=arf_exposure)

    constant = 2.3
    mdl = Const1D('flat')
    mdl.c0 = constant

    # Turn off integration on this model, so that it is not integrated
    # across the bin width.
    #
    mdl.integrate = False

    nchans = rdata.e_min.size
    channels = np.arange(1, nchans + 1, dtype=np.int16)
    counts = np.ones(nchans, dtype=np.int16)
    pha = DataPHA('test-pha', channel=channels, counts=counts,
                  exposure=pha_exposure)

    pha.set_arf(adata)

    rsp = Response1D(pha)
    wrapped = rsp(mdl)

    assert isinstance(wrapped, ArithmeticModel)

    expname = 'apply_arf({})'.format(mdl_label)
    assert wrapped.name == expname

    expected = exposure * constant * specresp

    pha.set_analysis('channel')
    out = wrapped([4, 5])
    assert_allclose(out, expected)