static LWLINE* lwline_from_gserialized_buffer(uint8_t *data_ptr, uint8_t g_flags, size_t *g_size)
{
	uint8_t *start_ptr = data_ptr;
	LWLINE *line;
	uint32_t npoints = 0;

	assert(data_ptr);

	line = (LWLINE*)lwalloc(sizeof(LWLINE));
	line->srid = SRID_UNKNOWN; /* Default */
	line->bbox = NULL;
	line->type = LINETYPE;
	line->flags = g_flags;

	data_ptr += 4; /* Skip past the type. */
	npoints = lw_get_uint32_t(data_ptr); /* Zero => empty geometry */
	data_ptr += 4; /* Skip past the npoints. */

	if ( npoints > 0 )
		line->points = ptarray_construct_reference_data(FLAGS_GET_Z(g_flags), FLAGS_GET_M(g_flags), npoints, data_ptr);
		
	else
		line->points = ptarray_construct(FLAGS_GET_Z(g_flags), FLAGS_GET_M(g_flags), 0); /* Empty linestring */

	data_ptr += FLAGS_NDIMS(g_flags) * npoints * sizeof(double);

	if ( g_size )
		*g_size = data_ptr - start_ptr;

	return line;
}

static LWCIRCSTRING* lwcircstring_from_gserialized_buffer(uint8_t *data_ptr, uint8_t g_flags, size_t *g_size)
{
	uint8_t *start_ptr = data_ptr;
	LWCIRCSTRING *circstring;
	uint32_t npoints = 0;

	assert(data_ptr);

	circstring = (LWCIRCSTRING*)lwalloc(sizeof(LWCIRCSTRING));
	circstring->srid = SRID_UNKNOWN; /* Default */
	circstring->bbox = NULL;
	circstring->type = CIRCSTRINGTYPE;
	circstring->flags = g_flags;

	data_ptr += 4; /* Skip past the circstringtype. */
	npoints = lw_get_uint32_t(data_ptr); /* Zero => empty geometry */
	data_ptr += 4; /* Skip past the npoints. */

	if ( npoints > 0 )
		circstring->points = ptarray_construct_reference_data(FLAGS_GET_Z(g_flags), FLAGS_GET_M(g_flags), npoints, data_ptr);		
	else
		circstring->points = ptarray_construct(FLAGS_GET_Z(g_flags), FLAGS_GET_M(g_flags), 0); /* Empty circularstring */

	data_ptr += FLAGS_NDIMS(g_flags) * npoints * sizeof(double);

	if ( g_size )
		*g_size = data_ptr - start_ptr;

	return circstring;
}

char* gbox_to_string(const GBOX *gbox)
{
	static int sz = 128;
	char *str = NULL;

	if ( ! gbox )
		return strdup("NULL POINTER");

	str = (char*)lwalloc(sz);

	if ( FLAGS_GET_GEODETIC(gbox->flags) )
	{
		snprintf(str, sz, "GBOX((%.8g,%.8g,%.8g),(%.8g,%.8g,%.8g))", gbox->xmin, gbox->ymin, gbox->zmin, gbox->xmax, gbox->ymax, gbox->zmax);
		return str;
	}
	if ( FLAGS_GET_Z(gbox->flags) && FLAGS_GET_M(gbox->flags) )
	{
		snprintf(str, sz, "GBOX((%.8g,%.8g,%.8g,%.8g),(%.8g,%.8g,%.8g,%.8g))", gbox->xmin, gbox->ymin, gbox->zmin, gbox->mmin, gbox->xmax, gbox->ymax, gbox->zmax, gbox->mmax);
		return str;
	}
	if ( FLAGS_GET_Z(gbox->flags) )
	{
		snprintf(str, sz, "GBOX((%.8g,%.8g,%.8g),(%.8g,%.8g,%.8g))", gbox->xmin, gbox->ymin, gbox->zmin, gbox->xmax, gbox->ymax, gbox->zmax);
		return str;
	}
	if ( FLAGS_GET_M(gbox->flags) )
	{
		snprintf(str, sz, "GBOX((%.8g,%.8g,%.8g),(%.8g,%.8g,%.8g))", gbox->xmin, gbox->ymin, gbox->mmin, gbox->xmax, gbox->ymax, gbox->mmax);
		return str;
	}
	snprintf(str, sz, "GBOX((%.8g,%.8g),(%.8g,%.8g))", gbox->xmin, gbox->ymin, gbox->xmax, gbox->ymax);
	return str;
}

int gbox_overlaps(const GBOX *g1, const GBOX *g2)
{

	/* Make sure our boxes are consistent */
	if ( FLAGS_GET_GEODETIC(g1->flags) != FLAGS_GET_GEODETIC(g2->flags) )
		lwerror("gbox_overlaps: cannot compare geodetic and non-geodetic boxes");

	/* Check X/Y first */
	if ( g1->xmax < g2->xmin || g1->ymax < g2->ymin ||
	     g1->xmin > g2->xmax || g1->ymin > g2->ymax )
		return LW_FALSE;
		
	/* If both geodetic or both have Z, check Z */
	if ( (FLAGS_GET_Z(g1->flags) && FLAGS_GET_Z(g2->flags)) || 
	     (FLAGS_GET_GEODETIC(g1->flags) && FLAGS_GET_GEODETIC(g2->flags)) )
	{
		if ( g1->zmax < g2->zmin || g1->zmin > g2->zmax )
			return LW_FALSE;
	}
	
	/* If both have M, check M */
	if ( FLAGS_GET_M(g1->flags) && FLAGS_GET_M(g2->flags) )
	{
		if ( g1->mmax < g2->mmin || g1->mmin > g2->mmax )
			return LW_FALSE;
	}
	
	return LW_TRUE;
}

static LWPOINT* lwpoint_from_gserialized_buffer(uint8_t *data_ptr, uint8_t g_flags, size_t *g_size)
{
	uint8_t *start_ptr = data_ptr;
	LWPOINT *point;
	uint32_t npoints = 0;

	assert(data_ptr);

	point = (LWPOINT*)lwalloc(sizeof(LWPOINT));
	point->srid = SRID_UNKNOWN; /* Default */
	point->bbox = NULL;
	point->type = POINTTYPE;
	point->flags = g_flags;

	data_ptr += 4; /* Skip past the type. */
	npoints = lw_get_uint32_t(data_ptr); /* Zero => empty geometry */
	data_ptr += 4; /* Skip past the npoints. */

	if ( npoints > 0 )
		point->point = ptarray_construct_reference_data(FLAGS_GET_Z(g_flags), FLAGS_GET_M(g_flags), 1, data_ptr);
	else
		point->point = ptarray_construct(FLAGS_GET_Z(g_flags), FLAGS_GET_M(g_flags), 0); /* Empty point */

	data_ptr += npoints * FLAGS_NDIMS(g_flags) * sizeof(double);

	if ( g_size )
		*g_size = data_ptr - start_ptr;

	return point;
}

Datum LWGEOM_m_point(PG_FUNCTION_ARGS)
{
	GSERIALIZED *geom;
	LWPOINT *point = NULL;
	LWGEOM *lwgeom;
	POINT3DM p;

	geom = (GSERIALIZED *)PG_DETOAST_DATUM(PG_GETARG_DATUM(0));

	if ( gserialized_get_type(geom) != POINTTYPE )
		lwerror("Argument to M() must be a point");

	lwgeom = lwgeom_from_gserialized(geom);
	point = lwgeom_as_lwpoint(lwgeom);
	
	if ( lwgeom_is_empty(lwgeom) )
		PG_RETURN_NULL();

	/* no M in input */
	if ( ! FLAGS_GET_M(point->flags) ) PG_RETURN_NULL();

	getPoint3dm_p(point->point, 0, &p);

	PG_FREE_IF_COPY(geom, 0);

	PG_RETURN_FLOAT8(p.m);
}

/**
 * @brief Merge two given POINTARRAY and returns a pointer
 * on the new aggregate one.
 * Warning: this function free the two inputs POINTARRAY
 * @return #POINTARRAY is newly allocated
 */
POINTARRAY *
ptarray_merge(POINTARRAY *pa1, POINTARRAY *pa2)
{
	POINTARRAY *pa;
	size_t ptsize = ptarray_point_size(pa1);

	if (FLAGS_GET_ZM(pa1->flags) != FLAGS_GET_ZM(pa2->flags))
		lwerror("ptarray_cat: Mixed dimension");

	pa = ptarray_construct( FLAGS_GET_Z(pa1->flags),
	                        FLAGS_GET_M(pa1->flags),
	                        pa1->npoints + pa2->npoints);

	memcpy(         getPoint_internal(pa, 0),
	                getPoint_internal(pa1, 0),
	                ptsize*(pa1->npoints));

	memcpy(         getPoint_internal(pa, pa1->npoints),
	                getPoint_internal(pa2, 0),
	                ptsize*(pa2->npoints));

	lwfree(pa1);
	lwfree(pa2);

	return pa;
}

static int lwcircstring_calculate_gbox_cartesian(LWCIRCSTRING *curve, GBOX *gbox)
{
	uint8_t flags = gflags(FLAGS_GET_Z(curve->flags), FLAGS_GET_M(curve->flags), 0);
	GBOX tmp;
	POINT4D p1, p2, p3;
	int i;

	if ( ! curve ) return LW_FAILURE;
	if ( curve->points->npoints < 3 ) return LW_FAILURE;

	tmp.flags = flags;

	/* Initialize */
	gbox->xmin = gbox->ymin = gbox->zmin = gbox->mmin = MAXFLOAT;
	gbox->xmax = gbox->ymax = gbox->zmax = gbox->mmax = -1 * MAXFLOAT;

	for ( i = 2; i < curve->points->npoints; i += 2 )
	{
		getPoint4d_p(curve->points, i-2, &p1);
		getPoint4d_p(curve->points, i-1, &p2);
		getPoint4d_p(curve->points, i, &p3);

		if (lw_arc_calculate_gbox_cartesian(&p1, &p2, &p3, &tmp) == LW_FAILURE)
			continue;

		gbox_merge(&tmp, gbox);
	}

	return LW_SUCCESS;
}

/**
* Create a new point. Null point array implies empty. Null dimensionality
* implies no specified dimensionality in the WKT.
*/
LWGEOM* wkt_parser_point_new(POINTARRAY *pa, char *dimensionality)
{
	uint8_t flags = wkt_dimensionality(dimensionality);
	LWDEBUG(4,"entered");
	
	/* No pointarray means it is empty */
	if( ! pa )
		return lwpoint_as_lwgeom(lwpoint_construct_empty(SRID_UNKNOWN, FLAGS_GET_Z(flags), FLAGS_GET_M(flags)));

	/* If the number of dimensions is not consistent, we have a problem. */
	if( wkt_pointarray_dimensionality(pa, flags) == LW_FALSE )
	{
		ptarray_free(pa);
		SET_PARSER_ERROR(PARSER_ERROR_MIXDIMS);
		return NULL;
	}

	/* Only one point allowed in our point array! */	
	if( pa->npoints != 1 )
	{
		ptarray_free(pa);
		SET_PARSER_ERROR(PARSER_ERROR_LESSPOINTS);
		return NULL;
	}		

	return lwpoint_as_lwgeom(lwpoint_construct(SRID_UNKNOWN, NULL, pa));
}

LWGEOM* wkt_parser_polygon_finalize(LWGEOM *poly, char *dimensionality)
{
	uint8_t flags = wkt_dimensionality(dimensionality);
	int flagdims = FLAGS_NDIMS(flags);
	LWDEBUG(4,"entered");
	
	/* Null input implies empty return */
	if( ! poly )
		return lwpoly_as_lwgeom(lwpoly_construct_empty(SRID_UNKNOWN, FLAGS_GET_Z(flags), FLAGS_GET_M(flags)));

	/* If the number of dimensions are not consistent, we have a problem. */
	if( flagdims > 2 )
	{
		if ( flagdims != FLAGS_NDIMS(poly->flags) )
		{
			lwgeom_free(poly);
			SET_PARSER_ERROR(PARSER_ERROR_MIXDIMS);
			return NULL;
		}
	
		/* Harmonize the flags in the sub-components with the wkt flags */
		if( LW_FAILURE == wkt_parser_set_dims(poly, flags) )
		{
			lwgeom_free(poly);
			SET_PARSER_ERROR(PARSER_ERROR_OTHER);
			return NULL;
		}
	}
	
	return poly;
}

/**
* Read the dimensionality from a flag, if provided. Then check that the
* dimensionality matches that of the pointarray. If the dimension counts
* match, ensure the pointarray is using the right "Z" or "M".
*/
static int wkt_pointarray_dimensionality(POINTARRAY *pa, uint8_t flags)
{	
	int hasz = FLAGS_GET_Z(flags);
	int hasm = FLAGS_GET_M(flags);
	int ndims = 2 + hasz + hasm;

	/* No dimensionality or array means we go with what we have */
	if( ! (flags && pa) )
		return LW_TRUE;
		
	LWDEBUGF(5,"dimensionality ndims == %d", ndims);
	LWDEBUGF(5,"FLAGS_NDIMS(pa->flags) == %d", FLAGS_NDIMS(pa->flags));
	
	/* 
	* ndims > 2 implies that the flags have something useful to add,
	* that there is a 'Z' or an 'M' or both.
	*/
	if( ndims > 2 )
	{
		/* Mismatch implies a problem */
		if ( FLAGS_NDIMS(pa->flags) != ndims )
			return LW_FALSE;
		/* Match means use the explicit dimensionality */
		else
		{
			FLAGS_SET_Z(pa->flags, hasz);
			FLAGS_SET_M(pa->flags, hasm);
		}
	}

	return LW_TRUE;
}

int
ptarray_append_point(POINTARRAY *pa, const POINT4D *pt, int repeated_points)
{

	/* Check for pathology */
	if( ! pa || ! pt ) 
	{
		lwerror("ptarray_append_point: null input");
		return LW_FAILURE;
	}

	/* Check for duplicate end point */
	if ( repeated_points == LW_FALSE && pa->npoints > 0 )
	{
		POINT4D tmp;
		getPoint4d_p(pa, pa->npoints-1, &tmp);
		LWDEBUGF(4,"checking for duplicate end point (pt = POINT(%g %g) pa->npoints-q = POINT(%g %g))",pt->x,pt->y,tmp.x,tmp.y);

		/* Return LW_SUCCESS and do nothing else if previous point in list is equal to this one */
		if ( (pt->x == tmp.x) && (pt->y == tmp.y) &&
		     (FLAGS_GET_Z(pa->flags) ? pt->z == tmp.z : 1) &&
		     (FLAGS_GET_M(pa->flags) ? pt->m == tmp.m : 1) )
		{
			return LW_SUCCESS;
		}
	}

	/* Append is just a special case of insert */
	return ptarray_insert_point(pa, pt, pa->npoints);
}

int gbox_is_valid(const GBOX *gbox)
{
	/* X */
	if ( ! isfinite(gbox->xmin) || isnan(gbox->xmin) ||
	     ! isfinite(gbox->xmax) || isnan(gbox->xmax) )
		return LW_FALSE;
		
	/* Y */
	if ( ! isfinite(gbox->ymin) || isnan(gbox->ymin) ||
	     ! isfinite(gbox->ymax) || isnan(gbox->ymax) )
		return LW_FALSE;
		
	/* Z */
	if ( FLAGS_GET_GEODETIC(gbox->flags) || FLAGS_GET_Z(gbox->flags) )
	{
		if ( ! isfinite(gbox->zmin) || isnan(gbox->zmin) ||
		     ! isfinite(gbox->zmax) || isnan(gbox->zmax) )
			return LW_FALSE;
	}

	/* M */
	if ( FLAGS_GET_M(gbox->flags) )
	{
		if ( ! isfinite(gbox->mmin) || isnan(gbox->mmin) ||
		     ! isfinite(gbox->mmax) || isnan(gbox->mmax) )
			return LW_FALSE;
	}
	
	return LW_TRUE;		
}

int gbox_merge(const GBOX *new_box, GBOX *merge_box)
{
	assert(merge_box);

	if ( FLAGS_GET_ZM(merge_box->flags) != FLAGS_GET_ZM(new_box->flags) )
		return LW_FAILURE;

	if ( new_box->xmin < merge_box->xmin) merge_box->xmin = new_box->xmin;
	if ( new_box->ymin < merge_box->ymin) merge_box->ymin = new_box->ymin;
	if ( new_box->xmax > merge_box->xmax) merge_box->xmax = new_box->xmax;
	if ( new_box->ymax > merge_box->ymax) merge_box->ymax = new_box->ymax;

	if ( FLAGS_GET_Z(merge_box->flags) || FLAGS_GET_GEODETIC(merge_box->flags) )
	{
		if ( new_box->zmin < merge_box->zmin) merge_box->zmin = new_box->zmin;
		if ( new_box->zmax > merge_box->zmax) merge_box->zmax = new_box->zmax;
	}
	if ( FLAGS_GET_M(merge_box->flags) )
	{
		if ( new_box->mmin < merge_box->mmin) merge_box->mmin = new_box->mmin;
		if ( new_box->mmax > merge_box->mmax) merge_box->mmax = new_box->mmax;
	}

	return LW_SUCCESS;
}

/**
* Create a new linestring. Null point array implies empty. Null dimensionality
* implies no specified dimensionality in the WKT. Check for numpoints >= 2 if
* requested.
*/
LWGEOM* wkt_parser_linestring_new(POINTARRAY *pa, char *dimensionality)
{
	uint8_t flags = wkt_dimensionality(dimensionality);
	LWDEBUG(4,"entered");

	/* No pointarray means it is empty */
	if( ! pa )
		return lwline_as_lwgeom(lwline_construct_empty(SRID_UNKNOWN, FLAGS_GET_Z(flags), FLAGS_GET_M(flags)));

	/* If the number of dimensions is not consistent, we have a problem. */
	if( wkt_pointarray_dimensionality(pa, flags) == LW_FALSE )
	{
		ptarray_free(pa);
		SET_PARSER_ERROR(PARSER_ERROR_MIXDIMS);
		return NULL;
	}
	
	/* Apply check for not enough points, if requested. */	
	if( (global_parser_result.parser_check_flags & LW_PARSER_CHECK_MINPOINTS) && (pa->npoints < 2) )
	{
		ptarray_free(pa);
		SET_PARSER_ERROR(PARSER_ERROR_MOREPOINTS);
		return NULL;
	}

	return lwline_as_lwgeom(lwline_construct(SRID_UNKNOWN, NULL, pa));
}

/**
* Peak into a #GSERIALIZED datum to find the bounding box. If the
* box is there, copy it out and return it. If not, calculate the box from the
* full object and return the box based on that. If no box is available,
* return #LW_FAILURE, otherwise #LW_SUCCESS.
*/
int 
gserialized_datum_get_gidx_p(Datum gsdatum, GIDX *gidx)
{
	GSERIALIZED *gpart;
	int result = LW_SUCCESS;

	POSTGIS_DEBUG(4, "entered function");

	/*
	** The most info we need is the 8 bytes of serialized header plus the 32 bytes
	** of floats necessary to hold the 8 floats of the largest XYZM index
	** bounding box, so 40 bytes.
	*/
	gpart = (GSERIALIZED*)PG_DETOAST_DATUM_SLICE(gsdatum, 0, 40);

	POSTGIS_DEBUGF(4, "got flags %d", gpart->flags);

	/* Do we even have a serialized bounding box? */
	if ( FLAGS_GET_BBOX(gpart->flags) )
	{
		/* Yes! Copy it out into the GIDX! */
		size_t size = gbox_serialized_size(gpart->flags);
		POSTGIS_DEBUG(4, "copying box out of serialization");
		memcpy(gidx->c, gpart->data, size);
		/* if M is present but Z is not, pad Z and shift M */
		if ( FLAGS_GET_M(gpart->flags) && ! FLAGS_GET_Z(gpart->flags) )
		{
			size += 2 * sizeof(float);
			GIDX_SET_MIN(gidx,3,GIDX_GET_MIN(gidx,2));
			GIDX_SET_MAX(gidx,3,GIDX_GET_MAX(gidx,2));
			GIDX_SET_MIN(gidx,2,-1*FLT_MAX);
			GIDX_SET_MAX(gidx,2,FLT_MAX);
		}
		SET_VARSIZE(gidx, VARHDRSZ + size);
		result = LW_SUCCESS;
	}
	else
	{
		/* No, we need to calculate it from the full object. */
		GSERIALIZED *g = (GSERIALIZED*)PG_DETOAST_DATUM(gsdatum);
		LWGEOM *lwgeom = lwgeom_from_gserialized(g);
		GBOX gbox;
		if ( lwgeom_calculate_gbox(lwgeom, &gbox) == LW_FAILURE )
		{
			POSTGIS_DEBUG(4, "could not calculate bbox, returning failure");
			lwgeom_free(lwgeom);
			return LW_FAILURE;
		}
		lwgeom_free(lwgeom);
		result = gidx_from_gbox_p(gbox, gidx);
	}
	
	if ( result == LW_SUCCESS )
	{
		POSTGIS_DEBUGF(4, "got gidx %s", gidx_to_string(gidx));
	}

	return result;
}

LWGEOM* wkt_parser_collection_finalize(int lwtype, LWGEOM *col, char *dimensionality) 
{
	uint8_t flags = wkt_dimensionality(dimensionality);
	int flagdims = FLAGS_NDIMS(flags);
	
	/* No geometry means it is empty */
	if( ! col )
	{
		return lwcollection_as_lwgeom(lwcollection_construct_empty(lwtype, SRID_UNKNOWN, FLAGS_GET_Z(flags), FLAGS_GET_M(flags)));
	}

	/* There are 'Z' or 'M' tokens in the signature */
	if ( flagdims > 2 )
	{
		/* If the number of dimensions are not consistent, we have a problem. */
		if( flagdims != FLAGS_NDIMS(col->flags) )
		{
			lwgeom_free(col);
			SET_PARSER_ERROR(PARSER_ERROR_MIXDIMS);
			return NULL;
		}

		/* For GEOMETRYCOLLECTION, the exact type of the dimensions must match too */
		if( lwtype == COLLECTIONTYPE &&
		    ( (FLAGS_GET_Z(flags) != FLAGS_GET_Z(col->flags)) ||
		      (FLAGS_GET_M(flags) != FLAGS_GET_M(col->flags)) ) )
		{
			lwgeom_free(col);
			SET_PARSER_ERROR(PARSER_ERROR_MIXDIMS);
			return NULL;
		}
		
		/* Harmonize the collection dimensionality */
		if( LW_FAILURE == wkt_parser_set_dims(col, flags) )
		{
			lwgeom_free(col);
			SET_PARSER_ERROR(PARSER_ERROR_OTHER);
			return NULL;
		}
	}
		
	/* Set the collection type */
	col->type=lwtype;
			
	return col;
}

/*
* Generic GEOMETRY
*/
static void lwgeom_to_wkt_sb(const LWGEOM *geom, stringbuffer_t *sb, int precision, uint8_t variant)
{
    LWDEBUGF(4, "lwgeom_to_wkt_sb: type %s, hasz %d, hasm %d",
             lwtype_name(geom->type), (geom->type),
             FLAGS_GET_Z(geom->flags)?1:0, FLAGS_GET_M(geom->flags)?1:0);

    switch (geom->type)
    {
    case POINTTYPE:
        lwpoint_to_wkt_sb((LWPOINT*)geom, sb, precision, variant);
        break;
    case LINETYPE:
        lwline_to_wkt_sb((LWLINE*)geom, sb, precision, variant);
        break;
    case POLYGONTYPE:
        lwpoly_to_wkt_sb((LWPOLY*)geom, sb, precision, variant);
        break;
    case MULTIPOINTTYPE:
        lwmpoint_to_wkt_sb((LWMPOINT*)geom, sb, precision, variant);
        break;
    case MULTILINETYPE:
        lwmline_to_wkt_sb((LWMLINE*)geom, sb, precision, variant);
        break;
    case MULTIPOLYGONTYPE:
        lwmpoly_to_wkt_sb((LWMPOLY*)geom, sb, precision, variant);
        break;
    case COLLECTIONTYPE:
        lwcollection_to_wkt_sb((LWCOLLECTION*)geom, sb, precision, variant);
        break;
    case CIRCSTRINGTYPE:
        lwcircstring_to_wkt_sb((LWCIRCSTRING*)geom, sb, precision, variant);
        break;
    case COMPOUNDTYPE:
        lwcompound_to_wkt_sb((LWCOMPOUND*)geom, sb, precision, variant);
        break;
    case CURVEPOLYTYPE:
        lwcurvepoly_to_wkt_sb((LWCURVEPOLY*)geom, sb, precision, variant);
        break;
    case MULTICURVETYPE:
        lwmcurve_to_wkt_sb((LWMCURVE*)geom, sb, precision, variant);
        break;
    case MULTISURFACETYPE:
        lwmsurface_to_wkt_sb((LWMSURFACE*)geom, sb, precision, variant);
        break;
    case TRIANGLETYPE:
        lwtriangle_to_wkt_sb((LWTRIANGLE*)geom, sb, precision, variant);
        break;
    case TINTYPE:
        lwtin_to_wkt_sb((LWTIN*)geom, sb, precision, variant);
        break;
    case POLYHEDRALSURFACETYPE:
        lwpsurface_to_wkt_sb((LWPSURFACE*)geom, sb, precision, variant);
        break;
    default:
        lwerror("lwgeom_to_wkt_sb: Type %d - %s unsupported.",
                geom->type, lwtype_name(geom->type));
    }
}

POINTARRAY *
ptarray_simplify(POINTARRAY *inpts, double epsilon, unsigned int minpts)
{
	int *stack;			/* recursion stack */
	int sp=-1;			/* recursion stack pointer */
	int p1, split;
	double dist;
	POINTARRAY *outpts;
	POINT4D pt;

	double eps_sqr = epsilon * epsilon;

	/* Allocate recursion stack */
	stack = lwalloc(sizeof(int)*inpts->npoints);

	p1 = 0;
	stack[++sp] = inpts->npoints-1;

	LWDEBUGF(2, "Input has %d pts and %d dims", inpts->npoints,
	                                            FLAGS_NDIMS(inpts->flags));

	/* Allocate output POINTARRAY, and add first point. */
	outpts = ptarray_construct_empty(FLAGS_GET_Z(inpts->flags), FLAGS_GET_M(inpts->flags), inpts->npoints);
	getPoint4d_p(inpts, 0, &pt);
	ptarray_append_point(outpts, &pt, LW_FALSE);

	LWDEBUG(3, "Added P0 to simplified point array (size 1)");

	do
	{

		ptarray_dp_findsplit(inpts, p1, stack[sp], &split, &dist);

		LWDEBUGF(3, "Farthest point from P%d-P%d is P%d (dist. %g)", p1, stack[sp], split, dist);

		if (dist > eps_sqr || ( outpts->npoints+sp+1 < minpts && dist >= 0 ) )
		{
			LWDEBUGF(4, "Added P%d to stack (outpts:%d)", split, sp);
			stack[++sp] = split;
		}
		else
		{
			getPoint4d_p(inpts, stack[sp], &pt);
			LWDEBUGF(4, "npoints , minpoints %d %d", outpts->npoints, minpts);
			ptarray_append_point(outpts, &pt, LW_FALSE);
			
			LWDEBUGF(4, "Added P%d to simplified point array (size: %d)", stack[sp], outpts->npoints);

			p1 = stack[sp--];
		}

		LWDEBUGF(4, "stack pointer = %d", sp);
	}
	while (! (sp<0) );

	lwfree(stack);
	return outpts;
}

/* Convert a double-based GBOX into a float-based GIDX,
   ensuring the float box is larger than the double box */
static int gidx_from_gbox_p(GBOX box, GIDX *a)
{
	int ndims;

	ndims = (FLAGS_GET_GEODETIC(box.flags) ? 3 : FLAGS_NDIMS(box.flags));
	SET_VARSIZE(a, VARHDRSZ + ndims * 2 * sizeof(float));

	GIDX_SET_MIN(a,0,nextDown_f(box.xmin));
	GIDX_SET_MAX(a,0,nextUp_f(box.xmax));
	GIDX_SET_MIN(a,1,nextDown_f(box.ymin));
	GIDX_SET_MAX(a,1,nextUp_f(box.ymax));

	/* Geodetic indexes are always 3d, geocentric x/y/z */
	if ( FLAGS_GET_GEODETIC(box.flags) )
	{
		GIDX_SET_MIN(a,2,nextDown_f(box.zmin));
		GIDX_SET_MAX(a,2,nextUp_f(box.zmax));
	}
	else
	{
		/* Cartesian with Z implies Z is third dimension */
		if ( FLAGS_GET_Z(box.flags) )
		{
			GIDX_SET_MIN(a,2,nextDown_f(box.zmin));
			GIDX_SET_MAX(a,2,nextUp_f(box.zmax));
			if ( FLAGS_GET_M(box.flags) )
			{
				GIDX_SET_MIN(a,3,nextDown_f(box.mmin));
				GIDX_SET_MAX(a,3,nextUp_f(box.mmax));
			}
		}
		/* Unless there's no Z, in which case M is third dimension */
		else if ( FLAGS_GET_M(box.flags) )
		{
			GIDX_SET_MIN(a,2,nextDown_f(box.mmin));
			GIDX_SET_MAX(a,2,nextUp_f(box.mmax));
		}
	}

	POSTGIS_DEBUGF(5, "converted %s to %s", gbox_to_string(&box), gidx_to_string(a));

	return G_SUCCESS;
}