I have rasters of Sentinel-2 data and want to split them im subrasters (1km*1km). This is working so far but takes a lot of time due to the for loops in Python.
ds = gdal.Open(d)
llx, xres, xskew, uly, yskew, yres = ds.GetGeoTransform()
lly = uly + (ds.RasterYSize * yres)
lrx = llx + (ds.RasterXSize * xres)
lly_s = int(str(lly)[0:4]+'000')
llx_s = int(str(llx)[0:3] + '000')
for i, j in itertools.product(range(0, math.ceil(lrx-llx_s), 1000), range(0, math.ceil(uly-lly_s), 1000)):
ulx_k = llx_s + i #upper left x == lower left x
lrx_k = ulx_k + 1000
lry_k = lly_s + j
uly_k = lry_k + 1000
translateOpt = gdal.TranslateOptions(format="Gtiff",projWin=(ulx_k, uly_k, lrx_k, lry_k),projWinSRS='epsg:25832', creationOptions = ["TILED=YES","COMPRESS=LZW"])
ds = gdal.Translate(destName = f'{d[:-4]}_1000/{ulx_k}_{lry_k}.tif', srcDS=d, options = translateOpt)
ds.FlushCache()
del ds
This is my code. It is working but really slow (kind of expected). Just wanted to know if someone knows a way to improve the performance.