Commit 54b7dd2e authored by Larkin Heintzman's avatar Larkin Heintzman

fixed terrain interpolation error

parent add8abf8
...@@ -184,6 +184,7 @@ def get_terrain_map(lat_lon = [0,0], sample_dist = 10, extent = 100, heading = 0 ...@@ -184,6 +184,7 @@ def get_terrain_map(lat_lon = [0,0], sample_dist = 10, extent = 100, heading = 0
for scx in range(sc): for scx in range(sc):
for scy in range(sc): # maybe the least efficient way to do this for scy in range(sc): # maybe the least efficient way to do this
idx = np.where(np.all(np.abs(qs - np.asarray([X[scx,scy],Y[scx,scy]])) <= 1e-9,axis = 1)) idx = np.where(np.all(np.abs(qs - np.asarray([X[scx,scy],Y[scx,scy]])) <= 1e-9,axis = 1))
# print(idx)
try: try:
es_square[scx,scy] = es[idx] es_square[scx,scy] = es[idx]
except ValueError as e: except ValueError as e:
...@@ -251,7 +252,8 @@ if __name__ == "__main__": ...@@ -251,7 +252,8 @@ if __name__ == "__main__":
scaled_extent = np.ceil(scale_factor*extent).astype(np.int) scaled_extent = np.ceil(scale_factor*extent).astype(np.int)
# interpolate terrain to match size/resolution of other layers # interpolate terrain to match size/resolution of other layers
f = interpolate.interp2d(np.arange(0, extent, sample_dist), np.arange(0, extent, sample_dist), e, kind='cubic') c = np.int(extent/sample_dist)
f = interpolate.interp2d(np.linspace(0, extent, c), np.linspace(0, extent, c), e, kind='cubic')
x_temp = np.linspace(0,extent,scaled_extent) # get correct size of terrain map x_temp = np.linspace(0,extent,scaled_extent) # get correct size of terrain map
y_temp = np.linspace(0,extent,scaled_extent) y_temp = np.linspace(0,extent,scaled_extent)
e_interp = f(x_temp, y_temp) e_interp = f(x_temp, y_temp)
...@@ -260,6 +262,10 @@ if __name__ == "__main__": ...@@ -260,6 +262,10 @@ if __name__ == "__main__":
plt.show() plt.show()
plt.imshow(e) plt.imshow(e)
plt.show() plt.show()
plt.imshow(x)
plt.show()
plt.imshow(y)
plt.show()
# elv_filename = "map_layers\\elv_data_"+file_id+".csv" # elv_filename = "map_layers\\elv_data_"+file_id+".csv"
# if save_files: # if save_files:
......
...@@ -38,10 +38,6 @@ def grab_features(anchor_point, extent, sample_dist = 10, heading = 0, save_file ...@@ -38,10 +38,6 @@ def grab_features(anchor_point, extent, sample_dist = 10, heading = 0, save_file
name_list = ['rivers_bdd', 'rivers', 'roads', 'lakes', 'powerlines', 'railroads', 'trails'] name_list = ['rivers_bdd', 'rivers', 'roads', 'lakes', 'powerlines', 'railroads', 'trails']
inac_layers = ['rivers_bdd', 'lakes'] inac_layers = ['rivers_bdd', 'lakes']
# url_list = [riverw_url, river_url]
# name_list = ['rivers_bdd', 'rivers','bleh','bleh']
# inac_layers = ['rivers_bdd']
gis = GIS(username="larkinheintzman",password="Meepp97#26640") # linked my arcgis pro account gis = GIS(username="larkinheintzman",password="Meepp97#26640") # linked my arcgis pro account
ap_meters = lat_lon2meters(anchor_point[0], anchor_point[1]) ap_meters = lat_lon2meters(anchor_point[0], anchor_point[1])
...@@ -79,6 +75,9 @@ def grab_features(anchor_point, extent, sample_dist = 10, heading = 0, save_file ...@@ -79,6 +75,9 @@ def grab_features(anchor_point, extent, sample_dist = 10, heading = 0, save_file
query_cnt = query_cnt + 1 query_cnt = query_cnt + 1
print("error on query: {}".format(e)) print("error on query: {}".format(e))
print("{} layer failed on query, trying again ...".format(name_list[i])) print("{} layer failed on query, trying again ...".format(name_list[i]))
gis = GIS(username="larkinheintzman",password="Meepp97#26640") # linked my arcgis pro account
lyr = FeatureLayer(url=url, gis=gis)
if query_cnt > 2 and not q: if query_cnt > 2 and not q:
print("{} layer failed too many times, leaving empty".format(name_list[i])) print("{} layer failed too many times, leaving empty".format(name_list[i]))
if save_files: if save_files:
...@@ -293,7 +292,8 @@ def grab_features(anchor_point, extent, sample_dist = 10, heading = 0, save_file ...@@ -293,7 +292,8 @@ def grab_features(anchor_point, extent, sample_dist = 10, heading = 0, save_file
e = np.flipud(e) e = np.flipud(e)
# interpolate terrain to match size/resolution of other layers # interpolate terrain to match size/resolution of other layers
f = interpolate.interp2d(np.arange(0, extent, sample_dist), np.arange(0, extent, sample_dist), e, kind='cubic') c = np.int(extent/sample_dist)
f = interpolate.interp2d(np.linspace(0, extent, c), np.linspace(0, extent, c), e, kind='cubic')
x_temp = np.linspace(0,extent,scaled_extent) # get correct size of terrain map x_temp = np.linspace(0,extent,scaled_extent) # get correct size of terrain map
y_temp = np.linspace(0,extent,scaled_extent) y_temp = np.linspace(0,extent,scaled_extent)
e_interp = f(x_temp, y_temp) e_interp = f(x_temp, y_temp)
......
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