Map of Birds Eye Coverage in EMEA & APAC (Part 3)

Step 2: Which tiles are covered by the shape at level 14?

First we use SQL Server 2008 to determine the minimum bounding rectangle, for example:

SELECT Shape.STEnvelope().STAsText() FROM FromBlom WHERE ID=24

In this example STEnvelope() calculated the minimum bounding rectangle and STAsText() provides the result as so called Well Known Text (WKT) which is one of the OGC-specifications we support. Besides that we support also Well Known Binary (WKB) and the Geographic Markup Language (GML).

The WKT for such a polygon looks like this:

POLYGON ((52.150335311889648 0.017314910888671875, 52.309449434280396 0.017314910888671875, 52.309449434280396 0.25343132019042969, 52.150335311889648 0.25343132019042969, 52.150335311889648 0.017314910888671875))

image

Once we have the bounding rectangle we can start to determine the tiles in it:

  1. We determine the locations of the nodes in the bounding rectangle
    SELECT Shape.STEnvelope().STPointN(1).STAsText() FROM FromBlom WHERE ID=24

    STPointN(i) gives me a simple feature of type POINT like ‘POINT (52.150335311889648 0.017314910888671875)’
    for a point of my choice. We will only need point 2 and 4 for our calculations.
    image

  2. Now we determine the latitude and longitude with the functions
    SELECT Shape.STEnvelope().STPointN(1).STX FROM FromBlom WHERE ID=24
    SELECT Shape.STEnvelope().STPointN(1).STY FROM FromBlom WHERE ID=24
    The result is of data type float e.g. 52.1503353118896 and 0.0173149108886719
  3. Well and then we use some mathematical functions which are well described in Joe Schwarz’s article on MSDN to determine the
    quadkey of these tiles. I implemented these functions in the following T-SQL script:

    DECLARE @ID INT;
    DECLARE @i INT;
    DECLARE @Lat FLOAT;
    DECLARE @Lon FLOAT;
    DECLARE @lvl INT;
    DECLARE @pixelX FLOAT;
    DECLARE @pixelY FLOAT;
    DECLARE @tileXdecMin INT;
    DECLARE @tileYdecMin INT;
    DECLARE @tileXdecMax INT;
    DECLARE @tileYdecMax INT;
    DECLARE @tileXdec INT;
    DECLARE @tileYdec INT;
    DECLARE @tileXdecDummy INT;
    DECLARE @tileYdecDummy INT;
    DECLARE @tileXbin VARCHAR(38);
    DECLARE @tileYbin VARCHAR(38);
    DECLARE @j INT;
    DECLARE @quadKey VARCHAR(19);
    
    SET @lvl = 14;
    
    DECLARE IdCursor CURSOR FOR SELECT ID FROM FromBlom WHERE Country='IE' ORDER BY ID;
    OPEN IdCursor;
    FETCH NEXT FROM IdCursor INTO @ID;
    WHILE @@FETCH_STATUS = 0
        BEGIN
            SET @i = 2;
            WHILE @i <= 4 BEGIN
                -- Get latitude and longitude of the bounding reactangle
                SET @Lat = (SELECT Shape.STEnvelope().STPointN(@i).STX FROM FromBlom WHERE ID=@ID);
                SET @Lon = (SELECT Shape.STEnvelope().STPointN(@i).STY FROM FromBlom WHERE ID=@ID);
                -- Calculate pixel X and Y of the CP
                SET @pixelX = (@Lon + 180) / 360 * 256 * POWER(2, @lvl);
                SET @pixelY = (0.5 - (LOG((1 + SIN(@Lat * PI() / 180)) / (1 - SIN(@Lat * PI() / 180))) 
    / 4 / PI())) * 256 * POWER(2, @lvl); -- Calculate decimal tile X and Y IF (@i=2) BEGIN SET @tileXdecMin = FLOOR(@pixelX / 256); SET @tileYdecMin = FLOOR(@pixelY / 256); END; ELSE BEGIN SET @tileXdecMax = FLOOR(@pixelX / 256); SET @tileYdecMax = FLOOR(@pixelY / 256); END; SET @i = @i+2; END; -- Now we loop through the tiles in this rectangle SET @tileYdec = @tileYdecMin; WHILE @tileYdec <= @tileYdecMax BEGIN SET @tileXdec = @tileXdecMin; WHILE @tileXdec <= @tileXdecMax BEGIN SET @tileYdecDummy = @tileYdec; SET @tileXdecDummy = @tileXdec; -- Calculate binary tile X and Y SET @tileXbin = ''; SET @tileYbin = ''; SET @j = @lvl-1; WHILE @j >= 0 BEGIN SET @tileXdecDummy = @tileXdecDummy - POWER(2, @j); IF (@tileXdecDummy < 0) BEGIN SET @tileXbin = @tileXbin + '0' SET @tileXdecDummy = @tileXdecDummy + POWER(2, @j); END ELSE BEGIN SET @tileXbin = @tileXbin + '1' END SET @tileYdecDummy = @tileYdecDummy - POWER(2, @j); IF (@tileYdecDummy < 0) BEGIN SET @tileYbin = @tileYbin + '0' SET @tileYdecDummy = @tileYdecDummy + POWER(2, @j); END ELSE BEGIN SET @tileYbin = @tileYbin + '1' END SET @j = @j-1; END; -- Calculate Quadkey SET @j = 1; SET @quadKey = ''; WHILE @j <= @lvl BEGIN SET @quadKey = @quadKey + CAST(2 * CAST(SUBSTRING(@tileYbin, @j, 1) AS INT) + 1 *
    CAST(SUBSTRING(@tileXbin, @j, 1) AS INT) AS VARCHAR(1)); SET @j = @j+1; END; INSERT INTO [VE-Tiles] (ID, LVL, Quadkey) VALUES (@ID, @lvl, @quadKey) SET @tileXdec = @tileXdec+1; END; SET @tileYdec = @tileYdec+1; END; PRINT ''; FETCH NEXT FROM IdCursor INTO @ID; END; CLOSE IdCursor; DEALLOCATE IdCursor; image

    If you had a closer look at my script you noticed that I inserted the quadkeys in a new table and in the next step I will create a
    shape which represents the bounding box of each tile:

    DECLARE @myQuadkey VARCHAR(19);
    DECLARE @myDummy VARCHAR(38);
    DECLARE @myLevel INT;
    DECLARE @myQuadKeyBin VARCHAR(38);
    DECLARE @i INT;
    DECLARE @TileXBin VARCHAR(19);
    DECLARE @TileYBin VARCHAR(19); 
    DECLARE @TileX INT;
    DECLARE @TileY INT;
    DECLARE @PixelXMin INT;
    DECLARE @PixelXMax INT;
    DECLARE @PixelYMin INT;
    DECLARE @PixelYMax INT;
    DECLARE @LatMin FLOAT;
    DECLARE @LatMax FLOAT;
    DECLARE @LongMin FLOAT;
    DECLARE @LongMax FLOAT;
    DECLARE Tile_Cursor CURSOR FOR SELECT Quadkey FROM [VE-Tiles] WITH (NOLOCK) WHERE BBox IS NULL;
    OPEN Tile_Cursor;
    FETCH NEXT FROM Tile_Cursor INTO @myQuadkey;
    WHILE @@FETCH_STATUS = 0
        BEGIN
    -------------------------------------------------------------------
    --        Convert Quadkey in binary Quadkey
    -------------------------------------------------------------------
            SET @myLevel = LEN(@myQuadkey);
            SET @myQuadKeyBin = '';
            SET @i = 0;
            WHILE @i < @myLevel BEGIN
                SET @myDummy = SUBSTRING(@myQuadkey, @myLevel-@i, 1);
                IF @myDummy = 3 SET @myQuadKeyBin = '11' + @myQuadKeyBin;
                IF @myDummy = 2 SET @myQuadKeyBin = '10' + @myQuadKeyBin;
                IF @myDummy = 1 SET @myQuadKeyBin = '01' + @myQuadKeyBin;
                IF @myDummy = 0 SET @myQuadKeyBin = '00' + @myQuadKeyBin;
                SET @i = @i+1
            END
    -------------------------------------------------------------------
    --        Break binary Quadtree in binary tileX and tileY
    -------------------------------------------------------------------
            SET @myDummy = @myQuadKeyBin;
            SET @TileXBin = '';
            SET @TileYBin = '';
            SET @i = 0;
            WHILE @i < LEN(@myQuadKeyBin) BEGIN
                SET @myDummy = SUBSTRING(@myQuadKeyBin, LEN(@myQuadKeyBin)-@i, 1)
                IF (CAST(@i AS FLOAT)/2 - FLOOR(@i/2) = 0) BEGIN
                    SET @TileXBin = @myDummy + @TileXBin;
                END
                ELSE BEGIN
                    SET @TileYBin = @myDummy + @TileYBin;
                END
                SET @i = @i+1
            END
    -------------------------------------------------------------------
    --        Convert binary tileX and tileY to decimal tileX and tileY
    -------------------------------------------------------------------
            SET @TileX = '';
            SET @TileY = '';
            SET @i = 0;
            WHILE @i < @myLevel BEGIN
                SET @myDummy = SUBSTRING(@TileXBin, LEN(@TileXBin)-@i, 1)
                SET @TileX = @TileX + @myDummy * POWER(2, @i)
                SET @i = @i+1
            END            
            SET @i = 0;
            WHILE @i < @myLevel BEGIN
                SET @myDummy = SUBSTRING(@TileYBin, LEN(@TileYBin)-@i, 1)
                SET @TileY = @TileY + @myDummy * POWER(2, @i)
                SET @i = @i+1
            END
    -------------------------------------------------------------------
    --        Determine PixelX and PixelY-Ranges
    -------------------------------------------------------------------
            SET @PixelXMin = @TileX * 256
            SET @PixelXMax = ((@TileX + 1) * 256) - 1
            SET @PixelYMin = @TileY * 256
            SET @PixelYMax = ((@TileY + 1) * 256) - 1
    -------------------------------------------------------------------
    --        Calculate Bounding-Box
    -------------------------------------------------------------------
            SET @LongMin = (CAST(@PixelXMin AS FLOAT) / 256 / POWER(2, CAST(@myLevel AS FLOAT)) * 360) - 180
            SET @LongMax = (CAST(@PixelXMax AS FLOAT) / 256 / POWER(2, CAST(@myLevel AS FLOAT)) * 360) - 180
            SET @LatMin = ASIN((EXP((0.5 - CAST(@PixelYMin AS FLOAT) / 256 / POWER(2, CAST(@myLevel AS FLOAT))) 
    * 4 * PI()) - 1) / (EXP((0.5 - CAST(@PixelYMin AS FLOAT) / 256
    / POWER(2, CAST(@myLevel AS FLOAT))) * 4 * PI()) + 1)) * 180 / PI() SET @LatMax = ASIN((EXP((0.5 - CAST(@PixelYMax AS FLOAT) / 256 / POWER(2, CAST(@myLevel AS FLOAT)))
    * 4 * PI()) - 1) / (EXP((0.5 - CAST(@PixelYMax AS FLOAT) / 256
    / POWER(2, CAST(@myLevel AS FLOAT))) * 4 * PI()) + 1)) * 180 / PI() ------------------------------------------------------------------- -- Update Geometry ------------------------------------------------------------------- EXEC('UPDATE [VE-Tiles] SET BBox=geometry::STGeomFromText(' + '''' + 'POLYGON((' + @LatMin + ' '
    + @LongMin + ', ' + @LatMin + ' ' + @LongMax + ', ' + @LatMax + ' ' + @LongMax + ', '
    + @LatMax + ' ' + @LongMin + ', ' + @LatMin + ' ' + @LongMin + '))' + ''''
    + ', 4326) WHERE QuadKey=' + '''' + @myQuadkey + ''''); FETCH NEXT FROM Tile_Cursor INTO @myQuadkey; END; CLOSE Tile_Cursor; DEALLOCATE Tile_Cursor;

    With a call like the one shown below we can retrieve the result for one shape
    SELECT BBox.STAsText() FROM [VE-Tiles] WHERE ID=24
    and here is how it looks like in Virtual Earth:
    image 

Advertisements
This entry was posted in SQL Server Spatial. Bookmark the permalink.

Leave a Reply

Fill in your details below or click an icon to log in:

WordPress.com Logo

You are commenting using your WordPress.com account. Log Out / Change )

Twitter picture

You are commenting using your Twitter account. Log Out / Change )

Facebook photo

You are commenting using your Facebook account. Log Out / Change )

Google+ photo

You are commenting using your Google+ account. Log Out / Change )

Connecting to %s