Fixed merge conflict in utils

hera_pspec/__init__.py

@@ -1,8 +1,7 @@
 """
 __init__.py file for hera_pspec
 """
-
-from hera_pspec import version, conversions, grouping, pspecbeam, plot, pstokes
+from hera_pspec import version, conversions, grouping, pspecbeam, plot, pstokes, testing
 from hera_pspec import uvpspec_utils as uvputils
 
 from hera_pspec.uvpspec import UVPSpec

hera_pspec/data/zen.even.xx.LST.1.28828.uvOCRSA/history

@@ -2,3 +2,4 @@
   Read/written with pyuvdata version: 1.2.1.  Git origin: https://github.com/HERA-Team/pyuvdata.git.  Git hash: c6c034802549acd0bc8381233b711a84ba08a3ac.  Git branch: master.  Git description: v1.2-186-gc6c0348. Combined data along baseline-time axis using pyuvdata. 20 zen.2458099.42378.xx.HH.uvOCRS-zen.2458099.43124.xx.HH.uvOCRS-zen.2458102.40886.xx.HH.uvOCRS-zen.2458102.41632.xx.HH.uvOCRS-zen.2458102.42378.xx.HH.uvOCRS-zen.2458104.40886.xx.HH.uvOCRS-zen.2458104.41631.xx.HH.uvOCRS-zen.2458104.42377.xx.HH.uvOCRS-zen.2458106.40141.xx.HH.uvOCRS-zen.2458106.40886.xx.HH.uvOCRS-zen.2458106.41632.xx.HH.uvOCRS-zen.2458108.39395.xx.HH.uvOCRS-zen.2458108.40140.xx.HH.uvOCRS-zen.2458108.40886.xx.HH.uvOCRS-zen.2458110.38650.xx.HH.uvOCRS-zen.2458110.39395.xx.HH.uvOCRS-zen.2458110.40141.xx.HH.uvOCRS-zen.2458112.38649.xx.HH.uvOCRS-zen.2458112.39395.xx.HH.uvOCRS-zen.2458112.40141.xx.HH.uvOCRS-zen.2458114.37559.xx.HH.uvOCRS-zen.2458114.38305.xx.HH.uvOCRS-zen.2458114.39051.xx.HH.uvOCRS-zen.2458116.37159.xx.HH.uvOCRS-zen.2458116.37904.xx.HH.uvOCRS-zen.2458116.38650.xx.HH.uvOCRS
   Downselected to specific times using pyuvdata.
   Read/written with pyuvdata version: 1.2.1.  Git origin: https://github.com/HERA-Team/pyuvdata.  Git hash: 7fd7a52cb9bdf1f80dc49d50eea821990c9803b0.  Git branch: partial_io_read_miriad.  Git description: v1.2-205-g7fd7a52.  Downselected to specific antennas, times, frequencies using pyuvdata.
+  Read/written with pyuvdata version: 1.2.1.  Git origin: https://github.com/HERA-Team/pyuvdata.  Git hash: 445e76acdcdbfd37d1468d4d4ed6dddd1725614c.  Git branch: utils_get_miriad_antpos.  Git description: v1.2-259-g445e76a.

hera_pspec/data/zen.odd.xx.LST.1.28828.uvOCRSA/history

@@ -1,2 +1,3 @@
  /users/heramgr/anaconda2/envs/hera/bin/lstbin_run.py --dlst 0.001565970008 --lst_start 5.0 --ntimes_per_file 60 --file_ext odd.{}.{}.{:7.5f}.uvOCRS --outdir ./even_odd --overwrite --output_file_select 20 --sig_clip --sigma 4.0 --min_N 5 --rephase /lustre/aoc/projects/hera/H1C_IDR2/IDR2_1/2458098/zen.*.*.xx.HH.uvOCRS /lustre/aoc/projects/hera/H1C_IDR2/IDR2_1/2458101/zen.*.*.xx.HH.uvOCRS /lustre/aoc/projects/hera/H1C_IDR2/IDR2_1/2458103/zen.*.*.xx.HH.uvOCRS /lustre/aoc/projects/hera/H1C_IDR2/IDR2_1/2458105/zen.*.*.xx.HH.uvOCRS /lustre/aoc/projects/hera/H1C_IDR2/IDR2_1/2458107/zen.*.*.xx.HH.uvOCRS /lustre/aoc/projects/hera/H1C_IDR2/IDR2_1/2458109/zen.*.*.xx.HH.uvOCRS /lustre/aoc/projects/hera/H1C_IDR2/IDR2_1/2458111/zen.*.*.xx.HH.uvOCRS /lustre/aoc/projects/hera/H1C_IDR2/IDR2_1/2458113/zen.*.*.xx.HH.uvOCRS /lustre/aoc/projects/hera/H1C_IDR2/IDR2_1/2458115/zen.*.*.xx.HH.uvOCRS Input files: zen.2458098.42378.xx.HH.uvOCRS-zen.2458098.43124.xx.HH.uvOCRS-zen.2458101.40887.xx.HH.uvOCRS-zen.2458101.41633.xx.HH.uvOCRS-zen.2458101.42378.xx.HH.uvOCRS-zen.2458103.40886.xx.HH.uvOCRS-zen.2458103.41632.xx.HH.uvOCRS-zen.2458103.42377.xx.HH.uvOCRS-zen.2458105.40140.xx.HH.uvOCRS-zen.2458105.40886.xx.HH.uvOCRS-zen.2458105.41632.xx.HH.uvOCRS-zen.2458107.39395.xx.HH.uvOCRS-zen.2458107.40141.xx.HH.uvOCRS-zen.2458107.40887.xx.HH.uvOCRS-zen.2458109.38649.xx.HH.uvOCRS-zen.2458109.39395.xx.HH.uvOCRS-zen.2458109.40141.xx.HH.uvOCRS-zen.2458111.38650.xx.HH.uvOCRS-zen.2458111.39396.xx.HH.uvOCRS-zen.2458111.40141.xx.HH.uvOCRS-zen.2458113.37903.xx.HH.uvOCRS-zen.2458113.38649.xx.HH.uvOCRS-zen.2458113.39395.xx.HH.uvOCRS-zen.2458115.37158.xx.HH.uvOCRS-zen.2458115.37904.xx.HH.uvOCRS-zen.2458115.38649.xx.HH.uvOCRS
   Read/written with pyuvdata version: 1.2.1.  Git origin: https://github.com/HERA-Team/pyuvdata.  Git hash: 7fd7a52cb9bdf1f80dc49d50eea821990c9803b0.  Git branch: partial_io_read_miriad.  Git description: v1.2-205-g7fd7a52.  Downselected to specific antennas, times, frequencies using pyuvdata.
+  Read/written with pyuvdata version: 1.2.1.  Git origin: https://github.com/HERA-Team/pyuvdata.  Git hash: 445e76acdcdbfd37d1468d4d4ed6dddd1725614c.  Git branch: utils_get_miriad_antpos.  Git description: v1.2-259-g445e76a.

hera_pspec/pstokes.py

@@ -269,3 +269,48 @@ def construct_pstokes(dset1, dset2, pstokes='pI', run_check=True, antenna_nums=N
         uvdS.check()
 
     return uvdS
+
+
+def filter_dset_on_stokes_pol(dsets, pstokes):
+    """
+    Given a list of UVData objects with dipole linear polarizations,
+    and a desired output pstokes, return the two UVData objects from
+    the input dsets that can be used in construct_pstokes to make 
+    the desired pseudo-Stokes visibility. If a single UVData object 
+    has multiple polarizations, this function only considers its first.
+
+    Parameters
+    ----------
+    dsets : list
+        List of UVData objects with linear dipole polarizations
+
+    pstokes : str or int
+        Pseudo-stokes polarization one wants to form out of input dsets.
+        Ex. 'pI', or 'pU', or 1, ...
+
+    Returns
+    -------
+    inp_dsets : list
+        List of two UVData objects from the input dsets that can be fed
+        to construct_pstokes to make the desired pseudo-Stokes visibility.
+    """
+    # type check
+    assert isinstance(dsets, list), "dsets must be fed as a list of UVData objects"
+    assert np.all(isinstance(d, UVData) for d in dsets), "dsets must be fed as a list of UVData objects"
+
+    # get polarization of each dset
+    pols = [d.polarization_array[0] for d in dsets]
+
+    # convert pstokes to integer if a string
+    if isinstance(pstokes, (str, np.str)):
+        pstokes = pyuvdata.utils.polstr2num(pstokes)
+    assert pstokes in [1, 2, 3, 4], "pstokes must be fed as a pseudo-Stokes parameter"
+
+    # get two necessary dipole pols given pstokes
+    desired_pols = pol_weights[pstokes].keys()
+    assert desired_pols[0] in pols and desired_pols[1] in pols, "necessary input pols {} and {} not found in dsets".format(*desired_pols)
+
+    inp_dsets = [dsets[pols.index(desired_pols[0])], dsets[pols.index(desired_pols[1])]]
+
+    return inp_dsets
+

hera_pspec/testing.py

@@ -0,0 +1,168 @@
+#!/usr/bin/env python2
+import numpy as np
+import copy, operator, itertools
+from collections import OrderedDict as odict
+from hera_pspec import uvpspec, pspecdata, conversions, pspecbeam, utils
+from pyuvdata import UVData
+
+
+def build_vanilla_uvpspec(beam=None):
+    """
+    Build an example vanilla UVPSpec object from scratch, with all necessary metadata.
+
+    Parameters
+    ----------
+    beam : PSpecBeamBase subclass
+
+    Returns
+    -------
+    uvp : UVPSpec object
+    """
+    uvp = uvpspec.UVPSpec()
+
+    Ntimes = 10
+    Nfreqs = 50
+    Ndlys = Nfreqs
+    Nspws = 1
+    Nspwdlys = Nspws * Nfreqs
+
+    # [((1, 2), (1, 2)), ((2, 3), (2, 3)), ((1, 3), (1, 3))]
+    blpairs = [1002001002, 2003002003, 1003001003]
+    bls = [1002, 2003, 1003]
+    Nbls = len(bls)
+    Nblpairs = len(blpairs)
+    Nblpairts = Nblpairs * Ntimes
+
+    blpair_array = np.tile(blpairs, Ntimes)
+    bl_array = np.array(bls)
+    bl_vecs = np.array([[  5.33391548e+00,  -1.35907816e+01,  -7.91624188e-09],
+                        [ -8.67982998e+00,   4.43554478e+00,  -1.08695203e+01],
+                        [ -3.34591450e+00,  -9.15523687e+00,  -1.08695203e+01]])
+    time_array = np.repeat(np.linspace(2458042.1, 2458042.2, Ntimes), Nblpairs)
+    time_1_array = time_array
+    time_2_array = time_array
+    lst_array = np.repeat(np.ones(Ntimes, dtype=np.float), Nblpairs)
+    lst_1_array = lst_array
+    lst_2_array = lst_array
+    time_avg_array = time_array
+    lst_avg_array = lst_array
+    spws = np.arange(Nspws)
+    spw_array = np.tile(spws, Ndlys)
+    freq_array = np.repeat(np.linspace(100e6, 105e6, Nfreqs, endpoint=False), Nspws)
+    dly_array = np.fft.fftshift(np.repeat(np.fft.fftfreq(Nfreqs, np.median(np.diff(freq_array))), Nspws))
+    pol_array = np.array([-5])
+    Npols = len(pol_array)
+    vis_units = 'unknown'
+    norm_units = 'Hz str'
+    weighting = 'identity'
+    channel_width = np.median(np.diff(freq_array))
+    history = 'example'
+    taper = "none"
+    norm = "I"
+    git_hash = "random"
+    scalar_array = np.ones((Nspws, Npols), np.float)
+    label1 = 'red'
+    label2 = 'blue'
+    labels = np.array([label1, label2])
+    label_1_array = np.ones((Nspws, Nblpairts, Npols), np.int) * 0
+    label_2_array = np.ones((Nspws, Nblpairts, Npols), np.int) * 1
+    if beam is not None:
+        OmegaP, OmegaPP = beam.get_Omegas(beam.primary_beam.polarization_array[0])
+        beam_freqs = beam.beam_freqs
+
+    # HERA coordinates in Karoo Desert, SA
+    telescope_location = np.array([5109325.85521063, 
+                                   2005235.09142983, 
+                                  -3239928.42475397])
+
+    cosmo = conversions.Cosmo_Conversions()
+
+    data_array, wgt_array, integration_array, nsample_array = {}, {}, {}, {}
+    for s in spws:
+        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)
+        integration_array[s] = np.ones((Nblpairts, Npols), dtype=np.float)
+        nsample_array[s] = np.ones((Nblpairts, Npols), dtype=np.float)
+
+    params = ['Ntimes', 'Nfreqs', 'Nspws', 'Nspwdlys', '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']
+
+    if beam is not None:
+        params += ['OmegaP', 'OmegaPP', 'beam_freqs']
+
+    # Set all parameters
+    for p in params:
+        setattr(uvp, p, locals()[p])
+
+    uvp.check()
+
+    return uvp, cosmo
+
+def uvpspec_from_data(data, bls, spw_ranges=None, beam=None, taper='none', cosmo=None, verbose=False):
+    """
+    Build an example UVPSpec object from a visibility file and PSpecData.
+
+    Parameters
+    ----------
+    data : UVData object or str
+        This can be a UVData object or a string filepath to a miriad file.
+
+    bls : list
+        This is a list of at least 2 baseline tuples.
+        Ex: [(24, 25), (37, 38), ...]
+
+    spw_ranges : list
+        List of spectral window tuples. See PSpecData.pspec docstring for details.
+
+    beam : PSpecBeamBase subclass or str
+        This can be a subclass of PSpecBeamBase of a string filepath to a 
+        UVBeam healpix map.
+
+    taper : string
+        Optional tapering applied to the data before OQE.
+
+    cosmo : Cosmo_Conversions object
+
+    verbose : bool
+        if True, report feedback to standard output
+
+    Returns
+    -------
+    uvp : UVPSpec object
+    """
+    # load data
+    if isinstance(data, str):
+        uvd = UVData()
+        uvd.read_miriad(data)
+    elif isinstance(data, UVData):
+        uvd = data
+
+    # get pol
+    pol = uvd.polarization_array[0]
+
+    # load beam
+    if isinstance(beam, str):
+        beam = pspecbeam.PSpecBeamUV(beam, cosmo=cosmo)
+    if beam is not None and cosmo is not None:
+        beam.cosmo = cosmo
+
+    # instantiate pspecdata
+    ds = pspecdata.PSpecData(dsets=[uvd, uvd], wgts=[None, None], labels=['d1', 'd2'], beam=beam)
+
+    # get red bls
+    bls1, bls2, _ = utils.construct_blpairs(bls, exclude_auto_bls=True)
+
+    # run pspec
+    uvp = ds.pspec(bls1, bls2, (0, 1), (pol, pol), input_data_weight='identity', spw_ranges=spw_ranges, 
+                   taper=taper, verbose=verbose)
+
+    return uvp
+

hera_pspec/tests/test_container.py

@@ -3,14 +3,14 @@
 import numpy as np
 import os, sys, copy
 from hera_pspec.data import DATA_PATH
-from hera_pspec import PSpecContainer, UVPSpec
-from test_uvpspec import build_example_uvpspec
+from hera_pspec import PSpecContainer, UVPSpec, testing
+
 
 class Test_PSpecContainer(unittest.TestCase):
 
     def setUp(self):
         self.fname = os.path.join(DATA_PATH, '_test_container.hdf5')
-        self.uvp, self.cosmo = build_example_uvpspec()
+        self.uvp, self.cosmo = testing.build_vanilla_uvpspec()
 
     def tearDown(self):
         # Remove HDF5 file

hera_pspec/tests/test_grouping.py

@@ -3,8 +3,7 @@
 import numpy as np
 import os
 from hera_pspec.data import DATA_PATH
-from test_uvpspec import build_example_uvpspec
-from hera_pspec import uvpspec, conversions, parameter, pspecbeam, pspecdata
+from hera_pspec import uvpspec, conversions, parameter, pspecbeam, pspecdata, testing
 from hera_pspec import uvpspec_utils as uvputils
 from hera_pspec import grouping
 
@@ -13,7 +12,7 @@ class Test_grouping(unittest.TestCase):
     def setUp(self):
         beamfile = os.path.join(DATA_PATH, 'NF_HERA_Beams.beamfits')
         self.beam = pspecbeam.PSpecBeamUV(beamfile)
-        uvp, cosmo = build_example_uvpspec(beam=self.beam)
+        uvp, cosmo = testing.build_vanilla_uvpspec(beam=self.beam)
         uvp.check()
         self.uvp = uvp
 

hera_pspec/tests/test_plot.py

@@ -3,7 +3,7 @@
 import matplotlib
 import numpy as np
 import os, copy, sys
-from hera_pspec import pspecdata, pspecbeam, conversions, plot
+from hera_pspec import pspecdata, pspecbeam, conversions, plot, utils
 from hera_pspec.data import DATA_PATH
 from pyuvdata import UVData
 
@@ -42,7 +42,7 @@ def setUp(self):
 
         # Construct list of baseline pairs to calculate power spectra for
         bls = [(24,25), (37,38), (38,39),]
-        self.bls1, self.bls2, blp = pspecdata.construct_blpairs(
+        self.bls1, self.bls2, blp = utils.construct_blpairs(
                         bls, exclude_permutations=False, exclude_auto_bls=True)
 
         # Calculate the power spectrum