Jump to content
Sign in to follow this  
Dana

function or math problem

Recommended Posts

Dana

I'm working with two scripts, trying to put the functionality of script 2 into script 1. Script 2 reads a .gpx (GPS track) file containing track points with latitude and longitude information, looks for stops, and compares the stop points with another database file containing known locations, finding the closest match. Both scripts contain multiple UDF's. When I wrote script 2 with this in mind, I wrote script 2 with nothing but UDFs; the first thing script 2 does (after declaring some global variables) is to call the first UDF, which is the main routine, which itself calls other UDFs.

I then cut and pasted the entire contents of script 2 into script 1, checking that there are no variable name conflicts. The only changes I made to the script 2 functions were to return back to the main script 1 GUI rather than exiting at the end Depending on user interaction, the UDF which was the main part of script 2 is called, which then calls the other UDFs that were in script 2 and are now in script 1. The problem is, it doesn't work.

I don't expect people to wade through all of both scripts, but the problem seems to originate in one of the smaller UDFs from script 2. This function takes two pairs of latitude and longitude, and calculates the distance between the two:

Func _LLDist($lat1, $lon1, $lat2, $lon2)
;using Sinnott's formula for great circle distance between two points
;first convert degrees to radians for trig functions, using ACos(-1) to derive pi
ConsoleWrite("1197 " & $lat1 & ", " & $lon1 & "; " & $lat2 & ", " & $lon2)
$lat1 = $lat1 * (ACos(-1) / 180)
$lon1 = $lon1 * (ACos(-1) / 180)
$lat2 = $lat2 * (ACos(-1) / 180)
$lon2 = $lon2 * (ACos(-1) / 180)
Local $R = 3440.064794816415 ; Earth radius in nautical miles (6371000 meters)
Local $dlon = $lon2 - $lon1
Local $dlat = $lat2 - $lat1
Local $a = (Sin($dlat / 2)) ^ 2 + Cos($lat1) * Cos($lat2) * (Sin($dlon / 2)) ^ 2
Local $c = 2 * ASin(Sqrt($a))
Local $d = $R * $c
ConsoleWrite(" 1208 dist - " & $d & @CRLF)
Return $d
EndFunc ;==>_LLDist

If I run script 2 by itself, I get correct console output:

1197 41.38123769313097, -72.50517829321325; 41.38124062679708, -72.50517854467034 1208 dist - 0.000176502556754481

If I run script 1 on the same data, I get:

1197 41.38123769313097, -72.50517829321325; 41.38124062679708, -72.50517854467034 1208 dist - -1.#IND

As you can see, the input data is the same but the output is not.

Additionally, after this happens, other data in some arrays seems to get corrupted, and, for example, an array value that's supposed to range from 1-20,000 or so shows values like -9223372036854762545. A divide by zero causing a memory corruption? Or something else?

Any ideas?

Edited by Dana

Share this post


Link to post
Share on other sites
Melba23

Dana,

Interesting. :huh2:

Have you tried seeing what the other mathematical manipulations used during the function return in the 2 cases as well as the final result? That might pin down the exact line on which the error occurs. :huh:

Just one other point - although I do not see why it should affect the result when using the same data - thisline:

Local $a = (Sin($dlat / 2)) ^ 2 + Cos($lat1) * Cos($lat2) * (Sin($dlon / 2)) ^ 2

could do with some more parentheses to reduce its dependence on operator precedence. ;)

M23


Any of my own code posted anywhere on the forum is available for use by others without any restriction of any kind._______My UDFs:

Spoiler

ArrayMultiColSort ---- Sort arrays on multiple columns
ChooseFileFolder ---- Single and multiple selections from specified path treeview listing
Date_Time_Convert -- Easily convert date/time formats, including the language used
ExtMsgBox --------- A highly customisable replacement for MsgBox
GUIExtender -------- Extend and retract multiple sections within a GUI
GUIFrame ---------- Subdivide GUIs into many adjustable frames
GUIListViewEx ------- Insert, delete, move, drag, sort, edit and colour ListView items
GUITreeViewEx ------ Check/clear parent and child checkboxes in a TreeView
Marquee ----------- Scrolling tickertape GUIs
NoFocusLines ------- Remove the dotted focus lines from buttons, sliders, radios and checkboxes
Notify ------------- Small notifications on the edge of the display
Scrollbars ----------Automatically sized scrollbars with a single command
StringSize ---------- Automatically size controls to fit text
Toast -------------- Small GUIs which pop out of the notification area

 

Share this post


Link to post
Share on other sites
jchd

That reminds me of ticket #2068

You should use:

Global Const $Pi = 3.1415926535897932385

only once at code init, then use $Pi to keep away from this bug. Perhaps your more complex code (2 UDFs) is invoking ACos(-1) more than 6 times and triggers the issue.

Then I tried a quick'n'dirty Mathematica session to compute the actual result we should get if things were close to perfect:

In[2]:= SetPrecision[{coords={lat1,lon1,lat2,lon2}=(# Pi/180) & /@{41.38123769313097,-72.50517829321325,41.38124062679708,-72.50517854467034},
{dlat=lat2-lat1,dlon=lon2-lon1},
{a=Sin[dlat/2]^2+Cos[lat1] Cos[lat2] Sin[dlon/2]^2,
c=2 ArcSin[Sqrt[a]]},
d=3440.064794816415 c},
20]
Out[2]= {
{41.381237693130970001,-72.505178293213248253,41.381240626797080040,-72.505178544670343399},{0.72223884629551826642,-1.2654540859620937443,0.72223889749765102231,-1.2654540903508479932},
{5.1202132755889806504*10^-8, -4.3887542489073894103*10^-9},
{6.5812556737731055044*10^-16,  5.1307916246026236769*10^-8},
0.00017650255637334405522}

As you can see I requested 20 digits of guaranteed precision for the whole block.

The fact that the points are very close does not seem to impact result accuracy too much, despite here the value of a and c (your $a and $c) are fairly close to 0. Thanks to the "haversine" trick and larger than average precision we don't experience numerical instability in this case, but it might show up silently for closer points or when using standard IEEE floating-point like AutoIt does.

To avoid spurious instability, you can consider points too close as points in a plane instead of a sphere, in which case the error made is unsignificant. I kindly live it to you to compute which conditions make machine-precision rounding error equivalent to sphere-to-plane mapping systematic error.

Also beware of [close to] antipodal points!

Edited by jchd

This wonderful site allows debugging and testing regular expressions (many flavors available). An absolute must have in your bookmarks.
Another excellent RegExp tutorial. Don't forget downloading your copy of up-to-date pcretest.exe and pcregrep.exe here
RegExp tutorial: enough to get started
PCRE v8.33 regexp documentation latest available release and currently implemented in AutoIt beta.

SQLitespeed is another feature-rich premier SQLite manager (includes import/export). Well worth a try.
SQLite Expert (freeware Personal Edition or payware Pro version) is a very useful SQLite database manager.
An excellent eBook covering almost every aspect of SQLite3: a must-read for anyone doing serious work.
SQL tutorial (covers "generic" SQL, but most of it applies to SQLite as well)
A work-in-progress SQLite3 tutorial. Don't miss other LxyzTHW pages!
SQLite official website with full documentation (may be newer than the SQLite library that comes standard with AutoIt)

Share this post


Link to post
Share on other sites
Spiff59

I posted a bugtracker a while back describing this:

http://www.autoitscript.com/trac/autoit/ticket/2068

The work-around might be replacing the following:

$lat1 = $lat1 * (ACos(-1) / 180)
$lon1 = $lon1 * (ACos(-1) / 180)
$lat2 = $lat2 * (ACos(-1) / 180)
$lon2 = $lon2 * (ACos(-1) / 180)

With a global at the top of your script like:

$iACos180 = ACos(-1) / 180

And then:

$lat1 = $lat1 * $iACos180
$lon1 = $lon1 * $iACos180
$lat2 = $lat2 * $iACos180
$lon2 = $lon2 * $iACos180

If you've other "ACos(-1)" statements without the division by 180 then maybe just make a global for "ACos(-1)".

But if you execute six "ACos(-1)" commands in row... Blam!

Edit: As Johnny Carson used to say "Weird stuff!"

Edit2: The ticket does claim it was fixed in 3.3.9.0... and testing shows it has been. So leave your code intact and run the beta if you like.

Edited by Spiff59

Share this post


Link to post
Share on other sites
Dana

Strange, I hadn't heard of the Acos(-1) bug. I'll check if I'm using the pre or post bugfix version. I had considered starting with $pi = Acos(-1) and then just using $pi thereafter, will try that tomorrow. I am calling the _LLDist function repeatedly, up to maybe a hundred times (I first sort the big database of 20K entries to quickly find points within 1 minute of latitude, usually 60 or so but it varies, then use the function check the actual distance for each of those points, within 1 mile vertially but could be thousands of miles away horizontally, and find the closest).

The odd thing is that the script 2 version uses exactly the same logic, same number of iterations, but doesn't show the problem. Script 1 invokes script 2 only once, after which it's the same.

Share this post


Link to post
Share on other sites
jchd

What you describe seems a good job for the Rtree extension of SQLite.

Don't forget the limited accuracy of FP I mentionned! For really close points, a small perturbation in lattitudes or longitudes can result in a large variation in the (meaningless in that case) result.

Edited by jchd

This wonderful site allows debugging and testing regular expressions (many flavors available). An absolute must have in your bookmarks.
Another excellent RegExp tutorial. Don't forget downloading your copy of up-to-date pcretest.exe and pcregrep.exe here
RegExp tutorial: enough to get started
PCRE v8.33 regexp documentation latest available release and currently implemented in AutoIt beta.

SQLitespeed is another feature-rich premier SQLite manager (includes import/export). Well worth a try.
SQLite Expert (freeware Personal Edition or payware Pro version) is a very useful SQLite database manager.
An excellent eBook covering almost every aspect of SQLite3: a must-read for anyone doing serious work.
SQL tutorial (covers "generic" SQL, but most of it applies to SQLite as well)
A work-in-progress SQLite3 tutorial. Don't miss other LxyzTHW pages!
SQLite official website with full documentation (may be newer than the SQLite library that comes standard with AutoIt)

Share this post


Link to post
Share on other sites
Dana

Well, glory be! Globally defining $pi = Acos(-1) at the beginning of the main scrlipt and using $pi in the calculation function fixed everything (except the zillions of ConsoleWrites scattered throughout the file from the debug attempt). Presumably upgrading from 3.3.6.1 might have fixed it, too, but 3.3.9 is still beta (for a year now?).

Still can't see why it behaved differently between the two implementations, unless it were a memory usage thing and I was already on the hairy edge (both of the orignal scripts have a fair amount of array manipulation).

jhcd, this is a brute force approach, but it's adequate for the task. The closest two points that I need to differentiate are some 500' apart, and I don't require an absolute distance, just which point from the database is closed to the GPS derived point. Actually, what takes longest is reading (3 seconds) the [up to] 20,000 line input text file (the file extension is .xls but it's really a .csv file, which I have no control over) and sorting it (_ArraySort, 4 seconds) by latitude for the first rough filter. Using a "top half, bottom half" iterative sort for points within 1 mile north/south and only using the _LLDist function on the 50-100 pre-sort results is quite quick.

Share this post


Link to post
Share on other sites

Create an account or sign in to comment

You need to be a member in order to leave a comment

Create an account

Sign up for a new account in our community. It's easy!

Register a new account

Sign in

Already have an account? Sign in here.

Sign In Now
Sign in to follow this  

×