def noise_power(self, freq):
"""Returns a function to calculate the noise PS at the given freq.
"""
z = freq_to_z(freq)
beam_size = self.beam_size(freq)
A_pix = beam_size**2
A_survey = self.f_sky * 4 * np.pi
tau = A_pix / A_survey * units.year * self.num_year * self.beam_num
# Calculate the comoving size of a frequency bin (at the given freq)
d = self.proper_distance
dxf = (d(freq_to_z(freq - self.freq_width)) - d(z))
# Define the window function in k-space for the parallel and perpendicular directions.
# Use a sinc function for parallel as it is the FT of a top-hat bin. This is probably a bad choice.
def window_par(kpar):
y = kpar * dxf / (4 * np.pi)
return np.sinc(y) * (np.abs(y) < 1.0)
# Azimuthally average over the X and Y window functions to produce an
# overall k_perp window function. Do this by averaging for a set number
# of points k values and then generating an interpolating function to
# appeoximate the full result.
def _int(phi, k):
# Integrand to average over
x = (3e2 * k * d(z)) / (freq * 2 * np.pi)
xx = x * np.cos(phi)
xy = x * np.sin(phi)
return (self.window_x(xx) * self.window_y(xy))**2
def _w_xy_average(k):
# Full averaged window function
return scipy.integrate.fixed_quad(_int, 0, 2 * np.pi, args=(k,), n=1024)[0]
# Generate a log interpolated approximation
k_val = np.linspace(0, self.kmax, 256)
int_val = np.array([_w_xy_average(k)**0.5 for k in k_val])
_w_perp_interp = scipy.interpolate.interp1d(k_val, np.log(int_val))
def window_perp(kperp):
return np.exp(_w_perp_interp(kperp))
# Calculate the comoving volume of a single pixel (beam)
V_pix = A_pix * d(z)**2 * dxf
# Receiver temperature contribution to instrumental Stokes I
T_recv_I = self.T_recv / 2**0.5
return inv_noise_ps_21cm(T_recv_I + self.T_sky(freq), tau, V_pix,
self.freq_width, window_par, window_perp)
评论列表
文章目录