Checking a list of rasters to see if they correspond to a list of points.

I have a list of large rasters that I would like to check systematically to see if they cover the same area as a list of points. Below is a program to walk the files and subdirectories of a folder and for each .tif, see if it contains at least 1 point. I used Google Colab to run this.

!pip install rasterio
!pip install geopandas
import rasterio
import rasterio.plot
import geopandas as gpd
from shapely.geometry import Point, Polygon


def plot_tif(tif_path, amuquall_df):
    tif_file =
    # print(tif_path)
    # print(tif_file.shape)
    # print(tif_file.bounds,tif_file.bounds[1])

    # make a polygon of the tif
    t_left = tif_file.bounds[0]
    t_bottom = tif_file.bounds[1]
    t_right = tif_file.bounds[2]
    t_top = tif_file.bounds[3]
    coords = [(t_left, t_bottom), (t_right, t_bottom), (t_right, t_top), (t_left, t_top)]    
    poly = Polygon(coords)
    print("raster crs: ",
    amuquall_df = amuquall_df.to_crs(
    print("points crs: ",
    tif_contains_points = False
    for pt in amuquall_df['geometry']:
            tif_contains_points = True

import os
usable_tifs = []
unusable_tifs = []
for dirpath, dirs, files in os.walk("/content/drive/MyDrive/crane_ane/Satellite_Data"):
    for filename in files:
        fname = os.path.join(dirpath,filename)
        if fname.endswith('.tif'):
            has_points = plot_tif(fname, amuquall_df)
            if(has_points == True):
            if(has_points == False):
print("*** usable tifs ***")
# print("*** UNusable tifs ***")
# print("\n".join(unusable_tifs))

Here is the output:

*** usable tifs *** /content/drive/MyDrive/crane_ane/Satellite_Data/SPOT/SPOT Images/Amuq.tif /content/drive/MyDrive/crane_ane/Satellite_Data/SPOT/SPOT Images/mosaic.tif /content/drive/MyDrive/crane_ane/Satellite_Data/SPOT/SPOT Images/nimadata1056733980.tif

… etc.

Leave a Reply

Your email address will not be published. Required fields are marked *