GIS - 地圖開發

--- tags: 前端 --- # GIS - 地圖開發 > GIS 地理資訊系統(Geographic Information System) ## 概念向量(vector) | 網格(raster) -|- 縮放不失真、儲存使用量小 | 渲染效率高，容易切割資源檔 ![](http://www.seekacreative.co.nz/uploads/3/0/6/7/30679021/vector-raster_orig.jpg) ### 座標參考系統CRS: * 台灣常用的：TWD97, WGS84, TWD67 * latitude(lat)：緯度 * longitude(lng、lon)：經度 * altitude：高度 * 投影法：麥卡托(Mercator)(方向準確) * 大圓距離[Haversine](https://www.baeldung.com/cs/haversine-formula) * 線上地圖服務WMS：openstreet、stamen、cartodb、esri、mapbox、google * 檔案格式 * shapefile：地理資訊系統始祖ESRI公司發明的 * KML：谷歌推出的XML結構 * geojson：開源通用JSON結構 * Point | MultiPoint * LineSting | MultiLineString * Polygon | MultiPolygon ```json { type:'FeatureCollection', bbox:[lng1,lat1,lng2,lat2], features:[ { type:'feature', geometry:{ type:'LineString', coordinates:[ [lng,lat], [lng,lat] ] }, properties: {} } ] } ``` ## WMS資源圖層|網址 -|- 特殊管制 | https://eghouse.hccg.gov.tw/arcgis/rest/services/Tiled3857/Nature_policy/MapServer/tile/{z}/{y}/{x} 水源水質 | https://eghouse.hccg.gov.tw/arcgis/rest/services/Tiled3857/Nature_water/MapServer/tile/{z}/{y}/{x} 地質敏感 | https://eghouse.hccg.gov.tw/arcgis/rest/services/Tiled3857/Nature_geo/MapServer/tile/{z}/{y}/{x} 環保與汙染 | https://eghouse.hccg.gov.tw/arcgis/rest/services/Tiled3857/Nature_environment/MapServer/tile/{z}/{y}/{x} 坡度示意圖 | https://eghouse.hccg.gov.tw/arcgis/rest/services/Tiled3857/SlopeTW3857_fix/MapServer/tile/{z}/{y}/{x} 國土利用 | https://wmts.nlsc.gov.tw/wmts/LUIMAP/default/EPSG:3857/{z}/{y}/{x} 公有土地 | https://wmts.nlsc.gov.tw/wmts/LAND_OPENDATA/default/EPSG:3857/{z}/{y}/{x} 學校 | https://wmts.nlsc.gov.tw/wmts/SCHOOL/default/EPSG:3857/{z}/{y}/{x} 村里界 | https://wmts.nlsc.gov.tw/wmts/Village/default/EPSG:3857/{z}/{y}/{x} 鄉市鎮界 | https://wmts.nlsc.gov.tw/wmts/TOWN/default/EPSG:3857/{z}/{y}/{x} 縣市界 | https://wmts.nlsc.gov.tw/wmts/CITY/default/EPSG:3857/{z}/{y}/{x} OSM | https://{a-c}.tile.openstreetmap.org/{z}/{x}/{y}.png 通用 | https://wmts.nlsc.gov.tw/wmts/EMAP5/default/EPSG:3857/{z}/{y}/{x} 灰階 | https://wmts.nlsc.gov.tw/wmts/EMAP01/default/EPSG:3857/{z}/{y}/{x} 航照 | https://wmts.nlsc.gov.tw/wmts/PHOTO2/default/EPSG:3857/{z}/{y}/{x} Google衛星 | https://mts1.google.com/vt/lyrs=s@186112443&hl=x-local&src=app&x={x}&y={y}&z={z}&s=Galile Google地形 | https://mts1.google.com/vt/lyrs=p@186112443&hl=x-local&src=app&x={x}&y={y}&z={z}&s=Galile ESRI衛星 | https://server.arcgisonline.com/ArcGIS/rest/services/World_Imagery/MapServer/tile/{z}/{y}/{x} 交通路網 | https://wmts.nlsc.gov.tw/wmts/EMAPX99/default/EPSG:3857/{z}/{y}/{x} ## Leaflet + 氣象地圖 * 可用python的folium產生網頁檔 * 開源近10年，擴充設計完善(熱力圖、叢集)，適用於視覺化 * API設計相當簡單，適合入門 * [Leaflet](https://leafletjs.com/) * [Leaflet.markercluster](https://github.com/Leaflet/Leaflet.markercluster) * [Leaflet.heat](https://github.com/Leaflet/Leaflet.heat) * [topojson](https://github.com/topojson/topojson) 引入套件 ```html <link rel="stylesheet" href="https://unpkg.com/leaflet@1.7.1/dist/leaflet.css" /> <script src="https://unpkg.com/leaflet@1.7.1/dist/leaflet.js"></script> <link rel="stylesheet" href="https://unpkg.com/leaflet.markercluster@1.4.1/dist/MarkerCluster.Default.css" /> <script src="https://unpkg.com/leaflet.markercluster@1.4.1/dist/leaflet.markercluster.js"></script> <script src="https://www.unpkg.com/leaflet.heat@0.2.0/dist/leaflet-heat.js"></script> <script src="https://unpkg.com/topojson@3.0.2/dist/topojson.min.js"></script> ``` 配置畫面 ```htmlembedded= <body> <div id="map"></div> </body> <style> html, body, #map { margin: 0; width: 100%; height: 100%; } </style> ``` 初始化地圖 ```javascript= const map = L.map('map').setView([24, 121], 10); L.tileLayer('https://{s}.tile.openstreetmap.org/{z}/{x}/{y}.png', { attribution: '氣象地圖' }).addTo(map); const markerCluster = L.markerClusterGroup().addTo(map); ``` 撰寫Leaflet的[擴充套件](https://gist.github.com/rclark/5779673) ```javascript= L.TopoJSON = L.GeoJSON.extend({ addData(jsonData) { if (jsonData.type === "Topology") { for (let key in jsonData.objects) { let geojson = topojson.feature(jsonData, jsonData.objects[key]); L.GeoJSON.prototype.addData.call(this, geojson); } } else L.GeoJSON.prototype.addData.call(this, jsonData); } }); L.topoJson = (data, options)=>new L.TopoJSON(data, options); ``` 查詢氣象API的[中英文名對照表](https://opendata.cwb.gov.tw/opendatadoc/DIV2/A0001-001.pdf) 選取多行`ctrl+alt+↓`、選取相同字`ctrl+d` ```javascript= const ZH_TW = { lat: '緯度', lon: '經度', locationName: '測站名稱', stationId: '測站ID', obsTime: '觀測時間', CITY: '縣市', TOWN: '地區', ELE: '高度', WDIR: '風向', WDSD: '風速', TEMP: '溫度', HUMD: '相對濕度', PRES: '測站氣壓', H_24R: '日累積雨量', H_FX: '小時最大陣風風速', H_XD: '小時最大陣風風向', H_FXT: '小時最大陣風時間', D_TX: '本日最高溫', D_TXT: '本日最高溫發生時間', D_TN: '本日最低溫', D_TNT: '本日最低溫發生時間', CITY: '縣市', CITY_SN: '縣市編號', TOWN: '鄉鎮', TOWN_SN: '鄉鎮編號', } const NAME = ['WDIR', 'WDSD', 'TEMP', 'HUMD', 'PRES', 'H_24R', 'H_FX', 'H_XD', 'D_TX', 'D_TN'] ``` 取得氣象API ```javascript= async function getWeather(markerCluster) { let json = await fetch('https://opendata.cwb.gov.tw/api/v1/rest/datastore/O-A0001-001?Authorization=rdec-key-123-45678-011121314').then(res => res.json()) let { success, result: { fields }, records: { location } } = json return location.map(loc => { let { lat, lon, locationName, parameter, stationId, time: { obsTime }, weatherElement } = loc let addr = Object.fromEntries(parameter.map(x => [x.parameterName, x.parameterValue])) let weather = Object.fromEntries(weatherElement.map(x => [x.elementName, x.elementValue])) let data = { lat, lon, locationName, stationId, ...addr, ...weather, obsTime } let table = document.createElement('table') table.className = 'table' for (let key in data) { let row = table.insertRow() row.insertCell().textContent = ZH_TW[key] row.insertCell().textContent = data[key] == -99 ? '無資料' : data[key] } L.marker([lat, lon]).bindPopup(table).addTo(markerCluster) return data }) } ``` 分層設色 ```javascript= function getChoroplenth(topoJson = {}, colorMap = {}, popupMap = {}) { let config = { style: (feature) => ({ color: "#000", opacity: 1, weight: 1, fillColor: colorMap[feature.properties.name], fillOpacity: 0.8 }) , onEachFeature: (feature, layer) => layer.bindPopup(feature.properties.name + ':' + popupMap[feature.properties.name]) } return L.topoJson(topoJson, config).addTo(map) } ``` 工具函式 ```javascript= // 色票轉換 https://stackoverflow.com/questions/36721830/convert-hsl-to-rgb-and-hex function hslToHex(h, s, l) { l /= 100; const a = s * Math.min(l, 1 - l) / 100; const f = n => { const k = (n + h / 30) % 12; const color = l - a * Math.max(Math.min(k - 3, 9 - k, 1), -1); return Math.round(255 * color).toString(16).padStart(2, '0'); }; return `#${f(0)}${f(8)}${f(4)}`; } function groupBy(rows = [], key = '') { let group = {} for (let row of rows) { let rKey = row[key] if (!(rKey in group)) group[rKey] = [] group[rKey].push(row) } return group } function mean(arr = []) { return arr.reduce((acc, cur) => acc + cur, 0) / arr.length } function normalize(arr = []) { let min = Math.min(...arr) let max = Math.max(...arr) return arr.map(x => (x - min) / (max - min)) } ``` 主程式 ```javascript= async function mounted() { let topo = await fetch('https://raw.githubusercontent.com/g0v/twgeojson/master/json/twCounty2010merge.topo.json').then(res => res.text()) topo = JSON.parse(topo.replace(/台/g, '臺').replace(/桃園縣/g, '桃園市')) let rows = await getWeather(markerCluster); const heat = L.heatLayer(rows.map(row => [row.lat, row.lon])).addTo(map); // 定義控制項 const checkboxes = { '標記叢集': markerCluster, '熱力圖': heat } const radios = { '無分層設色': L.layerGroup() } // 群組化 let cityGroup = groupBy(rows, 'CITY') for (let name of NAME) { let meanMap = {} let colorMap = {} let popupMap = {} for (let key in cityGroup) { let arr = cityGroup[key].map(row => parseFloat(row[name])).filter(x => x != -99 && !isNaN(x)) meanMap[key] = mean(arr) } let means = Object.values(meanMap) let norms = normalize(means) Object.keys(meanMap).map((key, i) => { colorMap[key] = hslToHex(norms[i] * 360, 50, 50) popupMap[key] = means[i].toFixed(2) }) radios[ZH_TW[name]] = getChoroplenth(topo, colorMap, popupMap); } map.fitBounds(Object.values(radios)[1].getBounds()); L.control.layers(radios, checkboxes).addTo(map); } mounted() ``` ## Google Map > 谷歌地圖開發者服務，適合用於開發多功能資訊整合 * 地圖JS * 街景 * 商家景點 * 地理編碼(geocoding) * 路徑規劃(direction) 引入套件 ```htmlembedded= <div id="map"></div> <script src="https://maps.googleapis.com/maps/api/js?key={API_KEY}&callback=initMap"></script> <script src="https://developers.google.com/maps/documentation/javascript/examples/markerclusterer/markerclusterer.js"></script> ``` 基本圖層 ```javascript= var map = new google.maps.Map(document.getElementById('map'),{center:{lat,lng},zoom:8,mapTypeId: 'terrain'}) // 標記叢集 var marker = new google.maps.Marker({position: {lat,lng}, map}) var imagePath = 'https://developers.google.com/maps/documentation/javascript/examples/markerclusterer/m' var markerCluster = new MarkerClusterer(map, [marker], {imagePath}); // 熱力圖 new google.maps.visualization.HeatmapLayer({data: [{lat,lng}],dissipating:false,map: map}) // 匯入資料 map.data.addGeoJson(GeoJson) // 事件監聽 map.addEventListener('click',e=>{e.latlng.lat()}) new google.maps.Polyline({ path: flightPlanCoordinates, geodesic: true, strokeColor: '#FF0000', strokeOpacity: 1.0, strokeWeight: 2 }).setMap(map) //加上setMap(null)移除、加上getLength取得長度(Polygon) new google.maps.Rectangle({ bounds: {north: 33.685, south: 33.671, east: -116.234,west: -116.251}, editable: true }).setMap(map); ``` 路徑規劃 ```javascript= var request = { origin: {lat,lng}, //傳入關鍵字(用+號分隔)、經緯度、place_id destination: end, travelMode: 'DRIVING' //DRIVING, WALKING, BICYCLING, TRANSIT } var directionsService = new google.maps.DirectionsService(); var directionsDisplay = new google.maps.DirectionsRenderer(); directionsDisplay.setMap(map) directionsService.route(request, function(result, status) { if (status == 'OK') { // 另一個狀態為"ZERO_RESULTS" directionsDisplay.setDirections(result); result[i].legs[j] // {start_address, end_address, distance{text}} var distance = result[0].legs.reduce((acc,cur)=>acc+cur.distance.value,0) } }) ``` ## OpenLaysers * 收錄了完整地圖功能模組，適合開發地圖編輯軟體 * 資料結構高達8層，開發相對不易 ```graphviz digraph{ map->layers; map->controls; map->interactions; map->overlays; map->view; view->projection, center, zoom, rotation overlays->element, offset, positioning layers->source; source->raster; source->cluster; source->vector; cluster->vector; raster->XYZ; vector->feature; feature->style; feature->geom; geom->coord; geom->type; type->Point, LineString, MultiLineString, Polygon, MultiPolygon } ``` ## 手刻地圖 - 核心數學 ```javascript= const tilePixel = { w: 256, h: 256 } const tileSize = 6378137 * Math.PI const World = { bbox: [-tileSize, -tileSize, tileSize, tileSize], w: 2 * tileSize, h: 2 * tileSize, } const View = { bbox: [], center: { x: 13437548.485305637, y: 2772337.9239074644 }, w: this.canvas.clientWidth, h: this.canvas.clientHeight, zoom: 9, minZoom: 0, maxZoom: 20, } class Point{ toLngLat(){ var d = 180 / Math.PI return [ point[0] * d / 6378137, (2 * Math.atan(Math.exp(point[1] / 6378137)) - (Math.PI / 2)) * d ] } getArea(aoa) { var radius = 6371008.8, area = 0 let x1 = aoa[aoa.length - 1][0] let y1 = aoa[aoa.length - 1][1] for (let i = 0; i < aoa.length; i++) { let x2 = aoa[i][0] let y2 = aoa[i][1] area += (x2 - x1) * Math.PI / 180 * (2 + Math.sin(y1 * Math.PI / 180) + Math.sin(y2 * Math.PI / 180)) x1 = x2 y1 = y2 } return Math.abs(area * radius * radius / 2.0) } client2point(client) { var view = this.view, W = view.bbox[2] - view.bbox[0], H = view.bbox[1] - view.bbox[3] return [ view.bbox[0] + client[0] * W / view.w, view.bbox[3] + client[1] * H / view.h ] } point2client(point) { var view = this.view, W = view.bbox[2] - view.bbox[0], H = view.bbox[1] - view.bbox[3] return [ (point[0] - view.bbox[0]) / W * view.w, (point[1] - view.bbox[3]) / H * view.h ] } } class LngLat{ toPoint(lnglat) { // SphericalMercator var d = Math.PI / 180, max = 85.0511287798, lat = Math.max(Math.min(max, lnglat[1]), -max), sin = Math.sin(lat * d); return [ 6378137 * lnglat[0] * d, 6378137 * Math.log((1 + sin) / (1 - sin)) / 2 ] } getLength(c1, c2){ var r = 6371008.8, lat1 = c1[1] * Math.PI / 180, lat2 = c2[1] * Math.PI / 180, dlat = (lat2 - lat1) / 2, dlon = (c2[0] - c1[0]) * Math.PI / 180 / 2, a = Math.sin(dlat) * Math.sin(dlat) + Math.cos(lat1) * Math.cos(lat2) * Math.sin(dlon) * Math.sin(dlon) return 2 * r * Math.atan2(Math.sqrt(a), Math.sqrt(1 - a)) } } class Bounds{ overlaps(bounds){ let x = (bounds.max.x > this.min.x) && (bounds.min.x < this.max.x) let y = (bounds.max.y > this.min.y) && (bounds.min.y < this.max.y) return x&&y } contains(point){ let x = (point.x > this.min.x) && (point.x < this.max.x) let y = (point.y > this.min.y) && (point.y < this.max.y) return x&&y } } function getBbox(aoa) { return [ Math.min(...aoa.map(c => c[0])), Math.min(...aoa.map(c => c[1])), Math.max(...aoa.map(c => c[0])), Math.max(...aoa.map(c => c[1])) ] } function drawSmoothPath(ctx, points, rate = 0.2) { const add = (p1, p2) => ([p1[0] + p2[0], p1[1] + p2[1]]) const subtract = (p1, p2) => ([p1[0] - p2[0], p1[1] - p2[1]]) const curveTo = (p1, p2, p3) => { ctx.bezierCurveTo(p1[0], p1[1], p2[0], p2[1], p3[0], p3[1]) } let cp = [] //control point for (let i = 0; i < points.length - 2; i++) { let p = subtract(points[i + 2], points[i]) cp.push([Math.round(p[0] * rate), Math.round(p[1] * rate)]) } points.map((p, i, ps) => { switch (i) { case 0: ctx.lineTo(p[0], p[1]); break; case 1: curveTo(ps[0], subtract(p, cp[i - 1]), p); break; case ps.length - 1: curveTo(add(ps[i - 1], cp[i - 2]), p, p); break; default: curveTo(add(ps[i - 1], cp[i - 2]), subtract(p, cp[i - 1]), p); break; } }) } const EasingFunctions = { // no easing, no acceleration linear: t => t, // accelerating from zero velocity easeInQuad: t => t * t, // decelerating to zero velocity easeOutQuad: t => t * (2 - t), // acceleration until halfway, then deceleration easeInOutQuad: t => t < .5 ? 2 * t * t : -1 + (4 - 2 * t) * t, // accelerating from zero velocity easeInCubic: t => t * t * t, // decelerating to zero velocity easeOutCubic: t => (--t) * t * t + 1, // acceleration until halfway, then deceleration easeInOutCubic: t => t < .5 ? 4 * t * t * t : (t - 1) * (2 * t - 2) * (2 * t - 2) + 1, // accelerating from zero velocity easeInQuart: t => t * t * t * t, // decelerating to zero velocity easeOutQuart: t => 1 - (--t) * t * t * t, // acceleration until halfway, then deceleration easeInOutQuart: t => t < .5 ? 8 * t * t * t * t : 1 - 8 * (--t) * t * t * t, // accelerating from zero velocity easeInQuint: t => t * t * t * t * t, // decelerating to zero velocity easeOutQuint: t => 1 + (--t) * t * t * t * t, // acceleration until halfway, then deceleration easeInOutQuint: t => t < .5 ? 16 * t * t * t * t * t : 1 + 16 * (--t) * t * t * t * t } function getShiftedCoords(coords, width) { var x2, y2, x42, y42 var shift_coord = [] for (let i = 0; i < coords.length; i++) { if (i == 0) { var y1 = coords[i][0] var x1 = coords[i][1] y2 = coords[i + 1][0] x2 = coords[i + 1][1] var s124 = width / Math.sqrt(Math.pow(x2 - x1, 2) + Math.pow(y2 - y1, 2)) x42 = (y2 - y1) * s124 y42 = -(x2 - x1) * s124 shift_coord.push([y1 + y42, x1 + x42]) } else if (i == coords.length - 1) { shift_coord.push([y2 + y42, x2 + x42]) } else { var y3 = coords[i + 1][0] var x3 = coords[i + 1][1] var s235 = width / Math.sqrt(Math.pow(x3 - x2, 2) + Math.pow(y3 - y2, 2)) var x52 = (y3 - y2) * s235 var y52 = -(x3 - x2) * s235 var radian = Math.atan2(x42 * y52 - y42 * x52, x42 * x52 + y42 * y52) var scalar = 1 / 2 / Math.pow(Math.cos(radian / 2), 2) var sx = (x42 + x52) * scalar var sy = (y42 + y52) * scalar shift_coord.push([y2 + sy, x2 + sx]) y2 = coords[i + 1][0] x2 = coords[i + 1][1] x42 = x52 y42 = y52 } } return shift_coord } function perpendicular(p, p1, p2) { var x = p1[0], y = p1[1], dx = p2[0] - p1[0], dy = p2[1] - p1[1] if (dx !== 0 || dy !== 0) { var t = ((p[0] - x) * dx + (p[1] - y) * dy) / (dx * dx + dy * dy) //scale of projection if (t > 1) { x = p2[0] y = p2[1] } else if (t > 0) { x += dx * t y += dy * t } } dx = p[0] - x dy = p[1] - y return { distance: Math.sqrt(dx * dx + dy * dy), closest: [x, y] } } function radialDistFilter(points, epsilon) { if (epsilon == 0) return points var epsilon_squ = epsilon * epsilon return points.filter((point, i, points) => { if (i == 0) return true var dx = point[0] - points[i - 1][0] var dy = point[1] - points[i - 1][1] return dx * dx + dy * dy > epsilon_squ }) } function DouglasPeucker(points, epsilon) { // points = this.radialDistFilter(points,epsilon) if (epsilon == 0) return points var dmax = 0, index = -1, firstPoint = points[0], lastPoint = points[points.length - 1] for (var i = 1; i < points.length; i++) { var d = this.perpendicular(points[i], firstPoint, lastPoint).distance if (d > dmax) { dmax = d index = i } } if (dmax >= epsilon) { return [ ...this.DouglasPeucker(points.slice(0, index + 1), epsilon).slice(0, -1), ...this.DouglasPeucker(points.slice(index), epsilon) ] } else { return [points[0], points[points.length - 1]] } } ``` ## Geolocation 網頁定位 API ```graphviz digraph{ navigator->geolocation; geolocation->getCurrentPosition,watchPosition,clearWatch; watchPosition->fn_success,fn_error,options; fn_success->position; position->coords; coords->latitude,longitude,accuracy,altitude,heading,speed,timestamp; options->enableHighAccuracy,timeout,maximumAge; } ```

Syntax	Example	Reference
# Header	Header	基本排版
- Unordered List	Unordered List
1. Ordered List	Ordered List
- [ ] Todo List	Todo List
> Blockquote	Blockquote
Bold font	Bold font
Italics font	Italics font
~~Strikethrough~~	~~Strikethrough~~
19^th^	19^th
H~2~O	H₂O
++Inserted text++	Inserted text
==Marked text==	Marked text
[link text](https:// "title")	Link
![image alt](https:// "title")	Image
`Code`	`Code`	在筆記中貼入程式碼
```javascript var i = 0; ```	`var i = 0;`
:smile:		Emoji list
{%youtube youtube_id %}	Externals
$L^aT_eX$	L^aT_eX
:::info This is a alert area. :::	This is a alert area.