Compare planet downloaded images with different resolution in Python

Question

I have downloaded satellite images using Planet API for an Area of Interest. (Image Details: PSScene4Band, analytic)The problem is Planet returns images with different number of pixels for the same subarea, with varying dates. I want to compare NDVI pixel by pixel of these images. Since the number of pixels differs, it's not possible.

Attached images I downloaded.Images link here

File 20180502_181704_1025_subarea.tif has size: 2027, 1193.

File 20170919_180930_0f22_subarea.tif has size: 2028, 1194

My code

import numpy, sys
from osgeo import gdal
from osgeo.gdalconst import * 
import numpy as np
import random
import matplotlib.pyplot as plt


# register all of the GDAL drivers 
gdal.AllRegister()

# open the image 1
inDs = gdal.Open("20170919_180930_0f22_subarea.tif")
if inDs is None:
   print('Could not open image file')
   sys.exit(1)
r=np.array(inDs.GetRasterBand(3).ReadAsArray(), dtype=float)
n=np.array(inDs.GetRasterBand(4).ReadAsArray(), dtype=float)

np.seterr(divide='ignore', invalid='ignore') 
#Ignore the divided by zero or 
Nan appears
ndvi1=(n-r)/(n+r) # The NDVI formula


# open the image 2
inDs2 = gdal.Open("sample/20180502_181704_1025_subarea.tif")
if inDs2 is None:
   print('Could not open image file')
   sys.exit(1)

r=np.array(inDs2.GetRasterBand(3).ReadAsArray(), dtype=float)
n=np.array(inDs2.GetRasterBand(4).ReadAsArray(), dtype=float)

np.seterr(divide='ignore', invalid='ignore') #Ignore the divided by zero or 
Nan appears
ndvi2=(n-r)/(n+r) # The NDVI formula

ndvi = ndvi2 - ndvi1

Problem occurs in the line ndvi = ndvi2 - ndvi1 since the numpy array dimensions differ, as the image resolution (rows x columns) differs. Note: I assume ignoring the pixels may result in the wrong result; as I have no clue which pixel is missing or extra.

I did make sure that both the TIFF images have same corner coordinates, but pixels differ.

Answer 1

You will have to resample the data. I recommend resampling with GDAL, see the following question for more details on resampling in Python using GDAL: Resampling a raster from python without using gdalwarp

As others have mentioned, you'll need to consider which resampling algorithm to use. Nearest-neighbor is one way to avoid changing the data values, and is commonly used for cases like this. Other methods can give nicer images, but may distort the data slightly.

Answer 2

If rescaling the images is not an option, you will have to forgo the the idea of comparing the images pixel by pixel.

You will need to do something like a histogram comparison, or correlation.

See if one of the methods described here suits you. Or here. Or here