Jak obliczyć odległość do obiektu za pomocą gdal_proximity?

28

Korzystam z gdal_proximity, aby znaleźć odległość do najbliższej głównej rzeki w USA (48 stanów niższych). Przewidziałem fluktuacje sieci NHD + do Conusa Albersa (epsg: 5070), wybrane rzeki z porządkiem strumienia> 5 i zrasteryzowane, płonące rzeki jako 255, bez rzeki jako 0. To w porządku, ale teraz muszę znaleźć odległość do najbliższej rzeki dla miejsc w promieniu 50 km. Plik wejściowy ma rozdzielczość 30m w skali kontynentalnej, więc jest bardzo duży, ale konwersja powinna być prostym poleceniem gdal_proximity:

gdal_proximity.bat -values 255 -distunits GEO -maxdist 50000 -nodata -999 infile.tif outfile.tif -co COMPRESS=DEFLATE -co BIGTIFF=YES -co TILED=YES

Wydaje się, że to prawie działa, ale generuje dziwny geometryczny wzór na wyjściu (patrz zdjęcie). Dane obecne w danych wyjściowych zostały poprawnie przetworzone. Czy ktoś może zasugerować, dlaczego brakuje tak dużej ilości danych wyjściowych?

Odległość od rzeki

Edycja: Aby sprawdzić, czy jest to spowodowane przez dowolny z opcjonalnych parametrów, ponownie uruchomiłem gdal_proximity w tej konfiguracji:

gdal_proximity.bat H:\data\tmp\NHDplus_network_flowline_SO6plus.tif H:/data/tmp/NHDplus_network_flowline_SO6plus_proximity.tif -values 255 -maxdist 50000 -of GTiff

Co przyniosło zasadniczo ten sam wynik:

Odległość od rzeki, brak opcjonalnych parametrów

Moją jedyną myślą jest to, że może to być związane z rozmiarem rastra (~ 100 gb nieskompresowany. O ile mi wiadomo, nie ma ograniczenia wielkości BigTiff, ale może istnieje ograniczenie tego, co gdal może analizować skutecznie?

qgis gdal qgis-processing R Rhodes
źródło
1
co się stanie, jeśli wyłączysz kafelek = TAK? Czy to działa, jeśli zmienisz GEO na PIXEL? (Wynik może być nieodpowiedni, ale może zawęzić problem)
Steven Kay,
Dziękuję za sugestię - dodałem odpowiedź na pierwotne pytanie.
R Rhodes,
W jakiej rozdzielczości jest twój plik infile.tif?
shahryar
2
Czy możesz spróbować odczytać dane za pomocą GDAL partiami (liniami) i sprawdzić, czy problemem są same dane, czy QGIS nie jest w stanie ich wizualizować? Pierwszym krokiem do znalezienia tego problemu jest zmniejszenie zasięgu przestrzennego do przykładowego AOI.
RutgerH

Odpowiedzi:

3

Podejrzewam, że osiągasz gdzieś limit pamięci, być może kiedy pamięć RAM jest wyczerpana, a system operacyjny zrzuca plik stronicowania. Monitoruj zasoby systemowe podczas procesu. Nie jest dla mnie jasne, dlaczego twoje wyniki występują w zakrzywionych pokosach, ale upewnij się, że wszystkie dane zostały rzutowane (zapisane) do tego samego układu współrzędnych.

Przyjrzyjmy się numerycznym typom danych , aby pomóc temu algorytmowi. Sieć zrasteryzowanego strumienia musi zawierać tylko wartości binarne, abyśmy mogli zaoszczędzić na zasobach, używając Bytetypu danych rastrowych. Wypal wartość 1 dla strumieni i 0 dla tła:

gdal_rasterize -l streams -burn 1 -tr 50 50 -a_nodata 0 -te -2339101 311625 2227004 3134200 -ot Byte -of GTiff streams.shp streams.tif

Następnie interesująca nas bliskość jest dodatnia i mniejsza lub równa 50 000 m. Odpowiednim typem danych jest 16-bitowa liczba całkowita bez znaku UInt16. Ponadto, jeśli ustawimy „brak danych” na maks. 65535, możemy zachować wartość 0 dla komórek strumienia.

Jeśli to konieczne, możesz również obniżyć do 8-bitowej liczby całkowitej bez znaku UInt8i nadal mieć dokładność zbliżoną do ~ 200 m.

gdal_proximity.bat -srcband 1 -distunits GEO -values 1 -maxdist 50000 -nodata 65535 -ot UInt16 -of GTiff streams.tif proximity.tif

* Uwaga: użyłem komórki o wielkości 50 m. Gdal_proximity zużył ~ 20 GB pamięci RAM i zajął ~ 5 minut na moim komputerze. Jeśli masz ograniczoną pamięć RAM, podziel raster wejściowy na możliwe do zarządzania rozmiary, jak zauważyli inni.

wyniki gdal_proximity

CrystallineEntity
źródło