Are transits correctly calculated for eccentric orbits?

[decaelus]

Dear PyAstronomy Team, Thanks for your excellent package. I'm using now for a port of an old code, and it's really coming in handy.

I have a question about the way in which transits are calculated when the orbit is allowed to be eccentric in PyAstronomy version 0.13.0b0. Inside the modelSuite/XTran/forTrans directory in the mandelAgol.py file on line 143, I see the following:

    if orbit == "keplerian":
      self["w"] = -90.0

If I'm reading the documentation correctly, that line is setting the argument of periapsis to -90 degrees. Why make that assumption? If one is modeling an orbit with a different orientation, doesn't this line make the transit calculation come out wrong?

Thanks, Brian

[sczesla]

Hi,

Thanks for scrutinizing the code. As far as I can see, this line defines a starting condition, which is adopted when an instance of that model is created. It can be changed at any point (see code example below) and may be a good or bad assumption for the actual data you have. I think the reason it was implemented like this, is that in the special case of a circular orbit (e=0), it makes the time of periastron passage (tau) the central transit time, which may be convenient.

Cheers

# Import some unrelated modules
import numpy as np
import matplotlib.pylab as plt
# ... and now the forTrans module
from PyAstronomy.modelSuite import forTrans as ft

# Create MandelAgolLC object with
# circular orbit and quadratic limb darkening
ma = ft.MandelAgolLC(orbit="keplerian", ld="quad")

# See the available parameters and their current values
ma.parameterSummary()

# Set parameters
ma["per"] = 1.0
ma["i"] = 90.0
ma["a"] = 6.5
ma["tau"] = 0.0
ma["p"] = 0.16
ma["linLimb"] = 0.6
ma["quadLimb"] = 0.0
ma["b"] = 0.

t = np.linspace(-0.5,0.5,500)
m1 = ma.evaluate(t)

# Now change w
ma["w"] = +123.7

# See the available parameters and their current values
ma.parameterSummary()

m2 = ma.evaluate(t)

plt.plot(t, m1, 'b.-')
plt.plot(t, m2, 'r.-')
plt.show()