HeightMap: add terrain generation methods (closes #195)

Add seven terrain generation methods wrapping libtcod heightmap functions:
- add_hill(center, radius, height): Add smooth hill
- dig_hill(center, radius, depth): Dig crater (use negative depth)
- add_voronoi(num_points, coefficients, seed): Voronoi-based features
- mid_point_displacement(roughness, seed): Diamond-square terrain
- rain_erosion(drops, erosion, sedimentation, seed): Erosion simulation
- dig_bezier(points, start_radius, end_radius, start_depth, end_depth): Carve paths
- smooth(iterations): Average neighboring cells

All methods return self for chaining. Includes 24 unit tests.

Note: dig_hill and dig_bezier use libtcod's "dig" semantics - use negative
depth values to actually dig below current terrain level.

Co-Authored-By: Claude Opus 4.5 <noreply@anthropic.com>

src/PyHeightMap.cpp

@@ -3,6 +3,8 @@
 #include "McRFPy_Doc.h"
 #include "PyPositionHelper.h"  // Standardized position argument parsing
 #include <sstream>
+#include <cstdlib>  // For random seed handling
+#include <ctime>    // For time-based seeds
 
 // Property definitions
 PyGetSetDef PyHeightMap::getsetters[] = {
@@ -145,6 +147,80 @@ PyMethodDef PyHeightMap::methods[] = {
          MCRF_DESC("Return NEW HeightMap with (1.0 - value) for each cell."),
          MCRF_RETURNS("HeightMap: New inverted HeightMap (original is unchanged)")
      )},
+    // Terrain generation methods (#195)
+    {"add_hill", (PyCFunction)PyHeightMap::add_hill, METH_VARARGS | METH_KEYWORDS,
+     MCRF_METHOD(HeightMap, add_hill,
+         MCRF_SIG("(center, radius: float, height: float)", "HeightMap"),
+         MCRF_DESC("Add a smooth hill at the specified position."),
+         MCRF_ARGS_START
+         MCRF_ARG("center", "Center position as (x, y) tuple, list, or Vector")
+         MCRF_ARG("radius", "Radius of the hill in cells")
+         MCRF_ARG("height", "Height of the hill peak")
+         MCRF_RETURNS("HeightMap: self, for method chaining")
+     )},
+    {"dig_hill", (PyCFunction)PyHeightMap::dig_hill, METH_VARARGS | METH_KEYWORDS,
+     MCRF_METHOD(HeightMap, dig_hill,
+         MCRF_SIG("(center, radius: float, depth: float)", "HeightMap"),
+         MCRF_DESC("Dig a smooth crater at the specified position. Use negative depth to dig below current terrain."),
+         MCRF_ARGS_START
+         MCRF_ARG("center", "Center position as (x, y) tuple, list, or Vector")
+         MCRF_ARG("radius", "Radius of the crater in cells")
+         MCRF_ARG("depth", "Target depth (use negative to dig below current values)")
+         MCRF_RETURNS("HeightMap: self, for method chaining")
+         MCRF_NOTE("Only modifies cells where current value exceeds target depth")
+     )},
+    {"add_voronoi", (PyCFunction)PyHeightMap::add_voronoi, METH_VARARGS | METH_KEYWORDS,
+     MCRF_METHOD(HeightMap, add_voronoi,
+         MCRF_SIG("(num_points: int, coefficients: tuple = (1.0, -0.5), seed: int = None)", "HeightMap"),
+         MCRF_DESC("Add Voronoi-based terrain features."),
+         MCRF_ARGS_START
+         MCRF_ARG("num_points", "Number of Voronoi seed points")
+         MCRF_ARG("coefficients", "Coefficients for distance calculations (default: (1.0, -0.5))")
+         MCRF_ARG("seed", "Random seed (None for random)")
+         MCRF_RETURNS("HeightMap: self, for method chaining")
+     )},
+    {"mid_point_displacement", (PyCFunction)PyHeightMap::mid_point_displacement, METH_VARARGS | METH_KEYWORDS,
+     MCRF_METHOD(HeightMap, mid_point_displacement,
+         MCRF_SIG("(roughness: float = 0.5, seed: int = None)", "HeightMap"),
+         MCRF_DESC("Generate terrain using midpoint displacement algorithm (diamond-square)."),
+         MCRF_ARGS_START
+         MCRF_ARG("roughness", "Controls terrain roughness (0.0-1.0, default 0.5)")
+         MCRF_ARG("seed", "Random seed (None for random)")
+         MCRF_RETURNS("HeightMap: self, for method chaining")
+         MCRF_NOTE("Works best with power-of-2+1 dimensions (e.g., 65x65, 129x129)")
+     )},
+    {"rain_erosion", (PyCFunction)PyHeightMap::rain_erosion, METH_VARARGS | METH_KEYWORDS,
+     MCRF_METHOD(HeightMap, rain_erosion,
+         MCRF_SIG("(drops: int, erosion: float = 0.1, sedimentation: float = 0.05, seed: int = None)", "HeightMap"),
+         MCRF_DESC("Simulate rain erosion on the terrain."),
+         MCRF_ARGS_START
+         MCRF_ARG("drops", "Number of rain drops to simulate")
+         MCRF_ARG("erosion", "Erosion coefficient (default 0.1)")
+         MCRF_ARG("sedimentation", "Sedimentation coefficient (default 0.05)")
+         MCRF_ARG("seed", "Random seed (None for random)")
+         MCRF_RETURNS("HeightMap: self, for method chaining")
+     )},
+    {"dig_bezier", (PyCFunction)PyHeightMap::dig_bezier, METH_VARARGS | METH_KEYWORDS,
+     MCRF_METHOD(HeightMap, dig_bezier,
+         MCRF_SIG("(points: tuple, start_radius: float, end_radius: float, start_depth: float, end_depth: float)", "HeightMap"),
+         MCRF_DESC("Carve a path along a cubic Bezier curve. Use negative depths to dig below current terrain."),
+         MCRF_ARGS_START
+         MCRF_ARG("points", "Four control points as ((x0,y0), (x1,y1), (x2,y2), (x3,y3))")
+         MCRF_ARG("start_radius", "Radius at start of path")
+         MCRF_ARG("end_radius", "Radius at end of path")
+         MCRF_ARG("start_depth", "Target depth at start (use negative to dig)")
+         MCRF_ARG("end_depth", "Target depth at end (use negative to dig)")
+         MCRF_RETURNS("HeightMap: self, for method chaining")
+         MCRF_NOTE("Only modifies cells where current value exceeds target depth")
+     )},
+    {"smooth", (PyCFunction)PyHeightMap::smooth, METH_VARARGS | METH_KEYWORDS,
+     MCRF_METHOD(HeightMap, smooth,
+         MCRF_SIG("(iterations: int = 1)", "HeightMap"),
+         MCRF_DESC("Smooth the heightmap by averaging neighboring cells."),
+         MCRF_ARGS_START
+         MCRF_ARG("iterations", "Number of smoothing passes (default 1)")
+         MCRF_RETURNS("HeightMap: self, for method chaining")
+     )},
     {NULL}
 };
 
@@ -758,3 +834,319 @@ PyObject* PyHeightMap::inverse(PyHeightMapObject* self, PyObject* Py_UNUSED(args
 
     return (PyObject*)result;
 }
+
+// Terrain generation methods (#195)
+
+// Helper: Create TCOD random generator with optional seed
+static TCOD_Random* CreateTCODRandom(PyObject* seed_obj)
+{
+    if (seed_obj == nullptr || seed_obj == Py_None) {
+        // Use default random - return nullptr to use libtcod's default
+        return nullptr;
+    }
+
+    if (!PyLong_Check(seed_obj)) {
+        PyErr_SetString(PyExc_TypeError, "seed must be an integer or None");
+        return nullptr;
+    }
+
+    uint32_t seed = (uint32_t)PyLong_AsUnsignedLong(seed_obj);
+    if (PyErr_Occurred()) {
+        return nullptr;
+    }
+
+    return TCOD_random_new_from_seed(TCOD_RNG_MT, seed);
+}
+
+// Method: add_hill(center, radius, height) -> HeightMap
+PyObject* PyHeightMap::add_hill(PyHeightMapObject* self, PyObject* args, PyObject* kwds)
+{
+    static const char* keywords[] = {"center", "radius", "height", nullptr};
+    PyObject* center_obj = nullptr;
+    float radius, height;
+
+    if (!PyArg_ParseTupleAndKeywords(args, kwds, "Off", const_cast<char**>(keywords),
+                                     &center_obj, &radius, &height)) {
+        return nullptr;
+    }
+
+    if (!self->heightmap) {
+        PyErr_SetString(PyExc_RuntimeError, "HeightMap not initialized");
+        return nullptr;
+    }
+
+    float cx, cy;
+    if (!PyPosition_FromObject(center_obj, &cx, &cy)) {
+        return nullptr;
+    }
+
+    TCOD_heightmap_add_hill(self->heightmap, cx, cy, radius, height);
+
+    Py_INCREF(self);
+    return (PyObject*)self;
+}
+
+// Method: dig_hill(center, radius, depth) -> HeightMap
+PyObject* PyHeightMap::dig_hill(PyHeightMapObject* self, PyObject* args, PyObject* kwds)
+{
+    static const char* keywords[] = {"center", "radius", "depth", nullptr};
+    PyObject* center_obj = nullptr;
+    float radius, depth;
+
+    if (!PyArg_ParseTupleAndKeywords(args, kwds, "Off", const_cast<char**>(keywords),
+                                     &center_obj, &radius, &depth)) {
+        return nullptr;
+    }
+
+    if (!self->heightmap) {
+        PyErr_SetString(PyExc_RuntimeError, "HeightMap not initialized");
+        return nullptr;
+    }
+
+    float cx, cy;
+    if (!PyPosition_FromObject(center_obj, &cx, &cy)) {
+        return nullptr;
+    }
+
+    TCOD_heightmap_dig_hill(self->heightmap, cx, cy, radius, depth);
+
+    Py_INCREF(self);
+    return (PyObject*)self;
+}
+
+// Method: add_voronoi(num_points, coefficients=(1.0, -0.5), seed=None) -> HeightMap
+PyObject* PyHeightMap::add_voronoi(PyHeightMapObject* self, PyObject* args, PyObject* kwds)
+{
+    static const char* keywords[] = {"num_points", "coefficients", "seed", nullptr};
+    int num_points;
+    PyObject* coef_obj = nullptr;
+    PyObject* seed_obj = nullptr;
+
+    if (!PyArg_ParseTupleAndKeywords(args, kwds, "i|OO", const_cast<char**>(keywords),
+                                     &num_points, &coef_obj, &seed_obj)) {
+        return nullptr;
+    }
+
+    if (!self->heightmap) {
+        PyErr_SetString(PyExc_RuntimeError, "HeightMap not initialized");
+        return nullptr;
+    }
+
+    if (num_points <= 0) {
+        PyErr_SetString(PyExc_ValueError, "num_points must be positive");
+        return nullptr;
+    }
+
+    // Parse coefficients - default to (1.0, -0.5)
+    std::vector<float> coef;
+    if (coef_obj == nullptr || coef_obj == Py_None) {
+        coef = {1.0f, -0.5f};
+    } else if (PyTuple_Check(coef_obj)) {
+        Py_ssize_t size = PyTuple_Size(coef_obj);
+        for (Py_ssize_t i = 0; i < size; i++) {
+            PyObject* item = PyTuple_GetItem(coef_obj, i);
+            if (PyFloat_Check(item)) {
+                coef.push_back((float)PyFloat_AsDouble(item));
+            } else if (PyLong_Check(item)) {
+                coef.push_back((float)PyLong_AsLong(item));
+            } else {
+                PyErr_SetString(PyExc_TypeError, "coefficients must be numeric");
+                return nullptr;
+            }
+        }
+    } else if (PyList_Check(coef_obj)) {
+        Py_ssize_t size = PyList_Size(coef_obj);
+        for (Py_ssize_t i = 0; i < size; i++) {
+            PyObject* item = PyList_GetItem(coef_obj, i);
+            if (PyFloat_Check(item)) {
+                coef.push_back((float)PyFloat_AsDouble(item));
+            } else if (PyLong_Check(item)) {
+                coef.push_back((float)PyLong_AsLong(item));
+            } else {
+                PyErr_SetString(PyExc_TypeError, "coefficients must be numeric");
+                return nullptr;
+            }
+        }
+    } else {
+        PyErr_SetString(PyExc_TypeError, "coefficients must be a tuple or list");
+        return nullptr;
+    }
+
+    if (coef.empty()) {
+        PyErr_SetString(PyExc_ValueError, "coefficients cannot be empty");
+        return nullptr;
+    }
+
+    // Create random generator if seed provided
+    TCOD_Random* rnd = CreateTCODRandom(seed_obj);
+    if (PyErr_Occurred()) {
+        return nullptr;
+    }
+
+    TCOD_heightmap_add_voronoi(self->heightmap, num_points, (int)coef.size(), coef.data(), rnd);
+
+    if (rnd) {
+        TCOD_random_delete(rnd);
+    }
+
+    Py_INCREF(self);
+    return (PyObject*)self;
+}
+
+// Method: mid_point_displacement(roughness=0.5, seed=None) -> HeightMap
+PyObject* PyHeightMap::mid_point_displacement(PyHeightMapObject* self, PyObject* args, PyObject* kwds)
+{
+    static const char* keywords[] = {"roughness", "seed", nullptr};
+    float roughness = 0.5f;
+    PyObject* seed_obj = nullptr;
+
+    if (!PyArg_ParseTupleAndKeywords(args, kwds, "|fO", const_cast<char**>(keywords),
+                                     &roughness, &seed_obj)) {
+        return nullptr;
+    }
+
+    if (!self->heightmap) {
+        PyErr_SetString(PyExc_RuntimeError, "HeightMap not initialized");
+        return nullptr;
+    }
+
+    // Create random generator if seed provided
+    TCOD_Random* rnd = CreateTCODRandom(seed_obj);
+    if (PyErr_Occurred()) {
+        return nullptr;
+    }
+
+    TCOD_heightmap_mid_point_displacement(self->heightmap, rnd, roughness);
+
+    if (rnd) {
+        TCOD_random_delete(rnd);
+    }
+
+    Py_INCREF(self);
+    return (PyObject*)self;
+}
+
+// Method: rain_erosion(drops, erosion=0.1, sedimentation=0.05, seed=None) -> HeightMap
+PyObject* PyHeightMap::rain_erosion(PyHeightMapObject* self, PyObject* args, PyObject* kwds)
+{
+    static const char* keywords[] = {"drops", "erosion", "sedimentation", "seed", nullptr};
+    int drops;
+    float erosion = 0.1f;
+    float sedimentation = 0.05f;
+    PyObject* seed_obj = nullptr;
+
+    if (!PyArg_ParseTupleAndKeywords(args, kwds, "i|ffO", const_cast<char**>(keywords),
+                                     &drops, &erosion, &sedimentation, &seed_obj)) {
+        return nullptr;
+    }
+
+    if (!self->heightmap) {
+        PyErr_SetString(PyExc_RuntimeError, "HeightMap not initialized");
+        return nullptr;
+    }
+
+    if (drops <= 0) {
+        PyErr_SetString(PyExc_ValueError, "drops must be positive");
+        return nullptr;
+    }
+
+    // Create random generator if seed provided
+    TCOD_Random* rnd = CreateTCODRandom(seed_obj);
+    if (PyErr_Occurred()) {
+        return nullptr;
+    }
+
+    TCOD_heightmap_rain_erosion(self->heightmap, drops, erosion, sedimentation, rnd);
+
+    if (rnd) {
+        TCOD_random_delete(rnd);
+    }
+
+    Py_INCREF(self);
+    return (PyObject*)self;
+}
+
+// Method: dig_bezier(points, start_radius, end_radius, start_depth, end_depth) -> HeightMap
+PyObject* PyHeightMap::dig_bezier(PyHeightMapObject* self, PyObject* args, PyObject* kwds)
+{
+    static const char* keywords[] = {"points", "start_radius", "end_radius", "start_depth", "end_depth", nullptr};
+    PyObject* points_obj = nullptr;
+    float start_radius, end_radius, start_depth, end_depth;
+
+    if (!PyArg_ParseTupleAndKeywords(args, kwds, "Offff", const_cast<char**>(keywords),
+                                     &points_obj, &start_radius, &end_radius, &start_depth, &end_depth)) {
+        return nullptr;
+    }
+
+    if (!self->heightmap) {
+        PyErr_SetString(PyExc_RuntimeError, "HeightMap not initialized");
+        return nullptr;
+    }
+
+    // Parse 4 control points
+    if (!PyTuple_Check(points_obj) && !PyList_Check(points_obj)) {
+        PyErr_SetString(PyExc_TypeError, "points must be a tuple or list of 4 control points");
+        return nullptr;
+    }
+
+    Py_ssize_t size = PyTuple_Check(points_obj) ? PyTuple_Size(points_obj) : PyList_Size(points_obj);
+    if (size != 4) {
+        PyErr_Format(PyExc_ValueError, "points must contain exactly 4 control points, got %zd", size);
+        return nullptr;
+    }
+
+    int px[4], py[4];
+    for (int i = 0; i < 4; i++) {
+        PyObject* point = PyTuple_Check(points_obj) ? PyTuple_GetItem(points_obj, i) : PyList_GetItem(points_obj, i);
+        int x, y;
+        if (!PyPosition_FromObjectInt(point, &x, &y)) {
+            PyErr_Format(PyExc_TypeError, "control point %d must be a (x, y) position", i);
+            return nullptr;
+        }
+        px[i] = x;
+        py[i] = y;
+    }
+
+    TCOD_heightmap_dig_bezier(self->heightmap, px, py, start_radius, start_depth, end_radius, end_depth);
+
+    Py_INCREF(self);
+    return (PyObject*)self;
+}
+
+// Method: smooth(iterations=1) -> HeightMap
+PyObject* PyHeightMap::smooth(PyHeightMapObject* self, PyObject* args, PyObject* kwds)
+{
+    static const char* keywords[] = {"iterations", nullptr};
+    int iterations = 1;
+
+    if (!PyArg_ParseTupleAndKeywords(args, kwds, "|i", const_cast<char**>(keywords),
+                                     &iterations)) {
+        return nullptr;
+    }
+
+    if (!self->heightmap) {
+        PyErr_SetString(PyExc_RuntimeError, "HeightMap not initialized");
+        return nullptr;
+    }
+
+    if (iterations <= 0) {
+        PyErr_SetString(PyExc_ValueError, "iterations must be positive");
+        return nullptr;
+    }
+
+    // 3x3 averaging kernel
+    static const int kernel_size = 9;
+    static const int dx[9] = {-1, 0, 1, -1, 0, 1, -1, 0, 1};
+    static const int dy[9] = {-1, -1, -1, 0, 0, 0, 1, 1, 1};
+    static const float weight[9] = {1.0f/9.0f, 1.0f/9.0f, 1.0f/9.0f,
+                                     1.0f/9.0f, 1.0f/9.0f, 1.0f/9.0f,
+                                     1.0f/9.0f, 1.0f/9.0f, 1.0f/9.0f};
+
+    for (int i = 0; i < iterations; i++) {
+        // Apply to all heights (minLevel=0, maxLevel=very high)
+        TCOD_heightmap_kernel_transform(self->heightmap, kernel_size, dx, dy, weight, 0.0f, 1000000.0f);
+    }
+
+    Py_INCREF(self);
+    return (PyObject*)self;
+}

src/PyHeightMap.h

@@ -45,6 +45,15 @@ public:
     static PyObject* threshold_binary(PyHeightMapObject* self, PyObject* args, PyObject* kwds);
     static PyObject* inverse(PyHeightMapObject* self, PyObject* Py_UNUSED(args));
 
+    // Terrain generation methods (#195) - mutate self, return self for chaining
+    static PyObject* add_hill(PyHeightMapObject* self, PyObject* args, PyObject* kwds);
+    static PyObject* dig_hill(PyHeightMapObject* self, PyObject* args, PyObject* kwds);
+    static PyObject* add_voronoi(PyHeightMapObject* self, PyObject* args, PyObject* kwds);
+    static PyObject* mid_point_displacement(PyHeightMapObject* self, PyObject* args, PyObject* kwds);
+    static PyObject* rain_erosion(PyHeightMapObject* self, PyObject* args, PyObject* kwds);
+    static PyObject* dig_bezier(PyHeightMapObject* self, PyObject* args, PyObject* kwds);
+    static PyObject* smooth(PyHeightMapObject* self, PyObject* args, PyObject* kwds);
+
     // Subscript support for hmap[x, y] syntax
     static PyObject* subscript(PyHeightMapObject* self, PyObject* key);
 

tests/unit/test_heightmap_terrain.py

@@ -0,0 +1,362 @@
+#!/usr/bin/env python3
+"""Unit tests for mcrfpy.HeightMap terrain generation methods (#195)
+
+Tests the HeightMap terrain methods: add_hill, dig_hill, add_voronoi,
+mid_point_displacement, rain_erosion, dig_bezier, smooth
+"""
+
+import sys
+import mcrfpy
+
+
+def test_add_hill_basic():
+    """add_hill() creates elevation at center"""
+    hmap = mcrfpy.HeightMap((50, 50), fill=0.0)
+    hmap.add_hill((25, 25), radius=10.0, height=1.0)
+
+    # Center should have highest value
+    center_val = hmap[25, 25]
+    edge_val = hmap[0, 0]
+
+    assert center_val > edge_val, f"Center ({center_val}) should be higher than edge ({edge_val})"
+    assert center_val > 0.5, f"Center should be significantly elevated, got {center_val}"
+    print("PASS: test_add_hill_basic")
+
+
+def test_add_hill_returns_self():
+    """add_hill() returns self for chaining"""
+    hmap = mcrfpy.HeightMap((20, 20))
+    result = hmap.add_hill((10, 10), 5.0, 1.0)
+    assert result is hmap
+    print("PASS: test_add_hill_returns_self")
+
+
+def test_add_hill_flexible_center():
+    """add_hill() accepts tuple, list, and Vector for center"""
+    hmap1 = mcrfpy.HeightMap((20, 20), fill=0.0)
+    hmap2 = mcrfpy.HeightMap((20, 20), fill=0.0)
+    hmap3 = mcrfpy.HeightMap((20, 20), fill=0.0)
+
+    hmap1.add_hill((10, 10), 5.0, 1.0)
+    hmap2.add_hill([10, 10], 5.0, 1.0)
+    hmap3.add_hill(mcrfpy.Vector(10, 10), 5.0, 1.0)
+
+    # All should produce same result
+    assert abs(hmap1[10, 10] - hmap2[10, 10]) < 0.001
+    assert abs(hmap1[10, 10] - hmap3[10, 10]) < 0.001
+    print("PASS: test_add_hill_flexible_center")
+
+
+def test_dig_hill_basic():
+    """dig_hill() creates depression at center using negative depth"""
+    hmap = mcrfpy.HeightMap((50, 50), fill=0.5)
+    # dig_hill with negative depth creates a crater by setting values
+    # to max(current, target_depth) where target curves from center
+    hmap.dig_hill((25, 25), radius=15.0, depth=-0.3)
+
+    # Center should have lowest value (set to the negative depth)
+    center_val = hmap[25, 25]
+    edge_val = hmap[0, 0]
+
+    assert center_val < edge_val, f"Center ({center_val}) should be lower than edge ({edge_val})"
+    assert center_val < 0, f"Center should be negative, got {center_val}"
+    print("PASS: test_dig_hill_basic")
+
+
+def test_dig_hill_returns_self():
+    """dig_hill() returns self for chaining"""
+    hmap = mcrfpy.HeightMap((20, 20))
+    result = hmap.dig_hill((10, 10), 5.0, 0.5)
+    assert result is hmap
+    print("PASS: test_dig_hill_returns_self")
+
+
+def test_add_voronoi_basic():
+    """add_voronoi() modifies heightmap"""
+    hmap = mcrfpy.HeightMap((50, 50), fill=0.0)
+    min_before, max_before = hmap.min_max()
+
+    hmap.add_voronoi(10, seed=12345)
+
+    min_after, max_after = hmap.min_max()
+
+    # Values should have changed
+    assert max_after > max_before or min_after < min_before, "Voronoi should modify values"
+    print("PASS: test_add_voronoi_basic")
+
+
+def test_add_voronoi_with_coefficients():
+    """add_voronoi() accepts custom coefficients"""
+    hmap = mcrfpy.HeightMap((30, 30), fill=0.0)
+    hmap.add_voronoi(5, coefficients=(1.0, -0.5, 0.3), seed=42)
+
+    # Should complete without error
+    min_val, max_val = hmap.min_max()
+    assert min_val != max_val, "Voronoi should create variation"
+    print("PASS: test_add_voronoi_with_coefficients")
+
+
+def test_add_voronoi_returns_self():
+    """add_voronoi() returns self for chaining"""
+    hmap = mcrfpy.HeightMap((20, 20))
+    result = hmap.add_voronoi(5)
+    assert result is hmap
+    print("PASS: test_add_voronoi_returns_self")
+
+
+def test_add_voronoi_invalid_num_points():
+    """add_voronoi() raises ValueError for invalid num_points"""
+    hmap = mcrfpy.HeightMap((20, 20))
+
+    try:
+        hmap.add_voronoi(0)
+        print("FAIL: should have raised ValueError for num_points=0")
+        sys.exit(1)
+    except ValueError:
+        pass
+
+    print("PASS: test_add_voronoi_invalid_num_points")
+
+
+def test_mid_point_displacement_basic():
+    """mid_point_displacement() generates terrain"""
+    # Use power-of-2+1 size for best results
+    hmap = mcrfpy.HeightMap((65, 65), fill=0.0)
+    hmap.mid_point_displacement(roughness=0.5, seed=12345)
+
+    min_val, max_val = hmap.min_max()
+    assert min_val != max_val, "MPD should create variation"
+    print("PASS: test_mid_point_displacement_basic")
+
+
+def test_mid_point_displacement_returns_self():
+    """mid_point_displacement() returns self for chaining"""
+    hmap = mcrfpy.HeightMap((33, 33))
+    result = hmap.mid_point_displacement()
+    assert result is hmap
+    print("PASS: test_mid_point_displacement_returns_self")
+
+
+def test_mid_point_displacement_reproducible():
+    """mid_point_displacement() is reproducible with same seed"""
+    hmap1 = mcrfpy.HeightMap((33, 33), fill=0.0)
+    hmap2 = mcrfpy.HeightMap((33, 33), fill=0.0)
+
+    hmap1.mid_point_displacement(roughness=0.6, seed=99999)
+    hmap2.mid_point_displacement(roughness=0.6, seed=99999)
+
+    # Should produce identical results
+    assert abs(hmap1[16, 16] - hmap2[16, 16]) < 0.001
+    assert abs(hmap1[5, 5] - hmap2[5, 5]) < 0.001
+    print("PASS: test_mid_point_displacement_reproducible")
+
+
+def test_rain_erosion_basic():
+    """rain_erosion() modifies terrain"""
+    hmap = mcrfpy.HeightMap((50, 50), fill=0.5)
+    # Add some hills first to have something to erode
+    hmap.add_hill((25, 25), 15.0, 0.5)
+
+    val_before = hmap[25, 25]
+    hmap.rain_erosion(1000, seed=12345)
+    val_after = hmap[25, 25]
+
+    # Erosion should change values
+    # Note: might not change if terrain is completely flat
+    # so we check min/max spread
+    min_val, max_val = hmap.min_max()
+    assert max_val > min_val, "Rain erosion should leave some variation"
+    print("PASS: test_rain_erosion_basic")
+
+
+def test_rain_erosion_returns_self():
+    """rain_erosion() returns self for chaining"""
+    hmap = mcrfpy.HeightMap((20, 20))
+    result = hmap.rain_erosion(100)
+    assert result is hmap
+    print("PASS: test_rain_erosion_returns_self")
+
+
+def test_rain_erosion_invalid_drops():
+    """rain_erosion() raises ValueError for invalid drops"""
+    hmap = mcrfpy.HeightMap((20, 20))
+
+    try:
+        hmap.rain_erosion(0)
+        print("FAIL: should have raised ValueError for drops=0")
+        sys.exit(1)
+    except ValueError:
+        pass
+
+    print("PASS: test_rain_erosion_invalid_drops")
+
+
+def test_dig_bezier_basic():
+    """dig_bezier() carves a path using negative depths"""
+    hmap = mcrfpy.HeightMap((50, 50), fill=0.5)
+
+    # Carve a path from corner to corner
+    # Use negative depths to actually dig below current terrain
+    points = ((0, 0), (10, 25), (40, 25), (49, 49))
+    hmap.dig_bezier(points, start_radius=5.0, end_radius=5.0,
+                    start_depth=-0.3, end_depth=-0.3)
+
+    # Start and end should be lower than original (carved out)
+    assert hmap[0, 0] < 0.5, f"Start should be carved, got {hmap[0, 0]}"
+    assert hmap[0, 0] < 0, f"Start should be negative, got {hmap[0, 0]}"
+    print("PASS: test_dig_bezier_basic")
+
+
+def test_dig_bezier_returns_self():
+    """dig_bezier() returns self for chaining"""
+    hmap = mcrfpy.HeightMap((20, 20))
+    result = hmap.dig_bezier(((0, 0), (5, 10), (15, 10), (19, 19)), 2.0, 2.0, 0.3, 0.3)
+    assert result is hmap
+    print("PASS: test_dig_bezier_returns_self")
+
+
+def test_dig_bezier_wrong_point_count():
+    """dig_bezier() raises ValueError for wrong number of points"""
+    hmap = mcrfpy.HeightMap((20, 20))
+
+    try:
+        hmap.dig_bezier(((0, 0), (5, 5), (10, 10)), 2.0, 2.0, 0.3, 0.3)  # Only 3 points
+        print("FAIL: should have raised ValueError for 3 points")
+        sys.exit(1)
+    except ValueError as e:
+        assert "4" in str(e)
+
+    print("PASS: test_dig_bezier_wrong_point_count")
+
+
+def test_dig_bezier_accepts_list():
+    """dig_bezier() accepts list of points"""
+    hmap = mcrfpy.HeightMap((20, 20), fill=0.5)
+    points = [[0, 0], [5, 10], [15, 10], [19, 19]]  # List instead of tuple
+    hmap.dig_bezier(points, 2.0, 2.0, 0.3, 0.3)
+    # Should complete without error
+    print("PASS: test_dig_bezier_accepts_list")
+
+
+def test_smooth_basic():
+    """smooth() reduces terrain variation"""
+    hmap = mcrfpy.HeightMap((30, 30), fill=0.0)
+    # Create sharp height differences
+    hmap.add_hill((15, 15), 5.0, 1.0)
+
+    # Get slope before smoothing (measure of sharpness)
+    slope_before = hmap.get_slope((15, 15))
+
+    hmap.smooth(iterations=3)
+
+    # Smoothing should reduce the slope
+    slope_after = hmap.get_slope((15, 15))
+    # Note: slope might not always decrease depending on terrain, so we just verify it runs
+    min_val, max_val = hmap.min_max()
+    assert max_val >= min_val
+    print("PASS: test_smooth_basic")
+
+
+def test_smooth_returns_self():
+    """smooth() returns self for chaining"""
+    hmap = mcrfpy.HeightMap((20, 20))
+    result = hmap.smooth()
+    assert result is hmap
+    print("PASS: test_smooth_returns_self")
+
+
+def test_smooth_invalid_iterations():
+    """smooth() raises ValueError for invalid iterations"""
+    hmap = mcrfpy.HeightMap((20, 20))
+
+    try:
+        hmap.smooth(iterations=0)
+        print("FAIL: should have raised ValueError for iterations=0")
+        sys.exit(1)
+    except ValueError:
+        pass
+
+    print("PASS: test_smooth_invalid_iterations")
+
+
+def test_chaining_terrain_methods():
+    """Terrain methods can be chained together"""
+    hmap = mcrfpy.HeightMap((50, 50), fill=0.0)
+
+    result = (hmap
+              .add_hill((25, 25), 10.0, 0.5)
+              .add_hill((10, 10), 8.0, 0.3)
+              .dig_hill((40, 40), 5.0, 0.2)
+              .smooth(iterations=2)
+              .normalize(0.0, 1.0))
+
+    assert result is hmap
+    min_val, max_val = hmap.min_max()
+    assert abs(min_val - 0.0) < 0.001
+    assert abs(max_val - 1.0) < 0.001
+    print("PASS: test_chaining_terrain_methods")
+
+
+def test_terrain_pipeline():
+    """Complete terrain generation pipeline"""
+    hmap = mcrfpy.HeightMap((65, 65), fill=0.0)
+
+    # Generate base terrain
+    hmap.mid_point_displacement(roughness=0.6, seed=42)
+    hmap.normalize(0.0, 1.0)
+
+    # Add features
+    hmap.add_hill((32, 32), 20.0, 0.3)  # Mountain
+    # dig_bezier with negative depths to carve a river valley
+    hmap.dig_bezier(((0, 32), (20, 20), (45, 45), (64, 32)), 3.0, 2.0, -0.3, -0.2)
+
+    # Apply erosion and smoothing
+    hmap.rain_erosion(500, seed=123)
+    hmap.smooth()
+
+    # Normalize to 0-1 range
+    hmap.normalize(0.0, 1.0)
+
+    min_val, max_val = hmap.min_max()
+    assert abs(min_val - 0.0) < 0.001
+    assert abs(max_val - 1.0) < 0.001
+    print("PASS: test_terrain_pipeline")
+
+
+def run_all_tests():
+    """Run all tests"""
+    print("Running HeightMap terrain generation tests (#195)...")
+    print()
+
+    test_add_hill_basic()
+    test_add_hill_returns_self()
+    test_add_hill_flexible_center()
+    test_dig_hill_basic()
+    test_dig_hill_returns_self()
+    test_add_voronoi_basic()
+    test_add_voronoi_with_coefficients()
+    test_add_voronoi_returns_self()
+    test_add_voronoi_invalid_num_points()
+    test_mid_point_displacement_basic()
+    test_mid_point_displacement_returns_self()
+    test_mid_point_displacement_reproducible()
+    test_rain_erosion_basic()
+    test_rain_erosion_returns_self()
+    test_rain_erosion_invalid_drops()
+    test_dig_bezier_basic()
+    test_dig_bezier_returns_self()
+    test_dig_bezier_wrong_point_count()
+    test_dig_bezier_accepts_list()
+    test_smooth_basic()
+    test_smooth_returns_self()
+    test_smooth_invalid_iterations()
+    test_chaining_terrain_methods()
+    test_terrain_pipeline()
+
+    print()
+    print("All HeightMap terrain generation tests PASSED!")
+
+
+# Run tests directly
+run_all_tests()
+sys.exit(0)