replaced callhorizons with astroquery.jplhorizons

_pp_conf.py

@@ -23,16 +23,25 @@
     print('Module numpy not found. Please install with: pip install numpy')
     sys.exit()
 
+try:
+    from astroquery.jplhorizons import Horizons
+except ImportError:
+    print('Module astroquery not found. Please install with: pip install '
+          'astroquery')
+    sys.exit()
+
+
 # import pipeline-specific modules
 from toolbox import *
 
 # only import if Python3 is used
 if sys.version_info > (3, 0):
     from past.builtins import execfile
 else:
-    warnings.warn(('This is 2018. You should really be using Python 3 by now. '
-                   'PP now has only limited support for Python 2.7. Please '
-                   'install Python 3 (preferentially Anaconda 3).'))
+    raise RuntimeError(
+        'This is 2018. You should really be using Python 3 by now. '
+        'PP now has only limited support for Python 2.7. Please '
+        'install Python 3 (preferentially Anaconda 3).')
 
 
 def setup_diagnostics():

doc/install.rst

@@ -25,9 +25,8 @@ Package Index`_ through `pip`_):
 * `numpy`_
 * `scipy`_
 * `astropy`_
-* `astroquery`_
+* `astroquery`_ (version >= 0.3.9)
 * `matplotlib`_
-* `callhorizons`_
 * `future`_    
 * `skimage`_
 
@@ -81,7 +80,7 @@ Unpack SCAMP and install using::
 
 Install Python modules::
 
-  pip install numpy scipy astropy astroquery matplotlib matplotlib callhorizons 
+  pip install --upgrade --user numpy scipy astropy astroquery matplotlib
 
 Add these lines to the ``.bashrc`` file in your home directory and
 replace ``<path>`` with the actual path to the PP directory::
@@ -120,7 +119,7 @@ Update pip::
 
 Install extra python modules in Anaconda::
 
-  pip install –-user numpy scipy matplotlib callhorizons future astroquery Pillow==2.6.1
+  pip install --upgrade –-user numpy scipy astropy astroquery matplotlib
   pip install astropy
   conda update astropy
 
@@ -249,7 +248,6 @@ and I will take care of implementing your telescope.
 .. _astropy: http://www.astropy.org/
 .. _astroquery: https://github.com/astropy/astroquery
 .. _matplotlib: http://matplotlib.org/
-.. _callhorizons: https://pypi.python.org/pypi/CALLHORIZONS
 .. _future: http://python-future.org/
 .. _skimage: https://scikit-image.org/
 .. _imagemagick: http://www.imagemagick.org/

pp_combine.py

@@ -27,14 +27,14 @@
 import os
 import sys
 import shutil
-import callhorizons
 import logging
 import subprocess
 import argparse
 import shlex
 import time
 from astropy.io import fits
 from past.utils import old_div
+from astroquery.jplhorizons import Horizons
 
 # pipeline-specific modules
 import _pp_conf
@@ -142,10 +142,11 @@ def combine(filenames, obsparam, comoving, targetname,
             # use ephemerides from Horizons if no manual rates are provided
             if manual_rates is None:
                 # call HORIZONS to get target coordinates
-                eph = callhorizons.query(targetname.replace('_', ' '))
-                eph.set_discreteepochs(date)
+                obj = Horizons(targetname.replace('_', ' '), epochs=date,
+                               location=str(obsparam['observatory_code']))
                 try:
-                    n = eph.get_ephemerides(str(obsparam['observatory_code']))
+                    eph = obj.ephemerides()
+                    n = len(eph)
                 except ValueError:
                     print('Target (%s) not an asteroid' % targetname)
                     logging.warning('Target (%s) not an asteroid' % targetname)

pp_distill.py

@@ -4,9 +4,6 @@
                  of select moving or fixed sources
     v1.0: 2016-01-24, mommermiscience@gmail.com
 """
-from __future__ import print_function
-from __future__ import division
-
 # Photometry Pipeline
 # Copyright (C) 2016-2018  Michael Mommert, mommermiscience@gmail.com
 
@@ -25,23 +22,14 @@
 # <http://www.gnu.org/licenses/>.
 
 
-from past.utils import old_div
 import numpy
 import os
 import sys
 import logging
 import argparse
-import time
 import sqlite3
-from astropy.io import fits
-import matplotlib
-matplotlib.use('Agg')
-import matplotlib.pylab as plt
-from scipy.optimize import minimize
-import callhorizons
-from astropy.table import Table
-from astropy.io import ascii
-import scipy.ndimage.interpolation
+from astroquery.jplhorizons import Horizons
+
 try:
     from astroquery.vizier import Vizier
 except ImportError:
@@ -88,7 +76,7 @@ def manual_positions(posfile, catalogs, display=True):
                                                 display=False))
 
         if display:
-            print(old_div(len(objects), len(catalogs)), 'object(s) found')
+            print(len(objects)/len(catalogs), 'object(s) found')
 
         return (list(numpy.hstack(objects)))
 
@@ -129,7 +117,7 @@ def manual_positions(posfile, catalogs, display=True):
                             'dec.deg':  positions[cat_idx]['dec']})
 
     if display:
-        print(old_div(len(objects), len(catalogs)), 'object(s) found')
+        print(len(objects)/len(catalogs), 'object(s) found')
 
     return objects
 
@@ -151,7 +139,7 @@ def pick_controlstar(catalogs, display=True):
                                        'ra.deg', 'dec.deg', 'FLAGS'],
                                    extract_other_catalog=[
                                        'ra.deg', 'dec.deg', 'FLAGS', 'MAG_APER'],
-                                   tolerance=old_div(1., 3600.))
+                                   tolerance=1/3600)
 
     objects = []
     if len(match[0][0]) > 0:
@@ -199,17 +187,19 @@ def moving_primary_target(catalogs, man_targetname, offset, is_asteroid=None,
         if man_targetname is not None:
             targetname = man_targetname.replace('_', ' ')
         for smallbody in [True, False]:
-            eph = callhorizons.query(targetname.replace('_', ' '),
-                                     smallbody=smallbody)
-            #eph = callhorizons.query(targetname, smallbody=False)
-            eph.set_discreteepochs(cat.obstime[0])
+            obj = Horizons(targetname.replace('_', ' '),
+                           id_type={True: 'smallbody',
+                                    False: 'majorbody'}[smallbody],
+                           epochs=cat.obstime[0],
+                           location=obsparam['observatory_code'])
             n = 0
             try:
-                n = eph.get_ephemerides(obsparam['observatory_code'])
+                eph = obj.ephemerides()
+                n = len(eph)
             except ValueError:
                 if display and smallbody is True:
-                    print("'%s' is not an asteroid" % targetname)
-                    logging.warning("'%s' is not an asteroid" %
+                    print("'%s' is not a small body" % targetname)
+                    logging.warning("'%s' is not a small body" %
                                     targetname)
                 if display and smallbody is False:
                     print("'%s' is not a Solar System object" % targetname)
@@ -232,13 +222,14 @@ def moving_primary_target(catalogs, man_targetname, offset, is_asteroid=None,
         if man_targetname is not None:
             targetname = man_targetname.replace('_', ' ')
             cat.obj = targetname
-        eph = callhorizons.query(targetname.replace('_', ' '),
-                                 smallbody=is_asteroid)
-        #eph = callhorizons.query(targetname, smallbody=False)
-        eph.set_discreteepochs(cat.obstime[0])
-
+        obj = Horizons(targetname.replace('_', ' '),
+                       id_type={True: 'smallbody',
+                                False: 'majorbody'}[is_asteroid],
+                       epochs=cat.obstime[0],
+                       location=obsparam['observatory_code'])
         try:
-            n = eph.get_ephemerides(obsparam['observatory_code'])
+            eph = obj.ephemerides()
+            n = len(eph)
         except ValueError:
             # if is_asteroid:
             #     if display and not message_shown:
@@ -257,7 +248,7 @@ def moving_primary_target(catalogs, man_targetname, offset, is_asteroid=None,
         if n is None or n == 0:
             logging.warning('WARNING: No position from Horizons! ' +
                             'Name (%s) correct?' % cat.obj.replace('_', ' '))
-            logging.warning('HORIZONS call: %s' % eph.url)
+            logging.warning('HORIZONS call: %s' % obj.uri)
             if display and not message_shown:
                 print('  no Horizons data for %s ' % cat.obj.replace('_', ' '))
                 message_shown = True
@@ -266,11 +257,11 @@ def moving_primary_target(catalogs, man_targetname, offset, is_asteroid=None,
             objects.append({'ident': eph[0]['targetname'].replace(" ", "_"),
                             'obsdate.jd': cat.obstime[0],
                             'cat_idx': cat_idx,
-                            'ra.deg': eph[0]['RA']-old_div(offset[0], 3600.),
-                            'dec.deg': eph[0]['DEC']-old_div(offset[1], 3600.)})
+                            'ra.deg': eph[0]['RA']-offset[0]/3600,
+                            'dec.deg': eph[0]['DEC']-offset[1]/3600})
             logging.info('Successfully grabbed Horizons position for %s ' %
                          cat.obj.replace('_', ' '))
-            logging.info('HORIZONS call: %s' % eph.url)
+            logging.info('HORIZONS call: %s' % obj.uri)
             if display and not message_shown:
                 print(cat.obj.replace('_', ' '), "identified")
                 message_shown = True
@@ -305,7 +296,7 @@ def fixed_targets(fixed_targets_file, catalogs, display=True):
                             'dec.deg': obj['dec']})
 
     if display:
-        print(old_div(len(objects), len(catalogs)), 'targets read')
+        print(len(objects)/len(catalogs), 'targets read')
 
     return objects
 
@@ -361,7 +352,7 @@ def serendipitous_variablestars(catalogs, display=True):
                             'dec.deg': star['DEJ2000']})
 
     if display:
-        print(old_div(len(objects), len(catalogs)), 'variable stars found')
+        print(len(objects)/len(catalogs), 'variable stars found')
 
     return objects
 
@@ -446,7 +437,7 @@ def serendipitous_asteroids(catalogs, display=True):
                       'target pool').format(obj['name']))
 
     if display:
-        print(old_div(len(objects), len(catalogs)), 'asteroids found')
+        print(len(objects)/len(catalogs), 'asteroids found')
 
     return objects
 
@@ -522,7 +513,7 @@ def distill(catalogs, man_targetname, offset, fixed_targets_file, posfile,
         print('#-----------------------')
 
     if display:
-        print(old_div(len(objects), len(catalogs)),
+        print(len(objects)/len(catalogs),
               'potential target(s) per frame identified.')
 
     # extract source data for identified targets

pp_photometry.py

@@ -5,9 +5,6 @@
 
     v1.0: 2015-12-30, mommermiscience@gmail.com
 """
-from __future__ import print_function
-from __future__ import division
-
 # Photometry Pipeline
 # Copyright (C) 2016-2018  Michael Mommert, mommermiscience@gmail.com
 
@@ -26,20 +23,15 @@
 # <http://www.gnu.org/licenses/>.
 
 
-from past.utils import old_div
 import numpy
 import os
 import sys
-import subprocess
 import logging
 import argparse
-import time
-from copy import deepcopy
 from astropy.io import fits
 import matplotlib
 matplotlib.use('Agg')
-import matplotlib.pylab as plt
-import callhorizons
+from astroquery.jplhorizons import Horizons
 
 # only import if Python3 is used
 if sys.version_info > (3, 0):
@@ -133,19 +125,21 @@ def curve_of_growth_analysis(filenames, parameters,
                 date = hdu[0].header['MIDTIMJD']
 
             # call HORIZONS to get target coordinates
-            eph = callhorizons.query(targetname.replace('_', ' '))
-            eph.set_discreteepochs(date)
+            obj = Horizons(targetname.replace('_', ' '),
+                           epochs=date,
+                           location=str(obsparam['observatory_code']))
             try:
-                n = eph.get_ephemerides(str(obsparam['observatory_code']))
+                eph = obj.ephemerides()
+                n = len(eph)
             except ValueError:
-                print('Target (%s) not an asteroid' % targetname)
-                logging.warning('Target (%s) not an asteroid' % targetname)
+                print('Target (%s) not a small body' % targetname)
+                logging.warning('Target (%s) not a small body' % targetname)
                 n = None
 
             if n is None or n == 0:
                 logging.warning('WARNING: No position from Horizons!' +
                                 'Name (%s) correct?' % targetname)
-                logging.warning('HORIZONS call: %s' % eph.url)
+                logging.warning('HORIZONS call: %s' % obj.uri)
                 logging.info('proceeding with background sources analysis')
                 parameters['background_only'] = True
             else:
@@ -164,22 +158,21 @@ def curve_of_growth_analysis(filenames, parameters,
             residuals = numpy.sqrt((data['ra.deg']-target_ra)**2 +
                                    (data['dec.deg']-target_dec)**2)
             target_idx = numpy.argmin(residuals)
-            if residuals[target_idx] > old_div(_pp_conf.pos_epsilon, 3600.):
+            if residuals[target_idx] > _pp_conf.pos_epsilon/3600:
                 logging.warning(('WARNING: frame %s, large residual to ' +
                                  'HORIZONS position of %s: %f arcsec; ' +
                                  'ignore this frame') %
                                 (filename, targetname,
                                  residuals[numpy.argmin(residuals)]*3600.))
             else:
-                target_flux.append(old_div(data[target_idx]['FLUX_' +
-                                                            _pp_conf.photmode],
-                                           max(data[target_idx]['FLUX_' +
-                                                                _pp_conf.photmode])))
+                target_flux.append(data[target_idx]['FLUX_'+_pp_conf.photmode] /
+                                   max(data[target_idx][
+                                       'FLUX_'+_pp_conf.photmode]))
                 target_snr.append(
                     data[target_idx]['FLUX_'+_pp_conf.photmode] /
                     data[target_idx]['FLUXERR_'+_pp_conf.photmode] /
-                    max(old_div(data[target_idx]['FLUX_'+_pp_conf.photmode],
-                                data[target_idx]['FLUXERR_'+_pp_conf.photmode])))
+                    max(data[target_idx]['FLUX_'+_pp_conf.photmode] /
+                        data[target_idx]['FLUXERR_'+_pp_conf.photmode]))
                 n_target_identified += 1
 
         # extract background source fluxes and snrs
@@ -194,13 +187,12 @@ def curve_of_growth_analysis(filenames, parameters,
                     continue
 
                 # create growth curve
-                background_flux.append(old_div(src['FLUX_'+_pp_conf.photmode],
-                                               max(src['FLUX_'+_pp_conf.photmode])))
+                background_flux.append(src['FLUX_'+_pp_conf.photmode] /
+                                       max(src['FLUX_'+_pp_conf.photmode]))
                 background_snr.append(src['FLUX_'+_pp_conf.photmode] /
                                       src['FLUXERR_'+_pp_conf.photmode] /
-                                      max(old_div(src['FLUX_' +
-                                                      _pp_conf.photmode],
-                                                  src['FLUXERR_'+_pp_conf.photmode])))
+                                      max(src['FLUX_'+_pp_conf.photmode] /
+                                          src['FLUXERR_'+_pp_conf.photmode]))
 
     # investigate curve-of-growth
 
@@ -211,7 +203,7 @@ def curve_of_growth_analysis(filenames, parameters,
     n_target = len(target_flux)
     if n_target > 0:
         target_flux = (numpy.median(target_flux, axis=0),
-                       old_div(numpy.std(target_flux, axis=0), numpy.sqrt(n_target)))
+                       numpy.std(target_flux, axis=0)/numpy.sqrt(n_target))
         target_snr = numpy.median(target_snr, axis=0)
     else:
         target_flux = (numpy.zeros(len(aprads)), numpy.zeros(len(aprads)))
@@ -220,8 +212,8 @@ def curve_of_growth_analysis(filenames, parameters,
     n_background = len(background_flux)
     if n_background > 0:
         background_flux = (numpy.median(background_flux, axis=0),
-                           old_div(numpy.std(background_flux, axis=0),
-                                   numpy.sqrt(n_background)))
+                           numpy.std(background_flux, axis=0) /
+                           numpy.sqrt(n_background))
         background_snr = numpy.median(background_snr, axis=0)
     else:
         background_flux = (numpy.zeros(len(aprads)), numpy.zeros(len(aprads)))