Menghasilkan poligon berukuran sama dengan PyQGIS?

Saya ingin membuat poligon di sepanjang baris untuk menggunakannya untuk AtlasCreator di langkah selanjutnya.

ArcMap memiliki alat yang disebut Fitur Strip Map Index .

Dengan alat ini saya dapat memilih tinggi dan lebar poligon saya (katakan 8km x 4km) dan menghasilkan / memutarnya sepanjang garis secara otomatis.

Salah satu atribut yang dihasilkan dari setiap poligon adalah sudut rotasi yang saya perlukan untuk memutar panah utara saya di Generator Atlas sesudahnya.

Adakah yang punya ide bagaimana menyelesaikan tugas ini di QGIS / dengan pyQGIS? Algoritma Grass- atau SAGA atau prossessing-toolbox-model yang dapat digunakan di dalam plugin kustom juga tidak masalah;) Edit1: Saya tidak hanya memerlukan luasan cetak tetapi juga poligon itu sendiri karena saya ingin mencetak peta dengan semua poligon / luasan sebagai semacam peta ikhtisar.

Sunting2: Saya menawarkan hadiah karena saya masih mencari solusi PyQGIS yang dapat digunakan dalam QGIS-Plugin tanpa perlu menginstal perangkat lunak selain dari QGIS (tidak ada RDBMS seperti PostGIS / Oracle)

Pertanyaan menarik! Ini adalah sesuatu yang ingin saya coba sendiri, jadi cobalah.

Anda dapat melakukan ini di PostGRES / POSTGIS dengan fungsi yang menghasilkan seperangkat poligon.

Dalam kasus saya, saya memiliki tabel dengan satu fitur (MULTILINESTRING) yang mewakili jalur kereta api. Perlu menggunakan CRS dalam meter, saya menggunakan osgb (27700). Saya telah melakukan 4km x 2km 'halaman'.

Di sini, Anda dapat melihat hasilnya ... hal-hal hijau adalah jaringan jalan, terpotong ke penyangga 1 km di sekitar rel, yang sesuai dengan ketinggian poligon dengan baik.

Inilah fungsinya ...

CREATE OR REPLACE FUNCTION getAllPages(wid float, hite float, srid integer, overlap float) RETURNS SETOF geometry AS
$BODY$
DECLARE
    page geometry; -- holds each page as it is generated
    myline geometry; -- holds the line geometry
    startpoint geometry;
    endpoint geometry;
    azimuth float; -- angle of rotation
    curs float := 0.0 ; -- how far along line left edge is
    step float;
    stepnudge float;
    currpoly geometry; -- used to make pages
    currline geometry;
    currangle float;
    numpages float;
BEGIN
    -- drop ST_LineMerge call if using LineString 
    -- replace this with your table.
    SELECT ST_LineMerge(geom) INTO myline from traced_osgb; 
    numpages := ST_Length(myline)/wid;

    step := 1.0/numpages;
    stepnudge := (1.0-overlap) * step; 
    FOR r in 1..cast (numpages as integer)
    LOOP
        -- work out current line segment

        startpoint :=  ST_SetSRID(ST_Line_Interpolate_Point(myline,curs),srid);
        endpoint :=  ST_SetSRID(ST_Line_Interpolate_Point(myline,curs+step),srid);
        currline := ST_SetSRID(ST_MakeLine(startpoint,endpoint),srid);

        -- make a polygon of appropriate size at origin of CRS
        currpoly := ST_SetSRID(ST_Extent(ST_MakeLine(ST_MakePoint(0.0,0.0),ST_MakePoint(wid,hite))),srid);

        -- then nudge downwards so the midline matches the current line segment
        currpoly := ST_Translate(currpoly,0.0,-hite/2.0);

        -- Rotate to match angle
        -- I have absolutely no idea how this bit works. 
        currangle := -ST_Azimuth(startpoint,endpoint) - (PI()/2.0) + PI();
        currpoly := ST_Rotate(currpoly, currangle);

        -- then move to start of current segment
        currpoly := ST_Translate(currpoly,ST_X(startpoint),ST_Y(startpoint));

        page := currpoly;

        RETURN NEXT page as geom; -- yield next result
        curs := curs + stepnudge;
    END LOOP;
    RETURN;
END
$BODY$
LANGUAGE 'plpgsql' ;

Menggunakan fungsi ini

Ini sebuah contoh; 4km x 2km halaman, epsg: 27700 dan 10% tumpang tindih

select st_asEwkt(getallpages) from getAllPages(4000.0, 2000.0, 27700, 0.1);

Setelah menjalankan ini, Anda dapat mengekspor dari PgAdminIII ke file csv. Anda dapat mengimpor ini ke QGIS, tetapi Anda mungkin perlu mengatur CRS secara manual untuk layer - QGIS tidak menggunakan SRID di EWKT untuk mengatur layer CRS untuk Anda: /

Menambahkan atribut bearing

Ini mungkin lebih mudah dilakukan dalam postgis, bisa dilakukan dalam ekspresi QGIS tetapi Anda harus menulis beberapa kode. Sesuatu seperti ini...

create table pages as (
    select getallpages from getAllPages(4000.0, 2000.0, 27700, 0.1)
);

alter table pages add column bearing float;

update pages set bearing=ST_Azimuth(ST_PointN(getallpages,1),ST_PointN(getallpages,2));

Peringatan

Agak sedikit diretas, dan hanya memiliki kesempatan untuk menguji pada satu dataset.

Tidak 100% yakin dua simpul mana yang harus Anda pilih pada pembaruan atribut bantalan itu query.. mungkin perlu bereksperimen.

Saya harus mengakui bahwa saya tidak tahu mengapa saya perlu melakukan formula berbelit-belit untuk memutar poligon agar sesuai dengan segmen garis saat ini. Saya pikir saya bisa menggunakan output dari ST_Azimuth () di ST_Rotate (), tetapi sepertinya tidak.

Ada beberapa solusi yang berbeda. Dan ini dapat bekerja dengan polyline sederhana dan beberapa entitas terpilih

diagram blok:

1. Parameter

   1. pilih orientasi untuk pembangkitan dan indeks baca (kiri-ke-kanan, utara-ke-selatan ...)
   2. atur ukuran objek

   shape = (4000,8000) # (<width>,<length>)

   3. define superposition coef (10% secara default?)
2. init
   1. Memesan polyline (bandingkan titik awal dan akhir) memesan tergantung pada pilihan orientasi Anda> membuat simpul memesan Featureeclass OrderNodes

3. loop pada OrderNodes

   1. buat Anda titik pertama sebagai jangkar

   2. untuk setiap titik menambahkannya pada dict x, y, id dan menghitung vektor

   3. menghasilkan poligon (lebih dari panjang dan orientasi vektor) dengan mengurangi superposisi (10% / 2)> poligon kiri 5% poligon kanan 5% dengan titik jangkar yang sama
   4. Berhenti ketika titik titik preseden keluar dari poligon atau jika vektor len> untuk membentuk panjang
   5. Buat poligon dengan solusi yang baik sebelumnya dan atur titik jangkar dengan posisi bagus terakhir
   6. Lakukan perulangan baru dan setel ulang dict x, y, id untuk menghasilkan objek poligon berikutnya.

Anda dapat mengubah proposisi ini jika tidak benar-benar jelas atau berkomentar.

Steven Kays menjawab dalam pyqgis. Cukup pilih garis di lapisan Anda sebelum menjalankan skrip. Script tidak mendukung linemerging sehingga tidak dapat bekerja pada layer dengan multilinestring

#!python
# coding: utf-8

# https://gis.stackexchange.com/questions/173127/generating-equal-sized-polygons-along-line-with-pyqgis
from qgis.core import QgsMapLayerRegistry, QgsGeometry, QgsField, QgsFeature, QgsPoint
from PyQt4.QtCore import QVariant


def getAllPages(layer, width, height, srid, overlap):
    for feature in layer.selectedFeatures():
        geom = feature.geometry()
        if geom.type() <> QGis.Line:
            print "Geometry type should be a LineString"
            return 2
        pages = QgsVectorLayer("Polygon?crs=epsg:"+str(srid), 
                      layer.name()+'_id_'+str(feature.id())+'_pages', 
                      "memory")
        fid = QgsField("fid", QVariant.Int, "int")
        angle = QgsField("angle", QVariant.Double, "double")
        attributes = [fid, angle]
        pages.startEditing()
        pagesProvider = pages.dataProvider()
        pagesProvider.addAttributes(attributes)
        curs = 0
        numpages = geom.length()/(width)
        step = 1.0/numpages
        stepnudge = (1.0-overlap) * step
        pageFeatures = []
        r = 1
        currangle = 0
        while curs <= 1:
            # print 'r =' + str(r)
            # print 'curs = ' + str(curs)
            startpoint =  geom.interpolate(curs*geom.length())
            endpoint = geom.interpolate((curs+step)*geom.length())
            x_start = startpoint.asPoint().x()
            y_start = startpoint.asPoint().y()
            x_end = endpoint.asPoint().x()
            y_end = endpoint.asPoint().y()
            # print 'x_start :' + str(x_start)
            # print 'y_start :' + str(y_start)
            currline = QgsGeometry().fromWkt('LINESTRING({} {}, {} {})'.format(x_start, y_start, x_end, y_end))
            currpoly = QgsGeometry().fromWkt(
                'POLYGON((0 0, 0 {height},{width} {height}, {width} 0, 0 0))'.format(height=height, width=width))
            currpoly.translate(0,-height/2)
            azimuth = startpoint.asPoint().azimuth(endpoint.asPoint())
            currangle = (startpoint.asPoint().azimuth(endpoint.asPoint())+270)%360
            # print 'azimuth :' + str(azimuth)
            # print 'currangle : ' +  str(currangle)

            currpoly.rotate(currangle, QgsPoint(0,0))
            currpoly.translate(x_start, y_start)
            currpoly.asPolygon()
            page = currpoly
            curs = curs + stepnudge
            feat = QgsFeature()
            feat.setAttributes([r, currangle])
            feat.setGeometry(page)
            pageFeatures.append(feat)
            r = r + 1

        pagesProvider.addFeatures(pageFeatures)
        pages.commitChanges()
        QgsMapLayerRegistry.instance().addMapLayer(pages)
    return 0

layer = iface.activeLayer()
getAllPages(layer, 500, 200, 2154, 0.4)

Dua jawaban (pada saat posting) adalah cerdik dan dijelaskan dengan baik. Namun, ada juga solusi SANGAT sederhana namun efektif yang memungkinkan untuk ini (dengan asumsi Anda akan menerima semua peta Anda sejajar dengan Utara di jalan tradisional, daripada arah Utara acak berdasarkan sungai). Jika Anda ingin rotasi, itu mungkin tetapi sedikit lebih rumit (lihat bagian bawah).

Pertama lihat posting saya di sini . Ini memberi Anda cara untuk membuat cakupan peta untuk Atlas. Metode yang Anda inginkan adalah adaptasi dari 'Workflow 2' di caranya. Membagi fitur linier Anda dengan simpul atau panjang dan buffer fitur dengan jumlah berapa pun. Jumlah yang Anda buffer sebagian akan menentukan tumpang tindih (tapi lihat di bawah) tetapi yang lebih penting, itu menciptakan fitur dengan suatu area. Anda dapat menggunakan sejumlah plugin untuk membagi baris tetapi GRASS v.split.length dan v.split.vert adalah opsi yang baik (tersedia di Processing Toolbox).

Setelah mengaktifkan Atlas Generation di Map Composer dan memilih layer buffered Anda, kembali ke tab item dan pilih objek peta Anda. Lihat 'Dikendalikan oleh Atlas', dan dalam kasus penggunaan Anda, saya akan memilih fitur Margin. Ini akan mengontrol tumpang tindih Anda di antara peta (atau Anda dapat memilih skala tetap).

Anda dapat mempratinjau Atlas Anda menggunakan Tombol Pratinjau Atlas di bilah alat atas penggubah dan melihat berapa banyak halaman yang akan dihasilkan. Catatan Anda dapat memilih untuk mengekspor semua halaman dalam satu PDF atau sebagai file terpisah.

Untuk membuat peta berputar di sepanjang garis, ada bidang rotasi di properti Item Komposer Peta. Anda perlu mengatur ekspresi (gunakan tombol kecil di sebelah kanan kotak rotasi). Pilih variabel sebagai opsi Anda dan kemudian Edit. Pembuat ekspresi akan muncul dan di sana Anda dapat mengakses geometri atau bidang fitur atlas. Anda kemudian dapat membangun ekspres untuk memutar peta sesuai dengan rotasi fitur (Anda dapat menghitung bantalan menggunakan titik awal dan akhir dari setiap segmen garis dan sedikit pemicu). Ulangi proses yang sama untuk memutar panah Utara Anda (menggunakan ekspresi yang sama atau variabel yang dihitung sebelumnya).