1 #include "line.h"
2 
3 /* Line functions */
4 
5 PG_FUNCTION_INFO_V1(sphereline_in);
6 PG_FUNCTION_INFO_V1(sphereline_from_point);
7 PG_FUNCTION_INFO_V1(sphereline_from_points);
8 PG_FUNCTION_INFO_V1(sphereline_from_trans);
9 PG_FUNCTION_INFO_V1(sphereline_meridian);
10 PG_FUNCTION_INFO_V1(sphereline_swap_beg_end);
11 PG_FUNCTION_INFO_V1(sphereline_turn);
12 PG_FUNCTION_INFO_V1(sphereline_begin);
13 PG_FUNCTION_INFO_V1(sphereline_end);
14 PG_FUNCTION_INFO_V1(sphereline_length);
15 PG_FUNCTION_INFO_V1(sphereline_cont_point);
16 PG_FUNCTION_INFO_V1(sphereline_cont_point_neg);
17 PG_FUNCTION_INFO_V1(sphereline_cont_point_com);
18 PG_FUNCTION_INFO_V1(sphereline_cont_point_com_neg);
19 PG_FUNCTION_INFO_V1(spherecircle_cont_line);
20 PG_FUNCTION_INFO_V1(spherecircle_cont_line_neg);
21 PG_FUNCTION_INFO_V1(spherecircle_cont_line_com);
22 PG_FUNCTION_INFO_V1(spherecircle_cont_line_com_neg);
23 PG_FUNCTION_INFO_V1(sphereline_overlap_circle);
24 PG_FUNCTION_INFO_V1(sphereline_overlap_circle_neg);
25 PG_FUNCTION_INFO_V1(sphereline_overlap_circle_com);
26 PG_FUNCTION_INFO_V1(sphereline_overlap_circle_com_neg);
27 PG_FUNCTION_INFO_V1(sphereline_equal);
28 PG_FUNCTION_INFO_V1(sphereline_equal_neg);
29 PG_FUNCTION_INFO_V1(sphereline_crosses);
30 PG_FUNCTION_INFO_V1(sphereline_crosses_neg);
31 PG_FUNCTION_INFO_V1(sphereline_overlap);
32 PG_FUNCTION_INFO_V1(sphereline_overlap_neg);
33 PG_FUNCTION_INFO_V1(spheretrans_from_line);
34 PG_FUNCTION_INFO_V1(spheretrans_line);
35 PG_FUNCTION_INFO_V1(spheretrans_line_inverse);
36 
37 /*
38  * Swaps the beginning and ending of the line.
39  */
40 static void
sline_swap_beg_end(SLine * out,const SLine * in)41 sline_swap_beg_end(SLine *out, const SLine *in)
42 {
43 	SLine	l;
44 	SEuler	se;
45 
46 	l.length = in->length;
47 	l.phi = -in->length;
48 	l.theta = PI;
49 	l.psi = 0.0;
50 	seuler_set_zxz(&se);
51 	se.phi = in->phi;
52 	se.theta = in->theta;
53 	se.psi = in->psi;
54 	euler_sline_trans(out, &l, &se);
55 }
56 
57 bool
sline_eq(const SLine * l1,const SLine * l2)58 sline_eq(const SLine *l1, const SLine *l2)
59 {
60 	if (FPne(l1->length, l2->length))
61 	{
62 		return false;
63 	}
64 	else
65 	{
66 		SEuler	e1, e2;
67 
68 		seuler_set_zxz(&e1);
69 		seuler_set_zxz(&e2);
70 		e1.phi = l1->phi;
71 		e1.theta = l1->theta;
72 		e1.psi = l1->psi;
73 		e2.phi = (FPeq(l2->length, PID)) ? (l1->phi) : (l2->phi);
74 		e2.theta = l2->theta;
75 		e2.psi = l2->psi;
76 		return (strans_eq(&e1, &e2));
77 	}
78 	return false;
79 }
80 
81 
82 bool
sline_from_points(SLine * sl,const SPoint * pbeg,const SPoint * pend)83 sline_from_points(SLine *sl, const SPoint *pbeg, const SPoint *pend)
84 {
85 	SEuler	se;
86 	float8	l;
87 
88 	l = spoint_dist(pbeg, pend);
89 
90 	if (FPeq(l, PI))
91 	{
92 		if (FPeq(pbeg->lng, pend->lng))
93 		{
94 			sline_meridian(sl, pbeg->lng);
95 			return true;
96 		}
97 		return false;
98 	}
99 
100 	if (spherevector_to_euler(&se, pbeg, pend))
101 	{
102 		sl->phi = se.phi;
103 		sl->theta = se.theta;
104 		sl->psi = se.psi;
105 		sl->length = l;
106 	}
107 	else
108 	{
109 		sl->phi = PIH;
110 		sl->theta = pbeg->lat;
111 		sl->psi = pbeg->lng - PIH;
112 		sl->length = 0.0;
113 	}
114 	return (true);
115 
116 }
117 
118 void
sline_meridian(SLine * sl,float8 lng)119 sline_meridian(SLine *sl, float8 lng)
120 {
121 	sl->phi = -PIH;
122 	sl->theta = PIH;
123 	sl->psi = (lng < 0.0) ? lng + PID : lng;
124 	sl->length = PI;
125 }
126 
127 void
sline_begin(SPoint * p,const SLine * l)128 sline_begin(SPoint *p, const SLine *l)
129 {
130 	const SPoint	tmp = {0.0, 0.0};
131 	SEuler			se;
132 
133 	sphereline_to_euler(&se, l);
134 	euler_spoint_trans(p, &tmp, &se);
135 }
136 
137 void
sline_end(SPoint * p,const SLine * l)138 sline_end(SPoint *p, const SLine *l)
139 {
140 	SPoint	tmp = {0.0, 0.0};
141 	SEuler	se;
142 
143 	tmp.lng = l->length;
144 	sphereline_to_euler(&se, l);
145 	euler_spoint_trans(p, &tmp, &se);
146 }
147 
148  /*
149   * Place begin of a line "l" as vector "v".
150   */
151 static void
sline_vector_begin(Vector3D * v,const SLine * l)152 sline_vector_begin(Vector3D *v, const SLine *l)
153 {
154 	const Vector3D	tmp = {1.0, 0.0, 0.0};
155 	SEuler			se;
156 
157 	sphereline_to_euler(&se, l);
158 	euler_vector_trans(v, &tmp, &se);
159 }
160 
161  /*
162   * Place end of a line "l" as vector "v".
163   */
164 static void
sline_vector_end(Vector3D * v,const SLine * l)165 sline_vector_end(Vector3D *v, const SLine *l)
166 {
167 	Vector3D	tmp = {0.0, 0.0, 0.0};
168 	SEuler		se;
169 
170 	tmp.x = cos(l->length);
171 	tmp.y = sin(l->length);
172 	sphereline_to_euler(&se, l);
173 	euler_vector_trans(v, &tmp, &se);
174 }
175 
176 void
sline_min_max_lat(const SLine * sl,float8 * minlat,float8 * maxlat)177 sline_min_max_lat(const SLine *sl, float8 *minlat, float8 *maxlat)
178 {
179 	float8		inc = sl->theta - floor(sl->theta / PID) * PID;
180 
181 	inc = (inc > PI) ? (PID - inc) : (inc);
182 
183 	if (FPzero(inc) || FPeq(inc, PI))
184 	{
185 		*minlat = *maxlat = 0.0;
186 		return;
187 	}
188 	else
189 	{
190 		SEuler		se;
191 		SLine		nl;
192 		SPoint		tp;
193 
194 		float8		lng;
195 
196 		seuler_set_zxz(&se);
197 		se.phi = -sl->psi;
198 		se.theta = (inc > PIH) ? (PI - 2 * inc) : (0.0);
199 		se.psi = 0.0;
200 		euler_sline_trans(&nl, sl, &se);
201 
202 		/* Now ascending node at (0,0), line ascending */
203 
204 		sline_begin(&tp, &nl);
205 		*minlat = *maxlat = tp.lat;
206 
207 		sline_end(&tp, &nl);
208 		*minlat = Min(tp.lat, *minlat);
209 		*maxlat = Max(tp.lat, *maxlat);
210 
211 		for (lng = PIH; lng < PID; lng += PI)
212 		{
213 			tp.lng = lng;
214 			tp.lat = (lng < PI) ? (inc) : (-inc);
215 			if (spoint_at_sline(&tp, &nl))
216 			{
217 				*minlat = Min(tp.lat, *minlat);
218 				*maxlat = Max(tp.lat, *maxlat);
219 			}
220 		}
221 	}
222 }
223 
224 int32
sphereline_latitude_points(const SLine * sl,float8 lat,SPoint * p1,SPoint * p2)225 sphereline_latitude_points(const SLine *sl, float8 lat, SPoint *p1, SPoint *p2)
226 {
227 	float8		inc = sl->theta - floor(sl->theta / PID) * PID;
228 	int32		ret = 0;
229 
230 	if (FPgt(lat, PIH))
231 		return 0;
232 	if (FPlt(lat, -PIH))
233 		return 0;
234 
235 	inc = (inc > PI) ? (PID - inc) : (inc);
236 
237 	if (FPzero(inc) || FPeq(inc, PI))
238 	{
239 
240 		if (FPzero(lat))
241 			return -1;
242 		else
243 			return 0;
244 	}
245 	else
246 	{
247 		SLine		nl;
248 		float8		rot = (inc > PIH) ? (sl->psi - PIH) : (sl->psi + PIH);
249 		bool		p1b,
250 					p2b;
251 
252 		/* Transform maximum latitude of full line to longitude 0.0 */
253 		memcpy((void *) &nl, (void *) sl, sizeof(SLine));
254 		nl.psi = (inc > PIH) ? (PIH) : (-PIH);
255 
256 		p1->lat = p2->lat = lat;
257 		p1->lng = p2->lng = 0.0;
258 
259 		if (FPeq(inc, PIH))
260 		{
261 			p1->lng = PIH;
262 			p2->lng = -PIH;
263 			ret = 2;
264 		}
265 		else
266 		{
267 			float8		ainc;
268 
269 			ainc = fabs(inc - ((inc > PIH) ? (PI) : (0.0)));
270 			if (FPgt(lat, ainc))
271 				return 0;
272 			else if (FPlt(lat, -ainc))
273 				return 0;
274 			else if (FPeq(lat, ainc))
275 			{
276 				p1->lng = 0.0;
277 				ret = 1;
278 			}
279 			else if (FPeq(lat, -ainc))
280 			{
281 				p1->lng = PI;
282 				ret = 1;
283 			}
284 			else
285 			{
286 				p1->lng = acos(sin(lat) * cos(ainc) / (sin(ainc) * cos(lat)));
287 				p2->lng = PID - p1->lng;
288 				ret = 2;
289 			}
290 		}
291 
292 		if (ret == 1)
293 		{
294 			p1b = spoint_at_sline(p1, &nl);
295 			if (!p1b)
296 			{
297 				ret = 0;
298 			}
299 		}
300 		else if (ret == 2)
301 		{
302 			p1b = spoint_at_sline(p1, &nl);
303 			p2b = spoint_at_sline(p2, &nl);
304 
305 			if (p1b && p2b)
306 			{
307 				ret = 2;
308 			}
309 			else if (!p1b && p2b)
310 			{
311 				ret = 1;
312 				memcpy((void *) p1, (void *) p2, sizeof(SPoint));
313 			}
314 			else if (p1b && !p2b)
315 			{
316 				ret = 1;
317 			}
318 			else
319 			{
320 				ret = 0;
321 			}
322 		}
323 
324 		if (ret > 0)
325 		{
326 			p1->lng += rot;
327 			p2->lng += rot;
328 			spoint_check(p1);
329 			spoint_check(p2);
330 		}
331 	}
332 	return ret;
333 }
334 
335 int8
sphereline_circle_pos(const SLine * sl,const SCIRCLE * sc)336 sphereline_circle_pos(const SLine *sl, const SCIRCLE *sc)
337 {
338 	float8			i, mi;
339 	const float8	step = (PI - 0.01);
340 	SPoint			p[2] = {{0.0, 0.0}, {0.0, 0.0}};
341 	SCIRCLE			c;
342 	bool			bbeg, bend;
343 	SEuler			se;
344 	int				contain;
345 
346 	if (FPzero(sl->length))
347 	{
348 		/* line is point */
349 		sline_begin(&p[0], sl);
350 		if (spoint_in_circle(&p[0], &c))
351 		{
352 			return PGS_CIRCLE_CONT_LINE;
353 		}
354 		else
355 		{
356 			return PGS_CIRCLE_LINE_AVOID;
357 		}
358 	}
359 
360 	contain = 0;
361 	sphereline_to_euler_inv(&se, sl);
362 	euler_scircle_trans(&c, sc, &se);
363 
364 	mi = sl->length / step;
365 
366 	/* split line in segments and check for each of it */
367 	for (i = 0.0; i < mi; i += 1.0)
368 	{
369 		p[0].lng = i * step;
370 		p[1].lng = (((i + 1.0) > mi) ? (sl->length) : ((i + 1.0) * step));
371 
372 		bbeg = spoint_in_circle(&p[0], &c);
373 		bend = spoint_in_circle(&p[1], &c);
374 
375 		if (bbeg && bend)
376 		{
377 			contain++;
378 		}
379 		else if (bbeg || bend)
380 		{
381 			return PGS_CIRCLE_LINE_OVER;
382 		}
383 		else if (FPle(((c.center.lat < 0) ? (-c.center.lat) : (c.center.lat)),
384 					  c.radius) &&
385 				 FPge(c.center.lng, p[0].lng) &&
386 				 FPle(c.center.lng, p[1].lng))
387 		{
388 			return PGS_CIRCLE_LINE_OVER;
389 		}
390 		else if (contain > 0)
391 		{
392 			return PGS_CIRCLE_LINE_OVER;
393 		}
394 	}
395 
396 	if (contain > 0 && contain == (floor(mi) + 1))
397 	{
398 		return PGS_CIRCLE_CONT_LINE;
399 	}
400 
401 	return PGS_CIRCLE_LINE_AVOID;
402 }
403 
404 bool
sline_circle_touch(const SLine * sl,const SCIRCLE * sc)405 sline_circle_touch(const SLine *sl, const SCIRCLE *sc)
406 {
407 	/* we assume here, line and circle overlap */
408 	SEuler	se;
409 	SCIRCLE	tc;
410 
411 	sphereline_to_euler_inv(&se, sl);
412 	euler_scircle_trans(&tc, sc, &se);
413 	if (FPge(tc.center.lng, 0.0) && FPle(tc.center.lng, sl->length))
414 	{
415 		if (FPeq(fabs(tc.center.lat), sc->radius))
416 		{
417 			return true;
418 		}
419 		return false;
420 	}
421 	else
422 	{
423 		SPoint		p;
424 
425 		p.lng = p.lat = 0.0;
426 		if (FPeq(spoint_dist(&p, &tc.center), sc->radius))
427 		{
428 			return true;
429 		}
430 		p.lng = sl->length;
431 		if (FPeq(spoint_dist(&p, &tc.center), sc->radius))
432 		{
433 			return true;
434 		}
435 		return false;
436 	}
437 }
438 
439 int8
sline_sline_pos(const SLine * l1,const SLine * l2)440 sline_sline_pos(const SLine *l1, const SLine *l2)
441 {
442 	const SLine	   *il1, *il2;
443 	Vector3D		v[2][2],
444 					vtmp;
445 	SEuler			se;
446 	SLine			sl1, sl2, lseg;
447 	SPoint			p[4];
448 	const float8	seg_length = (PI - 0.1);
449 	float8			seg_begin;
450 
451 	if (sline_eq(l1, l2))
452 	{
453 		return PGS_LINE_EQUAL;
454 	}
455 	sline_swap_beg_end(&sl1, l1);
456 	if (sline_eq(&sl1, l2))
457 	{
458 		return PGS_LINE_CONT_LINE;
459 	}
460 
461 	/* transform the larger line into equator ( begin at (0,0) ) */
462 	if (FPge(l1->length, l2->length))
463 	{
464 		il1 = l1;
465 		il2 = l2;
466 	}
467 	else
468 	{
469 		il1 = l2;
470 		il2 = l1;
471 	}
472 
473 	if (FPzero(il1->length))
474 	{
475 		/* both are points */
476 		return PGS_LINE_AVOID;
477 	}
478 
479 	sl1.phi    = sl1.theta = sl1.psi = 0.0;
480 	p[0].lat   = p[0].lng  = p[1].lat = 0.0;
481 	sl1.length = p[1].lng  = il1->length;
482 	v[0][0].x  = 1.0;
483 	v[0][0].y  = 0.0;
484 	v[0][0].z  = 0.0;
485 	v[0][1].x  = cos(il1->length);
486 	v[0][1].y  = sin(il1->length);
487 	v[0][1].z  = 0.0;
488 
489 	sphereline_to_euler_inv(&se, il1);
490 	sline_vector_begin(&vtmp, il2);
491 	euler_vector_trans(&v[1][0], &vtmp, &se);
492 	sline_vector_end(&vtmp, il2);
493 	euler_vector_trans(&v[1][1], &vtmp, &se);
494 	vector3d_spoint(&p[2], &v[1][0]);
495 	vector3d_spoint(&p[3], &v[1][1]);
496 
497 	/* check connected lines */
498 	if (FPgt(il2->length, 0.0) && (vector3d_eq(&v[0][0], &v[1][0]) ||
499 								   vector3d_eq(&v[0][0], &v[1][1]) ||
500 								   vector3d_eq(&v[0][1], &v[1][0]) ||
501 								   vector3d_eq(&v[0][1], &v[1][1])))
502 	{
503 		return PGS_LINE_CONNECT;
504 	}
505 
506 	/* Check, sl2 is at equator */
507 	if (FPzero(p[2].lat) && FPzero(p[3].lat))
508 	{
509 		bool	a1 = spoint_at_sline(&p[2], &sl1);
510 		bool	a2 = spoint_at_sline(&p[3], &sl1);
511 
512 		if (a1 && a2)
513 		{
514 			if (il1 == l2)
515 			{
516 				return PGS_LINE_OVER;
517 			}
518 			else
519 			{
520 				return PGS_LINE_CONT_LINE;
521 			}
522 		}
523 		else if (a1 || a2)
524 		{
525 			return PGS_LINE_OVER;
526 		}
527 		return PGS_LINE_AVOID;
528 	}
529 
530 	/* Now sl2 is not at equator */
531 
532 	if (FPle(il2->length, seg_length))
533 	{
534 		bool		a1 = (FPge(p[2].lat, 0.0) && FPle(p[3].lat, 0.0));
535 		/* sl2 crosses equator desc. */
536 		bool		a2 = (FPle(p[2].lat, 0.0) && FPge(p[3].lat, 0.0));
537 		/* sl2 crosses equator asc. */
538 
539 		if (a1 || a2)
540 		{
541 			SPoint		sp;
542 
543 			euler_sline_trans(&sl2, il2, &se);
544 			sphereline_to_euler_inv(&se, &sl2);
545 			sp.lng = ((a1) ? (PI) : (0.0)) - se.phi;
546 			/* node */
547 			sp.lat = 0;
548 			spoint_check(&sp);
549 			if (FPge(sp.lng, 0.0) && FPle(sp.lng, p[1].lng))
550 			{
551 				return PGS_LINE_CROSS;
552 			}
553 		}
554 		return PGS_LINE_AVOID;
555 	}
556 
557 	/*
558 	 * We split the smaller line in segments with length less than 180 deg
559 	 */
560 	euler_sline_trans(&sl2, il2, &se);
561 	for (seg_begin = 0.0; seg_begin < il2->length; seg_begin += seg_length)
562 	{
563 
564 		lseg.length = ((seg_begin + seg_length) > il2->length) ?
565 					  (il2->length - seg_begin) : seg_length;
566 		lseg.phi    = sl2.phi + seg_begin;
567 		lseg.theta  = sl2.theta;
568 		lseg.psi    = sl2.psi;
569 
570 		if (sline_sline_pos(&sl1, &lseg) != PGS_LINE_AVOID)
571 		{
572 			return PGS_LINE_CROSS;
573 		}
574 	}
575 	return PGS_LINE_AVOID;
576 }
577 
578 void
sphereline_to_euler_inv(SEuler * se,const SLine * sl)579 sphereline_to_euler_inv(SEuler *se, const SLine *sl)
580 {
581 	sphereline_to_euler(se, sl);
582 	spheretrans_inv(se);
583 }
584 
585 void
sphereline_to_euler(SEuler * se,const SLine * sl)586 sphereline_to_euler(SEuler *se, const SLine *sl)
587 {
588 	seuler_set_zxz(se);
589 	se->phi = sl->phi;
590 	se->theta = sl->theta;
591 	se->psi = sl->psi;
592 }
593 
594 void
euler_sline_trans(SLine * out,const SLine * in,const SEuler * se)595 euler_sline_trans(SLine *out, const SLine *in, const SEuler *se)
596 {
597 	SEuler	stmp[2];
598 
599 	sphereline_to_euler(&stmp[0], in);
600 	seuler_trans_zxz(&stmp[1], &stmp[0], se);
601 	out->phi = stmp[1].phi;
602 	out->theta = stmp[1].theta;
603 	out->psi = stmp[1].psi;
604 	out->length = in->length;
605 }
606 
607 bool
spoint_at_sline(const SPoint * p,const SLine * sl)608 spoint_at_sline(const SPoint *p, const SLine *sl)
609 {
610 	SEuler	se;
611 	SPoint	sp;
612 
613 	sphereline_to_euler_inv(&se, sl);
614 	euler_spoint_trans(&sp, p, &se);
615 
616 	if (FPzero(sp.lat))
617 	{
618 		if (FPge(sp.lng, 0.0) && FPle(sp.lng, sl->length))
619 		{
620 			return true;
621 		}
622 		else
623 		{
624 			return false;
625 		}
626 	}
627 	else
628 	{
629 		return false;
630 	}
631 }
632 
633 void
sline_center(SPoint * c,const SLine * sl)634 sline_center(SPoint *c, const SLine *sl)
635 {
636 	SEuler	se;
637 	SPoint	p;
638 
639 	p.lng = sl->length / 2.0;
640 	p.lat = 0.0;
641 	sphereline_to_euler(&se, sl);
642 	euler_spoint_trans(c, &p, &se);
643 }
644 
645 Datum
sphereline_in(PG_FUNCTION_ARGS)646 sphereline_in(PG_FUNCTION_ARGS)
647 {
648 	SLine		   *sl = (SLine *) palloc(sizeof(SLine));
649 	char		   *c = PG_GETARG_CSTRING(0);
650 	unsigned char	etype[3];
651 	float8			eang[3],
652 					length;
653 	SEuler			se,
654 					stmp,
655 					so;
656 	int				i;
657 
658 	void		sphere_yyparse(void);
659 
660 	init_buffer(c);
661 	sphere_yyparse();
662 	if (get_line(&eang[0], &eang[1], &eang[2], etype, &length))
663 	{
664 
665 		for (i = 0; i < 3; i++)
666 		{
667 			switch (i)
668 			{
669 				case 0:
670 					se.phi_a = etype[i];
671 					break;
672 				case 1:
673 					se.theta_a = etype[i];
674 					break;
675 				case 2:
676 					se.psi_a = etype[i];
677 					break;
678 			}
679 		}
680 		se.phi = eang[0];
681 		se.theta = eang[1];
682 		se.psi = eang[2];
683 
684 		stmp.phi = stmp.theta = stmp.psi = 0.0;
685 		stmp.phi_a = stmp.theta_a = stmp.psi_a = EULER_AXIS_Z;
686 		seuler_trans_zxz(&so, &se, &stmp);
687 
688 		sl->phi = so.phi;
689 		sl->theta = so.theta;
690 		sl->psi = so.psi;
691 
692 		if (FPge(length, PID))
693 		{
694 			length = PID;
695 		}
696 		sl->length = length;
697 
698 	}
699 	else
700 	{
701 		reset_buffer();
702 		pfree(sl);
703 		sl = NULL;
704 		elog(ERROR, "sphereline_in: parse error");
705 	}
706 	reset_buffer();
707 
708 	PG_RETURN_POINTER(sl);
709 
710 }
711 
712 Datum
sphereline_equal(PG_FUNCTION_ARGS)713 sphereline_equal(PG_FUNCTION_ARGS)
714 {
715 	SLine	   *l1 = (SLine *) PG_GETARG_POINTER(0);
716 	SLine	   *l2 = (SLine *) PG_GETARG_POINTER(1);
717 
718 	PG_RETURN_BOOL(sline_eq(l1, l2));
719 }
720 
721 Datum
sphereline_equal_neg(PG_FUNCTION_ARGS)722 sphereline_equal_neg(PG_FUNCTION_ARGS)
723 {
724 	SLine	   *l1 = (SLine *) PG_GETARG_POINTER(0);
725 	SLine	   *l2 = (SLine *) PG_GETARG_POINTER(1);
726 
727 	PG_RETURN_BOOL(!sline_eq(l1, l2));
728 }
729 
730 Datum
sphereline_from_point(PG_FUNCTION_ARGS)731 sphereline_from_point(PG_FUNCTION_ARGS)
732 {
733 	SLine	   *sl = (SLine *) palloc(sizeof(SLine));
734 	SPoint	   *p = (SPoint *) PG_GETARG_POINTER(0);
735 
736 	sline_from_points(sl, p, p);
737 	PG_RETURN_POINTER(sl);
738 }
739 
740 Datum
sphereline_from_points(PG_FUNCTION_ARGS)741 sphereline_from_points(PG_FUNCTION_ARGS)
742 {
743 	SLine	   *sl = (SLine *) palloc(sizeof(SLine));
744 	SPoint	   *beg = (SPoint *) PG_GETARG_POINTER(0);
745 	SPoint	   *end = (SPoint *) PG_GETARG_POINTER(1);
746 
747 	if (!sline_from_points(sl, beg, end))
748 	{
749 		pfree(sl);
750 		sl = NULL;
751 		elog(ERROR, "sphereline_from_points: length of line must be not equal 180 degrees");
752 	}
753 
754 	PG_RETURN_POINTER(sl);
755 }
756 
757 Datum
sphereline_from_trans(PG_FUNCTION_ARGS)758 sphereline_from_trans(PG_FUNCTION_ARGS)
759 {
760 	SLine	   *sl = (SLine *) palloc(sizeof(SLine));
761 	SEuler	   *se = (SEuler *) PG_GETARG_POINTER(0);
762 	float8		l = PG_GETARG_FLOAT8(1);
763 
764 	if (FPlt(l, 0.0))
765 	{
766 		pfree(sl);
767 		elog(ERROR, "sphereline_from_trans: length of line must be >= 0");
768 		PG_RETURN_NULL();
769 	}
770 	else
771 	{
772 		SEuler	tmp;
773 
774 		if (FPgt(l, PID))
775 		{
776 			l = PID;
777 		}
778 		strans_zxz(&tmp, se);
779 		sl->phi = tmp.phi;
780 		sl->theta = tmp.theta;
781 		sl->psi = tmp.psi;
782 		sl->length = l;
783 	}
784 
785 	PG_RETURN_POINTER(sl);
786 }
787 
788 Datum
sphereline_meridian(PG_FUNCTION_ARGS)789 sphereline_meridian(PG_FUNCTION_ARGS)
790 {
791 	SLine	   *out = (SLine *) palloc(sizeof(SLine));
792 	float8		lng = PG_GETARG_FLOAT8(0);
793 
794 	sline_meridian(out, lng);
795 	PG_RETURN_POINTER(out);
796 }
797 
798 Datum
sphereline_swap_beg_end(PG_FUNCTION_ARGS)799 sphereline_swap_beg_end(PG_FUNCTION_ARGS)
800 {
801 	SLine	   *in = (SLine *) PG_GETARG_POINTER(0);
802 	SLine	   *out = (SLine *) palloc(sizeof(SLine));
803 
804 	sline_swap_beg_end(out, in);
805 	PG_RETURN_POINTER(out);
806 }
807 
808 Datum
sphereline_turn(PG_FUNCTION_ARGS)809 sphereline_turn(PG_FUNCTION_ARGS)
810 {
811 	SLine	   *in = (SLine *) PG_GETARG_POINTER(0);
812 
813 	if (FPzero(in->length))
814 	{
815 		PG_RETURN_NULL();
816 	}
817 	else
818 	{
819 		SLine	   *out = (SLine *) palloc(sizeof(SLine));
820 		SEuler		se;
821 		SLine		tmp;
822 
823 		tmp.phi = 0.0;
824 		tmp.theta = PI;
825 		tmp.psi = 0.0;
826 		tmp.length = PID - in->length;
827 
828 		sphereline_to_euler(&se, in);
829 		euler_sline_trans(out, &tmp, &se);
830 
831 		PG_RETURN_POINTER(out);
832 	}
833 	PG_RETURN_NULL();
834 }
835 
836 Datum
sphereline_begin(PG_FUNCTION_ARGS)837 sphereline_begin(PG_FUNCTION_ARGS)
838 {
839 	SLine	   *sl = (SLine *) PG_GETARG_POINTER(0);
840 	SPoint	   *sp = (SPoint *) palloc(sizeof(SPoint));
841 
842 	sline_begin(sp, sl);
843 	PG_RETURN_POINTER(sp);
844 }
845 
846 Datum
sphereline_end(PG_FUNCTION_ARGS)847 sphereline_end(PG_FUNCTION_ARGS)
848 {
849 	SLine	   *sl = (SLine *) PG_GETARG_POINTER(0);
850 	SPoint	   *sp = (SPoint *) palloc(sizeof(SPoint));
851 
852 	sline_end(sp, sl);
853 	PG_RETURN_POINTER(sp);
854 }
855 
856 Datum
sphereline_length(PG_FUNCTION_ARGS)857 sphereline_length(PG_FUNCTION_ARGS)
858 {
859 	SLine	   *sl = (SLine *) PG_GETARG_POINTER(0);
860 
861 	PG_RETURN_FLOAT8(sl->length);
862 }
863 
864 
865 Datum
spherecircle_cont_line(PG_FUNCTION_ARGS)866 spherecircle_cont_line(PG_FUNCTION_ARGS)
867 {
868 	SCIRCLE    *c = (SCIRCLE *) PG_GETARG_POINTER(0);
869 	SLine	   *l = (SLine *) PG_GETARG_POINTER(1);
870 	int8		b = sphereline_circle_pos(l, c);
871 
872 	PG_RETURN_BOOL(b == PGS_CIRCLE_CONT_LINE);
873 }
874 
875 Datum
spherecircle_cont_line_neg(PG_FUNCTION_ARGS)876 spherecircle_cont_line_neg(PG_FUNCTION_ARGS)
877 {
878 	SCIRCLE    *c = (SCIRCLE *) PG_GETARG_POINTER(0);
879 	SLine	   *l = (SLine *) PG_GETARG_POINTER(1);
880 	int8		b = sphereline_circle_pos(l, c);
881 
882 	PG_RETURN_BOOL(b != PGS_CIRCLE_CONT_LINE);
883 }
884 
885 Datum
spherecircle_cont_line_com(PG_FUNCTION_ARGS)886 spherecircle_cont_line_com(PG_FUNCTION_ARGS)
887 {
888 	SCIRCLE    *c = (SCIRCLE *) PG_GETARG_POINTER(1);
889 	SLine	   *l = (SLine *) PG_GETARG_POINTER(0);
890 	int8		b = sphereline_circle_pos(l, c);
891 
892 	PG_RETURN_BOOL(b == PGS_CIRCLE_CONT_LINE);
893 }
894 
895 Datum
spherecircle_cont_line_com_neg(PG_FUNCTION_ARGS)896 spherecircle_cont_line_com_neg(PG_FUNCTION_ARGS)
897 {
898 	SCIRCLE    *c = (SCIRCLE *) PG_GETARG_POINTER(1);
899 	SLine	   *l = (SLine *) PG_GETARG_POINTER(0);
900 	int8		b = sphereline_circle_pos(l, c);
901 
902 	PG_RETURN_BOOL(b != PGS_CIRCLE_CONT_LINE);
903 }
904 
905 Datum
sphereline_overlap_circle(PG_FUNCTION_ARGS)906 sphereline_overlap_circle(PG_FUNCTION_ARGS)
907 {
908 	SLine	   *l = (SLine *) PG_GETARG_POINTER(0);
909 	SCIRCLE    *c = (SCIRCLE *) PG_GETARG_POINTER(1);
910 	int8		b = sphereline_circle_pos(l, c);
911 
912 	PG_RETURN_BOOL(b > PGS_CIRCLE_LINE_AVOID);
913 }
914 
915 Datum
sphereline_overlap_circle_neg(PG_FUNCTION_ARGS)916 sphereline_overlap_circle_neg(PG_FUNCTION_ARGS)
917 {
918 	SLine	   *l = (SLine *) PG_GETARG_POINTER(0);
919 	SCIRCLE    *c = (SCIRCLE *) PG_GETARG_POINTER(1);
920 	int8		b = sphereline_circle_pos(l, c);
921 
922 	PG_RETURN_BOOL(b == PGS_CIRCLE_LINE_AVOID);
923 }
924 
925 Datum
sphereline_overlap_circle_com(PG_FUNCTION_ARGS)926 sphereline_overlap_circle_com(PG_FUNCTION_ARGS)
927 {
928 	SLine	   *l = (SLine *) PG_GETARG_POINTER(1);
929 	SCIRCLE    *c = (SCIRCLE *) PG_GETARG_POINTER(0);
930 	int8		b = sphereline_circle_pos(l, c);
931 
932 	PG_RETURN_BOOL(b > PGS_CIRCLE_LINE_AVOID);
933 }
934 
935 Datum
sphereline_overlap_circle_com_neg(PG_FUNCTION_ARGS)936 sphereline_overlap_circle_com_neg(PG_FUNCTION_ARGS)
937 {
938 	SLine	   *l = (SLine *) PG_GETARG_POINTER(1);
939 	SCIRCLE    *c = (SCIRCLE *) PG_GETARG_POINTER(0);
940 	int8		b = sphereline_circle_pos(l, c);
941 
942 	PG_RETURN_BOOL(b == PGS_CIRCLE_LINE_AVOID);
943 }
944 
945 Datum
sphereline_crosses(PG_FUNCTION_ARGS)946 sphereline_crosses(PG_FUNCTION_ARGS)
947 {
948 	SLine	   *l1 = (SLine *) PG_GETARG_POINTER(0);
949 	SLine	   *l2 = (SLine *) PG_GETARG_POINTER(1);
950 	int8		r = sline_sline_pos(l1, l2);
951 
952 	PG_RETURN_BOOL(r == PGS_LINE_CROSS);
953 }
954 
955 Datum
sphereline_crosses_neg(PG_FUNCTION_ARGS)956 sphereline_crosses_neg(PG_FUNCTION_ARGS)
957 {
958 	SLine	   *l1 = (SLine *) PG_GETARG_POINTER(0);
959 	SLine	   *l2 = (SLine *) PG_GETARG_POINTER(1);
960 	int8		r = sline_sline_pos(l1, l2);
961 
962 	PG_RETURN_BOOL(r != PGS_LINE_CROSS);
963 }
964 
965 Datum
sphereline_overlap(PG_FUNCTION_ARGS)966 sphereline_overlap(PG_FUNCTION_ARGS)
967 {
968 	SLine	   *l1 = (SLine *) PG_GETARG_POINTER(0);
969 	SLine	   *l2 = (SLine *) PG_GETARG_POINTER(1);
970 
971 	PG_RETURN_BOOL(sline_sline_pos(l1, l2) > PGS_LINE_AVOID);
972 }
973 
974 Datum
sphereline_overlap_neg(PG_FUNCTION_ARGS)975 sphereline_overlap_neg(PG_FUNCTION_ARGS)
976 {
977 	SLine	   *l1 = (SLine *) PG_GETARG_POINTER(0);
978 	SLine	   *l2 = (SLine *) PG_GETARG_POINTER(1);
979 
980 	PG_RETURN_BOOL(sline_sline_pos(l1, l2) == PGS_LINE_AVOID);
981 }
982 
983 Datum
sphereline_cont_point(PG_FUNCTION_ARGS)984 sphereline_cont_point(PG_FUNCTION_ARGS)
985 {
986 	SLine	   *l = (SLine *) PG_GETARG_POINTER(0);
987 	SPoint	   *p = (SPoint *) PG_GETARG_POINTER(1);
988 
989 	PG_RETURN_BOOL(spoint_at_sline(p, l));
990 }
991 
992 Datum
sphereline_cont_point_neg(PG_FUNCTION_ARGS)993 sphereline_cont_point_neg(PG_FUNCTION_ARGS)
994 {
995 	SLine	   *l = (SLine *) PG_GETARG_POINTER(0);
996 	SPoint	   *p = (SPoint *) PG_GETARG_POINTER(1);
997 
998 	PG_RETURN_BOOL(!spoint_at_sline(p, l));
999 }
1000 
1001 Datum
sphereline_cont_point_com(PG_FUNCTION_ARGS)1002 sphereline_cont_point_com(PG_FUNCTION_ARGS)
1003 {
1004 	SLine	   *l = (SLine *) PG_GETARG_POINTER(1);
1005 	SPoint	   *p = (SPoint *) PG_GETARG_POINTER(0);
1006 
1007 	PG_RETURN_BOOL(spoint_at_sline(p, l));
1008 }
1009 
1010 Datum
sphereline_cont_point_com_neg(PG_FUNCTION_ARGS)1011 sphereline_cont_point_com_neg(PG_FUNCTION_ARGS)
1012 {
1013 	SLine	   *l = (SLine *) PG_GETARG_POINTER(1);
1014 	SPoint	   *p = (SPoint *) PG_GETARG_POINTER(0);
1015 
1016 	PG_RETURN_BOOL(!spoint_at_sline(p, l));
1017 }
1018 
1019 Datum
spheretrans_line(PG_FUNCTION_ARGS)1020 spheretrans_line(PG_FUNCTION_ARGS)
1021 {
1022 	SLine	   *sl = (SLine *) PG_GETARG_POINTER(0);
1023 	SEuler	   *se = (SEuler *) PG_GETARG_POINTER(1);
1024 	SLine	   *ret = (SLine *) palloc(sizeof(SLine));
1025 
1026 	euler_sline_trans(ret, sl, se);
1027 	PG_RETURN_POINTER(ret);
1028 }
1029 
1030 Datum
spheretrans_line_inverse(PG_FUNCTION_ARGS)1031 spheretrans_line_inverse(PG_FUNCTION_ARGS)
1032 {
1033 	Datum		sl = PG_GETARG_DATUM(0);
1034 	SEuler	   *se = (SEuler *) PG_GETARG_POINTER(1);
1035 	SEuler		tmp;
1036 	Datum		ret;
1037 
1038 	spheretrans_inverse(&tmp, se);
1039 	ret = DirectFunctionCall2(spheretrans_line,
1040 							  sl, PointerGetDatum(&tmp));
1041 	PG_RETURN_DATUM(ret);
1042 
1043 }
1044 
1045 Datum
spheretrans_from_line(PG_FUNCTION_ARGS)1046 spheretrans_from_line(PG_FUNCTION_ARGS)
1047 {
1048 	SLine	   *l = (SLine *) PG_GETARG_POINTER(0);
1049 	SEuler	   *e = (SEuler *) palloc(sizeof(SEuler));
1050 
1051 	sphereline_to_euler(e, l);
1052 	PG_RETURN_POINTER(e);
1053 }
1054