Geo functions PR

CustomFunctions/GeoBoxBounds.fmfn

@@ -0,0 +1,104 @@
+// GeoBoxBounds ( latitude ; longitude ; radiusInMeters )
+// version 2012-09-27, Jeremy Bante, http://scr.im/jbante 
+
+/**
+ * =====================================
+ * GeoBoxBounds ( latitude ; longitude ; radiusInMeters )
+ *
+ * PURPOSE:
+ *		Calculates the approximate coordinate bounds for an approximately square region
+ *		around a point on the surface of the Earth. This can be useful for finding data
+ *		within a certain distance of a reference point.
+ *
+ * RETURNS:
+ *		A return-delimited list of 4 numbers:
+ *		1. The south latitude boundary for the calculated region
+ *		2. The north latitude boundary for the calculated region
+ *		3. The west longitude boundary for the calculated region
+ *		4. The east longitude boundary for the calculated region
+ *
+ * PARAMETERS:
+ *		latitude: Degrees between -90 and 90 for center of geobox
+ *		longitude: Degrees between -180 and 180 for center of geobox
+ *		radiusInMeters: The distance in meters from the center of the region to the
+ *		closest points on the edges of the square bounding region.
+ *
+ * DEPENDENCIES: none
+ *
+ * EXAMPLE: Find records within a search radius
+ *		Set Variable [$geoBox; Value:GeoBoxBounds ( $latitude ; $longitude ; $radius)]
+ *		Enter Find Mode []
+ *		Set Field [Table::latitude; Value:GetValue ( $geoBox ; 1 ) & "..." & GetValue ( $geoBox ; 2 )]
+ *		Set Field [Table::longitude; Value:GetValue ( $geoBox ; 3 ) & "..." & GetValue ( $geoBox ; 4 )]
+ *		Perform Find []
+ *
+ * NOTES:
+ *		Use FileMaker's GetValue function to retrieve each value from the result.
+ *
+ *		When the geobox includes a pole, the function will return the appropriate values
+ *		to perform a find for coordinates enclosed by the lower latitude. For example:
+ *		[90, 88, -180, 180] for a region around the north pole.
+ *
+ *		When the geobox crosses the antimeridian (180º/-180º longitude), the function will
+ *		return a west boundary coordinate greater than the east boundary coordinate. In
+ *		this situation, a coordinate search within the geobox should be broken into
+ *		separate searches: one search from the west bound to 180º, and another from -180º
+ *		to the east bound.
+ *
+ *		The bounding region's approximation to a square shape is less accurate as
+ *		radiusInMeters grows large and as the region approaches either pole. However, the
+ *		region will always contain the circular region of the same radius.
+ *
+ *		Square bounding regions are more computationally convenient for finding data.
+ *		Developers needing an exact circular region can use a square region to constrain
+ *		the available data to a smaller found set, then use the GeoDistance function to
+ *		omit points between the circular and square neighborhoods of the reference point.
+ *
+ * HISTORY:
+ *		MODIFIED on 2012-09-27 by Jeremy Bante <http://scr.im/jbante> to account for
+ *		geoboxes that may span a pole or meridian.
+ *		CREATED on 2012-09-26 by Jeremy Bante <http://scr.im/jbante>.
+ *
+ * REFERENCES:
+ *		Radius of Earth:	http://www.wolframalpha.com/input/?i=radius+of+earth
+ * =====================================
+ */
+
+Let ( [
+	~earthRadius = 6367500;	// meters
+
+	// calculate latitude bounds
+	~latitudeOffset = Degrees ( radiusInMeters / ~earthRadius );
+	~north = latitude + ~latitudeOffset;
+	~south = latitude - ~latitudeOffset;
+	~wrapNorth = ~north > 90;
+	~wrapSouth = ~south < -90;
+	~northBound =
+		Case (
+			~wrapNorth ; 90;
+			~wrapSouth ; Max ( ~north ; -180 - ~south );
+			~north
+		);
+	~southBound =
+		Case (
+			~wrapSouth ; -90;
+			~wrapNorth ; Min ( ~south ; 180 - ~north );
+			~south
+		);
+
+	// calculate longitude bounds
+	~longitudeOffset = ~latitudeOffset * Cos ( Radians ( latitude ) );
+	~eastBound =
+		If ( ~wrapNorth or ~wrapSouth ; 180 ; /* Else */ longitude + ~longitudeOffset );
+	~eastBound = If ( ~eastBound > 180 ; ~eastBound - 360 ; /* Else */ ~eastBound );
+	~westBound =
+		If ( ~wrapNorth or ~wrapSouth ; -180 ; /* Else */ longitude - ~longitudeOffset );
+	~westBound = If ( ~westBound < -180 ; ~westBound + 360 ; /* Else */ ~westBound )
+];
+	List (
+		~southBound;
+		~northBound;
+		~westBound;
+		~eastBound
+	)
+)

CustomFunctions/GeoDistance.fmfn

@@ -0,0 +1,72 @@
+// GeoDistance ( latitude1 ; longitude1 ; latitude2 ; longitude2 )
+// version 2019-01-04, Jeremy Bante, http://scr.im/jbante
+
+/*******************************************************************************
+ * LocationDistance ( latitude1 ; longitude1 ; latitude2 ; longitude2 )
+ * Calculates the distance between two points on the surface of the Earth. This
+ * function assumes that coordinates are in the WGS84 coordinate system.
+ *
+ * @parameter latitude1: Latitude of the first point
+ * @parameter longitude1: Longitude of the first point
+ * @parameter latitude2: Latitude of the second point
+ * @parameter longitude2: Longitude of the second point
+ *
+ * @return The distance between the two points, in meters.
+ *
+ * @history 2012-09-26 - Jeremy Bante <http://scr.im/fugue> - Created
+ * @history 2014-09-15 - Jeremy Bante <http://scr.im/fugue> - Using the
+ * Andoyer-Lambert-Thomas formula approximating distances on an ellipsoid
+ * instead of a sphere to improve accuracy.
+ * @history 2017-03-16 - Erik Shagdar - updated function to use the
+ * specified parameters.
+ * @history 2019-01-04 - Erik Shagdar - Return 0 if the same point is used.
+ *
+ * @see http://www.dtic.mil/dtic/tr/fulltext/u2/703541.pdf
+ * @see http://en.wikipedia.org/wiki/World_Geodetic_System
+ ******************************************************************************/
+
+If ( latitude1 = latitude2 and longitude1 = longitude2 ;
+	0 ;
+	Let ( [
+		// WGS84
+		_equatorialRadius = 6378137 ;   // meters
+		_flattening = .0033528106647475 ;   // 1 / 298.257223563
+
+		_θ1 =
+			Atan (
+				( 1 - _flattening ) * Tan ( Radians ( latitude1 ) )
+			) ;
+		_θ2 =
+			Atan ( ( 1 - _flattening ) * Tan ( Radians ( latitude2 ) ) ) ;
+		_θm = ( _θ1 + _θ2 ) / 2 ;
+		_sinθm = Sin ( _θm ) ;
+		_cosθm = Cos ( _θm ) ;
+		_Δθm = ( _θ2 - _θ1 ) / 2 ;
+		_sinΔθm = Sin ( _Δθm ) ;
+		_cosΔθm = Cos ( _Δθm ) ;
+		_Δλ = Radians ( longitude2 - longitude1 ) ;
+		_h = _cosΔθm ^ 2 - _sinθm ^ 2 ;
+		_l = _sinΔθm ^ 2 + _h * Sin ( _Δλ / 2 ) ^ 2 ;
+		_cos_d = 1 - 2 * _l ;
+		_d = Acos ( _cos_d ) ;
+		_sin_d = Sin ( _d ) ;
+		_uu = 2 * _sinθm ^ 2 * _cosΔθm ^ 2 / ( 1 - _l ) ;
+		_vv = 2 * _sinΔθm ^ 2 * _cosθm ^ 2 / _l ;
+		_x = _uu + _vv ;
+		_t = _d / _sin_d ;
+		_e = 2 * _cos_d ;
+		_y = _uu - _vv ;
+		_dd = 4 * _t ^ 2 ;
+		_b = 2 * _dd ;
+		_a = _dd * _e ;
+		_c = _t - .5 * ( _a - _e ) ;
+		_n1 = _x * ( _a + _c * _x ) ;
+		_n2 = _y * ( _b + _e * _y ) ;
+		_n3 = _dd * _x * _y ;
+		_δ1d = .25 * _flattening * ( _t * _x - _y ) ;
+		_δ2d = ( _n1 - _n2 + _n3 ) * _flattening ^ 2 / 64 ;
+		_distance = _equatorialRadius * _sin_d * ( _t - _δ1d + _δ2d )
+	] ;
+	_distance
+	)
+)