# ------------------------------------------------------------------- # 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 # ------------------------------------------------------------------- # xp_scenery.py # This class builds an elevation mesh from provided DEM data, and # generates all required files that are needed to provide a usable # scenery package for X-Plane. # ------------------------------------------------------------------- import os import math import urllib.request import numpy from defines import * from log import * from PIL import Image, ImageFilter, ImageEnhance class mstr_xp_scenery: # Set required variables def __init__(self, lat, lng, mlat, mlng, vstep, latlngfld): mstr_msg("xp_scenery", "[X-Plane] Scenery generator instantiated") self._lat = lat self._lng = lng self._mlat = mlat self._mlng = mlng self._vstep = vstep self._latlngfld = latlngfld self._demfn = self.build_dem_filename() self._dsfstring = "" self._demdata = None # To be populated when the mesh is built self._demcoord = None # Also to be populated when mesh is built # Build the correct file name for the elevation model def build_dem_filename(self, xes=False): fn = "" if self._lat > 0: fn = fn + "N" else: fn = fn + "S" if abs(self._lat) < 10: fn = fn + "0" + str(self._lat) if abs(self._lat) >= 10 and abs(self._lat) <= 90: fn = fn + str(self._lat) if self._lng > 0: fn = fn + "E" else: fn = fn + "W" if abs(self._lng) < 10: fn = fn + "00" + str(self._lng) if abs(self._lng) >= 10 and abs(self._lng) <= 99: fn = fn + "0" + str(self._lng) if abs(self._lng) >= 100 : fn = fn + str(self._lng) if xes == False: fn = fn + ".hgt" if xes == True: fn = fn + ".xes" mstr_msg("xp_scenery", "[X-Plane] DEM file name constructed: " + fn) return fn # Build the DSF for the ortho photo overlay def build_and_convert_dsf(self): end = self.find_earthnavdata_number() llf = self.xplane_latlng_folder(end) meshtxt = mstr_datafolder + "_cache/mesh_"+self._latlngfld+".txt" scr = mstr_datafolder + "z_orthographic/data/" + self._latlngfld + "/meshscript.txt" cmd = mstr_xp_dsftool + " --text2dsf " + meshtxt + " '" + mstr_datafolder + "z_orthographic/Earth nav data/" + llf + "/" + self._latlngfld + ".dsf'" os.system(cmd) # Find exact with of longitude def find_width_of_longitude(self, lat): dm = math.cos(math.radians(lat)) * 111.321 # <- 1 deg width at equator in km return round(dm * 1000, 3) # 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._lng / 10) * 10) earthnavdata.append(lat) earthnavdata.append(lng) return earthnavdata # Construct an X-Plane compatible folder name for latitude and longitude def xplane_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 # Acquires the elevation data as needed - # you can either acquire the older STRM set with lower resolution, # or a modern LIDAR scan with higher resolution. However, the # LIDAR files are much larger and generates a more taxing mesh. # If you are testing, acquire the low resolution file. # Both versions are coming from my repository at marstr.online . def acquire_elevation_data(self): mstr_msg("xp_scenery", "[X-Plane] Acquiring DEM model data") url = "https://marstr.online/dem/" url = url + self._demfn urllib.request.urlretrieve(url, mstr_datafolder + "_cache/" + self._demfn) mstr_msg("xp_scenery", "[X-Plane] DEM data acquired") # Download the X-Plane definition file from my server def acquire_xes_data(self): mstr_msg("xp_scenery", "[X-Plane] Acquiring XES file") url = "https://marstr.online/xes/" xesfn = self.build_dem_filename(True) xp_folder = self.xplane_latlng_folder([self._lat, self._lng]) url = url + xp_folder + ".xes" urllib.request.urlretrieve(url, mstr_datafolder + "_cache/" + xesfn) mstr_msg("xp_scenery", "[X-Plane] XES data acquired") # This generates all .ter files def build_ter_files(self): mstr_msg("xp_scenery", "[X-Plane] Generating and writing terrain (.ter) files") cur_lat = self._lat cur_lng = self._lng xp_folder = self.xplane_latlng_folder([self._lat, self._lng]) for lat in range(1, self._mlat+1): for lng in range(1, self._mlng+1): wdt = self.find_width_of_longitude(cur_lat) dmt = wdt * mstr_zl_18 cnt_x = cur_lat + (self._vstep/2) cnt_y = cur_lng + (mstr_zl_18/2) ddsf_water = mstr_datafolder + "z_orthographic/orthos/" + self._latlngfld + "/" + str(lat) + "_" + str(lng) + "_water.png" terstr = "" terstr = terstr + "A\r\n" terstr = terstr + "800\r\n" terstr = terstr + "TERRAIN\r\n" terstr = terstr + "\r\n" terstr = terstr + "LOAD_CENTER " + str(cnt_x) + " " + str(cnt_y) + " " + str(dmt) + " 2048\r\n" terstr = terstr + "TEXTURE_NOWRAP ../../orthos/" + self._latlngfld + "/" + str(lat)+"_"+str(lng)+".dds\r\n" if mstr_xp_scn_normalmaps == True: terstr = terstr + "NORMAL_TEX 1.0 ../../normals/" + self._latlngfld + "/" + str(lat)+"_"+str(lng)+".png\r\n" if os.path.isfile(ddsf_water) == True: terstr = terstr + "BORDER_TEX ../../orthos/" + self._latlngfld + "/" + str(lat) + "_" + str(lng) + "_water.png\r\n" terfln = mstr_datafolder + "z_orthographic/terrain/" + self._latlngfld + "/" + str(lat)+"_"+str(lng)+".ter" with open(terfln, 'w') as textfile: textfile.write(terstr) cur_lng = round(cur_lng + mstr_zl_18, 6) cur_lng = self._lng cur_lat = round(cur_lat + self._vstep, 6) mstr_msg("xp_scenery", "[X-Plane] Terrain files written") # This generates the entire terrain mesh def generate_terrain_mesh(self): # Get the DEM model file name, and acquire important info about the data meshfn = mstr_datafolder + "_cache/" + self.build_dem_filename() siz = os.path.getsize(meshfn) dim = int(math.sqrt(siz/2)) assert dim*dim*2 == siz, 'Invalid file size' self._demdata = numpy.fromfile(meshfn, numpy.dtype('>i2'), dim*dim).reshape((dim, dim)) self._demdata = self._demdata[::-1] # Invert order so that we can start from bottom left # We want to achieve perfect stepping for each data point in the DEM. demstep = round( 1 / len(self._demdata), 6) # Generate an array which contains only the coordinates self._demcoord = [] for r in range(0, len(self._demdata)): row = [] for c in range(0, len(self._demdata)): lat = round(self._lat + r * demstep, 6) lng = round(self._lng + c * demstep, 6) crd = [ lat, lng, self._demdata[r][c]] #crd = [ lat, lng ] row.append(crd) self._demcoord.append(row) mstr_msg("xp_scenery", "[X-Plane] Populating DSF information file") # The complete string to write into the DSF txt file dsf_str = "" dsf_str = dsf_str + "PROPERTY sim/west " + str(int(self._lng)) + "\r\n" dsf_str = dsf_str + "PROPERTY sim/east " + str((int(self._lng) + 1)) + "\r\n" dsf_str = dsf_str + "PROPERTY sim/south " + str(int(self._lat)) + "\r\n" dsf_str = dsf_str + "PROPERTY sim/north " + str((int(self._lat) + 1)) + "\r\n" dsf_str = dsf_str + "PROPERTY sim/require_object 0/6\r\n" dsf_str = dsf_str + "PROPERTY planet earth\r\n" dsf_str = dsf_str + "PROPERTY sim/creation_agent Orthographic\r\n" #dsf_str = dsf_str + "TERRAIN_DEF terrain_Water\r\n" # The file to be converted into DSF later meshtxt = mstr_datafolder + "_cache/mesh_"+self._latlngfld+".txt" with open(meshtxt, 'w') as textfile: textfile.write(dsf_str) for lat in range(1, self._mlat+1): for lng in range(1, self._mlng+1): # Write the line only if an ortho exists of course. ddsf = mstr_datafolder + "z_orthographic/orthos/" + self._latlngfld + "/" + str(lat) + "_" + str(lng) + ".dds" ddsf_water = mstr_datafolder + "z_orthographic/orthos/" + self._latlngfld + "/" + str(lat) + "_" + str(lng) + "_water.png" if os.path.isfile(ddsf) == True: dsfstr = "TERRAIN_DEF terrain/" + self._latlngfld + "/" + str(lat) + "_" + str(lng) + ".ter\r\n" # Let's check if this tile needs water beneath if os.path.isfile(ddsf_water) == True: dsfstr = dsfstr + "TERRAIN_DEF terrain_Water\r\n" with open(meshtxt, 'a') as textfile: textfile.write(dsfstr) # OK. So. Let's build the mesh. """ # First, the ground water mesh with open(meshtxt, 'a') as textfile: textfile.write("BEGIN_PATCH 0 0.000000 -1.000000 1 5\r\n") # Vertical row (Latitude Row) for lat_r in range(0, len(self._demcoord)-2): # Horizontal row (Longitude Column) for lng_c in range(0, len(self._demcoord)-2): # Lat/lng coordinate lat_crd = self._demcoord[lat_r][lng_c][0] lat_crd_t = self._demcoord[lat_r+1][lng_c][0] lng_crd = self._demcoord[lat_r][lng_c][1] lng_crd_r = self._demcoord[lat_r][lng_c+1][1] # Coords of triangle vertices # 0 - Longitude # 1 - Latitude # 2 - Height in m t1_v1 = [ lng_crd_r, lat_crd, 0 ] t1_v2 = [ lng_crd, lat_crd_t, 0 ] t1_v3 = [ lng_crd_r, lat_crd_t, 0 ] t2_v1 = [ lng_crd, lat_crd_t, 0 ] t2_v2 = [ lng_crd_r, lat_crd, 0 ] t2_v3 = [ lng_crd, lat_crd, 0 ] t1_v1 = [ lng_crd_r, lat_crd, self._demcoord[lat_r][lng_c+1][2] ] t1_v2 = [ lng_crd, lat_crd_t, self._demcoord[lat_r+1][lng_c][2] ] t1_v3 = [ lng_crd_r, lat_crd_t, self._demcoord[lat_r+1][lng_c+1][2] ] t2_v1 = [ lng_crd, lat_crd_t, self._demcoord[lat_r+1][lng_c][2] ] t2_v2 = [ lng_crd_r, lat_crd, self._demcoord[lat_r][lng_c+1][2] ] t2_v3 = [ lng_crd, lat_crd, self._demcoord[lat_r][lng_c][2] ] # Write down the two triangles t_str = "" t_str = t_str + "BEGIN_PRIMITIVE 0\r\n" t_str = t_str + "PATCH_VERTEX " + str(t1_v1[0]) + " " + str(t1_v1[1]) + " " + str(t1_v1[2]) + " 0.000015 0.000015\r\n" t_str = t_str + "PATCH_VERTEX " + str(t1_v2[0]) + " " + str(t1_v2[1]) + " " + str(t1_v2[2]) + " 0.000015 0.000015\r\n" t_str = t_str + "PATCH_VERTEX " + str(t1_v3[0]) + " " + str(t1_v3[1]) + " " + str(t1_v3[2]) + " 0.000015 0.000015\r\n" t_str = t_str + "END_PRIMITIVE 0\r\n" t_str = t_str + "BEGIN_PRIMITIVE 0\r\n" t_str = t_str + "PATCH_VERTEX " + str(t2_v1[0]) + " " + str(t2_v1[1]) + " " + str(t2_v1[2]) + " 0.000015 0.000015\r\n" t_str = t_str + "PATCH_VERTEX " + str(t2_v2[0]) + " " + str(t2_v2[1]) + " " + str(t2_v2[2]) + " 0.000015 0.000015\r\n" t_str = t_str + "PATCH_VERTEX " + str(t2_v3[0]) + " " + str(t2_v3[1]) + " " + str(t2_v3[2]) + " 0.000015 0.000015\r\n" t_str = t_str + "END_PRIMITIVE 0\r\n" # Send to the file with open(meshtxt, 'a') as textfile: textfile.write(t_str) t_str = "" # Water mesh ends with open(meshtxt, 'a') as textfile: textfile.write("END PATCH\r\n") """ # Current patch curpatch = 0 for lat in range(1, self._mlat+1): for lng in range(1, self._mlng+1): # Create the patch only if the matching ortho exists. # This way we make sure that we hit the same order as the .ter files. # We can also detect which lat and lng coord we are on. #ddsf = mstr_datafolder + "z_orthographic/orthos/" + self._latlngfld + "/1_1.dds" ddsf = mstr_datafolder + "z_orthographic/orthos/" + self._latlngfld + "/" + str(lat) + "_" + str(lng) + ".dds" ddsf_water = mstr_datafolder + "z_orthographic/orthos/" + self._latlngfld + "/" + str(lat) + "_" + str(lng) + "_water.png" if os.path.isfile(ddsf) == True: scangrid = self.find_height_scan_start_end_points([ self._lat+((lat-1)*self._vstep), self._lng+((lng-1)*mstr_zl_18) ]) #sloped = self.build_sloped_scangrid(scangrid) # Base coords for this ortho base_lat = self._lat + ((lat-1) * self._vstep) base_lng = self._lng + ((lng-1) * mstr_zl_18) # Begin a new patch mstr_msg("xp_scenery", "[X-Plane] Processing ortho patch " + str(curpatch)) with open(meshtxt, 'a') as textfile: textfile.write("BEGIN_PATCH " + str(curpatch) + " 0.000000 -1.000000 1 7\r\n") # Step for each ortho vertex odiv = 4 latstep = self._vstep/odiv lngstep = mstr_zl_18 /odiv uv_step = 1 / odiv # Height values hgt_bl = 0 hgt_br = 0 hgt_tr = 0 hgt_tl = 0 # Generate the ortho tile for y in range(0,odiv): for x in range(0,odiv): # Coordinates lat_b = round(base_lat + (y*latstep), 6) lat_t = round(base_lat + ((y+1)*latstep), 6) lng_l = round(base_lng + (x*lngstep), 6) lng_r = round(base_lng + ((x+1)*lngstep), 6) # Minimal adjustment if x == 0: lng_l = base_lng if y == 0: lat_b = base_lat if y == 3: lat_t = base_lat + self._vstep if x == 3: lng_r = base_lng + mstr_zl_18 # Corrections, just in case if lat_b > self._lat + 1: lat_b = self._lat+1 if lat_t > self._lat + 1: lat_t = self._lat+1 if lng_l > self._lng + 1: lng_l = self._lng+1 if lng_r > self._lng + 1: lng_r = self._lng+1 # Height indexes hgt_bl_idx = self.find_height_for_coord([lat_b, lng_l]) hgt_br_idx = self.find_height_for_coord([lat_b, lng_r]) hgt_tr_idx = self.find_height_for_coord([lat_t, lng_r]) hgt_tl_idx = self.find_height_for_coord([lat_t, lng_l]) hgt_bl = round(self._demcoord[ hgt_bl_idx[0] ][ hgt_bl_idx[1] ][2], 6) hgt_br = round(self._demcoord[ hgt_br_idx[0] ][ hgt_br_idx[1] ][2], 6) hgt_tr = round(self._demcoord[ hgt_tr_idx[0] ][ hgt_tr_idx[1] ][2], 6) hgt_tl = round(self._demcoord[ hgt_tl_idx[0] ][ hgt_tl_idx[1] ][2], 6) # Coords of triangle vertices # 0 - Longitude # 1 - Latitude # 2 - Height in m t1_v1 = [ lng_r, lat_b, hgt_br ] t1_v2 = [ lng_l, lat_t, hgt_tl ] t1_v3 = [ lng_r, lat_t, hgt_tr ] t2_v1 = [ lng_l, lat_t, hgt_tl ] t2_v2 = [ lng_r, lat_b, hgt_br ] t2_v3 = [ lng_l, lat_b, hgt_bl ] # Write down the two triangles t_str = "" t_str = t_str + "BEGIN_PRIMITIVE 0\r\n" t_str = t_str + "PATCH_VERTEX " + str(t1_v1[0]) + " " + str(t1_v1[1]) + " " + str(t1_v1[2]) + " 0.000015 0.000015 " + str((x+1) * uv_step) + " " + str(y*uv_step) + "\r\n" t_str = t_str + "PATCH_VERTEX " + str(t1_v2[0]) + " " + str(t1_v2[1]) + " " + str(t1_v2[2]) + " 0.000015 0.000015 " + str(x * uv_step) + " " + str((y+1)*uv_step) + "\r\n" t_str = t_str + "PATCH_VERTEX " + str(t1_v3[0]) + " " + str(t1_v3[1]) + " " + str(t1_v3[2]) + " 0.000015 0.000015 " + str((x+1) * uv_step) + " " + str((y+1)*uv_step) + "\r\n" t_str = t_str + "END_PRIMITIVE 0\r\n" t_str = t_str + "BEGIN_PRIMITIVE 0\r\n" t_str = t_str + "PATCH_VERTEX " + str(t2_v1[0]) + " " + str(t2_v1[1]) + " " + str(t2_v1[2]) + " 0.000015 0.000015 " + str(x * uv_step) + " " + str((y+1)*uv_step) + "\r\n" t_str = t_str + "PATCH_VERTEX " + str(t2_v2[0]) + " " + str(t2_v2[1]) + " " + str(t2_v2[2]) + " 0.000015 0.000015 " + str((x+1) * uv_step) + " " + str(y*uv_step) + "\r\n" t_str = t_str + "PATCH_VERTEX " + str(t2_v3[0]) + " " + str(t2_v3[1]) + " " + str(t2_v3[2]) + " 0.000015 0.000015 " + str(x * uv_step) + " " + str(y*uv_step) + "\r\n" t_str = t_str + "END_PRIMITIVE 0\r\n" # Send to the file with open(meshtxt, 'a') as textfile: textfile.write(t_str) t_str = "" # Height value: # hgtindex = self.find_height_for_coord([lat, lng]) # height = self._demcoord[hgtindex[0]][hgtindex[1]][2] # End this patch with open(meshtxt, 'a') as textfile: textfile.write("END PATCH\r\n") # Increase patch number curpatch = curpatch + 1 # Let's check if this tile needs water beneath needs_water = False if os.path.isfile(ddsf_water) == True: needs_water = True if needs_water == True: # Begin a new patch with open(meshtxt, 'a') as textfile: textfile.write("BEGIN_PATCH " + str(curpatch) + " 0.000000 -1.000000 1 5\r\n") # Generate the ortho tile for y in range(0,odiv): for x in range(0,odiv): # Coordinates lat_b = round(base_lat + (y*latstep), 6) lat_t = round(base_lat + ((y+1)*latstep), 6) lng_l = round(base_lng + (x*lngstep), 6) lng_r = round(base_lng + ((x+1)*lngstep), 6) # Minimal adjustment if x == 0: lng_l = base_lng if y == 0: lat_b = base_lat if y == 3: lat_t = base_lat + self._vstep if x == 3: lng_r = base_lng + mstr_zl_18 # Corrections, just in case if lat_b > self._lat + 1: lat_b = self._lat+1 if lat_t > self._lat + 1: lat_t = self._lat+1 if lng_l > self._lng + 1: lng_l = self._lng+1 if lng_r > self._lng + 1: lng_r = self._lng+1 hgt_bl_idx = self.find_height_for_coord([lat_b, lng_l]) hgt_br_idx = self.find_height_for_coord([lat_b, lng_r]) hgt_tr_idx = self.find_height_for_coord([lat_t, lng_r]) hgt_tl_idx = self.find_height_for_coord([lat_t, lng_l]) hgt_bl = round(self._demcoord[ hgt_bl_idx[0] ][ hgt_bl_idx[1] ][2] - .1, 6) hgt_br = round(self._demcoord[ hgt_br_idx[0] ][ hgt_br_idx[1] ][2] - .1, 6) hgt_tr = round(self._demcoord[ hgt_tr_idx[0] ][ hgt_tr_idx[1] ][2] - .1, 6) hgt_tl = round(self._demcoord[ hgt_tl_idx[0] ][ hgt_tl_idx[1] ][2] - .1, 6) # Coords of triangle vertices # 0 - Longitude # 1 - Latitude # 2 - Height in m t1_v1 = [ lng_r, lat_b, hgt_br ] t1_v2 = [ lng_l, lat_t, hgt_tl ] t1_v3 = [ lng_r, lat_t, hgt_tr ] t2_v1 = [ lng_l, lat_t, hgt_tl ] t2_v2 = [ lng_r, lat_b, hgt_br ] t2_v3 = [ lng_l, lat_b, hgt_bl ] # Write down the two triangles t_str = "" t_str = t_str + "BEGIN_PRIMITIVE 0\r\n" t_str = t_str + "PATCH_VERTEX " + str(t1_v1[0]) + " " + str(t1_v1[1]) + " " + str(t1_v1[2]) + " 0.000015 0.000015\r\n" t_str = t_str + "PATCH_VERTEX " + str(t1_v2[0]) + " " + str(t1_v2[1]) + " " + str(t1_v2[2]) + " 0.000015 0.000015\r\n" t_str = t_str + "PATCH_VERTEX " + str(t1_v3[0]) + " " + str(t1_v3[1]) + " " + str(t1_v3[2]) + " 0.000015 0.000015\r\n" t_str = t_str + "END_PRIMITIVE 0\r\n" t_str = t_str + "BEGIN_PRIMITIVE 0\r\n" t_str = t_str + "PATCH_VERTEX " + str(t2_v1[0]) + " " + str(t2_v1[1]) + " " + str(t2_v1[2]) + " 0.000015 0.000015\r\n" t_str = t_str + "PATCH_VERTEX " + str(t2_v2[0]) + " " + str(t2_v2[1]) + " " + str(t2_v2[2]) + " 0.000015 0.000015\r\n" t_str = t_str + "PATCH_VERTEX " + str(t2_v3[0]) + " " + str(t2_v3[1]) + " " + str(t2_v3[2]) + " 0.000015 0.000015\r\n" t_str = t_str + "END_PRIMITIVE 0\r\n" # Send to the file with open(meshtxt, 'a') as textfile: textfile.write(t_str) # End this patch with open(meshtxt, 'a') as textfile: textfile.write("END PATCH\r\n") # Increase patch number curpatch = curpatch + 1 # Find the next best matching height for a point def find_height_for_coord(self, coord): idx = [0,0] dst = 99999 ste = self.find_height_scan_start_end_points(coord) for r in range(ste[0], ste[1]+1): for d in range(ste[2], ste[3]+1): dist = math.dist(coord, [self._demcoord[r][d][0], self._demcoord[r][d][1]]) if dist < dst: dst = dist idx = [r,d] return idx # Find the starting and end points to scan for heights in the DEM grid def find_height_scan_start_end_points(self, stc): startend = [0,0,0,0] stp = 1 / len(self._demdata) # Bottom lt = self._lat while lt < stc[0]: lt = lt + stp startend[0] = startend[0] + 1 # Top lt = self._lat while lt < stc[0]+self._vstep: lt = lt+stp startend[1] = startend[1] + 1 # Left ln = self._lng while ln < stc[1]: ln = ln + stp startend[2] = startend[2] + 1 # Right ln = self._lng while ln < stc[1]+mstr_zl_18: ln = ln + stp startend[3] = startend[3] + 1 # Make sure we have everything startend[0] = startend[0]-1 startend[1] = startend[1]+1 startend[2] = startend[2]-1 startend[3] = startend[3]+1 # Some corrections if startend[0] < 0: startend[0] = 0 if startend[1] > len(self._demdata)-1: startend[1] = startend[1] = len(self._demdata)-1 if startend[2] < 0: startend[2] = 0 if startend[3] > len(self._demdata)-1: startend[3] = startend[3] = len(self._demdata)-1 """ t = self._lat while t < startcoord[0]: t = t + self._vstep startend[0] = startend[0]+1 t = self._lat while t < startcoord[0]+self._vstep: t = t + self._vstep startend[1] = startend[1]+1 t = self._lng while t < startcoord[1]: t = t + mstr_zl_18 startend[2] = startend[2]+1 t = self._lng while t < startcoord[1]+mstr_zl_18: t = t + mstr_zl_18 startend[3] = startend[3]+1 # Some corrections startend[0] = startend[0]-1 if startend[0] < 0: startend[0] = 0 startend[1] = startend[1]+1 if startend[1] > len(self._demdata)-1: startend[1] = startend[1] = len(self._demdata)-1 startend[2] = startend[2]-1 if startend[2] < 0: startend[2] = 0 startend[3] = startend[3]+1 if startend[3] > len(self._demdata)-1: startend[3] = startend[3] = len(self._demdata)-1 """ return startend # Function to subdivide between two vectors def subdivide_vectors(self, v1, v2, subdivisions): #return np.linspace(v1, v2, subdivisions + 2, axis=0) # +2 to include endpoints return numpy.linspace(v1, v2, subdivisions + 2, axis=0) # +2 to include endpoints # This build a scangrid with increased resolution, extrapolated from existing points. # With this we can accurately depict the height of a point a long the slope of the ground mesh. def build_sloped_scangrid(self, grid): # Contains the data as defined by the grid passed in tmp_dem = [] # Acquire original grid data for l in range(grid[2], grid[3]+1): row = [] for c in range(grid[0], grid[1]+1): row.append(self._demcoord[l][c]) tmp_dem.append(row) # Subdivide the array subdivisions = 10 result = [] for i in range(len(tmp_dem) - 1): for j in range(len(tmp_dem[i]) - 1): subdivided = self.subdivide_vectors(tmp_dem[i][j], tmp_dem[i][j + 1], subdivisions) if i > 0: # Avoid duplicating start point subdivided = subdivided[1:] result.append(subdivided) # Combine all subdivisions into one array result = numpy.vstack(result) return result