Saturday, September 11, 2010

mysql function, well, procedure

Well this is fantastic. I have been looking for a routine that determines whether a point is in a polygon and i stumbled upon this at: http://forums.mysql.com/read.php?23,201913,201913

I have quoted it at the bottom of this post. it refers to a c function (something i am not familiar with) that performs the "contains" task but the author (Adam Smith) has converted it to a MySql procedure and, as far as I can tell from my experiments, its works! I also tested it on a donut polygon setup and it works a treat there too. As I am going to have to iterate 1000's of records and update them accordingly this is a great bit of code to add to my library.

DELIMITER //
CREATE FUNCTION minimum(x DECIMAL(18,15),y DECIMAL(18,15)) RETURNS DECIMAL(18,15) DETERMINISTIC
BEGIN
DECLARE result DECIMAL(18,15);
IF x<y THEN
  SET result = x;
ELSE
  SET result = y;
END IF;
RETURN result;
END;
//

CREATE FUNCTION maximum(x DECIMAL(18,15),y DECIMAL(18,15)) RETURNS DECIMAL(18,15) DETERMINISTIC
BEGIN
DECLARE result DECIMAL(18,15);
IF x>y THEN
  SET result = x;
ELSE
  SET result = y;
END IF;
RETURN result;
END;
//

CREATE FUNCTION point_inside_polygon(x DECIMAL(18,15),y DECIMAL(18,15),p TEXT) RETURNS INT(1) DETERMINISTIC
BEGIN
DECLARE counter INT DEFAULT 0;
DECLARE result INT(1) DEFAULT 0;
DECLARE n INT DEFAULT 0;
DECLARE str TEXT;
DECLARE str2 TEXT;
DECLARE pos INT;
DECLARE coords VARCHAR(50);
DECLARE coords_pos INT;
DECLARE px DECIMAL(18,15);
DECLARE py DECIMAL(18,15);
DECLARE p1x DECIMAL(18,15);
DECLARE p1y DECIMAL(18,15);
DECLARE p2x DECIMAL(18,15);
DECLARE p2y DECIMAL(18,15);
DECLARE m DECIMAL(18,15);
DECLARE i INT;
DECLARE j INT;
DECLARE modulus INT;
DECLARE xinters DECIMAL(18,15);
SET str = REPLACE(REPLACE(p,'POLYGON((',''),'))','');
SET str2 = CONCAT(str,','); 
SET pos = INSTR(str,',');
SET coords = SUBSTRING(str,1,pos-1);
SET coords_pos = INSTR(coords,' ');
SET p1x = SUBSTRING(coords,1,coords_pos-1);
SET p1y = SUBSTRING(coords,coords_pos+1);
WHILE pos>0 DO
  SET n = n + 1;
  SET str = SUBSTRING(str,pos+1);
  SET pos = INSTR(str,',');
END WHILE;
SET i = 0;
SET n = n + 1;
WHILE i<=n DO
  SET modulus = i % n;
  SET j = 0;
  SET str = str2;
  WHILE j<=modulus DO
    SET j = j + 1;
    SET pos = INSTR(str,',');
    SET coords = SUBSTRING(str,1,pos-1);
    SET coords_pos = INSTR(coords,' ');
    SET px = SUBSTRING(coords,1,coords_pos-1);
    SET py = SUBSTRING(coords,coords_pos+1);
    SET str = SUBSTRING(str,pos+1);
  END WHILE;
  SET p2x = px;
  SET p2y = py;
  SELECT minimum(p1y,p2y) INTO m;
  IF y>m THEN
    SELECT maximum(p1y,p2y) INTO m;
    IF y<=m THEN
      SELECT maximum(p1x,p2x) INTO m;
      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;
END WHILE;
IF counter%2!=0 THEN
  SET result = 1;
END IF;
RETURN result;
END;
//

DELIMITER ;

-- usage example:
CREATE TABLE `polygons` ( 
`id` int(11) NOT NULL auto_increment, 
`polygon_data` polygon NOT NULL, 
PRIMARY KEY (`id`), 
SPATIAL KEY `polygon_data` (`polygon_data`(32)) 
) ENGINE=MyISAM DEFAULT CHARSET=latin1;

SET @x = -79.12345678901234;
SET @y = 43.12345678901234;
SET @point = CONCAT('POINT(',@x,' ',@y,')');
SELECT id FROM polygons WHERE MBRContains(polygon_data,GeomFromText(@point)) AND point_inside_polygon(@x,@y,ASTEXT(polygon_data));

No comments:

Post a Comment