fogleman / hmm

Heightmap meshing utility.

hmm

hmm is a heightmap meshing utility.

If you've done any 3D game development, 3D printing, or other such things, you've likely wanted to convert a grayscale heightmap image into a 3D mesh. The naive way is pretty simple but generates huge meshes with millions of triangles. After hacking my way through various solutions over the years, I finally decided I needed to write a good tool for this purpose.

hmm is a modern implementation of a nice algorithm from the 1995 paper Fast Polygonal Approximation of Terrains and Height Fields by Garland and Heckbert. The meshes produced by hmm satisfy the Delaunay condition and can satisfy a specified maximal error or maximal number of triangles or vertices. It's also very fast.

Dependencies

• C++11 or higher
• glm

Installation

brew install glm # on macOS
sudo apt-get install libglm-dev # on Ubuntu / Debian

git clone https://github.com/fogleman/hmm.git
cd hmm
make
make install

Usage

heightmap meshing utility
usage: hmm --zscale=float [options] ... infile outfile.stl
options:
  -z, --zscale           z scale relative to x & y (float)
  -x, --zexagg           z exaggeration (float [=1])
  -e, --error            maximum triangulation error (float [=0.001])
  -t, --triangles        maximum number of triangles (int [=0])
  -p, --points           maximum number of vertices (int [=0])
  -b, --base             solid base height (float [=0])
      --level            auto level input to full grayscale range
      --invert           invert heightmap
      --blur             gaussian blur sigma (int [=0])
      --gamma            gamma curve exponent (float [=0])
      --border-size      border size in pixels (int [=0])
      --border-height    border z height (float [=1])
      --normal-map       path to write normal map png (string [=])
      --shade-path       path to write hillshade png (string [=])
      --shade-alt        hillshade light altitude (float [=45])
      --shade-az         hillshade light azimuth (float [=0])
  -q, --quiet            suppress console output
  -?, --help             print this message

hmm supports a variety of file formats like PNG, JPG, etc. for the input heightmap. The output is always a binary STL file. The only other required parameter is -z, which specifies how much to scale the Z axis in the output mesh.

$ hmm input.png output.stl -z ZSCALE

You can also provide a maximal allowed error, number of triangles, or number of vertices. (If multiple are specified, the first one reached is used.)

$ hmm input.png output.stl -z 100 -e 0.001 -t 1000000

Visual Guide

Click on the image below to see examples of various command line arguments. You can try these examples yourself with this heightmap: gale.png.

Z Scale

The required -z parameter defines the distance between a fully black pixel and a fully white pixel in the vertical Z axis, with units equal to one pixel width or height. For example, if each pixel in the heightmap represented a 1x1 meter square area, and the vertical range of the heightmap was 100 meters, then -z 100 should be used.

Z Exaggeration

The -x parameter is simply an extra multiplier on top of the provided Z scale. It is provided as a convenience so you don't have to do multiplication in your head just to exaggerate by, e.g. 2x, since Z scales are often derived from real world data and can have strange values like 142.2378.

Max Error

The -e parameter defines the maximum allowed error in the output mesh, as a percentage of the total mesh height. For example, if -e 0.01 is used, then no pixel will have an error of more than 1% of the distance between a fully black pixel and a fully white pixel. This means that for an 8-bit input image, an error of e = 1 / 256 ~= 0.0039 will ensure that no pixel has an error greater than one full grayscale unit. (It may still be desirable to use a lower value like 0.5 / 256.)

Base Height

When the -b option is used to create a solid mesh, it defines the height of the base before the lowest part of the heightmesh appears, as a percentage of the heightmap's height. For example, if -z 100 -b 0.5 were used, then the final mesh would be about 150 units tall (if a fully white pixel exists in the input).

Border

A border can be added to the mesh with the --border-size and --border-height flags. The heightmap will be padded by border-size pixels before triangulating. The (pre-scaled) Z value of the border can be set with border-height which defaults to 1.

Filters

A Gaussian blur can be applied with the --blur flag. This is particularly useful for noisy images.

The heightmap can be inverted with the --invert flag. This is useful for lithophanes.

The heightmap can be auto-leveled with the --level flag. This will stretch the grayscale values to use the entire black => white range.

A gamma curve can be applied to the heightmap with the --gamma flag. This applies x = x ^ gamma to each pixel, where x is in [0, 1].

Normal Maps

A full resolution normal map can be generated with the --normal-map argument. This will save a normal map as an RGB PNG to the specified path. This is useful for rendering higher resolution bumps and details while using a lower resolution triangle mesh.

Hillshade Images

A grayscale hillshade image can be generated with the --shade-path argument. The altitude and azimuth of the light source can be changed with the --shade-alt and --shade-az arguments, which default to 45 degrees in altitude and 0 degrees from north (up).

Performance

Performance depends a lot on the amount of detail in the heightmap, but here are some figures for an example heightmap of a 40x40 kilometer area centered on Mount Everest. Various heightmap resolutions and permitted max errors are shown. Times computed on a 2018 13" MacBook Pro (2.7 GHz Intel Core i7).

Runtime in Seconds

Image Size / Error	e=0.01	e=0.001	e=0.0005	e=0.0001
9490 x 9490 px (90.0 MP)	6.535	13.102	19.394	58.949
4745 x 4745 px (22.5 MP)	1.867	4.903	8.886	33.327
2373 x 2373 px (5.6 MP)	0.559	2.353	4.930	14.243
1187 x 1187 px (1.4 MP)	0.168	1.021	1.961	3.709

Number of Triangles Output

Image Size / Error	e=0.01	e=0.001	e=0.0005	e=0.0001
9490 x 9490 px (90.0 MP)	33,869	1,084,972	2,467,831	14,488,022
4745 x 4745 px (22.5 MP)	33,148	1,032,263	2,323,772	11,719,491
2373 x 2373 px (5.6 MP)	31,724	935,787	1,979,227	6,561,070
1187 x 1187 px (1.4 MP)	27,275	629,352	1,160,079	2,347,713

TODO

• reconstruct grayscale image?
• better error handling, especially for file I/O
• better overflow handling - what's the largest supported heightmap?
• mesh validation?

hmm's Issues

visual studio and cmake build

https://github.com/shelllee/hmm

fix: cmdline.h error on windows
add: cmake config
update: gitignone

there will be one pull request when fully tested

Faster error queue

Another small optimization that you could port back: mapbox/delatin@b69b2f6

Since queue comparisons happen much more often than other operations like swap/push/pop, we can optimize them by keeping all errors in the same order as the queue (instead of a lookup table like queue indices). This gave me a modest ~3-5% improvement on overall running time, and also makes the code slightly easier to understand.

support float32 images

would it be possible to support float32 input files, like geotiffs, without conversion into a 8bit image format first? Totally understandable if additional dependencies are not desired. I am sure a short python/bash script could automate the conversion and correct zscale...

why is zscale required?

Why is zscale a required argument? Is it nonlinear - if you scale before vs after does it make a difference?

I ask because I postprocess the mesh to scale x/y to something independent of heightmap resolution, and at that time I have to consider the hmm z scaling. It'd be easier for my use case to either pass an x/y scaling as well as z to hmm, or to not pass a z scaling to hmm and do it myself afterwards alongside x/y. But if I have to pass an accurate z scaling to hmm for the algorithm to work well, then that's a good reason for it to be the way it is.

Windows build

It looks like for windows you need to include <algorithm> in order to use e.g. std::min.

GLM_ENABLE_EXPERIMENTAL seems to be required

i get the following errors when trying to make

08:01 AM noon ∈ hmm (master) make                                                        2 ↵
Creating directories
Beginning release build
Compiling: src/base.cpp -> build/release/base.o
In file included from src/base.cpp:3:0:
/usr/include/glm/gtx/hash.hpp:16:3: error: #error "GLM: GLM_GTX_hash is an experimental extension and may change in the future. Use #define GLM_ENABLE_EXPERIMENTAL before including it, if you really want to use it."
 # error "GLM: GLM_GTX_hash is an experimental extension and may change in the future. Use #define GLM_ENABLE_EXPERIMENTAL before including it, if you really want to use it."
   ^~~~~
In file included from /usr/include/glm/gtx/hash.hpp:27:0,
                 from src/base.cpp:3:
/usr/include/glm/gtx/dual_quaternion.hpp:24:3: error: #error "GLM: GLM_GTX_dual_quaternion is an experimental extension and may change in the future. Use #define GLM_ENABLE_EXPERIMENTAL before including it, if you really want to use it."
 # error "GLM: GLM_GTX_dual_quaternion is an experimental extension and may change in the future. Use #define GLM_ENABLE_EXPERIMENTAL before including it, if you really want to use it."
   ^~~~~
Makefile:176: recipe for target 'build/release/base.o' failed
make[1]: *** [build/release/base.o] Error 1
Makefile:116: recipe for target 'release' failed
make: *** [release] Error 2

Adding a function to find the triangle containing a point

Based on what I have seen in your code I was wondering about the possibility of inserting the following functions into triangulator.cpp and *.h in order to extend your code to be able to get the location of the triangle containing a point. I do no have much of C++ experience hence why I am not making a pull request.

Here are the functions,

// Private Function to check if a point is inside a triangle using barycentric coordinates
bool PointInTriangle(const glm::ivec2 p, const glm::ivec2 a, const glm::ivec2 b, const glm::ivec2 c) const {
    auto sign = [](const glm::ivec2 p1, const glm::ivec2 p2, const glm::ivec2 p3) {
        return (p1.x - p3.x) * (p2.y - p3.y) - (p2.x - p3.x) * (p1.y - p3.y);
    };

    bool d1 = sign(p, a, b) < 0.0f;
    bool d2 = sign(p, b, c) < 0.0f;
    bool d3 = sign(p, c, a) < 0.0f;

    return ((d1 == d2) && (d2 == d3));
}

// Public  function to triangle contianing point pt
int locatePoint(const glm::ivec2 pt) const {
    for (int i = 0; i < m_Queue.size(); ++i) {
        const int t = m_Queue[i];
        const int e0 = t * 3;
        const int e1 = e0 + 1;
        const int e2 = e0 + 2;

        const int p0 = m_Triangles[e0];
        const int p1 = m_Triangles[e1];
        const int p2 = m_Triangles[e2];

        const glm::ivec2 a = m_Points[p0];
        const glm::ivec2 b = m_Points[p1];
        const glm::ivec2 c = m_Points[p2];

        if (PointInTriangle(pt, a, b, c)) {
            return t;
        }
    }
    return -1;  // Point is not inside any triangle
}

Thank you for considering these.

conda-forge package

FYI hmm is now available as a conda-forge package: https://github.com/conda-forge/cpp-hmm-feedstock
To install it:

$ conda install cpp-hmm -c conda-forge

@fogleman please let me know if you want me to add you as a maintainer of the feedstock.
Thanks for this great library!

Is there any plan to add mask support?

Hi, thanks for this wonderful utility.

And is there any plan to add mask support? This can be used to generate models of areas with non-rectangular boundaries.

Thanks!

Homebrew priority queue -> std::set?

I see that std::priority_queue doesn't have the required Remove method. Googling the stackoverflow to learn how people live with that in C++ land got this: people sometimes use std::set for that, if queue can't keep element multiple times. It has erase method to remove something. Underlying structure is then a binary tree instead of heap, allowing faster removal without repacking everything.

Hope this helps, feel free to ignore, any comment why is it bad idea if it is would be appreciated.

hmm/src/triangulator.cpp

Line 303 in 0e0fdaa

// priority queue functions

Round edges

This is a great tool that I was seeking for several days now!

What I want to do is to create lithophanes that fit into a frame. The frames have round edges. My first idea was to add black corners in the input image, but that comes out very messy.

Is there a chance that you would implement the possibility to define a radius for the edges of the picture itself, for the base and for the frame or give me some advice where to start implementing it? I am not really familiar with the creating of meshes.

License is missing

:)

Work with curved height maps?

Am I right that this only works with "flat" height maps: F(x, y) = z where x, y, z are cartesian coordinates? For example, there is SRTM DEM set where each tile represents a part of Earth surface, so it's curved: F(f, l) = h where f and l are latitude and longitude (angles). Is there any way this can be used for such height maps?

Optimize triangle data footprint

@fogleman Just something I noticed when porting this to JS — currently, triangle data is always appended to triangles / halfedges / candidates / errors / queueIndexes arrays, leaving a ton of redundant old triangles that flipped on legalization.

An easy fix is to add an optional argument to addTriangle which is an existing halfedge index, and use that instead of a new one if provided. Then you can pass that to the first addTriangle after every QueueRemove or QueuePop.

This reduces the size of the arrays above at the end of a run by ~4x (they should exactly match the number of output triangles), and also modestly improves overall running time because of better cache locality / less allocation.

Windows Build without GCC: cmdline.h depends on `#include <cxxabi.h>`

Awesome utility program and thanks for sharing it on Github!

Currently cannot build with MSVC due to "fatal error C1083: Cannot open include file: 'cxxabi.h': No such file or directory".

The root cause is the dependency on this abi::__cxa_demangle() in this code section:

static inline std::string demangle(const std::string &name)
{
  int status=0;
  char *p=abi::__cxa_demangle(name.c_str(), 0, 0, &status);
  std::string ret(p);
  free(p);
  return ret;
}

Also see: PDAL/PDAL#335 - similar issue, might offer similar solution.