코드 2020의 출현: 12일째 QGIS와 Python 사용

46670 단어 gispythonadventofcode
오늘의 도전은 나룻배 운전을 포함한다.좌표를 처리하기 때문에 QGIS에서 재미있는 가시화를 할 수 있습니다
이 글에서 언급한 것은 지리정보시스템(GIS), 변환 매트릭스, 모방 변환

1부 도전


Link to challenge on Advent of Code 2020 website
도전은 또 다른 도전이다. 명령을 읽고 그것을 실행하라.이런 상황에서 당신은 나룻배로 지시에 따라 항행한다.모든 명령에는 명령과 숫자가 포함되어 있다.
그 중 네 개의 명령은 나룻배를 나침반 방향으로 이동시킨다: N, E, S, W;그 중 두 명령은 나룻배의 항로를 바꾸었다: L, R: 나룻배는 나침반이 아닌 항로를 따라 이동하라는 명령도 있었다. F.
이것은 우리의 나룻배가 내부 상태가 있다는 것을 의미한다. 그것은 우주에서의 위치와 항행 방향이다.

QGIS


이외에 나는 이 시리즈의 블로그에서 파이톤의 다양한 응용 프로그램을 토론하고 싶지 않다. 나는 QGIS, 원본을 개설한 GIS (지리정보시스템) 프로그램을 사용해서 파이톤으로 스크립트를 작성할 수 있다.
GIS는 완전한 업계, 기능 집합, 도구 집합이며,python은 나쁜 선택이 아니라 이미지와 데이터 처리 라이브러리 때문이다.GIS 작업에 사용되는python 상용 라이브러리는 GDAL 귀속과shapely를 포함합니다.그러나 오늘 우리는 이것들을 사용하지 않고 내장된 파이톤QGIS을 사용합니다. 이것은 자신의SDK/API를 가지고 있습니다.

나는 위성과 무인기 이미지를 처리할 때 QGIS를 자주 사용하지만, 주로 스크립트를 작성하는 것이 아니라 데이터를 불러오고 검사하는 데 쓰인다.따라서 내가 자주 사용하지만 아직 깊이 연구하지 않은 프로그램의python 스크립트 기능을 탐색하는 것은 매우 재미있다.

비밀 번호


실제 Python 코드는 매우 간단합니다. 데이터를 가져오고, 그 위에서 순환하며, 상태를 업데이트하고, 나룻배를 이동해서, 그 행방을 추적합니다.
이를 위해, 명령과 선박의 항로 (선박의 유일한 고유 상태) 를 읽고, 운동 벡터와 새로운 항로를 출력하며, 외부 순환은 선박의 외적 상태를 추적하는 데 사용될 함수를 작성할 것이다.
상태를 줄이려면 이 함수를 통과해야 하는 수량을 제외하고는 상태를 내재와 외재로 나눌 충분한 이유가 없다.
from math import cos, sin, radians

vectors = {
    "N": (0, 1, 0, 0),
    "E": (1, 0, 0, 0),
    "S": (0, -1, 0, 0),
    "W": (-1, 0, 0, 0),
    "L": (0, 0, -1, 0),
    "R": (0, 0, 1, 0),
    "F": (0, 0, 0, 1)
}

def process_command(command, heading):
    code = command[0]
    amount = int(command[1:])

    easting, northing, rotate, forward = vectors[code]
    north_delta = amount*(northing + forward*cos(radians(heading)))
    east_delta = amount*(easting + forward*sin(radians(heading)))
    heading += rotate * amount

    return east_delta, north_delta, heading 
이곳의 작은 기교는 대량의if/elif 순환을 피하고 각 사례를 각각 처리하는 것이다.반대로 우리는 코드의 지점을 줄이고 규칙을 일종의 전환으로 응용하는 유사한 선형 대수 방법을 사용했을 뿐이다. (실제로 나는 어제의 해결 방안을 나의 생각을 이 방향으로 밀어붙인 탓으로 돌렸다.)
이것이야말로 코드의 핵심이다.하지만 더 흥미로운 것은 QGIS를 시각화한다는 점이다.

형상화


이를 위해 QGIS 기능에 이러한 데이터 포인트를 삽입하는 등 전체 코드는 다음과 같습니다.
from math import cos, sin, radians
from PyQt5.QtWidgets import QFileDialog
from qgis.PyQt.QtCore import QVariant
from qgis.core import (
    Qgis,
    QgsField,
    QgsProject,
    QgsVectorLayer,
    QgsFeature,
    QgsGeometry,
    QgsPoint,
    QgsPointXY,
    QgsCoordinateReferenceSystem,
    QgsMessageLog
)

vectors = {
    "N": (0, 1, 0, 0),
    "E": (1, 0, 0, 0),
    "S": (0, -1, 0, 0),
    "W": (-1, 0, 0, 0),
    "L": (0, 0, -1, 0),
    "R": (0, 0, 1, 0),
    "F": (0, 0, 0, 1)
}

def process_command(command, heading):
    code = command[0]
    amount = int(command[1:])

    easting, northing, rotate, forward = vectors[code]
    north_delta = amount*(northing + forward*cos(radians(heading)))
    east_delta = amount*(easting + forward*sin(radians(heading)))
    heading += rotate * amount

    return east_delta, north_delta, heading


def main():
    data_file = QFileDialog.getOpenFileName(None, "Open data file")
    if not data_file:
        return

    QgsMessageLog.logMessage(f"Opening file {data_file}", "AoC", Qgis.Info)

    # create vector layers
    vl_points = QgsVectorLayer("Point", "AoC Points",  "memory")
    pr_points = vl_points.dataProvider()
    pr_points.addAttributes([
        QgsField("idx", QVariant.Int),
        QgsField("cmd", QVariant.String),
        QgsField("heading", QVariant.Double),
    ])
    vl_points.updateFields()

    vl_lines = QgsVectorLayer("LineString", "AoC Lines",  "memory")
    pr_lines = vl_lines.dataProvider()

    # read data
    data = open(data_file[0]).readlines()
    QgsMessageLog.logMessage(f"{data}", "AoC", Qgis.Info)

    heading = 90 # ship starts facign east
    northing = 0
    easting = 0
    line = [QgsPoint(0, 0)]
    for idx, entry in enumerate(data):
        east_delta, north_delta, heading = process_command(entry, heading)
        northing += north_delta
        easting += east_delta
        QgsMessageLog.logMessage(f"cmd[{idx}]:{entry} E{easting:.0f}, N{northing:.0f}, H{heading}", "AoC", Qgis.Info)

        # add feature
        line.append(QgsPoint(easting, northing))

        feat = QgsFeature()
        feat.setGeometry(QgsGeometry.fromPointXY(QgsPointXY(easting, northing)))
        feat.setAttributes([idx, entry, heading])
        pr_points.addFeatures([feat])
    vl_points.updateExtents()

    # construct linestring
    feat = QgsFeature()
    feat.setGeometry(QgsGeometry.fromPolyline(line))
    pr_lines.addFeatures([feat])
    vl_lines.updateExtents()

    QgsProject.instance().addMapLayer(vl_points)
    QgsProject.instance().addMapLayer(vl_lines)

    # disable crs (we need flat world)
    crs = QgsCoordinateReferenceSystem(None)
    QgsProject.instance().setCrs(crs)

main()
이것은 우리 프로젝트에 두 개의 QGIS 벡터층을 만들고 CRS를 평탄한 지구를 원하기 때문에 '없음' 으로 설정할 것이다.그리고 우리는 기본 설정을 받아들이지 않고 기호를 약간 조정하여 코드에 설정된 heading 속성을 기호로 회전시키는 것을 포함하여 예쁜 색깔과 모양을 제공할 수 있다.

우리는 이 귀여운 나룻배 운동도를 만들었다.

화살표가 나룻배의 항행 방향을 표시하기;텍스트 레이블은 Viridis 색상 체계를 사용하여 짙은 남색에서 노란색으로 변경하는 명령입니다(Viridis는 다른 각도보다 빨간색과 녹색 블라인드에 적합합니다).
코드에서 시각화를 사용자 정의할 수 있지만 이 연습의 목적에 따라 UI를 사용하여 나머지 시각화를 구성하는 것이 더 쉽습니다.

도전의 두 번째


도전의 두 번째 부분에서 일이 좀 커졌다.지금은 선박을 직접 이동하는 것이 아니라 가상 항로점을 이동한 후 선박이 이 항로점으로 이동한다.이곳의 항로점은 사실상 운동 벡터를 정의하는 또 다른 내재 상태이다.
현재, 네 개의 명령이 나침반 방향에서 항로점을 이동한다. N, E, S, W;그 중 두 개의 명령은 그 원점을 둘러싸고 항로점을 회전시킨다. L, R: 그리고 한 명령은 나룻배가 항로점에 정의된 벡터에 따라 이동하는 것을 알려준다.
따라서 나는 코드를 조정하여 벡터에서'앞으로'분량을 삭제하는데 다음과 같은 주요 변화가 있다.

vectors = {
    "N": (0, 1, 0),
    "E": (1, 0, 0),
    "S": (0, -1, 0),
    "W": (-1, 0, 0),
    "L": (0, 0, -1),
    "R": (0, 0, 1),
}

def process_command(command, ship, waypoint):
    code = command[0]
    amount = int(command[1:])

    if code == "F":
        ship[0] += amount*waypoint[0]
        ship[1] += amount*waypoint[1]
    else:
        easting, northing, rotate = vectors[code]
        angle = radians(-amount*rotate)
        old_e, old_n = waypoint
        waypoint[0] = amount*easting + old_e*cos(angle) - old_n*sin(angle)
        waypoint[1] = amount*northing + old_e*sin(angle) + old_n*cos(angle)
여기서 나는 함수를 갱신하여 전체 선박과 항로점 상태를 받아들이고'F'상태를 분할하여 항로점 벡터를 사용하여 선박의 위치를 갱신했다.나머지 코드는 현재 항로점 상태를 업데이트하는 데만 사용됩니다.코드 자체는 회전과 평이를 실현하기 위해 모방 변환 행렬을 효과적으로 실현했다.
코드의 나머지 부분은 시각화를 위해 점 저장 방식을 업데이트합니다.
from math import cos, sin, radians, degrees, atan2
from PyQt5.QtWidgets import QFileDialog
from qgis.PyQt.QtCore import QVariant
from qgis.core import (
    Qgis,
    QgsField,
    QgsProject,
    QgsVectorLayer,
    QgsFeature,
    QgsGeometry,
    QgsPoint,
    QgsPointXY,
    QgsCoordinateReferenceSystem,
    QgsMessageLog
)

vectors = {
    "N": (0, 1, 0),
    "E": (1, 0, 0),
    "S": (0, -1, 0),
    "W": (-1, 0, 0),
    "L": (0, 0, -1),
    "R": (0, 0, 1),
}

def process_command(command, ship, waypoint):
    code = command[0]
    amount = int(command[1:])

    if code == "F":
        ship[0] += amount*waypoint[0]
        ship[1] += amount*waypoint[1]
    else:
        easting, northing, rotate = vectors[code]
        angle = radians(-amount*rotate)
        old_e, old_n = waypoint
        waypoint[0] = amount*easting + old_e*cos(angle) - old_n*sin(angle)
        waypoint[1] = amount*northing + old_e*sin(angle) + old_n*cos(angle)

def main():
    data_file = QFileDialog.getOpenFileName(None, "Open data file")
    if not data_file:
        return

    QgsMessageLog.logMessage(f"Opening file {data_file}", "AoC", Qgis.Info)

    # create vector layers
    vl_points = QgsVectorLayer("Point", "AoC Points",  "memory")
    pr_points = vl_points.dataProvider()
    pr_points.addAttributes([
        QgsField("idx", QVariant.Int),
        QgsField("cmd", QVariant.String),
        QgsField("heading", QVariant.Double),
    ])
    vl_points.updateFields()

    vl_lines = QgsVectorLayer("LineString", "AoC Lines",  "memory")
    pr_lines = vl_lines.dataProvider()

    # read data
    data = open(data_file[0]).readlines()
    QgsMessageLog.logMessage(f"{data}", "AoC", Qgis.Info)

    heading = 90 # ship starts facign east
    ship = [0, 0]
    waypoint = [10, 1]
    line = [QgsPoint(ship[0], ship[1])]

    for idx, entry in enumerate(data):
        # plot feature
        if entry[0] == "F":
            line.append(QgsPoint(ship[0], ship[1]))

            feat = QgsFeature()
            feat.setGeometry(QgsGeometry.fromPointXY(QgsPointXY(ship[0], ship[1])))
            heading = degrees(atan2(waypoint[0], waypoint[1]))
            feat.setAttributes([idx, entry, heading])
            pr_points.addFeatures([feat])

        # add feature
        process_command(entry, ship, waypoint)
        QgsMessageLog.logMessage(f"cmd[{idx}]:{entry} E{ship[0]:.0f}, N{ship[1]:.0f}, WE{waypoint[0]:.0f}, WN{waypoint[1]:.0f}", "AoC", Qgis.Info)

    # last point
    line.append(QgsPoint(ship[0], ship[1]))

    feat = QgsFeature()
    feat.setGeometry(QgsGeometry.fromPointXY(QgsPointXY(ship[0], ship[1])))
    heading = degrees(atan2(waypoint[0], waypoint[1]))
    feat.setAttributes([idx, entry, heading])
    pr_points.addFeatures([feat])

    vl_points.updateExtents()

    # construct linestring
    feat = QgsFeature()
    feat.setGeometry(QgsGeometry.fromPolyline(line))
    pr_lines.addFeatures([feat])
    vl_lines.updateExtents()

    QgsProject.instance().addMapLayer(vl_points)
    QgsProject.instance().addMapLayer(vl_lines)

    # disable crs (we need flat world)
    crs = QgsCoordinateReferenceSystem(None)
    QgsProject.instance().setCrs(crs)

main()
이렇게 하면 다음과 같은 시각 형상이 생성됩니다.

화살표가 현재 항로점 벡터를 가리키기;텍스트 레이블은 명령 번호일 뿐입니다.
도전에 필요한 답안은 로그 메시지에서 출력된 배의 최종 위치를 간단하게 찾을 수 있다.
앞으로!

좋은 웹페이지 즐겨찾기