Skip to content

Commit

Permalink
changed all references of "stokes" to "pol" across repo
Browse files Browse the repository at this point in the history
  • Loading branch information
@nkern
nkern committed Apr 27, 2018
1 parent 3aa945a commit c069177
Show file tree
Hide file tree
Showing 5 changed files with 48 additions and 60 deletions.
12 changes: 6 additions & 6 deletions hera_pspec/noise.py
View file
Original file line number Diff line number Diff line change
Expand Up @@ -67,7 +67,7 @@ def add_beam(self, beam):

self.beam = beam

def calc_scalar(self, freqs, stokes, num_steps=5000, little_h=True):
def calc_scalar(self, freqs, pol, num_steps=5000, little_h=True):
"""
Calculate noise power spectrum prefactor from Eqn. (1) of Pober et al. 2014, ApJ 782, 66,
equal to
Expand All @@ -78,7 +78,7 @@ def calc_scalar(self, freqs, stokes, num_steps=5000, little_h=True):
----------
freqs : float ndarray, holds frequency bins of spectral window in Hz
stokes : str, specification of (pseudo) Stokes polarization, or linear dipole polarization
pol : str, specification of polarization to calculate scalar for
See pyuvdata.utils.polstr2num for options.
num_steps : number of frequency bins to use in numerical integration of scalar
Expand All @@ -93,13 +93,13 @@ def calc_scalar(self, freqs, stokes, num_steps=5000, little_h=True):
self.subband : float ndarray, frequencies in spectral window used to calculate self.scalar
self.stokes : str, stokes polarization used to calculate self.scalar
self.pol : str, polarization used to calculate self.scalar
"""
# parse stokes
# compute scalar
self.scalar = self.beam.compute_pspec_scalar(freqs.min(), freqs.max(), len(freqs), num_steps=num_steps,
stokes=stokes, little_h=little_h, noise_scalar=True)
pol=pol, little_h=little_h, noise_scalar=True)
self.subband = freqs
self.stokes = stokes
self.pol = pol

def calc_P_N(self, Tsys, t_int, Ncoherent=1, Nincoherent=None, form='Pk', k=None, little_h=True):
"""
Expand Down
64 changes: 27 additions & 37 deletions hera_pspec/pspecbeam.py
View file
Original file line number Diff line number Diff line change
Expand Up @@ -101,15 +101,15 @@ def __init__(self, cosmo=None):
else:
self.cosmo = conversions.Cosmo_Conversions()

def compute_pspec_scalar(self, lower_freq, upper_freq, num_freqs, num_steps=5000, stokes='I',
def compute_pspec_scalar(self, lower_freq, upper_freq, num_freqs, num_steps=5000, pol='I',
taper='none', little_h=True, noise_scalar=False):
"""
Computes the scalar function to convert a power spectrum estimate
in "telescope units" to cosmological units
See arxiv:1304.4991 and HERA memo #27 for details.
Currently, only the Stokes "I", "XX" and "YY" beam are supported.
Currently, only the "I", "XX" and "YY" polarization beams are supported.
See Equations 4 and 5 of Moore et al. (2017) ApJ 836, 154
or arxiv:1502.05072 for details.
Expand All @@ -131,10 +131,9 @@ def compute_pspec_scalar(self, lower_freq, upper_freq, num_freqs, num_steps=5000
numerical integral
Default: 10000
stokes: str, optional
Which Stokes parameter's to compute the beam scalar for.
pol: str, optional
Which polarization to compute the beam scalar for.
'I', 'Q', 'U', 'V', 'XX', 'YY', 'XY', 'YX', although
currently only 'I' and linear dipole pol (e.g. 'XX') are implemented.
Default: 'I'
taper : str, optional
Expand Down Expand Up @@ -162,7 +161,7 @@ def compute_pspec_scalar(self, lower_freq, upper_freq, num_freqs, num_steps=5000
pspec_freqs = np.linspace(lower_freq, upper_freq, num_freqs, endpoint=False)

# Get omega_ratio
omega_ratio = self.power_beam_sq_int(stokes) / self.power_beam_int(stokes)**2
omega_ratio = self.power_beam_sq_int(pol) / self.power_beam_int(pol)**2

# Get scalar
scalar = _compute_pspec_scalar(self.cosmo, self.beam_freqs, omega_ratio, pspec_freqs,
Expand All @@ -171,7 +170,7 @@ def compute_pspec_scalar(self, lower_freq, upper_freq, num_freqs, num_steps=5000

return scalar

def Jy_to_mK(self, freqs, stokes='I'):
def Jy_to_mK(self, freqs, pol='I'):
"""
Return the multiplicative factor, M [mK / Jy], to convert a visibility from Jy -> mK,
Expand All @@ -186,10 +185,9 @@ def Jy_to_mK(self, freqs, stokes='I'):
----------
freqs : float ndarray, contains frequencies to evaluate conversion factor [Hz]
stokes: str, optional
Which Stokes parameter's to compute the beam scalar for.
pol: str, optional
Which polarization to compute the beam scalar for.
'I', 'Q', 'U', 'V', 'XX', 'YY', 'XY', 'YX', although
currently only 'I' and linear dipole pol (e.g. 'XX') are implemented.
Default: 'I'
Returns
Expand All @@ -207,7 +205,7 @@ def Jy_to_mK(self, freqs, stokes='I'):
if np.max(freqs) > self.beam_freqs.max():
raise ValueError("Warning: max freq {} > self.beam_freqs.max(), extrapolating...".format(np.max(freqs)))

Op = interp1d(self.beam_freqs/1e6, self.power_beam_int(stokes=stokes), kind='quadratic', fill_value='extrapolate')(freqs/1e6)
Op = interp1d(self.beam_freqs/1e6, self.power_beam_int(pol=pol), kind='quadratic', fill_value='extrapolate')(freqs/1e6)

return 1e-20 * conversions.cgs_units.c**2 / (2 * conversions.cgs_units.kb * freqs**2 * Op)

Expand All @@ -233,7 +231,7 @@ def __init__(self, fwhm, beam_freqs, cosmo=None):
else:
self.cosmo = conversions.Cosmo_Conversions()

def power_beam_int(self, stokes='I'):
def power_beam_int(self, pol='I'):
"""
Computes the integral of the beam over solid angle to give
a beam area (in sr). Uses analytic formula that the answer
Expand All @@ -242,16 +240,14 @@ def power_beam_int(self, stokes='I'):
Trivially this returns an array (i.e., a function of frequency),
but the results are frequency independent.
Currently, only the Stokes "I", "XX" and "YY" beam are supported.
See Equations 4 and 5 of Moore et al. (2017) ApJ 836, 154
or arxiv:1502.05072 for details.
Parameters
----------
stokes: str, optional
Which Stokes parameter's to compute the beam scalar for.
'I', 'Q', 'U', 'V', 'XX', 'YY', 'XY', 'YX', although
currently only 'I' and linear dipole pol (e.g. 'XX') are implemented.
pol: str, optional
Which polarization to compute the beam scalar for.
'I', 'Q', 'U', 'V', 'XX', 'YY', 'XY', 'YX'
Default: 'I'
Returns
Expand All @@ -260,7 +256,7 @@ def power_beam_int(self, stokes='I'):
"""
return np.ones_like(self.beam_freqs) * 2. * np.pi * self.fwhm**2 / (8. * np.log(2.))

def power_beam_sq_int(self, stokes='I'):
def power_beam_sq_int(self, pol='I'):
"""
Computes the integral of the beam**2 over solid angle to give
a beam area (in str). Uses analytic formula that the answer
Expand All @@ -269,16 +265,14 @@ def power_beam_sq_int(self, stokes='I'):
Trivially this returns an array (i.e., a function of frequency),
but the results are frequency independent.
Currently, only the Stokes "I", "XX" and "YY" beam are supported.
See Equations 4 and 5 of Moore et al. (2017) ApJ 836, 154
or arxiv:1502.05072 for details.
Parameters
----------
stokes: str, optional
Which Stokes parameter's to compute the beam scalar for.
'I', 'Q', 'U', 'V', 'XX', 'YY', 'XY', 'YX', although
currently only 'I' and linear dipole pol (e.g. 'XX') are implemented.
pol: str, optional
Which polarization to compute the beam scalar for.
'I', 'Q', 'U', 'V', 'XX', 'YY', 'XY', 'YX'
Default: 'I'
Returns
Expand Down Expand Up @@ -310,56 +304,52 @@ def __init__(self, beam_fname, cosmo=None):
else:
self.cosmo = conversions.Cosmo_Conversions()

def power_beam_int(self, stokes='I'):
def power_beam_int(self, pol='I'):
"""
Computes the integral of the beam over solid angle to give
a beam area (in str) as a function of frequency. Uses function
in pyuvdata.
Currently, only the Stokes "I", "XX" and "YY" beam are supported.
See Equations 4 and 5 of Moore et al. (2017) ApJ 836, 154
or arxiv:1502.05072 for details.
Parameters
----------
stokes: str, optional
Which Stokes parameter's to compute the beam scalar for.
'I', 'Q', 'U', 'V', 'XX', 'YY', 'XY', 'YX', although
currently only 'I' and linear dipole pol (e.g. 'XX') are implemented.
pol: str, optional
Which polarization to compute the beam scalar for.
'I', 'Q', 'U', 'V', 'XX', 'YY', 'XY', 'YX'
Default: 'I'
Returns
-------
primary beam area: float, array-like
"""
if hasattr(self.primary_beam, 'get_beam_area'):
return self.primary_beam.get_beam_area(stokes)
return self.primary_beam.get_beam_area(pol)
else:
raise NotImplementedError("Outdated version of pyuvdata.")

def power_beam_sq_int(self, stokes='I'):
def power_beam_sq_int(self, pol='I'):
"""
Computes the integral of the beam**2 over solid angle to give
a beam**2 area (in str) as a function of frequency. Uses function
in pyuvdata.
Currently, only the Stokes "I", "XX" and "YY" beam are supported.
See Equations 4 and 5 of Moore et al. (2017) ApJ 836, 154
or arxiv:1502.05072 for details.
Parameters
----------
stokes: str, optional
Which Stokes parameter's to compute the beam scalar for.
'I', 'Q', 'U', 'V', 'XX', 'YY', 'XY', 'YX', although
currently only 'I' and linear dipole pol (e.g. 'XX') are implemented.
pol: str, optional
Which polarization to compute the beam scalar for.
'I', 'Q', 'U', 'V', 'XX', 'YY', 'XY', 'YX'
Default: 'I'
Returns
-------
primary beam area: float, array-like
"""
if hasattr(self.primary_beam, 'get_beam_area'):
return self.primary_beam.get_beam_sq_area(stokes)
return self.primary_beam.get_beam_sq_area(pol)
else:
raise NotImplementedError("Outdated version of pyuvdata.")
14 changes: 6 additions & 8 deletions hera_pspec/pspecdata.py
View file
Original file line number Diff line number Diff line change
Expand Up @@ -759,20 +759,18 @@ def delays(self):
return utils.get_delays(self.freqs[self.spw_range[0]:self.spw_range[1]]) * 1e9 # convert to ns


def scalar(self, stokes='I', taper='none', little_h=True, num_steps=2000, beam=None):
def scalar(self, pol='I', taper='none', little_h=True, num_steps=2000, beam=None):
"""
Computes the scalar function to convert a power spectrum estimate
in "telescope units" to cosmological units, using self.spw_range to set spectral window.
See arxiv:1304.4991 and HERA memo #27 for details.
Currently this is only for Stokes I.
Parameters
----------
stokes: str, optional
Which Stokes parameter's beam to compute the scalar for.
'I', 'Q', 'U', 'V', although currently only 'I' is implemented
pol: str, optional
Which polarization to compute the scalar for.
e.g. 'I', 'Q', 'U', 'V', 'XX', 'YY'...
Default: 'I'
taper : str, optional
Expand Down Expand Up @@ -806,12 +804,12 @@ def scalar(self, stokes='I', taper='none', little_h=True, num_steps=2000, beam=N
# calculate scalar
if beam is None:
scalar = self.primary_beam.compute_pspec_scalar(
start, end, len(freqs), stokes=stokes,
start, end, len(freqs), pol=pol,
taper=taper, little_h=little_h,
num_steps=num_steps)
else:
scalar = beam.compute_pspec_scalar(start, end, len(freqs),
stokes=stokes, taper=taper,
pol=pol, taper=taper,
little_h=little_h,
num_steps=num_steps)
return scalar
Expand Down
2 changes: 1 addition & 1 deletion hera_pspec/tests/test_noise.py
View file
Original file line number Diff line number Diff line change
Expand Up @@ -43,7 +43,7 @@ def test_scalar(self):
freqs = np.linspace(150e6, 160e6, 100, endpoint=False)
self.sense.calc_scalar(freqs, 'I', num_steps=5000, little_h=True)
nt.assert_true(np.isclose(freqs, self.sense.subband).all())
nt.assert_true(self.sense.stokes, 'I')
nt.assert_true(self.sense.pol, 'I')

def test_calc_P_N(self):
# calculate scalar
Expand Down
16 changes: 8 additions & 8 deletions hera_pspec/tests/test_pspecbeam.py
View file
Original file line number Diff line number Diff line change
Expand Up @@ -44,7 +44,7 @@ def test_UVbeam(self):
lower_freq = 120.*10**6
upper_freq = 128.*10**6
num_freqs = 20
scalar = self.bm.compute_pspec_scalar(lower_freq, upper_freq, num_freqs, stokes='I', num_steps=2000)
scalar = self.bm.compute_pspec_scalar(lower_freq, upper_freq, num_freqs, pol='I', num_steps=2000)

# Check that user-defined cosmology can be specified
bm = pspecbeam.PSpecBeamUV(self.beamfile,
Expand All @@ -55,11 +55,11 @@ def test_UVbeam(self):
self.assertEqual(Om_pp.ndim, 1)

# Check that errors are raised for other Stokes parameters
for stokes in ['Q', 'U', 'V', 'Z']:
nt.assert_raises(NotImplementedError, self.bm.power_beam_int, stokes=stokes)
nt.assert_raises(NotImplementedError, self.bm.power_beam_sq_int, stokes=stokes)
for pol in ['Q', 'U', 'V', 'Z']:
nt.assert_raises(NotImplementedError, self.bm.power_beam_int, pol=pol)
nt.assert_raises(NotImplementedError, self.bm.power_beam_sq_int, pol=pol)
nt.assert_raises(NotImplementedError, self.bm.compute_pspec_scalar,
lower_freq, upper_freq, num_freqs, stokes=stokes)
lower_freq, upper_freq, num_freqs, pol=pol)

self.assertAlmostEqual(Om_p[0], 0.078694909518866998)
self.assertAlmostEqual(Om_p[18], 0.065472512282419112)
Expand All @@ -72,7 +72,7 @@ def test_UVbeam(self):
self.assertAlmostEqual(scalar/567871703.75268996, 1.0, delta=1e-4)

# convergence of integral
scalar_large_Nsteps = self.bm.compute_pspec_scalar(lower_freq, upper_freq, num_freqs, stokes='I', num_steps=10000)
scalar_large_Nsteps = self.bm.compute_pspec_scalar(lower_freq, upper_freq, num_freqs, pol='I', num_steps=10000)
self.assertAlmostEqual(scalar / scalar_large_Nsteps, 1.0, delta=1e-5)

# test taper execution
Expand All @@ -93,7 +93,7 @@ def test_UVbeam(self):
nt.assert_raises(TypeError, self.bm.Jy_to_mK, np.array([1]))

# test noise scalar
sclr = self.bm.compute_pspec_scalar(lower_freq, upper_freq, num_freqs, stokes='I', num_steps=2000, noise_scalar=True)
sclr = self.bm.compute_pspec_scalar(lower_freq, upper_freq, num_freqs, pol='I', num_steps=2000, noise_scalar=True)
nt.assert_almost_equal(sclr, 70.983962969086235)


Expand All @@ -103,7 +103,7 @@ def test_Gaussbeam(self):
lower_freq = 120.*10**6
upper_freq = 128.*10**6
num_freqs = 20
scalar = self.gauss.compute_pspec_scalar(lower_freq, upper_freq, num_freqs, stokes='I', num_steps=2000)
scalar = self.gauss.compute_pspec_scalar(lower_freq, upper_freq, num_freqs, pol='I', num_steps=2000)

# Check that user-defined cosmology can be specified
bgauss = pspecbeam.PSpecBeamGauss(0.8,
Expand Down

0 comments on commit c069177

Please sign in to comment.