Berechnung der Weltraumhintergrund-Röntgenstrahlung (CXB) mit Python durch Angabe des Flussbereichs

Hintergrund

Bei der Schätzung der Röntgenemission im Weltraumhintergrund (CXB) möchten wir nicht nur einen Teil der Gesamthelligkeit wissen, sondern auch, wie viel bei der Abbildung von Beobachtungen, ausgenommen der durch die Bildgebung erfassten Hellpunktquellen. Bleibt ohne räumliche Zersetzung? Muss geschätzt werden.

Hier zeige ich Ihnen, wie Sie mit Python gemäß https://iopscience.iop.org/article/10.1086/374335 weich (1-2 keV) und hart (2-10 keV) berechnen.

Quellcode

Informationen zum Lesen des Codes finden Sie unter Colab of calc_cxb.

python

#!/bin/evn python 

import matplotlib.pyplot as plt
import numpy as np

# THE RESOLVED FRACTION OF THE COSMIC X-RAY BACKGROUND A. Moretti, 
# The Astrophysical Journal, Volume 588, Number 2 (2003)
# https://iopscience.iop.org/article/10.1086/374335

smax = 1e-1 # maxinum of flux range (erg cm^-2 c^-1)
smin = 1e-17 # minimum of flux range (erg cm^-2 c^-1)
n = int(1e5) # number of grid in flux 
fcut = 5e-14 # user-specified maxinum flux (erg cm^-2 c^-1)

def calc_cxb(s,fcut=5e-14, ftotal_walker=2.18e-11, plot=True, soft=True):

	if soft: # 1-2 keV
		print("..... soft X-ray 1-2 keV is chosen")
		a1=1.82; a2=0.6; s0=1.48e-14; ns=6150
		foutname="soft"
	else: # hard 2-10 keV 
		print("..... hard X-ray 2-10 keV is chosen")
		a1=1.57; a2=0.44; s0=4.5e-15; ns=5300
		foutname="hard"


	n_larger_s = ns * 2.0e-15 ** a1 / (s**a1 + s0**(a1-a2)*s**a2) # Moretti et al. 2013, eq(2)

	if plot:
		F = plt.figure(figsize=(6,6))
		ax = plt.subplot(1,1,1)
		plt.plot(s, n_larger_s, ".", label="test")
		plt.xscale('log')
		plt.ylabel(r"$N(>S)/deg^2$")
		plt.xlabel(r"$S (erg cm^{-2} s^{-1})$")
		plt.yscale('log')
		plt.savefig(foutname + "_ns.png ")
		plt.show()

	diff_n_larger_s = np.abs(np.diff(n_larger_s))
	diff_s = np.diff(s)
	div_ns = np.abs(diff_n_larger_s/diff_s)

	s_lastcut = s[:-1]
	totalflux = np.sum(div_ns*s_lastcut*diff_s) # Moretti et al. 2013, eq(4)

	fluxcut = np.where(s > fcut)[0][:-1]
	particalflux = np.sum(div_ns[fluxcut]*s_lastcut[fluxcut]*diff_s[fluxcut])
	fdiff = totalflux - particalflux 
	fcxb = ftotal_walker - particalflux # Walker et al. 2016, eq(1)

	print("totalflux    = ", "%.4e"%totalflux, " [erg cm^-2 c^-1]", " from ", "%.4e"%smin, " to ",  "%.4e"%smax)
	print("particalflux = ", "%.4e"%particalflux, " [erg cm^-2 c^-1]", " from ", "%.4e"%fcut, " to ",  "%.4e"%smax)
	print("fdiff        = totalflux - particalflux     = ", "%.4e"%fdiff, " [erg cm^-2 c^-1]", " from ", "%.4e"%smin, " to ",  "%.4e"%fcut)	
	print("fcxb(2-10keV)= ftotal_walker - particalflux = ", "%.4e"%fcxb, " [erg cm^-2 c^-1]")

	if plot:
		F = plt.figure(figsize=(6,6))
		ax = plt.subplot(1,1,1)
		plt.plot(s[:-1], div_ns, ".", label="test")
		plt.xscale('log')
		plt.ylabel(r"$dN/dS/deg^2$")
		plt.xlabel(r"$S (erg cm^{-2} s^{-1})$")
		plt.yscale('log')
		plt.savefig(foutname + "_nsdiff.png ")		
		plt.show()

	return s, n_larger_s, diff_n_larger_s, div_ns, diff_s

#s = np.linspace(sexcl,smax,n)
s = np.logspace(np.log10(smin),np.log10(smax),n)

# plot soft band
s, n_larger_s, diff_n_larger_s, div_ns, diff_s = calc_cxb(s,fcut=fcut)
# plot hard band
s, n_larger_s, diff_n_larger_s, div_ns, diff_s = calc_cxb(s,fcut=fcut,soft=False)

Ausführungsergebnis

1-2 keV

2-10 keV

Berechnungsmethode

Wenn es linear ist, ist erwartungsgemäß nicht genügend Speicher vorhanden, um es zu berechnen. Im Protokoll des Raums s (erg cm ^ -2 s ^ -1) wurde für die Integration nur eine rechteckige Näherung verwendet.

Diese Funktion wird berechnet durch n_larger_s = ns * 2.0e-15 ** a1 / (s ** a1 + s0 ** (a1-a2) * s ** a2).

Diese Integration wird berechnet durch totalflux = np.sum (div_ns * s_lastcut * diff_s).

fcxb = ftotal_walker --particalflux # Walker et al. 2016, Gleichung (1) ist ein Bonus, und der CXB von 2-10 keV in diesem Artikel wird geschätzt und die Berechnung basiert darauf.