将纬度/经度点转换为麦卡托投影法上的像素(x，y)

我正在尝试将一个 lat/long 点转换成2 d 点，这样我就可以将它显示在世界图像上——这是一个麦卡托投影法。

我已经看到了各种各样的方法来实现这一点，并且还有一些关于堆栈溢出的问题——我已经尝试了不同的代码片段，尽管我得到了正确的经度到像素，但是纬度总是关闭的——看起来似乎变得更加合理了。

我需要的公式，以考虑到图像的大小，宽度等。

我试过这段代码:

double minLat = -85.05112878;
double minLong = -180;
double maxLat = 85.05112878;
double maxLong = 180;


// Map image size (in points)
double mapHeight = 768.0;
double mapWidth = 991.0;


// Determine the map scale (points per degree)
double xScale = mapWidth/ (maxLong - minLong);
double yScale = mapHeight / (maxLat - minLat);


// position of map image for point
double x = (lon - minLong) * xScale;
double y = - (lat + minLat) * yScale;


System.out.println("final coords: " + x + " " + y);

在我尝试的例子中，纬度似乎偏离了大约30像素。有什么帮助或建议吗？

更新

基于这个问题: 从 Lat/lon 到 xy

我已经尝试使用提供的代码，但我仍然有一些问题与纬度转换，经度是罚款。

int mapWidth = 991;
int mapHeight = 768;


double mapLonLeft = -180;
double mapLonRight = 180;
double mapLonDelta = mapLonRight - mapLonLeft;


double mapLatBottom = -85.05112878;
double mapLatBottomDegree = mapLatBottom * Math.PI / 180;
double worldMapWidth = ((mapWidth / mapLonDelta) * 360) / (2 * Math.PI);
double mapOffsetY = (worldMapWidth / 2 * Math.log((1 + Math.sin(mapLatBottomDegree)) / (1 - Math.sin(mapLatBottomDegree))));


double x = (lon - mapLonLeft) * (mapWidth / mapLonDelta);
double y = 0.1;
if (lat < 0) {
lat = lat * Math.PI / 180;
y = mapHeight - ((worldMapWidth / 2 * Math.log((1 + Math.sin(lat)) / (1 - Math.sin(lat)))) - mapOffsetY);
} else if (lat > 0) {
lat = lat * Math.PI / 180;
lat = lat * -1;
y = mapHeight - ((worldMapWidth / 2 * Math.log((1 + Math.sin(lat)) / (1 - Math.sin(lat)))) - mapOffsetY);
System.out.println("y before minus: " + y);
y = mapHeight - y;
} else {
y = mapHeight / 2;
}
System.out.println(x);
System.out.println(y);

当使用原始代码时，如果纬度值是正值，那么它返回一个负值，所以我稍微修改了一下它，并使用极端纬度进行了测试——应该是0点和766点，它工作得很好。然而，当我尝试使用不同的纬度值 ex: 58.07(就在英国北部)时，它显示为西班牙北部。

107728 次浏览
小开

你不能仅仅像那样把经度/纬度转换成 x/y，因为世界不是平的。你看过这个帖子了吗？将经度/纬度转换为 X/Y 坐标

更新-1/18/13

我决定尝试一下，我是这样做的:-

public class MapService {
// CHANGE THIS: the output path of the image to be created
private static final String IMAGE_FILE_PATH = "/some/user/path/map.png";


// CHANGE THIS: image width in pixel
private static final int IMAGE_WIDTH_IN_PX = 300;


// CHANGE THIS: image height in pixel
private static final int IMAGE_HEIGHT_IN_PX = 500;


// CHANGE THIS: minimum padding in pixel
private static final int MINIMUM_IMAGE_PADDING_IN_PX = 50;


// formula for quarter PI
private final static double QUARTERPI = Math.PI / 4.0;


// some service that provides the county boundaries data in longitude and latitude
private CountyService countyService;


public void run() throws Exception {
// configuring the buffered image and graphics to draw the map
BufferedImage bufferedImage = new BufferedImage(IMAGE_WIDTH_IN_PX,
IMAGE_HEIGHT_IN_PX,
BufferedImage.TYPE_INT_RGB);


Graphics2D g = bufferedImage.createGraphics();
Map<RenderingHints.Key, Object> map = new HashMap<RenderingHints.Key, Object>();
map.put(RenderingHints.KEY_INTERPOLATION, RenderingHints.VALUE_INTERPOLATION_BICUBIC);
map.put(RenderingHints.KEY_RENDERING, RenderingHints.VALUE_RENDER_QUALITY);
map.put(RenderingHints.KEY_ANTIALIASING, RenderingHints.VALUE_ANTIALIAS_ON);
RenderingHints renderHints = new RenderingHints(map);
g.setRenderingHints(renderHints);


// min and max coordinates, used in the computation below
Point2D.Double minXY = new Point2D.Double(-1, -1);
Point2D.Double maxXY = new Point2D.Double(-1, -1);


// a list of counties where each county contains a list of coordinates that form the county boundary
Collection<Collection<Point2D.Double>> countyBoundaries = new ArrayList<Collection<Point2D.Double>>();


// for every county, convert the longitude/latitude to X/Y using Mercator projection formula
for (County county : countyService.getAllCounties()) {
Collection<Point2D.Double> lonLat = new ArrayList<Point2D.Double>();


for (CountyBoundary countyBoundary : county.getCountyBoundaries()) {
// convert to radian
double longitude = countyBoundary.getLongitude() * Math.PI / 180;
double latitude = countyBoundary.getLatitude() * Math.PI / 180;


Point2D.Double xy = new Point2D.Double();
xy.x = longitude;
xy.y = Math.log(Math.tan(QUARTERPI + 0.5 * latitude));


// The reason we need to determine the min X and Y values is because in order to draw the map,
// we need to offset the position so that there will be no negative X and Y values
minXY.x = (minXY.x == -1) ? xy.x : Math.min(minXY.x, xy.x);
minXY.y = (minXY.y == -1) ? xy.y : Math.min(minXY.y, xy.y);


lonLat.add(xy);
}


countyBoundaries.add(lonLat);
}


// readjust coordinate to ensure there are no negative values
for (Collection<Point2D.Double> points : countyBoundaries) {
for (Point2D.Double point : points) {
point.x = point.x - minXY.x;
point.y = point.y - minXY.y;


// now, we need to keep track the max X and Y values
maxXY.x = (maxXY.x == -1) ? point.x : Math.max(maxXY.x, point.x);
maxXY.y = (maxXY.y == -1) ? point.y : Math.max(maxXY.y, point.y);
}
}


int paddingBothSides = MINIMUM_IMAGE_PADDING_IN_PX * 2;


// the actual drawing space for the map on the image
int mapWidth = IMAGE_WIDTH_IN_PX - paddingBothSides;
int mapHeight = IMAGE_HEIGHT_IN_PX - paddingBothSides;


// determine the width and height ratio because we need to magnify the map to fit into the given image dimension
double mapWidthRatio = mapWidth / maxXY.x;
double mapHeightRatio = mapHeight / maxXY.y;


// using different ratios for width and height will cause the map to be stretched. So, we have to determine
// the global ratio that will perfectly fit into the given image dimension
double globalRatio = Math.min(mapWidthRatio, mapHeightRatio);


// now we need to readjust the padding to ensure the map is always drawn on the center of the given image dimension
double heightPadding = (IMAGE_HEIGHT_IN_PX - (globalRatio * maxXY.y)) / 2;
double widthPadding = (IMAGE_WIDTH_IN_PX - (globalRatio * maxXY.x)) / 2;


// for each country, draw the boundary using polygon
for (Collection<Point2D.Double> points : countyBoundaries) {
Polygon polygon = new Polygon();


for (Point2D.Double point : points) {
int adjustedX = (int) (widthPadding + (point.getX() * globalRatio));


// need to invert the Y since 0,0 starts at top left
int adjustedY = (int) (IMAGE_HEIGHT_IN_PX - heightPadding - (point.getY() * globalRatio));


polygon.addPoint(adjustedX, adjustedY);
}


g.drawPolygon(polygon);
}


// create the image file
ImageIO.write(bufferedImage, "PNG", new File(IMAGE_FILE_PATH));
}
}

结果: 图像宽度 = 600px，图像高度 = 600px，图像填充 = 50px

enter image description here

结果: 图像宽度 = 300px，图像高度 = 500px，图像填充 = 50px

enter image description here

小开

墨卡托映射投影是朗伯圆锥共形映射投影的一种特殊极限情形 赤道作为单一的标准平行线。所有其他纬度的平行线都是直线和子午线 也是与赤道成直角的直线，间距相等 投影的倾斜形式。它很少用于土地制图的目的，但几乎普遍用于 导航图。除了保形外，它还有一个特殊的性质，就是在它上面画的直线 因此，导航员可以从直线的角度推导出他们的航向 与经脉相连

Mercator projection

由球面纬度 φ 和经度 λ 导出投影东北坐标的公式 是:

E = FE + R (λ – λₒ)
N = FN + R ln[tan(π/4 + φ/2)]

其中 λ哦是自然起源的经度，FE 和 FN 是假东移和假北移。 在球面墨卡托实际上没有使用这些值，所以您可以将公式简化为

derivation of the mercator projection (wikipedia)

伪代码的例子，所以这可以适用于每一种编程语言。

latitude    = 41.145556; // (φ)
longitude   = -73.995;   // (λ)


mapWidth    = 200;
mapHeight   = 100;


// get x value
x = (longitude+180)*(mapWidth/360)


// convert from degrees to radians
latRad = latitude*PI/180;


// get y value
mercN = ln(tan((PI/4)+(latRad/2)));
y     = (mapHeight/2)-(mapWidth*mercN/(2*PI));

资料来源:

  1. OGP 地理信息委员会，指导说明7，第2部分: 坐标转换和转换
  2. 麦卡托投影法的推导
  3. 国家地图集: 地图预测
  4. 墨卡托地图投影

剪辑 用 PHP 创建了一个工作示例(因为我的 Java 很烂)

Https://github.com/mfeldheim/mapstuff.git

编辑2

麦卡托投影法的动画不错 Https://amp-reddit-com.cdn.ampproject.org/v/s/amp.reddit.com/r/educationalgifs/comments/5lhk8y/how_the_mercator_projection_distorts_the_poles/?usqp=mq331aqjcaeoavgbgaeb&amp_js_v=0.1

小开

Java 版本的原始 谷歌地图 JavaScript API v3 Java 脚本代码如下，它工作起来没有问题

public final class GoogleMapsProjection2
{
private final int TILE_SIZE = 256;
private PointF _pixelOrigin;
private double _pixelsPerLonDegree;
private double _pixelsPerLonRadian;


public GoogleMapsProjection2()
{
this._pixelOrigin = new PointF(TILE_SIZE / 2.0,TILE_SIZE / 2.0);
this._pixelsPerLonDegree = TILE_SIZE / 360.0;
this._pixelsPerLonRadian = TILE_SIZE / (2 * Math.PI);
}


double bound(double val, double valMin, double valMax)
{
double res;
res = Math.max(val, valMin);
res = Math.min(res, valMax);
return res;
}


double degreesToRadians(double deg)
{
return deg * (Math.PI / 180);
}


double radiansToDegrees(double rad)
{
return rad / (Math.PI / 180);
}


PointF fromLatLngToPoint(double lat, double lng, int zoom)
{
PointF point = new PointF(0, 0);


point.x = _pixelOrigin.x + lng * _pixelsPerLonDegree;


// Truncating to 0.9999 effectively limits latitude to 89.189. This is
// about a third of a tile past the edge of the world tile.
double siny = bound(Math.sin(degreesToRadians(lat)), -0.9999,0.9999);
point.y = _pixelOrigin.y + 0.5 * Math.log((1 + siny) / (1 - siny)) *- _pixelsPerLonRadian;


int numTiles = 1 << zoom;
point.x = point.x * numTiles;
point.y = point.y * numTiles;
return point;
}


PointF fromPointToLatLng(PointF point, int zoom)
{
int numTiles = 1 << zoom;
point.x = point.x / numTiles;
point.y = point.y / numTiles;


double lng = (point.x - _pixelOrigin.x) / _pixelsPerLonDegree;
double latRadians = (point.y - _pixelOrigin.y) / - _pixelsPerLonRadian;
double lat = radiansToDegrees(2 * Math.atan(Math.exp(latRadians)) - Math.PI / 2);
return new PointF(lat, lng);
}


public static void main(String []args)
{
GoogleMapsProjection2 gmap2 = new GoogleMapsProjection2();


PointF point1 = gmap2.fromLatLngToPoint(41.850033, -87.6500523, 15);
System.out.println(point1.x+"   "+point1.y);
PointF point2 = gmap2.fromPointToLatLng(point1,15);
System.out.println(point2.x+"   "+point2.y);
}
}


public final class PointF
{
public double x;
public double y;


public PointF(double x, double y)
{
this.x = x;
this.y = y;
}
}
小开 
 public static String getTileNumber(final double lat, final double lon, final int zoom) {
int xtile = (int)Math.floor( (lon + 180) / 360 * (1<<zoom) ) ;
int ytile = (int)Math.floor( (1 - Math.log(Math.tan(Math.toRadians(lat)) + 1 /  Math.cos(Math.toRadians(lat))) / Math.PI) / 2 * (1<<zoom) ) ;
if (xtile < 0)
xtile=0;
if (xtile >= (1<<zoom))
xtile=((1<<zoom)-1);
if (ytile < 0)
ytile=0;
if (ytile >= (1<<zoom))
ytile=((1<<zoom)-1);
return("" + zoom + "/" + xtile + "/" + ytile);
}
}
小开

只有 JAVA？

这里是 Python 代码! 请参考 将纬度/经度点转换为麦卡托投影法上的像素(x，y)

import math
from numpy import log as ln


# Define the size of map
mapWidth    = 200
mapHeight   = 100




def convert(latitude, longitude):
# get x value
x = (longitude + 180) * (mapWidth / 360)


# convert from degrees to radians
latRad = (latitude * math.pi) / 180


# get y value
mercN = ln(math.tan((math.pi / 4) + (latRad / 2)))
y     = (mapHeight / 2) - (mapWidth * mercN / (2 * math.pi))
    

return x, y


print(convert(41.145556, 121.2322))

回答:

(167.35122222222225, 24.877939817552335)
小开

我是新来的，只是来写作，因为我已经关注这个社区好几年了。我很高兴能有所贡献。

我花了几乎一天的时间来寻找你的问题鼓励我继续寻找。

我得到了下面的函数，它可以工作! 本文学分: https://towardsdatascience.com/geotiff-coordinate-querying-with-javascript-5e6caaaf88cf

var bbox = [minLong, minLat, maxLong, maxLat];
var pixelWidth = mapWidth;
var pixelHeight = mapHeight;
var bboxWidth = bbox[2] - bbox[0];
var bboxHeight = bbox[3] - bbox[1];


var convertToXY = function(latitude, longitude) {
var widthPct = ( longitude - bbox[0] ) / bboxWidth;
var heightPct = ( latitude - bbox[1] ) / bboxHeight;
var x = Math.floor( pixelWidth * widthPct );
var y = Math.floor( pixelHeight * ( 1 - heightPct ) );
return { x, y };
}

提交答案