[Gdal-dev] GDAL Capabilities

[Chapman, Martin]

Here is the code to do it yourself:

int CIntersectDlgDlg::pointnpolygon(int npol, double *xp, double *yp, double x, double y) {

int i, j, c = 0;

for ( i = 0, j = npol - 1; i < npol; j = i++ ) 

{

if ((((yp[i] <= y) && (y < yp[j])) || ((yp[j] <= y) && (y < yp[i]))) && 

(x < (xp[j] - xp[i]) * (y - yp[i]) / (yp[j] - yp[i]) + xp[i]))

{

c = !c;

}

}

return c;

}

int CIntersectDlgDlg::linenpolygon(line l1, line l2)

{

return ((ccw(l1.p1, l1.p2, l2.p1) != ccw(l1.p1, l1.p2, l2.p2)) && (ccw(l2.p1, l2.p2, l1.p1) != ccw(l2.p1, l2.p2, l1.p2))); }

int CIntersectDlgDlg::ccw(point p1, point p2, point p3)

{

double dx1, dx2, dy1, dy2;

dx1 = p2.x - p1.x; 

dy1 = p2.y - p1.y;

dx2 = p3.x - p2.x; 

dy2 = p3.y - p2.y;

if (dy1 * dx2 < dy2 * dx1) 

return 1;

else 

return 0;

}

Let me know if you have questions or if you can make it better...

Martin

	-----Original Message----- 
	From: Frank Warmerdam [mailto:warmerdam at pobox.com] 
	Sent: Wed 1/21/2004 10:04 AM 
	To: gdal-dev at remotesensing.org; PRouse at weather.com 
	Cc: 
	Subject: Re: [Gdal-dev] GDAL Capabilities
	
	

	Patrick Rouse wrote:
	> Hi,
	>
	> I am looking for an opensource solution to a simple gis/topology problem.  I am
	> an Arcinfo guy, so I know how to do this within the realm of aml and visual
	> basic, but I have been tasked with creating a simple opensource (c++), portable
	> program to generate a dynamic polygon based on lat/lon points input by the user
	> and determine which polygons on a shapefile (us counties) fall within or
	> intersect this dynamic polygon.  A simple check to see if the county polygons in
	> the shapefile fall within or intersect  the dynamic polygon.  Browsing the
	> documentation it seems that this library is capable of this.  I am not much of a
	> c++ programmer, but I will be assisting a real c++ programmer in this project
	> and I am looking for possible solutions to present to him.  Let me know if this
	> library can handle this operation and if possible links to any sample code used
	> for opening shapefiles and performing topological tests.
	
	Patrick,
	
	GDAL's OGR API can be used to read the shapefiles for the counties, but
	GDAL/OGR does not include functions to overlay polygons on polygons, so there
	is no rigerous mechanism to determine which (if any) county polygons
	intersect your dynamic polygon.
	
	You can get an approximate answer by computing the extents of the dynamic
	polygon, and setting this rectangular extent as the SpatialFilter being
	reading back counties.  This way you will only get counties back whose
	extent rectangle intersects the extent rectangle of your dynamic polygon.
	
	However, to do the job properly you need overlay software.  I would suggest
	the GEOS library as a component to do this.
	
	   http://geos.refractions.net/
	
	It's possible that GEOS is already integrated with a shapefile reader in
	which case you don't need GDAL/OGR at all.  At some point I intend to
	integrate GEOS into OGR so that various spatial operators on OGRGeometry
	objects can be done by internally calling GEOS.   However, there is no
	schedule for that work, so there isn't any point in your waiting for it.
	
	An overall sweeter solution to your problem, might be to load your county
	shapefiles into PostGIS, and then use the spatial operators in PostGIS
	to find which ones intersect your dynamic polygon.  This could be done from
	SQL.  PostGIS is a spatial extension to the Postgres RDBMS.  The PostGIS
	spatial predicates are built on GEOS.
	
	Best regards,
	--
	---------------------------------------+--------------------------------------
	I set the clouds in motion - turn up   | Frank Warmerdam, warmerdam at pobox.com
	light and sound - activate the windows | http://pobox.com/~warmerdam
	and watch the world go round - Rush    | Geospatial Programmer for Rent
	
	
	_______________________________________________
	Gdal-dev mailing list
	Gdal-dev at remotesensing.org
	http://remotesensing.org/mailman/listinfo/gdal-dev