本文转自 http://www.cnblogs.com/Penrod/archive/2009/12/16/1625658.html

1，假定地图块为256X256像素大小。从第一级整个世界一个地图块，按照比例尺倍增的方式放大。第N级的地图比例尺计算方法如下：

　　　　1 : (20037508.3427892*2*100)/(256/96*2.54) /(2^(N-1))= 1: 591658710.90912992125984251968504/(2^(N-1)) 

2，转换成地图图片(Map Tile)，第N级的总像素数计算方法如下:

 　　256*(2^(N-1)) 

3，原始数据(经纬度)按墨卡托投影(把地球视为正球体)。投影文件如下：

PROJCS["Google_Mercator",GEOGCS["GCS_WGS_1984",DATUM["D_WGS_1984",SPHEROID["WGS_1984",6378137.0,0.0]],

PRIMEM["Greenwich",0.0],UNIT["Degree",0.0174532925199433]],PROJECTION["Mercator"],PARAMETER["False_Easting",0.0],

PARAMETER["False_Northing",0.0],PARAMETER["Central_Meridian",0.0],PARAMETER["Standard_Parallel_1",0.0],UNIT["Meter",1.0]]

4，投影后世界地图东西向范围(-20037508.3427892,20037508.3427892)(米),南北向忽略变形大的高纬度地区取东西向相同范围。

经纬度范围东西[-180，180]，南北[-85.05112877980659，85.05112877980659]。 

转换成地图图片后，第17级像素坐标和经纬度的对应关系可以达到精确度，本文后续的像素坐标指在该级别下的坐标。

经纬度(dLng,dLat)和第17级像素坐标(X,Y)转换。

dLng = (X * 180.0/16777216-180.0) ;

dLat = (atan(exp((Y-16777216)/(-5340353.8065374867)))*2-(PI/2))/PI*180;  

 X = (int)((dLng + 180)*16777216/180 + 0.5);

 Y = (int)((log(tan((dLat/180*PI+(PI/2))/2)) *(-5340353.8065374867)) +16777216 + 0.5); 

第N级的像素坐标(nX,nY)和第17级像素坐标转换(X,Y)转换。

nX = X / (2^(17-N));

nY = Y / (2^(17-N)); 

X = nX*(2^(17-N)); 

Y = nY*(2^(17-N)); 

5，地图块编号规则，按照从"左"到"右"递增列号，从"上"到"下"递增行号，列号和行号从0开始。

6，第N级像素坐标(nX,nY)和该点所在列号行号(col,row)以及该点在该块上的偏移坐标(dX,dY)的转换。

nX  = col * 256 + dX;

nY  = row * 256 + dY;

col = nX / 256;

dX = nX % 256;

row = nY / 256;

dY = nY % 256;