Page 1 of 1

Verkettung von GDAL-Funktionen

Posted: 23 Dec 2024, 12:54
by Anonymous
Ich versuche, drei Einzelbandbilder zu lesen und sie mit gdal/python zu einem 3-Bandbild zu kombinieren.
Ich habe eine funktionierende Lösung gefunden, aber sie scheint nicht optimal zu sein:
from osgeo import gdal
import numpy as np

red_path = "F:/B4.TIF"
green_path = "F:/B2.TIF"
blue_path = "F:/B1.TIF"

red = gdal.Open(red_path)
gt = red.GetGeoTransform()
prj = red.GetProjection()

red_img = gdal.Open(red_path)
red_band = red_img.GetRasterBand(1)
red_array = red_band.ReadAsArray().astype(np.float32)

green_img = gdal.Open(green_path)
green_band = green_img.GetRasterBand(1)
green_array = green_band.ReadAsArray().astype(np.float32)

blue_img = gdal.Open(blue_path)
blue_band = blue_img.GetRasterBand(1)
blue_array = blue_band.ReadAsArray().astype(np.float32)

driver = gdal.GetDriverByName("GTiff")
driver.Register()

final_img = driver.Create("combined_colors.tiff",
xsize=red_array.shape[1],
ysize=red_array.shape[0],
bands=3,
eType=gdal.GDT_Float32)
final_img.SetGeoTransform(gt)
final_img.SetProjection(prj)

first_band = final_img.GetRasterBand(1)
first_band.WriteArray(red_array)
first_band.SetNoDataValue(np.nan)
first_band.FlushCache()
first_band = None

second_band = final_img.GetRasterBand(2)
second_band.WriteArray(green_array)
second_band.SetNoDataValue(np.nan)
second_band.FlushCache()
second_band = None

third_band = final_img.GetRasterBand(3)
third_band.WriteArray(blue_array)
third_band.SetNoDataValue(np.nan)
third_band.FlushCache()
third_band = None

final_img = None

Jetzt möchte ich die gleichen Schritte ausführen, aber in einer saubereren Syntax. Deshalb habe ich versucht, mehrere Vorgänge in einer Zeile auszuführen, aber es funktioniert nicht mehr.
from osgeo import gdal
import numpy as np

red_path = "F:/B4.TIF"
green_path = "F:/B2.TIF"
blue_path = "F:/B1.TIF"

red = gdal.Open(red_path)
gt = red.GetGeoTransform()
prj = red.GetProjection()

red = gdal.Open(red_path).GetRasterBand(1).ReadAsArray().astype(np.float32)
green = gdal.Open(green_path).GetRasterBand(1).ReadAsArray().astype(np.float32)
blue = gdal.Open(blue_path).GetRasterBand(1).ReadAsArray().astype(np.float32)

driver = gdal.GetDriverByName("GTiff")
driver.Register()

final_img = driver.Create("combined_colors.tiff",
xsize=red.shape[1],
ysize=red.shape[0],
bands=3,
eType=gdal.GDT_Float32)
final_img.SetGeoTransform(gt)
final_img.SetProjection(prj)

final_img.GetRasterBand(1).WriteArray(red)
final_img.GetRasterBand(2).WriteArray(green)
final_img.GetRasterBand(3).WriteArray(blue)
final_img.FlushCache()

final_img = None

Weiß jemand, warum der zweite Code nicht funktioniert?
Vielen Dank im Voraus