CCS C Software and Maintenance Offers
FAQFAQ   FAQForum Help   FAQOfficial CCS Support   SearchSearch  RegisterRegister 

ProfileProfile   Log in to check your private messagesLog in to check your private messages   Log inLog in 

CCS does not monitor this forum on a regular basis.

Please do not post bug reports on this forum. Send them to CCS Technical Support

using float64 and cos; sin; acos and asin
Goto page Previous  1, 2
 
Post new topic   Reply to topic    CCS Forum Index -> General CCS C Discussion
View previous topic :: View next topic  
Author Message
Ttelmah



Joined: 11 Mar 2010
Posts: 19559

View user's profile Send private message

PostPosted: Fri May 04, 2012 2:40 am     Reply with quote

It is important to understand 'limitations' of everything.

The simple flat space distance between points solution, is quick and easy _provided the points are close_. However it goes severely wrong, when the distances get greater. The advantage of the Haversine formula is that it works pretty well in both case.
All the formulae go wrong somewhere because of simple things like the fact that the Earth is not a sphere.

However just to show how good the Haversine solution is, I've coded it with a couple of tweaks (one thing that improves accuracy a _lot_ is not using the pow function to square the values - simply multiply them together):

Code:

#define torad (0.01745329251994)
#define toraddiv2 (0.00872664625997)

float32 CalculateDistance(float32 Lat1,float32 Lon1, float32 Lat2, float32 Lon2) {
    int32 Diameter = 12742018; // Earth's mean diameter in metres
    float32 DeltaLat, DeltaLon, temp1, temp2;
    DeltaLat= (Lat2 - Lat1) * toraddiv2;
    DeltaLon = (Lon2 - Lon1) * toraddiv2;
    temp1 = sin(DeltaLat);
    temp1*=temp1;
    temp2 = sin(DeltaLon);
    temp2*=temp2;
    temp2*=cos(lat2*torad);
    temp2*=cos(lat1*torad);
    temp1+=temp2;
    temp2=sqrt(temp1);
    temp1=asin(temp2);
    return(temp1*Diameter); // Return calculated distance in metres
}


This gives the distance between two coordinates in metres.

For the original test using 0.1 degrees apart, it gives 11118.4385m, which given that this is based upon a mean Earth diameter of 12742018, which would imply a circumference of 40030230m, implying a distance for 0.1 degrees of 1/3600th this - which gives 11119.508m, is pretty remarkably accurate. Even better the same formula using two points 90degrees apart, gives 10007560m, which is equally good.

I'd suggest that this is probably as close as it is possible to get on the PIC. It might be worth seeing if it improves at all with float64 (I haven't spent the time to do this), butt it can only improve by a meter or two, since it is good already!.... Smile

Best Wishes
longtu



Joined: 15 Jun 2009
Posts: 6

View user's profile Send private message

PostPosted: Fri May 04, 2012 3:10 am     Reply with quote

Dear Ttemah,

Could you please test your code with distance less than 50m. Example:
- Point 1: (0.0002,0)
- Point 2: (0.0000,0)

I want to detect with between 10 - 50m only.

Thanks you advance.
_________________
----------------------------------------------------------
[email protected]
longtu



Joined: 15 Jun 2009
Posts: 6

View user's profile Send private message

PostPosted: Fri May 04, 2012 3:21 am     Reply with quote

Dear Ttemah,

I have just test your code(using CCS V4.114). The result is WRONG.

_________________
----------------------------------------------------------
[email protected]
Ttelmah



Joined: 11 Mar 2010
Posts: 19559

View user's profile Send private message

PostPosted: Fri May 04, 2012 4:07 am     Reply with quote

What PIC.....

I am using a 18F4550. It gives correct results on this.

I have a 'nasty' feeling something is totally screwed in the trig functions on the PIC24 family...

For <50m, I'd go to simple flat space geometry ever time. The difference between a 'flat' solution, and a spherical solution, only starts to become even remotely noticeable for distances over a few km, and the maths is much more work. I wouldn't expect to have a reasonable hope of working meaningfully at scales of a few metres. when you think that sin(0.0002 degrees), only differs from 'zero', by 3, in the sixth decimal place....
You wouldn't start using spherical geometry, to measure the distance across the floor of your living room, from the lengths of the walls, would you?.

One suggestion, would be to test the two 'delta' values, using FABS, and if they are both less than (say) 1, perform a solution using the simple Pythagorean approach. Save arithmetic for short distances, and probably give better overall results.

However it gives a figure of 20.9573m for 0.0002, which is still only 1.2m wrong. Hopefully quite impressive. Smile

Best Wishes
jeremiah



Joined: 20 Jul 2010
Posts: 1358

View user's profile Send private message

PostPosted: Fri May 04, 2012 5:27 am     Reply with quote

I haven't checked the most recent version of the compiler, but some of the 4.12x set of versions had a bug with atan2, or specifically with the setup prior to calling atan2 (which haversine uses). I found that they clobbered one of the working registers used to setup for the call. I worked around this by going into math.h and copying the atan2 function, pasting it in my file and renaming it to my_atan2 and calling that. I don't know if that is your issue, but that is what I faced. I reported it to CCS so they may have fixed that and it is some other issue.
longtu



Joined: 15 Jun 2009
Posts: 6

View user's profile Send private message

PostPosted: Fri May 04, 2012 5:45 am     Reply with quote

I am using PIC18F26K22. Do you think my problem cause by PIC18F26K22?
_________________
----------------------------------------------------------
[email protected]
Ttelmah



Joined: 11 Mar 2010
Posts: 19559

View user's profile Send private message

PostPosted: Fri May 04, 2012 7:53 am     Reply with quote

jeremiah wrote:
I haven't checked the most recent version of the compiler, but some of the 4.12x set of versions had a bug with atan2, or specifically with the setup prior to calling atan2 (which haversine uses). I found that they clobbered one of the working registers used to setup for the call. I worked around this by going into math.h and copying the atan2 function, pasting it in my file and renaming it to my_atan2 and calling that. I don't know if that is your issue, but that is what I faced. I reported it to CCS so they may have fixed that and it is some other issue.


If you look at the version I posted, I don't use atan2. Mine is based on the 2r arcsin(sqrt(h)) version instead, 'pre-solving' 2r as an integer number of meters. However I'd not be surprised if the fault you have found applies to other trig functions....

Best Wishes
Display posts from previous:   
Post new topic   Reply to topic    CCS Forum Index -> General CCS C Discussion All times are GMT - 6 Hours
Goto page Previous  1, 2
Page 2 of 2

 
Jump to:  
You cannot post new topics in this forum
You cannot reply to topics in this forum
You cannot edit your posts in this forum
You cannot delete your posts in this forum
You cannot vote in polls in this forum


Powered by phpBB © 2001, 2005 phpBB Group