Select to view content in your preferred language

Loop through raster names in two folders, perform Con (Spatial Analyst tool) if the dates match

733
5
Jump to solution
09-06-2022 07:41 AM
VasudhaChaturvedi
New Contributor II

I cannot perform the Con tool over all the raster images with same date in two folders. What I am trying to achieve is producing an image where the pixel values of 1 are converted to 0 and the pixel values of 0 are replaced with the pixel values of another image using the images with the same dates. 

 

This is the script I have worked upon so far. 

 

import arcpy, os
from arcpy.sa import *

# Set the current workspace
location_a = "D:/datasets/regression_interpolation_rasters/airtemp_snowmask_regression_rasters"
location_b = "D:/datasets/regression_interpolation_rasters/regression_rasters"
location = "D:/datasets/regression_interpolation_rasters/snowmask_rasters"

# get the names of rasters
arcpy.env.workspace = location_a
snow_raster = arcpy.ListRasters("", "TIF")
for snow_ras in snow_raster:
    snow = Raster(snow_ras)
    #print(type(snow))
    snow_name = snow_ras.split("_")
    snow_values = snow_name[2]
    #print(snow_values) 
    
# get the names of rasters 
arcpy.env.workspace = location_b
airtemp_raster = arcpy.ListRasters("", "TIF")
#print(airtemp_raster)

for airtemp_ras in airtemp_raster:
    airtemp = Raster(airtemp_ras)
    airtemp_name = airtemp_ras.split("_")
    airtemp_values = airtemp_name[0]
    #print(airtemp_values)
    if airtemp_values == snow_values:
        output_raster = Con(snow_raster==1, 0, airtemp_raster)  
        output = os.path.join(location,airtemp_ras.split("_")[0] +'_airtemp'+ '.tif')
        output_raster.save(output)
            #print("yes")
    else:
        print("action cannot take place")

 

 

All it is showing the else statement even though the dates match.The placement of the dates vary but the format of dates are similar.These are the images in location_b

VasudhaChaturvedi_1-1662474912616.png

 

These are the images in location_a

VasudhaChaturvedi_2-1662474932205.png

 

0 Kudos
1 Solution

Accepted Solutions
JohannesLindner
MVP Frequent Contributor
  • You don't save your snow_values anywhere, so snow_values is just the date part of the last raster in location_a, which is why the comparison in line 28 returns False.
  • You have location_a and location_b the wrong way around (or you mislabeled your screenshots).
  • Your variable names are really misleading.

Here's what you want to do (untested):

# build a dictionary {date: raster} for the snow rasters
# I'm going by your screenshots: snow rasters in location_b
snow_raster_dict = dict()
arcpy.env.workspace = location_b
snow_raster_list = arcpy.ListRasters("", "TIF")
for snow_name in snow_raster_list:
    snow_raster = Raster(snow_name)
    #print(type(snow_raster))
    snow_name_parts = snow_name.split("_")
    snow_date = snow_name_parts[2]
    snow_raster_dict[snow_date] = snow_raster
#print(snow_raster_dict)    

# get the names of airtemp rasters 
# again, going by your screenshots -> location_a
arcpy.env.workspace = location_a
airtemp_raster_list = arcpy.ListRasters("", "TIF")
#print(airtemp_raster_list)

for airtemp_name in airtemp_raster_list:
    airtemp_raster = Raster(airtemp_name)
    airtemp_name_parts = airtemp_name.split("_")
    airtemp_date = airtemp_name_parts[0]
    # try to get the corresponding snow raster
    try:
        snow_raster = snow_raster_dict[airtemp_date]
        # do the calculation
        output_raster = Con(snow_raster==1, 0, airtemp_raster)  
        output = os.path.join(location,airtemp_name.split("_")[0] +'_airtemp'+ '.tif')
        output_raster.save(output)
        #print(f"created new airtemp raster for {airtemp_date}")
    # catch KeyErrors (date not found in snow_raster_dict)
    except KeyError:
        print(f"no snow raster found for {airtemp_date}")

 


Have a great day!
Johannes

View solution in original post

5 Replies
by Anonymous User
Not applicable

Your 'snow_val' in your conditional is the last item assigned to it in the first loop, so its not iterating as intended.  You need to update your loops to nest the

for airtemp_ras in airtemp_raster:

logic so that checks the current snow_val to each air_temp item.

 

 

0 Kudos
JohannesLindner
MVP Frequent Contributor
  • You don't save your snow_values anywhere, so snow_values is just the date part of the last raster in location_a, which is why the comparison in line 28 returns False.
  • You have location_a and location_b the wrong way around (or you mislabeled your screenshots).
  • Your variable names are really misleading.

Here's what you want to do (untested):

# build a dictionary {date: raster} for the snow rasters
# I'm going by your screenshots: snow rasters in location_b
snow_raster_dict = dict()
arcpy.env.workspace = location_b
snow_raster_list = arcpy.ListRasters("", "TIF")
for snow_name in snow_raster_list:
    snow_raster = Raster(snow_name)
    #print(type(snow_raster))
    snow_name_parts = snow_name.split("_")
    snow_date = snow_name_parts[2]
    snow_raster_dict[snow_date] = snow_raster
#print(snow_raster_dict)    

# get the names of airtemp rasters 
# again, going by your screenshots -> location_a
arcpy.env.workspace = location_a
airtemp_raster_list = arcpy.ListRasters("", "TIF")
#print(airtemp_raster_list)

for airtemp_name in airtemp_raster_list:
    airtemp_raster = Raster(airtemp_name)
    airtemp_name_parts = airtemp_name.split("_")
    airtemp_date = airtemp_name_parts[0]
    # try to get the corresponding snow raster
    try:
        snow_raster = snow_raster_dict[airtemp_date]
        # do the calculation
        output_raster = Con(snow_raster==1, 0, airtemp_raster)  
        output = os.path.join(location,airtemp_name.split("_")[0] +'_airtemp'+ '.tif')
        output_raster.save(output)
        #print(f"created new airtemp raster for {airtemp_date}")
    # catch KeyErrors (date not found in snow_raster_dict)
    except KeyError:
        print(f"no snow raster found for {airtemp_date}")

 


Have a great day!
Johannes
VasudhaChaturvedi
New Contributor II

I corrected the labeling of the raster images in the screenshot. The screenshot labelling was wrong, which has been corrected now

0 Kudos
VasudhaChaturvedi
New Contributor II

Thanks @JohannesLindner for your suggestion. I updated the code as per your suggestion, the images whose dates match are being stored but somehow all the images have the same pixel value of 1 replaced by 0 which is incorrect. 

For example this is one of the output stored where the pixel value of 1 is replaced with pixel value 0, and pixel value of 0 is replaced with the pixel value of another raster.

VasudhaChaturvedi_0-1662534896210.png

But the binary image with 1(grey color) and 0 (black color) pixel values looks like this:

VasudhaChaturvedi_1-1662534993193.png

And the image with which the pixel values of 0 in binary image will be replaced is:

VasudhaChaturvedi_2-1662535097814.png 

 

 

 

0 Kudos
VasudhaChaturvedi
New Contributor II

Your code actually worked @JohannesLindner. Thank you for your help. It was just a minor mistake that was causing a problem. Instead of single = I used double == signs for the try statement. So the final code to solve my query is here.

import arcpy, os
from arcpy.sa import *
from arcpy import env

# Set the current workspace
location_a = "D:/datasets/planetscope/binary_snow_images_converted"
location_b = "D:/datasets/regression_interpolation_rasters/linear_regression_interpolation_rasters_without_snowmask"
location = "D:/datasets/regression_interpolation_rasters/regression_rasters_after_snowmask"
arcpy.env.extent = "D:/datasets/DEM ALTA_VALTELLINA LOMBARDY REGION/clip_dtm.tif"
arcpy.env.snapRaster = "D:/datasets/DEM ALTA_VALTELLINA LOMBARDY REGION/clip_dtm.tif"
arcpy.env.cellSize = "D:/datasets/regression_interpolation_rasters/linear_regression_interpolation_rasters_without_snowmask/2017-10-01_airtemp.tif"

# get the names of rasters
snow_raster_dict = dict()
all_snow_rasters = []
snow_raster_list = arcpy.ListRasters("", "TIF")

for snow_name in snow_raster_list:
    snow_raster = Raster(snow_name)
    snow_name_parts = snow_name.split("_")
    snow_date = snow_name_parts[2]
    snow_raster_dict[snow_date] = snow_raster
    
# get the names of rasters 
arcpy.env.workspace = location_b
airtemp_raster_list = arcpy.ListRasters("", "TIF")

for airtemp_name in airtemp_raster_list:
    airtemp_raster = Raster(airtemp_name)
    airtemp_name_parts = airtemp_name.split("_")
    airtemp_date = airtemp_name_parts[0]

    try:
        snow_raster = snow_raster_dict[airtemp_date]
        output_raster = Con(snow_raster==1, 0, airtemp_raster)  
        output = os.path.join(location,airtemp_name.split("_")[0] +'_airtemp'+ '.tif')
        output_raster.save(output)
    except KeyError:
        print(f"no snow raster found for {airtemp_date}")
0 Kudos