# Galactic Equator zenith transit

I need to find the two times of zenith passage of the Galactic Equator for a given Earth Lat/Lon on any day between 3500 BC and 1500 AD. Since the Galactic Equator is a line and not a point, I try to use PyEphem to create a zenith point for my geographic location (see code below), then find the Galactic coordinates of that point. When the Galactic latitude of this zenith point is zero, that's the time of zenith passage of the Galactic Equator.

However, I think I need my zenith point to be in Equatorial coordinates (not Ecliptic) to accomodate the geographic location and the date/time, but I can't figure out how to make that work.

Here is my Python script based on PyEphem:

import ephem
from math import *

D = (input("Enter date (yyyy/mm/dd): "))
d = ephem.Date(D)
print(d)

location = ephem.Observer()
location.lon = '-87.0'
location.lat = '20.0'
location.elevation = 12
location.pressure = 1013.0
import csv
my_list = []
with open(r'C:\Users\Mark\Desktop\PYTHON\gzt.csv', 'w', newline='') as f:
for T in range(1,1440):
dm = ephem.Date(d + ephem.minute*T)
location.date = dm
ec = ephem.Ecliptic(ra, dec, epoch = d)
eq = ephem.Equatorial(ec)
#print(eq.ra, eq.dec)
z = ephem.FixedBody()
z._ra = eq.ra
z._dec =eq.dec
#z._epoch = d
z.compute()
#print(z.ra, z.dec)
zg = ephem.Galactic(z.a_ra, z.a_dec)
dmt = (('%s') % (dm))
#print(dmt, 'Zlat = ', zg.lat)
AZ = (('%s  %s') % (dmt, zg.lat))
my_list = [AZ]
writer = csv.writer(f, dialect='excel-tab')
writer.writerow(my_list)

• I cannot help you with the programming, but my suggestion would be to find the declination of the zenith (it equals the latitude). Then find the Right Ascension (RA) of the galactic equator where it cross the declination. Then find when that RA is at the zenith, that is, when the local sidereal time equals the RA. That will give you the 0, 1, 2, or 3 times in a day when the galactic equator is at the zenith. – JohnHoltz Apr 10 '20 at 12:38
• Thank you! James K (I think) for fixing my post. – Mark Apr 10 '20 at 16:39
• Thanks JohnHoltz. It's more involved than I imagined. I see now I must at the least pass the hour angle in Sidereal Time from each time step to the Equatorial coordinates. So I will try to adjust my code to do that part now. I'm not an Astronomer or a Programmer, but I need both of those fields for my research, so I teach myself what I can. But there is a "Swiss Cheese" effect to that. – Mark Apr 10 '20 at 17:09

If I eliminate ec and z, and do this:
    eq = ephem.Equatorial(ra, dec, epoch = d)

the crossings I find in gzt.csv agree with Stellarium within a minute of time.
If you're open to another package by the same developer, skyfield.almanac can help with such calculations. In this case you would implement a function which takes a time argument and returns an integer indicating which side of the galactic equator the zenith is on. Then you would pass this function as an argument to find_discrete.