Repository: RadioAstronomySoftwareGroup/pyuvdata (python)

Issue:
get_antenna_redundancies fails to properly group North-South oriented baselines.
This example demonstrates the issue:

```python
import numpy as np
import pyuvdata
import matplotlib.pyplot as plt

eps = 1e-10 # the size of the error induced by pyuvdata ENU <-> ECEF transformations

# antenna positions at corners of a square with two small perturbations
ant_pos = np.array([
    [-eps,1.,0.],
    [1.,1.,0.],
    [eps,0.,0.],
    [1.,0.,0.]
])

ant_nums = np.arange(4)

red_grps, _, _ = pyuvdata.utils.get_antenna_redundancies(ant_nums, ant_pos, tol=0.1)

for grp in red_grps:
    bl_num = grp[0] # get one baseline from each redundant group
    ai, aj = pyuvdata.utils.baseline_to_antnums(bl_num, 4)
    bij = ant_pos[aj] - ant_pos[ai]
    
    plt.plot(bij[0], bij[1], '.k', ms=10)
plt.ylim(-1.1,1.1)
plt.xlim(-1.1,1.1)
plt.show()
```
produces this figure of baseline positions for this demo array

![redundancy_group_failure](https://user-images.githubusercontent.com/13594524/62745625-e7f43200-ba19-11e9-9787-c3a7420425cc.png)

`red_grps` has 5 groups when there are really only 4. The points at (0, 1) and (0, -1) should be in the same redundant group. This looks like a problem with the initial attempt in `get_antenna_redundancies` to select only half the points.

Pull request:
Redundancy branch cut fix
<!--- Provide a general summary of your changes in the Title above -->

## Description
<!--- Describe your changes in detail -->
Fixes a bug in `utils.get_antenna_redundancies` for North-South baselines. This method flipped baselines according to the u-positive convention. Under this convention, North-South baselines (u ~ 0) would be flipped only if u is exactly zero and v is negative. This means that a baseline with u = -ε and v < 0 would be flipped, but another with u = +ε w and v<0 would not. Redundancies among North-South baselines to the given tolerance are not identified.

This fixes that defining ` u = 0` as `np.isclose(u, 0, atol=tol)` (ie., the specified redundancy tolerance, and `u<0` as `u < -tol` (likewise for v and w). This means that very nearly north-south baselines will be flipped if v<0 (ie., pointing south), ensuring that North-south redundancies are identified.

The fix is tested using the example code from issue #650 .

## Motivation and Context
<!--- Why is this change required? What problem does it solve? -->
This fixes a significant problem for the redundancy finder when handling realistic baselines.
<!--- If it fixes an open issue, please link to the issue here. If this PR closes an issue, put the word 'closes' before the issue link to auto-close the issue when the PR is merged. -->
closes #650 

## Types of changes
<!--- What types of changes does your code introduce? Put an `x` in all the boxes that apply: -->
- [x] Bug fix (non-breaking change which fixes an issue)
- [ ] New feature (non-breaking change which adds functionality)
- [ ] Breaking change (fix or feature that would cause existing functionality to change)
- [ ] Documentation change (documentation changes only)
- [ ] Version change
- [ ] Build or continuous integration change


## Checklist:
<!--- You may remove the checklists that don't apply to your change type(s) or just leave them empty -->
<!--- Go over all the following points, and replace the space with an `x` in all the boxes that apply. -->
<!--- If you're unsure about any of these, don't hesitate to ask. We're here to help! -->
- [x] I have read the [contribution guide](https://github.com/RadioAstronomySoftwareGroup/pyuvdata/blob/master/.github/CONTRIBUTING.md).
- [x] My code follows the code style of this project.

Bug fix checklist:
- [x] My fix includes a new test that breaks as a result of the bug (if possible).
- [x] All new and existing tests pass in both python 2 and python 3.
- [x] I have updated the [CHANGELOG](https://github.com/RadioAstronomySoftwareGroup/pyuvdata/blob/master/CHANGELOG.md).

Patch:
diff --git a/CHANGELOG.md b/CHANGELOG.md
index 00215d68..53a71fa8 100644
--- a/CHANGELOG.md
+++ b/CHANGELOG.md
@@ -4,6 +4,9 @@ All notable changes to this project will be documented in this file.
 ## [Unreleased]
 
 ### Added
+- `utils.uvcalibrate(.., undo=True)` kwarg for undo-ing a calibration.
+- `utils.uvcalibrate` updates `UVData.vis_units` if `UVCal.gain_scale` is set.
+- `UVCal.gain_scale` non-required attribute
 - UVData.get_redundancies method to replace old get_baseline_redundancies and get_antenna_redundancies.
 - option for `UVBeam.interp` to return a new beam object.
 - `UVFlag` information on Read The Docs
@@ -14,6 +17,7 @@ All notable changes to this project will be documented in this file.
 - `UVFlag` inherits from `UVBase` object.
 
 ### Fixed
+- A bug in get_antenna_redundancies for nearly N-S baselines.
 - A bug where `conj_pol` could not handle cardinal direction polarizations.
 - A bug that gave the wrong error message when calling `UVData.phase_to_time` without an Astropy Time object.
 
diff --git a/docs/tutorial.rst b/docs/tutorial.rst
index 7fb060f0..2c77b8ec 100644
--- a/docs/tutorial.rst
+++ b/docs/tutorial.rst
@@ -959,22 +959,29 @@ considered redundant if they are within a specified tolerance distance (default
 The :func:`utils.get_baseline_redundancies` function accepts an array of baseline indices
 and an array of baseline vectors (ie, uvw coordinates) as input, and finds
 redundancies among the vectors as given. If the ``with_conjugates`` option is
-selected, it will flip baselines such that ``u > 0``, or ``v > 0 if u = 0``, or
-``w > 0 if u = v = 0``. In this case, a list of ``conjugates`` is returned as well,
+selected, it will include baselines that are redundant when reversed in the same group.
+In this case, a list of ``conjugates`` is returned as well,
 which contains indices for the baselines that were flipped for the redundant groups.
 In either mode of operation, this will only return baseline indices that are in the list passed in.
 
 The :func:`utils.get_antenna_redundancies` function accepts an array of
 antenna indices and an array of antenna positions as input, defines baseline vectors
-and indices under the positive-u condition described above, and runs
-:func:`utils.get_baseline_redundancies` to find redundant baselines.
-This is similar to running get_baseline_redundancies with the `with_conjugates`
-option, except that the baseline indices are returned such that none need to be
-flipped to be redundant. This is more like what's defined in the `hera_cal` package.
+and indices in the convention that ``ant1<ant2``, and runs
+:func:`utils.get_baseline_redundancies` to find redundant baselines. It will then apply the conjugates
+list to the groups it finds.
+
+There is a subtle difference between the purposes of the two functions. `utils.get_antenna_redundancies`
+gives you all redundant baselines from the antenna positions, and does not necessarily reflect the baselines
+in a file. This is similar to what is written in the `hera_cal` package.
+Alternatively, `utils.get_baseline_redundancies` may be given the actual
+baseline vectors in a file and it will search for redundancies among those.
 
 The method :meth:`~pyuvdata.UVData.get_redundancies` is provided as a convenience. If
 run with the `use_antpos` option, it will mimic the behavior of `utils.get_antenna_redundancies`.
-Otherwise it will return redundancies in the existing data using `utils.get_baesline_redundancies`.
+Otherwise it will return redundancies in the existing data using `utils.get_baseline_redundancies`.
+If run with `use_antpos` and the `conjugate_bls` option, it will also adjust the data_array and baseline_array
+so that the baselines in the returned groups correspond with the baselines listed on the object (i.e., except for
+antenna pairs with no associated data).
 
 ::
 
@@ -1039,7 +1046,8 @@ baselines from the antenna positions and fills in the full data array based on r
 
     # Note -- Compressing and inflating changes the baseline order, reorder before comparing.
     >>> uv0.inflate_by_redundancy(tol=tol)
-    >>> uv_backup.reorder_blts()
+    >>> uv_backup.reorder_blts(conj_convention="u>0", uvw_tol=tol)
+    >>> uv0.reorder_blts()
     >>> np.all(uv0.baseline_array == uv_backup.baseline_array)
     True
 
@@ -1151,6 +1159,9 @@ a) Calibration of UVData by UVCal
   >>> UVC.read_calfits('pyuvdata/data/zen.2458116.30448.HH.flagged_abs.calfits')
   >>> UV_calibrated = utils.uvcalibrate(UV, UVC, inplace=False)
 
+  # We can also un-calibrate using the same UVCal
+  >>> UV_uncalibrated = utils.uvcalibrate(UV_calibrated, UVC, inplace=False, undo=True)
+
 UVCal: Selecting data
 -----------------------
 The select method lets you select specific antennas (by number or name),
diff --git a/pyuvdata/calfits.py b/pyuvdata/calfits.py
index 00d7b5d5..8b92d74f 100644
--- a/pyuvdata/calfits.py
+++ b/pyuvdata/calfits.py
@@ -134,6 +134,8 @@ class CALFITS(UVCal):
             prihdr['BL_RANGE'] = '[' + ', '.join([str(b) for b in self.baseline_range]) + ']'
         if self.diffuse_model is not None:
             prihdr['DIFFUSE'] = self.diffuse_model
+        if self.gain_scale is not None:
+            prihdr['GNSCALE'] = self.gain_scale
         prihdr['INTTIME'] = self.integration_time
         prihdr['CHWIDTH'] = self.channel_width
         prihdr['XORIENT'] = self.x_orientation
@@ -403,6 +405,7 @@ class CALFITS(UVCal):
                 hdr.remove('HISTORY')
             self.time_range = list(map(float, hdr.pop('TMERANGE').split(',')))
             self.gain_convention = hdr.pop('GNCONVEN')
+            self.gain_scale = hdr.pop("GNSCALE", None)
             self.x_orientation = hdr.pop('XORIENT')
             self.cal_type = hdr.pop('CALTYPE')
             if self.cal_type == 'delay':
diff --git a/pyuvdata/utils.py b/pyuvdata/utils.py
index d9231cb6..15a5e985 100644
--- a/pyuvdata/utils.py
+++ b/pyuvdata/utils.py
@@ -394,7 +394,7 @@ def unphase_uvw(ra, dec, uvw):
 
 
 def uvcalibrate(uvdata, uvcal, inplace=True, prop_flags=True, flag_missing=True,
-                Dterm_cal=False, delay_convention='minus'):
+                Dterm_cal=False, delay_convention='minus', undo=False):
     """
     Calibrate a UVData object with a UVCal object.
 
@@ -417,6 +417,9 @@ def uvcalibrate(uvdata, uvcal, inplace=True, prop_flags=True, flag_missing=True,
     delay_convention : str, optional
         Exponent sign to use in conversion of 'delay' to 'gain' cal_type
         if the input uvcal is not inherently 'gain' cal_type. Default to 'minus'.
+    undo : bool, optional
+        If True, undo the provided calibration. i.e. apply the calibration with
+        flipped gain_convention. Flag propagation rules apply the same.
 
     Returns
     -------
@@ -466,13 +469,21 @@ def uvcalibrate(uvdata, uvcal, inplace=True, prop_flags=True, flag_missing=True,
                 uvdata.flag_array[blt_inds, 0, :, pol_ind] += mask
 
             # apply to data
-            if uvcal.gain_convention == 'multiply':
+            mult_gains = uvcal.gain_convention == 'multiply'
+            if undo:
+                mult_gains = not mult_gains
+            if mult_gains:
                 uvdata.data_array[blt_inds, 0, :, pol_ind] *= gain
-            elif uvcal.gain_convention == 'divide':
+            else:
                 uvdata.data_array[blt_inds, 0, :, pol_ind] /= gain
 
-    # update history
+    # update attributes
     uvdata.history += "\nCalibrated with pyuvdata.utils.uvcalibrate."
+    if undo:
+        uvdata.vis_units = 'UNCALIB'
+    else:
+        if uvcal.gain_scale is not None:
+            uvdata.vis_units = uvcal.gain_scale
 
     if not inplace:
         return uvdata
@@ -1062,6 +1073,14 @@ def antnums_to_baseline(ant1, ant2, Nants_telescope, attempt256=False):
         return np.int64(baseline)
 
 
+def baseline_index_flip(baseline, Nants_telescope):
+    """
+    Change baseline number to reverse antenna order.
+    """
+    ant1, ant2 = baseline_to_antnums(baseline, Nants_telescope)
+    return antnums_to_baseline(ant2, ant1, Nants_telescope)
+
+
 def get_baseline_redundancies(baselines, baseline_vecs, tol=1.0, with_conjugates=False):
     """
     Find redundant baseline groups
@@ -1098,14 +1117,14 @@ def get_baseline_redundancies(baselines, baseline_vecs, tol=1.0, with_conjugates
 
     if with_conjugates:
         conjugates = []
-        for bl in baseline_vecs:
-            if bl[0] == 0:
-                if bl[1] == 0:
-                    conjugates.append(bl[2] < 0)
-                else:
-                    conjugates.append(bl[1] < 0)
-            else:
-                conjugates.append(bl[0] < 0)
+        for bv in baseline_vecs:
+            uneg = bv[0] < -tol
+            uzer = np.isclose(bv[0], 0.0, atol=tol)
+            vneg = bv[1] < -tol
+            vzer = np.isclose(bv[1], 0.0, atol=tol)
+            wneg = bv[2] < -tol
+            conjugates.append(uneg or (uzer and vneg) or (uzer and vzer and wneg))
+
         conjugates = np.array(conjugates, dtype=bool)
         baseline_vecs[conjugates] *= (-1)
         baseline_ind_conj = baselines[conjugates]
@@ -1164,8 +1183,6 @@ def get_antenna_redundancies(antenna_numbers, antenna_positions, tol=1.0, includ
     """
     Find redundant baseline groups based on antenna positions.
 
-    Include all possible redundant baselines based on antenna positions.
-
     Parameters
     ----------
     antenna_numbers : array_like of int
@@ -1185,8 +1202,24 @@ def get_antenna_redundancies(antenna_numbers, antenna_positions, tol=1.0, includ
         List of vectors describing redundant group centers
     lengths : list of float
         List of redundant group baseline lengths in meters
+
+    Notes
+    -----
+
+    The baseline numbers refer to antenna pairs (a1, a2) such that
+    the baseline vector formed from ENU antenna positions,
+        blvec = enu[a1] - enu[a2]
+    is close to the other baselines in the group.
+
+    This is achieved by putting baselines in a form of the u>0
+    convention, but with a tolerance in defining the signs of
+    vector components.
+
+    To guarantee that the same baseline numbers are present in a UVData
+    object, ``UVData.conjugate_bls('u>0', uvw_tol=tol)``, where `tol` is
+    the tolerance used here.
     """
-    Nants = antenna_numbers.shape[0]
+    Nants = antenna_numbers.size
 
     bls = []
     bl_vecs = []
@@ -1197,18 +1230,20 @@ def get_antenna_redundancies(antenna_numbers, antenna_positions, tol=1.0, includ
             mini = aj
         for ai in range(mini, Nants):
             anti, antj = antenna_numbers[ai], antenna_numbers[aj]
-            bidx = antnums_to_baseline(anti, antj, Nants)
-            bv = antenna_positions[aj] - antenna_positions[ai]
-            # Enforce u-positive orientation
-            if (bv[0] < 0 or ((bv[0] == 0) and bv[1] < 0)
-                    or ((bv[0] == 0) and (bv[1] == 0) and bv[2] < 0)):
-                bv *= (-1)
-                bidx = antnums_to_baseline(antj, anti, Nants)
+            bidx = antnums_to_baseline(antj, anti, Nants)
+            bv = antenna_positions[ai] - antenna_positions[aj]
             bl_vecs.append(bv)
             bls.append(bidx)
     bls = np.array(bls)
     bl_vecs = np.array(bl_vecs)
-    return get_baseline_redundancies(bls, bl_vecs, tol=tol, with_conjugates=False)
+    gps, vecs, lens, conjs = get_baseline_redundancies(bls, bl_vecs, tol=tol, with_conjugates=True)
+    # Flip the baselines in the groups.
+    for gi, gp in enumerate(gps):
+        for bi, bl in enumerate(gp):
+            if bl in conjs:
+                gps[gi][bi] = baseline_index_flip(bl, Nants)
+
+    return gps, vecs, lens
 
 
 def _reraise_context(fmt, *args):
diff --git a/pyuvdata/uvcal.py b/pyuvdata/uvcal.py
index 6f86fa83..0f0a0a6f 100644
--- a/pyuvdata/uvcal.py
+++ b/pyuvdata/uvcal.py
@@ -232,6 +232,12 @@ class UVCal(UVBase):
                                               expected_type=str, description=desc)
 
         # --- truly optional parameters ---
+        desc = ('The gain scale of the calibration, which indicates the units of the '
+                'calibrated visibilities. For example, Jy or K str.')
+        self._gain_scale = uvp.UVParameter('gain_scale', form='str',
+                                           expected_type=str,
+                                           description=desc, required=False)
+
         desc = ('Array of input flags, True is flagged. shape: (Nants_data, Nspws, '
                 'Nfreqs, Ntimes, Njones), type = bool.')
         self._input_flag_array = uvp.UVParameter('input_flag_array',
diff --git a/pyuvdata/uvdata.py b/pyuvdata/uvdata.py
index 11fc38d5..85b26381 100644
--- a/pyuvdata/uvdata.py
+++ b/pyuvdata/uvdata.py
@@ -1033,7 +1033,7 @@ class UVData(UVBase):
             self.phase(phase_center_ra, phase_center_dec, phase_center_epoch,
                        phase_frame=output_phase_frame)
 
-    def conjugate_bls(self, convention='ant1<ant2', use_enu=True):
+    def conjugate_bls(self, convention='ant1<ant2', use_enu=True, uvw_tol=0.0):
         """
         Conjugate baselines according to one of the supported conventions.
 
@@ -1050,6 +1050,9 @@ class UVData(UVBase):
             Use true antenna positions to determine uv location (as opposed to
             uvw array). Only applies if `convention` is 'u<0', 'u>0', 'v<0', 'v>0'.
             Set to False to use uvw array values.
+        uvw_tol : float
+            Defines a tolerance on uvw coordinates for setting the
+            u>0, u<0, v>0, or v<0 conventions. Defaults to 0m.
 
         Raises
         ------
@@ -1087,29 +1090,29 @@ class UVData(UVBase):
             elif convention == 'ant2<ant1':
                 index_array = np.asarray(self.ant_2_array > self.ant_1_array).nonzero()
             elif convention == 'u<0':
-                index_array = np.asarray((uvw_array_use[:, 0] > 0)
-                                         | (uvw_array_use[:, 1] < 0) & (uvw_array_use[:, 0] == 0)
-                                         | ((uvw_array_use[:, 2] < 0)
-                                            & (uvw_array_use[:, 0] == 0)
-                                            & (uvw_array_use[:, 1] == 0))).nonzero()
+                index_array = np.asarray((uvw_array_use[:, 0] > uvw_tol)
+                                         | (uvw_array_use[:, 1] > uvw_tol) & np.isclose(uvw_array_use[:, 0], 0, atol=uvw_tol)
+                                         | (uvw_array_use[:, 2] > uvw_tol)
+                                         & np.isclose(uvw_array_use[:, 0], 0, atol=uvw_tol)
+                                         & np.isclose(uvw_array_use[:, 1], 0, atol=uvw_tol)).nonzero()
             elif convention == 'u>0':
-                index_array = np.asarray((uvw_array_use[:, 0] < 0)
-                                         | (uvw_array_use[:, 1] < 0) & (uvw_array_use[:, 0] == 0)
-                                         | ((uvw_array_use[:, 2] < 0)
-                                            & (uvw_array_use[:, 0] == 0)
-                                            & (uvw_array_use[:, 1] == 0))).nonzero()
+                index_array = np.asarray((uvw_array_use[:, 0] < -uvw_tol)
+                                         | ((uvw_array_use[:, 1] < -uvw_tol) & np.isclose(uvw_array_use[:, 0], 0, atol=uvw_tol))
+                                         | ((uvw_array_use[:, 2] < -uvw_tol)
+                                         & np.isclose(uvw_array_use[:, 0], 0, atol=uvw_tol)
+                                         & np.isclose(uvw_array_use[:, 1], 0, atol=uvw_tol))).nonzero()
             elif convention == 'v<0':
-                index_array = np.asarray((uvw_array_use[:, 1] > 0)
-                                         | (uvw_array_use[:, 0] < 0) & (uvw_array_use[:, 1] == 0)
-                                         | ((uvw_array_use[:, 2] < 0)
-                                            & (uvw_array_use[:, 0] == 0)
-                                            & (uvw_array_use[:, 1] == 0))).nonzero()
+                index_array = np.asarray((uvw_array_use[:, 1] > uvw_tol)
+                                         | (uvw_array_use[:, 0] > uvw_tol) & np.isclose(uvw_array_use[:, 1], 0, atol=uvw_tol)
+                                         | (uvw_array_use[:, 2] > uvw_tol)
+                                         & np.isclose(uvw_array_use[:, 0], 0, atol=uvw_tol)
+                                         & np.isclose(uvw_array_use[:, 1], 0, atol=uvw_tol)).nonzero()
             elif convention == 'v>0':
-                index_array = np.asarray((uvw_array_use[:, 1] < 0)
-                                         | (uvw_array_use[:, 0] < 0) & (uvw_array_use[:, 1] == 0)
-                                         | ((uvw_array_use[:, 2] < 0)
-                                            & (uvw_array_use[:, 0] == 0)
-                                            & (uvw_array_use[:, 1] == 0))).nonzero()
+                index_array = np.asarray((uvw_array_use[:, 1] < -uvw_tol)
+                                         | (uvw_array_use[:, 0] < -uvw_tol) & np.isclose(uvw_array_use[:, 1], 0, atol=uvw_tol)
+                                         | (uvw_array_use[:, 2] < -uvw_tol)
+                                         & np.isclose(uvw_array_use[:, 0], 0, atol=uvw_tol)
+                                         & np.isclose(uvw_array_use[:, 1], 0, atol=uvw_tol)).nonzero()
         else:
             index_array = convention
 
@@ -1212,7 +1215,7 @@ class UVData(UVBase):
                       'reorder_pols in version 1.5', DeprecationWarning)
         self.reorder_pols(order=order)
 
-    def reorder_blts(self, order='time', minor_order=None, conj_convention=None,
+    def reorder_blts(self, order='time', minor_order=None, conj_convention=None, uvw_tol=0.0,
                      conj_convention_use_enu=True, run_check=True, check_extra=True,
                      run_check_acceptability=True):
         """
@@ -1233,6 +1236,10 @@ class UVData(UVBase):
         conj_convention : str or array_like of int
             Optionally conjugate baselines to make the baselines have the
             desired orientation. See conjugate_bls for allowed values and details.
+        uvw_tol : float
+            If conjugating baselines, sets a tolerance for determining the signs
+            of u,v, and w, and whether or not they are zero.
+            See conjugate_bls for details.
         conj_convention_use_enu: bool
             If `conj_convention` is set, this is passed to conjugate_bls, see that
             method for details.
@@ -1289,7 +1296,7 @@ class UVData(UVBase):
 
         if conj_convention is not None:
             self.conjugate_bls(convention=conj_convention,
-                               use_enu=conj_convention_use_enu)
+                               use_enu=conj_convention_use_enu, uvw_tol=uvw_tol)
 
         if isinstance(order, str):
             if minor_order is None:
@@ -4424,7 +4431,8 @@ class UVData(UVBase):
         return np.diff(np.sort(list(set(self.time_array))))[0] * 86400
 
     def get_redundancies(self, tol=1.0, use_antpos=False,
-                         include_conjugates=True, include_autos=True):
+                         include_conjugates=False, include_autos=True,
+                         conjugate_bls=False):
         """
         Get redundant baselines to a given tolerance. This can be used to identify
         redundant baselines present in the data, or find all possible redundant
@@ -4444,6 +4452,9 @@ class UVData(UVBase):
         include_autos : bool
             Option to include autocorrelations in the full redundancy list.
             Only used if use_antpos is True.
+        conjugate_bls : bool
+            If using antenna positions, this will conjugate baselines on this
+            object to correspond with those in the returned groups.
 
         Returns
         -------
@@ -4453,7 +4464,7 @@ class UVData(UVBase):
             List of vectors describing redundant group uvw centers
         lengths : list of float
             List of redundant group baseline lengths in meters
-        conjugate_bls : list of int, or None, optional
+        conjugates : list of int, or None, optional
             List of indices for baselines that must be conjugated to fit into their
             redundant groups.
             Will return None if use_antpos is True and include_conjugates is True
@@ -4469,6 +4480,9 @@ class UVData(UVBase):
             antpos, numbers = self.get_ENU_antpos(center=False)
             result = uvutils.get_antenna_redundancies(numbers, antpos, tol=tol,
                                                       include_autos=include_autos)
+            if conjugate_bls:
+                self.conjugate_bls(convention='u>0', uvw_tol=tol)
+
             if include_conjugates:
                 result = result + (None,)
             return result
@@ -4489,10 +4503,8 @@ class UVData(UVBase):
         warnings.warn("UVData.get_antenna_redundancies has been replaced with get_redundancies,"
                       "and will be removed in version 1.6.",
                       DeprecationWarning)
-        if kwargs.pop('conjugate_bls', False):
-            self.conjugate_bls('u>0')
         kwargs['use_antpos'] = True
-        red_gps, blvecs, lens, conjs = self.get_redundancies(*args, **kwargs)
+        red_gps, blvecs, lens = self.get_redundancies(*args, **kwargs)
         return red_gps, blvecs, lens
 
     def get_baseline_redundancies(self, *args, **kwargs):
@@ -4535,7 +4547,7 @@ class UVData(UVBase):
             if inplace is False, return the compressed UVData object
         """
 
-        red_gps, centers, lengths, conjugates = self.get_redundancies(tol)
+        red_gps, centers, lengths, conjugates = self.get_redundancies(tol, include_conjugates=True)
 
         bl_ants = [self.baseline_to_antnums(gp[0]) for gp in red_gps]
         return self.select(bls=bl_ants, inplace=inplace, metadata_only=metadata_only,
@@ -4559,7 +4571,7 @@ class UVData(UVBase):
         """
 
         self.conjugate_bls(convention='u>0')
-        red_gps, centers, lengths, _, = self.get_redundancies(tol=tol, use_antpos=True)
+        red_gps, centers, lengths = self.get_redundancies(tol=tol, use_antpos=True, conjugate_bls=True)
 
         # Stack redundant groups into one array.
         group_index, bl_array_full = zip(*[(i, bl) for i, gp in enumerate(red_gps) for bl in gp])


Test patch:
diff --git a/pyuvdata/tests/test_utils.py b/pyuvdata/tests/test_utils.py
index 79946d72..b5660b38 100644
--- a/pyuvdata/tests/test_utils.py
+++ b/pyuvdata/tests/test_utils.py
@@ -475,14 +475,13 @@ def test_redundancy_finder():
     uvd.read_uvfits(os.path.join(DATA_PATH, 'fewant_randsrc_airybeam_Nsrc100_10MHz.uvfits'))
 
     uvd.select(times=uvd.time_array[0])
-    uvd.unphase_to_drift()   # uvw_array is now equivalent to baseline positions
-    uvtest.checkWarnings(uvd.conjugate_bls, func_kwargs={'convention': 'u>0', 'use_enu': True},
-                         message=['The default for the `center`'],
-                         nwarnings=1, category=DeprecationWarning)
+    uvd.unphase_to_drift(use_ant_pos=True)   # uvw_array is now equivalent to baseline positions
+    uvd.conjugate_bls(convention='ant1<ant2', use_enu=True)
 
     tol = 0.05  # meters
 
     bl_positions = uvd.uvw_array
+    bl_pos_backup = copy.deepcopy(uvd.uvw_array)
 
     pytest.raises(ValueError, uvutils.get_baseline_redundancies,
                   uvd.baseline_array, bl_positions[0:2, 0:1])
@@ -545,31 +544,24 @@ def test_redundancy_finder():
     assert len(baseline_groups_ants) == 19
 
     # Check with conjugated baseline redundancies returned
-    u16_0 = bl_positions[16, 0]
     # Ensure at least one baseline has u==0 and v!=0 (for coverage of this case)
     bl_positions[16, 0] = 0
     baseline_groups, vec_bin_centers, lens, conjugates = uvutils.get_baseline_redundancies(
         uvd.baseline_array, bl_positions, tol=tol, with_conjugates=True)
 
     # restore baseline (16,0) and repeat to get correct groups
-    bl_positions[16, 0] = u16_0
+    bl_positions = bl_pos_backup
     baseline_groups, vec_bin_centers, lens, conjugates = uvutils.get_baseline_redundancies(
         uvd.baseline_array, bl_positions, tol=tol, with_conjugates=True)
 
-    # Should get the same groups as with the antenna method:
-    baseline_groups_flipped = []
-    for bgp in baseline_groups:
-        bgp_new = []
-        for bl in bgp:
-            ai, aj = uvutils.baseline_to_antnums(bl, uvd.Nants_telescope)
+    # Apply flips to compare with get_antenna_redundancies().
+    bl_gps_unconj = copy.deepcopy(baseline_groups)
+    for gi, gp in enumerate(bl_gps_unconj):
+        for bi, bl in enumerate(gp):
             if bl in conjugates:
-                bgp_new.append(uvutils.antnums_to_baseline(aj, ai, uvd.Nants_telescope))
-            else:
-                bgp_new.append(uvutils.antnums_to_baseline(ai, aj, uvd.Nants_telescope))
-        bgp_new.sort()
-        baseline_groups_flipped.append(bgp_new)
-    baseline_groups = [sorted(bgp) for bgp in baseline_groups]
-    assert np.all(sorted(baseline_groups_ants) == sorted(baseline_groups_flipped))
+                bl_gps_unconj[gi][bi] = uvutils.baseline_index_flip(bl, len(antnums))
+    bl_gps_unconj = [sorted(bgp) for bgp in bl_gps_unconj]
+    assert np.all(sorted(baseline_groups_ants) == sorted(bl_gps_unconj))
     for gi, gp in enumerate(baseline_groups):
         for bl in gp:
             bl_ind = np.where(uvd.baseline_array == bl)
@@ -597,7 +589,12 @@ def test_redundancy_conjugates():
 
     expected_conjugates = []
     for i, (u, v, w) in enumerate(bl_vecs):
-        if (u < 0) or (v < 0 and u == 0) or (w < 0 and u == v == 0):
+        uneg = u < -tol
+        uzer = np.isclose(u, 0.0, atol=tol)
+        vneg = v < -tol
+        vzer = np.isclose(v, 0.0, atol=tol)
+        wneg = w < -tol
+        if uneg or (uzer and vneg) or (uzer and vzer and wneg):
             expected_conjugates.append(bl_inds[i])
     bl_gps, vecs, lens, conjugates = uvutils.get_baseline_redundancies(
         bl_inds, bl_vecs, tol=tol, with_conjugates=True)
@@ -941,6 +938,8 @@ def test_uvcalibrate_apply_gains():
     uvd.read(os.path.join(DATA_PATH, 'zen.2457698.40355.xx.HH.uvcAA'))
     uvc = pyuvdata.UVCal()
     uvc.read_calfits(os.path.join(DATA_PATH, 'zen.2457698.40355.xx.gain.calfits'))
+    # assign gain scale manually
+    uvc.gain_scale = 'Jy'
     # downselect to match each other
     uvd.select(frequencies=uvd.freq_array[0, :10])
     uvc.select(times=uvc.time_array[:3])
@@ -952,11 +951,18 @@ def test_uvcalibrate_apply_gains():
     uvc.gain_convention = 'divide'
     uvdcal = uvutils.uvcalibrate(uvd, uvc, prop_flags=True, flag_missing=False, inplace=False)
     np.testing.assert_array_almost_equal(uvdcal.get_data(key), uvd.get_data(key) / (uvc.get_gains(ant1) * uvc.get_gains(ant2).conj()).T)
+    assert uvdcal.vis_units == 'Jy'
+
+    # test undo
+    uvdcal = uvutils.uvcalibrate(uvdcal, uvc, prop_flags=True, flag_missing=False, inplace=False, undo=True)
+    np.testing.assert_array_almost_equal(uvd.get_data(key), uvdcal.get_data(key))
+    assert uvdcal.vis_units == 'UNCALIB'
 
     # multiplication calibrate
     uvc.gain_convention = 'multiply'
     uvdcal = uvutils.uvcalibrate(uvd, uvc, prop_flags=False, flag_missing=False, inplace=False)
     np.testing.assert_array_almost_equal(uvdcal.get_data(key), uvd.get_data(key) * (uvc.get_gains(ant1) * uvc.get_gains(ant2).conj()).T)
+    assert uvdcal.vis_units == 'Jy'
 
     # test delay conversion runs through
     uvc.read_calfits(os.path.join(DATA_PATH, 'zen.2457698.40355.xx.delay.calfits'))
@@ -999,3 +1005,27 @@ def test_uvcalibrate_flag_propagation():
     assert not uvdcal.get_flags(20, 72, 'xx').max()  # assert no flags exist
     uvdcal = uvutils.uvcalibrate(uvd, uvc_sub, prop_flags=True, flag_missing=True, inplace=False)
     assert uvdcal.get_flags(20, 72, 'xx').min()  # assert completely flagged
+
+
+def test_upos_tol_reds():
+    # Checks that the u-positive convention in get_antenna_redundancies
+    # is enforced to the specificed tolerance.
+
+    # Make a layout with two NS baselines, one with u ~ -2*eps, and another with u == 0
+    # This would previously cause one to be flipped, when they should be redundant.
+
+    eps = 1e-5
+    tol = 3 * eps
+
+    ant_pos = np.array([
+        [-eps, 1., 0.],
+        [1., 1., 0.],
+        [eps, 0., 0.],
+        [1., 0., 0.]
+    ])
+
+    ant_nums = np.arange(4)
+
+    red_grps, _, _ = uvutils.get_antenna_redundancies(ant_nums, ant_pos, tol=tol)
+
+    assert len(red_grps) == 4
diff --git a/pyuvdata/tests/test_uvcal.py b/pyuvdata/tests/test_uvcal.py
index c4b12ee3..012e72f6 100644
--- a/pyuvdata/tests/test_uvcal.py
+++ b/pyuvdata/tests/test_uvcal.py
@@ -51,7 +51,7 @@ def uvcal_data():
                         '_input_flag_array', '_freq_range',
                         '_observer', '_git_origin_cal',
                         '_git_hash_cal', '_total_quality_array',
-                        '_extra_keywords']
+                        '_extra_keywords', '_gain_scale']
 
     extra_properties = ['gain_array', 'delay_array', 'sky_field',
                         'sky_catalog', 'ref_antenna_name', 'Nsources',
@@ -59,7 +59,7 @@ def uvcal_data():
                         'input_flag_array', 'freq_range',
                         'observer', 'git_origin_cal',
                         'git_hash_cal', 'total_quality_array',
-                        'extra_keywords']
+                        'extra_keywords', 'gain_scale']
 
     other_properties = ['pyuvdata_version_str']
 
@@ -1482,6 +1482,27 @@ def test_uvcal_get_methods():
     pytest.raises(ValueError, uvc.get_gains, 10)
 
 
+def test_write_read_optional_attrs():
+    # read a test file
+    cal_in = UVCal()
+    testfile = os.path.join(DATA_PATH, 'zen.2457698.40355.xx.gain.calfits')
+    cal_in.read_calfits(testfile)
+
+    # set some optional parameters
+    cal_in.gain_scale = 'Jy'
+    cal_in.sky_field = 'GLEAM'
+
+    # write
+    write_file_calfits = os.path.join(DATA_PATH, 'test/test.calfits')
+    cal_in.write_calfits(write_file_calfits, clobber=True)
+
+    # read and compare
+    cal_in2 = UVCal()
+    cal_in2.read_calfits(write_file_calfits)
+    assert cal_in == cal_in2
+    os.remove(write_file_calfits)
+
+
 def test_deprecated_x_orientation():
     cal_in = UVCal()
     testfile = os.path.join(DATA_PATH, 'zen.2457698.40355.xx.gain.calfits')
diff --git a/pyuvdata/tests/test_uvdata.py b/pyuvdata/tests/test_uvdata.py
index 3e86d484..ac107443 100644
--- a/pyuvdata/tests/test_uvdata.py
+++ b/pyuvdata/tests/test_uvdata.py
@@ -3499,21 +3499,23 @@ def test_deprecated_redundancy_funcs():
                                                           category=DeprecationWarning,
                                                           message='UVData.get_baseline_redundancies has been replaced')
 
-    red_gps_new, _, _, _ = uv0.get_redundancies(include_autos=False, use_antpos=True)
+    red_gps_new, _, _, = uv0.get_redundancies(include_autos=False, use_antpos=True)
     assert red_gps_new == redant_gps
 
 
 def test_get_antenna_redundancies():
     uv0 = UVData()
     uv0.read_uvfits(os.path.join(DATA_PATH, 'fewant_randsrc_airybeam_Nsrc100_10MHz.uvfits'))
-    uv0.conjugate_bls('u>0')
-    red_gps, centers, lengths, _ = uv0.get_redundancies(use_antpos=True, include_autos=False)
+
+    old_bl_array = np.copy(uv0.baseline_array)
+    red_gps, centers, lengths = uv0.get_redundancies(use_antpos=True, include_autos=False, conjugate_bls=True)
+    # new and old baseline Numbers are not the same (different conjugation)
+    assert not np.allclose(uv0.baseline_array, old_bl_array)
 
     # assert all baselines are in the data (because it's conjugated to match)
     for i, gp in enumerate(red_gps):
         for bl in gp:
             assert bl in uv0.baseline_array
-    old_bl_array = np.copy(uv0.baseline_array)
 
     # conjugate data differently
     uv0.conjugate_bls(convention='ant1<ant2')
@@ -3523,8 +3525,9 @@ def test_get_antenna_redundancies():
 
     assert conjs is None
 
-    # new and old baseline Numbers are not the same (different conjugation)
-    assert not np.allclose(uv0.baseline_array, old_bl_array)
+    apos, anums = uv0.get_ENU_antpos()
+    new_red_gps, new_centers, new_lengths = uvutils.get_antenna_redundancies(
+        anums, apos, include_autos=False)
 
     # all redundancy info is the same
     assert red_gps == new_red_gps
@@ -3539,11 +3542,10 @@ def test_redundancy_contract_expand():
     uv0 = UVData()
     uv0.read_uvfits(os.path.join(DATA_PATH, 'fewant_randsrc_airybeam_Nsrc100_10MHz.uvfits'))
 
-    tol = 0.02   # Fails at lower precision because some baselines falling into multiple redundant groups
+    tol = 0.02   # Fails at lower precision because some baselines fall into multiple redundant groups
 
     # Assign identical data to each redundant group:
-    uv0.conjugate_bls('u>0')
-    red_gps, centers, lengths, _ = uv0.get_redundancies(tol=tol, use_antpos=True)
+    red_gps, centers, lengths = uv0.get_redundancies(tol=tol, use_antpos=True, conjugate_bls=True)
     for i, gp in enumerate(red_gps):
         for bl in gp:
             inds = np.where(bl == uv0.baseline_array)
@@ -3561,31 +3563,31 @@ def test_redundancy_contract_expand():
     uvtest.checkWarnings(
         uv2.inflate_by_redundancy,
         [tol],
-        nwarnings=2,
-        category=[DeprecationWarning, UserWarning],
-        message=['The default for the `center` keyword has changed.',
+        nwarnings=3,
+        category=[DeprecationWarning, DeprecationWarning, UserWarning],
+        message=['The default for the `center` keyword', 'The default for the `center` keyword',
                  'Missing some redundant groups. Filling in available data.']
     )
     uv2.history = uv0.history
     # Inflation changes the baseline ordering into the order of the redundant groups.
     # reorder bls for comparison
-    uv0.reorder_blts()
-    uv2.reorder_blts()
+    uv0.reorder_blts(conj_convention='u>0')
+    uv2.reorder_blts(conj_convention='u>0')
     uv2._uvw_array.tols = [0, tol]
-
     assert uv2 == uv0
 
     uv3 = uv2.compress_by_redundancy(tol=tol, inplace=False)
     uvtest.checkWarnings(
         uv3.inflate_by_redundancy,
         [tol],
-        nwarnings=2,
-        category=[DeprecationWarning, UserWarning],
-        message=['The default for the `center` keyword has changed.',
+        nwarnings=3,
+        category=[DeprecationWarning, DeprecationWarning, UserWarning],
+        message=['The default for the `center` keyword', 'The default for the `center` keyword',
                  'Missing some redundant groups. Filling in available data.']
     )
     # Confirm that we get the same result looping inflate -> compress -> inflate.
-    uv3.reorder_blts()
+    uv3.reorder_blts(conj_convention='u>0')
+    uv2.reorder_blts(conj_convention='u>0')
 
     uv2.history = uv3.history
     assert uv2 == uv3
@@ -3602,8 +3604,7 @@ def test_redundancy_contract_expand_nblts_not_nbls_times_ntimes():
     tol = 1.0
 
     # Assign identical data to each redundant group:
-    uv0.conjugate_bls('u>0')
-    red_gps, centers, lengths, _ = uv0.get_redundancies(tol=tol, use_antpos=True)
+    red_gps, centers, lengths = uv0.get_redundancies(tol=tol, use_antpos=True, conjugate_bls=True)
     for i, gp in enumerate(red_gps):
         for bl in gp:
             inds = np.where(bl == uv0.baseline_array)
@@ -3614,9 +3615,9 @@ def test_redundancy_contract_expand_nblts_not_nbls_times_ntimes():
 
     # check inflating gets back to the original
     uvtest.checkWarnings(uv2.inflate_by_redundancy, {tol: tol},
-                         nwarnings=2, category=[DeprecationWarning, UserWarning],
-                         message=['The default for the `center` keyword has changed.',
-                                  'Missing some redundant groups. Filling in available data.'])
+                         nwarnings=3, category=[DeprecationWarning, DeprecationWarning, UserWarning],
+                         message=['The default for the `center` keyword'] * 2
+                         + ['Missing some redundant groups. Filling in available data.'])
 
     uv2.history = uv0.history
     # Inflation changes the baseline ordering into the order of the redundant groups.
@@ -3652,8 +3653,7 @@ def test_compress_redundancy_metadata_only():
     tol = 0.01
 
     # Assign identical data to each redundant group:
-    uv0.conjugate_bls('u>0')
-    red_gps, centers, lengths, _ = uv0.get_redundancies(tol=tol, use_antpos=True)
+    red_gps, centers, lengths = uv0.get_redundancies(tol=tol, use_antpos=True, conjugate_bls=True)
     for i, gp in enumerate(red_gps):
         for bl in gp:
             inds = np.where(bl == uv0.baseline_array)
@@ -3709,9 +3709,9 @@ def test_redundancy_missing_groups():
     uvtest.checkWarnings(
         uv1.inflate_by_redundancy,
         [tol],
-        nwarnings=2,
-        category=[DeprecationWarning, UserWarning],
-        message=['The default for the `center` keyword has changed.',
+        nwarnings=3,
+        category=[DeprecationWarning, DeprecationWarning, UserWarning],
+        message=['The default for the `center` keyword', 'The default for the `center` keyword',
                  'Missing some redundant groups. Filling in available data.']
     )
 
@@ -3754,7 +3754,7 @@ def test_quick_redundant_vs_redundant_test_array():
     groups = [[bl for bl in grp if bl != -1] for grp in groups]
     groups.sort(key=len)
 
-    redundant_groups, centers, lengths, conj_inds = uv.get_redundancies(tol=tol)
+    redundant_groups, centers, lengths, conj_inds = uv.get_redundancies(tol=tol, include_conjugates=True)
     redundant_groups.sort(key=len)
     assert groups == redundant_groups
 
@@ -3795,7 +3795,7 @@ def test_redundancy_finder_when_nblts_not_nbls_times_ntimes():
     groups = [[bl for bl in grp if bl != -1] for grp in groups]
     groups.sort(key=len)
 
-    redundant_groups, centers, lengths, conj_inds = uv.get_redundancies(tol=tol)
+    redundant_groups, centers, lengths, conj_inds = uv.get_redundancies(tol=tol, include_conjugates=True)
     redundant_groups.sort(key=len)
     assert groups == redundant_groups