Skip to content
New issue

Have a question about this project? Sign up for a free GitHub account to open an issue and contact its maintainers and the community.

By clicking “Sign up for GitHub”, you agree to our terms of service and privacy statement. We’ll occasionally send you account related emails.

Already on GitHub? Sign in to your account

replace the .sh scripts with py scriprts #72

Open
Pomax opened this issue Feb 5, 2023 · 5 comments
Open

replace the .sh scripts with py scriprts #72

Pomax opened this issue Feb 5, 2023 · 5 comments

Comments

@Pomax
Copy link

Pomax commented Feb 5, 2023

Since we already need to have python installed as scripting language for a local (non-docker) installation, we shouldn't strictly speaking need bash (or zshell, etc) to run the create_dataset operation, we can implement that in Python instead and suddenly things will even work on Windows. The requirements would need to have requests added, but then we can use a meshes.py:

meshes = [
    ('https://srtm.csi.cgiar.org/wp-content/uploads/files/250m/SRTM_NE_250m_TIF.rar',
     'some.SRTM_NE_250m_TIF.rar', 'SRTM_NE_250m.tif', 10, 10,),
    ('https://srtm.csi.cgiar.org/wp-content/uploads/files/250m/SRTM_SE_250m_TIF.rar',
     'some.SRTM_SE_250m_TIF.rar', 'SRTM_SE_250m.tif', 10, 10),
    ('https://srtm.csi.cgiar.org/wp-content/uploads/files/250m/SRTM_W_250m_TIF.rar',
     'some.SRTM_W_250m_TIF.rar', 'SRTM_W_250m.tif', 10, 20),
]

With as download script:

import os
import requests
import meshes

for mesh in meshes.meshes:
    url, filename, raster, xtiles, ytiles = mesh
    if os.path.exists(filename):
        print(f'Skipping download for {filename} (file already exists)')
        continue
    print(f'Downloading {filename}...')
    with requests.get(url, stream=True) as r:
        r.raise_for_status()
        with open(filename, 'wb') as f:
            for chunk in r.iter_content(chunk_size=8192):
                f.write(chunk)
            f.close()

print(f'\nRemember to unrar these files before running create_tiles.')

And then with gdal installed (which for Windows users means installing from .whl file as counterpart to the current instructions for non-Windows users), the tile creation script could be:

from osgeo import gdal
import meshes

xtiles = 10
ytiles = 10
tilesets = [raster for __url, __filename, raster in meshes.meshes]

for raster in tilesets:
    raster_name = raster.replace('.tif','')
    print(f'Converting {raster_name} to tiles')

    data = gdal.Open(raster)
    ulx, xres, xskew, uly, yskew, yres  = data.GetGeoTransform()
    lrx = ulx + (data.RasterXSize * xres)
    lry = uly + (data.RasterYSize * yres)
    xmin = ulx
    xsize = lrx - xmin
    ysize = uly - lry
    xdif = xsize / xtiles

    for x in range(0, xtiles-1):
        xmax = xmin + xdif
        ymax = uly
        ydif = ysize / ytiles

        for y in range(0,ytiles-1):
            ymin = ymax - ydif
            output_name = f'{raster_name}_{x}_{y}.tif'
            print(f'Writing {output_name}')
            gdal.Translate(output_name, data, projWin = [xmin, ymax, xmax, ymin])
            ymax = ymin

        xmin= xmax

print('\nFinished conversion')

The user would still need to unrar the archives themselves, but there are a million different utils for this, so guessing whether someone has 7zip, or winrar, or any number of other utils installed when they can just double click the file and have it extract using whatever is installed is probably unnecessary.

As a final blocker, of course, gunicorn does not work on Windows, but looking at the server.py, while gunicorn is convenient, it's not really being used for anything beyond simple route handing, which is pretty easily rewritten to a class that extends BaseHTTPRequestHandler.

Filed as issue rather than as a PR, because while "you're already going to need one scripting language, no need for any other ones" is a powerful enabler, not everyone cares about platform independence as much as some users do =)

@Pomax
Copy link
Author

Pomax commented Feb 5, 2023

(To preempt "Why not just use docker?", mostly it's because docker is a monster to install on Windows. Any project that you can run without needing Docker is a sigh of relief)

@Pomax
Copy link
Author

Pomax commented Feb 5, 2023

Example rewrite without gunicorn (obviously not polished, but good enough for my needs):

import configparser
import json
import os
from http.server import BaseHTTPRequestHandler, HTTPServer
from urllib.parse import parse_qs, urlparse

from gdal_interfaces import GDALTileInterface


class InternalException(ValueError):
    """
    Utility exception class to handle errors internally and return error codes to the client
    """
    pass


print('Reading config file ...')
parser = configparser.ConfigParser()
parser.read('config.ini')

HOST = '127.0.0.1'
PORT = 9000
DATA_FOLDER = parser.get('server', 'data-folder')
OPEN_INTERFACES_SIZE = parser.getint('server', 'open-interfaces-size')


interface = GDALTileInterface(
    DATA_FOLDER, '%s/summary.json' % DATA_FOLDER, OPEN_INTERFACES_SIZE)
interface.create_summary_json()


def get_elevation(lat, lng):
    """
    Get the elevation at point (lat,lng) using the currently opened interface
    :param lat:
    :param lng:
    :return:
    """
    try:
        elevation = interface.lookup(lat, lng)
    except:
        return {
            'latitude': lat,
            'longitude': lng,
            'error': 'No such coordinate (%s, %s)' % (lat, lng)
        }

    return {
        'latitude': lat,
        'longitude': lng,
        'elevation': elevation
    }


def lat_lng_from_location(location_with_comma):
    """
    Parse the latitude and longitude of a location in the format "xx.xxx,yy.yyy" (which we accept as a query string)
    :param location_with_comma:
    :return:
    """
    try:
        lat, lng = [float(i) for i in location_with_comma.split(',')]
        return lat, lng
    except:
        raise InternalException(json.dumps(
            {'error': 'Bad parameter format "%s".' % location_with_comma}))


def query_to_locations(query):
    """
    Grab a list of locations from the query and turn them into [(lat,lng),(lat,lng),...]
    :return:
    """
    locations = query['locations'][0] if 'locations' in query else None
    if not locations:
        raise InternalException(json.dumps(
            {'error': '"Locations" is required.'}))

    return [lat_lng_from_location(l) for l in locations.split('|')]


def body_to_locations(payload):
    """
    Grab a list of locations from the body and turn them into [(lat,lng),(lat,lng),...]
    :return:
    """
    try:
        locations = payload.get('locations', None)
    except Exception:
        raise InternalException(json.dumps({'error': 'Invalid JSON.'}))

    if not locations:
        raise InternalException(json.dumps(
            {'error': '"Locations" is required in the body.'}))

    latlng = []
    for l in locations:
        try:
            latlng += [(l['latitude'], l['longitude'])]
        except KeyError:
            raise InternalException(json.dumps(
                {'error': '"%s" is not in a valid format.' % l}))

    return latlng


def do_lookup(locations):
    """
    Generic method which gets the locations in [(lat,lng),(lat,lng),...] format by calling get_locations_func
    and returns an answer ready to go to the client.
    :return:
    """
    return {'results': [get_elevation(lat, lng) for (lat, lng) in locations]}


class OpenElevationServer(BaseHTTPRequestHandler):
    def set_headers(self):
        self.send_response(200)
        self.send_header('Access-Control-Allow-Origin', '*')
        self.send_header('Cache-Control', 'no-cache')
        self.send_header('Content-type', 'application/json')
        self.end_headers()

    def log_request(self, code='-', size='-'):
        # Don't log regular requests, only error
        return

    def do_GET(self):
        """
        GET method. Uses query_to_locations.
        :return:
        """
        if self.path == '/favicon.ico':
            return self.send_response(404)

        query = parse_qs(urlparse(self.path).query)

        try:
            data = do_lookup(query_to_locations(query))
            response = json.dumps(data).encode('utf-8')
            self.set_headers()
            self.wfile.write(response)
        except :
            self.send_response(400)

    def do_OPTIONS(self):
        self.send_response(200)
        self.send_header('Access-Control-Allow-Credentials', 'true')
        self.send_header('Access-Control-Allow-Origin', '*')
        self.send_header('Access-Control-Allow-Methods', 'GET, POST, OPTIONS')
        self.send_header("Access-Control-Allow-Headers", "X-Requested-With, Content-type")
        self.end_headers()

    def do_POST(self):
        """
        GET method. Uses body_to_locations.
        :return:
        """
        payload = json.loads(self.rfile.read(
            int(self.headers['Content-Length'])))

        try:
            data = do_lookup(body_to_locations(payload))
            response = json.dumps(data).encode('utf-8')
            self.set_headers()
            self.wfile.write(response)
        except:
            self.send_response(400)


def run():
    try:
        webServer = HTTPServer((HOST, PORT), OpenElevationServer)
        print(f'Server started http://{HOST}:{PORT}')
        webServer.serve_forever()
    except KeyboardInterrupt:
        webServer.server_close()
        print('Server stopped')
        os._exit(1)


if __name__ == "__main__":
    run()

@Pomax
Copy link
Author

Pomax commented Feb 6, 2023

Mind you, when I run this, things mostly work fine, but not when I'm working with coordinates in New Zealand. Calling http://localhost:9000/?locations=-44.9856891,168.321971 yields an error index 14551 is out of bounds for axis 1 with size 10081 in the terminal.

Is there a simple or clever way to check which tiles cover which part of the planet, to see if there's any gaps?

@aliasfoxkde
Copy link

Perfect, this was very helpful for what I'm doing. Thanks!

@davidetedesco1
Copy link

Hi where i can get meshes.py?

Sign up for free to join this conversation on GitHub. Already have an account? Sign in to comment
Labels
None yet
Projects
None yet
Development

No branches or pull requests

3 participants