PyQGIS: KOMPSAT-3 위성영상 시각화 코드 예제

안녕하세요? 이번 글은 PyQGIS를 통해 KOMPSAT-3 위성영상을 시각화하는 코드 예제를 제공해 봅니다. QGIS 메뉴 바에서 "플러그인 > 파이썬 콘솔"에서 실행해보시면 되겠습니다.

 

먼저, 필요한 라이브러리를 호출하겠습니다.

import os
from qgis.core import (
    QgsRasterLayer,
    QgsProject,
    QgsContrastEnhancement,
    QgsRasterMinMaxOrigin,
    QgsRectangle,
    QgsRasterTransparency,
)
from qgis.utils import iface
from osgeo import gdal

 

KOMPSAT-3 위성영상 폴더를 정의하고, 파일 명명 규칙에 따라 해당 폴더에 위치한 4장의 다중분광 밴드(Multispectral)의 경로를 설정합니다.

출처: https://en.wikipedia.org/wiki/KOMPSAT-3

# 경로 설정
folder_name = "K3_20240114043957_62211_09341264_L1G"
raster_dir = "D:GEODATA/{folder_name}")

b_band_path = os.path.join(raster_dir, f"{folder_name}_B.tif")
g_band_path = os.path.join(raster_dir, f"{folder_name}_G.tif")
r_band_path = os.path.join(raster_dir, f"{folder_name}_R.tif")
n_band_path = os.path.join(raster_dir, f"{folder_name}_N.tif")

 

다중분광 밴드들을 하나의 VRT 파일로 생성합니다.

# VRT 파일 생성
vrt_path = os.path.join(raster_dir, "k3_composite.vrt")
vrt_options = gdal.BuildVRTOptions(separate=True)
vrt_result = gdal.BuildVRT(
    vrt_path,
    [b_band_path, g_band_path, r_band_path, n_band_path],
    options=vrt_options
)
print("VRT file created at:", vrt_path)

 

VRT 파일을 QgsRasterLayer로 불러와서 QGIS 프로젝트에 추가해 보겠습니다.

# VRT 파일을 QgsRasterLayer로 불러오기
vrt_layer = QgsRasterLayer(vrt_path, "K3 Composite", "gdal")

# 레이어 유효성 검사 및 QGIS 프로젝트에 추가
if not vrt_layer.isValid():
    print(f"Failed to load layer from {vrt_path}")
else:
    QgsProject.instance().addMapLayer(vrt_layer)
    print(f"K3 Composite Layer added and styled successfully!")

 

레이어 추가 결과는 다음과 같습니다.

 

현재 레이어는 RGB Composite가 1:2:3 밴드로 설정되어 있습니다.

 

이는 "설정 > 옵션 > 래스터 탭 > 밴드 및 리샘플링"에서 RGB 밴드 선택이 적색 밴드 1, 녹색 밴드 2, 청색 밴드 3으로 기본 설정되어 있기 때문입니다. RGB 밴드 선택을 3:2:1 밴드로 설정해서 Natural Color(자연색)를 구현해 보겠습니다.

 

다음과 같이 코드를 추가했습니다.

# 레이어 유효성 검사 및 QGIS 프로젝트에 추가
if not vrt_layer.isValid():
    print(f"Failed to load layer from {vrt_path}")
else:
    # RGB 밴드 선택
    renderer = vrt_layer.renderer()
    renderer.setRedBand(3)  # Red band를 Band 3으로 설정
    renderer.setGreenBand(2)  # Green band를 Band 2로 설정
    renderer.setBlueBand(1)  # Blue band를 Band 1로 설정

    QgsProject.instance().addMapLayer(vrt_layer)
    print(f"K3 Composite Layer added and styled successfully!")

 

결과는 아래와 같습니다. RGB 밴드 선택은 의도대로 변경되었는데, 최소 / 최대 값 설정을 사용자 정의에서 다른 옵션으로 바꿔주는 것이 좋겠습니다. 여기서는 "누적 카운트 삭감" 옵션을 적용해 보겠습니다.

 

아래와 같이 대비 향상 설정 적용을 위한 코드를 추가합니다.

# 레이어 유효성 검사 및 QGIS 프로젝트에 추가
if not vrt_layer.isValid():
    print(f"Failed to load layer from {vrt_path}")
else:

    # RGB 밴드 선택
    renderer = vrt_layer.renderer()
    renderer.setRedBand(3)  # Red band를 Band 3으로 설정
    renderer.setGreenBand(2)  # Green band를 Band 2로 설정
    renderer.setBlueBand(1)  # Blue band를 Band 1로 설정

    # 대비 향상 알고리즘: 최소/최대값까지 균일화
    contrast_algorithm = QgsContrastEnhancement.StretchToMinimumMaximum

    # 최소/최대값 설정: 누적 카운트 삭감
    limits = QgsRasterMinMaxOrigin.CumulativeCut

    # 대비 향상 설정 적용
    vrt_layer.setContrastEnhancement(contrast_algorithm, limits, QgsRectangle(), 0, True)

    QgsProject.instance().addMapLayer(vrt_layer)
    print(f"K3 Composite Layer added and styled successfully!")

 

이제 누적 카운트 삭감 옵션까지 적용되었습니다. 이번에는 RGB 밴드 값이 0인 경우에 투명도를 적용해 보겠습니다.

 

아래와 같이 투명 픽셀 설정 코드를 추가해 봅니다. 이렇게 하면 "투명도 > 사용자 정의 투명 옵션 > 투명 픽셀 목록"에 적색, 녹색, 청색이 모두 0인 경우 투명도가 100으로 적용되는 조건이 추가됩니다.

    # 투명 픽셀 설정
    transparency  = vrt_layer.renderer().rasterTransparency()
    pixel = QgsRasterTransparency.TransparentThreeValuePixel()
    pixel.blue = 0
    pixel.green = 0
    pixel.red = 0
    pixel.percentTransparent = 100
    transparency.setTransparentThreeValuePixelList([pixel])
    vrt_layer.triggerRepaint()

 

결과는 아래와 같습니다.

 

전체 코드 중 최소/최대값 설정은 평균 +/- 표준편차로, RGB 밴드 선택을 False Color(위색)으로 설정하고자 한다면, 아래와 같이 일부 코드를 변경해 주시면 됩니다.

import os
from qgis.core import (
    QgsRasterLayer,
    QgsProject,
    QgsContrastEnhancement,
    QgsRasterMinMaxOrigin,
    QgsRectangle,
    QgsRasterTransparency,
)
from qgis.utils import iface
from osgeo import gdal

# 경로 설정
folder_name = "K3_20240114043957_62211_09341264_L1G"
raster_dir = "D:GEODATA/{folder_name}")

b_band_path = os.path.join(raster_dir, f"{folder_name}_B.tif")
g_band_path = os.path.join(raster_dir, f"{folder_name}_G.tif")
r_band_path = os.path.join(raster_dir, f"{folder_name}_R.tif")
n_band_path = os.path.join(raster_dir, f"{folder_name}_N.tif")

# VRT 파일 생성
vrt_path = os.path.join(raster_dir, "k3_composite.vrt")
vrt_options = gdal.BuildVRTOptions(separate=True)
vrt_result = gdal.BuildVRT(
    vrt_path,
    [b_band_path, g_band_path, r_band_path, n_band_path],
    options=vrt_options
)
print("VRT file created at:", vrt_path)

# VRT 파일을 QgsRasterLayer로 불러오기
vrt_layer = QgsRasterLayer(vrt_path, "K3 Composite", "gdal")

# 레이어 유효성 검사 및 QGIS 프로젝트에 추가
if not vrt_layer.isValid():
    print(f"Failed to load layer from {vrt_path}")
else:

    # RGB 밴드 선택
    renderer = vrt_layer.renderer()
    renderer.setRedBand(4)  # Red band 설정
    renderer.setGreenBand(3)  # Green band 설정
    renderer.setBlueBand(2)  # Blue band 설정

    # 대비 향상 알고리즘: 최소/최대값까지 균일화
    contrast_algorithm = QgsContrastEnhancement.StretchToMinimumMaximum

    # 최소/최대값 설정: 평균 +/- 표준편차
    limits = QgsRasterMinMaxOrigin.StdDev

    # 대비 향상 설정 적용
    vrt_layer.setContrastEnhancement(contrast_algorithm, limits, QgsRectangle(), 0, True)

    QgsProject.instance().addMapLayer(vrt_layer)
    print(f"K3 Composite Layer added and styled successfully!")

    # 투명 픽셀 설정
    transparency  = vrt_layer.renderer().rasterTransparency()
    pixel = QgsRasterTransparency.TransparentThreeValuePixel()
    pixel.blue = 0
    pixel.green = 0
    pixel.red = 0
    pixel.percentTransparent = 100
    transparency.setTransparentThreeValuePixelList([pixel])
    vrt_layer.triggerRepaint()