Równoległe operacje GIS w PyQGIS?

Częstym wymogiem w GIS jest zastosowanie narzędzia przetwarzania do wielu plików lub zastosowanie procesu dla wielu funkcji w jednym pliku do innego pliku.

Wiele z tych operacji jest krępująco równoległych, ponieważ wyniki obliczeń w żaden sposób nie wpływają na żadną inną operację w pętli. Nie tylko to, ale często pliki wejściowe są różne.

Klasycznym przykładem jest kafelkowanie plików kształtów w stosunku do plików zawierających wielokąty w celu przycięcia ich.

Oto (przetestowana) klasyczna metoda proceduralna umożliwiająca osiągnięcie tego w skrypcie Python dla QGIS. (fyi wyjściowe pliki pamięci tymczasowej do prawdziwych plików ponad połowę czasu na przetwarzanie moich plików testowych)

import processing
import os
input_file="/path/to/input_file.shp"
clip_polygons_file="/path/to/polygon_file.shp"
output_folder="/tmp/test/"
input_layer = QgsVectorLayer(input_file, "input file", "ogr")
QgsMapLayerRegistry.instance().addMapLayer(input_layer)
tile_layer  = QgsVectorLayer(clip_polygons_file, "clip_polys", "ogr")
QgsMapLayerRegistry.instance().addMapLayer(tile_layer)
tile_layer_dp=input_layer.dataProvider()
EPSG_code=int(tile_layer_dp.crs().authid().split(":")[1])
tile_no=0
clipping_polygons = tile_layer.getFeatures()
for clipping_polygon in clipping_polygons:
    print "Tile no: "+str(tile_no)
    tile_no+=1
    geom = clipping_polygon.geometry()
    clip_layer=QgsVectorLayer("Polygon?crs=epsg:"+str(EPSG_code)+\
    "&field=id:integer&index=yes","clip_polygon", "memory")
    clip_layer_dp = clip_layer.dataProvider()
    clip_layer.startEditing()
    clip_layer_feature = QgsFeature()
    clip_layer_feature.setGeometry(geom)
    (res, outFeats) = clip_layer_dp.addFeatures([clip_layer_feature])
    clip_layer.commitChanges()
    clip_file = os.path.join(output_folder,"tile_"+str(tile_no)+".shp")
    write_error = QgsVectorFileWriter.writeAsVectorFormat(clip_layer, \
    clip_file, "system", \
    QgsCoordinateReferenceSystem(EPSG_code), "ESRI Shapefile")
    QgsMapLayerRegistry.instance().addMapLayer(clip_layer)
    output_file = os.path.join(output_folder,str(tile_no)+".shp")
    processing.runalg("qgis:clip", input_file, clip_file, output_file)
    QgsMapLayerRegistry.instance().removeMapLayer(clip_layer.id())

Byłoby dobrze, gdyby mój plik wejściowy miał 2 GB, a plik przycinający wielokąt zawiera ponad 400 wielokątów. Wynikowy proces trwa ponad tydzień na mojej czterordzeniowej maszynie. Przez cały czas trzy rdzenie są na biegu jałowym.

Rozwiązaniem, które mam w głowie, jest wyeksportowanie procesu do plików skryptów i uruchomienie ich asynchronicznie, na przykład przy użyciu równoległego GNU. Jednak szkoda, że ​​trzeba zrezygnować z QGIS w rozwiązaniu specyficznym dla systemu operacyjnego, zamiast używać czegoś natywnego dla Pythona QGIS. Więc moje pytanie brzmi:

Czy mogę równolegle zawstydzać równoległe operacje geograficzne natywnie w Pythonie QGIS?

Jeśli nie, to może ktoś już ma kod, aby wysłać tego rodzaju pracę do asynchronicznych skryptów powłoki?

Jeśli zmienisz program tak, aby czytał nazwę pliku z wiersza poleceń i podzieliłeś plik wejściowy na mniejsze części, możesz zrobić coś takiego za pomocą GNU Parallel:

parallel my_processing.py {} /path/to/polygon_file.shp ::: input_files*.shp

Spowoduje to uruchomienie 1 zadania na rdzeń.

Wszystkie nowe komputery mają wiele rdzeni, ale większość programów ma charakter szeregowy i dlatego nie będzie używać wielu rdzeni. Jednak wiele zadań można bardzo zrównoleglać:

• Uruchom ten sam program na wielu plikach
• Uruchom ten sam program dla każdej linii w pliku
• Uruchom ten sam program dla każdego bloku w pliku

GNU Parallel to ogólny paralelizator, który ułatwia równoległe uruchamianie zadań na tej samej maszynie lub na wielu maszynach, do których masz dostęp ssh.

Jeśli masz 32 różne zadania, które chcesz uruchomić na 4 procesorach, prostym sposobem na zrównoleglenie jest uruchomienie 8 zadań na każdym procesorze:

Zamiast tego GNU Parallel odradza nowy proces po zakończeniu - utrzymując procesory aktywne, a tym samym oszczędzając czas:

Instalacja

Jeśli GNU Parallel nie jest spakowany dla twojej dystrybucji, możesz wykonać osobistą instalację, która nie wymaga dostępu do roota. Można to zrobić w 10 sekund, wykonując następujące czynności:

(wget -O - pi.dk/3 || curl pi.dk/3/ || fetch -o - http://pi.dk/3) | bash

Inne opcje instalacji można znaleźć na stronie http://git.savannah.gnu.org/cgit/parallel.git/tree/README

Ucz się więcej

Zobacz więcej przykładów: http://www.gnu.org/software/parallel/man.html

Obejrzyj filmy wprowadzające: https://www.youtube.com/playlist?list=PL284C9FF2488BC6D1

Przejrzyj samouczek: http://www.gnu.org/software/parallel/parallel_tutorial.html

Zarejestruj się, aby otrzymać listę e-mail, aby uzyskać pomoc: https://lists.gnu.org/mailman/listinfo/parallel

Zamiast korzystać z metody GNU Parallel, możesz użyć modułu mutliprocess Pythona , aby utworzyć pulę zadań i wykonać je. Nie mam dostępu do konfiguracji QGIS, aby ją przetestować, ale w Pythonie 2.6 dodano proces wieloprocesowy, więc jeśli używasz wersji 2.6 lub nowszej, powinna być dostępna. Istnieje wiele przykładów online dotyczących korzystania z tego modułu.

Oto równoległe rozwiązanie GNU. Z pewną ostrożnością można wprowadzić najbardziej zawstydzająco równoległe algorytmy ogr lub saga oparte na systemie Linux do pracy z nim w instalacji QGIS.

Oczywiście to rozwiązanie wymaga instalacji GNU równolegle. Aby na przykład zainstalować GNU równolegle w Ubuntu, przejdź do terminala i wpisz

sudo apt-get -y install parallel

NB: Nie mogłem uruchomić polecenia powłoki równoległej w Popen lub podprocesie, co wolałbym, więc zhakowałem eksport do skryptu bash i zamiast tego uruchomiłem to z Popenem.

Oto specyficzne polecenie powłoki przy użyciu równoległego, które zawinąłem w python

parallel ogr2ogr -skipfailures -clipsrc tile_{1}.shp output_{1}.shp input.shp ::: {1..400}

Każde {1} zostaje zamienione na liczbę z zakresu {1..400}, a następnie czterysta poleceń powłoki jest zarządzanych przez GNU równolegle, aby jednocześnie używać wszystkich rdzeni mojego i7 :).

Oto kod napisany przeze mnie, aby rozwiązać przykładowy problem, który opublikowałem. Można go wkleić bezpośrednio po końcu kodu w pytaniu.

import stat
from subprocess import Popen
from subprocess import PIPE
feature_count=tile_layer.dataProvider().featureCount()
subprocess_args=["parallel", \
"ogr2ogr","-skipfailures","-clipsrc",\
os.path.join(output_folder,"tile_"+"{1}"+".shp"),\
os.path.join(output_folder,"output_"+"{1}"+".shp"),\
input_file,\
" ::: ","{1.."+str(feature_count)+"}"]
#Hacky part where I write the shell command to a script file
temp_script=os.path.join(output_folder,"parallelclip.sh")
f = open(temp_script,'w')
f.write("#!/bin/bash\n")
f.write(" ".join(subprocess_args)+'\n')
f.close()
st = os.stat(temp_script)
os.chmod(temp_script, st.st_mode | stat.S_IEXEC)
#End of hacky bash script export
p = Popen([os.path.join(output_folder,"parallelclip.sh")],\
stdin=PIPE, stdout=PIPE, stderr=PIPE)
#Below is the commented out Popen line I couldn't get to work
#p = Popen(subprocess_args, stdin=PIPE, stdout=PIPE, stderr=PIPE)
output, err = p.communicate(b"input data that is passed to subprocess' stdin")
rc = p.returncode
print output
print err

#Delete script and old clip files
os.remove(os.path.join(output_folder,"parallelclip.sh"))
for i in range(feature_count):
    delete_file = os.path.join(output_folder,"tile_"+str(i+1)+".shp")
    nosuff=os.path.splitext(delete_file)[0]
    suffix_list=[]
    suffix_list.append('.shx')
    suffix_list.append('.dbf')
    suffix_list.append('.qpj')
    suffix_list.append('.prj')
    suffix_list.append('.shp')
    suffix_list.append('.cpg')
    for suffix in suffix_list:
        try:
            os.remove(nosuff+suffix)
        except:
            pass

Pozwól mi powiedzieć, że to naprawdę coś, gdy zobaczysz, że wszystkie rdzenie zapalają się do pełnego hałasu :). Specjalne podziękowania dla Ole i zespołu, który zbudował Gnu Parallel.

Byłoby miło mieć rozwiązanie wieloplatformowe i byłoby miło, gdybym mógł wymyślić wieloprocesowy moduł python dla osadzonego pytona qgis, ale niestety tak nie było.

Niezależnie od tego to rozwiązanie będzie dla mnie przydatne, a może i dla Ciebie.