Commit 1ccfab24 authored by Bryson Howell's avatar Bryson Howell

bringing in changes from unity_test branch on original repo.

parent be716c01
......@@ -104,9 +104,11 @@ def get_terrain_map(lat_lon = [0,0], sample_dist = 10, extent = 100, heading = 0
# gis = GIS("pro")
# gis = GIS(url="http://virginiatech.maps.arcgis.com", client_id="rluxzSWjZS6TfeXs", username="hlarkin3_virginiatech", password="arcgisheintzman97#26640", verify_cert=False)
gis = GIS(username="larkinheintzman",password="Meepp97#26640")
#gis = GIS(username="larkinheintzman",password="Meepp97#26640")
gis = GIS(api_key="AAPKe7cb1ab4f2ba44748cf53ac4f30d0caavu5t_-uMdWP8SRjtoea3s66-cRyFZmMTNx4rqy2w5sjwegU_tbJyBesd0LCMmdtV")
if verbosity:
print("Successfully logged in as: " + gis.properties.user.username)
#print("Successfully logged in as: " + gis.properties.user.username)
print("Logged into GIS")
elv_map = gis.content.get('58a541efc59545e6b7137f961d7de883')
elv_layer = elv_map.layers[0]
......@@ -212,18 +214,21 @@ def get_terrain_map(lat_lon = [0,0], sample_dist = 10, extent = 100, heading = 0
# flip elevation data up to down to match other layers
e = np.flipud(e)
# interpolate terrain to match size/resolution of other layers
c = np.int(extent/sample_dist)
c = np.int32(extent/sample_dist)
scale_factor = 3/20 # factor to get 6.66667m mapping from 1m mapping (1/6.6667)
scaled_extent = np.ceil(scale_factor*extent).astype(np.int)
factor = scaled_extent/e.shape[0]
e_interp = ndimage.zoom(e, factor, order = 3)
print(show_plot)
if show_plot:
# Attaching 3D axis to the figure
grid_x, grid_y = np.mgrid[min(x):max(x):100j, min(y):max(y):100j]
grid_x, grid_y = np.mgrid[np.min(x):np.max(x):100j, np.min(y):np.max(y):100j]
points = np.array([x, y])
grid_z = griddata(points.transpose(), z, (grid_x, grid_y), method='cubic', fill_value=np.mean(z))
print(np.shape(points))
#grid_z = griddata(points.transpose(), e, (grid_x, grid_y), method='cubic', fill_value=np.mean(e))
grid_z = griddata(points.transpose(), e, (grid_x, grid_y), method='cubic', fill_value=np.mean(e))
print(grid_z)
fig = plt.figure()
ax = fig.gca(projection='3d')
......@@ -255,7 +260,9 @@ if __name__ == "__main__":
[e,e_interp,x,y,data,ll_pt] = get_terrain_map(lat_lon=anchor_point,
sample_dist = sample_dist,
extent = extent,
heading = -heading) # because flipping
heading = -heading,
show_plot = True,
verbosity = True) # because flipping
plt.imshow(e_interp)
plt.title('e interp')
......@@ -296,4 +303,3 @@ if __name__ == "__main__":
# elv_filename = "map_layers\\elv_data_"+file_id+".csv"
# if save_files:
# np.savetxt(elv_filename,e_interp,delimiter=",", fmt='%f')
......@@ -46,12 +46,12 @@ def grab_features(anchor_point, extent, sample_dist = 10, case_name = 'blah', he
name_list = ['rivers_bdd', 'rivers', 'roads', 'lakes', 'powerlines', 'railroads', 'trails']
inac_layers = ['rivers_bdd', 'lakes']
gis = GIS(username="larkinheintzman",password="Meepp97#26640") # linked my arcgis pro account
gis = GIS(api_key="AAPKe7cb1ab4f2ba44748cf53ac4f30d0caavu5t_-uMdWP8SRjtoea3s66-cRyFZmMTNx4rqy2w5sjwegU_tbJyBesd0LCMmdtV") # linked Bryson's pro account
ap_meters = lat_lon2meters(anchor_point[0], anchor_point[1])
scale_factor = 3/20 # factor to get 6.66667m mapping from 1m mapping (1/6.6667)
scaled_extent = np.ceil(scale_factor*extent).astype(np.int)
scaled_extent = np.ceil(scale_factor*extent).astype(np.int32)
viz_cnt = 0
viz_map = np.zeros([scaled_extent,scaled_extent,len(name_list)+len(inac_layers)])
......@@ -90,7 +90,7 @@ def grab_features(anchor_point, extent, sample_dist = 10, case_name = 'blah', he
print("error on query: {}".format(e))
print("{} layer failed on query, trying again ...".format(name_list[i]))
# gis.
gis = GIS(username="larkinheintzman",password="Meepp97#26640") # linked my arcgis pro account
gis = GIS(api_key="AAPKe7cb1ab4f2ba44748cf53ac4f30d0caavu5t_-uMdWP8SRjtoea3s66-cRyFZmMTNx4rqy2w5sjwegU_tbJyBesd0LCMmdtV") # linked my arcgis pro account
lyr = FeatureLayer(url=url, gis=gis)
print("query time {}".format(query_endtime - query_starttime))
......@@ -99,16 +99,16 @@ def grab_features(anchor_point, extent, sample_dist = 10, case_name = 'blah', he
print("{} layer failed too many times, leaving empty".format(name_list[i]))
if save_files:
if save_to_folder:
fn = "map_layers\\" + case_name + "\\"+name_list[i]+"_data_"+file_id+".csv"
fn = "./map_layers/" + case_name + "/"+name_list[i]+"_data_"+file_id+".csv"
np.savetxt(fn,bin_map,delimiter=",", fmt='%f')
if name_list[i] in inac_layers:
fn = "map_layers\\" + case_name + "\\" + name_list[i] + "_inac_data_" + file_id + ".csv"
fn = "./map_layers/" + case_name + "/" + name_list[i] + "_inac_data_" + file_id + ".csv"
np.savetxt(fn, bin_map, delimiter=",", fmt='%f')
else:
fn = "map_layers\\" + name_list[i]+"_data_"+file_id+".csv"
fn = "./map_layers/" + name_list[i]+"_data_"+file_id+".csv"
np.savetxt(fn,bin_map,delimiter=",", fmt='%f')
if name_list[i] in inac_layers:
fn = "map_layers\\" + name_list[i] + "_inac_data_" + file_id + ".csv"
fn = "./map_layers/" + name_list[i] + "_inac_data_" + file_id + ".csv"
np.savetxt(fn, bin_map, delimiter=",", fmt='%f')
continue
print("{} layer sucessfully queried".format(name_list[i]))
......@@ -213,7 +213,7 @@ def grab_features(anchor_point, extent, sample_dist = 10, case_name = 'blah', he
pts_inac = np.delete(pts_inac, np.where(np.equal(pt,pts_inac).all(1)), axis = 0)
# binarization step
pts_inac = np.round(pts_inac).astype(np.int)
pts_inac = np.round(pts_inac).astype(np.int32)
# flip y axis
pts_inac[:,1] = inac_bin_map.shape[1] - pts_inac[:,1]
# remove any points outside limits of binary map (fixes round versus ceil issues)
......@@ -226,8 +226,8 @@ def grab_features(anchor_point, extent, sample_dist = 10, case_name = 'blah', he
# trim features to scaled extent
[x_pts, y_pts] = trim_extent(x_pts, y_pts, scaled_extent)
# binarization step
x_pts_idx = np.round(x_pts).astype(np.int)
y_pts_idx = np.round(y_pts).astype(np.int)
x_pts_idx = np.round(x_pts).astype(np.int32)
y_pts_idx = np.round(y_pts).astype(np.int32)
# flip y axis because indices are fliped
y_pts_idx = bin_map.shape[1] - y_pts_idx
......@@ -246,10 +246,10 @@ def grab_features(anchor_point, extent, sample_dist = 10, case_name = 'blah', he
viz_cnt = viz_cnt + 1
if save_files:
if save_to_folder:
fn = "map_layers\\" + case_name + "\\" + name_list[i] + "_inac_data_" + file_id + ".csv"
fn = "./map_layers/" + case_name + "/" + name_list[i] + "_inac_data_" + file_id + ".csv"
np.savetxt(fn, inac_bin_map, delimiter=",", fmt='%f')
else:
fn = "map_layers\\" + name_list[i] + "_inac_data_" + file_id + ".csv"
fn = "./map_layers/" + name_list[i] + "_inac_data_" + file_id + ".csv"
np.savetxt(fn, inac_bin_map, delimiter=",", fmt='%f')
......@@ -257,10 +257,10 @@ def grab_features(anchor_point, extent, sample_dist = 10, case_name = 'blah', he
viz_cnt = viz_cnt + 1
if save_files:
if save_to_folder:
fn = "map_layers\\" + case_name + "\\" + name_list[i]+"_data_"+file_id+".csv"
fn = "./map_layers/" + case_name + "/" + name_list[i]+"_data_"+file_id+".csv"
np.savetxt(fn,bin_map,delimiter=",", fmt='%f')
else:
fn = "map_layers\\" + name_list[i]+"_data_"+file_id+".csv"
fn = "./map_layers/" + name_list[i]+"_data_"+file_id+".csv"
np.savetxt(fn,bin_map,delimiter=",", fmt='%f')
......@@ -272,10 +272,10 @@ def grab_features(anchor_point, extent, sample_dist = 10, case_name = 'blah', he
if save_files:
if save_to_folder:
elv_filename = "map_layers\\" + case_name + "\\elv_data_" + file_id + ".csv"
elv_filename = "./map_layers/" + case_name + "/elv_data_" + file_id + ".csv"
np.savetxt(elv_filename,e_interp,delimiter=",", fmt='%f')
else:
elv_filename = "map_layers\\elv_data_" + file_id + ".csv"
elv_filename = "./map_layers/elv_data_" + file_id + ".csv"
np.savetxt(elv_filename,e_interp,delimiter=",", fmt='%f')
if plot_data:
......@@ -406,7 +406,8 @@ if __name__ == "__main__":
# [48.64606, -122.4247,'[48.64606, -122.4247]']
]
base_dir = 'C:/Users/Larkin/ags_grabber'
#You should be running this script from the base dir
base_dir = './'
# for item in ics:
# try:
......@@ -427,14 +428,14 @@ if __name__ == "__main__":
anchor_point = [float(ics_pt[0]), float(ics_pt[1])]
extent = 10e3
save_flag = True
plot_flag = False
plot_flag = True
file_extension = 'temp'
sample_dist = int(extent/100)
heading = 0
start_time = time.time()
grab_features(anchor_point = anchor_point, extent = extent, sample_dist = sample_dist, case_name = str(ics_pt[2]),
heading = heading, save_files = save_flag, save_to_folder = False, file_id = file_extension, plot_data = plot_flag)
heading = heading, save_files = save_flag, save_to_folder = True, file_id = file_extension, plot_data = plot_flag)
time.sleep(1) # wait for... files to settle?
# run matlab
......@@ -442,4 +443,4 @@ if __name__ == "__main__":
res = eng.importmap_py(str(ics_pt[2]), base_dir)
print("------- total time = {} seconds, iteration {}/{} ------".format(time.time() - start_time,i,len(ics)))
eng.quit()
\ No newline at end of file
eng.quit()
import numpy as np
import matplotlib.pyplot as plt
from scipy.io import loadmat
from PIL import Image
from PIL import features
#import libtiff
import cv2
#Loads matlab matrices created by feature_set.py as numpy arrays
def main():
#Load the matlab matrix
dir = './matlab_data/'
savedir = './unity_data/'
filename = 'BW_LFandInac_Zelev_kentland.mat'
matdict = loadmat(dir+filename, appendmat=False)
#Testing stuff
k = matdict.keys()
print(k)
#These are square numpy arrays of equal size
elev = matdict['sZelev'] #Values are height of grid cell (in meters)
linfeat = matdict['BWLF'] #Values are 1 for linear feature present, 0 for none
obstacles = matdict['BWInac'] #Values are 1 for obstacle, 0 for clear
check_elev(elev)
make_featmap(linfeat,savedir)
make_heightmap(elev,savedir)
return
#Returns the difference between the min and max elevation.
#To be used as terrain height in Unity.
def check_elev(elev):
diff = np.max(elev)-np.min(elev)
print(diff)
print("Height of environment is %f\n" % diff)
return
#Saves image of linear features so it can be used as a mesh in Unity
def make_featmap(linfeat,savedir):
name = savedir + "kentland_features"
#Convert to png format image
print(np.shape(linfeat))
feat_int8 = (linfeat*255).astype(np.uint8)
img = Image.fromarray(feat_int8)
img.save('./unity_data/kentland_features.png')
return
#Converts the terrain height data in the .mat files into a grayscale heightmap
#that can be loaded into Unity. Heightmaps are 16-bit grayscale images saved in
#the RAW format.
def make_heightmap(elev,savedir):
name = savedir + "kentland_heightmap.raw"
#Makes sure libtiff is working right
print(features.check('libtiff'))
#fig = plt.figure()
#plt.imshow(elev)
#plt.show()
#Normalize values to be between 0-1 and preview as grayscale image
elev_gray = np.float16((elev-np.min(elev))/(np.max(elev)-np.min(elev)))
#elev_gray = np.ushort(elev)
#elev_gray = np.uint16(elev)
#From array expects uint8
fig = plt.figure()
#plt.imshow(elev_gray)
#plt.savefig('./unity_data/test.png',pad_inches=0)
#plt.show()
#f = open(name,'w')
#Convert to raw format image
elev_int8 = (elev_gray*255).astype(np.uint8)
img = Image.fromarray(elev_int8)
img.save('./unity_data/kentland_gimptest.png')
img = img.convert('L')
img.save('./unity_data/kentland_gray.raw', 'TIFF', compression='raw')
#img = img.convert('1')
#img = img.convert('I;16')
img = img.convert('1')
img.save(name, 'TIFF', compression='tiff_raw_16') #want tiff_raw_16
img.show()
#cv2.imshow('Kentland Heightmap',elev_gray)
#cv2.waitKey(0)
#cv2.imwrite('kentland_heightmap.raw',elev_gray)
return
if __name__ == '__main__':
main()
Markdown is supported
0% or
You are about to add 0 people to the discussion. Proceed with caution.
Finish editing this message first!
Please register or to comment