added RUWE warning

logan-pearce · Nov 26, 2019 · addbe94 · addbe94
1 parent e246324
commit addbe94
Show file tree

Hide file tree

Showing 8 changed files with 113 additions and 14 deletions.
diff --git a/LOFTI_Tutorial.ipynb b/LOFTI_Tutorial.ipynb
@@ -3778,14 +3778,87 @@
    "cell_type": "markdown",
    "metadata": {},
    "source": [
-    "Parameter    Mean    Std    Mode    68% Min Cred Int    95% Min Cred Int\n",
-    "$a \\; (AU)$    440.273860033696    778.7828994666231    167.91369902260953    (120.78742643128001, 335.5722381365086)    (120.78742643128001, 1440.6607692140242)\n",
-    "$e$    0.7363533729151298    0.25417521870515686    0.9354707541002492    (0.7079602799311652, 0.9944789694374682)    (0.11403236837364339, 0.9944789694374682)\n",
-    "$ i \\; (deg)$    102.21552571902514    12.315962182778556    96.93037897540756    (93.85757612975858, 104.21407925155971)    (91.20949059969848, 129.6667081849308)\n",
-    "$ \\omega \\; (deg)$    170.14767463219783    99.64182945955947    110.92869998371395    (13.407759062840405, 220.04147746284139)    (5.060697734275732, 336.04000373968086)\n",
-    "$\\Omega \\; (deg)$    37.72898415688386    94.40123009172831    -29.405473241926238    (-45.72234936546789, 143.51214034948345)    (-119.84686434942489, 179.66687380258452)\n",
-    "$T_0 \\; (yr)$    -5266.072718030055    41405.54729054067    1390.9519132615283    (-387.5345674226005, 1637.8029937252832)    (-13979.8862767336, 1939.049329542796)\n",
-    "$a\\,(1-e) \\; (AU)$    109.41920890604749    219.0675705854163    26.869332294301532    (0.9391792185192334, 117.30940253868981)    (0.9391792185192334, 263.2353692706762)"
+    "Parameter    Mean    Std    Mode    68% Min Cred Int    95% Min Cred Int <br>\n",
+    "$a \\; (AU)$    440.273860033696    778.7828994666231    167.91369902260953    (120.78742643128001, 335.5722381365086)    (120.78742643128001, 1440.6607692140242) <br>\n",
+    "$e$    0.7363533729151298    0.25417521870515686    0.9354707541002492    (0.7079602799311652, 0.9944789694374682)    (0.11403236837364339, 0.9944789694374682) <br>\n",
+    "$ i \\; (deg)$    102.21552571902514    12.315962182778556    96.93037897540756    (93.85757612975858, 104.21407925155971)    (91.20949059969848, 129.6667081849308)<br>\n",
+    "$ \\omega \\; (deg)$    170.14767463219783    99.64182945955947    110.92869998371395    (13.407759062840405, 220.04147746284139)    (5.060697734275732, 336.04000373968086)<br>\n",
+    "$\\Omega \\; (deg)$    37.72898415688386    94.40123009172831    -29.405473241926238    (-45.72234936546789, 143.51214034948345)    (-119.84686434942489, 179.66687380258452)<br>\n",
+    "$T_0 \\; (yr)$    -5266.072718030055    41405.54729054067    1390.9519132615283    (-387.5345674226005, 1637.8029937252832)    (-13979.8862767336, 1939.049329542796)<br>\n",
+    "$a\\,(1-e) \\; (AU)$    109.41920890604749    219.0675705854163    26.869332294301532    (0.9391792185192334, 117.30940253868981)    (0.9391792185192334, 263.2353692706762)<br>"
+   ]
+  },
+  {
+   "cell_type": "markdown",
+   "metadata": {},
+   "source": [
+    "# GIGO: Poor Gaia solutions will give poor orbit fits!"
+   ]
+  },
+  {
+   "cell_type": "markdown",
+   "metadata": {},
+   "source": [
+    "The Gaia astrometric solution assumes that all point sources are single stars, and that their motion was linear (no acceleration) during the time series images used to arrive at the astrometric solution.  If a source is actually an unresolved multiple system, or if the binary is resolved but orbital acceleration was large enough to affect the astrometric solution, the solution will not be reliable.  This is best quantified in the Renormalized Unit Wight Error (RUWE) parameter (see https://gea.esac.esa.int/archive/documentation/GDR2/Gaia_archive/chap_datamodel/sec_dm_main_tables/ssec_dm_ruwe.html).  And RUWE >~ 1.2 indicates an unreliable solution. <br><br>\n",
+    "LOFTI will look at the RUWE and raise a warning if it's greater than 1.2, and ask if you really want to do this.<br><br>\n",
+    "Kepler 444 is a triple system (Dupuy et al. 2016).  Kepler 444 A is resolved in Gaia, but Kepler 444 B and C aren't resolved. "
+   ]
+  },
+  {
+   "cell_type": "code",
+   "execution_count": 3,
+   "metadata": {},
+   "outputs": [
+    {
+     "name": "stdout",
+     "output_type": "stream",
+     "text": [
+      "Computing constraints.\n",
+      "Created TAP+ (v1.0.1) - Connection:\n",
+      "\tHost: gea.esac.esa.int\n",
+      "\tUse HTTPS: False\n",
+      "\tPort: 80\n",
+      "\tSSL Port: 443\n",
+      "Finished computing constraints:\n",
+      "Delta RA, err in mas: 1754.9243128632852 0.3189517009128414\n",
+      "Delta Dec, err in mas: 544.4546301538144 0.5744290942303749\n",
+      "\n",
+      "pmRA, err in km/s: -1.4105342203867228 0.16512004559558907\n",
+      "pmDec, err in km/s: 1.9312076961988849 0.175587192939648\n",
+      "deltaRV, err im km/s (pos towards observer): 0.0 0.0\n",
+      "\n",
+      "Total relative velocity [km/s]: 2.391478612227127 +/- 0.24103006406237912\n",
+      "Total plane-of-sky relative velocity [mas/yr]: 15.03161860446753 +/- 1.5086728000824803\n",
+      "\n",
+      "sep,err [mas] 1837.4412911195748 0.34862183330431934 pa,err [deg]: 72.76381454087225 0.01736614347812121\n",
+      "sep [AU] 61.666994140223395\n",
+      "sep, err [km] (9225251015.846796, 0.0) (1750327.4457696672, 0.0)\n",
+      "D_star 33.561341218498995 +\\- 0.28704379047656037\n",
+      "Delta Gmag 3.6348228\n",
+      "RUWE source 1: 16.686954\n",
+      "RUWE source 2: 1.0002892\n",
+      "\n",
+      "Does this look good? Hit enter to start the fit, n to exit: \n",
+      "Yeehaw let's go\n",
+      "WARNING: RUWE for one or more of your solutions is greater than 1.2. This indicates \n",
+      "            that the source might be an unresolved binary or experiencing acceleration \n",
+      "            during the observation.  Orbit fit results may not be trustworthy.  Do you \n",
+      "            wish to continue?\n",
+      "            Hit enter to proceed, n to exit: n\n"
+     ]
+    }
+   ],
+   "source": [
+    "Kepler444BC = 2101486923382009472\n",
+    "Kepler444A = 2101486923385239808\n",
+    "fitorbit(Kepler444BC, Kepler444A, verbose=True)"
+   ]
+  },
+  {
+   "cell_type": "markdown",
+   "metadata": {},
+   "source": [
+    "The RUWE for Kepler 444 BC is much too large to give a reliable orbit fit, and the script raised a warning to let us know this.  GL 896 is another case where the RUWE for one source is too large, this time because the orbit exhibited acceleration during the Gaia observations.  See Pearce et al. 2019 for further discussion."
    ]
   },
   {
@@ -3812,7 +3885,7 @@
    "name": "python",
    "nbconvert_exporter": "python",
    "pygments_lexer": "ipython3",
-   "version": "3.7.0"
+   "version": "3.7.3"
   }
  },
  "nbformat": 4,

diff --git a/LOFTI_Tutorial.pdf b/LOFTI_Tutorial.pdf
diff --git a/dist/lofti_gaia-0.4.tar.gz b/dist/lofti_gaia-0.4.tar.gz
diff --git a/lofti_gaia.egg-info/PKG-INFO b/lofti_gaia.egg-info/PKG-INFO
@@ -1,6 +1,6 @@
 Metadata-Version: 1.0
 Name: lofti-gaia
-Version: 0.3
+Version: 0.4
 Summary: Orbit fitting with Gaia astrometry
 Home-page: https://github.com/logan-pearce/lofti_gaia
 Author: Logan Pearce

diff --git a/lofti_gaia/.DS_Store b/lofti_gaia/.DS_Store
diff --git a/lofti_gaia/__pycache__/lofti.cpython-37.pyc b/lofti_gaia/__pycache__/lofti.cpython-37.pyc
diff --git a/lofti_gaia/lofti.py b/lofti_gaia/lofti.py
@@ -23,6 +23,9 @@ def prepareconstraints(source_id1, source_id2):
           pa [deg] (tuple, flt): position angle[0] and error[1] in degrees east of north
           delta_mag (flt): contrast in magnitudes from primary to secondary star
           d_star [pc] (flt): distance to system in parsecs
+          ruwe (tuple, flt): renormilized unit weight error.  RUWE <~ 1.2 indicates reliable astrometric
+              solution.  
+              See https://gea.esac.esa.int/archive/documentation/GDR2/Gaia_archive/chap_datamodel/sec_dm_main_tables/ssec_dm_ruwe.html
     """
 
     from astroquery.gaia import Gaia
@@ -42,6 +45,13 @@ def prepareconstraints(source_id1, source_id2):
 
     job = Gaia.launch_job("SELECT * FROM gaiadr2.gaia_source WHERE source_id = "+str(source_id2))
     k = job.get_results()
+
+    # Retrieve RUWE for both sources
+    job = Gaia.launch_job("SELECT * FROM gaiadr2.ruwe WHERE source_id = "+str(source_id1))
+    jruwe = job.get_results()
+
+    job = Gaia.launch_job("SELECT * FROM gaiadr2.ruwe WHERE source_id = "+str(source_id2))
+    kruwe = job.get_results()
 
     # Parallaxes:
     plxa, plxaerr = j[0]['parallax'], j[0]['parallax_error']
@@ -122,7 +132,7 @@ def prepareconstraints(source_id1, source_id2):
     return [deltaRA, deltaRA_err], [deltaDec, deltaDec_err], [pmRA_kms, pmRA_err_kms], \
            [pmDec_kms, pmDec_err_kms], [deltarv, deltarverr], [total_pos_velocity, total_pos_velocity_error], \
            [total_velocity_kms, total_velocity_error_kms], [rho, rhoerr], [pa, paerr], \
-           delta_mag, [d_star,d_star_err]
+           delta_mag, [d_star,d_star_err], [jruwe['ruwe'],kruwe['ruwe']]
 
 
 
@@ -206,11 +216,13 @@ def fitorbit(source_id1, source_id2,
     import time as tm
     from lofti_gaia.loftifittingtools import draw_priors, calc_OFTI, to_si, update_progress
     import pickle
+    import warnings
+    warnings.filterwarnings("ignore")
 
     print('Computing constraints.')
     # Compute constraints:
     deltaRA, deltaDec, pmRA_kms, pmDec_kms, deltarv, total_pos_velocity, total_velocity_kms, \
-    rho, pa, delta_mag, d_star = prepareconstraints(source_id1, source_id2)
+    rho, pa, delta_mag, d_star, ruwe = prepareconstraints(source_id1, source_id2)
 
     if verbose == True:
         print('Finished computing constraints:')
@@ -230,12 +242,26 @@ def fitorbit(source_id1, source_id2,
         print('sep, err [km]',to_si(rho[0],0,d_star[0]),to_si(rho[1],0,d_star[0]))
         print('D_star',d_star[0],'+\-',d_star[1])
         print('Delta Gmag',delta_mag)
+        print('RUWE source 1:', ruwe[0][0])
+        print('RUWE source 2:', ruwe[1][0])
         print()
-        yn = input('Does this look good? Hit enter to start the fit, n to exit')
+        yn = input('Does this look good? Hit enter to start the fit, n to exit: ')
         if yn == 'n':
             return None
         else:
             print("Yeehaw let's go")
+
+    if ruwe[0]>1.2 or ruwe[1]>1.2:
+        yn = input('''WARNING: RUWE for one or more of your solutions is greater than 1.2. This indicates 
+            that the source might be an unresolved binary or experiencing acceleration 
+            during the observation.  Orbit fit results may not be trustworthy.  Do you 
+            wish to continue?
+            Hit enter to proceed, n to exit: ''')
+        if yn == 'n':
+            return None
+        else:
+            pass
+
 
     #################### Begin the fit: ####################
 

diff --git a/setup.py b/setup.py
@@ -1,7 +1,7 @@
 from setuptools import setup
 
 setup(name='lofti_gaia',
-      version='0.3',
+      version='0.5',
       description='Orbit fitting with Gaia astrometry',
       url='https://github.com/logan-pearce/lofti_gaia',
       author='Logan Pearce',