On 7/20/06, Steve Lime <[EMAIL PROTECTED]> wrote:
Puneet: Why is speed such a big deal or will you be doing this often?
I will be doing this a few times a year (so, often, although not
necessarily often enough), but this portion is only a part of a longer
process. So, I have to optimize each step as much as I can for overall
better rate of return.
Anyway, if it were me I'd just write a quick MapScript script using
perl. Create a new shapefileObj, a new XBase table, open the file and
start looping. The loop itself should have like 3 lines of code. The
whole script should be about 15...
yup, here is what I have right now (using Geo::ShapeFile) --
foreach point
foreach polygon
pointIsInPoly($point, $poly)
pointIsInPoly {
my ($x, $y) = @$p; # point being tested
my @xy = @$poly; # the poly being tested against
my $n = @xy / 2; # Number of points in polygon.
my @i = map {2 * $_} 0 .. (@xy/2); # The even indices of @xy.
my @x = map {$xy[$_ ]} @i; # Even indices: x-coordinates.
my @y = map {$xy[$_ + 1]} @i; # Odd indices: y-coordinates.
my ($i, $j); # Indices.
my $side = 0; # 0 = outside, 1 = inside.
for ($i = 0, $j = $n - 1 ; $i < $n; $j = $i++) {
if (
(
# If the y is between the (y-) borders ...
(($y[$i] <= $y) && ($y < $y[$j])) ||
(($y[$j] <= $y) && ($y < $y[$i]))
)
and
# ...the (x,y) to infinity line crosses the edge
# from the ith point to the jth point...
($x
<
($x[$j] - $x[$i] ) *
($y - $y[$i]) / ($y[$j] - $y[$i]) + $x[$i] )) {
$side = not $side; # Jump the fence.
}
}
return $side ? 1 : 0;
}
>>> P Kishor <[EMAIL PROTECTED]> 7/20/2006 7:19:06 AM >>>
Folks,
I have seen this question asked in the past, and have a similar
problem. I want to convert many lon-lat to a shapefile of points
quickly. The lon-lat are in a csv text file (values.txt) like so
foo,bar,baz,long,lat,qux
1, bar1, baz1, -87.796341, 41.907504, s10
...
...
So, I installed FW_Tools 1.05, created a DSN out of the text file
(this is and has to be on a Windows box), and ogrinfo happily gave me
information that I wanted, but with a twist...
OGRFeature(values.txt):1056
foo (Integer) = 1
bar(String) = ba1
baz(String) = baz1
long (Real) = -87.796341
lat (Real) = 41.907504
qux(String) = 10.0000
So, what's with interpreting qux? It converted my "s10" into a string
of value "10.0000".
Anyway, then I created an OVF file (values.ovf) like so
<OGRVRTDataSource>
<OGRVRTLayer name="values">
<SrcDataSource>ODBC:values</SrcDataSource>
<SrcLayer>values.txt</SrcLayer>
<GeometryField encoding="PointFromColumns" x="Long" y="Lat"/>
<GeometryType>wkbPoint</GeometryType>
</OGRVRTLayer>
</OGRVRTDataSource>
well,
ogr2ogr -f "ESRI Shapefile" . values.ovf values
but while that seemed to do something without any error, I got no
output.
So, is what I am doing the right way of going about this? Is there a
better or another way to try out and benchmark? Oh, did I mention? --
The values.txt file has about 5.25 million rows, so speed is
important.
--
Puneet Kishor http://punkish.eidesis.org/
Nelson Inst. for Env. Studies, UW-Madison http://www.ies.wisc.edu/
Open Source Geospatial Foundation https://edu.osgeo.org/
