orthographic/orthographic.py

544 lines
23 KiB
Python

# -------------------------------------------------------------------
# ORTHOGRAPHIC
# Your personal aerial satellite. Always on. At any altitude.*
# Developed by MarStrMind
# License: Open Software License 3.0
# Up to date version always on marstr.online
# -------------------------------------------------------------------
# orthographic.py
# Main class which handles the generation of the ortho tile.
# -------------------------------------------------------------------
import math
import os
import glob
import threading
from multiprocessing import Process
from defines import *
from log import *
from osmxml import *
from maskgen import *
from layergen import *
from photogen import *
from tileprep import *
from xp_scenery import *
# The main class which handles the rest
class mstr_orthographic:
# Constructor of class. Takes longitude and latitude.
def __init__(self, lat, lng, outfolder, pwd, prep=False):
self._lat = lat
self._long = lng
self._output = outfolder
self._pwd = pwd
self._vstep = self._findVerticalStepping()
self._latlngfld = self.latlng_folder([lat,lng])
self._prep = prep
mstr_msg("orthographic", "Initiated with LAT: " + str(lat) + ", LNG: " + str(lng))
# It did happen that the generation of photos crashed as, for some reason,
# a file in _cache was apparently used by another process (hint: it was
# not). I therefore need this test before deleting a file in _cache, so
# that generation of the orthos can move forward.
def _isFileAccessibleWin(self, src):
a = False
if os.path.isfile(src) == True:
try:
os.rename(src, src)
a = True
except OSError as e:
a = False
return a
# Need a same call for POSIX
def _isFileAccessiblePosix(self, src):
a = True
if os.access(src, os.W_OK) == False:
a = False
return a
# This will determine the vertical stepping in degrees in order to generate
# masks with a 1:1 square ratio. This is important as X-Plane textures for
# orthos can only be a power of 2, such as 2048x2048
def _findVerticalStepping(self):
scale = 1 / math.cos(math.radians(self._lat))
maxlat = (1 / scale) * mstr_zl_18
return maxlat
# To write down X-Plane .ter files, we will need to know the exact size
# of the particular longitude we are in, as this value varies depending
# on where you are on a sphere.
# Returned values is in meters.
# The current latitude is needed.
def _findWidthOfLongitude(self, lat):
dm = math.cos(math.radians(lat)) * 111.321 # <- 1 deg width at equator in km
return round(dm * 1000, 3)
# In this case we only want to acquire PBF for a latitude and longitude. Normally
# not needed for standard ortho generation.
def _generateData(self):
# The tile is constructed of many smaller parts. We walk through the
# smallest possible, from which the bigger ones are later built.
bb_lat = self._lat
bb_lng = self._long
bb_lat_edge = self._lat+self._vstep
bb_lng_edge = self._long+mstr_zl_18
cur_tile_x = 1
cur_tile_y = 1
osmxml = mstr_osmxml(0,0)
mstr_msg("orthographic", "Set initial coordinates and bounding box for OSM acquisition")
# The highest encountered tile numbers
# This is needed to produce the zoom level 16 images
top_lat = 1
top_lng = 1
# We need to know the highest possible latitude and longitude tile numbers,
# in case we render at the edge
mlat = 1
mlng = 1
while bb_lat < self._lat + 1:
bb_lat = bb_lat + self._vstep
mlat = mlat+1
while bb_lng < self._long + 1:
bb_lng = bb_lng + mstr_zl_18
mlng = mlng+1
mstr_msg("orthographic", "Max lat tile: " + str(mlat) + " - max lng tile: " + str(mlng))
maxlatlng = [ mlat, mlng ]
# Reset these two
bb_lat = self._lat
bb_lng = self._long
# Previously, I downloaded all XML files in one go - but to ease the
# stress on OSM servers and my server, we will do acquire the data
# only for the current processed part of the tile.
for lat_grid in range(1, maxlatlng[0]+1):
for lng_grid in range(1, maxlatlng[1]+1):
# Adjust bounding box
osmxml.adjust_bbox(bb_lat, bb_lng, bb_lat_edge, bb_lng_edge)
osmxml.generate_osm(cur_tile_y, cur_tile_x) # <- This acquires current OSM info
# Adjust longitude coordinates
cur_tile_x = cur_tile_x+1
bb_lng = bb_lng + mstr_zl_18
bb_lng_edge = bb_lng_edge + mstr_zl_18
mstr_msg("orthographic", "Adjustment of longitude performed")
# Adjust peak longitude tile number
if cur_tile_x > top_lng:
top_lng = cur_tile_x
# Adjust latitude and all other values when we get here
cur_tile_y = cur_tile_y+1
cur_tile_x = 1
bb_lng = self._long
bb_lng_edge = self._long + mstr_zl_18
bb_lat = bb_lat + self._vstep
bb_lat_edge = bb_lat_edge + self._vstep
mstr_msg("orthographic", "Adjustment of latitude performed")
# Adjust peak latitude number
if cur_tile_y > top_lat:
top_lat = cur_tile_y
# Start the multi-threaded build of all orthos
# amtsmt = AmountSimultaneous - so how many orthos you want to
# generate at the same time. You may need to fine tune this value
# so that you don't overload your machine.
def _generateOrthos_mt(self, amtsmt):
# Need to know maximum values first
bb_lat = self._lat
bb_lng = self._long
bb_lat_edge = self._lat+self._vstep
bb_lng_edge = self._long+mstr_zl_18
mlat = 1
mlng = 1
while bb_lat < self._lat + 1:
bb_lat = bb_lat + self._vstep
mlat = mlat+1
while bb_lng < self._long + 1:
bb_lng = bb_lng + mstr_zl_18
mlng = mlng+1
mstr_msg("orthographic", "Max lat tile: " + str(mlat) + " - max lng tile: " + str(mlng))
maxlatlng = [ mlat, mlng ]
procs = []
for p in range(1, amtsmt+1):
proc = Process(target=self._buildOrtho, args=[1, p, amtsmt])
procs.append(proc)
proc.start()
mstr_msg("orthographic", "Ortho threads started")
# Starts a threading loop to build orthos, with the defined starting point in
# the lat-lng grid. You will also need to provide the horizontal stepping so
# that the thread keeps running.
def _buildOrtho(self, v, h, step):
# Starting point
grid_lat = v
grid_lng = h
# The tile is constructed of many smaller parts. We walk through the
# smallest possible, from which the bigger ones are later built.
bb_lat = self._lat
bb_lng = self._long
bb_lat_edge = self._lat+self._vstep
bb_lng_edge = self._long+mstr_zl_18
# We need to know the highest possible latitude and longitude tile numbers,
# in case we render at the edge
mlat = 1
mlng = 1
while bb_lat < self._lat + 1:
bb_lat = bb_lat + self._vstep
mlat = mlat+1
while bb_lng < self._long + 1:
bb_lng = bb_lng + mstr_zl_18
mlng = mlng+1
mstr_msg("orthographic", "Max lat tile: " + str(mlat) + " - max lng tile: " + str(mlng))
maxlatlng = [ mlat, mlng ]
while grid_lat <= maxlatlng[0]:
ddsf = mstr_datafolder + "z_orthographic/orthos/" + self._latlngfld + "/" + str(grid_lat) + "_" + str(grid_lng) + ".dds"
if os.path.isfile(ddsf) == False:
# Reset these two
bb_lat = self._lat + ((grid_lat-1)*self._vstep)
bb_lng = self._long + ((grid_lng-1)*mstr_zl_18)
bb_lat_edge = self._lat + ((grid_lat-1)*self._vstep) + self._vstep
bb_lng_edge = self._long + ((grid_lng-1)*mstr_zl_18) + mstr_zl_18
osmxml = mstr_osmxml()
osmxml.adjust_bbox(bb_lat, bb_lng, bb_lat_edge, bb_lng_edge)
osmxml.acquire_osm(grid_lat, grid_lng)
# Let the user know
mstr_msg("orthographic", "Generating missing orthophoto " + str(grid_lat) + "-" + str(grid_lng))
# Check for work to be done
layers = self.determineLayerWork(osmxml)
# We need to walk through the array of layers,
# in their z-order.
# For each layer, we will generate the mask, the layer image
# itself, and finally, compose the ortho photo.
mstr_msg("orthographic", "Beginning generation of layers")
# In here we store the layers
photolayers = []
waterlayers = []
# The masks are handed to layergen in sequence. The layers are then
# in turn handed to photogen.
curlyr = 1
wtr_info = False
for layer in layers:
# Let the user know
mstr_msg("orthographic", "Processing layer " + str(curlyr) + " of " + str(len(layers)))
# Generate the mask
mg = mstr_maskgen( [self._lat, grid_lat, self._long, grid_lng], self._vstep, layer[0], layer[1], layer[2])
if layer[0] == "building":
mg.set_tile_width(self._findWidthOfLongitude(bb_lat))
mg.set_latlng_numbers(self._lat, grid_lat, self._long, grid_lng)
mask = mg._build_mask(osmxml)
# Generate the layer
lg = mstr_layergen(layer[0], layer[1], self._lat, grid_lat, self._long, grid_lng, layer[2])
lg.set_max_latlng_tile(maxlatlng)
lg.set_latlng_folder(self._latlngfld)
#lg.open_db()
lg.open_tile_info()
lyr = lg.genlayer(mask, osmxml)
photolayers.append(lyr)
if (layer[0] == "natural" and layer[1] == "water") or (layer[0] == "water" and layer[1] == "lake") or (layer[0] == "water" and layer[1] == "pond") or (layer[0] == "water" and layer[1] == "river") or (layer[0] == "waterway" and layer[1] == "river"):
waterlayers.append(lyr)
if wtr_info == False:
wtr_info = True
wtrfile = mstr_datafolder + "z_orthographic/data/" + self._latlngfld + "/wtrfile"
with open(wtrfile, 'a') as textfile:
textfile.write(str(grid_lat) + " " + str(grid_lng) + "\r\n")
curlyr = curlyr+1
mstr_msg("orthographic", "All layers created")
# We should have all layers now.
# Snap a photo with our satellite :)
mstr_msg("orthographic", "Generating ortho photo")
pg = mstr_photogen(self._lat, self._long, grid_lat, grid_lng, maxlatlng[0], maxlatlng[1])
pg.genphoto(photolayers, waterlayers)
mstr_msg("orthographic", " -- Ortho photo generated -- ")
print("")
print("")
# Perform adjustment of grid position
n_lng = grid_lng + step
if n_lng > maxlatlng[1]:
np = n_lng - maxlatlng[1]
grid_lng = np
grid_lat = grid_lat+1
else:
grid_lng = n_lng
# Prepares the entire tile
def _prepareTile(self):
mstr_msg("orthographic", "Beginning construction of tile")
# We need to know which platform we are on
os_platform = os.name
# Create the _cache folder, should it not exist.
# Temporary images for the ortho tile generation go here
if not os.path.exists(self._output + "/_cache"):
os.makedirs(self._output + "/_cache")
mstr_msg("orthographic", "Created _cache folder.")
# Generate the Tiles/lat-lng folder for the finished tile
if not os.path.exists(self._output + "/z_orthographic"):
os.makedirs(self._output + "/z_orthographic")
mstr_msg("orthographic", "Created z_orthographic folder")
# Generate the orthos folder
if not os.path.exists(self._output + "/z_orthographic/orthos"):
os.makedirs(self._output + "/z_orthographic/orthos")
mstr_msg("orthographic", "Created tile orthos folder")
if not os.path.exists(self._output + "/z_orthographic/orthos" + self._latlngfld):
os.makedirs(self._output + "/z_orthographic/orthos/" + self._latlngfld, exist_ok=True)
# Generate the database folder
if not os.path.exists(self._output + "/z_orthographic/data"):
os.makedirs(self._output + "/z_orthographic/data")
mstr_msg("orthographic", "Created tile database folder")
if not os.path.exists(self._output + "/z_orthographic/data/" + self._latlngfld):
os.makedirs(self._output + "/z_orthographic/data/" + self._latlngfld)
# X-Plane specific
if mstr_xp_genscenery == True:
btnum = self.find_earthnavdata_number()
btstr = self.latlng_folder(btnum)
if not os.path.exists(self._output + "/z_orthographic/terrain"):
os.makedirs(self._output + "/z_orthographic/terrain")
mstr_msg("orthographic", "[X-Plane] Created terrain files folder")
if not os.path.exists(self._output + "/z_orthographic/terrain/" + self._latlngfld):
os.makedirs(self._output + "/z_orthographic/terrain/" + self._latlngfld)
if not os.path.exists(self._output + "/z_orthographic/Earth nav data"):
os.makedirs(self._output + "/z_orthographic/Earth nav data")
mstr_msg("orthographic", "[X-Plane] Created Earth nav folder")
if not os.path.exists(self._output + "/z_orthographic/Earth nav data/" + btstr):
os.makedirs(self._output + "/z_orthographic/Earth nav data/" + btstr)
if mstr_xp_scn_normalmaps == True:
if not os.path.exists(self._output + "/z_orthographic/normals"):
os.makedirs(self._output + "/z_orthographic/normals")
mstr_msg("orthographic", "[X-Plane] created tile normal maps folder")
if not os.path.exists(self._output + "/z_orthographic/normals/" + self._latlngfld):
os.makedirs(self._output + "/z_orthographic/normals/" + self._latlngfld)
# The tile is constructed of many smaller parts. We walk through the
# smallest possible, from which the bigger ones are later built.
bb_lat = self._lat
bb_lng = self._long
bb_lat_edge = self._lat+self._vstep
bb_lng_edge = self._long+mstr_zl_18
cur_tile_x = 1
cur_tile_y = 1
#osmxml = mstr_osmxml(0,0)
mstr_msg("orthographic", "Set initial coordinates and bounding box for OSM acquisition")
# The highest encountered tile numbers
# This is needed to produce the zoom level 16 images
top_lat = 1
top_lng = 1
# We need to know the highest possible latitude and longitude tile numbers,
# in case we render at the edge
mlat = 1
mlng = 1
while bb_lat < self._lat + 1:
bb_lat = bb_lat + self._vstep
mlat = mlat+1
while bb_lng < self._long + 1:
bb_lng = bb_lng + mstr_zl_18
mlng = mlng+1
mstr_msg("orthographic", "Max lat tile: " + str(mlat) + " - max lng tile: " + str(mlng))
maxlatlng = [ mlat, mlng ]
# Reset these two
bb_lat = self._lat
bb_lng = self._long
# We will now prepare the graphic tile generation. We do this by only generating
# the masks and determine which sources to use in the actual images.
# Previously, I downloaded all XML files in one go - but to ease the
# stress on OSM servers and my server, we will do acquire the data
# only for the current processed part of the tile.
for lat_grid in range(1, maxlatlng[0]+1):
for lng_grid in range(1, maxlatlng[1]+1):
# Adjust bounding box
osmxml = mstr_osmxml()
osmxml.adjust_bbox(bb_lat, bb_lng, bb_lat_edge, bb_lng_edge)
osmxml.acquire_osm(lat_grid, lng_grid)
mstr_msg("orthographic", "Adjusted bounding box for XML object")
# Check for work to be done
layers = self.determineLayerWork(osmxml)
curlyr = 1
for layer in layers:
if layer[2] == False and layer[0] != "building":
# Let the user know
mstr_msg("orthographic", "Processing layer " + str(curlyr) + " of " + str(len(layers)))
# Generate the mask
mg = mstr_maskgen( [self._lat, cur_tile_y, self._long, cur_tile_x], self._vstep, layer[0], layer[1], layer[2])
mask = mg._build_mask(osmxml, is_prep=True) # We need an object here
tp = mstr_tileprep(self._lat, self._long, lat_grid, lng_grid, layer[0], layer[1], mask, False)
tp._prepareTile()
curlyr = curlyr+1
# Adjust longitude coordinates
cur_tile_x = cur_tile_x+1
bb_lng = bb_lng + mstr_zl_18
bb_lng_edge = bb_lng_edge + mstr_zl_18
mstr_msg("orthographic", "Adjustment of longitude performed")
# Adjust peak longitude tile number
if cur_tile_x > top_lng:
top_lng = cur_tile_x
# Adjust latitude and all other values when we get here
cur_tile_y = cur_tile_y+1
cur_tile_x = 1
bb_lng = self._long
bb_lng_edge = self._long + mstr_zl_18
bb_lat = bb_lat + self._vstep
bb_lat_edge = bb_lat_edge + self._vstep
mstr_msg("orthographic", "Adjustment of latitude performed")
# Adjust peak latitude number
if cur_tile_y > top_lat:
top_lat = cur_tile_y
# Generates X-Plane 11/12 scenery with
# - the finished orthos
# - a current LIDAR scan of the terrain
def generate_xp_scenery(self):
mstr_msg("orthographic", "[X-Plane] Generation of scenery started")
# This call appears quite often... surely this can be done better
mlat = 1
mlng = 1
bb_lat = self._lat
bb_lng = self._long
bb_lat_edge = self._lat+self._vstep
bb_lng_edge = self._long+mstr_zl_18
while bb_lat < self._lat + 1:
bb_lat = bb_lat + self._vstep
mlat = mlat+1
while bb_lng < self._long + 1:
bb_lng = bb_lng + mstr_zl_18
mlng = mlng+1
mstr_msg("orthographic", "Max lat tile: " + str(mlat) + " - max lng tile: " + str(mlng))
maxlatlng = [ mlat, mlng ]
# The object that handles it all
xpscn = mstr_xp_scenery(self._lat, self._long, maxlatlng[0], maxlatlng[1], self._vstep, self._latlngfld)
mstr_msg("orthographic", "[X-Plane] Scenery object instantiated")
# Download LIDAR scan from our endpoint
xpscn.acquire_elevation_data()
mstr_msg("orthographic", "[X-Plane] Elevation data acquired")
# Generate the .ter files
xpscn.build_ter_files()
mstr_msg("orthographic", "[X-Plane] Terrain files (.ter) generated and written")
# And lastly, generate the mesh
xpscn.generate_terrain_mesh()
mstr_msg("orthographic", "[X-Plane] Scenery mesh constructed")
# Convert the DSF
xpscn.build_and_convert_dsf()
mstr_msg("orthographic", "[X-Plane] DSF generated")
# Checks which layers need to be generated, and what kind of layer it is
def determineLayerWork(self, xmlobj):
mstr_msg("orthographic", "Checking for work to be performed")
layers = []
#tilexml = mstr_datafolder + "_cache/tile.xml"
#xml = mstr_osmxml(0,0)
way = xmlobj.acquire_waypoint_data()
rls = xmlobj.acquire_relations()
for l in mstr_ortho_layers:
# Check if there is anything to render
has_way = False
has_rls = False
for w in way:
if w[2] == l[0] and w[3] == l[1]:
has_way = True
break
for r in rls:
if l[0] in r[1] and l[1] in r[1]:
has_rls = True
break
if has_way == True or has_rls == True:
mstr_msg("orthographic", "Adding: " + l[0]+":"+l[1])
is_line = False
for s in mstr_ortho_layers:
if s[0] == l[0] and s[1] == l[1]:
if isinstance(s[2], int) == False:
is_line = False
break
if isinstance(s[2], int) == True:
is_line = True
break
ly = (l[0], l[1], is_line)
layers.append(ly)
mstr_msg("orthographic", "A total of " + str(len(layers)) + " layers were found")
return layers
# Construct a folder name for latitude and longitude
def latlng_folder(self, numbers):
fstr = ""
if numbers[0] >= 0: fstr = "+"
if numbers[0] < 0: fstr = "-"
if abs(numbers[0]) < 10: fstr = fstr + "0" + str(numbers[0])
if abs(numbers[0]) >= 10 and numbers[0] <= 90: fstr = fstr + str(numbers[0])
if numbers[1] >= 0: fstr = fstr + "+"
if numbers[1] < 0: fstr = fstr + "-"
if abs(numbers[1]) < 10: fstr = fstr + "00" + str(numbers[1])
if abs(numbers[1]) >= 10 and numbers[0] <= 99: fstr = fstr + "0" + str(numbers[1])
if abs(numbers[1]) >= 100 : fstr = fstr + str(numbers[1])
return fstr
# Find the next "by-ten" numbers for the current latitude and longitude
def find_earthnavdata_number(self):
earthnavdata = []
lat = abs(int(self._lat / 10) * 10)
lng = abs(int(self._long / 10) * 10)
earthnavdata.append(lat)
earthnavdata.append(lng)
return earthnavdata