背景
现在我有一个提供mvt矢量切片的服务、一个基于mapboxjs的gis展示页面。
想要做一个功能：在浏览地图时，用户鼠标点击地图上某个位置后，发起一个查询请求，该请求需要在地图服务端查询是否命中某些个feature？如果命中，则返回命中的位于最上层的一个feature，否则返回空。
思路
- 最开始的思路是在请求中传入 x,y经纬度参数，服务端接收后作为一个point去判断这个point与某些feature的包含关系。
结果：经过测试，不可取。
原因：这样忽略了用户当前所处 zoom（地图缩放等级） 和 硬件 pixel（像素值）的因素，会导致你传过去的这个point始终都是在 最大zoom比如 22下的一个判断，这样当用户将 zoom 缩扩到 10+ 或其他值下，本该一个鼠标点击就可以命中那个feature的，而实际却没有。
假设地图上有一个 “惠” 字，假设我在不同zoom下总能用鼠标命中同一个位置(x,y)，那么在不同zoom下请看效果：
-
zoom=22时：明显无法命中该字
image.png -
zoom=16时：明显可以命中该字
image.png
- 要解决思路1的问题，需要 引入 zoom 和 pixSize 共同将请求中的point 换算扩展成适当的 bbox ，再去判断空间相交关系即可。
解决：
已知：赤道处：y=0°纬度值
地球周长：40075016 米
切片规格：512*512，即 (29)2
计算 y=0,zoom=0处的“米每像素”的值，设为 metersPerPixel？
公式为：
var metersPerPixel = function(latitude, zoomLevel) {
var earthCircumference = 40075017;
var latitudeRadians = latitude * (Math.PI/180);
return earthCircumference * Math.cos(latitudeRadians) / Math.pow(2, zoomLevel + 9);
};
已知屏幕分辨率19201080，鼠标在浏览器上手抓的大小固定为20pix20pix大小，所以这里 pixValue = 20;
计算 y=0,zoom=7处的“鼠标抓取范围边长”的值？
在 y=0,zoom=7时，metersPerPixel = 611，所以 边长为 20 * 611 = 12220 米；
计算 y=34.578915,zoom=7处的“鼠标抓取范围边长”的值？
在 y=34.578915,zoom=7时，metersPerPixel = 503.654，所以 边长为 20 * 503.654 = 10073.1 米；
image.png
上图赤道处，米每像素的值，可以参考 https://docs.mapbox.com/help/glossary/zoom-level/
可以试着计算非赤道处的该值；
可以试着计算 256256 切片规格下的该值，256256，即 (28)2；
如果在手机上浏览，可以试着将pixValue改为手指触碰后抓取范围，然后计算看看。
java代码：
proj4j
<dependency>
<groupId>org.osgeo</groupId>
<artifactId>proj4j</artifactId>
<version>0.1.0</version>
</dependency>
Proj4jUtil
import org.osgeo.proj4j.*;
/**
* @ClassName Proj4jUtil
* @Author ckk
* @Date 2024/1/9 13:38
*/
public class Proj4jUtil {
public static CRSFactory factory = new CRSFactory();
public static CoordinateReferenceSystem sourceCRS4326 = factory.createFromName("EPSG:4326"); // WGS84坐标系统
public static CoordinateReferenceSystem targetCRS3857 = factory.createFromName("EPSG:3857"); // EPSG:3857坐标系统
public static CoordinateTransformFactory ctFactory = new CoordinateTransformFactory();
public static CoordinateTransform tranform4326To3857 = ctFactory.createTransform(sourceCRS4326, targetCRS3857);
public static CoordinateTransform tranform3857To4326 = ctFactory.createTransform(targetCRS3857, sourceCRS4326);
public static ProjCoordinate convertTo3857(double longitude, double latitude) {
ProjCoordinate transform = tranform4326To3857.transform(new ProjCoordinate(longitude, latitude), new ProjCoordinate());
return transform;
}
public static ProjCoordinate convertTo4326(double x, double y) {
ProjCoordinate transform = tranform3857To4326.transform(new ProjCoordinate(x, y), new ProjCoordinate());
return transform;
}
// public static void main(String[] args) {
// ProjCoordinate projCoordinate = Proj4jUtil.convertTo3857(117.10150553692836, 34.3724595227);
// ProjCoordinate projCoordinate1 = Proj4jUtil.convertTo4326(13035679.96749657, 4078924.561071809);
// System.out.println("11");
// }
}
MapBoxVectorTileUtil
import org.osgeo.proj4j.ProjCoordinate;
import java.util.HashMap;
import java.util.Map;
/**
* 缩放级别决定了在地图上可以看到多少世界。Mapbox提供23个缩放级别的地图，其中0是最低的缩放级别（完全缩小），22
* 是最高的（完全放大）。在低缩放级别时，一小套地图瓷砖覆盖了大片地理区域。在更高的缩放级别，更多的瓷砖覆盖了较小的地理区域。
* @ClassName a
* @Author ckk
* @Date 2023/12/22 13:36
*/
public class MapBoxVectorTileUtil {
public final static int earthCircumference = 40075016;
public final static int tileSize = 512;
/**
* 计算点击处的抓取范围，默认返回 3857 投影下的范围
* 默认支持 512*512 切片尺寸的gis画布，如需扩展要更改 tileSize 参数
* @param x 经度
* @param y 纬度
* @param zoom 缩放级别
* @param pixSize 像素大小（抓取工具的像素大小）
* @param output4326 是否返回 4326 投影下的范围，默认为 false
* @return
*/
public static String logLat2pixelBbox(double x, double y, double zoom, int pixSize, boolean output4326){
// 计算 米每像素
double metersPerPixel = metersPerPixel(y, zoom);
// 4326转3857
ProjCoordinate projCoordinate = Proj4jUtil.convertTo3857(x, y);
// 计算外扩
double extent = metersPerPixel * pixSize / 2;
// 计算鼠标或触摸的抓取范围
double minX = -extent + projCoordinate.x;
double minY = -extent + projCoordinate.y;
double maxX = extent + projCoordinate.x;
double maxY = extent + projCoordinate.y;
StringBuilder bbox = new StringBuilder();
if(!output4326){
// 构造bbox
bbox.append("SRID=3857;POLYGON((");
bbox.append(minX).append(" ").append(minY).append(", ");
bbox.append(minX).append(" ").append(maxY).append(", ");
bbox.append(maxX).append(" ").append(maxY).append(", ");
bbox.append(maxX).append(" ").append(minY).append(", ");
bbox.append(minX).append(" ").append(minY);
bbox.append("))");
}else{
// 3857转4326
ProjCoordinate leftBottom = Proj4jUtil.convertTo4326(minX, minY);
ProjCoordinate leftTop = Proj4jUtil.convertTo4326(minX, maxY);
ProjCoordinate rightTop = Proj4jUtil.convertTo4326(maxX, maxY);
ProjCoordinate rightBottom = Proj4jUtil.convertTo4326(maxX, minY);
// 构造bbox
bbox.append("SRID=4326;POLYGON((");
bbox.append(leftBottom.x).append(" ").append(leftBottom.y).append(", ");
bbox.append(leftTop.x).append(" ").append(leftTop.y).append(", ");
bbox.append(rightTop.x).append(" ").append(rightTop.y).append(", ");
bbox.append(rightBottom.x).append(" ").append(rightBottom.y).append(", ");
bbox.append(leftBottom.x).append(" ").append(leftBottom.y);
bbox.append("))");
}
return bbox.toString();
}
/**
* 计算y纬度和zoom下，米每像素的值
* @param y 纬度
* @param zoom 缩放级别
* @return 米每像素
*/
public static double metersPerPixel(double y, double zoom){
double latitudeRadians = y * Math.PI / 180;
double metersPerPixel = earthCircumference * Math.cos(latitudeRadians) / (Math.pow(2, zoom) * tileSize);
return metersPerPixel;
}
// public static void main(String[] args) {
// // 测试 19.137666927404464 下能够命中 惠河
// String bbox = logLat2pixelBbox(117.10702183509386, 34.371835846018826, 19.137666927404464, 20, false);
// System.out.println(bbox);
//
// // 测试 20.482818424170375 下不能够命中 惠河
// String bbox1 = logLat2pixelBbox(117.10702183509386, 34.371835846018826, 20.482818424170375, 20, false);
// System.out.println(bbox1);
// }
}
