#include <iostream>
#include <vector>

using namespace std;

#include <geos.h>

#include <DistanceOp.h> 
#include <EnhancedPrecisionOp.h> 

void main()
{
	GeometryFactory GF;

	geos::geom:: Polygon *PG1,*PG2;

	LinearRing *LR;

	vector<Geometry*> holes;

	CoordinateArraySequence *pts;

	pts=new CoordinateArraySequence();

	pts->add(Coordinate(0,0,0));
	pts->add(Coordinate(0,1,0));
	pts->add(Coordinate(1,1,0));
	pts->add(Coordinate(1,0,0));
	pts->add(Coordinate(0,0,0));

	LR=GF.createLinearRing(*pts);

	delete pts;

	PG1=GF.createPolygon(*LR,holes); 

	GF.destroyGeometry(LR);

	pts=new CoordinateArraySequence();

	pts->add(Coordinate(1.25,0.25,0));
	pts->add(Coordinate(1.25,0.75,0));
	pts->add(Coordinate(1.75,0.75,0));
	pts->add(Coordinate(1.75,0.25,0));
	pts->add(Coordinate(1.25,0.25,0));

	LR=GF.createLinearRing(*pts);

	delete pts;

	PG2=GF.createPolygon(*LR,holes); 

	GF.destroyGeometry(LR);

	CoordinateSequence *cs;

	cs=PG1->getCoordinates();

	for(size_t i=0;i<cs->size();i++) cout<<cs->getX(i)<<" "<<cs->getY(i)<<endl;

	cout<<endl;

	cs=PG2->getCoordinates();

	for(size_t i=0;i<cs->size();i++) cout<<cs->getX(i)<<" "<<cs->getY(i)<<endl;

	cout<<endl;

	getchar();

	cout<<PG1->distance(PG2)<<endl;

	cout<<PG2->distance(PG1)<<endl;

	typedef geos::operation::distance::DistanceOp DOP;

	DOP *Dis = new DOP(PG1, PG2);

	cout<<"in"<<endl;
	cs=Dis->closestPoints();
	cout<<"out"<<endl;

	cout<<cs->getSize()<<endl;

	for(size_t i=0;i<cs->size();i++) cout<<cs->getX(i)<<" "<<cs->getY(i)<<endl;

	delete Dis;

	GF.destroyGeometry(PG1);	

	GF.destroyGeometry(PG2);	

	return;
};