Skip to content

Commit 5ce87c6

Browse files
derickrderickr
committed
Added the fraction_along_gc_line() function.
This is the code behind http://drck.me/gc-a98.
1 parent 801eaec commit 5ce87c6

File tree

5 files changed

+243
-0
lines changed

5 files changed

+243
-0
lines changed

‎geospatial.c

Lines changed: 105 additions & 0 deletions
Original file line numberDiff line numberDiff line change
@@ -35,12 +35,20 @@ ZEND_BEGIN_ARG_INFO_EX(haversine_args, 0, 0, 4)
3535
ZEND_ARG_INFO(0, radius)
3636
ZEND_END_ARG_INFO()
3737

38+
ZEND_BEGIN_ARG_INFO_EX(fraction_along_gc_line_args, 0, 0, 3)
39+
ZEND_ARG_INFO(0, geoJsonPoint1)
40+
ZEND_ARG_INFO(0, geoJsonPoint2)
41+
ZEND_ARG_INFO(0, fraction)
42+
ZEND_ARG_INFO(0, radius)
43+
ZEND_END_ARG_INFO()
44+
3845
/* {{{ geospatial_functions[]
3946
*
4047
* Every user visible function must have an entry in geospatial_functions[].
4148
*/
4249
const zend_function_entry geospatial_functions[] = {
4350
PHP_FE(haversine, haversine_args)
51+
PHP_FE(fraction_along_gc_line, fraction_along_gc_line_args)
4452
{ NULL, NULL, NULL }
4553
};
4654
/* }}} */
@@ -97,6 +105,56 @@ PHP_MINFO_FUNCTION(geospatial)
97105
}
98106
/* }}} */
99107

108+
/* {{{ Helpers */
109+
void retval_point_from_coordinates(zval *return_value, double lon, double lat)
110+
{
111+
zval *coordinates;
112+
113+
array_init(return_value);
114+
MAKE_STD_ZVAL(coordinates);
115+
array_init(coordinates);
116+
add_assoc_string_ex(return_value, "type", sizeof("type"), "Point", 1);
117+
add_next_index_double(coordinates, lon);
118+
add_next_index_double(coordinates, lat);
119+
add_assoc_zval_ex(return_value, "coordinates", sizeof("coordinates"), coordinates);
120+
}
121+
122+
int geojson_point_to_lon_lat(zval *point, double *lon, double *lat)
123+
{
124+
zval **type, **coordinates, **z_lon, **z_lat;
125+
HashTable *coords;
126+
127+
if (zend_hash_find(HASH_OF(point), "type", sizeof("type"), (void**) &type) != SUCCESS) {
128+
return 0;
129+
}
130+
if (Z_TYPE_PP(type) != IS_STRING || strcmp(Z_STRVAL_PP(type), "Point") != 0) {
131+
return 0;
132+
}
133+
if (zend_hash_find(HASH_OF(point), "coordinates", sizeof("coordinates"), (void**) &coordinates) != SUCCESS) {
134+
return 0;
135+
}
136+
if (Z_TYPE_PP(coordinates) != IS_ARRAY) {
137+
return 0;
138+
}
139+
coords = HASH_OF(*coordinates);
140+
if (coords->nNumOfElements != 2) {
141+
return 0;
142+
}
143+
if (zend_hash_index_find(coords, 0, (void**) &z_lon) != SUCCESS) {
144+
return 0;
145+
}
146+
if (zend_hash_index_find(coords, 1, (void**) &z_lat) != SUCCESS) {
147+
return 0;
148+
}
149+
convert_to_double_ex(z_lon);
150+
convert_to_double_ex(z_lat);
151+
*lon = Z_DVAL_PP(z_lon);
152+
*lat = Z_DVAL_PP(z_lat);
153+
return 1;
154+
}
155+
156+
/* }}} */
157+
100158
double php_geo_haversine(double from_lat, double from_long, double to_lat, double to_long)
101159
{
102160
double delta_lat, delta_long;
@@ -131,6 +189,53 @@ PHP_FUNCTION(haversine)
131189
}
132190
/* }}} */
133191

192+
void php_geo_fraction_along_gc_line(double from_lat, double from_long, double to_lat, double to_long, double fraction, double radius, double *res_lat, double *res_long)
193+
{
194+
double distance;
195+
double a, b, x, y, z;
196+
197+
/* First we calculate the distance */
198+
distance = php_geo_haversine(from_lat, from_long, to_lat, to_long);
199+
200+
a = sin((1 - fraction) * distance) / sin(distance);
201+
b = sin(fraction * distance) / sin(distance);
202+
x = a * cos(from_lat) * cos(from_long) + b * cos(to_lat) * cos(to_long);
203+
y = a * cos(from_lat) * sin(from_long) + b * cos(to_lat) * sin(to_long);
204+
z = a * sin(from_lat) + b * sin(to_lat);
205+
206+
*res_lat = atan2(z, sqrt(x * x + y * y));
207+
*res_long = atan2(y, x);
208+
}
209+
210+
/* {{{ proto GeoJSON fraction_along_gc_line(GeoJSONPoint from, GeoJSONPoint to, double fraction [, double radius ])
211+
* Calculates a lat/long pair at a fraction (0-1) of the distance along a GC line */
212+
PHP_FUNCTION(fraction_along_gc_line)
213+
{
214+
zval *from_geojson, *to_geojson;
215+
double from_lat, from_long, to_lat, to_long, fraction;
216+
double radius = GEO_EARTH_RADIUS;
217+
double res_lat, res_long;
218+
219+
if (zend_parse_parameters(ZEND_NUM_ARGS() TSRMLS_CC, "aad|d", &from_geojson, &to_geojson, &fraction, &radius) == FAILURE) {
220+
return;
221+
}
222+
223+
geojson_point_to_lon_lat(from_geojson, &from_long, &from_lat);
224+
geojson_point_to_lon_lat(to_geojson, &to_long, &to_lat);
225+
226+
php_geo_fraction_along_gc_line(
227+
from_lat * GEO_DEG_TO_RAD,
228+
from_long * GEO_DEG_TO_RAD,
229+
to_lat * GEO_DEG_TO_RAD,
230+
to_long * GEO_DEG_TO_RAD,
231+
fraction, radius,
232+
&res_lat, &res_long
233+
);
234+
235+
retval_point_from_coordinates(return_value, res_long / GEO_DEG_TO_RAD, res_lat / GEO_DEG_TO_RAD);
236+
}
237+
/* }}} */
238+
134239
/*
135240
* Local variables:
136241
* tab-width: 4

‎php_geospatial.h

Lines changed: 1 addition & 0 deletions
Original file line numberDiff line numberDiff line change
@@ -47,6 +47,7 @@ PHP_MINIT_FUNCTION(geospatial);
4747
PHP_MINFO_FUNCTION(geospatial);
4848

4949
PHP_FUNCTION(haversine);
50+
PHP_FUNCTION(fraction_along_gc_line);
5051

5152
#endif /* PHP_GEOSPATIAL_H */
5253

‎tests/fraction_along-001.phpt

Lines changed: 81 additions & 0 deletions
Original file line numberDiff line numberDiff line change
@@ -0,0 +1,81 @@
1+
--TEST--
2+
Test for "fraction_along_gc_line" #1
3+
--FILE--
4+
<?php
5+
$point1 = [ 'type' => 'Point', 'coordinates' => [ 5, 10 ] ];
6+
$point2 = [ 'type' => 'Point', 'coordinates' => [ 15, 10 ] ];
7+
8+
var_dump(fraction_along_gc_line($point1, $point2, 0));
9+
var_dump(fraction_along_gc_line($point1, $point2, 0.2));
10+
var_dump(fraction_along_gc_line($point1, $point2, 0.4));
11+
var_dump(fraction_along_gc_line($point1, $point2, 0.6));
12+
var_dump(fraction_along_gc_line($point1, $point2, 0.8));
13+
var_dump(fraction_along_gc_line($point1, $point2, 1));
14+
?>
15+
--EXPECT--
16+
array(2) {
17+
["type"]=>
18+
string(5) "Point"
19+
["coordinates"]=>
20+
array(2) {
21+
[0]=>
22+
float(5)
23+
[1]=>
24+
float(10)
25+
}
26+
}
27+
array(2) {
28+
["type"]=>
29+
string(5) "Point"
30+
["coordinates"]=>
31+
array(2) {
32+
[0]=>
33+
float(6.9998522347268)
34+
[1]=>
35+
float(10.023944943799)
36+
}
37+
}
38+
array(2) {
39+
["type"]=>
40+
string(5) "Point"
41+
["coordinates"]=>
42+
array(2) {
43+
[0]=>
44+
float(8.9999260791276)
45+
[1]=>
46+
float(10.035925156339)
47+
}
48+
}
49+
array(2) {
50+
["type"]=>
51+
string(5) "Point"
52+
["coordinates"]=>
53+
array(2) {
54+
[0]=>
55+
float(11.000073920872)
56+
[1]=>
57+
float(10.035925156339)
58+
}
59+
}
60+
array(2) {
61+
["type"]=>
62+
string(5) "Point"
63+
["coordinates"]=>
64+
array(2) {
65+
[0]=>
66+
float(13.000147765273)
67+
[1]=>
68+
float(10.023944943799)
69+
}
70+
}
71+
array(2) {
72+
["type"]=>
73+
string(5) "Point"
74+
["coordinates"]=>
75+
array(2) {
76+
[0]=>
77+
float(15)
78+
[1]=>
79+
float(10)
80+
}
81+
}

‎tests/fraction_along-002.phpt

Lines changed: 21 additions & 0 deletions
Original file line numberDiff line numberDiff line change
@@ -0,0 +1,21 @@
1+
--TEST--
2+
Test for "fraction_along_gc_line" #2
3+
--FILE--
4+
<?php
5+
$point1 = [ 'type' => 'Point', 'coordinates' => [ 0.1062, 51.5171 ] ];
6+
$point2 = [ 'type' => 'Point', 'coordinates' => [ 3.2200, 55.9500 ] ];
7+
8+
var_dump(fraction_along_gc_line($point1, $point2, 0.5));
9+
?>
10+
--EXPECT--
11+
array(2) {
12+
["type"]=>
13+
string(5) "Point"
14+
["coordinates"]=>
15+
array(2) {
16+
[0]=>
17+
float(1.5809481271999)
18+
[1]=>
19+
float(53.743611334154)
20+
}
21+
}

‎tests/fraction_along-003.phpt

Lines changed: 35 additions & 0 deletions
Original file line numberDiff line numberDiff line change
@@ -0,0 +1,35 @@
1+
--TEST--
2+
Test for "fraction_along_gc_line" #3
3+
--FILE--
4+
<?php
5+
$point1 = [ 'type' => 'Point', 'coordinates' => [ 0, 90 ] ];
6+
$point2 = [ 'type' => 'Point', 'coordinates' => [ 180, 0 ] ];
7+
var_dump(fraction_along_gc_line($point1, $point2, 0.5));
8+
9+
$point1 = [ 'type' => 'Point', 'coordinates' => [ 0, 90 ] ];
10+
$point2 = [ 'type' => 'Point', 'coordinates' => [ 3, -90 ] ];
11+
var_dump(fraction_along_gc_line($point1, $point2, 0.5));
12+
?>
13+
--EXPECT--
14+
array(2) {
15+
["type"]=>
16+
string(5) "Point"
17+
["coordinates"]=>
18+
array(2) {
19+
[0]=>
20+
float(180)
21+
[1]=>
22+
float(45)
23+
}
24+
}
25+
array(2) {
26+
["type"]=>
27+
string(5) "Point"
28+
["coordinates"]=>
29+
array(2) {
30+
[0]=>
31+
float(1.5)
32+
[1]=>
33+
float(0)
34+
}
35+
}

0 commit comments

Comments
 (0)