Quadkey to Lat/Long in SQL

Question

I need to transform quadkey column to latitude and longitude using SQL (Postgres to be more specific). It seems like a common problem, but somehow I cannot find a solution.

Example value:

Quadkey: 12302313230320012

I'd like to get lat/lon of the center of the tile.

Any help will be appreciated.

Answer 1

OK, I have solved the problem. It can be wrapped up into one function, but I am posting two steps in case someone will have slightly different problem:

First we can use this function to convert Quadkey to TileXYZ:

create or replace function quadkey_to_xyz(
    quadkey text,
    out x_tile int,
    out y_tile int,
    out z_tile int) as
$$
declare
    lenQk   int;
    idx     int = 0;
    revQk   text[];
    qkArray text[];
    digit   char;
    mask    int;
BEGIN
    x_tile = 0;
    y_tile = 0;
    z_tile = 0;

    quadkey = TRIM(quadkey);
    lenQk = char_length(quadkey);
    if lenQk = 0 then
        return;
    end if;
    z_tile = lenQk;

    qkArray = regexp_split_to_array(quadkey, '');
    revQk = ARRAY(select qkArray[i] from generate_subscripts(qkArray, 1) as s(i) order by i desc);

    idx = 0;
    foreach digit in array revQk
        loop
            mask = 1 << idx;
            if digit = '1' then
                x_tile = x_tile | mask;
            elseif digit = '2' then
                y_tile = y_tile | mask;
            elseif digit = '3' then
                x_tile = x_tile | mask;
                y_tile = y_tile | mask;
            elseif digit != '0' then
                raise exception 'Unexpected quadkey digit: %', digit;
            end if;
            idx = idx + 1;
        end loop;

    return;
    end;
$$ language plpgsql immutable parallel safe returns null on null input ;

And then we can use these formulas to calculate lat and long

Long = 256*x_tile * 360.0 / (256*pow(2, z_tile))-180;
Lat =  asin((exp((0.5-y_tile/pow(2,z_tile))*4*pi())-1)/(exp((0.5 - y_tile / pow(2, z_tile)) * 4 * pi()) + 1)) * 180 / pi();