Minetest world maker using UK LIDAR data

lidar2blocks.ipynb

{
 "cells": [
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "# Turing LIDAR data into voxels"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 1,
   "metadata": {},
   "outputs": [],
   "source": [
    "import pylas\n",
    "import numpy as np\n",
    "import random"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "The notebook creates a file `blocks.csv` which can be used with the [csv2terrain](https://github.com/quantumjim/csv2terrain/blob/master/README.md) Minetest mod. The file describes the terrain around a given point, created from LIDAR data.\n",
    "\n",
    "To do this, you'll need to download the data for the relvant 5km x 5km tile from the 2019 'National LIDAR Data Programme Point Cloud' [here](https://environment.data.gov.uk/DefraDataDownload/?Mode=survey). Then set the path to the `.laz.` file below."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 2,
   "metadata": {},
   "outputs": [],
   "source": [
    "# for example the tile TN47NW\n",
    "file = 'P_10772/TL4075_P_10772_20191211_20191211.laz'\n",
    "# For Ely Cathedral use 'P_10786/TL5080_P_10786_20181113_20181113.laz'\n",
    "# For King's College Chapel use 'P_10771/TL4055_P_10771_20191211_20191211.laz'"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "You'll also need the coordinate of the place that you want to be the center of the map. Use the X (Easting), Y (Northing) values given [here](https://gridreferencefinder.com)."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 3,
   "metadata": {},
   "outputs": [],
   "source": [
    "# for example, the church of Sutton-in-the-Isle\n",
    "coord = (544837,278978)\n",
    "# For Ely Cathedral use (554053,280271)\n",
    "# For King's College Chapel use (544733 , 258395)"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "The last thing needed to specify the map is its size. The `size` you specify is the length and width in terms of blocks. In the real world, this corresponds to... Actually I don't really know. On the order of 1m."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 4,
   "metadata": {},
   "outputs": [],
   "source": [
    "# I dare go no higher than 256\n",
    "size = 256"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "The following function opens the file and extracts the data required for. In `heights` we store the z coordinate given the x and y of each point in the point cloud. In `classification` we similarly store the 'classification', which is explained [here](https://www.arcgis.com/apps/MapJournal/index.html?appid=c6cef6cc642a48838d38e722ea8ccfee)."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 5,
   "metadata": {},
   "outputs": [],
   "source": [
    "def get_points(x0,y0,file,size):\n",
    "    \n",
    "    # read in all the data for the tile\n",
    "    cloud = pylas.read(file)\n",
    "\n",
    "    # make a mask for the points needed\n",
    "    x, y = cloud.x.copy(), cloud.y.copy()\n",
    "    mask = (x>=x0-size/2) & (x<x0+size/2) & (y>=y0-size/2) & (y<y0+size/2)\n",
    "\n",
    "    # extract height (which is stored at index 2) and classification (which seems to be at 5)\n",
    "    # also divide x and y by 100, because that is required for some reason\n",
    "    height = {}\n",
    "    classification = {}\n",
    "    for point in cloud.points[mask]:\n",
    "        xx = int(point[0]/100-x0+size/2)\n",
    "        yy = int(size-1-(point[1]/100-y0+size/2))\n",
    "        height[xx,yy] = point[2]\n",
    "        classification[xx,yy] = point[5]\n",
    "    \n",
    "    # rejig the heights so the minimum is 0, and also do the mysterious 'division by 100' trick\n",
    "    min_h = min(height.values())\n",
    "    for point in height:\n",
    "        height[point] = int((height[point]-min_h)/100)\n",
    "\n",
    "    return height,classification"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "This function is use by `get_blocks` which takes in the three peices of info from earlier and writes the `blocks.csv` file to disk."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 6,
   "metadata": {},
   "outputs": [],
   "source": [
    "def get_blocks(file,coord,size):\n",
    "    \n",
    "    # first, get the info we need\n",
    "    height,classification = get_points(coord[0],coord[1],file,size)\n",
    "\n",
    "    # the following list specifies what values for the classification correspond to in terms of Minetest blocks\n",
    "    # https://wiki.minetest.net/Itemstrings/default\n",
    "    block_type = [ 'water_source',\n",
    "                   'cobble',\n",
    "                   'dirt_with_grass',\n",
    "                   'mossycobble',\n",
    "                   'stone',\n",
    "                   'leaves',\n",
    "                   'stonebrick',\n",
    "                   'mossycobble',\n",
    "                   'mossycobble']\n",
    "\n",
    "    # now we go through each position and assign a block type\n",
    "    # for x,y in the map, grass goes at the lowest height\n",
    "    blocks = {}\n",
    "    for x in range(size):\n",
    "        for y in range(size):\n",
    "            blocks[x,y,1] = 'dirt_with_grass'\n",
    "    # for each x,y with a point in the point cloud, the corresponding block is used up to the corresponding height\n",
    "    # this will lead to some overlaps if two points have the same x,y, but oh well\n",
    "    for x in range(size):\n",
    "        for y in range(size):\n",
    "            if (x,y) in height:\n",
    "                z = height[x,y]\n",
    "                for zz in range(1,z+2):\n",
    "                    blocks[x,y,zz] = block_type[classification[x,y]]\n",
    "            else:\n",
    "                # if there was no point, some water is added at the second lowest height\n",
    "                blocks[x,y,1] = block_type[0]\n",
    "    \n",
    "    # now we write all the required info in the file\n",
    "    # the player will start at the middle of the map\n",
    "    max_height = max(height.values())+2\n",
    "    x,y = int(size/2),int(size/2)\n",
    "    if (x,y) in height:\n",
    "        player = x,height[x,y]+2,y\n",
    "    else:\n",
    "        player = x,int(max_height/2),y\n",
    "    with open('blocks.csv','w') as file:\n",
    "        file.write( '0,0,0,min,\\n' )\n",
    "        file.write( str(size)+','+str(max_height)+','+str(size)+','+'max'+',\\n' )\n",
    "        file.write( str(player[0])+','+str(player[1])+','+str(player[2])+','+'player'+',\\n' )\n",
    "        for (x,y,z) in blocks:\n",
    "            file.write( str(x)+','+str(z)+','+str(y)+','+blocks[x,y,z]+',\\n' )"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "All that's left is to actually run the function."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 7,
   "metadata": {},
   "outputs": [],
   "source": [
    "get_blocks(file,coord,size)"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "For details on how to create the map in Minetest, see the details of the [csv2terrain](https://github.com/quantumjim/csv2terrain/blob/master/README.md) mod."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": []
  }
 ],
 "metadata": {
  "kernelspec": {
   "display_name": "Python 3",
   "language": "python",
   "name": "python3"
  },
  "language_info": {
   "codemirror_mode": {
    "name": "ipython",
    "version": 3
   },
   "file_extension": ".py",
   "mimetype": "text/x-python",
   "name": "python",
   "nbconvert_exporter": "python",
   "pygments_lexer": "ipython3",
   "version": "3.8.5"
  }
 },
 "nbformat": 4,
 "nbformat_minor": 4
}