Skip to content

Commit

Permalink
precommit linting/ruff formatter
Browse files Browse the repository at this point in the history
  • Loading branch information
changgoo committed Oct 27, 2023
1 parent 5e78b2e commit f641d09
Show file tree
Hide file tree
Showing 8 changed files with 1,496 additions and 1,009 deletions.
4 changes: 3 additions & 1 deletion astro_tigress/__init__.py
Original file line number Diff line number Diff line change
Expand Up @@ -3,4 +3,6 @@
__license__ = "MIT"
__description__ = "Python packages for the TIGRESS public data release"

from astro_tigress.tigress_read import *
from astro_tigress.tigress_read import Model, DataMHD, DataRad, DataChem, DataCO

__all__ = ["Model", "DataMHD", "DataRad", "DataChem", "DataCO"]
188 changes: 101 additions & 87 deletions astro_tigress/const.py
Original file line number Diff line number Diff line change
Expand Up @@ -11,225 +11,239 @@
# CGS PHYSICAL CONSTANTS
# %&%&%&%&%&%&%&%&%&%&%&%&%&%&%&%&%&%&%&%&%&%&%&%&%&%&%&%&

c = 2.99792458e10 # speed of light CGS
h = 6.6260755e-27 # Planck's constant CGS
g = 6.67259e-8 # Grav const CGS
kb = 1.380658e-16 # Boltzmann's const CGS
a = 7.56591e-15 # Radiation constant CGS
sb = 5.67051e-5 # sigma (stefan-boltzmann const) CGS
qe = 4.803206e-10 # Charge of electron CGS
ev = 1.60217733e-12 # Electron volt CGS
na = 6.0221367e23 # Avagadro's Number
me = 9.1093897e-28 # electron mass CGS
mp = 1.6726231e-24 # proton mass CGS
mn = 1.674929e-24 # neutron mass CGS
mh = 1.673534e-24 # hydrogen mass CGS
amu = 1.6605402e-24 # atomic mass unit CGS
c = 2.99792458e10 # speed of light CGS
h = 6.6260755e-27 # Planck's constant CGS
g = 6.67259e-8 # Grav const CGS
kb = 1.380658e-16 # Boltzmann's const CGS
a = 7.56591e-15 # Radiation constant CGS
sb = 5.67051e-5 # sigma (stefan-boltzmann const) CGS
qe = 4.803206e-10 # Charge of electron CGS
ev = 1.60217733e-12 # Electron volt CGS
na = 6.0221367e23 # Avagadro's Number
me = 9.1093897e-28 # electron mass CGS
mp = 1.6726231e-24 # proton mass CGS
mn = 1.674929e-24 # neutron mass CGS
mh = 1.673534e-24 # hydrogen mass CGS
amu = 1.6605402e-24 # atomic mass unit CGS

# %&%&%&%&%&%&%&%&%&%&%&%&%&%&%&%&%&%&%&%&%&%&%&%&%&%&%&%&
# ASTRONOMICAL CONSTANTS
# %&%&%&%&%&%&%&%&%&%&%&%&%&%&%&%&%&%&%&%&%&%&%&%&%&%&%&%&

# GENERAL
au = 1.496e13 # astronomical unit CGS
pc = 3.0857e18 # parsec CGS
yr = 3.155815e7 # sidereal year CGS
Myr = 1.0e6 * yr # million years CGS
ms = 1.98900e+33 # solar mass CGS
mj = 1.8986e30 # jupiter mass CGS
aj = 5.202545 * au # jupiter semi-major axis CGS
rs = 6.9599e10 # sun's radius CGS
ls = 3.839e33 # sun's luminosity CGS
mm = 7.35000e+25 # moon mass CGS
mer = 5.97400e+27 # earth mass CGS
rer = 6.378e8 # earth's radius CGS
mmsn_er = 1700.0 # g/cm^2, MMSN at 1AU
medd = 3.60271e+34 # Eddington mass CGS
km = 1.0e5 # km in cm
au = 1.496e13 # astronomical unit CGS
pc = 3.0857e18 # parsec CGS
yr = 3.155815e7 # sidereal year CGS
Myr = 1.0e6 * yr # million years CGS
ms = 1.98900e33 # solar mass CGS
mj = 1.8986e30 # jupiter mass CGS
aj = 5.202545 * au # jupiter semi-major axis CGS
rs = 6.9599e10 # sun's radius CGS
ls = 3.839e33 # sun's luminosity CGS
mm = 7.35000e25 # moon mass CGS
mer = 5.97400e27 # earth mass CGS
rer = 6.378e8 # earth's radius CGS
mmsn_er = 1700.0 # g/cm^2, MMSN at 1AU
medd = 3.60271e34 # Eddington mass CGS
km = 1.0e5 # km in cm

# RADIO SPECIFIC
jy = 1.e-23 # Jansky CGS
restfreq_hi = 1420405751.786 # 21cm transition (Hz)
restfreq_co = 115271201800. # CO J=1-0 (Hz)
cm2perkkms_hi = 1.823e18 # HI column per intensity (thin)
jy = 1.0e-23 # Jansky CGS
restfreq_hi = 1420405751.786 # 21cm transition (Hz)
restfreq_co = 115271201800.0 # CO J=1-0 (Hz)
cm2perkkms_hi = 1.823e18 # HI column per intensity (thin)

# OTHER WAVELENGTHS
ksun = 3.28 # abs K mag of sun (Binney & Merrifield)
ksun = 3.28 # abs K mag of sun (Binney & Merrifield)

# %&%&%&%&%&%&%&%&%&%&%&%&%&%&%&%&%&%&%&%&%&%&%&%&%&%&%&%&
# GEOMETRY
# %&%&%&%&%&%&%&%&%&%&%&%&%&%&%&%&%&%&%&%&%&%&%&%&%&%&%&%&

srdeg = 3283.0 # Degrees squared in a steradian
dtor = 0.0174532925199 # Degrees per radian
pi = 3.14159265359 # Pi
srdeg = 3283.0 # Degrees squared in a steradian
dtor = 0.0174532925199 # Degrees per radian
pi = 3.14159265359 # Pi

# %&%&%&%&%&%&%&%&%&%&%&%&%&%&%&%&%&%&%&%&%&%&%&%&%&%&%&%&
# TELESCOPE STUFF
# %&%&%&%&%&%&%&%&%&%&%&%&%&%&%&%&%&%&%&%&%&%&%&%&%&%&%&%&


def hpbw(lam_cm=None, freq_hz=None, diam_m=None, units=True):
"""
Return the half-power beam width for a telescope.
"""
if (lam_cm == None and freq_hz == None) or \
diam_m == None:
if (lam_cm is None and freq_hz is None) or diam_m is None:
hpbw = None
else:
if lam_cm == None:
lam_cm = c/freq_hz
hpbw = 1.2 * lam_cm / (diam_m*1e2)
if lam_cm is None:
lam_cm = c / freq_hz
hpbw = 1.2 * lam_cm / (diam_m * 1e2)

if units == True:
return {"units":"sr",
"value":hpbw}
if units:
return {"units": "sr", "value": hpbw}
else:
return hpbw


# %&%&%&%&%&%&%&%&%&%&%&%&%&%&%&%&%&%&%&%&%&%&%&%&%&%&%&%&
# PLANCK FUNCTION CALCULATIONS
# %&%&%&%&%&%&%&%&%&%&%&%&%&%&%&%&%&%&%&%&%&%&%&%&%&%&%&%&


def bb_inu(temp_k=None, lam_cm=None, units=True):
"""
Planck function. Returns f_nu given T and lambda.
"""
nu = c / lam_cm
fnu = 2.0*h*nu**3/c**2 /(np.exp(h*nu/(kb*temp_k))-1.0)
fnu = 2.0 * h * nu**3 / c**2 / (np.exp(h * nu / (kb * temp_k)) - 1.0)

if units == True:
return {"units":"erg/s/cm^2/sr/Hz",
"value":fnu}
if units:
return {"units": "erg/s/cm^2/sr/Hz", "value": fnu}
else:
return fnu


def bb_ilam(temp_k=None, lam_cm=None, units=True):
"""
Planck function. Returns f_lambda given T and lambda.
"""
flam = 2.0*h*c**2/lam_cm**5 / (np.exp(h*c/(lam_cm*kb*temp_k))-1.0)
flam = 2.0 * h * c**2 / lam_cm**5 / (np.exp(h * c / (lam_cm * kb * temp_k)) - 1.0)

if units == True:
return {"units":"erg/s/cm^2/sr/cm",
"value":flam}
if units:
return {"units": "erg/s/cm^2/sr/cm", "value": flam}
else:
return flam


### NOT FINISHED - FIX ###


def peak_inu(temp_k=None, units=True):
"""
Wavelength for peak of F_nu given T. NOT FUNCTIONAL YET.
"""
peak_lam = 2.821439372/kb/temp_k/h
peak_nu = c / peak_lam
peak_lam = 2.821439372 / kb / temp_k / h
# peak_nu = c / peak_lam

if units == True:
return {"units":"cm",
"value":peak_lam}
if units:
return {"units": "cm", "value": peak_lam}
else:
return peak_lam


### NOT FINISHED - FIX ###


def peak_ilam(temp_k=None, units=True):
"""
Wavelength for peak of F_lambda given T. NOT FUNCTIONAL YET.
"""
peak_lam = h*c/kb/temp_k/4.965114231
if units == True:
return {"units":"cm",
"value":peak_lam}
peak_lam = h * c / kb / temp_k / 4.965114231
if units:
return {"units": "cm", "value": peak_lam}
else:
return peak_lam


# %&%&%&%&%&%&%&%&%&%&%&%&%&%&%&%&%&%&%&%&%&%&%&%&%&%&%&%&
# TRANSLATE HMS/DMS STRINGS BAND AND FORTH FROM DEC. DEG
# %&%&%&%&%&%&%&%&%&%&%&%&%&%&%&%&%&%&%&%&%&%&%&%&%&%&%&%&

# CASA's QA tool should do this but sucks


def hms(val):
"""
Decimal value to 3-element tuple of hour, min, sec. Just calls dms/15
"""
return dms(val/15.0)
return dms(val / 15.0)


def dms(val):
"""
Decimal value to 3-element tuple of deg, min, sec.
"""
if type(val) == type(np.array([1,2,3])):
if isinstance(val, type(np.array([1, 2, 3]))):
val_abs = np.abs(val)
d = np.floor(val_abs)
m = np.floor((val_abs - d)*60.0)
s = ((val_abs - d - m/60.0))*3600.0
d *= +2.0*(val >= 0.0) - 1.0
m = np.floor((val_abs - d) * 60.0)
s = (val_abs - d - m / 60.0) * 3600.0
d *= +2.0 * (val >= 0.0) - 1.0
else:
val_abs = abs(val)
d = math.floor(val_abs)
m = math.floor((val_abs - d)*60.0)
s = ((val_abs - d - m/60.0))*3600.0
d *= +2.0*(val >= 0.0) - 1.0
m = math.floor((val_abs - d) * 60.0)
s = (val_abs - d - m / 60.0) * 3600.0
d *= +2.0 * (val >= 0.0) - 1.0
return (d, m, s)


def ten(vec, hms=False):
"""
Convert 3-element vector to decimal degrees. Use HMS for hours.
"""
if type(val) == type(np.array([1,2,3])):
dec = (np.abs(vec[0])+vec[1]/60.+vec[2]/3600.)
dec *= +2.0*(vec[0] >= 0.0) - 1.0
if isinstance(vec, type(np.array([1, 2, 3]))):
dec = np.abs(vec[0]) + vec[1] / 60.0 + vec[2] / 3600.0
dec *= +2.0 * (vec[0] >= 0.0) - 1.0
else:
dec = (abs(vec[0])+vec[1]/60.+vec[2]/3600.)
dec *= +2.0*(vec[0] >= 0.0) - 1.0
if hms == True:
return 15.0*dec
dec = abs(vec[0]) + vec[1] / 60.0 + vec[2] / 3600.0
dec *= +2.0 * (vec[0] >= 0.0) - 1.0
if hms:
return 15.0 * dec
else:
return dec

def string_sixty(inp_str,hms=False):

def string_sixty(inp_str, hms=False):
"""
Convert an HMS/DMS style string to a numeric vector. Again, QA
should do this.
"""
if hms==True:
if hms:
pass
else:
pass

def sixty_string(inp_val,hms=False,colons=False):

def sixty_string(inp_val, hms=False, colons=False):
"""
Convert a numeric vector to a string. Again, QA should do this.
"""
if hms==True:
if hms:
out_str = "%02dh%02dm%05.3fs" % (inp_val[0], inp_val[1], inp_val[2])
elif colons==True:
elif colons:
out_str = "%+02d:%02d:%05.3f" % (inp_val[0], inp_val[1], inp_val[2])
else:
out_str = "%+02dd%02dm%05.3fs" % (inp_val[0], inp_val[1], inp_val[2])
return out_str


# %&%&%&%&%&%&%&%&%&%&%&%&%&%&%&%&%&%&%&%&%&%&%&%&%&%&%&%&
# LAPLACE COEFFICIENT AND ITS DIRIVATIVES CALCULATIONS
# %&%&%&%&%&%&%&%&%&%&%&%&%&%&%&%&%&%&%&%&%&%&%&%&%&%&%&%&


def binte(phi, n, s, alpha):
"""the integration of Laplace coefficient"""
return np.cos(n*phi) / (1 - 2*alpha*np.cos(phi) + alpha**2)**s
return np.cos(n * phi) / (1 - 2 * alpha * np.cos(phi) + alpha**2) ** s


def laplace_b(n, s, alpha):
"""the Laplace coefficient b^n_s(alpha)"""
return scipy.integrate.quad(binte, 0., 2.*math.pi, (n, s, alpha)
)[0] / math.pi
return scipy.integrate.quad(binte, 0.0, 2.0 * math.pi, (n, s, alpha))[0] / math.pi


def dlaplace_b_dalpha(n, s, alpha):
def deri(alpha0, n0, s0):
return laplace_b(n0, s0, alpha0)
return scipy.misc.derivative(deri, alpha, dx=1.e-5, args=(n, s))

return scipy.misc.derivative(deri, alpha, dx=1.0e-5, args=(n, s))


# %&%&%&%&%&%&%&%&%&%&%&%&%&%&%&%&%&%&%&%&%&%&%&%&%&%&%&%&
# KEPLERIAN ORBIT COEFFICIENTS IN CGS
# %&%&%&%&%&%&%&%&%&%&%&%&%&%&%&%&%&%&%&%&%&%&%&%&%&%&%&%&
def nkep(a, mstar=ms):
"""mean motion of orbit at semi-major axis a in CGS"""
return math.sqrt(g * mstar / float(a)**3)
return math.sqrt(g * mstar / float(a) ** 3)
Loading

0 comments on commit f641d09

Please sign in to comment.