Created
August 6, 2025 10:40
-
-
Save zabop/fa593bec0da2331de5995602000b3bfb to your computer and use it in GitHub Desktop.
This file contains hidden or bidirectional Unicode text that may be interpreted or compiled differently than what appears below. To review, open the file in an editor that reveals hidden Unicode characters.
Learn more about bidirectional Unicode characters
| { | |
| "cells": [ | |
| { | |
| "cell_type": "code", | |
| "execution_count": 1, | |
| "id": "2ed81173", | |
| "metadata": {}, | |
| "outputs": [ | |
| { | |
| "name": "stdout", | |
| "output_type": "stream", | |
| "text": [ | |
| "Input file size is 4, 5\n", | |
| "0...10...20...30...40...50...60...70...80...90...100 - done.\n" | |
| ] | |
| } | |
| ], | |
| "source": [ | |
| "from PIL import Image\n", | |
| "import numpy as np\n", | |
| "\n", | |
| "arr = np.array([\n", | |
| " [1, 0.2, 0, 1],\n", | |
| " [0.2, 0.4, 0.6, 0.4], \n", | |
| " [0.4, 0.6, 0.8, 0.8], \n", | |
| " [0.6, 0.8, 0.6, 0.4], \n", | |
| " [0.8, 0.6, 1, 1], \n", | |
| " ])\n", | |
| "arr *= 255\n", | |
| "\n", | |
| "im = Image.fromarray(arr.astype(np.uint8))\n", | |
| "im.save(\"example.png\", format=\"PNG\", compress_level=0)\n", | |
| "\n", | |
| "ulx, uly = 935764.3,6668575.4\n", | |
| "\n", | |
| "pixelsize = 250\n", | |
| "lrx, lry = ulx + pixelsize * arr.shape[1], uly - pixelsize * arr.shape[0]\n", | |
| "\n", | |
| "!rm -f georeferenced_mock_small.tif\n", | |
| "!gdal_translate \\\n", | |
| " -a_srs EPSG:25831 \\\n", | |
| " -a_ullr {ulx} {uly} {lrx} {lry} \\\n", | |
| " -a_nodata 255 \\\n", | |
| " \"example.png\" \"example.tif\"" | |
| ] | |
| }, | |
| { | |
| "cell_type": "code", | |
| "execution_count": 2, | |
| "id": "e85050f4", | |
| "metadata": {}, | |
| "outputs": [], | |
| "source": [ | |
| "from typing import Literal, Callable\n", | |
| "import shapely.geometry\n", | |
| "import geopandas as gpd\n", | |
| "import numpy as np\n", | |
| "\n", | |
| "maxc = 20037508.342789244\n", | |
| "sidelength: Callable[[int], float] = lambda z: 2 * maxc / (2**z)\n", | |
| "\n", | |
| "def get_ulx_uly_lrx_lry(\n", | |
| " z: int, x: int, y: int\n", | |
| ") -> tuple[\n", | |
| " float, float, float, float\n", | |
| "]: # upper left x, upper left y, lower right x, lower right y\n", | |
| "\n", | |
| " ulx = -maxc + sidelength(z) * x\n", | |
| " uly = maxc - sidelength(z) * y\n", | |
| "\n", | |
| " lrx = -maxc + sidelength(z) * (x + 1)\n", | |
| " lry = maxc - sidelength(z) * (y + 1)\n", | |
| "\n", | |
| " return ulx, uly, lrx, lry\n", | |
| "\n", | |
| "z,x,y = 11,1085,595\n", | |
| "\n", | |
| "pixelcount = 39\n", | |
| "ulx, uly, lrx, lry = get_ulx_uly_lrx_lry(z,x,y)\n", | |
| "\n", | |
| "s = shapely.geometry.box(ulx, uly, lrx, lry)\n", | |
| "gpd.GeoSeries(s).set_crs(3857).to_file(\"tile.geojson\")\n", | |
| "\n", | |
| "s = sidelength(z) / pixelcount\n", | |
| "\n", | |
| "pixels : list[shapely.geometry.Polygon] = []\n", | |
| "\n", | |
| "for col_index in range(pixelcount):\n", | |
| " for row_index in range(pixelcount):\n", | |
| " corners: list[np.ndarray[tuple[Literal[2]], np.dtype[np.float64]]] = []\n", | |
| " for each in [[1, 1], [1, -1], [-1, -1], [-1, 1], [1, 1]]:\n", | |
| " tc = np.array([ulx + s * col_index + s/2, uly - s * row_index - s / 2])\n", | |
| " v = 0.5 * np.array([each]) * np.array([s, s])\n", | |
| " [corner] = tc + v\n", | |
| " corners.append(corner)\n", | |
| " p = shapely.geometry.Polygon(corners)\n", | |
| " pixels.append(p)\n", | |
| "\n", | |
| "gpd.GeoSeries(pixels).set_crs(3857).to_file(\"pixels.gpkg\")" | |
| ] | |
| }, | |
| { | |
| "cell_type": "code", | |
| "execution_count": 3, | |
| "id": "07b3d9e8", | |
| "metadata": {}, | |
| "outputs": [], | |
| "source": [ | |
| "!\\\n", | |
| "gdalwarp -te 1193640.633701 8375052.315150 1213208.512942 8394620.194391 \\\n", | |
| " -dstnodata 99 \\\n", | |
| " -t_srs EPSG:3857 \\\n", | |
| " -ts 39 39 \\\n", | |
| " -q \\\n", | |
| " -overwrite \\\n", | |
| " -r average \\\n", | |
| " -wo sample_steps=all -wo source_extra=10 -wo init_dest=99 \\\n", | |
| " example.tif warped_example.tif\n" | |
| ] | |
| }, | |
| { | |
| "cell_type": "code", | |
| "execution_count": 5, | |
| "id": "609f10d4", | |
| "metadata": {}, | |
| "outputs": [], | |
| "source": [ | |
| "!open ." | |
| ] | |
| }, | |
| { | |
| "cell_type": "markdown", | |
| "id": "3ab9b2d1", | |
| "metadata": {}, | |
| "source": [ | |
| "Non-nodata pixel shows up where a nodata pixel should be." | |
| ] | |
| } | |
| ], | |
| "metadata": { | |
| "kernelspec": { | |
| "display_name": "venv (3.13.2)", | |
| "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.13.2" | |
| } | |
| }, | |
| "nbformat": 4, | |
| "nbformat_minor": 5 | |
| } |
Sign up for free
to join this conversation on GitHub.
Already have an account?
Sign in to comment