Implementazione MySQL dell'algoritmo di ray-casting?

Question

Dobbiamo trovare un metodo rapido e abbastanza accurato per il point-in-poligon per valori lat/long e poligoni su google maps. Dopo alcune ricerche - sono imbattuto in alcuni post su MySQL estensioni geometriche, e ha implementare anche questo -

SELECT id, Contains(PolyFromText('POLYGON(".$polygonpath.")') , PointFromText(concat(\"POINT(\", latitude, \" \", longitude, \")\"))) AS 
      CONTAINS 
FROM tbl_points 

che non ha però funzionato con poligoni composti da un gran numero di punti :(

Dopo un po 'di più ricerca - si è imbattuto in un algoritmo standard chiamato algoritmo Ray-casting, ma prima di provare a sviluppare una query per questo in MySQL, volevo rischiare se qualcuno lo avesse già fatto o si fosse imbattuto in un link utile che mostra come implementare l'algoritmo in Server MySQL/SQL

Quindi, tagliando breve - la domanda è:

Qualcuno può fornire l'implementazione di algoritmo Ray casting per MySQL/SQL?

dettaglio aggiuntive:

• poligoni sono o concava, convessa o complesso.
• Targeting per un'esecuzione rapida con una precisione del 100%.

Answer 1

La seguente funzione (versione di MySQL dell'algoritmo raycasting) ha scosso il mio mondo:

CREATE FUNCTION myWithin(p POINT, poly POLYGON) RETURNS INT(1) DETERMINISTIC 
BEGIN 
DECLARE n INT DEFAULT 0; 
DECLARE pX DECIMAL(9,6); 
DECLARE pY DECIMAL(9,6); 
DECLARE ls LINESTRING; 
DECLARE poly1 POINT; 
DECLARE poly1X DECIMAL(9,6); 
DECLARE poly1Y DECIMAL(9,6); 
DECLARE poly2 POINT; 
DECLARE poly2X DECIMAL(9,6); 
DECLARE poly2Y DECIMAL(9,6); 
DECLARE i INT DEFAULT 0; 
DECLARE result INT(1) DEFAULT 0; 
SET pX = X(p); 
SET pY = Y(p); 
SET ls = ExteriorRing(poly); 
SET poly2 = EndPoint(ls); 
SET poly2X = X(poly2); 
SET poly2Y = Y(poly2); 
SET n = NumPoints(ls); 
WHILE i<n DO 
SET poly1 = PointN(ls, (i+1)); 
SET poly1X = X(poly1); 
SET poly1Y = Y(poly1); 
IF ((((poly1X <= pX) && (pX < poly2X)) || ((poly2X <= pX) && (pX < poly1X))) && (pY > (poly2Y - poly1Y) * (pX - poly1X)/(poly2X - poly1X) + poly1Y)) THEN 
SET result = !result; 
END IF; 
SET poly2X = poly1X; 
SET poly2Y = poly1Y; 
SET i = i + 1; 
END WHILE; 
RETURN result; 
End; 

Aggiungi

DELIMITER ;; 

prima della funzione, come richiesto. L'utilizzo della funzione è:

SELECT myWithin(point, polygon) as result; 

dove

point = Point(lat,lng) 
polygon = Polygon(lat1 lng1, lat2 lng2, lat3 lng3, .... latn lngn, lat1 lng1) 

Si prega di notare che il poligono deve essere chiuso (normalmente è chiuso se si desidera recuperare un dato KML o GoogleMap standard, ma solo assicurarsi che è - nota LAT1 set lng1 è ripetuta alla fine)

non avevo punti e poligoni nel mio database come campi geometrici, così ho dovuto fare qualcosa di simile:

01.235.164,106 mila

Select myWithin(PointFromText(concat("POINT(", latitude, " ", longitude, ")")),PolyFromText('POLYGON((lat1 lng1, ..... latn lngn, lat1 lng1))')) as result 

Spero che questo possa aiutare qualcuno.

Answer 2

vorrei scrivere un'UDF personalizzata che implementa l'algoritmo ray-fusione in C o Delphi o qualsiasi strumento di alto livello si utilizza:

collegamenti per la scrittura di un'UDF
Ecco codice sorgente di un GIS MySQL implementazione che cerca il punto su una sfera (usala come modello per vedere come interagire con i tipi di dati spaziali in MySQL).
http://www.lenzg.net/archives/220-New-UDF-for-MySQL-5.1-provides-GIS-functions-distance_sphere-and-distance_spheroid.html

Dal manuale di MySQL:
http://dev.mysql.com/doc/refman/5.0/en/adding-functions.html

UDF tutorial per MS Visual C++
http://rpbouman.blogspot.com/2007/09/creating-mysql-udfs-with-microsoft.html

UDF tutorial Delphi:
Creating a UDF for MySQL in Delphi

Source-Codice in materia di algoritmo di ray-casting
pseudo-codice: http://rosettacode.org/wiki/Ray-casting_algorithm
articolo drDobbs (notare il link al codice nella parte superiore di questo articolo): http://drdobbs.com/cpp/184409586
Delphi (in realtà FreePascal): http://www.cabiatl.com/mricro/raycast/

Answer 3

Volevo usare la procedura sopra memorizzata in una tabella di poligoni, quindi ecco i comandi per farlo.

Dopo aver importato uno shapefile contenente poligoni in MySQL usando ogr2ogr come segue

ogr2ogr -f "mysql" MYSQL:"mydbname,host=localhost,user=root,password=mypassword,port=3306" -nln "mytablename" -a_srs "EPSG:4326" /path/to/shapefile.shp 

è quindi possibile utilizzare MBRwithin prefiltrare vostro tavolo e mywithin di finire come segue

DROP TEMPORARY TABLE IF EXISTS POSSIBLE_POLYS; 
CREATE TEMPORARY TABLE POSSIBLE_POLYS(OGR_FID INT,SHAPE POLYGON); 
INSERT INTO POSSIBLE_POLYS (OGR_FID, SHAPE) SELECT mytablename.OGR_FID,mytablename.SHAPE FROM mytablename WHERE MBRwithin(@testpoint,mytablename.SHAPE); 

DROP TEMPORARY TABLE IF EXISTS DEFINITE_POLY; 
CREATE TEMPORARY TABLE DEFINITE_POLY(OGR_FID INT,SHAPE POLYGON); 
INSERT INTO DEFINITE_POLY (OGR_FID, SHAPE) SELECT POSSIBLE_POLYS.OGR_FID,POSSIBLE_POLYS.SHAPE FROM POSSIBLE_POLYS WHERE mywithin(@testpoint,POSSIBLE_POLYS.SHAPE); 

in cui viene creato @testpoint , ad esempio, da

SET @longitude=120; 
SET @latitude=-30; 
SET @testpoint =(PointFromText(concat("POINT(", @longitude, " ", @latitude, ")"))); 

Answer 4

Per ogni evenienza, una funzione MySQL che accetta MULTIPOLYGON come ingresso: http://forums.mysql.com/read.php?23,286574,286574

DELIMITER $$ 

CREATE DEFINER=`root`@`localhost` FUNCTION `GISWithin`(pt POINT, mp MULTIPOLYGON) RETURNS int(1) 
    DETERMINISTIC 
BEGIN 

DECLARE str, xy TEXT; 
DECLARE x, y, p1x, p1y, p2x, p2y, m, xinters DECIMAL(16, 13) DEFAULT 0; 
DECLARE counter INT DEFAULT 0; 
DECLARE p, pb, pe INT DEFAULT 0; 

SELECT MBRWithin(pt, mp) INTO p; 

IF p != 1 OR ISNULL(p) THEN 
    RETURN p; 
END IF; 

SELECT X(pt), Y(pt), ASTEXT(mp) INTO x, y, str; 
SET str = REPLACE(str, 'POLYGON((',''); 
SET str = REPLACE(str, '))', ''); 
SET str = CONCAT(str, ','); 

SET pb = 1; 
SET pe = LOCATE(',', str); 
SET xy = SUBSTRING(str, pb, pe - pb); 
SET p = INSTR(xy, ' '); 
SET p1x = SUBSTRING(xy, 1, p - 1); 
SET p1y = SUBSTRING(xy, p + 1); 
SET str = CONCAT(str, xy, ','); 

WHILE pe > 0 DO 
    SET xy = SUBSTRING(str, pb, pe - pb); 
    SET p = INSTR(xy, ' '); 
    SET p2x = SUBSTRING(xy, 1, p - 1); 
    SET p2y = SUBSTRING(xy, p + 1); 

    IF p1y < p2y THEN SET m = p1y; ELSE SET m = p2y; END IF; 

    IF y > m THEN 
     IF p1y > p2y THEN SET m = p1y; ELSE SET m = p2y; END IF; 
     IF y <= m THEN 
      IF p1x > p2x THEN SET m = p1x; ELSE SET m = p2x; END IF; 
      IF x <= m THEN 
       IF p1y != p2y THEN 
        SET xinters = (y - p1y) * (p2x - p1x)/(p2y - p1y) + p1x; 
       END IF; 
       IF p1x = p2x OR x <= xinters THEN 
        SET counter = counter + 1; 
       END IF; 
      END IF; 
     END IF; 
    END IF; 

    SET p1x = p2x; 
    SET p1y = p2y; 
    SET pb = pe + 1; 
    SET pe = LOCATE(',', str, pb); 
END WHILE; 

RETURN counter % 2; 

END 

Answer 5

Ora è un'estensione spaziale al MySQL5.6.1 e sopra. Vedi function_st-contains a Docs.

Answer 6

In risposta alla funzione zarun per trovare lat/long all'interno di un poligono.

Avevo una tabella delle proprietà con informazioni lat/long. Quindi ho dovuto ottenere i record il cui lat/long si trova all'interno di polygon lats/longs (che ho ottenuto dalle API di Google). All'inizio ero stupido su come usare la funzione Zarun. Quindi ecco la query di soluzione per questo.

• Tabella: proprietà
• campi: id, latitudine, longitudine, letti ecc ...
• Query:

SELECT id 
FROM properties 
WHERE myWithin(
    PointFromText(concat("POINT(", latitude, " ", longitude, ")")), 
    PolyFromText('POLYGON((37.628134 -77.458334,37.629867 -77.449021,37.62324 -77.445416,37.622424 -77.457819,37.628134 -77.458334))') 
) = 1 limit 0,50; 

spero che farà risparmiare tempo per scemi come me;)

Answer 7

Ecco una versione che funziona con MULTIPOLYGONS (un adattamento di quello di Zarun che funziona solo per POLYGONS).

CREATE FUNCTION GISWithin(p POINT, multipoly MULTIPOLYGON) RETURNS INT(1) DETERMINISTIC 
BEGIN 
DECLARE n,i,m,x INT DEFAULT 0; 
DECLARE pX,pY,poly1X,poly1Y,poly2X,poly2Y DECIMAL(13,10); 
DECLARE ls LINESTRING; 
DECLARE poly MULTIPOLYGON; 
DECLARE poly1,poly2 POINT; 
DECLARE result INT(1) DEFAULT 0; 
SET pX = X(p); 
SET pY = Y(p); 
SET m = NumGeometries(multipoly); 
WHILE x<m DO 
SET poly = GeometryN(multipoly,x); 
SET ls = ExteriorRing(poly); 
SET poly2 = EndPoint(ls); 
SET poly2X = X(poly2); 
SET poly2Y = Y(poly2); 
SET n = NumPoints(ls); 
WHILE i<n DO 
SET poly1 = PointN(ls, (i+1)); 
SET poly1X = X(poly1); 
SET poly1Y = Y(poly1); 
IF ((((poly1X <= pX) && (pX < poly2X)) || ((poly2X <= pX) && (pX < poly1X))) && (pY > (poly2Y - poly1Y) * (pX - poly1X)/(poly2X - poly1X) + poly1Y)) THEN 
SET result = !result; 
END IF; 
SET poly2X = poly1X; 
SET poly2Y = poly1Y; 
SET i = i + 1; 
END WHILE; 
SET x = x + 1; 
END WHILE; 
RETURN result; 
End;