HeightmapTessellator.js 20 KB

123456789101112131415161718192021222324252627282930313233343536373839404142434445464748495051525354555657585960616263646566676869707172737475767778798081828384858687888990919293949596979899100101102103104105106107108109110111112113114115116117118119120121122123124125126127128129130131132133134135136137138139140141142143144145146147148149150151152153154155156157158159160161162163164165166167168169170171172173174175176177178179180181182183184185186187188189190191192193194195196197198199200201202203204205206207208209210211212213214215216217218219220221222223224225226227228229230231232233234235236237238239240241242243244245246247248249250251252253254255256257258259260261262263264265266267268269270271272273274275276277278279280281282283284285286287288289290291292293294295296297298299300301302303304305306307308309310311312313314315316317318319320321322323324325326327328329330331332333334335336337338339340341342343344345346347348349350351352353354355356357358359360361362363364365366367368369370371372373374375376377378379380381382383384385386387388389390391392393394395396397398399400401402403404405406407408409410411412413414415416417418419420421422423424425426427428429430431432433434435436437438439440441442443444445446447448449450451452453454455456457458459460461462463464465466467468469470471472473474475476477478479480481482483484485486487488489490491492493494495496497498499500501502503504505506507508509510511512513514515516517518519520521522523524525526527528529530531532533534535
  1. import AxisAlignedBoundingBox from "./AxisAlignedBoundingBox.js";
  2. import BoundingSphere from "./BoundingSphere.js";
  3. import Cartesian2 from "./Cartesian2.js";
  4. import Cartesian3 from "./Cartesian3.js";
  5. import defaultValue from "./defaultValue.js";
  6. import defined from "./defined.js";
  7. import DeveloperError from "./DeveloperError.js";
  8. import Ellipsoid from "./Ellipsoid.js";
  9. import EllipsoidalOccluder from "./EllipsoidalOccluder.js";
  10. import CesiumMath from "./Math.js";
  11. import Matrix4 from "./Matrix4.js";
  12. import OrientedBoundingBox from "./OrientedBoundingBox.js";
  13. import Rectangle from "./Rectangle.js";
  14. import TerrainEncoding from "./TerrainEncoding.js";
  15. import Transforms from "./Transforms.js";
  16. import WebMercatorProjection from "./WebMercatorProjection.js";
  17. /**
  18. * Contains functions to create a mesh from a heightmap image.
  19. *
  20. * @namespace HeightmapTessellator
  21. *
  22. * @private
  23. */
  24. const HeightmapTessellator = {};
  25. /**
  26. * The default structure of a heightmap, as given to {@link HeightmapTessellator.computeVertices}.
  27. *
  28. * @constant
  29. */
  30. HeightmapTessellator.DEFAULT_STRUCTURE = Object.freeze({
  31. heightScale: 1.0,
  32. heightOffset: 0.0,
  33. elementsPerHeight: 1,
  34. stride: 1,
  35. elementMultiplier: 256.0,
  36. isBigEndian: false,
  37. });
  38. const cartesian3Scratch = new Cartesian3();
  39. const matrix4Scratch = new Matrix4();
  40. const minimumScratch = new Cartesian3();
  41. const maximumScratch = new Cartesian3();
  42. /**
  43. * Fills an array of vertices from a heightmap image.
  44. *
  45. * @param {Object} options Object with the following properties:
  46. * @param {Int8Array|Uint8Array|Int16Array|Uint16Array|Int32Array|Uint32Array|Float32Array|Float64Array} options.heightmap The heightmap to tessellate.
  47. * @param {Number} options.width The width of the heightmap, in height samples.
  48. * @param {Number} options.height The height of the heightmap, in height samples.
  49. * @param {Number} options.skirtHeight The height of skirts to drape at the edges of the heightmap.
  50. * @param {Rectangle} options.nativeRectangle A rectangle in the native coordinates of the heightmap's projection. For
  51. * a heightmap with a geographic projection, this is degrees. For the web mercator
  52. * projection, this is meters.
  53. * @param {Number} [options.exaggeration=1.0] The scale used to exaggerate the terrain.
  54. * @param {Number} [options.exaggerationRelativeHeight=0.0] The height from which terrain is exaggerated.
  55. * @param {Rectangle} [options.rectangle] The rectangle covered by the heightmap, in geodetic coordinates with north, south, east and
  56. * west properties in radians. Either rectangle or nativeRectangle must be provided. If both
  57. * are provided, they're assumed to be consistent.
  58. * @param {Boolean} [options.isGeographic=true] True if the heightmap uses a {@link GeographicProjection}, or false if it uses
  59. * a {@link WebMercatorProjection}.
  60. * @param {Cartesian3} [options.relativeToCenter=Cartesian3.ZERO] The positions will be computed as <code>Cartesian3.subtract(worldPosition, relativeToCenter)</code>.
  61. * @param {Ellipsoid} [options.ellipsoid=Ellipsoid.WGS84] The ellipsoid to which the heightmap applies.
  62. * @param {Object} [options.structure] An object describing the structure of the height data.
  63. * @param {Number} [options.structure.heightScale=1.0] The factor by which to multiply height samples in order to obtain
  64. * the height above the heightOffset, in meters. The heightOffset is added to the resulting
  65. * height after multiplying by the scale.
  66. * @param {Number} [options.structure.heightOffset=0.0] The offset to add to the scaled height to obtain the final
  67. * height in meters. The offset is added after the height sample is multiplied by the
  68. * heightScale.
  69. * @param {Number} [options.structure.elementsPerHeight=1] The number of elements in the buffer that make up a single height
  70. * sample. This is usually 1, indicating that each element is a separate height sample. If
  71. * it is greater than 1, that number of elements together form the height sample, which is
  72. * computed according to the structure.elementMultiplier and structure.isBigEndian properties.
  73. * @param {Number} [options.structure.stride=1] The number of elements to skip to get from the first element of
  74. * one height to the first element of the next height.
  75. * @param {Number} [options.structure.elementMultiplier=256.0] The multiplier used to compute the height value when the
  76. * stride property is greater than 1. For example, if the stride is 4 and the strideMultiplier
  77. * is 256, the height is computed as follows:
  78. * `height = buffer[index] + buffer[index + 1] * 256 + buffer[index + 2] * 256 * 256 + buffer[index + 3] * 256 * 256 * 256`
  79. * This is assuming that the isBigEndian property is false. If it is true, the order of the
  80. * elements is reversed.
  81. * @param {Number} [options.structure.lowestEncodedHeight] The lowest value that can be stored in the height buffer. Any heights that are lower
  82. * than this value after encoding with the `heightScale` and `heightOffset` are clamped to this value. For example, if the height
  83. * buffer is a `Uint16Array`, this value should be 0 because a `Uint16Array` cannot store negative numbers. If this parameter is
  84. * not specified, no minimum value is enforced.
  85. * @param {Number} [options.structure.highestEncodedHeight] The highest value that can be stored in the height buffer. Any heights that are higher
  86. * than this value after encoding with the `heightScale` and `heightOffset` are clamped to this value. For example, if the height
  87. * buffer is a `Uint16Array`, this value should be `256 * 256 - 1` or 65535 because a `Uint16Array` cannot store numbers larger
  88. * than 65535. If this parameter is not specified, no maximum value is enforced.
  89. * @param {Boolean} [options.structure.isBigEndian=false] Indicates endianness of the elements in the buffer when the
  90. * stride property is greater than 1. If this property is false, the first element is the
  91. * low-order element. If it is true, the first element is the high-order element.
  92. *
  93. * @example
  94. * const width = 5;
  95. * const height = 5;
  96. * const statistics = Cesium.HeightmapTessellator.computeVertices({
  97. * heightmap : [1.0, 2.0, 3.0, 4.0, 5.0, 6.0, 7.0, 8.0, 9.0],
  98. * width : width,
  99. * height : height,
  100. * skirtHeight : 0.0,
  101. * nativeRectangle : {
  102. * west : 10.0,
  103. * east : 20.0,
  104. * south : 30.0,
  105. * north : 40.0
  106. * }
  107. * });
  108. *
  109. * const encoding = statistics.encoding;
  110. * const position = encoding.decodePosition(statistics.vertices, index);
  111. */
  112. HeightmapTessellator.computeVertices = function (options) {
  113. //>>includeStart('debug', pragmas.debug);
  114. if (!defined(options) || !defined(options.heightmap)) {
  115. throw new DeveloperError("options.heightmap is required.");
  116. }
  117. if (!defined(options.width) || !defined(options.height)) {
  118. throw new DeveloperError("options.width and options.height are required.");
  119. }
  120. if (!defined(options.nativeRectangle)) {
  121. throw new DeveloperError("options.nativeRectangle is required.");
  122. }
  123. if (!defined(options.skirtHeight)) {
  124. throw new DeveloperError("options.skirtHeight is required.");
  125. }
  126. //>>includeEnd('debug');
  127. // This function tends to be a performance hotspot for terrain rendering,
  128. // so it employs a lot of inlining and unrolling as an optimization.
  129. // In particular, the functionality of Ellipsoid.cartographicToCartesian
  130. // is inlined.
  131. const cos = Math.cos;
  132. const sin = Math.sin;
  133. const sqrt = Math.sqrt;
  134. const atan = Math.atan;
  135. const exp = Math.exp;
  136. const piOverTwo = CesiumMath.PI_OVER_TWO;
  137. const toRadians = CesiumMath.toRadians;
  138. const heightmap = options.heightmap;
  139. const width = options.width;
  140. const height = options.height;
  141. const skirtHeight = options.skirtHeight;
  142. const hasSkirts = skirtHeight > 0.0;
  143. const isGeographic = defaultValue(options.isGeographic, true);
  144. const ellipsoid = defaultValue(options.ellipsoid, Ellipsoid.WGS84);
  145. const oneOverGlobeSemimajorAxis = 1.0 / ellipsoid.maximumRadius;
  146. const nativeRectangle = Rectangle.clone(options.nativeRectangle);
  147. const rectangle = Rectangle.clone(options.rectangle);
  148. let geographicWest;
  149. let geographicSouth;
  150. let geographicEast;
  151. let geographicNorth;
  152. if (!defined(rectangle)) {
  153. if (isGeographic) {
  154. geographicWest = toRadians(nativeRectangle.west);
  155. geographicSouth = toRadians(nativeRectangle.south);
  156. geographicEast = toRadians(nativeRectangle.east);
  157. geographicNorth = toRadians(nativeRectangle.north);
  158. } else {
  159. geographicWest = nativeRectangle.west * oneOverGlobeSemimajorAxis;
  160. geographicSouth =
  161. piOverTwo -
  162. 2.0 * atan(exp(-nativeRectangle.south * oneOverGlobeSemimajorAxis));
  163. geographicEast = nativeRectangle.east * oneOverGlobeSemimajorAxis;
  164. geographicNorth =
  165. piOverTwo -
  166. 2.0 * atan(exp(-nativeRectangle.north * oneOverGlobeSemimajorAxis));
  167. }
  168. } else {
  169. geographicWest = rectangle.west;
  170. geographicSouth = rectangle.south;
  171. geographicEast = rectangle.east;
  172. geographicNorth = rectangle.north;
  173. }
  174. let relativeToCenter = options.relativeToCenter;
  175. const hasRelativeToCenter = defined(relativeToCenter);
  176. relativeToCenter = hasRelativeToCenter ? relativeToCenter : Cartesian3.ZERO;
  177. const includeWebMercatorT = defaultValue(options.includeWebMercatorT, false);
  178. const exaggeration = defaultValue(options.exaggeration, 1.0);
  179. const exaggerationRelativeHeight = defaultValue(
  180. options.exaggerationRelativeHeight,
  181. 0.0
  182. );
  183. const hasExaggeration = exaggeration !== 1.0;
  184. const includeGeodeticSurfaceNormals = hasExaggeration;
  185. const structure = defaultValue(
  186. options.structure,
  187. HeightmapTessellator.DEFAULT_STRUCTURE
  188. );
  189. const heightScale = defaultValue(
  190. structure.heightScale,
  191. HeightmapTessellator.DEFAULT_STRUCTURE.heightScale
  192. );
  193. const heightOffset = defaultValue(
  194. structure.heightOffset,
  195. HeightmapTessellator.DEFAULT_STRUCTURE.heightOffset
  196. );
  197. const elementsPerHeight = defaultValue(
  198. structure.elementsPerHeight,
  199. HeightmapTessellator.DEFAULT_STRUCTURE.elementsPerHeight
  200. );
  201. const stride = defaultValue(
  202. structure.stride,
  203. HeightmapTessellator.DEFAULT_STRUCTURE.stride
  204. );
  205. const elementMultiplier = defaultValue(
  206. structure.elementMultiplier,
  207. HeightmapTessellator.DEFAULT_STRUCTURE.elementMultiplier
  208. );
  209. const isBigEndian = defaultValue(
  210. structure.isBigEndian,
  211. HeightmapTessellator.DEFAULT_STRUCTURE.isBigEndian
  212. );
  213. let rectangleWidth = Rectangle.computeWidth(nativeRectangle);
  214. let rectangleHeight = Rectangle.computeHeight(nativeRectangle);
  215. const granularityX = rectangleWidth / (width - 1);
  216. const granularityY = rectangleHeight / (height - 1);
  217. if (!isGeographic) {
  218. rectangleWidth *= oneOverGlobeSemimajorAxis;
  219. rectangleHeight *= oneOverGlobeSemimajorAxis;
  220. }
  221. const radiiSquared = ellipsoid.radiiSquared;
  222. const radiiSquaredX = radiiSquared.x;
  223. const radiiSquaredY = radiiSquared.y;
  224. const radiiSquaredZ = radiiSquared.z;
  225. let minimumHeight = 65536.0;
  226. let maximumHeight = -65536.0;
  227. const fromENU = Transforms.eastNorthUpToFixedFrame(
  228. relativeToCenter,
  229. ellipsoid
  230. );
  231. const toENU = Matrix4.inverseTransformation(fromENU, matrix4Scratch);
  232. let southMercatorY;
  233. let oneOverMercatorHeight;
  234. if (includeWebMercatorT) {
  235. southMercatorY = WebMercatorProjection.geodeticLatitudeToMercatorAngle(
  236. geographicSouth
  237. );
  238. oneOverMercatorHeight =
  239. 1.0 /
  240. (WebMercatorProjection.geodeticLatitudeToMercatorAngle(geographicNorth) -
  241. southMercatorY);
  242. }
  243. const minimum = minimumScratch;
  244. minimum.x = Number.POSITIVE_INFINITY;
  245. minimum.y = Number.POSITIVE_INFINITY;
  246. minimum.z = Number.POSITIVE_INFINITY;
  247. const maximum = maximumScratch;
  248. maximum.x = Number.NEGATIVE_INFINITY;
  249. maximum.y = Number.NEGATIVE_INFINITY;
  250. maximum.z = Number.NEGATIVE_INFINITY;
  251. let hMin = Number.POSITIVE_INFINITY;
  252. const gridVertexCount = width * height;
  253. const edgeVertexCount = skirtHeight > 0.0 ? width * 2 + height * 2 : 0;
  254. const vertexCount = gridVertexCount + edgeVertexCount;
  255. const positions = new Array(vertexCount);
  256. const heights = new Array(vertexCount);
  257. const uvs = new Array(vertexCount);
  258. const webMercatorTs = includeWebMercatorT ? new Array(vertexCount) : [];
  259. const geodeticSurfaceNormals = includeGeodeticSurfaceNormals
  260. ? new Array(vertexCount)
  261. : [];
  262. let startRow = 0;
  263. let endRow = height;
  264. let startCol = 0;
  265. let endCol = width;
  266. if (hasSkirts) {
  267. --startRow;
  268. ++endRow;
  269. --startCol;
  270. ++endCol;
  271. }
  272. const skirtOffsetPercentage = 0.00001;
  273. for (let rowIndex = startRow; rowIndex < endRow; ++rowIndex) {
  274. let row = rowIndex;
  275. if (row < 0) {
  276. row = 0;
  277. }
  278. if (row >= height) {
  279. row = height - 1;
  280. }
  281. let latitude = nativeRectangle.north - granularityY * row;
  282. if (!isGeographic) {
  283. latitude =
  284. piOverTwo - 2.0 * atan(exp(-latitude * oneOverGlobeSemimajorAxis));
  285. } else {
  286. latitude = toRadians(latitude);
  287. }
  288. let v = (latitude - geographicSouth) / (geographicNorth - geographicSouth);
  289. v = CesiumMath.clamp(v, 0.0, 1.0);
  290. const isNorthEdge = rowIndex === startRow;
  291. const isSouthEdge = rowIndex === endRow - 1;
  292. if (skirtHeight > 0.0) {
  293. if (isNorthEdge) {
  294. latitude += skirtOffsetPercentage * rectangleHeight;
  295. } else if (isSouthEdge) {
  296. latitude -= skirtOffsetPercentage * rectangleHeight;
  297. }
  298. }
  299. const cosLatitude = cos(latitude);
  300. const nZ = sin(latitude);
  301. const kZ = radiiSquaredZ * nZ;
  302. let webMercatorT;
  303. if (includeWebMercatorT) {
  304. webMercatorT =
  305. (WebMercatorProjection.geodeticLatitudeToMercatorAngle(latitude) -
  306. southMercatorY) *
  307. oneOverMercatorHeight;
  308. }
  309. for (let colIndex = startCol; colIndex < endCol; ++colIndex) {
  310. let col = colIndex;
  311. if (col < 0) {
  312. col = 0;
  313. }
  314. if (col >= width) {
  315. col = width - 1;
  316. }
  317. const terrainOffset = row * (width * stride) + col * stride;
  318. let heightSample;
  319. if (elementsPerHeight === 1) {
  320. heightSample = heightmap[terrainOffset];
  321. } else {
  322. heightSample = 0;
  323. let elementOffset;
  324. if (isBigEndian) {
  325. for (
  326. elementOffset = 0;
  327. elementOffset < elementsPerHeight;
  328. ++elementOffset
  329. ) {
  330. heightSample =
  331. heightSample * elementMultiplier +
  332. heightmap[terrainOffset + elementOffset];
  333. }
  334. } else {
  335. for (
  336. elementOffset = elementsPerHeight - 1;
  337. elementOffset >= 0;
  338. --elementOffset
  339. ) {
  340. heightSample =
  341. heightSample * elementMultiplier +
  342. heightmap[terrainOffset + elementOffset];
  343. }
  344. }
  345. }
  346. heightSample = heightSample * heightScale + heightOffset;
  347. maximumHeight = Math.max(maximumHeight, heightSample);
  348. minimumHeight = Math.min(minimumHeight, heightSample);
  349. let longitude = nativeRectangle.west + granularityX * col;
  350. if (!isGeographic) {
  351. longitude = longitude * oneOverGlobeSemimajorAxis;
  352. } else {
  353. longitude = toRadians(longitude);
  354. }
  355. let u = (longitude - geographicWest) / (geographicEast - geographicWest);
  356. u = CesiumMath.clamp(u, 0.0, 1.0);
  357. let index = row * width + col;
  358. if (skirtHeight > 0.0) {
  359. const isWestEdge = colIndex === startCol;
  360. const isEastEdge = colIndex === endCol - 1;
  361. const isEdge = isNorthEdge || isSouthEdge || isWestEdge || isEastEdge;
  362. const isCorner =
  363. (isNorthEdge || isSouthEdge) && (isWestEdge || isEastEdge);
  364. if (isCorner) {
  365. // Don't generate skirts on the corners.
  366. continue;
  367. } else if (isEdge) {
  368. heightSample -= skirtHeight;
  369. if (isWestEdge) {
  370. // The outer loop iterates north to south but the indices are ordered south to north, hence the index flip below
  371. index = gridVertexCount + (height - row - 1);
  372. longitude -= skirtOffsetPercentage * rectangleWidth;
  373. } else if (isSouthEdge) {
  374. // Add after west indices. South indices are ordered east to west.
  375. index = gridVertexCount + height + (width - col - 1);
  376. } else if (isEastEdge) {
  377. // Add after west and south indices. East indices are ordered north to south. The index is flipped like above.
  378. index = gridVertexCount + height + width + row;
  379. longitude += skirtOffsetPercentage * rectangleWidth;
  380. } else if (isNorthEdge) {
  381. // Add after west, south, and east indices. North indices are ordered west to east.
  382. index = gridVertexCount + height + width + height + col;
  383. }
  384. }
  385. }
  386. const nX = cosLatitude * cos(longitude);
  387. const nY = cosLatitude * sin(longitude);
  388. const kX = radiiSquaredX * nX;
  389. const kY = radiiSquaredY * nY;
  390. const gamma = sqrt(kX * nX + kY * nY + kZ * nZ);
  391. const oneOverGamma = 1.0 / gamma;
  392. const rSurfaceX = kX * oneOverGamma;
  393. const rSurfaceY = kY * oneOverGamma;
  394. const rSurfaceZ = kZ * oneOverGamma;
  395. const position = new Cartesian3();
  396. position.x = rSurfaceX + nX * heightSample;
  397. position.y = rSurfaceY + nY * heightSample;
  398. position.z = rSurfaceZ + nZ * heightSample;
  399. Matrix4.multiplyByPoint(toENU, position, cartesian3Scratch);
  400. Cartesian3.minimumByComponent(cartesian3Scratch, minimum, minimum);
  401. Cartesian3.maximumByComponent(cartesian3Scratch, maximum, maximum);
  402. hMin = Math.min(hMin, heightSample);
  403. positions[index] = position;
  404. uvs[index] = new Cartesian2(u, v);
  405. heights[index] = heightSample;
  406. if (includeWebMercatorT) {
  407. webMercatorTs[index] = webMercatorT;
  408. }
  409. if (includeGeodeticSurfaceNormals) {
  410. geodeticSurfaceNormals[index] = ellipsoid.geodeticSurfaceNormal(
  411. position
  412. );
  413. }
  414. }
  415. }
  416. const boundingSphere3D = BoundingSphere.fromPoints(positions);
  417. let orientedBoundingBox;
  418. if (defined(rectangle)) {
  419. orientedBoundingBox = OrientedBoundingBox.fromRectangle(
  420. rectangle,
  421. minimumHeight,
  422. maximumHeight,
  423. ellipsoid
  424. );
  425. }
  426. let occludeePointInScaledSpace;
  427. if (hasRelativeToCenter) {
  428. const occluder = new EllipsoidalOccluder(ellipsoid);
  429. occludeePointInScaledSpace = occluder.computeHorizonCullingPointPossiblyUnderEllipsoid(
  430. relativeToCenter,
  431. positions,
  432. minimumHeight
  433. );
  434. }
  435. const aaBox = new AxisAlignedBoundingBox(minimum, maximum, relativeToCenter);
  436. const encoding = new TerrainEncoding(
  437. relativeToCenter,
  438. aaBox,
  439. hMin,
  440. maximumHeight,
  441. fromENU,
  442. false,
  443. includeWebMercatorT,
  444. includeGeodeticSurfaceNormals,
  445. exaggeration,
  446. exaggerationRelativeHeight
  447. );
  448. const vertices = new Float32Array(vertexCount * encoding.stride);
  449. let bufferIndex = 0;
  450. for (let j = 0; j < vertexCount; ++j) {
  451. bufferIndex = encoding.encode(
  452. vertices,
  453. bufferIndex,
  454. positions[j],
  455. uvs[j],
  456. heights[j],
  457. undefined,
  458. webMercatorTs[j],
  459. geodeticSurfaceNormals[j]
  460. );
  461. }
  462. return {
  463. vertices: vertices,
  464. maximumHeight: maximumHeight,
  465. minimumHeight: minimumHeight,
  466. encoding: encoding,
  467. boundingSphere3D: boundingSphere3D,
  468. orientedBoundingBox: orientedBoundingBox,
  469. occludeePointInScaledSpace: occludeePointInScaledSpace,
  470. };
  471. };
  472. export default HeightmapTessellator;