diff --git a/orbitize/hipparcos.py b/orbitize/hipparcos.py
index d5ec2f31..cb770712 100644
--- a/orbitize/hipparcos.py
+++ b/orbitize/hipparcos.py
@@ -201,6 +201,7 @@ def __init__(
         self.alphadec0_epoch = alphadec0_epoch
 
         self.hip_data = hip_data
+        self.brandt_correction = brandt_correction
 
         # dvd files don't contain the Hipparcos astrometric solution, so
         # we need to look it up
@@ -371,6 +372,7 @@ def _save(self, hf):
         if self.hip_data is not None:
             hip_data_group = hf.create_group("hip_data")
             hip_data_group.update(self.hip_data)
+        hf.attrs["brandt_correction"] = self.brandt_correction
 
     def compute_lnlike(self, raoff_model, deoff_model, samples, param_idx):
         """
diff --git a/orbitize/results.py b/orbitize/results.py
index 4d597246..1c68f62e 100644
--- a/orbitize/results.py
+++ b/orbitize/results.py
@@ -8,6 +8,7 @@
 import orbitize.basis
 import orbitize.plot
 import orbitize.gaia, orbitize.hipparcos
+from orbitize import DATADIR
 
 class Results(object):
     """
@@ -315,6 +316,7 @@ def load_results(self, filename, append=False):
                 hip_data = {key : hip_data_raw[key][()] for key in hip_data_raw.keys()}
             else:
                 hip_data = None
+            brandt_correction = bool(hf.attrs.get("brandt_correction", False))
 
             hipparcos_IAD = orbitize.hipparcos.HipparcosLogProb(
                 tmpfile,
@@ -322,7 +324,8 @@ def load_results(self, filename, append=False):
                 num_secondary_bodies,
                 alphadec0_epoch,
                 renormalize_errors,
-                hip_data,
+                hip_data=hip_data,
+                brandt_correction=brandt_correction
             )
 
             # load Gaia data
@@ -339,6 +342,7 @@ def load_results(self, filename, append=False):
 
             # alternatively load HGCA. Note requires hipparcos_IAD attribute
             gaiagost_data = hf.get("Gaia_GOST")
+            hgca_filepath = None
             if gaiagost_data is not None:
                 
                 tmpfile = 'tmpfile_gost_OkToDeleteAfterFitFinishes'
@@ -346,11 +350,16 @@ def load_results(self, filename, append=False):
                 tmptbl.write(tmpfile, format="ascii", overwrite=True)
 
                 hgca_entry_raw = hf.get("HGCA_entry")
+                hgca_entry = None
                 if hgca_entry_raw is None:
-                    raise KeyError("HGCA Entry not found in {} for initializing HGCALogProb".format(filename))
-                hgca_entry = {key : hgca_entry_raw[key][()] for key in hgca_entry_raw.keys()}
-
-                gaia = orbitize.gaia.HGCALogProb(int(hip_num), hipparcos_IAD, tmpfile, hgca_entry=hgca_entry)
+                    if os.path.exists(os.path.join(DATADIR, "HGCA_vEDR3.fits")):
+                        hgca_filepath = os.path.join(DATADIR, "HGCA_vEDR3.fits")
+                    else:
+                        raise KeyError("HGCA Entry not found in {} for initializing HGCALogProb".format(filename))
+                else:
+                    hgca_entry = {key : hgca_entry_raw[key][()] for key in hgca_entry_raw.keys()}
+
+                gaia = orbitize.gaia.HGCALogProb(int(hip_num), hipparcos_IAD, tmpfile, hgca_filepath=hgca_filepath, hgca_entry=hgca_entry)
                 hipparcos_IAD = None # HGCA handles hipparocs, so don't want to pass Hipparcos also into the system
 
 
diff --git a/tests/test_multiplanet.py b/tests/test_multiplanet.py
index 0ba452c6..a76173ef 100644
--- a/tests/test_multiplanet.py
+++ b/tests/test_multiplanet.py
@@ -231,12 +231,12 @@ def test_fit_selfconsist():
     samp.curr_pos[0, :, 2] = np.random.normal(c_params[0], 0.01, n_walkers)  # sma
     samp.curr_pos[0, :, 3] = np.random.normal(c_params[-1], 0.01, n_walkers)  # tau
     # we will make a fairly broad mass starting position
-    samp.curr_pos[0, :, 4] = np.random.uniform(mass_b * 0.25, mass_b * 4, n_walkers)
-    samp.curr_pos[0, :, 5] = np.random.uniform(mass_c * 0.25, mass_c * 4, n_walkers)
+    samp.curr_pos[0, :, 4] = np.random.uniform(mass_b * 0.5, mass_b * 2, n_walkers)
+    samp.curr_pos[0, :, 5] = np.random.uniform(mass_c * 0.5, mass_c * 2, n_walkers)
     samp.curr_pos[0, 0, 4] = mass_b
     samp.curr_pos[0, 0, 5] = mass_c
 
-    samp.run_sampler(n_walkers * 200, burn_steps=200)
+    samp.run_sampler(n_walkers * 1000, burn_steps=1000)
 
     res = samp.results