python - Retrieve pixel value with geographic coordinate as input with gdal

I have a list of X,Y coordinates in UTM called coords. I also have a .tif of a digital terrain model (DTM) also referenced in UTM. I would like to use the Python wrapper for gdal to return the pixel values (i.e. elevation values) for each coordinate in coords.

In my search for an answer so far I have only found this answer, which is related, but does not just use a simple list of coordinates and is seemingly a bit more complex than the simple task I am trying to achieve.

My concern is not so much the particular script that needs to be written, but the methods or functions within gdal that can perform this type of function. Is there a simple function within gdal that takes a geographic coordinate and returns a pixel value?

Answer

Here is the function I came up with, using a function I found in another stack post (that I unfortunately cannot remember the title of). It was originally written to be used with a point vector file instead of manually inputting the points like I am doing. Below is the simplified function, using affine and gdal, where data_source is an opened gdal object of a GeoTIFF and coord is a tuple of a geo-coordinate. This tuple must be in the same coordinate system as the GeoTIFF.

def retrieve_pixel_value(geo_coord, data_source):
    """Return floating-point value that corresponds to given point."""

    x, y = geo_coord[0], geo_coord[1]
    forward_transform =  \
        affine.Affine.from_gdal(*data_source.GetGeoTransform())
    reverse_transform = ~forward_transform
    px, py = reverse_transform * (x, y)
    px, py = int(px + 0.5), int(py + 0.5)
    pixel_coord = px, py

    data_array = np.array(data_source.GetRasterBand(1).ReadAsArray())
    return data_array[pixel_coord[0]][pixel_coord[1]]