Openlayers raster reprojection 栅格重投影原理分析

背景

OpenLayers 3 新版本能够以不同于服务器提供的坐标系显示来自WMS,WMTS,静态图像和许多其他来源的栅格数据。图像的地图投影的转换直接在Web浏览器中进行。任何Proj4js支持的坐标参考系统中的视图都是可能的,并且以前不兼容的图层现在可以组合和覆盖。

Raster Reproject意义

可以在浏览器端实现不同资源在不同投影下的转换,不再依赖于服务端处理,比如在Geoserver中指定投影类型,还是相当给力的

效果图

image

还可以看以下效果对比图


4326
3857

3857使用web Mercator 中国区域是方形,Reproject到4326时,就有些扁了,在相同的中心点,缩放级别下,地图的视察范围也不太一致,4326经过压扁后,能显示更多的范围。

再来一张放大后的,变形明显。

image

分析过程

reprojection-image为例说明,使用的是三角仿射变换triangle affine transformation。

三角形动态计算,以正好铺满当前视口

不同的缩放级别1


不同的缩放级别1

不同的缩放级别2


不同的缩放级别2

不同的缩放级别3


不同的缩放级别3

不同的缩放级别4


不同的缩放级别4

不同的比例尺下会更新三角格网的大小,以恰好平分地图视口,至少两个三角

每个三角形单独进行Reproject,不同三角之间互不影响

对三角格网渲染进行Degbu过滤后示意

渲染其中一部分三角


渲染偶数三角

这个图能方便明白Reproject原理

  • 首先划分三角格网
  • 根据位置,计算每个三角网对应的原始raster位置和区域,内部使用数学方法,求出转换参数,对每个三角网进行转换,这里其实并没有对每个像素进行Reproject,而是以三角为最小单位,进行Reproject, 理论上只要三角足够小,通过积分是可以达到同样的效果的,相比像素角度的处理效率也会高很多
  • 将多个三角进行无缘拼接

核心实现

三角Reproject

源码ol/reproj.js, render funcction

/**
 * Renders the source data into new canvas based on the triangulation.
 *
 * @param {number} width Width of the canvas.
 * @param {number} height Height of the canvas.
 * @param {number} pixelRatio Pixel ratio.
 * @param {number} sourceResolution Source resolution.
 * @param {import("./extent.js").Extent} sourceExtent Extent of the data source.
 * @param {number} targetResolution Target resolution.
 * @param {import("./extent.js").Extent} targetExtent Target extent.
 * @param {import("./reproj/Triangulation.js").default} triangulation
 * Calculated triangulation.
 * @param {Array<{extent: import("./extent.js").Extent,
 *                 image: (HTMLCanvasElement|HTMLImageElement|HTMLVideoElement)}>} sources
 * Array of sources.
 * @param {number} gutter Gutter of the sources.
 * @param {boolean=} opt_renderEdges Render reprojection edges.
 * @return {HTMLCanvasElement} Canvas with reprojected data.
 */
export function render(width, height, pixelRatio,
  sourceResolution, sourceExtent, targetResolution, targetExtent,
  triangulation, sources, gutter, opt_renderEdges) {

  var context = createCanvasContext2D(Math.round(pixelRatio * width),
    Math.round(pixelRatio * height));

  if (sources.length === 0) {
    return context.canvas;
  }

  context.scale(pixelRatio, pixelRatio);

  var sourceDataExtent = createEmpty();
  sources.forEach(function(src, i, arr) {
    extend(sourceDataExtent, src.extent);
  });

  var canvasWidthInUnits = getWidth(sourceDataExtent);
  var canvasHeightInUnits = getHeight(sourceDataExtent);
  
  // 加载image的内部canvas
  var stitchContext = createCanvasContext2D(
    Math.round(pixelRatio * canvasWidthInUnits / sourceResolution),
    Math.round(pixelRatio * canvasHeightInUnits / sourceResolution));

  var stitchScale = pixelRatio / sourceResolution;

  sources.forEach(function(src, i, arr) {
    var xPos = src.extent[0] - sourceDataExtent[0];
    var yPos = -(src.extent[3] - sourceDataExtent[3]);
    var srcWidth = getWidth(src.extent);
    var srcHeight = getHeight(src.extent);

    stitchContext.drawImage(
      src.image,
      gutter, gutter,
      src.image.width - 2 * gutter, src.image.height - 2 * gutter,
      xPos * stitchScale, yPos * stitchScale,
      srcWidth * stitchScale, srcHeight * stitchScale);
  });

  var targetTopLeft = getTopLeft(targetExtent);

  // 对三角进行循环处理
  triangulation.getTriangles().forEach(function(triangle, i, arr) {

    /* Calculate affine transform (src -> dst)
     * Resulting matrix can be used to transform coordinate
     * from `sourceProjection` to destination pixels.
     *
     * To optimize number of context calls and increase numerical stability,
     * we also do the following operations:
     * trans(-topLeftExtentCorner), scale(1 / targetResolution), scale(1, -1)
     * here before solving the linear system so [ui, vi] are pixel coordinates.
     *
     * Src points: xi, yi
     * Dst points: ui, vi
     * Affine coefficients: aij
     *
     * | x0 y0 1  0  0 0 |   |a00|   |u0|
     * | x1 y1 1  0  0 0 |   |a01|   |u1|
     * | x2 y2 1  0  0 0 | x |a02| = |u2|
     * |  0  0 0 x0 y0 1 |   |a10|   |v0|
     * |  0  0 0 x1 y1 1 |   |a11|   |v1|
     * |  0  0 0 x2 y2 1 |   |a12|   |v2|
     */
    var source = triangle.source;
    var target = triangle.target;
    var x0 = source[0][0], y0 = source[0][1];
    var x1 = source[1][0], y1 = source[1][1];
    var x2 = source[2][0], y2 = source[2][1];
    var u0 = (target[0][0] - targetTopLeft[0]) / targetResolution;
    var v0 = -(target[0][1] - targetTopLeft[1]) / targetResolution;
    var u1 = (target[1][0] - targetTopLeft[0]) / targetResolution;
    var v1 = -(target[1][1] - targetTopLeft[1]) / targetResolution;
    var u2 = (target[2][0] - targetTopLeft[0]) / targetResolution;
    var v2 = -(target[2][1] - targetTopLeft[1]) / targetResolution;

    // Shift all the source points to improve numerical stability
    // of all the subsequent calculations. The [x0, y0] is used here.
    // This is also used to simplify the linear system.
    var sourceNumericalShiftX = x0;
    var sourceNumericalShiftY = y0;
    x0 = 0;
    y0 = 0;
    x1 -= sourceNumericalShiftX;
    y1 -= sourceNumericalShiftY;
    x2 -= sourceNumericalShiftX;
    y2 -= sourceNumericalShiftY;

    var augmentedMatrix = [
      [x1, y1, 0, 0, u1 - u0],
      [x2, y2, 0, 0, u2 - u0],
      [0, 0, x1, y1, v1 - v0],
      [0, 0, x2, y2, v2 - v0]
    ];
    var affineCoefs = solveLinearSystem(augmentedMatrix);
    if (!affineCoefs) {
      return;
    }
    context.save();
    context.beginPath();
    var centroidX = (u0 + u1 + u2) / 3;
    var centroidY = (v0 + v1 + v2) / 3;
    var p0 = enlargeClipPoint(centroidX, centroidY, u0, v0);
    var p1 = enlargeClipPoint(centroidX, centroidY, u1, v1);
    var p2 = enlargeClipPoint(centroidX, centroidY, u2, v2);

    // 设置三角的切割范围,只显示这个范围内的图片,多个三角拼接起来就是整体图像
    context.moveTo(p1[0], p1[1]);
    context.lineTo(p0[0], p0[1]);
    context.lineTo(p2[0], p2[1]);
    context.clip();

    // 最关键的三个方法,通过一系列转换,保证在相应的比例尺下在当前范围内
    // 显示正确的Reproject图像
    // 直接设置下面3个参数,不好容易能达到想要的效果
    // 可以通过将相同的参数输出在本地,进行调试,观察三角内的图像
    context.transform(
      affineCoefs[0], affineCoefs[2], affineCoefs[1], affineCoefs[3], u0, v0);

    context.translate(sourceDataExtent[0] - sourceNumericalShiftX,
      sourceDataExtent[3] - sourceNumericalShiftY);

    context.scale(sourceResolution / pixelRatio,
      -sourceResolution / pixelRatio);

    context.drawImage(stitchContext.canvas, 0, 0);
    context.restore();
  });

  // 调试使用,是否显示三角
  if (opt_renderEdges) {
    context.save();

    context.strokeStyle = 'black';
    context.lineWidth = 1;

    triangulation.getTriangles().forEach(function(triangle, i, arr) {
      var target = triangle.target;
      var u0 = (target[0][0] - targetTopLeft[0]) / targetResolution;
      var v0 = -(target[0][1] - targetTopLeft[1]) / targetResolution;
      var u1 = (target[1][0] - targetTopLeft[0]) / targetResolution;
      var v1 = -(target[1][1] - targetTopLeft[1]) / targetResolution;
      var u2 = (target[2][0] - targetTopLeft[0]) / targetResolution;
      var v2 = -(target[2][1] - targetTopLeft[1]) / targetResolution;

      context.beginPath();
      context.moveTo(u1, v1);
      context.lineTo(u0, v0);
      context.lineTo(u2, v2);
      context.closePath();
      context.stroke();
    });

    context.restore();
  }
  return context.canvas;
}

我以其中一个三角的计算参数为示例,可以看到参数还是比较复杂的

context.moveTo(1136, 199);
context.lineTo(850, 200);
context.lineTo(1136, 401);
context.clip();

context.transform(
  0.011129369548469296, 0.00006403305449769778, 0.00009456173415458175, -0.0111338055978892, 851.2500000000001, 200.00000000000063);

context.translate(-413988.7802610396,
  835088.5497620346);

context.scale(351.73160173160176,
  -351.73160173160176);

context.drawImage(stitchContext.canvas, 0, 0);
context.restore();
image

其中,进行仿射变换的转换过程及数学原理,还需要进一步研究。

参考

©著作权归作者所有,转载或内容合作请联系作者
  • 序言:七十年代末,一起剥皮案震惊了整个滨河市,随后出现的几起案子,更是在滨河造成了极大的恐慌,老刑警刘岩,带你破解...
    沈念sama阅读 212,080评论 6 493
  • 序言:滨河连续发生了三起死亡事件,死亡现场离奇诡异,居然都是意外死亡,警方通过查阅死者的电脑和手机,发现死者居然都...
    沈念sama阅读 90,422评论 3 385
  • 文/潘晓璐 我一进店门,熙熙楼的掌柜王于贵愁眉苦脸地迎上来,“玉大人,你说我怎么就摊上这事。” “怎么了?”我有些...
    开封第一讲书人阅读 157,630评论 0 348
  • 文/不坏的土叔 我叫张陵,是天一观的道长。 经常有香客问我,道长,这世上最难降的妖魔是什么? 我笑而不...
    开封第一讲书人阅读 56,554评论 1 284
  • 正文 为了忘掉前任,我火速办了婚礼,结果婚礼上,老公的妹妹穿的比我还像新娘。我一直安慰自己,他们只是感情好,可当我...
    茶点故事阅读 65,662评论 6 386
  • 文/花漫 我一把揭开白布。 她就那样静静地躺着,像睡着了一般。 火红的嫁衣衬着肌肤如雪。 梳的纹丝不乱的头发上,一...
    开封第一讲书人阅读 49,856评论 1 290
  • 那天,我揣着相机与录音,去河边找鬼。 笑死,一个胖子当着我的面吹牛,可吹牛的内容都是我干的。 我是一名探鬼主播,决...
    沈念sama阅读 39,014评论 3 408
  • 文/苍兰香墨 我猛地睁开眼,长吁一口气:“原来是场噩梦啊……” “哼!你这毒妇竟也来了?” 一声冷哼从身侧响起,我...
    开封第一讲书人阅读 37,752评论 0 268
  • 序言:老挝万荣一对情侣失踪,失踪者是张志新(化名)和其女友刘颖,没想到半个月后,有当地人在树林里发现了一具尸体,经...
    沈念sama阅读 44,212评论 1 303
  • 正文 独居荒郊野岭守林人离奇死亡,尸身上长有42处带血的脓包…… 初始之章·张勋 以下内容为张勋视角 年9月15日...
    茶点故事阅读 36,541评论 2 327
  • 正文 我和宋清朗相恋三年,在试婚纱的时候发现自己被绿了。 大学时的朋友给我发了我未婚夫和他白月光在一起吃饭的照片。...
    茶点故事阅读 38,687评论 1 341
  • 序言:一个原本活蹦乱跳的男人离奇死亡,死状恐怖,灵堂内的尸体忽然破棺而出,到底是诈尸还是另有隐情,我是刑警宁泽,带...
    沈念sama阅读 34,347评论 4 331
  • 正文 年R本政府宣布,位于F岛的核电站,受9级特大地震影响,放射性物质发生泄漏。R本人自食恶果不足惜,却给世界环境...
    茶点故事阅读 39,973评论 3 315
  • 文/蒙蒙 一、第九天 我趴在偏房一处隐蔽的房顶上张望。 院中可真热闹,春花似锦、人声如沸。这庄子的主人今日做“春日...
    开封第一讲书人阅读 30,777评论 0 21
  • 文/苍兰香墨 我抬头看了看天上的太阳。三九已至,却和暖如春,着一层夹袄步出监牢的瞬间,已是汗流浃背。 一阵脚步声响...
    开封第一讲书人阅读 32,006评论 1 266
  • 我被黑心中介骗来泰国打工, 没想到刚下飞机就差点儿被人妖公主榨干…… 1. 我叫王不留,地道东北人。 一个月前我还...
    沈念sama阅读 46,406评论 2 360
  • 正文 我出身青楼,却偏偏与公主长得像,于是被迫代替她去往敌国和亲。 传闻我的和亲对象是个残疾皇子,可洞房花烛夜当晚...
    茶点故事阅读 43,576评论 2 349

推荐阅读更多精彩内容