Raster/Elevation Layer Query Samples

import os
import globalmapper as gm
import numpy as np
import matplotlib.pyplot as plt

layer_path = "Z:/GM-20069/raster_elevation_to_numpy_sample_data/Depth_to_Water.tif"
export_dir = os.path.join(os.path.dirname(os.path.dirname(__file__)), "../misc/raster_elevation_to_numpy_output/")
os.makedirs(export_dir, exist_ok=True)

# This script places raster elevation values (in meters), and (x,y) pixel coordinate values into a python list from a
# defined area. The (x,y) pixel cooridnates are converted to (x,y) ground coorindate values based on the current project
# of the map. Ground coordinate values and elevation values are placed into a numpy array for use with Matplotlib.
# As an example of data visualization with external packages, a heat grid of the retrived grid of values is exported to a
# specified directory.

# Closes all layers that may be open
layers_are_open = gm.GetLoadedLayerList()
if layers_are_open:
    arr_ptr, arr_size = layers_are_open
    initial_layers = gm.GM_LayerHandle_array_frompointer(arr_ptr)
    for i in range(arr_size):
        gm.CloseLayer(initial_layers[i])

# Force load a raster layer with unknown interpretation as a grid layer
err_LoadLayerList, arr_ptr, arr_size = gm.LoadLayerList(
    layer_path, gm.GM_LoadFlags_LoadUnknownAsGrid
)
loaded_array = gm.GM_LayerHandle_array_frompointer(arr_ptr)

### START SAMPLE
# Displays information about the raster layer
info = gm.GetLayerInfo(loaded_array[0])
print("Layer %d, 0x%x:" % (0, loaded_array[0]))
print("    Description: %s" % info.mDescription)
print("    Type: %s" % info.mTypeName)
print("    Raster Data: %s" % info.mHasRasterData)
print("    Pixel Height: %s" % info.mGlobalPixelHeight)
print("    Pixel Width: %s" % info.mGlobalPixelWidth)
print("    Max Elevation: %s" % info.mMaxElevation)
print("    Min Elevation: %s" % info.mMinElevation)

# Collects elevation values and x,y pixel coordinates converted to ground
# coordinate values

# Python lists of elevation and x,y ground coordinate values
elevation_values = []
x_coord = []
y_coord = []
# i is defined as a starting y pixel coordinate
i = 3073
# i, j are definded for 10x10 pixel coordinate grid
while i < 3083:
    j = 6066
    while j < 6076:
        # Collects elevation values as specific x, y pixel coordinates
        err_GetPixelElevation, elevation = gm.GetPixelElevation(loaded_array[0], j, i)
        # X,Y pixel coordinates are converted to ground coordinates
        err_ConvertCoordLayerPixelToGround, coordinateX, coordinateY = (
            gm.ConvertCoordLayerPixelToGround(loaded_array[0], j, i)
        )
        # All values are added to their respective lists
        elevation_values.append(elevation)
        x_coord.append(coordinateX)
        y_coord.append(coordinateY)
        j += 1
    i += 1

# Numpy arrays are created for elevation and x, y ground coordinates
# Elevation is changed from a 1-d to 2-d array with 10 values per row
elevation = np.array(elevation_values).reshape(-1, 10)
# Duplicates are removed from x and y ground coordinate arrays
x_values = np.unique(np.array(x_coord))
y_values = np.unique(np.array(y_coord))
# Rounds ground coordinates to 3 decimal points
x = np.round(x_values, 3)
y = np.round(y_values, 3)
# Displays a simple heat grid with Matplotlib
fig, ax = plt.subplots()
# Creates image for heat grid
im = ax.imshow(elevation)

# Sets labels for x and y coordinate values
ax.set_xticks(np.arange(len(x)), labels=x, fontsize=10)
ax.set_yticks(np.arange(len(y)), labels=y, fontsize=10)
ax.set_ylabel("Y Ground Coordinate Values")
ax.set_xlabel("X Ground Coordinate Values")

# Rotate the tick labels and set their alignment
plt.setp(ax.get_xticklabels(), rotation=45, ha="right", rotation_mode="anchor")

# Loop over data dimensions and create text annotations of elevation values
for i in range(len(y)):
    for j in range(len(x)):
        text = ax.text(
            j, i, elevation[i, j], ha="center", va="center", color="w", fontsize=6
        )

# Displays title for grid
ax.set_title("Change of Elevation (meters) of a Selected Area")
fig.tight_layout()
plt.savefig(export_dir + "heat_grid.jpg")

### END SAMPLE

assert (
    gm.GM_Error_None == err_LoadLayerList
), f"LoadLayerList: Expected {gm.GM_Error_None}, got {err_LoadLayerList}"
assert (
    gm.GM_Error_None == err_GetPixelElevation
), f"GetPixelElevation: Expected {gm.GM_Error_None}, got {err_GetPixelElevation}"
assert (
    gm.GM_Error_None == err_ConvertCoordLayerPixelToGround
), f"ConvertCoordLayerPixelToGround: Expected {gm.GM_Error_None}, got {err_ConvertCoordLayerPixelToGround}"

assert os.path.exists(layer_path), f"Sample data doesn't exist: {layer_path}"

assert os.path.isdir(export_dir), f"Error creating the export directory: {export_dir}"

# Expected path for file export
exported_filepath = os.path.join(export_dir, "heat_grid.jpg")
# Assert that the final expected output file exists and is not empty
assert os.path.exists(
    exported_filepath
), f"Expected output file does not exist: {exported_filepath}"
assert (
    os.path.getsize(exported_filepath) > 0
), f"Exported file is empty: {exported_filepath}"

import os
import globalmapper as gm

import_dir = "Z:/GM-20069/single_function_sample_data/"

# The pixel elevations for a specified row are collected, and statistical data is displayed, including the row's 
# average elevation and the number of points falling within defined elevation ranges.

# Close all layers that may be open
layers_are_open = gm.GetLoadedLayerList()
if layers_are_open:
    arr_ptr, arr_size = layers_are_open
    initial_layers = gm.GM_LayerHandle_array_frompointer(arr_ptr)
    for i in range(arr_size):
        gm.CloseLayer(initial_layers[i])

err_LoadLayerList, array_pointer, array_size_param = gm.LoadLayerList(
    os.path.join(import_dir, "Depth_to_Water.tif"),
    gm.GM_LoadFlags_HideAllPrompts,
)

# Get handle for the loaded layers
arr_ptr, arr_size = gm.GetLoadedLayerList()
layers_array = gm.GM_LayerHandle_array_frompointer(arr_ptr)
layers_list = gm.carray_to_list(layers_array, arr_size)

### START SAMPLE
# Get the metadata for the loaded image layer
info = gm.GetLayerInfo(layers_list[0])
metadata_array = gm.GM_AttrValue_array_frompointer(info.mMetadataList)
metadata_list = gm.carray_to_list(metadata_array, info.mMetadataListSize)
num_pxs = 0

# Find the attribute that holds the number of pixel columns
for attr in metadata_list:
    if attr.mName == 'NUM COLUMNS':
        num_pxs = int(attr.mVal)

if num_pxs == 0:
    print("Error: Attribute for the number of columns not found.")

# The number of pixel columns denotes the maximum number
# of pixels that can be held in each row of the loaded imagery
elevation_array = gm.float_array(num_pxs)
# Gets the elevations from the specified pixel row
err_GetPixelElevationRow = gm.GetPixelElevationRow(
    layers_list[0], 5858, elevation_array, -9999
)

# Display elevation results
elevation_list = gm.carray_to_list(elevation_array, num_pxs)
elev_range_0_50 = []
elev_range_50_100 = []
elev_range_100_150 = []
elev_range_150_On = []
valid = 0
total = 0
for i in range(18674):
    if elevation_list[i] != -9999.0:
        valid += 1
        total += elevation_list[i]
        if elevation_list[i] <= 50:
            elev_range_0_50.append(elevation_list[i])
        elif elevation_list[i] <= 100:
            elev_range_50_100.append(elevation_list[i])
        elif elevation_list[i] <= 150:
            elev_range_100_150.append(elevation_list[i])
        elif elevation_list[i] > 150:
            elev_range_150_On.append(elevation_list[i])

print("\nTotal valid points in pixel row:", valid)
print("Average elevation:", total / valid)
print("Pixels with elevation 0-50m:", len(elev_range_0_50))
print("Pixels with elevation 50-100m:", len(elev_range_50_100))
print("Pixels with elevation 100-150m:", len(elev_range_100_150))
print("Pixels with elevation >150m:", len(elev_range_150_On), "\n")

### END SAMPLE

assert (
    gm.GM_Error_None == err_LoadLayerList
), f"LoadLayerList: Expected {gm.GM_Error_None}, got {err_LoadLayerList}"
assert (
    gm.GM_Error_None == err_GetPixelElevationRow
), f"GetPixelElevationRow: Expected {gm.GM_Error_None}, got {err_GetPixelElevationRow}"

assert os.path.isdir(import_dir), f"Import directory doesn't exist: {import_dir}"

assert os.path.exists(
    os.path.join(import_dir, "Depth_to_Water.tif")
), "Sample data doesn't exist"

# These are the values they had when the sample was created
assert valid == 8960, "Incorrect data reading"
assert (total / valid) == 60.99821428571428, "Incorrect data reading"
assert len(elev_range_0_50) == 5309, "Incorrect data reading"
assert len(elev_range_50_100) == 1670, "Incorrect data reading"
assert len(elev_range_100_150) == 1287, "Incorrect data reading"
assert len(elev_range_150_On) == 694, "Incorrect data reading"

import os
import globalmapper as gm

import_dir = "Z:/GM-20069/single_function_sample_data/"

# Elevations are extracted from the pixels corresponding to the points in the point list, across all layers, starting 
# from the topmost layer. The collected elevations are then displayed.

# Close all layers that may be open
layers_are_open = gm.GetLoadedLayerList()
if layers_are_open:
    arr_ptr, arr_size = layers_are_open
    initial_layers = gm.GM_LayerHandle_array_frompointer(arr_ptr)
    for i in range(arr_size):
        gm.CloseLayer(initial_layers[i])

err_LoadLayerList1, array_pointer, array_size_param = gm.LoadLayerList(
    os.path.join(import_dir, "depth_img.jpg"),
    gm.GM_LoadFlags_HideAllPrompts,
)
err_LoadLayerList2, array_pointer, array_size_param = gm.LoadLayerList(
    os.path.join(import_dir, "Depth_to_Water.tif"),
    gm.GM_LoadFlags_HideAllPrompts,
)

# Get handle for the loaded layers
arr_ptr, arr_size = gm.GetLoadedLayerList()
layers_array = gm.GM_LayerHandle_array_frompointer(arr_ptr)
layers_list = gm.carray_to_list(layers_array, arr_size)

### START SAMPLE
point1 = gm.GM_Point_t(1553978, 326704)
point2 = gm.GM_Point_t(1093280, 256782)
point_list = [point1, point2]
size = len(point_list) * len(layers_list)
elevation_array = gm.float_array(size)
# Gets the elevation at each location from a layer/layer list (starting with the top-most layer)
err_GetLocationElevationList = gm.GetLocationElevationList(
    layers_list, point_list, elevation_array, 0
)

# Display elevations at the specified points (2 layers were loaded for this sample)
elevation_list = gm.carray_to_list(elevation_array, size)
print("\nLayer 1 data:", elevation_list[0], elevation_list[1])
print("Layer 0 data:", elevation_list[2], elevation_list[3], "\n")

### END SAMPLE

assert (
    gm.GM_Error_None == err_LoadLayerList1
), f"LoadLayerList1: Expected {gm.GM_Error_None}, got {err_LoadLayerList1}"
assert (
    gm.GM_Error_None == err_LoadLayerList2
), f"LoadLayerList2: Expected {gm.GM_Error_None}, got {err_LoadLayerList2}"
assert (
    gm.GM_Error_None == err_GetLocationElevationList
), f"GetLocationElevationList: Expected {gm.GM_Error_None}, got {err_GetLocationElevationList}"