Calculates solar distance (centre Sun-centre Earth) given a JDE. Based on low-precision Jean Meeus algorithm (precision maybe 1000’s of km)
"""Solar distance. A low-accuracy algorithm. Calculate solar distance given JDE. Results in AU km or light seconds. """ __version__="0.1" crid="Astronomy:Solar distance V"+__version__+""" © 2023 SteveG1CMZ""" #Jean Meeus low accuracy algorithm from math import * #from aeclipse import * #just for example dates showme=True #select edition1 or edition2 edition1=None #not implemented edition2=1 #24=25 Solar Coordinates #low accuracy (5sf au, 1000's km, 1/10s) rads=pi/180 aum = 149597870700 #au to m (exact) kms = 3.33564e-6 #km to light seconds def solar_distau(JDE): #=25.1 T centuries T=(JDE-2451545.0)/36525 #=25.2 mean long of Sun (degrees) L0=(edition2*280.46646) + 36000.76983*T + 0.0003032*T**2 #=25.3 mean anomaly of Sun (degrees) M=edition2*(357.52911 + 35999.05029*T - 0.0001537*T**2) #=25.4 eccentricity (Earth orbit) ecc=edition2*(0.016708634 - 0.000042037*T -0.0000001267*T**2) Mr=M*rads #Suns equation of centre (degrees) C= edition2*( (1.914602 - 0.004817*T - 0.000014*T**2)*sin(Mr) +(0.019993 - 0.000101*T)*sin(2*Mr) +0.000289*sin(3*Mr)) #Sun true Suntruelong=L0+C #degrees ??? Suntrueanomaly=M+C #degrees vr=Suntrueanomaly*rads #print("T L0 M ecc C long anom") #print(T,L0,M,ecc) #print(C,Suntruelong,Suntrueanomaly) #25.5 Sun-Earth distance (centre-centre) den=(1+ecc*cos(vr)) R=(1.000001018*(1-ecc**2))/den #au return round(R,5) #au def solar_distkm(JDE): Rau=solar_distau(JDE) return round(Rau*aum/1000) #km def solar_dists(JDE): Rkm=solar_distkm(JDE) return round(Rkm*kms,1) #light seconds def examples(): JDE=2448908.5 print(JDE) au=solar_distau(JDE) print(au,"au") km=solar_distkm(JDE) print(km,"km") print(solar_dists(JDE),"s") #today today=2460246 #when I wrote this line au=solar_distau(today) print(au,"au") km=solar_distkm(today) print(km,"km") print(solar_dists(today),"s") if showme or __name__=="__main__": print(crid) print("distance Earth-Sun (centre-centre)") #print("au: given a JDE") #print("see aeclipse to change Date to JDE") #25.a examples()