import numpy as np

import pandas as pd
from astropy.table import QTable

from astropy import units as u
from astropy import constants as const
from astropy.units import imperial
imperial.enable()

"""
Explanation: Units in Python
The Astropy package includes a powerful framework that allows users to attach units to scalars and arrays, and manipulate/combine these, keeping track of the units.
Astropy has a lot of built-in units. A complete list can be found here.
End of explanation
"""


u.m    # The unit of meters

u.s    # The unit of seconds

u.m / u.s    # combine them into a composite unit

"""
Explanation: Note: because we imported the units package as u, you cannot use u as a variable name.
You can use the units by themselves:
End of explanation
"""


u.m.find_equivalent_units()

"""
Explanation: For any unit you can find all of the built-in units that are equivalent:
End of explanation
"""


time_1 = 0.25 * u.s

time_1

position = np.arange(1,6,1) * u.m    # np.arange(x,y,z) - create an array of numbers between x and y in steps z

position

velocity = position / time_1

velocity

"""
Explanation: The units package is much more useful when you combine it with scalars or arrays to create Quantities
End of explanation
"""


velocity.value

velocity.unit

"""
Explanation: You can access the number and unit part of the Quantity separately:
End of explanation
"""


"The velocity of the first particle is {0:.1f}".format(velocity[0])

"The velocity of the first particle is {0.value:.1f} in the units of {0.unit:s}.".format(velocity[0])

"""
Explanation: This is useful in formatting output:
End of explanation
"""


distance = 100 * u.km

time_2 = distance / velocity

time_2

"""
Explanation: At example problem: How long would it take to go 100 km at velocities in velocity?
End of explanation
"""


time_2.decompose()

"""
Explanation: Notice that the units are a bit strange. We can simplify this using .decompose()
End of explanation
"""


time_2.to(u.h)

velocity.to(u.cm / u.h)

velocity.to(imperial.mi / u.h)

velocity.si                     # quick conversion to SI units

velocity.cgs                    # quick conversion to CGS units

density = 3000 * (u.kg / u.m**3)     # From last week's homework

density.to(u.kg / u.km**3)

"""
Explanation: Unit conversion is really easy!
End of explanation
"""


u.s.find_equivalent_units()

time_3 = 29 * u.day + 7 * u.h + 56 * u.min + 12 * u.s
time_4 = 2 * u.fortnight + 1.33 * u.day

time_3 - time_4

(time_3 - time_4).to(u.min)

"""
Explanation: Units make math with Time easy!
End of explanation
"""


ringo = u.def_unit('Ringos', 3.712 * imperial.yd)

position.to(ringo)

velocity.to(ringo / u.s)

"""
Explanation: You can define your own units
End of explanation
"""


dimless_y = (1 * u.m) / (1 * u.km)

dimless_y

dimless_y.unit

dimless_y.decompose()   # returns the scale of the dimentionless quanity

"""
Explanation: Dimentionless Units
End of explanation
"""


np.log(2 * u.m)

np.log((2 * u.km) / (1 * u.m))

np.log10((2 * u.km) / (1 * u.m))

"""
Explanation: Some math functions only make sense with dimentionless quanities
End of explanation
"""


np.sin(2 * u.m)

np.sin(2 * u.deg)

"""
Explanation: Or they expect the correct type of unit!
End of explanation
"""


np.sin(90)               # not really what I expected

np.sin(90 * u.deg)       # better

"""
Explanation: Using units can save you headaches.

All of the trig functions expect all angles to be in radians. 
If you forget this, it can lead to problems that are hard to debug
End of explanation
"""


planet_table = QTable.read('Planets.csv', format='ascii.csv')

planet_table[0:3]

planet_table['a'].unit = u.AU
planet_table[0:3]

"""
Explanation: QTables vs. DataFrames

Astropy tables are constructed to use units.
The units of the column will be saved with the table

Part II - Advantage Astropy
End of explanation
"""


planet_table['a'].to(u.km)

planet_table['a'].to(u.km)[2]

"""
Explanation: The unit is now part of the structure of the table, and will be saved when the table is saved.
End of explanation
"""


planet_table2 = pd.read_csv('Planets.csv')

planet_table2[0:3]

planet_table3 = QTable.from_pandas(planet_table2)

planet_table3['a'].unit = u.AU

planet_table3[0:3]

planet_table3['a']

planet_table3['a'].to(u.km)

"""
Explanation: You can use  pandas to read in tables and then convert them to QTables
End of explanation
"""


const.G

const.M_sun

"""
Explanation: Constants
The Astropy package also includes a whole bunch of built-in constants to make your life easier.
A complete list of the built-in constants can be found here.
End of explanation
"""


distance = planet_table['a'][0:3]  # Mercury, Venus, Earth

orbit_v = np.sqrt(const.G * const.M_sun / distance)

orbit_v

orbit_v.decompose()

orbit_v.to(u.km/u.s)

orbit_v.to(ringo/u.ms)

"""
Explanation: An Example: The velocity of an object in circular orbit around the Sun is
$$ v=\sqrt{GM_{\odot}\over d} $$
End of explanation
"""


def find_diameter(H,A):
    result = (1329 / np.sqrt(A)) * (10 ** (-0.2 * H))
    return result * u.km

H = 3.34
A = 0.09

asteroid_diameter = find_diameter(H,A)
asteroid_diameter

def find_mass(D):
    p = 3000 * (u.kg / u.m**3)
    result = p * (1/6) * np.pi * D**3
    return result

asteroid_mass = find_mass(asteroid_diameter)
asteroid_mass

asteroid_mass.decompose()

moon_mass = u.def_unit('Lunar\ Masses', 7.34767309e22 * u.kg)

asteroid_mass.to(moon_mass)

"""
Explanation: Functions and Units - (Last week's homework)
End of explanation
"""


