Google map tile chart algorithm analysis

 

From an article I saw last yearArticleHttp://www.earthplayer.com/bbs/viewthread.php? Tid = 321 is no longer accessible, and the address of the tile server is changed. Update it for backup.

Google map tile chartAlgorithmAnalysis
Map tile Encoding

Google map uses two algorithms to encode the location of tile.

For Google Ditu, the tile URL is similar to: http://mt2.google.cn/mt? V = cn1.6 & HL = ZH-CN & X = 6678 & Y = 3557 & Z = 13 & s = galile use X and Y to set tile coordinates and amplification factors. The amplification factor ranges from 17 (completely reduced) to 0 (maximum ratio ). When the amplification factor is 17, the whole earth is displayed in a tile. At this time, x = 0, y = 0; when the amplification factor is 16, the Earth is divided into 2x2 parts, then 0 <= x <= 1 and 0 <= Y <= 1. Each time you zoom in, each tile is divided into four parts. Therefore, when the amplification factor is Z, the number of horizontal and vertical tile displayed is 2 ^ (17-z ).

The theoretical algorithm will look like the following:

// Correct the latitude to go from 0 (North) to 180 (south ),
// Instead of 90 (North) to-90 (South)
Latitude = 90-latitude;

// Correct the longpolling to go from 0 to 360
Longpolling = 180 + longpolling;

// Find tile size from zoom level
Double lattilesize = 180/(POW (2, (17-zoom )));
Double longtilesize = 360/(POW (2, (17-zoom )));

// Find the tile coordinates
Int tilex = (INT) (longpolling/longtilesize );
Int tiley = (INT) (latitude/lattilesize );

The preceding figure assumes that the projection is a plane.

In fact, GoogleMap and 51ditu use the mocato projection during display. Therefore, the above algorithm needs to be modified.
In mocato projection, the distance between the two weft is not necessarily equal, so Titley corresponding to the latitude is based on its vertical position.

BelowCodeThe number of verticals of tile calculated from the latitude position of tile.

Collapse
/** <Summary> get the vertical tile number from a latitude using Mercator ptrojection formula </Summary> */

/** <Summary> This function obtains the sequence number of the vertical axis of the tile corresponding to the corresponding latitude. The parameter is latitude. </Summary> */

Private int getmercatorlatitude (double lati)
{
Double maxlat = math. Pi;

Double lat = lati;

If (LAT> 90) lat = lat-180;
If (LAT <-90) lat = lat + 180;

// Conversion DEGRE => radians
// Converts degrees to radians
Double Phi = math. Pi * LAT/180;

double res;
// other posts on the Internet. This is the plus sign.
// double temp = math. tan (math. PI/4 + PHI/2);
// res = math. log (temp);
// The following statement is the preceding merge
res = 0.5 * Math. log (1 + math. sin (PHI)/(1-math. sin (PHI);
double maxtiley = math. pow (2, zoom);
int result = (INT) (1-Res/maxlat)/2) * (maxtiley);

Return (result );
}

 

Coverage area:
In theory, the latitude range should be-90 degrees to 90 degrees, but in fact, because the mokto projection makes two levels of infinity, the coverage area is less than-90 to 90.

 

Add the python verification code.

 

Code
# -*-Coding: UTF-8 -*-


Import Math

Tileserver =   " Http://mt2.google.cn/mt? V = cn1.6 & HL = ZH-CN "

Def Getmercatorlatitude (lati, zoom ):

Maxlat = Math. Pi
Lat = Lati

If (Lat >   90 ):
Lat = Lat -   180
If (Lat <   - 90 ):
Lat = Lat +   180

# Converts degrees to radians.
Phi = Math. Pi * Lat /   180.0 ;

# Temp = math. Tan (math. PI/4.0 + PHI/2.0 );
# Res = math. Log (temp );
# The following sentence is the preceding merge.
Res =   0.5   * Math. Log (( 1   + Math. Sin (PHI )) / ( 1   - Math. Sin (PHI )))
Maxtiley = Math. Pow ( 2 , 17 - Zoom)
Result = (INT )((( 1   - Res / Maxlat) /   2 ) * (Maxtiley ))

Return Result;


Def Gettile (longpolling, latitude, zoom ):

Longpolling = 180 + Longpolling

Longtilesize = 360.0 / (POW ( 2 ,( 17 - Zoom )))

Tilex = Longpolling / Longtilesize

Tiley = Getmercatorlatitude (latitude, zoom)

Tilex = INT (math. Floor (tilex ))

Tiley = INT (math. Floor (tiley ))

# Return [tilex, tiley]

Return Tileserver +   " & Amp; X = "   + STR (tilex) +   " & Y = "   + STR (tiley) +   " & Z = "   + STR ( 17 - Zoom) +   " & S = galile "

Print Gettile ( 114.28 , 30.55 , 5 )

# Test result http://mt2.google.cn/mt? V = cn1.6 & HL = ZH-CN & X = 3348 & Y = 1682 & Z = 12 & s = galile