|
HELP!! I am not the type to give up but I am way past wits end! This is for work ....
The goal was simple...
Given the starting and ending latitude and longitude calculate the direction in compass degrees of travel as the crow flies to get to the destination location.
For example if I wanted to travel from say Sacramento to San Diego CA, I might have to travel I want to know how many degrees and what heading to get to Denver.
But ... Every time I search for a solution, the answer evades me and my hopes are crushed as I run into a combination of the following problems over and over again.
1. The site launches into a diatribe dozens of pages long written by some author who seems to delight in hearing or reading his incomprehensible diatribe till my eyes glaze over and I am over-whelmed.
2. The site starts to answer the question but then half way through it changes direction and ends up answering an entirely different question. Example: Half way through they stop and answer a different question like finding the distance of travel between the two points or the longitude and latitude of the mid-point and then just stops right there and never gets to answering the problem.
3. The week before last this start out saying they will show the formula for this and stating the question completely but when they get to the formula and example page, their page is missing the controls for entering the destination latitude, and also the formula they use leave this out too, and of course the answer is not forthcoming.
So finally last Friday I find a site that says if you want the formula for calculating this than look no further and they give me the code.
Here it is below and it works between say Paradise CA, and Denver Colorado, and the other way around but when I asked what was the direction from Paradise, CA, to San Diego CA, instead of say 220 Degrees (SW Generally), it comes out and says it's 44 Degrees.
I am not the type to quit, which was what some friends said when I mention that I have failed and am defeated, I can't solve this and I've spent many times over the amount of time it should have taken to solve this but I am defeated completely.
here's the code. Here's your chance any smart folks out there to solve this, but for me I am the retard, and I can't figure it out and I give up. Hopefully someone is up to the challenge. It's over for me. I will have to hunker down to pencil and paper and make some rube-goldberg piece of code because I cannot find an elegant solution out there and I am beyond pissed about it.
Here is the defines for it ....
#define MILES 3963.1
#define NauticalMiles 3443.9
#define KILOMETERS 6378
#define ParadiseLat 39.743
#define ParadiseLon -121.605
#define SanDiegoLat 32.723
#define SanDiegoLon -117.165
#define ColoSpringsLat 38.835
#define ColoSpringsLon -104.823
#define PI 3.14159265358979323846
Here is the function the main function below calls...
double degreesToRadian(double x)
{
return (M_PI * (x) / 180.);
}
Here is the main function below, it works for say Colorado Springs to Paradise CA and vice versa but falls flat on Paradise CA, to San Diego. It returns a rediculous answer of 44 degrees. Good luck, if you can find a solution you are far smarter than I am which right now doesn't sound like much of a challenge.
double latitudeLongitudeLocationsToAngleOfTravelInDegrees(double latRad1, double lonRad1, double latRad2, double lonRad2)
{
double theRadians = atan2(cos(latRad1)*sin(latRad2)-sin(latRad1)*cos(latRad2)*cos(lonRad2-lonRad1), sin(lonRad2-lonRad1)*cos(latRad2));
double theDegreesInRad = radianToDegrees(theRadians);
theDegreesInRad = (theDegreesInRad * PI / 180) + 360;
int theDegrees = radianToDegrees(theDegreesInRad);
theDegrees = theDegrees % 360;
return theDegrees;
}
Been there, done that, forgot why!
|
|
|
|
|
Hi Garry,
I am not familiar with long/lat but I do remember some high school math and the internet is helpful too;
this link[^] looks trustworthy, scroll down to the "bearing" chapter, IMO the first formula giving theta is what you need.
I would suggest you handle everything in radians (that's the unit sin, cos, atan2, ... will accept and
return anyway), only as a last step you could convert back to degrees.
Remark: this is the C# forum, your code is C code, that is something entirely different!
PS: yes, some sites first calculate distance, then derive bearing from it. Looks like a detour to me.
Luc Pattyn [Forum Guidelines] [My Articles]
I use ListBoxes for line-oriented text, and PictureBoxes for pictures, not drawings.
modified on Friday, June 10, 2011 12:28 PM
|
|
|
|
|
It could be C#. Just yesterday while working on my sieve I played with some # defines... which I later removed.
|
|
|
|
|
Hi PIEBALD,
C/C++ preprocessing supports defines (or macros) as key-value pairs,
whereas C# only accepts only key names, the value is fixed at "true"
hence it couldn't possibly be C/C++ code.
PIEBALDconsult wrote: on my sieve I played with some # defines... which I later removed
sure, sieves are intended for removing things. 
|
|
|
|
|
I pass my C# through the C-Preprocessor... just to prove a point.
|
|
|
|
|
PIEBALDconsult wrote: I pass my C# through the C-Preprocessor...
One could try and do that, I did consider it when I switched from Java to C#, but never actually did it.
|
|
|
|
|
Thanks kindly. I'll have a look at the link. The reason I posted it here is that CodeProject does not have C as a choice in the dropdown. The closest other choice for forum location was C++/MFC. Perhaps it should be called C,C++/MFC
I forgot to say I have this sick feeling in my gut as I realized the answer is not so simple if the user goes north or south over the pole, their north south heading will differ from Start to Destination.
I finally decided to drop down to saying X distance from West,East and North,South or vice versa but now it looks like I will have a four part answer.
I wish this weren't for work, but such is not the case, but like so many things about work, if the demand of the employer is not met however difficult or impossible the demand is, it can be the employment line for me.
Been there, done that, forgot why!
|
|
|
|
|
Garry Freemyer wrote: looks like I will have a four part answer.
No, I don't expect so; if you have two points on a circle in a 2D plane, you can apply sines and consines and everything no matter where the points are, it always works out; it suffices the coordinates (whether x/y or long/lat) uniquely define the location.
Garry Freemyer wrote: does not have C as a choice
I know, but C stuff is handled under C++/MFC over here.
Luc Pattyn [Forum Guidelines] [My Articles]
I use ListBoxes for line-oriented text, and PictureBoxes for pictures, not drawings.
modified on Friday, June 10, 2011 12:29 PM
|
|
|
|
|
That site you included is the first one I tried but the function below that is modeled after the formula there always returns zero no matter what I do when passed in the values. Since I do not understand the formula I cannot understand why it's not working.
Input values are as follows ....
double lat1 = degreesToRadian(ParadiseLat);
double lon1 = degreesToRadian(ParadiseLon);
double lat2 = degreesToRadian(SanDiegoLat);
double lon2 = degreesToRadian(SanDiegoLon);
double degreesToRadian(double x)
{
return (M_PI * (x) / 180.);
}
double radianToDegrees(double x)
{
return x * 57.2957795;
}
The function below always returns zero no matter when I expect it to be returning something close to south south west, absolutely not north east. Expected at least something more than 220 degrees, certainly not 44 degrees. I noticed that it appeared to be a mirror image of what I wanted, but I don't know how to fix it.
I know some are wondering why I don't understand this "Simple High School Trig", and the answer is I never got to take Trigonometry. The reason for that is a tale of Gross Misdiagnosis and Gross Negligence by school officials and Medical personal. Its a drama I prefer not to get into here.
Someone said San Diego was South East of Paradise CA .. well maybe so but when I looked at the map it looked South West.
*You will note that the return value is truncated to an integer in the final step just before the modulus % operator. This is because when I try this with anything but an integer, be it signed or unsigned, I get a compiler error. Apparently % only works for whole numbers.
int latitudeLongitudeLocationsToAngleOfTravelInDegrees(double latRad1, double lonRad1, double latRad2, double lonRad2)
{
//double theRadians = atan2(cos(latRad1)*sin(latRad2)-sin(latRad1)*cos(latRad2)*cos(lonRad2-lonRad1), sin(lonRad2-lonRad1)*cos(latRad2));
//double theDegreesInRad = radianToDegrees(theRadians);
//theDegreesInRad = (theDegreesInRad * PI / 180) + 360;
//int theDegrees = radianToDegrees(theDegreesInRad);
//theDegrees = theDegrees % 360;
double theAngleInRad = atan2(cos(latRad1)*sin(latRad2)-sin(latRad1)*cos(latRad2)*cos(lonRad2-lonRad1), sin(lonRad2-lonRad1)*cos(latRad2));
int theDegrees = radianToDegrees(theAngleInRad); // *See comment above.
return (theDegrees + 360) / 360;
// return theDegrees;
}
Been there, done that, forgot why!
modified on Tuesday, January 13, 2009 12:10 PM
|
|
|
|
|
Garry Freemyer wrote: return (theDegrees + 360) / 360;
can't be right, you want modulus, hence %
BTW: I don't like your 57.2957795 constant, what you want there is 180./M_PI so radianToDegrees inverts whatever degreesToRadian does as good as it possibly can.
Luc Pattyn [Forum Guidelines] [My Articles]
I use ListBoxes for line-oriented text, and PictureBoxes for pictures, not drawings.
modified on Friday, June 10, 2011 12:29 PM
|
|
|
|
|
Good catch. I can't believe I keep replacing the % with /.
I did the changes, and when I fed in longitude and latitude for Paradise CA to Colorado Springs Colorado, I got a wildly unexpected result.
Given the code fragments below the angle I was expecting was something like 57 degrees came back and said it was 1 degree! NORTH!!
I get the formula working for Paradise CA to Colorado Springs and it breaks down when I try for a Destination of San Diego, and if I fix San Diego, it breaks for Colorado Springs. I might as well be passing these figures into the random function. I would scream but you are trying to help me and it would not be polite.
...
#define ParadiseLat 39.743
#define ParadiseLon -121.605
#define SanDiegoLat 32.723
#define SanDiegoLon -117.165
#define ColoSpringsLat 38.835
#define ColoSpringsLon -104.823
...
double lat1 = degreesToRadian(ParadiseLat);
double lon1 = degreesToRadian(ParadiseLon);
double lat2 = degreesToRadian(ColoSpringsLat);
double lon2 = degreesToRadian(ColoSpringsLon);
....
double radianToDegrees(double x)
{
return x * 180./M_PI;
}
...
int latitudeLongitudeLocationsToAngleOfTravelInDegrees(double latRad1, double lonRad1, double latRad2, double lonRad2)
{
//double theRadians = atan2(cos(latRad1)*sin(latRad2)-sin(latRad1)*cos(latRad2)*cos(lonRad2-lonRad1), sin(lonRad2-lonRad1)*cos(latRad2));
//double theDegreesInRad = radianToDegrees(theRadians);
//theDegreesInRad = (theDegreesInRad * PI / 180) + 360;
//int theDegrees = radianToDegrees(theDegreesInRad);
//theDegrees = theDegrees % 360;
double theAngleInRad = atan2(cos(latRad1)*sin(latRad2)-sin(latRad1)*cos(latRad2)*cos(lonRad2-lonRad1), sin(lonRad2-lonRad1)*cos(latRad2));
int theDegrees = radianToDegrees(theAngleInRad);
return (theDegrees + 360) % 360;
// return theDegrees;
}
Been there, done that, forgot why!
|
|
|
|
|
Hi,
I would prefer you show actual code, "..." would not compile, function degreesToRadian now has vanished;
also please make life easier to yourself and us CPians, by dropping the statements you commented out, and by putting the code in PRE tags (use the CODE BLOCK button) that would make it more readable, like so:
#define ParadiseLat 39.743
#define ParadiseLon -121.605
#define SanDiegoLat 32.723
#define SanDiegoLon -117.165
Garry Freemyer wrote: I get the formula working for Paradise CA to Colorado Springs and it breaks down when I try for a Destination of San Diego, and if I fix San Diego, it breaks for Colorado Springs.
I can't imagine what that could possibly mean; the formulae and the code either are correct or they are not, how could you "fix" them for a specific city?
If you need to debug, I recommend you turn the atan2 line into three lines (and insert some Console.WriteLine statements), so you actually see the two values that you are feeding into atan2.
Luc Pattyn [Forum Guidelines] [My Articles]
I use ListBoxes for line-oriented text, and PictureBoxes for pictures, not drawings.
modified on Friday, June 10, 2011 12:30 PM
|
|
|
|
|
The error is I get a 1 or a zero back from the function calls.
I do not understand where degreestoRadian has gone. I specifically included it.
I agree with the idea of removing commented out code and to format it in a code block.
I am out of time. I cannot solve the problem the way my boss would prefer and so now I am going to say the destination is X units West, or East and Y units North or South of the starting point.
Some day when I am ready I will probably tackle it again, but this problem has become like gambling, I keep thinking I am going to hit the jackpot if I pull the lever one more time, but the trend says its way past time to give up for now.
Been there, done that, forgot why!
|
|
|
|
|
Fixing means messing with it without knowing what I am doing. I certainly do not know what any of these functions do. I'm not a trigonomitrist, I am just some programmer that if given a formula, has a fairly decent chance of translating it to a computer language, but that doesn't mean I know the ins and outs.
Been there, done that, forgot why!
|
|
|
|
|
I had a fortune from a cookie that said something like: "When solving a problem, it's best to know the answer beforehand."
What is the answer you expect?
|
|
|
|
|
Boss was expecting something like San Diego is about 220 degrees (Approximately SW) of Paradise CA. I am afraid the answer is not going to be so simple. Its got to be Starting Heading in the form of ~0 .. 360, Ending heading something like ~0 .. 360 degrees).
Been there, done that, forgot why!
|
|
|
|
|
According to Google maps, San Diego, CA is South East of Paradise, CA
Clickety[^]
The course from the origin is 152 degrees. From magnetic north it's 136 based on this:
Clickety[^]
|
|
|
|
|
Yeah, right, next you'll say Reno is west of Los Angeles.
|
|
|
|
|
This is far from "simple geography compass problem", in fact there are a number of answers and ways to calculate distance/bearing to a point. The earth being roughly sperical means that there are a few ways to calculate a general "route" between two places on the surface.
If you want to deal with this the very simple way its Euclidean distance[^] that you want.
I think this[^] is more what you're after though!
|
|
|
|
|
You can simplify the problem by imagining the earth to be a rectangular array, with the x coordinate going from 0 to 360 and the y coordinate going from -90 to 90.
The two cities' coordinates map to two points on this rectangular array. If the difference between the two cities' longitudes is more than 180, add 360 to the lower longitude, extending the map.
Now you just have to get the angle between the line segment joining these two points, and a horizontal ray starting from the start city, and going to the right.
Many sites have the formula for the angle between two line segments, like http://www.vb-helper.com/howto_find_angles.html[^]
|
|
|
|
|
I had some spare time today, so I decided to tackle your geodesic problem. Here is my C# code
(I tend to use C#, the essence is in method Calc(); I trust you can do the same thing in C):
using System;
using System.Text;
namespace LongLat {
class Program {
class City {
public string Name;
public double LongitudeDegrees;
public double LatitudeDegrees;
public double LongitudeRadians;
public double LatitudeRadians;
public City(string name, double longitude, double latitude) {
this.Name=name;
this.LongitudeDegrees=longitude;
this.LatitudeDegrees=latitude;
this.LongitudeRadians=longitude*Math.PI/180.0;
this.LatitudeRadians=latitude*Math.PI/180.0;
}
public override string ToString() {
return Name+"("+LongitudeDegrees+", "+LatitudeDegrees+")";
}
}
static void Main(string[] args) {
City Brussels=new City("Brussels", 4.35, 50.85);
City Antwerp=new City("Antwerp", 4.40, 51.21);
City AntwerpEast=new City("East-of-Antwerp", 5.40, 51.21);
City AntwerpWest=new City("West-of-Antwerp", 3.40, 51.21);
City AntwerpNorth=new City("North-of-Antwerp", 4.40, 52.21);
City AntwerpSouth=new City("South-of-Antwerp", 4.40, 50.21);
City Paris=new City("Paris", 2.34, 48.85);
City Paradise=new City("Paradise", -121.605, 39.743);
City SanDiego=new City("San-Diego", -117.165, 32.723);
City ColoStrings=new City("Colorado-Springs", -104.823, 38.835);
Calc(Antwerp, AntwerpEast);
Calc(Antwerp, AntwerpWest);
Calc(Antwerp, AntwerpNorth);
Calc(Antwerp, AntwerpSouth);
Calc(Brussels, Antwerp);
Calc(Brussels, Paris);
Calc(Paradise, SanDiego);
Calc(Paradise, ColoStrings);
Console.ReadLine();
}
static void log(string s) {
Console.WriteLine(s);
}
static double Calc(City city1, City city2) {
double deltaLong=city2.LongitudeRadians-city1.LongitudeRadians;
double y=Math.Cos(city1.LatitudeRadians)*Math.Sin(city2.LatitudeRadians)-
Math.Sin(city1.LatitudeRadians)*Math.Cos(city2.LatitudeRadians)*Math.Cos(deltaLong);
double x=Math.Sin(deltaLong)*Math.Cos(city2.LatitudeRadians);
double bearing=Math.Atan2(y, x);
bearing=180.0*bearing/Math.PI;
log("Bearing from "+city1+" to "+city2+" = "+bearing.ToString("N2"));
return bearing;
}
}
}
and this is the output:
Bearing from Antwerp(4.4, 51.21) to East-of-Antwerp(5.4, 51.21) = 0.39
Bearing from Antwerp(4.4, 51.21) to West-of-Antwerp(3.4, 51.21) = 179.61
Bearing from Antwerp(4.4, 51.21) to North-of-Antwerp(4.4, 52.21) = 90.00
Bearing from Antwerp(4.4, 51.21) to South-of-Antwerp(4.4, 50.21) = -90.00
Bearing from Brussels(4.35, 50.85) to Antwerp(4.4, 51.21) = 85.03
Bearing from Brussels(4.35, 50.85) to Paris(2.34, 48.85) = -123.72
Bearing from Paradise(-121.605, 39.743) to San-Diego(-117.165, 32.723) = -61.63
Bearing from Paradise(-121.605, 39.743) to Colorado-Springs(-104.823, 38.835) = 1.37
The above Calc() method is basically identical to your original code, so I do not know how you got wrong results in the first place.
Hope this helps.
|
|
|
|
|
Wow that's a nice bit of work!
I got the same results as you, but since Colorado Springs is east and a bit south of Paradise and since east is 90 degrees, I was expecting something like maybe 107 degrees not 1.37. Your answer seems to be in the form of -180 ... +180 and I can't seem to get the results normalized to the range of 0 .. 360 where 0 and 360 are north.
If I take -180 and +180 to be north than I guess 1.37 is roughly south, but Colorado Springs is East not south so I am confused.
I guess what I do not understand is how to arrive at a range of 00 ... 360 as described above.
Been there, done that, forgot why!
|
|
|
|
|
Hi Garry,
no no, the results use zero for East, 90 for North, etc, same as regular two-dimensional geometry.
Check with the AntwerpEast etc. cities I included for testing.
If you want North=0, you should subtract 90 degrees from the end result (then west=90, ...).
If you want North=0 and East=90, you need to invert the sign.
If you want 0...360 rather than -180...180, you should do bearing=(bearing+360)%360 as you have done earlier yourself.
Luc Pattyn [Forum Guidelines] [My Articles]
I use ListBoxes for line-oriented text (not TextBoxes), and PictureBoxes for pictures (not drawings).
modified on Friday, June 10, 2011 12:30 PM
|
|
|
|
|
Thanks much!
I tried the third option on this result already and the problem I am having is that % will not work with floats, doubles but only with integers I get a compile error when I try it.
x % y will not compile unless the types for x and y are integer only.
So I tried multiplication to get the number rounded up to an integer, and then multiply the 360 by the same factor and then do that math and I am trying to work with that.
I suspect I may have to invert the signs first and then do the third step or I might not get the expected results I need....
0 .. 360 where 0 and 360 are north and 90 is East, and West is 270.
Been there, done that, forgot why!
|
|
|
|
|
Garry Freemyer wrote: x % y will not compile
Sorry, my mistake.
For 0=North, 90=East, 180=South, 270=West, this is what I would add after bearing=180.0*bearing/Math.PI;:
bearing=90.0-bearing;
if (bearing<0.0) bearing+=360.0;
Luc Pattyn [Forum Guidelines] [My Articles]
I use ListBoxes for line-oriented text (not TextBoxes), and PictureBoxes for pictures (not drawings).
modified on Friday, June 10, 2011 12:31 PM
|
|
|
|
|
General 
News 
Suggestion 
Question 
Bug 
Answer 
Joke 
Praise 
Rant 
Admin 
Submit an article or tip
