compile_opt IDL2

;Constants in c.g.s. units.
G=DOUBLE(6.67408e-8)
m_solar=DOUBLE(1.98855e33)
AU=DOUBLE(1.49598e13)
parsec=DOUBLE(3.08568e18)
year=DOUBLE(365*24*60*60)

;Defining separation and orbital period using Kepler's laws.
m_star=1.0*m_solar
m_planet=9.543e-4*m_solar                       ;Can define any mass,using 1 Jupiter.
a=(FINDGEN(30)+1)*0.5*AU                        ;Separation.
T=SQRT(((a^3)*4.0*!Pi^2)/(G*(m_star+m_planet))) ;Period.

;Defining distance from Earth and angular separation ('a' and 'd' must be in same units)
separation=5*AU                    ;Separation for planet-star, using 5 AU.
d=(FINDGEN(100)+10)*parsec         ;Creates an array of 10 to 109 parsec distances.
arc=(3600*180*separation)/(d*!Pi)  ;Basically converting physical separation into angular one.


;Orbital time versus planet separation
orbit=PLOT(a/AU,t/year,XTITLE='a (AU)', $
  ytitle='Orbital Time(Years)',yrange=[0,10],xrange=[0,5])

;[CHECK]Linear relation between 'a' cubed and 'T' squared to check if relation is correct
orbit_relation=PLOT((a/AU)^3,(t/year)^2,XTITLE='a^3(AU^3)', $
    ytitle='Orbital Time Squared(Years^2)',yrange=[0,10],xrange=[0,10])

;Distance versus angular separation (depending on separation of star-planet)
angular=PLOT(d/parsec,arc,XTITLE='Distance (parsec)', $
  ytitle='Angular Seperation (arcseconds)',yrange=[0,1])


;Center of gravity for a star-planet system.
center=(a*m_planet)/(m_star+m_planet)
p= plot(center/AU,a/AU,YTITLE='a (AU)', $
    xtitle='Center of Mass (AU)')


end