// CurvatureOfTheEarth.cpp : Defines the entry point for the application.
//
#include "pch.h"
#include "framework.h"
#include "CurvatureOfTheEarth.h"
#include "GreatCircleDistance.h"
#include "Vincenty.h"
#define MAX_LOADSTRING 100
// Global Variables:
HINSTANCE hInst; // current instance
WCHAR szTitle[MAX_LOADSTRING]; // The title bar text
WCHAR szWindowClass[MAX_LOADSTRING]; // the main window class name
WCHAR line[128]; // general purpose buffer
std::wstring outputText; // output wide character text
// Forward declarations of functions included in this code module:
ATOM MyRegisterClass(HINSTANCE hInstance);
BOOL InitInstance(HINSTANCE, int);
LRESULT CALLBACK WndProc(HWND, UINT, WPARAM, LPARAM);
INT_PTR CALLBACK About(HWND, UINT, WPARAM, LPARAM);
INT_PTR CALLBACK TableDialog(HWND, UINT, WPARAM, LPARAM);
int APIENTRY wWinMain(_In_ HINSTANCE hInstance,
_In_opt_ HINSTANCE hPrevInstance,
_In_ LPWSTR lpCmdLine,
_In_ int nCmdShow)
{
UNREFERENCED_PARAMETER(hPrevInstance);
UNREFERENCED_PARAMETER(lpCmdLine);
// TODO: Place code here.
// Initialize global strings
LoadStringW(hInstance, IDS_APP_TITLE, szTitle, MAX_LOADSTRING);
LoadStringW(hInstance, IDC_CURVATUREOFTHEEARTH, szWindowClass, MAX_LOADSTRING);
MyRegisterClass(hInstance);
// Perform application initialization:
if (!InitInstance (hInstance, nCmdShow))
{
return FALSE;
}
HACCEL hAccelTable = LoadAccelerators(hInstance, MAKEINTRESOURCE(IDC_CURVATUREOFTHEEARTH));
MSG msg;
// Main message loop:
while (GetMessage(&msg, nullptr, 0, 0))
{
if (!TranslateAccelerator(msg.hwnd, hAccelTable, &msg))
{
TranslateMessage(&msg);
DispatchMessage(&msg);
}
}
return (int) msg.wParam;
}
//
// FUNCTION: MyRegisterClass()
//
// PURPOSE: Registers the window class.
//
ATOM MyRegisterClass(HINSTANCE hInstance)
{
WNDCLASSEXW wcex = { 0 };
wcex.cbSize = sizeof(WNDCLASSEX);
wcex.style = CS_HREDRAW | CS_VREDRAW;
wcex.lpfnWndProc = WndProc;
wcex.cbClsExtra = 0;
wcex.cbWndExtra = 0;
wcex.hInstance = hInstance;
wcex.hIcon = LoadIcon(hInstance, MAKEINTRESOURCE(IDI_CURVATUREOFTHEEARTH));
wcex.hCursor = LoadCursor(nullptr, IDC_ARROW);
wcex.hbrBackground = (HBRUSH)(COLOR_WINDOW+1);
wcex.lpszMenuName = MAKEINTRESOURCEW(IDC_CURVATUREOFTHEEARTH);
wcex.lpszClassName = szWindowClass;
wcex.hIconSm = LoadIcon(wcex.hInstance, MAKEINTRESOURCE(IDI_SMALL));
return RegisterClassExW(&wcex);
}
//
// FUNCTION: InitInstance(HINSTANCE, int)
//
// PURPOSE: Saves instance handle and creates main window
//
// COMMENTS:
//
// In this function, we save the instance handle in a global variable and
// create and display the main program window.
//
BOOL InitInstance(HINSTANCE hInstance, int nCmdShow)
{
hInst = hInstance; // Store instance handle in our global variable
HWND hWnd = CreateWindowW(szWindowClass, szTitle, WS_OVERLAPPEDWINDOW,
CW_USEDEFAULT, 0, CW_USEDEFAULT, 0, nullptr, nullptr, hInstance, nullptr);
if (!hWnd)
{
return FALSE;
}
ShowWindow(hWnd, nCmdShow);
UpdateWindow(hWnd);
return TRUE;
}
//
// FUNCTION: WndProc(HWND, UINT, WPARAM, LPARAM)
//
// PURPOSE: Processes messages for the main window.
//
// WM_COMMAND - process the application menu
// WM_PAINT - Paint the main window
// WM_DESTROY - post a quit message and return
//
//
LRESULT CALLBACK WndProc(HWND hWnd, UINT message, WPARAM wParam, LPARAM lParam)
{
switch (message)
{
case WM_CREATE:
DialogBox(hInst, MAKEINTRESOURCE(IDD_TABLE_DIALOG), hWnd, TableDialog);
break;
case WM_COMMAND:
{
int wmId = LOWORD(wParam);
// Parse the menu selections:
switch (wmId)
{
case IDM_ABOUT:
DialogBox(hInst, MAKEINTRESOURCE(IDD_ABOUTBOX), hWnd, About);
break;
case IDM_EXIT:
DestroyWindow(hWnd);
break;
default:
return DefWindowProc(hWnd, message, wParam, lParam);
}
}
break;
case WM_PAINT:
{
PAINTSTRUCT ps;
HDC hdc = BeginPaint(hWnd, &ps);
// TODO: Add any drawing code that uses hdc here...
EndPaint(hWnd, &ps);
}
break;
case WM_DESTROY:
PostQuitMessage(0);
break;
default:
return DefWindowProc(hWnd, message, wParam, lParam);
}
return 0;
}
// Message handler for about box.
INT_PTR CALLBACK About(HWND hDlg, UINT message, WPARAM wParam, LPARAM lParam)
{
UNREFERENCED_PARAMETER(lParam);
switch (message)
{
case WM_INITDIALOG:
return (INT_PTR)TRUE;
case WM_COMMAND:
if (LOWORD(wParam) == IDOK || LOWORD(wParam) == IDCANCEL)
{
EndDialog(hDlg, LOWORD(wParam));
return (INT_PTR)TRUE;
}
break;
}
return (INT_PTR)FALSE;
}
INT_PTR CALLBACK TableDialog(HWND hDlg, UINT message, WPARAM wParam, LPARAM lParam)
{
UNREFERENCED_PARAMETER(lParam);
static double deltaTime = 0;
static int height = 20, width = 80;
static size_t delta = 0;
static HFONT hFont = NULL;
static HWND hEditMultiline = NULL;
static GreatCircleDistance gcd;
static Vincenty vincenty;
switch (message)
{
case WM_INITDIALOG:
hFont = CreateFont(
-MulDiv(7, GetDeviceCaps(GetDC(hDlg), LOGPIXELSY), 72),
0, 0, 0, FW_BOLD, FALSE, FALSE, FALSE,
DEFAULT_CHARSET, OUT_DEFAULT_PRECIS,
CLIP_DEFAULT_PRECIS, DEFAULT_QUALITY,
FIXED_PITCH | FF_MODERN,
TEXT("Courier New")
);
hEditMultiline = CreateWindowEx(
WS_EX_CLIENTEDGE, // Extended style for sunken border
TEXT("EDIT"), // Class name
TEXT(""), // Initial text (can be blank)
WS_CHILD | WS_VISIBLE | WS_VSCROLL | ES_LEFT | ES_MULTILINE | ES_READONLY,
10, 0, 10 * width, 30 * height, // Position and size
hDlg, // Parent window handle
(HMENU)IDC_EDIT_MULTILINE, // Unique control ID
hInst, // Application instance
NULL // Extra parameter
);
SendMessage(hEditMultiline, WM_SETFONT, (WPARAM)hFont, TRUE);
return (INT_PTR)TRUE;
case WM_COMMAND:
if (LOWORD(wParam) == IDCANCEL)
{
EndDialog(hDlg, LOWORD(wParam));
return (INT_PTR)TRUE;
}
if (LOWORD(wParam) == IDC_BUTTON_COMPUTE)
{
outputText += L"Range in Yards Versus Curvature of the Earth in Feet\r\n\r\n";
outputText += L"Yards\t 0\t 100\t 200\t 300\t 400\t 500\t 600\t 700\t 800\t 900\r\n";
for (double row = 1000.0; row <= 40000.0; row += 1000.0)
{
swprintf_s(line, L"%6.1lf\t", row);
outputText += line;
for (double col = 0.0; col <= 900.0; col += 100.0)
{
double r = row + col;
double distance = 3.0 * gcd.Interpolate(r);
swprintf_s(line, L"%4.1lf\t", distance);
outputText += line;
}
outputText += L"\r\n";
}
outputText += L"Yards\t 0\t 100\t 200\t 300\t 400\t 500\t 600\t 700\t 800\t 900\r\n";
SetWindowText(hEditMultiline, outputText.c_str());
return (INT_PTR)TRUE;
}
}
return (INT_PTR)FALSE;
}
#pragma once
#include "Vincenty.h"
struct PointRaR
{
double ra; // curvature of the Earth correction in yards
double r; // flat Earth distance (chord) in yards
PointRaR(double ra, double r)
{
this->ra = ra;
this->r = r;
}
};
class GreatCircleDistance
{
private:
Vincenty vincenty;
std::vector<PointRaR> pts;
bool binarySearch(double x, int& lt, int& rt);
public:
friend int compare(
const PointRaR& lt,
const PointRaR& rt)
{
if (lt.r < rt.r)
return -1;
if (lt.r > rt.r)
return +1;
return 0;
};
// R is in yards, returns in yards
// Curvature of the Earth correction
double Interpolate(double R);
// construction of the curvature of the Earth table
GreatCircleDistance();
};
#pragma once
class Vincenty
{
public:
static double Re; // Radius of Earth in meters
double deltaSigma(
double phi1, double lambda1,
double phi2, double lambda2);
double distance(
double phi1, double lambda1,
double phi2, double lambda2);
double distanceConstantLambda(
double phi1, double phi2);
double distanceConstantXY(double z);
double x(double phi, double lambda);
double y(double phi, double lambda);
double z(double phi);
double phi(double z);
double lambda(double x, double y);
};
#include "pch.h"
#include "GreatCircleDistance.h"
bool GreatCircleDistance::binarySearch(
double x, int& lt, int& rt)
{
int n = static_cast<int>(pts.size()), L = 0, M, R = n - 1;
Label10:
if (x == pts[L].r)
{
lt = rt = L;
return true;
}
if (x == pts[R].r)
{
lt = rt = R;
return true;
}
if (x > pts[L].r && x < pts[R].r && R - L == 1)
{
lt = L;
rt = R;
return true;
}
if (x > pts[L].r && x < pts[R].r)
{
M = (L + R) / 2;
if (x > pts[M].r)
{
L = M;
goto Label10;
}
if (x < pts[M].r)
{
R = M;
goto Label10;
}
}
lt = rt = -1;
return false;
}
double GreatCircleDistance::Interpolate(double R)
{
int lt, rt;
if (binarySearch(R, lt, rt))
{
double x0 = pts[lt].ra, x1 = pts[rt].ra;
double y0 = pts[lt].r, y1 = pts[rt].r;
double deltaX = x1 - x0, deltaY = y1 - y0;
double ra = deltaX * (R - y0) / deltaY + x0;
return ra;
}
return -1.0;
}
GreatCircleDistance::GreatCircleDistance()
{
double deltaPhi = 0.000001, phi1 = 0.0, phi2 = deltaPhi, delta;
double deltaRa, d0 = vincenty.z(phi1), d1, r, ra;
int cnt = 0;
pts.push_back(PointRaR(0.0, 0.0));
while (cnt < 10000)
{
d1 = vincenty.z(phi2);
if (d0 >= d1)
{
delta = d0 - d1;
deltaRa = d0 * d0 - d1 * d1;
}
else
{
delta = d1 - d0;
deltaRa = d1 * d1 - d0 * d0;
}
r = delta;
ra = sqrt(deltaRa);
ra = r >= ra ? r - ra : ra - r;
pts.push_back(PointRaR(1.0936 * r, 1.0936 * ra));
phi2 += deltaPhi;
cnt++;
}
}
#include "pch.h"
#include "Vincenty.h"
double Vincenty::Re = 6378137.0; // radius of Earth in meters
double Vincenty::deltaSigma(
double phi1, double lambda1,
double phi2, double lambda2)
{
double deltaPhi = phi1 - phi2, deltaLambda = lambda1 - lambda2;
double cosPhi1 = cos(phi1), cosPhi2 = cos(phi2);
double sinPhi1 = sin(phi1), sinPhi2 = sin(phi2);
double cosDeltaLambda = cos(deltaLambda), sinDeltaLambda = sin(deltaLambda);
double numer1 = cosPhi2 * sinDeltaLambda;
double numer2 = cosPhi1 * sinPhi2 - sinPhi1 * cosPhi2 * cosDeltaLambda;
double numer = sqrt(numer1 * numer1 + numer2 * numer2);
double denom = sinPhi1 * sinPhi2 + cosPhi1 * cosPhi2 * cosDeltaLambda;
return atan2(numer, denom);
}
double Vincenty::distance(
double phi1, double lambda1,
double phi2, double lambda2)
{
return Re * deltaSigma(phi1, lambda1, phi2, lambda2);
}
double Vincenty::distanceConstantLambda(
double phi1, double phi2)
{
return Re * distance(phi1, 0.0, phi2, 0.0);
}
double Vincenty::distanceConstantXY(double z)
{
return Re * distance(phi(0.0), 0.0, phi(z), 0.0);
}
double Vincenty::x(double phi, double lambda)
{
return Re * sin(phi) * cos(lambda);
}
double Vincenty::y(double phi, double lambda)
{
return Re * sin(phi) * sin(lambda);
}
double Vincenty::z(double phi)
{
return Re * cos(phi);
}
double Vincenty::phi(double z)
{
return acos(z / sqrt(Re));
}
double Vincenty::lambda(double x, double y)
{
return asin(y / sqrt(x * x + y * y));
}
Numerical Explorations 