PostGIS 3.0.6dev-r@@SVN_REVISION@@
Loading...
Searching...
No Matches

◆ pta_unstroke()

LWGEOM * pta_unstroke ( const POINTARRAY points,
int32_t  srid 
)

Definition at line 975 of file lwstroke.c.

976{
977 int i = 0, j, k;
978 POINT4D a1, a2, a3, b;
979 POINT4D first, center;
980 char *edges_in_arcs;
981 int found_arc = LW_FALSE;
982 int current_arc = 1;
983 int num_edges;
984 int edge_type; /* non-zero if edge is part of an arc */
985 int start, end;
986 LWCOLLECTION *outcol;
987 /* Minimum number of edges, per quadrant, required to define an arc */
988 const unsigned int min_quad_edges = 2;
989
990 /* Die on null input */
991 if ( ! points )
992 lwerror("pta_unstroke called with null pointarray");
993
994 /* Null on empty input? */
995 if ( points->npoints == 0 )
996 return NULL;
997
998 /* We can't desegmentize anything shorter than four points */
999 if ( points->npoints < 4 )
1000 {
1001 /* Return a linestring here*/
1002 lwerror("pta_unstroke needs implementation for npoints < 4");
1003 }
1004
1005 /* Allocate our result array of vertices that are part of arcs */
1006 num_edges = points->npoints - 1;
1007 edges_in_arcs = lwalloc(num_edges + 1);
1008 memset(edges_in_arcs, 0, num_edges + 1);
1009
1010 /* We make a candidate arc of the first two edges, */
1011 /* And then see if the next edge follows it */
1012 while( i < num_edges-2 )
1013 {
1014 unsigned int arc_edges;
1015 double num_quadrants;
1016 double angle;
1017
1018 found_arc = LW_FALSE;
1019 /* Make candidate arc */
1020 getPoint4d_p(points, i , &a1);
1021 getPoint4d_p(points, i+1, &a2);
1022 getPoint4d_p(points, i+2, &a3);
1023 memcpy(&first, &a1, sizeof(POINT4D));
1024
1025 for( j = i+3; j < num_edges+1; j++ )
1026 {
1027 LWDEBUGF(4, "i=%d, j=%d", i, j);
1028 getPoint4d_p(points, j, &b);
1029 /* Does this point fall on our candidate arc? */
1030 if ( pt_continues_arc(&a1, &a2, &a3, &b) )
1031 {
1032 /* Yes. Mark this edge and the two preceding it as arc components */
1033 LWDEBUGF(4, "pt_continues_arc #%d", current_arc);
1034 found_arc = LW_TRUE;
1035 for ( k = j-1; k > j-4; k-- )
1036 edges_in_arcs[k] = current_arc;
1037 }
1038 else
1039 {
1040 /* No. So we're done with this candidate arc */
1041 LWDEBUG(4, "pt_continues_arc = false");
1042 current_arc++;
1043 break;
1044 }
1045
1046 memcpy(&a1, &a2, sizeof(POINT4D));
1047 memcpy(&a2, &a3, sizeof(POINT4D));
1048 memcpy(&a3, &b, sizeof(POINT4D));
1049 }
1050 /* Jump past all the edges that were added to the arc */
1051 if ( found_arc )
1052 {
1053 /* Check if an arc was composed by enough edges to be
1054 * really considered an arc
1055 * See http://trac.osgeo.org/postgis/ticket/2420
1056 */
1057 arc_edges = j - 1 - i;
1058 LWDEBUGF(4, "arc defined by %d edges found", arc_edges);
1059 if ( first.x == b.x && first.y == b.y ) {
1060 LWDEBUG(4, "arc is a circle");
1061 num_quadrants = 4;
1062 }
1063 else {
1064 lw_arc_center((POINT2D*)&first, (POINT2D*)&b, (POINT2D*)&a1, (POINT2D*)&center);
1065 angle = lw_arc_angle((POINT2D*)&first, (POINT2D*)&center, (POINT2D*)&b);
1066 int p2_side = lw_segment_side((POINT2D*)&first, (POINT2D*)&a1, (POINT2D*)&b);
1067 if ( p2_side >= 0 ) angle = -angle;
1068
1069 if ( angle < 0 ) angle = 2 * M_PI + angle;
1070 num_quadrants = ( 4 * angle ) / ( 2 * M_PI );
1071 LWDEBUGF(4, "arc angle (%g %g, %g %g, %g %g) is %g (side is %d), quadrants:%g", first.x, first.y, center.x, center.y, b.x, b.y, angle, p2_side, num_quadrants);
1072 }
1073 /* a1 is first point, b is last point */
1074 if ( arc_edges < min_quad_edges * num_quadrants ) {
1075 LWDEBUGF(4, "Not enough edges for a %g quadrants arc, %g needed", num_quadrants, min_quad_edges * num_quadrants);
1076 for ( k = j-1; k >= i; k-- )
1077 edges_in_arcs[k] = 0;
1078 }
1079
1080 i = j-1;
1081 }
1082 else
1083 {
1084 /* Mark this edge as a linear edge */
1085 edges_in_arcs[i] = 0;
1086 i = i+1;
1087 }
1088 }
1089
1090#if POSTGIS_DEBUG_LEVEL > 3
1091 {
1092 char *edgestr = lwalloc(num_edges+1);
1093 for ( i = 0; i < num_edges; i++ )
1094 {
1095 if ( edges_in_arcs[i] )
1096 edgestr[i] = 48 + edges_in_arcs[i];
1097 else
1098 edgestr[i] = '.';
1099 }
1100 edgestr[num_edges] = 0;
1101 LWDEBUGF(3, "edge pattern %s", edgestr);
1102 lwfree(edgestr);
1103 }
1104#endif
1105
1106 start = 0;
1107 edge_type = edges_in_arcs[0];
1108 outcol = lwcollection_construct_empty(COMPOUNDTYPE, srid, ptarray_has_z(points), ptarray_has_m(points));
1109 for( i = 1; i < num_edges; i++ )
1110 {
1111 if( edge_type != edges_in_arcs[i] )
1112 {
1113 end = i - 1;
1114 lwcollection_add_lwgeom(outcol, geom_from_pa(points, srid, edge_type, start, end));
1115 start = i;
1116 edge_type = edges_in_arcs[i];
1117 }
1118 }
1119 lwfree(edges_in_arcs); /* not needed anymore */
1120
1121 /* Roll out last item */
1122 end = num_edges - 1;
1123 lwcollection_add_lwgeom(outcol, geom_from_pa(points, srid, edge_type, start, end));
1124
1125 /* Strip down to singleton if only one entry */
1126 if ( outcol->ngeoms == 1 )
1127 {
1128 LWGEOM *outgeom = outcol->geoms[0];
1129 outcol->ngeoms = 0; lwcollection_free(outcol);
1130 return outgeom;
1131 }
1132 return lwcollection_as_lwgeom(outcol);
1133}
#define LW_FALSE
Definition liblwgeom.h:108
#define COMPOUNDTYPE
Definition liblwgeom.h:124
void * lwalloc(size_t size)
Definition lwutil.c:227
void lwfree(void *mem)
Definition lwutil.c:242
void lwcollection_free(LWCOLLECTION *col)
int getPoint4d_p(const POINTARRAY *pa, uint32_t n, POINT4D *point)
Definition lwgeom_api.c:125
LWCOLLECTION * lwcollection_construct_empty(uint8_t type, int32_t srid, char hasz, char hasm)
#define LW_TRUE
Return types for functions with status returns.
Definition liblwgeom.h:107
LWCOLLECTION * lwcollection_add_lwgeom(LWCOLLECTION *col, const LWGEOM *geom)
Appends geom to the collection managed by col.
LWGEOM * lwcollection_as_lwgeom(const LWCOLLECTION *obj)
Definition lwgeom.c:291
double lw_arc_center(const POINT2D *p1, const POINT2D *p2, const POINT2D *p3, POINT2D *result)
Determines the center of the circle defined by the three given points.
int ptarray_has_z(const POINTARRAY *pa)
Definition ptarray.c:37
int lw_segment_side(const POINT2D *p1, const POINT2D *p2, const POINT2D *q)
lw_segment_side()
Definition lwalgorithm.c:65
int ptarray_has_m(const POINTARRAY *pa)
Definition ptarray.c:44
#define LWDEBUG(level, msg)
Definition lwgeom_log.h:83
#define LWDEBUGF(level, msg,...)
Definition lwgeom_log.h:88
void lwerror(const char *fmt,...)
Write a notice out to the error handler.
Definition lwutil.c:190
static LWGEOM * geom_from_pa(const POINTARRAY *pa, int32_t srid, int is_arc, int start, int end)
Definition lwstroke.c:965
static double lw_arc_angle(const POINT2D *a, const POINT2D *b, const POINT2D *c)
Return ABC angle in radians TODO: move to lwalgorithm.
Definition lwstroke.c:869
static int pt_continues_arc(const POINT4D *a1, const POINT4D *a2, const POINT4D *a3, const POINT4D *b)
Returns LW_TRUE if b is on the arc formed by a1/a2/a3, but not within that portion already described ...
Definition lwstroke.c:891
uint32_t ngeoms
Definition liblwgeom.h:566
LWGEOM ** geoms
Definition liblwgeom.h:561
double y
Definition liblwgeom.h:376
double x
Definition liblwgeom.h:376
double x
Definition liblwgeom.h:400
double y
Definition liblwgeom.h:400
uint32_t npoints
Definition liblwgeom.h:413

References COMPOUNDTYPE, geom_from_pa(), LWCOLLECTION::geoms, getPoint4d_p(), lw_arc_angle(), lw_arc_center(), LW_FALSE, lw_segment_side(), LW_TRUE, lwalloc(), lwcollection_add_lwgeom(), lwcollection_as_lwgeom(), lwcollection_construct_empty(), lwcollection_free(), LWDEBUG, LWDEBUGF, lwerror(), lwfree(), LWCOLLECTION::ngeoms, POINTARRAY::npoints, pt_continues_arc(), ptarray_has_m(), ptarray_has_z(), POINT2D::x, POINT4D::x, POINT2D::y, and POINT4D::y.

Referenced by lwline_unstroke(), and lwpolygon_unstroke().

Here is the call graph for this function:
Here is the caller graph for this function: