#include <windows.h>
#include <stdio.h>
#include <math.h>
#include <time.h> // to see if daylight savings time

//#define AU 92955807.2730255 // 1 astronomical unit in miles
#define AU 149597870.700 // 1 astronomical unit in kilometers (officially designated)
//#define MILESPERKILOMETER 0.62137119224
#define PI 3.141592653589793
#define DEG2RAD 0.017453292519943295 // PI/180.0
#define RAD2DEG 57.295779513082323 // 180.0/PI
#define SunsEccentricity 0.016708
#define SunsEclipticLongitudeAtEpoch 280.466069
#define SunsEclipticLongitudeAtPerigeeAtEpoch 282.938346
#define MoonsEclipticLongitudeAtEpoch 218.316433
#define MoonsEclipticLongitudeAtPerigeeAtEpoch 83.353451
#define InclinationofMoonsOrbitVsEcliptic 5.1453964
#define EclipticObliquityOfEpoch2000 23.439292
#define EccentricityOfMoonsOrbit 0.0549
#define MoonsSemiMajorAxis 384400 // km
#define EccentricityOfSunsOrbit 0.016708
#define SunsSemiMajorAxis 149598500 // km

//double SunRadius;
//double MoonRadius, VenusRadius, EarthRadius, MarsRadius, JupiterRadius, SaturnRadius;
double rad, Rotation = 360.0, D, LonSign, dDay, dHour, dMinute, dSecond, DecimalHour, DecimalDay, JD, JD0, Days, dT, dR, dB, dT0, UT, GST, Adjust, LST, Zoom = 50.0;
double Epoch2000, EpochDecimalDay, DaysFromEpoch2000, SunsMeanAnomaly, SunsEquationOfTheCenter, SunsTrueAnomaly, SunsEclipticLongitude;
double MoonsMeanEclipticLongitude, MoonsMeanEclipticLongitudeAscendingNodeAtEpoch, MoonsMeanEclipticLongitudeAscendingNode, MoonsMeanAnomaly;
double AnnualEquation, Evection, MeanAnomalyCorrection, MoonsTrueAnomaly, VariationCorrection, MoonsEclipticLongitude, Y, X, T, MoonsEclipticLatitude;
double MoonsAge, PhaseAngle, PercentIllumination, DistanceToMoon, EarthDistanceToSun, SunAltitude, SunAzimuth, MoonAltitude, MoonAzimuth, VenusAltitude, VenusAzimuth;
double DistanceFromCenter, Declination, RightAscension, Azimuth, Altitude, HourAngle, T0, T1, T2, ObserversLatitude = 37.407367, ObserversLongitude = -122.116302;
double EarthMeanAnomoly, EarthTrueAnomoly, EarthEquationOfTheCenter, EarthHeliocentricEclipticLongitude, EarthEclipticLongitudeAtPerihelion, EarthHeliocentricEclipticLatitude, EarthRadiusVectorLength;
double MeanAnomoly[5], EquationOfTheCenter[5], TrueAnomoly[5], HeliocentricEclipticLongitude[5], HeliocentricEclipticLatitude[5], RadiusVectorLength[5], Adjustment[5], GeocentricEclipticLongitude[5], GeocentricEclipticLatitude[5], HeliocentricEclipticLatitude[5];
double b[5], PlanetAltitude[5], PlanetAzimuth[5];
double DistanceToEarth[5], DistanceToSun[5];
double Moons_a = 0.00257;
double Moons_e = 0.0549;
double MoonsRadiusVectorLength;
int vertMoon, horizMoon;
int vertPlanet[5], horizPlanet[5];
int xPosDown, yPosDown, xPosUp, yPosUp, chars;
int ScreenWidth, ScreenHeight, xCenter, yCenter, UToffset, xCenterSun, yCenterSun, vertSun, horizSun, vertMoon, horizMoon, vertVenus, horizVenus, chars, Radius;
int LonDeg, x, y, z, m2, Year, Month, Day, OrigDay, Hour, H, Minute, Second, A = 0, B = 0, EpochYear, EpochMonth; // 0 because Year is never < 0 and never < Oct 15, 1582;

// planets:         (0=Venus,   1=Earth    2=Mars,    3=Jupiter, 4=Saturn):
double Tau[5] =     {0.615197,  1.000017,  1.880848,  11.862615, 29.447498};  // orbital period in tropical years
double m[5] =       {0.814998,  1.0,       0.107447,  317.828133,95.160904};  // mass relative to Earth
double r[5] =       {6051.8,    6378.14,   3389.5,    69911,     58232};      // radius in km
double e[5] =       {0.0067767, 0.0167112, 0.093394,  0.048393,  0.053862};   // orbital eccentricity
double a[5] =       {0.723336,  1.000003,  1.523710,  5.202887,  9.536676};   // length of the orbital semi-major axis in AUs
// b[z] = a[z]*sqrt(1-(e[z]*e[z])); // b is semi-minor length
double theta[5] =   {16.92,                9.36,      196.74,	 165.60};     // angular diameter at 1 AU
double V[5] =       {-4.40,                -1.52,     -9.40,     -8.88};      // visual magnitude at 1 AU
double mu[5] =      {324.860,   398.600,   42828,     126687000, 37931000};   // standard gravitational parameter in km /s
double iota[5] =    {3.394676,  -0.000015, 1.849691,  1.3043975, 2.485992};   // inclination of the orbital plane w.r.t. the eclipse
double epsilon[5] = {181.979100,100.464572,-4.553432, 34.396441, 49.954244};  // ecliptic longitude at the epoch
double omega[5] =   {131.602467,102.937682,-23.943630,14.728480, 92.598878};  // ecliptic longitude at perihelion
double Omega[5] =   {76.679843,            49.559539, 100.473909,113.662424}; // ecliptic longitude of the ascending node at epoch

DaysInYearMonths[] = {0, 31, 28, 31, 30, 31, 30, 31, 31, 30, 31, 30, 31};
DaysInLeapYearMonths[] = {0, 31, 29, 31, 30, 31, 30, 31, 31, 30, 31, 30, 31};
DWORD fileSize, dwBytesRead;
BOOL first = TRUE, showsun, showmoon, showPlanet[5];
char String[64];
char cLatitude[32];
char cLongitude[32];
char szAppName[] = "SimpleSolarSystem";
char LatLon[] = "C:\\Users\\Public\\LatLon.txt";
char Planet[5][8] = {"Venus","Earth","Mars","Jupiter","Saturn"};
char String[64];

time_t rawtime; // for timeinfo
struct tm *timeinfo; // for daylight savings time
SYSTEMTIME SystemTime;
HANDLE hFile;
HWND hwnd;
HINSTANCE hInstance;
HBRUSH hBrush, hDaySkyBrush, hNightSkyBrush, hSunBrush, hPlanetBrush, hMoonBrush, hCenterBrush;
HPEN hPen, hSunPen, hMoonPen, hPlanetPen, hOrbitPen;
HGDIOBJ hOldPen, hOldBrush;
HDC hdc, hdcMem;
HBITMAP hBitmap, hMemBitmap;
PAINTSTRUCT ps;
RECT rect;

LRESULT CALLBACK WndProc(HWND, UINT, WPARAM, LPARAM);

void format_commas(int n, char *out)
{
    int c;
    char buf[20];
    char *p;

    sprintf(buf, "%d", n);
    c = 2 - strlen(buf) % 3;
    for (p = buf; *p != 0; p++) {
       *out++ = *p;
       if (c == 1)
           *out++ = ',';
       c = (c + 1) % 3;
    }
    *--out = 0;
}

double CorrectTo360(double d)
{
	x = (int)d;
	if (x < 0) {
		while (x < 0) {
			x += 360;
			d += 360.0;
		};
	}
	else if (x >= 360) {
		while (x >= 360) {
			x -= 360;
			d -= 360.0;
		};
	}
	return d;
}

void ConvertEclipticToEquatorial(double MoonsEclipticLatitude, double MoonsEclipticLongitude)
{ // eclipticlongitude is lamda AND eclipticlatitude is Beta AND obliquity of the ecliptic is epsilon
	Declination = RAD2DEG*asin(sin(MoonsEclipticLatitude*DEG2RAD)*cos(EclipticObliquityOfEpoch2000*DEG2RAD) + cos(MoonsEclipticLatitude*DEG2RAD)*sin(EclipticObliquityOfEpoch2000*DEG2RAD)*sin(MoonsEclipticLongitude*DEG2RAD));
	Y = ((sin(MoonsEclipticLongitude*DEG2RAD))*cos(EclipticObliquityOfEpoch2000*DEG2RAD)) - (tan(MoonsEclipticLatitude*DEG2RAD)*sin(EclipticObliquityOfEpoch2000*DEG2RAD));
	X = cos(MoonsEclipticLongitude*DEG2RAD);
	RightAscension = RAD2DEG*atan(Y/X);
	if ((X < 0) && (Y > 0)) RightAscension += 180.0;
	else if ((X < 0.0) && (Y < 0.0)) RightAscension += 180.0;
	else if ((X > 0.0) && (Y < 0.0)) RightAscension += 360.0;
	RightAscension /= 15.0; // convert to hours
}

void ConvertEquatorialToHorizon(double Declination, double HourAngle)
{
	T0 = (sin(Declination*DEG2RAD)*sin(ObserversLatitude*DEG2RAD)) + (cos(Declination*DEG2RAD)*cos(ObserversLatitude*DEG2RAD)*cos(HourAngle*DEG2RAD));
	Altitude = RAD2DEG*asin(T0);
	T1 = (sin(Declination*DEG2RAD) - (sin(ObserversLatitude*DEG2RAD)*sin(Altitude*DEG2RAD)));
	T2 = T1 / (cos(ObserversLatitude*DEG2RAD)*cos(Altitude*DEG2RAD));
	Azimuth = RAD2DEG*acos(T2);
	D = sin(HourAngle*DEG2RAD);
	if (D > 0.0)
		Azimuth = 360.0 - Azimuth;
}

double GetJulianDay(double DecimalDay, int Year, int Month)
{
	if (Month > 2) {
		y = Year;
		m2 = Month;
	}
	else {
		y = Year - 1;
		m2 = Month + 12;
	}
	A = y / 100;
	B = 2 - A + (A / 4);
	D = (double)(B + (int)(365.25*(double)y) + (int)(30.6001*((double)(m2+1))) + DecimalDay + 1720994.5);
	return D;
}

double UTtoGST(double UT, int Year, int Month)
{
	double jd;

	jd = GetJulianDay(dDay, Year, Month);
	JD0 = GetJulianDay(0.0, Year, 1);
	Days = jd-JD0;
	dT = (JD0-2415020.0)/36525.0;
	dR = 6.6460656+(2400.051262*dT)+(0.00002581*dT*dT);
	dB = 24.0-dR+(double)(24*(Year-1900));
	dT0 = (0.0657098*Days)-dB;
	return (dT0+(1.002738*UT));
}

void GetAltitudeAzimuth(void)
{
	dDay = (double)Day;
	if (timeinfo->tm_isdst) {
		H = Hour-1;// H is NOT daylight savings time -- H is used in calculations -- Hour is displayed
		if (H == -1) {
			H = 23;
			dDay--;
		}
		if (H == 24) {
			H = 0;
			dDay++;
		}
	}
	else
		H = Hour;
	dHour = (double)H;
	dMinute = (double)Minute;
	dSecond = (double)Second;
	DecimalHour = dHour + (dMinute/60.0) + (dSecond/3600.0);
	LonDeg = abs((int)ObserversLongitude);
	if (ObserversLongitude < 0.0) // West
		UToffset = (int)(LonDeg / 15.0); // doesn't account for daylight savings time, so don't send daylight savings time!
	else // East
		UToffset = (int)((360.0 - LonDeg) / 15.0);
	DecimalDay = dDay + (DecimalHour/24.0) + (UToffset/24.0);
	UT = DecimalHour + UToffset;
	if (UT >= 24.0)
		UT -= 24.0;
//	UT += 63.8/3600.0; // rough correction to seconds since Epoch Year 2000
	JD = GetJulianDay(DecimalDay, Year, Month);
	GST = UTtoGST(UT, Year, Month);
	if (GST < 1.0) {
		D = 24.0;
		GST = D + GST;
	}
	Adjust = ObserversLongitude / 15.0;
	LST = GST + Adjust;
	if (LST < 1.0) LST = LST + 24.0;
	else if (LST > 24.0) LST = LST - 24.0;
	EpochDecimalDay = 1.5;
	EpochMonth = 1;
	EpochYear = 2000;
	Epoch2000 = GetJulianDay(EpochDecimalDay, EpochYear, EpochMonth);
	DaysFromEpoch2000 = JD - Epoch2000;

//	MOON
	MoonsMeanEclipticLongitude = (13.176339686*DaysFromEpoch2000) + MoonsEclipticLongitudeAtEpoch;
	MoonsMeanEclipticLongitude = CorrectTo360(MoonsMeanEclipticLongitude);
	MoonsMeanEclipticLongitudeAscendingNodeAtEpoch = 125.044522; // Table 7.1
	MoonsMeanEclipticLongitudeAscendingNode = MoonsMeanEclipticLongitudeAscendingNodeAtEpoch - (0.0529539*DaysFromEpoch2000);
	MoonsMeanEclipticLongitudeAscendingNode = CorrectTo360(MoonsMeanEclipticLongitudeAscendingNode);
	MoonsMeanAnomaly = MoonsMeanEclipticLongitude - (0.1114041*DaysFromEpoch2000) - MoonsEclipticLongitudeAtPerigeeAtEpoch;
	MoonsMeanAnomaly = CorrectTo360(MoonsMeanAnomaly); // 12.
	AnnualEquation = 0.1858*sin(SunsMeanAnomaly*DEG2RAD); // 13.
	Evection = 1.2739*sin(((2*(MoonsMeanEclipticLongitude-SunsEclipticLongitude))-MoonsMeanAnomaly)*DEG2RAD);
	MeanAnomalyCorrection = MoonsMeanAnomaly + Evection - AnnualEquation - (0.37*sin(SunsMeanAnomaly*DEG2RAD));
	MoonsTrueAnomaly = (6.2886*sin(MeanAnomalyCorrection*DEG2RAD)) + (0.214*sin(2*MeanAnomalyCorrection*DEG2RAD));
	MoonsMeanEclipticLongitude = (MoonsMeanEclipticLongitude + MoonsTrueAnomaly + Evection - AnnualEquation); // 7.3.9 lamda prime
	VariationCorrection = 0.6583*sin(2.0*(MoonsMeanEclipticLongitude-SunsEclipticLongitude)*DEG2RAD);
	MoonsEclipticLongitude = MoonsMeanEclipticLongitude + VariationCorrection;
	MoonsMeanEclipticLongitudeAscendingNode -= 0.16*sin(SunsMeanAnomaly*DEG2RAD);
	MoonsMeanEclipticLongitudeAscendingNode = CorrectTo360(MoonsMeanEclipticLongitudeAscendingNode);
	Y = sin((MoonsEclipticLongitude-MoonsMeanEclipticLongitudeAscendingNode)*DEG2RAD)*cos(InclinationofMoonsOrbitVsEcliptic*DEG2RAD);
	X = cos((MoonsEclipticLongitude-MoonsMeanEclipticLongitudeAscendingNode)*DEG2RAD);
	T = RAD2DEG*atan(Y/X);
	if ((X < 0) && (Y > 0)) T = T + 180.0;
	else if ((X < 0.0) && (Y < 0.0)) T = T + 180.0;
	else if ((X > 0.0) && (Y < 0.0)) T = T + 360.0;
	MoonsEclipticLongitude = MoonsMeanEclipticLongitudeAscendingNode + T;
	MoonsEclipticLongitude = CorrectTo360(MoonsEclipticLongitude); // Lmoon
	MoonsRadiusVectorLength = (Moons_a*(1-(Moons_e*Moons_e))) / (1+(Moons_e*cos(MoonsTrueAnomaly*DEG2RAD))); // R (in Astronomical Units)

// PLANETS
	for (z = 0; z < 5; z++) {
 		MeanAnomoly[z] = ((360.0*DaysFromEpoch2000)/(365.242191*Tau[z]))+epsilon[z]-omega[z];
		if (MeanAnomoly[z] >= 360.0) {
			do { MeanAnomoly[z] -= 360.0; } while (MeanAnomoly[z] >= 360.0);
		}
		EquationOfTheCenter[z] = (360.0/PI)*e[z]*sin(MeanAnomoly[z]*DEG2RAD); // E
		TrueAnomoly[z] = MeanAnomoly[z] + EquationOfTheCenter[z]; // upsilon (like small v)
		HeliocentricEclipticLongitude[z] = TrueAnomoly[z] + omega[z]; // L (omega = ecliptic longitude at perihelion)
		HeliocentricEclipticLatitude[z] = RAD2DEG*asin(sin((HeliocentricEclipticLongitude[z]-Omega[z])*DEG2RAD)*sin(iota[z]*DEG2RAD)); // upside-down V (Lamda)
		RadiusVectorLength[z] = (a[z]*(1-(e[z]*e[z]))) / (1+(e[z]*cos(TrueAnomoly[z]*DEG2RAD))); // R (in Astronomical Units)
	} // end of for (z = 0; z < 5; z++)
}

int WINAPI WinMain(HINSTANCE hInstance, HINSTANCE hPrevInstance, PSTR szCmdLine, int iCmdShow)
{
	MSG          msg;
	WNDCLASS     wndclass;

	wndclass.style         = CS_HREDRAW|CS_VREDRAW;
	wndclass.lpfnWndProc   = WndProc;
	wndclass.cbClsExtra    = 0;
	wndclass.cbWndExtra    = 0;
	wndclass.hInstance     = hInstance;
	wndclass.hIcon         = LoadIcon(NULL, IDI_APPLICATION);
	wndclass.hCursor       = LoadCursor(NULL, IDC_ARROW);
	wndclass.hbrBackground = (HBRUSH)GetStockObject(BLACK_BRUSH);
	wndclass.lpszMenuName  = NULL;
	wndclass.lpszClassName = szAppName;

	if (!RegisterClass(&wndclass))
		return 0;

	hwnd = CreateWindow(szAppName, szAppName,
		WS_OVERLAPPEDWINDOW,
		0, 0, 0, 0,
		NULL, NULL, hInstance, NULL);

	ShowWindow(hwnd, SW_SHOWMAXIMIZED);
	UpdateWindow(hwnd);

	while (GetMessage(&msg, NULL, 0, 0))
	{
		TranslateMessage(&msg);
		DispatchMessage(&msg);
	}
	return msg.wParam;
}

LRESULT CALLBACK WndProc(HWND hwnd, UINT message, WPARAM wParam, LPARAM lParam)
{
	switch(message)
	{
	case WM_CREATE:
		hDaySkyBrush = CreateSolidBrush(0xFFDFBF);
		hNightSkyBrush = CreateSolidBrush(0x606060);
		hPen = CreatePen(PS_SOLID, 1, 0);
		hSunBrush = CreateSolidBrush(0x00E0F0);
		hSunPen = CreatePen(PS_SOLID, 1, 0x00E0F0);
		hPlanetBrush = CreateSolidBrush(0xEEEEEE);
		hPlanetPen = CreatePen(PS_SOLID, 1, 0xEEEEEE);
		hMoonBrush = CreateSolidBrush(0xFF0000);
		hMoonPen = CreatePen(PS_SOLID, 1, 0xFF0000);
		hOrbitPen = CreatePen(PS_SOLID, 1, 0x404040);
		hCenterBrush = CreateSolidBrush(0);
		ScreenWidth = GetSystemMetrics(SM_CXFULLSCREEN);
		ScreenHeight = GetSystemMetrics(SM_CYFULLSCREEN);
		rect.bottom = ScreenHeight;
		rect.right = ScreenWidth;
		rect.left = rect.top = 0;
		xCenter = ScreenWidth / 2;
		yCenter = ScreenHeight / 2;
		if (INVALID_HANDLE_VALUE != (hFile = CreateFile(LatLon, GENERIC_READ, 0, NULL, OPEN_EXISTING, FILE_ATTRIBUTE_NORMAL, NULL))) {
			if (fileSize = GetFileSize(hFile, NULL)) {
				ReadFile(hFile, String, fileSize, &dwBytesRead, NULL); // current engine name
				CloseHandle(hFile);
				for (x = 0, y = 0; (x < (int)fileSize) && (String[x] != ','); x++, y++)
					cLatitude[y] = String[x];
				cLatitude[y] = 0;
				for (x += 2, y = 0; x < (int)fileSize; x++, y++)
					cLongitude[y] = String[x];
				cLongitude[y] = 0;
				ObserversLatitude = atof(cLatitude);
				ObserversLongitude = atof(cLongitude);
				LonDeg = abs((int)ObserversLongitude);
				if (ObserversLongitude <= 0.0) // West
					UToffset = (int)(LonDeg / 15.0); // doesn't account for daylight savings time, so don't send daylight savings time!
				else // East
					UToffset = (int)((360.0 - LonDeg) / 15.0);
			}
		}
		GetLocalTime(&SystemTime);
		Year = SystemTime.wYear;
		Month = SystemTime.wMonth;
		Day = SystemTime.wDay;
		OrigDay = Day;
		Hour = SystemTime.wHour;
		Minute = SystemTime.wMinute;
		Second = SystemTime.wSecond;
		time(&rawtime);
		timeinfo = localtime(&rawtime);
		break;

	case WM_KEYDOWN:
		switch (wParam)
		{
		case VK_ESCAPE:
			DestroyWindow(hwnd);
			break;
		case 187: // '+':
			Zoom += 50.0;
			InvalidateRect(hwnd, &rect, FALSE);
			break;
		case 189: // '-':
			if (Zoom > 50.0) {
				Zoom -= 50.0;
				InvalidateRect(hwnd, &rect, FALSE);
			}
			break;
		case 'R':
			if (GetKeyState(VK_SHIFT) & 0x8000)
				Rotation += 10.0;
			else
				Rotation -= 10.0;
			InvalidateRect(hwnd, &rect, FALSE);
			break;
/*
		if (wParam == 'D') {
			if (goback) {
				if (Day > 1)
					Day--;
				else {
					if (Year % 4)
						Day = DaysInLeapYearMonths[Month];
					else
						Day = DaysInYearMonths[Month];
					Month--;
					if (Month == 0) {
						Day = 31;
						Month = 12;
						Year--;
					}
				}
			}
			else {
				if ((Year % 4) && (Day < DaysInYearMonths[Month]))
					Day++;
				else if ((0 == (Year % 4)) && ((Day) < DaysInLeapYearMonths[Month]))
					Day++;
				else {
					Day = 1;
					Month++;
					if (Month == 13) {
						Day = 1;
						Month = 1;
						Year++;
					}
				}
			}
			InvalidateRect(hwnd, &rect, FALSE);
*/
		case 'D':
			if (GetKeyState(VK_SHIFT) & 0x8000) {
				if (Day > 1)
					Day--;
				else {
					if (Year % 4)
						Day = DaysInLeapYearMonths[Month];
					else
						Day = DaysInYearMonths[Month];
					Month--;
					if (Month == 0) {
						Month = 12;
						Year--;
					}
				}
			}
			else {
				if ((Year % 4) && (Day < DaysInLeapYearMonths[Month]))
					Day++;
				else if ((0 == (Year % 4)) && ((Day) < DaysInYearMonths[Month]))
					Day++;
				else {
					Day = 1;
					Month++;
					if (Month == 13) {
						Month = 1;
						Year++;
					}
				}
			}
			SendMessage(hwnd, WM_PAINT, 0, 0);
			InvalidateRect(hwnd, &rect, FALSE);
			break;
		case VK_UP:
			yCenter += 50;
			InvalidateRect(hwnd, &rect, FALSE);
			break;
		case VK_DOWN:
			yCenter -= 50;
			InvalidateRect(hwnd, &rect, FALSE);
			break;
		case VK_LEFT:
			xCenter += 50;
			InvalidateRect(hwnd, &rect, FALSE);
			break;
		case VK_RIGHT:
			xCenter -= 50;
			InvalidateRect(hwnd, &rect, FALSE);
			break;
		}
		break;

	case WM_MOUSEWHEEL:
		if (wParam == 0xFF880000)
			SendMessage(hwnd, WM_KEYDOWN, 189, 0);
		else if (wParam == 0x00780000)
			SendMessage(hwnd, WM_KEYDOWN, 187, 0);
		break;

	case WM_PAINT:
		hdc = BeginPaint(hwnd, &ps);
		if (first) {
			first = FALSE;
			if (!hdcMem) {
				hdcMem = CreateCompatibleDC(hdc);
				hMemBitmap = CreateCompatibleBitmap(hdc, rect.right, rect.bottom);
				SelectObject(hdcMem, hMemBitmap);
			}
		}
//		SunRadius = (695700.0*Zoom) / AU; // to show if greater than 5 pixel radius
//		MoonRadius = (1738.1*Zoom) / AU; // to show if greater than 2 pixel radius
//		VenusRadius = (r[0]*Zoom) / AU;
//		EarthRadius = (r[1]*Zoom) / AU;
//		MarsRadius = (r[2]*Zoom) / AU;
//		JupiterRadius = (r[3]*Zoom) / AU;
//		SaturnRadius = (r[4]*Zoom) / AU;
		FillRect(hdcMem, &rect, hCenterBrush);
		SetBkMode(hdcMem, TRANSPARENT);
		SetTextColor(hdcMem, 0xD0D0D0);
		//////////////////
		GetAltitudeAzimuth();
		//////////////////
		hOldPen = SelectObject(hdcMem, hOrbitPen);
		hOldBrush = SelectObject(hdcMem, hCenterBrush);
		for (z = 4; z >= 0; z--) { // planet circular orbits (eccentric but circular at these resolutions)
			Radius = (int)(RadiusVectorLength[z]*Zoom); // Astronomical Units
			Ellipse(hdcMem, xCenter-Radius, yCenter-Radius, xCenter+Radius, yCenter+Radius);
		}

		SelectObject(hdcMem, hOldBrush);
		SelectObject(hdcMem, hOldPen);

		hOldPen = SelectObject(hdcMem, hPlanetPen);
		hOldBrush = SelectObject(hdcMem, hPlanetBrush);
		for (z = 0; z < 5; z++) {
			rad = (Rotation-HeliocentricEclipticLongitude[z])*DEG2RAD; // default Rotation is 360.0- to make planets go counter-clockwise around sun (viewing from North Star)
			vertPlanet[z] = (int)(RadiusVectorLength[z]*Zoom*sin(rad)); // Astronomical Units
			horizPlanet[z] = (int)(RadiusVectorLength[z]*Zoom*cos(rad)); // Astronomical Units
			Ellipse(hdcMem, xCenter-2+horizPlanet[z], yCenter-2+vertPlanet[z], xCenter+2+horizPlanet[z], yCenter+2+vertPlanet[z]);
			TextOut(hdcMem, xCenter+5+horizPlanet[z], yCenter-8+vertPlanet[z],Planet[z], 8);
		}
		SelectObject(hdcMem, hOldBrush);
		SelectObject(hdcMem, hOldPen);

		rad = (Rotation-MoonsEclipticLongitude)*DEG2RAD;
		vertMoon = (int)(10.0*MoonsRadiusVectorLength*Zoom*sin(rad)); // 10.0 makes the moon 10 times further from Earth than proper proportion
		horizMoon = (int)(10.0*MoonsRadiusVectorLength*Zoom*cos(rad));
		hOldPen = SelectObject(hdcMem, hMoonPen);
		hOldBrush = SelectObject(hdcMem, hMoonBrush);
		Ellipse(hdcMem, xCenter+horizPlanet[1]-2+horizMoon, yCenter+vertPlanet[1]-2+vertMoon, xCenter+horizPlanet[1]+2+horizMoon, yCenter+vertPlanet[1]+2+vertMoon);
		SelectObject(hdcMem, hOldBrush);
		SelectObject(hdcMem, hOldPen);

		hOldPen = SelectObject(hdcMem, hSunPen);
		hOldBrush = SelectObject(hdcMem, hSunBrush);
//		if (SunRadius > 5.0)
//			Ellipse(hdcMem, xCenter-(int)(SunRadius+0.5), yCenter-(int)(SunRadius+0.5), xCenter+(int)(SunRadius+0.5), yCenter+(int)(SunRadius+0.5));
//		else
			Ellipse(hdcMem, xCenter-5, yCenter-5, xCenter+5, yCenter+5);
		SelectObject(hdcMem, hOldBrush);
		SelectObject(hdcMem, hOldPen);

		TextOut(hdcMem, 20, 100, "Press Esc to exit", 17);
		TextOut(hdcMem, 20, 120, "Press + or - OR scroll the mouse wheel to zoom in or out", 56);
		TextOut(hdcMem, 20, 140, "Press D to go forward a day, or Shift-D to go back a day", 56);
		TextOut(hdcMem, 20, 160, "Press R or Shift-R to rotate the solar system", 45);
		TextOut(hdcMem, 20, 180, "Press an Arrow key to move the solar system", 43);
//		chars = sprintf_s(String, 64, "Latitude, Longitude: %.6f, %.6f", ObserversLatitude, ObserversLongitude);
//		TextOut(hdcMem, 20, 220, String, chars);
		chars = sprintf_s(String, 64, "%i/%i/%i", Month, Day, Year);
		TextOut(hdcMem, 20,220, String, chars);
		chars = sprintf_s(String, 64, "%2i:%2i", Hour, Minute);
		TextOut(hdcMem, 20, 240, String, chars);
//		TextOut(hdcMem, 20, 800,  "Latitude and Longitude makes a slight difference in where the moon is.", 70);
//		TextOut(hdcMem, 20, 820,  "The Latitude and Longitude selected is mine, and you can change it.", 67);
//		TextOut(hdcMem, 20, 840,  "In Google Maps, right-click on the desired location,", 52);
//		TextOut(hdcMem, 20, 860,  "and select 'What's Here'.", 25);
//		TextOut(hdcMem, 20, 880,  "you'll see the latitude, longitude at the bottom of the map.", 60);
//		TextOut(hdcMem, 20, 900, "Left-click on the latitude, longitude to see it again", 53);
//		TextOut(hdcMem, 20, 920, "in another box, which you can copy and then paste", 49);
//		TextOut(hdcMem, 20, 940, "into a file in C:\\Users\\Public named LatLon.txt", 48);
//		TextOut(hdcMem, 20, 960, "(if that file doesn't exist, first create and save it in Notepad).", 66);

		BitBlt(hdc, 0, 0, rect.right, rect.bottom, hdcMem, 0, 0, SRCCOPY);
		EndPaint(hwnd, &ps);
		break;

	case WM_DESTROY:
		PostQuitMessage(0);
		break;
	}
	return DefWindowProc(hwnd, message, wParam, lParam);
}