JT-CIM-SDK/packages/Widgets/SpatialAnalysis/ViewShed.js

/* 引入Cesium */
// import * as Cesium from 'Cesium';
/* 引入克里金插值 */
import kriging from '@sakitam-gis/kriging';
/* 引入算法 */
import * as turf from "@turf/turf";

import CreateRemindertip from "../common/ReminderTip.js";

/**
 * 视域分析类
 */
class ViewShed {
	/**
	 * 默认初始化
	 * @param {Object} viewer 三维场景
	 */
	constructor(viewer) {
		if (!viewer) throw new Cesium.DeveloperError('no viewer object!');

		// let canvas = document.getElementById("canvasMap");
		// if (canvas == null) {
		// 	//创建画布
		// 	canvas = document.createElement('canvas');
		// 	canvas.id = "canvasMap";
		// 	canvas.style.display = "none";
		// 	canvas.style.top = "0px";
		// 	canvas.style.position = "absolute";
		// 	/* 加入到页面 */
		// 	document.body.append(canvas);
		// }
		// this.canvasEle = canvas;

		this.viewer = viewer;
		this.handler = undefined;
		this.lightCamera;
		this.pyramid;
		this.frustumPrimitive;
		this.viewershedPolygon;
	}

	/**
	 * 初始化handler
	 * @ignore 忽略注释，注释不生成Doc
	 */
	initHandler() {
		if (this.handler) {
			this.handler.destroy();
			this.handler = undefined;
		}

		let canvas = document.getElementById("canvasMap");
		if (canvas == null) {
			//创建画布
			canvas = document.createElement('canvas');
			canvas.id = "canvasMap";
			canvas.style.display = "none";
			canvas.style.top = "0px";
			canvas.style.position = "absolute";
			/* 加入到页面 */
			document.body.append(canvas);
		}
		this.canvasEle = canvas;
	}

	/**
	 * @ignore 忽略注释，注释不生成Doc
	 * @param {Object} pos0
	 * @param {Object} pos1
	 */
	ReturnDistance(pos0, pos1) {
		let distance = 0;
		let point1cartographic = Cesium.Cartographic.fromCartesian(pos0);
		let point2cartographic = Cesium.Cartographic.fromCartesian(pos1);
		let geodesic = new Cesium.EllipsoidGeodesic();
		geodesic.setEndPoints(point1cartographic, point2cartographic);
		let s = geodesic.surfaceDistance;
		return s;
	}

	/**
	 * @ignore 忽略注释，注释不生成Doc
	 * @param {Object} x
	 * @param {Object} y
	 * @param {Object} objectsToExclude
	 */
	getHeight(x, y, objectsToExclude) {
		let endCartographic = Cesium.Cartographic.fromDegrees(x, y);
		let endHeight = this.viewer.scene.sampleHeight(
			endCartographic,
			objectsToExclude
		);
		return endHeight;
	}

	/**
	 * @ignore 忽略注释，注释不生成Doc
	 * @param {Object} Cartesian3
	 */
	cartesian3ToDegree(Cartesian3) {
		let _ellipsoid = this.viewer.scene.globe.ellipsoid;
		let _cartographic = _ellipsoid.cartesianToCartographic(Cartesian3);
		let _lat = Cesium.Math.toDegrees(_cartographic.latitude);
		let _lng = Cesium.Math.toDegrees(_cartographic.longitude);
		let _alt = _cartographic.height;
		return [_lng, _lat, _alt];
	}

	/**
	 * @ignore 忽略注释，注释不生成Doc
	 * @param {Object} lng1
	 * @param {Object} lat1
	 * @param {Object} lng2
	 * @param {Object} lat2
	 */
	getAngle(lng1, lat1, lng2, lat2) {
		let dRotateAngle = Math.atan2(Math.abs(lng1 - lng2), Math.abs(lat1 - lat2));
		if (lng2 >= lng1) {
			dRotateAngle = lat2 < lat1 ? Math.PI - dRotateAngle : dRotateAngle;
		} else {
			dRotateAngle =
				lat2 >= lat1 ? 2 * Math.PI - dRotateAngle : Math.PI + dRotateAngle;
		}
		dRotateAngle = (dRotateAngle * 180) / Math.PI;
		return dRotateAngle;
	}

	/**
	 * @ignore 忽略注释，注释不生成Doc
	 * @param {Object} pointA
	 * @param {Object} pointB
	 */
	getPitch(pointA, pointB) {
		let transfrom = Cesium.Transforms.eastNorthUpToFixedFrame(pointA);
		const vector = Cesium.Cartesian3.subtract(
			pointB,
			pointA,
			new Cesium.Cartesian3()
		);
		let direction = Cesium.Matrix4.multiplyByPointAsVector(
			Cesium.Matrix4.inverse(transfrom, transfrom),
			vector,
			vector
		);
		Cesium.Cartesian3.normalize(direction, direction);
		return Cesium.Math.PI_OVER_TWO - Cesium.Math.acosClamped(direction.z);
	}

	/**
	 * @ignore 忽略注释，注释不生成Doc
	 */
	updateViewShed() {
		this.clear();
		this.setLightCamera();
		this.addVisualPyramid();
		this.createFrustum();
	}

	/**
	 * @ignore 忽略注释，注释不生成Doc
	 */
	clear() {
		if (this.pyramid) {
			this.viewer.entities.removeById(this.pyramid.id);
			this.pyramid = undefined;
		}
		if (this.frustumPrimitive) {
			this.viewer.scene.primitives.remove(this.frustumPrimitive);
			this.frustumPrimitive = undefined;
		}
		if (this.debugModelMatrixPrimitive) {
			this.viewer.scene.primitives.remove(this.debugModelMatrixPrimitive);
			this.debugModelMatrixPrimitive = undefined;
		}
	}

	/**
	 * @ignore 忽略注释，注释不生成Doc
	 */
	addVisualPyramid() {
		let options = this.ViewShedOptions;
		let position = options.viewPosition;
		let visualRange = Number(options.visualRange);
		let transform = Cesium.Transforms.eastNorthUpToFixedFrame(position);
		this.debugModelMatrixPrimitive = this.viewer.scene.primitives.add(
			new Cesium.DebugModelMatrixPrimitive({
				modelMatrix: transform,
				length: 5.0,
			})
		);
		const halfClock = options.horizontalViewAngle / 2;
		const halfCone = options.verticalViewAngle / 2;
		const pitch = Cesium.Math.toDegrees(options.pitch);
		const ellipsoid = new Cesium.EllipsoidGraphics({
			radii: new Cesium.Cartesian3(visualRange, visualRange, visualRange),
			minimumClock: Cesium.Math.toRadians(90 - options.direction - halfClock),
			maximumClock: Cesium.Math.toRadians(90 - options.direction + halfClock),
			minimumCone: Cesium.Math.toRadians(90 - pitch - halfCone),
			maximumCone: Cesium.Math.toRadians(90 - pitch + halfCone),
			fill: false,
			outline: true,
			subdivisions: 256,
			stackPartitions: 64,
			slicePartitions: 64,
			outlineColor: Cesium.Color.YELLOWGREEN.withAlpha(0.5),
		});
		const pyramidEntity = new Cesium.Entity({
			position: position,
			ellipsoid,
		});
		this.pyramid = this.viewer.entities.add(pyramidEntity);
	}

	/**
	 * @ignore 忽略注释，注释不生成Doc
	 */
	setLightCamera() {
		if (!this.lightCamera) {
			this.lightCamera = new Cesium.Camera(this.viewer.scene);
		}
		let options = this.ViewShedOptions;
		let visualRange = Number(options.visualRange);
		this.lightCamera.position = options.viewPosition;
		this.lightCamera.frustum.near = 0.1;
		this.lightCamera.frustum.far = visualRange;
		const hr = Cesium.Math.toRadians(options.horizontalViewAngle);
		const vr = Cesium.Math.toRadians(options.verticalViewAngle);
		this.lightCamera.frustum.aspectRatio =
			(visualRange * Math.tan(hr / 2) * 2) /
			(visualRange * Math.tan(vr / 2) * 2);
		this.lightCamera.frustum.fov = hr > vr ? hr : vr;
		this.lightCamera.setView({
			destination: options.viewPosition,
			orientation: {
				heading: Cesium.Math.toRadians(options.direction || 0),
				pitch: options.pitch || 0,
				roll: 0,
			},
		});
	}

	/**
	 * @ignore 忽略注释，注释不生成Doc
	 */
	createFrustum() {
		const scratchRight = new Cesium.Cartesian3();
		const scratchRotation = new Cesium.Matrix3();
		const scratchOrientation = new Cesium.Quaternion();
		const direction = this.lightCamera.directionWC;
		const up = this.lightCamera.upWC;
		let right = this.lightCamera.rightWC;
		right = Cesium.Cartesian3.negate(right, scratchRight);
		let rotation = scratchRotation;
		Cesium.Matrix3.setColumn(rotation, 0, right, rotation);
		Cesium.Matrix3.setColumn(rotation, 1, up, rotation);
		Cesium.Matrix3.setColumn(rotation, 2, direction, rotation);
		let orientation = Cesium.Quaternion.fromRotationMatrix(
			rotation,
			scratchOrientation
		);
		let instanceOutline = new Cesium.GeometryInstance({
			geometry: new Cesium.FrustumOutlineGeometry({
				frustum: this.lightCamera.frustum,
				origin: this.ViewShedOptions.viewPosition,
				orientation: orientation,
			}),
			id: "视椎体轮廓线" + Math.random().toString(36).substr(2),
			attributes: {
				color: Cesium.ColorGeometryInstanceAttribute.fromColor(
					new Cesium.Color(0.0, 1.0, 0.0, 1.0)
				),
				show: new Cesium.ShowGeometryInstanceAttribute(true),
			},
		});
		this.frustumPrimitive = this.viewer.scene.primitives.add(
			new Cesium.Primitive({
				geometryInstances: instanceOutline,
				appearance: new Cesium.PerInstanceColorAppearance({
					flat: true,
					translucent: false,
					closed: true,
				}),
			})
		);
	}

	/**
	 * @ignore 忽略注释，注释不生成Doc
	 */
	addViewershedPolygon(positionArr) {
		let polygon = new Cesium.PolygonGeometry({
			polygonHierarchy: new Cesium.PolygonHierarchy(
				Cesium.Cartesian3.fromDegreesArray(positionArr)
			),
			height: 0.0,
			extrudedHeight: 0.0,
			vertexFormat: Cesium.PerInstanceColorAppearance.VERTEX_FORMAT,
			stRotation: 0.0,
			ellipsoid: Cesium.Ellipsoid.WGS84,
			granularity: Cesium.Math.RADIANS_PER_DEGREE,
			perPositionHeight: false,
			closeTop: true,
			closeBottom: true,
			arcType: Cesium.ArcType.GEODESIC,
		});

		let polygonInstance = new Cesium.GeometryInstance({
			geometry: polygon,
			name: "ViewershedPolygon",
			attributes: {
				color: Cesium.ColorGeometryInstanceAttribute.fromColor(
					Cesium.Color.BLUE.withAlpha(0.6)
				),
				show: new Cesium.ShowGeometryInstanceAttribute(true),
			},
		});
		this.viewershedPolygon = this.viewer.scene.primitives.add(
			new Cesium.GroundPrimitive({
				geometryInstances: polygonInstance,
				appearance: new Cesium.EllipsoidSurfaceAppearance({
					aboveGround: true,
					material: new Cesium.Material({
						fabric: {
							type: "Image",
							uniforms: {
								image: this.returnImgae(),
							},
						},
					}),
				}),
			})
		);
	}

	/**
	 * @ignore 忽略注释，注释不生成Doc
	 */
	drawViewershed(precision) {
		const pos = this.cartesian3ToDegree(this.ViewShedOptions.viewPosition);
		const radius = this.ViewShedOptions.visualRange;
		const direction = this.ViewShedOptions.direction;
		let boundary = this.computeBoundaryOptions(pos, radius, direction);
		const bbox = boundary.bbox;
		let mask = turf.polygon([boundary.boundaryPoints]);
		const dis = this.ViewShedOptions.visualRange / (precision * 1000);
		let gridPoints = turf.pointGrid(bbox, dis, {
			mask: mask
		});
		console.log(this.ViewShedOptions.visualRange, precision, dis);

		let pointsResult = this.createTargetPoints(gridPoints, dis, pos);
		let variogram = kriging.train(
			pointsResult.values,
			pointsResult.lngs,
			pointsResult.lats,
			"exponential",
			0,
			100
		);
		let grid = kriging.grid([boundary.boundaryPoints], variogram, dis / 1000);
		const colors = [
			"#ff000080",
			"#ff000080",
			"#ff000080",
			"#ff000080",
			"#ff000080",
			"#ff000080",
			"#00ff0080",
			"#00ff0080",
			"#00ff0080",
			"#00ff0080",
			"#00ff0080",
			"#00ff0080",
		];

		this.canvasEle.width = 3840;
		this.canvasEle.height = 2160;
		kriging.plot(
			this.canvasEle,
			grid,
			[bbox[0], bbox[2]],
			[bbox[1], bbox[3]],
			colors
		);
		this.addViewershedPolygon(boundary.positionArr);
	}

	/**
	 * @ignore 忽略注释，注释不生成Doc
	 */
	computeBoundaryOptions(pos, radius, angle) {
		let Ea = 6378137; // 赤道半径
		let Eb = 6356725; // 极半径
		const lng = pos[0],
			lat = pos[1];
		const bbox = [lng, lat, lng, lat];
		let positionArr = [];
		let boundaryPoints = [];
		positionArr.push(lng, lat);
		boundaryPoints.push([lng, lat]);
		//正北是0°
		let start = angle + 45 > 360 ? angle - 45 - 360 : angle - 45;
		let end = start + 90;
		for (let i = start; i <= end; i++) {
			let dx = radius * Math.sin((i * Math.PI) / 180.0);
			let dy = radius * Math.cos((i * Math.PI) / 180.0);
			let ec = Eb + ((Ea - Eb) * (90.0 - lat)) / 90.0;
			let ed = ec * Math.cos((lat * Math.PI) / 180);
			let BJD = lng + ((dx / ed) * 180.0) / Math.PI;
			let BWD = lat + ((dy / ec) * 180.0) / Math.PI;
			positionArr.push(BJD, BWD);
			boundaryPoints.push([BJD, BWD]);
			this.refreshBBox(bbox, BJD, BWD);
		}
		boundaryPoints.push([lng, lat]);
		return {
			positionArr,
			boundaryPoints,
			bbox,
		};
	}

	/**
	 * 更新外围矩形 Bbox
	 * @ignore 忽略注释，注释不生成Doc
	 * @param {Array} result 外围矩形Bbox-[minX, minY, maxX, maxY]
	 * @param {Number} x 经度
	 * @param {Number} y 纬度
	 */
	refreshBBox(result, x, y) {
		result[0] = x < result[0] ? x : result[0];
		result[1] = y < result[1] ? y : result[1];
		result[2] = x > result[2] ? x : result[2];
		result[3] = y > result[3] ? y : result[3];
	}

	/**
	 * 插值点用射线判断通视性
	 * @ignore 忽略注释，注释不生成Doc
	 * @param {*} gridPoints 网格点
	 * @param {*} step 步长，可以理解成是精度
	 * @param {*} sourcePos 视域分析起点
	 * @returns kriging插值所需的参数对象{ values:[], lngs:[], lats:[]}
	 */
	createTargetPoints(gridPoints, step, sourcePos) {
		let positionArr = [];
		let objectsToExclude = [
			this.frustumPrimitive,
			this.pyramid,
			this.debugModelMatrixPrimitive,
		];
		let values = [],
			lngs = [],
			lats = [];
		let height = this.getHeight(sourcePos[0], sourcePos[1], objectsToExclude);
		positionArr.push({
			x: sourcePos[0],
			y: sourcePos[1],
			z: height,
		});
		let viewPoint = this.ViewShedOptions.viewPosition;
		for (let index = 0; index < gridPoints.features.length; index++) {
			const feature = gridPoints.features[index];
			const coords = feature.geometry.coordinates;
			const x = coords[0],
				y = coords[1];
			let h = this.getHeight(x, y, objectsToExclude);
			let endPoint = Cesium.Cartesian3.fromDegrees(x, y, h);
			let direction = Cesium.Cartesian3.normalize(
				Cesium.Cartesian3.subtract(
					endPoint,
					viewPoint,
					new Cesium.Cartesian3()
				),
				new Cesium.Cartesian3()
			);
			// 建立射线
			let ray = new Cesium.Ray(viewPoint, direction);
			let result = this.viewer.scene.pickFromRay(ray, objectsToExclude); // 计算交互点，返回第一个
			if (result) {
				let buffer = this.ReturnDistance(endPoint, result.position);
				// let M_color = Cesium.Color.GREEN;
				if (buffer > step) {
					// M_color = Cesium.Color.RED;
					values.push(0);
				} else {
					values.push(1);
				}
				lngs.push(x);
				lats.push(y);
			}
		}
		return {
			values,
			lngs,
			lats,
		};
	}

	/**
	 * @ignore 忽略注释，注释不生成Doc
	 * @returns base64图片
	 */
	returnImgae() {
		return this.canvasEle.toDataURL("image/png");
	}
}

/**
 * 通用对外公开函数 
 */
Object.assign(ViewShed.prototype, /** @lends ViewShed.prototype */ {
	/**
	 * 开始执行视域分析
	 * @param {number} precision 精度，值越大创建耗时越长，建议在10~20之间
	 */
	createViewshed: function(precision) {
		let _self = this;
		let scene = _self.viewer.scene;
		_self.initHandler();
		_self.clearAll();

		let count = 0;
		let toolTip = "左键点击创建起点";
		_self.handler = new Cesium.ScreenSpaceEventHandler(_self.viewer.canvas);
		_self.handler.setInputAction((event) => {
			count++;

			if (count === 1) {
				toolTip = "左键点击创建终点";
				let earthPosition = scene.pickPosition(event.position);
				let pos = _self.cartesian3ToDegree(earthPosition);
				_self.handler.setInputAction(function(move) {
					CreateRemindertip(toolTip, move.endPosition, true);
					let newPosition = scene.pickPosition(move.endPosition);
					if (Cesium.defined(newPosition)) {
						let pos1 = _self.cartesian3ToDegree(newPosition);
						let distance = Cesium.Cartesian3.distance(newPosition, earthPosition);
						let angle = _self.getAngle(pos[0], pos[1], pos1[0], pos1[1]);
						let pitch = _self.getPitch(earthPosition, newPosition);
						_self.ViewShedOptions = {
							viewPosition: earthPosition,
							endPosition: newPosition,
							direction: angle,
							pitch: pitch,
							horizontalViewAngle: 90,
							verticalViewAngle: 60,
							visibleAreaColor: Cesium.Color.GREEN,
							invisibleAreaColor: Cesium.Color.RED,
							visualRange: distance,
						};
						_self.updateViewShed();
					}
				}, Cesium.ScreenSpaceEventType.MOUSE_MOVE);
			}

			if (count === 2) {
				_self.initHandler();
				_self.drawViewershed(precision);
				CreateRemindertip(toolTip, event.endPosition, false);
			}

		}, Cesium.ScreenSpaceEventType.LEFT_CLICK);

		_self.handler.setInputAction(function(move) {
			CreateRemindertip(toolTip, move.endPosition, true);
		}, Cesium.ScreenSpaceEventType.MOUSE_MOVE);
	},

	/**
	 * 清除视域分析
	 */
	clearAll: function() {
		this.clear();
		if (this.viewershedPolygon) {
			this.viewer.scene.primitives.remove(this.viewershedPolygon);
			this.viewershedPolygon = undefined;
		}
	},
});

export default ViewShed;