Skip to content

Commit

Permalink
Fix issue with testing data being generated with wrong wgts_array dim…
Browse files Browse the repository at this point in the history 
…ensions; minor whitespace changes
  • Loading branch information
@philbull
philbull committed Sep 12, 2018
1 parent f0d3401 commit ee35efc
Show file tree
Hide file tree
Showing 2 changed files with 58 additions and 34 deletions.
50 changes: 30 additions & 20 deletions hera_pspec/testing.py
View file
Original file line number Diff line number Diff line change
Expand Up @@ -9,7 +9,8 @@

def build_vanilla_uvpspec(beam=None):
"""
Build an example vanilla UVPSpec object from scratch, with all necessary metadata.
Build an example vanilla UVPSpec object from scratch, with all necessary
metadata.
Parameters
----------
Expand Down Expand Up @@ -54,7 +55,8 @@ def build_vanilla_uvpspec(beam=None):
spw_freq_array = np.tile(np.arange(Nspws), Nfreqs)
spw_dly_array = np.tile(np.arange(Nspws), Ndlys)
spw_array = np.arange(Nspws)
freq_array = np.repeat(np.linspace(100e6, 105e6, Nfreqs, endpoint=False), Nspws)
freq_array = np.repeat(np.linspace(100e6, 105e6, Nfreqs, endpoint=False),
Nspws)
dly_array = np.repeat(utils.get_delays(freq_array, n_dlys=Ndlys), Nspws)
pol_array = np.array([-5])
Npols = len(pol_array)
Expand Down Expand Up @@ -84,26 +86,31 @@ def build_vanilla_uvpspec(beam=None):
store_cov=True
cosmo = conversions.Cosmo_Conversions()

data_array, wgt_array, integration_array, nsample_array, cov_array = {}, {}, {}, {}, {}
data_array, wgt_array = {}, {}
integration_array, nsample_array, cov_array = {}, {}, {}
for s in spw_array:
data_array[s] = np.ones((Nblpairts, Ndlys, Npols), dtype=np.complex) \
* blpair_array[:, None, None] / 1e9
wgt_array[s] = np.ones((Nblpairts, Ndlys, 2, Npols), dtype=np.float)
wgt_array[s] = np.ones((Nblpairts, Nfreqs, 2, Npols), dtype=np.float)
# NB: The wgt_array has dimensions Nfreqs rather than Ndlys; it has the
# dimensions of the input visibilities, not the output delay spectra
integration_array[s] = np.ones((Nblpairts, Npols), dtype=np.float)
nsample_array[s] = np.ones((Nblpairts, Npols), dtype=np.float)
cov_array[s] =np.moveaxis(np.array([[np.identity(Ndlys,dtype=np.complex)\
for m in range(Nblpairts)] for n in range(Npols)]),0,-1)
for m in range(Nblpairts)]
for n in range(Npols)]), 0, -1)


params = ['Ntimes', 'Nfreqs', 'Nspws', 'Nspwdlys', 'Nspwfreqs', 'Nspws', 'Nblpairs', 'Nblpairts',
'Npols', 'Ndlys', 'Nbls', 'blpair_array', 'time_1_array',
'time_2_array', 'lst_1_array', 'lst_2_array', 'spw_array',
params = ['Ntimes', 'Nfreqs', 'Nspws', 'Nspwdlys', 'Nspwfreqs', 'Nspws',
'Nblpairs', 'Nblpairts', 'Npols', 'Ndlys', 'Nbls',
'blpair_array', 'time_1_array', 'time_2_array',
'lst_1_array', 'lst_2_array', 'spw_array',
'dly_array', 'freq_array', 'pol_array', 'data_array', 'wgt_array',
'integration_array', 'bl_array', 'bl_vecs', 'telescope_location',
'vis_units', 'channel_width', 'weighting', 'history', 'taper', 'norm',
'git_hash', 'nsample_array', 'time_avg_array', 'lst_avg_array',
'cosmo', 'scalar_array', 'labels', 'norm_units', 'labels', 'label_1_array',
'label_2_array','store_cov','cov_array', 'spw_dly_array', 'spw_freq_array']
'vis_units', 'channel_width', 'weighting', 'history', 'taper',
'norm', 'git_hash', 'nsample_array', 'time_avg_array',
'lst_avg_array', 'cosmo', 'scalar_array', 'labels', 'norm_units',
'labels', 'label_1_array', 'label_2_array', 'store_cov',
'cov_array', 'spw_dly_array', 'spw_freq_array']

if beam is not None:
params += ['OmegaP', 'OmegaPP', 'beam_freqs']
Expand Down Expand Up @@ -192,21 +199,24 @@ def uvpspec_from_data(data, bl_grps, data_std=None, spw_ranges=None,

# instantiate pspecdata
ds = pspecdata.PSpecData(dsets=[uvd, uvd], dsets_std=[uvd_std, uvd_std],
wgts=[None, None], labels=['d1', 'd2'], beam=beam,
n_dlys=n_dlys)
wgts=[None, None], labels=['d1', 'd2'], beam=beam)

# get blpair groups
assert isinstance(bl_grps, list), "bl_grps must be a list"
if not isinstance(bl_grps[0], list): bl_grps = [bl_grps]
assert np.all([isinstance(blgrp, list) for blgrp in bl_grps]), "bl_grps must be fed as a list of lists"
assert np.all([isinstance(blgrp[0], tuple) for blgrp in bl_grps]), "bl_grps must be fed as a list of lists of tuples"
assert np.all([isinstance(blgrp, list) for blgrp in bl_grps]), \
"bl_grps must be fed as a list of lists"
assert np.all([isinstance(blgrp[0], tuple) for blgrp in bl_grps]), \
"bl_grps must be fed as a list of lists of tuples"
bls1, bls2 = [], []
for blgrp in bl_grps:
_bls1, _bls2, _ = utils.construct_blpairs(blgrp, exclude_auto_bls=True, exclude_permutations=True)
_bls1, _bls2, _ = utils.construct_blpairs(blgrp, exclude_auto_bls=True,
exclude_permutations=True)
bls1.extend(_bls1)
bls2.extend(_bls2)

# run pspec
uvp = ds.pspec(bls1, bls2, (0, 1), (pol, pol), input_data_weight='identity', spw_ranges=spw_ranges,
taper=taper, verbose=verbose,store_cov=store_cov)
uvp = ds.pspec(bls1, bls2, (0, 1), (pol, pol), input_data_weight='identity',
spw_ranges=spw_ranges, taper=taper, verbose=verbose,
store_cov=store_cov, n_dlys=n_dlys)
return uvp
42 changes: 28 additions & 14 deletions hera_pspec/tests/test_uvpspec.py
View file
Original file line number Diff line number Diff line change
Expand Up @@ -45,7 +45,7 @@ def test_get_funcs(self):
nt.assert_almost_equal(d[0,0], (101.1021011020000001+0j))
# get_wgts
w = self.uvp.get_wgts((0, ((1, 2), (1, 2)), 'xx'))
nt.assert_equal(w.shape, (10, 30, 2))
nt.assert_equal(w.shape, (10, 50, 2)) # should have Nfreq dim, not Ndlys
nt.assert_true(w.dtype == np.float)
nt.assert_equal(w[0,0,0], 1.0)
# get_integrations
Expand Down Expand Up @@ -506,8 +506,9 @@ def test_combine_uvpspec(self):
nt.assert_equal(out.Npols, 3)

# Test whether n_dlys != Nfreqs works
uvp4 = testing.uvpspec_from_data(uvd, bls, spw_ranges=[(20, 30), (60, 90)],
beam=beam, n_dlys=uvd.Nfreqs//2)
uvp4 = testing.uvpspec_from_data(uvd, bls, beam=beam,
spw_ranges=[(20, 30), (60, 90)],
n_dlys=[5, 15])
uvp4b = copy.deepcopy(uvp4)
uvp4b.pol_array[0] = -6
out = uvpspec.combine_uvpspec([uvp4, uvp4b], verbose=False)
Expand All @@ -518,10 +519,14 @@ def test_combine_uvpspec_std(self):
uvd = UVData()
uvd_std = UVData()
uvd.read_miriad(os.path.join(DATA_PATH, 'zen.even.xx.LST.1.28828.uvOCRSA'))
uvd_std.read_miriad(os.path.join(DATA_PATH,'zen.even.xx.LST.1.28828.uvOCRSA'))
beam = pspecbeam.PSpecBeamUV(os.path.join(DATA_PATH, "HERA_NF_dipole_power.beamfits"))
uvd_std.read_miriad(
os.path.join(DATA_PATH,'zen.even.xx.LST.1.28828.uvOCRSA'))
beam = pspecbeam.PSpecBeamUV(
os.path.join(DATA_PATH, "HERA_NF_dipole_power.beamfits"))
bls = [(37, 38), (38, 39), (52, 53)]
uvp1 = testing.uvpspec_from_data(uvd, bls, data_std=uvd_std, spw_ranges=[(20, 24), (64, 68)], beam=beam)
uvp1 = testing.uvpspec_from_data(uvd, bls, data_std=uvd_std,
spw_ranges=[(20, 24), (64, 68)],
beam=beam)
# test failure due to overlapping data
uvp2 = copy.deepcopy(uvp1)
nt.assert_raises(AssertionError, uvpspec.combine_uvpspec, [uvp1, uvp2])
Expand All @@ -535,24 +540,32 @@ def test_combine_uvpspec_std(self):
nt.assert_raises(AssertionError, uvpspec.combine_uvpspec, [uvp1, uvp2])

# test partial data overlap failure
uvp2 = testing.uvpspec_from_data(uvd, [(37, 38), (38, 39), (53, 54)], data_std=uvd_std, spw_ranges=[(20, 24), (64, 68)], beam=beam)
uvp2 = testing.uvpspec_from_data(uvd, [(37, 38), (38, 39), (53, 54)],
data_std=uvd_std,
spw_ranges=[(20, 24), (64, 68)],
beam=beam)
nt.assert_raises(AssertionError, uvpspec.combine_uvpspec, [uvp1, uvp2])
uvp2 = testing.uvpspec_from_data(uvd, bls, spw_ranges=[(20, 24), (64, 68)], data_std=uvd_std, beam=beam)
uvp2 = testing.uvpspec_from_data(uvd, bls,
spw_ranges=[(20, 24), (64, 68)],
data_std=uvd_std, beam=beam)
nt.assert_raises(AssertionError, uvpspec.combine_uvpspec, [uvp1, uvp2])
uvp2 = copy.deepcopy(uvp1)
uvp2.pol_array[0] = -6
uvp2 = uvpspec.combine_uvpspec([uvp1, uvp2], verbose=False)
nt.assert_raises(AssertionError, uvpspec.combine_uvpspec, [uvp1, uvp2])

# test concat across spw
uvp2 = testing.uvpspec_from_data(uvd, bls, spw_ranges=[(85, 91)], data_std=uvd_std, beam=beam)
uvp2 = testing.uvpspec_from_data(uvd, bls, spw_ranges=[(85, 91)],
data_std=uvd_std, beam=beam)
out = uvpspec.combine_uvpspec([uvp1, uvp2], verbose=False)
nt.assert_equal(out.Nspws, 3)
nt.assert_equal(out.Nfreqs, 14)
nt.assert_equal(out.Nspwdlys, 14)

# test concat across blpairts
uvp2 = testing.uvpspec_from_data(uvd, [(53, 54), (67, 68)], spw_ranges=[(20, 24), (64, 68)], data_std=uvd_std, beam=beam)
uvp2 = testing.uvpspec_from_data(uvd, [(53, 54), (67, 68)],
spw_ranges=[(20, 24), (64, 68)],
data_std=uvd_std, beam=beam)
out = uvpspec.combine_uvpspec([uvp1, uvp2], verbose=False)
nt.assert_equal(out.Nblpairs, 4)
nt.assert_equal(out.Nbls, 5)
Expand All @@ -566,17 +579,17 @@ def test_combine_uvpspec_std(self):
nt.assert_raises(AssertionError, uvpspec.combine_uvpspec, [uvp1, uvp2])

# test feed as strings
uvp1 = testing.uvpspec_from_data(uvd, bls, spw_ranges=[(20, 30)], data_std=uvd_std, beam=beam)
uvp1 = testing.uvpspec_from_data(uvd, bls, spw_ranges=[(20, 30)],
data_std=uvd_std, beam=beam)
uvp2 = copy.deepcopy(uvp1)
uvp2.pol_array[0] = -6
uvp1.write_hdf5('uvp1.hdf5', overwrite=True)
uvp2.write_hdf5('uvp2.hdf5', overwrite=True)
out = uvpspec.combine_uvpspec(['uvp1.hdf5', 'uvp2.hdf5'], verbose=False)
nt.assert_true(out.Npols, 2)
for ff in ['uvp1.hdf5', 'uvp2.hdf5']:
if os.path.exists(ff):
os.remove(ff)

if os.path.exists(ff): os.remove(ff)

# test UVPSpec __add__
uvp3 = copy.deepcopy(uvp1)
uvp3.pol_array[0] = -7
Expand All @@ -599,3 +612,4 @@ def test_conj_blpair():
blpair = uvputils._conj_blpair(101102103104, which='both')
nt.assert_equal(blpair, 102101104103)
nt.assert_raises(ValueError, uvputils._conj_blpair, 102101103104, which='foo')

0 comments on commit ee35efc

Please sign in to comment.