<commit_before>import unittest
import datetime
import os

import numpy as np

import pytz

import astropy.time
import astropy.coordinates
import astropy.io
import astropy.utils.data
import astropy.units as u

from desisurvey import utils, config


class TestUtils(unittest.TestCase):

    def setUp(self):
        # Configure a CSV reader for the Horizons output format.
        csv_reader = astropy.io.ascii.Csv()
        csv_reader.header.comment = r'[^ ]'
        csv_reader.data.start_line = 35
        csv_reader.data.end_line = 203
        # Read moon ephemerides for the first week of 2020.
        path = astropy.utils.data._find_pkg_data_path(
            os.path.join('data', 'horizons_2020_week1_moon.csv'),
            package='desisurvey')
        self.table = csv_reader.read(path)
        # Horizons CSV file has a trailing comma on each line.
        self.table.remove_column('col10')
        # Use more convenient column names.
        names = ('date', 'jd', 'sun', 'moon', 'ra', 'dec',
                 'az', 'alt', 'lst', 'frac')
        for old_name, new_name in zip(self.table.colnames, names):
            self.table[old_name].name = new_name

    def test_get_overhead(self):
        """Sanity checks on overhead time calculations"""
        c = config.Configuration()
        tro = c.readout_time()
        p0 = astropy.coordinates.SkyCoord(ra=300 * u.deg, dec=10 * u.deg)
        # Move with no readout and no slew only has focus overhead.
        self.assertEqual(utils.get_overhead_time(None, p0, tro),
                         c.focus_time())
        self.assertEqual(utils.get_overhead_time(p0, p0, 2 * tro),
                         c.focus_time())
        # Move with readout and no slew has focus and overhead in parallel.
        self.assertEqual(utils.get_overhead_time(None, p0),
                         max(c.focus_time(), c.readout_time()))
        self.assertEqual(utils.get_overhead_time(p0, p0, 0.2 * tro),
                         max(c.focus_time(), 0.8 * c.readout_time()))
        # Overhead with slew same when dRA == dDEC.
        for delta in (1, 45, 70):
            ra = p0.ra + [+delta, -delta, 0, 0] * u.deg
            dec = p0.dec + [0, 0, +delta, -delta] * u.deg
            p1 = astropy.coordinates.SkyCoord(ra=ra, dec=dec)
            dt = utils.get_overhead_time(p0, p1)
            self.assertTrue(dt.shape == (4,))
            self.assertEqual(dt[0], dt[1])
            self.assertEqual(dt[0], dt[2])
            self.assertEqual(dt[0], dt[3])

    def test_get_observer_to_sky(self):
        """Check (alt,az) -> (ra,dec) against JPL Horizons"""
        ra, dec = self.table['ra'], self.table['dec']
        alt, az = self.table['alt'], self.table['az']
        when = astropy.time.Time(self.table['jd'], format='jd')
        obs = utils.get_observer(when, alt * u.deg, az * u.deg)
        obs_sky = obs.transform_to(astropy.coordinates.ICRS)
        true_sky = astropy.coordinates.ICRS(ra=ra * u.deg, dec=dec * u.deg)
        sep = true_sky.separation(obs_sky).to(u.arcsec).value
        # Opening angle between true and observed (ra,dec) unit vectors
        # must be within 30 arcsec.
        self.assertTrue(np.max(np.fabs(sep)) < 30)

    def test_get_observer_from_sky(self):
        """Check (alt,az) -> (ra,dec) against JPL Horizons"""
        ra, dec = self.table['ra'], self.table['dec']
        alt, az = self.table['alt'], self.table['az']
        when = astropy.time.Time(self.table['jd'], format='jd')
        sky = astropy.coordinates.ICRS(ra=ra * u.deg, dec=dec * u.deg)
        true_altaz = utils.get_observer(when, alt * u.deg, az * u.deg)
        obs_altaz = sky.transform_to(utils.get_observer(when))
        sep = true_altaz.separation(obs_altaz).to(u.arcsec).value
        # Opening angle between true and observed (alt,az) unit vectors
        # must be within 30 arcsec.
        self.assertTrue(np.max(np.fabs(sep)) < 30)

    def test_get_observer_args(self):
        """Must provide both alt and az to get_observer()"""
        t = astropy.time.Time('2020-01-01')
        with self.assertRaises(ValueError):
            utils.get_observer(t, alt=0 * u.deg)
        with self.assertRaises(ValueError):
            utils.get_observer(t, az=0 * u.deg)

    def test_get_observer_units(self):
        """Alt,az must have angular units for get_observer()"""
        t = astropy.time.Time('2020-01-01')
        with self.assertRaises(TypeError):
            utils.get_observer(t, alt=1, az=1)
        with self.assertRaises(u.UnitsError):
            utils.get_observer(t, alt=1 * u.m, az=1 * u.m)

    def test_get_observer_time(self):
        """Must pass arg convertible to astropy Time to get_observer()"""
        utils.get_observer('2020-01-01')
        utils.get_observer(datetime.datetime(2020, 1, 1, 0, 0))
        with self.assertRaises(ValueError):
            utils.get_observer(datetime.date(2020, 1, 1))
        with self.assertRaises(ValueError):
            utils.get_observer(50000.)

    def test_zenith_airmass(self):
        """Airmass values monotically increase with zenith angle"""
        Z = np.arange(90) * np.pi / 180.
        X = utils.zenith_angle_to_airmass(Z)
        self.assertEqual(Z.shape, X.shape)
        self.assertTrue(np.all(np.diff(X) > 0))

    def test_zenith_airmass(self):
        """Zenith angles can have units"""
        Z1 = np.arange(90) * u.deg
        Z2 = Z1.to(u.rad)
        X1 = utils.zenith_angle_to_airmass(Z1)
        X2 = utils.zenith_angle_to_airmass(Z2)
        self.assertTrue(np.allclose(X1, X2))

    def test_airmass_scalar(self):
        """Scalar input returns scalar output"""
        X = utils.zenith_angle_to_airmass(0.)
        self.assertEqual(X.shape, ())

    def test_airmass_clip(self):
        """Zenith angles > 90 deg are clipped"""
        self.assertAlmostEqual(
            utils.zenith_angle_to_airmass(90 * u.deg),
            utils.zenith_angle_to_airmass(100 * u.deg))

    def test_get_airmass_lowest(self):
        """The lowest airmass occurs when dec=latitude"""
        t = astropy.time.Time('2020-01-01')
        ra = np.arange(360) * u.deg
        dec = utils.get_location().latitude
        Xinv = 1 / utils.get_airmass(t, ra, dec)
        self.assertTrue(np.max(Xinv) > 0.999)
        dec += 30 * u.deg
        Xinv = 1 / utils.get_airmass(t, ra, dec)
        self.assertTrue(np.max(Xinv) < 0.9)

    def test_get_airmass_always_visible(self):
        """An object at DEC=70d is always visible"""
        t = astropy.time.Time('2020-01-01')
        ra = np.arange(360) * u.deg
        dec = 70 * u.deg
        Xinv = 1 / utils.get_airmass(t, ra, dec)
        self.assertTrue(np.all(0.2 < Xinv) and np.all(Xinv < 0.8))

    def test_get_location(self):<commit_msg>Attemp to fix astropy 1.2 incompatibility<commit_after>import unittest
import datetime
import os

import numpy as np

import pytz

import astropy.time
import astropy.coordinates
import astropy.io
import astropy.utils.data
import astropy.units as u

from desisurvey import utils, config


class TestUtils(unittest.TestCase):

    def setUp(self):
        # Configure a CSV reader for the Horizons output format.
        csv_reader = astropy.io.ascii.Csv()
        csv_reader.header.comment = r'[^ ]'
        csv_reader.data.start_line = 35
        csv_reader.data.end_line = 203
        # Read moon ephemerides for the first week of 2020.
        path = astropy.utils.data._find_pkg_data_path(
            os.path.join('data', 'horizons_2020_week1_moon.csv'),
            package='desisurvey')
        self.table = csv_reader.read(path)
        # Horizons CSV file has a trailing comma on each line.
        self.table.remove_column('col10')
        # Use more convenient column names.
        names = ('date', 'jd', 'sun', 'moon', 'ra', 'dec',
                 'az', 'alt', 'lst', 'frac')
        for old_name, new_name in zip(self.table.colnames, names):
            self.table[old_name].name = new_name

    def test_get_overhead(self):
        """Sanity checks on overhead time calculations"""
        c = config.Configuration()
        tro = c.readout_time()
        p0 = astropy.coordinates.SkyCoord(ra=300 * u.deg, dec=10 * u.deg)
        # Move with no readout and no slew only has focus overhead.
        self.assertEqual(utils.get_overhead_time(None, p0, tro),
                         c.focus_time())
        self.assertEqual(utils.get_overhead_time(p0, p0, 2 * tro),
                         c.focus_time())
        # Move with readout and no slew has focus and overhead in parallel.
        self.assertEqual(utils.get_overhead_time(None, p0),
                         max(c.focus_time(), c.readout_time()))
        self.assertEqual(utils.get_overhead_time(p0, p0, 0.2 * tro),
                         max(c.focus_time(), 0.8 * c.readout_time()))
        # Overhead with slew same when dRA == dDEC.
        for delta in (1, 45, 70):
            ra = p0.ra + [+delta, -delta, 0, 0] * u.deg
            dec = p0.dec + [0, 0, +delta, -delta] * u.deg
            p1 = astropy.coordinates.SkyCoord(ra=ra, dec=dec)
            dt = utils.get_overhead_time(p0, p1)
            self.assertTrue(dt.shape == (4,))
            self.assertEqual(dt[0], dt[1])
            self.assertEqual(dt[0], dt[2])
            self.assertEqual(dt[0], dt[3])

    def test_get_observer_to_sky(self):
        """Check (alt,az) -> (ra,dec) against JPL Horizons"""
        ra, dec = self.table['ra'], self.table['dec']
        alt, az = self.table['alt'], self.table['az']
        when = astropy.time.Time(self.table['jd'], format='jd')
        obs = utils.get_observer(when, alt * u.deg, az * u.deg)
        obs_sky = obs.transform_to(astropy.coordinates.ICRS)
        true_sky = astropy.coordinates.ICRS(ra=ra * u.deg, dec=dec * u.deg)
        sep = true_sky.separation(obs_sky).to(u.arcsec).value
        # Opening angle between true and observed (ra,dec) unit vectors
        # must be within 30 arcsec.
        self.assertTrue(np.max(np.fabs(sep)) < 30)

    def test_get_observer_from_sky(self):
        """Check (alt,az) -> (ra,dec) against JPL Horizons"""
        ra, dec = self.table['ra'], self.table['dec']
        alt, az = self.table['alt'], self.table['az']
        when = astropy.time.Time(self.table['jd'], format='jd')
        sky = astropy.coordinates.ICRS(ra=ra * u.deg, dec=dec * u.deg)
        true_altaz = utils.get_observer(when, alt * u.deg, az * u.deg)
        obs_altaz = sky.transform_to(utils.get_observer(when))
        sep = true_altaz.separation(obs_altaz).to(u.arcsec).value
        # Opening angle between true and observed (alt,az) unit vectors
        # must be within 30 arcsec.
        self.assertTrue(np.max(np.fabs(sep)) < 30)

    def test_get_observer_args(self):
        """Must provide both alt and az to get_observer()"""
        t = astropy.time.Time('2020-01-01')
        with self.assertRaises(ValueError):
            utils.get_observer(t, alt=0 * u.deg)
        with self.assertRaises(ValueError):
            utils.get_observer(t, az=0 * u.deg)

    def test_get_observer_units(self):
        """Alt,az must have angular units for get_observer()"""
        t = astropy.time.Time('2020-01-01')
        with self.assertRaises(TypeError):
            utils.get_observer(t, alt=1, az=1)
        with self.assertRaises(u.UnitsError):
            utils.get_observer(t, alt=1 * u.m, az=1 * u.m)

    def test_get_observer_time(self):
        """Must pass arg convertible to astropy Time to get_observer()"""
        utils.get_observer('2020-01-01')
        utils.get_observer(datetime.datetime(2020, 1, 1, 0, 0))
        with self.assertRaises(ValueError):
            utils.get_observer(datetime.date(2020, 1, 1))
        with self.assertRaises(ValueError):
            utils.get_observer(50000.)

    def test_zenith_airmass(self):
        """Airmass values monotically increase with zenith angle"""
        Z = np.arange(90) * np.pi / 180.
        X = utils.zenith_angle_to_airmass(Z)
        self.assertEqual(Z.shape, X.shape)
        self.assertTrue(np.all(np.diff(X) > 0))

    def test_zenith_airmass(self):
        """Zenith angles can have units"""
        Z1 = np.arange(90) * u.deg
        Z2 = Z1.to(u.rad)
        X1 = utils.zenith_angle_to_airmass(Z1)
        X2 = utils.zenith_angle_to_airmass(Z2)
        self.assertTrue(np.allclose(X1, X2))

    def test_airmass_scalar(self):
        """Scalar input returns scalar output"""
        X = utils.zenith_angle_to_airmass(0.)
        self.assertEqual(X.shape, ())

    def test_airmass_clip(self):
        """Zenith angles > 90 deg are clipped"""
        self.assertAlmostEqual(
            utils.zenith_angle_to_airmass(90 * u.deg),
            utils.zenith_angle_to_airmass(100 * u.deg))

    def test_get_airmass_lowest(self):
        """The lowest airmass occurs when dec=latitude"""
        t = astropy.time.Time('2020-01-01')
        ra = np.arange(360) * u.deg
        dec = utils.get_location().latitude
        Xinv = 1 / utils.get_airmass(t, ra, dec)
        self.assertTrue(np.max(Xinv) > 0.999)
        dec = dec + 30 * u.deg
        Xinv = 1 / utils.get_airmass(t, ra, dec)
        self.assertTrue(np.max(Xinv) < 0.9)

    def test_get_airmass_always_visible(self):
        """An object at DEC=70d is always visible"""
        t = astropy.time.Time('2020-01-01')
        ra = np.arange(360) * u.deg
        dec = 70 * u.deg
        Xinv = 1 / utils.get_airmass(t, ra, dec)
        self.assertTrue(np.all(0.2 < Xinv) and np.all(Xinv < 0.8))

    def test_get_location(self):
{'fraction_non_alphanumeric': 0.1313355530891602, 'fraction_numerical': 0.02634110247872734, 'mean_word_length': 2.4330708661417324, 'pattern_counts': {'":': 0, '<': 13, '<?xml version=': 0, '>': 13, 'https://': 0, 'lorem ipsum': 0, 'www.': 0, 'xml': 0}, 'pii_count': 0, 'substrings_counts': 0, 'word_list_counts': {'cursed_substrings.json': 0, 'profanity_word_list.json': 0, 'sexual_word_list.json': 0, 'zh_pornsignals.json': 0}}
{'dir': 'git-commits-cleaned', 'id': '2992684', 'n_tokens_mistral': 4801, 'n_tokens_neox': 4521, 'n_words': 991}