Daha çox

PostGIS-də poliqon örtüklərinin zonal statistikası

PostGIS-də poliqon örtüklərinin zonal statistikası


Görünən sadə bir tapşırıq yerinə yetirməyə çalışıram - PostGIS-də bir poliqon örtükünə əsaslanan bir rasterdən ortalama piksel dəyərlərini çıxarın. Bu blog, bu blog və bu yazı da daxil olmaqla bir neçə mənbəyə müraciət etdim. Verilənlər bazasına Landsat kafelini (in-db) və üç poliqondan (ID 1,2,3) ibarət bir çoxbucaqlı örtük yüklədim. Landsat kafel, 15 m məkan qətnaməsinə və 4 banta sahib bir pansharped şəkildir.

Kömək sənədləri və yuxarıda göstərilən yazı / suallar əsasında aşağıdakı PostGIS sorğusuna gəldim:

SEÇİN id, (SUM ((ST_SummaryStats (a.rast, 1, true)). Sum) / SUM ((ST_SummaryStats (a.rast, 1, true)). Count)) FROM imagery.l8_2015_09_12 AS a, analyz_results .zonal_stats_test kimi b WHERE ST_Intersects (b.geom, a.rast) GROUP BY id;

Nəticələr:

1: 9902.49 2: 10079.68 3: 12355.90

Lakin ArcGIS-in "zona statistikası cədvəl" alətinə qarşı nəticələri sınadığım zaman aşağıdakı nəticələr əldə edirəm:

1: 8089.61 2: 7527.62 3: 12290.05

Yaxın, lakin eyni deyil. Demək olar ki, varam, amma kimsə mənə hər iki texnika arasındakı fərqin mənbəyini anlamağa kömək edə biləcəyini düşünürdüm. ST_summaryStats metodumda səhv bir şey olduğuna şübhə edirəm.


1) Rasteri statistikalarını hesablamadan əvvəl ST_Clip (rast, geom) və ya ST_Intersection (rast, geom) ilə "kəsmədiyinizə" təəccüb edirəm. İstədiyiniz nə çoxbucaqlıları kəsişən piksel dəyərlərindəki statistikalardır? Fərqli bir nəticə ilə nəticələnməyinizin bir səbəbi bu ola bilər, ancaq nəticənin yalnız kiçik hissələrinin kəsilməsi olardı, çünki nəticəniz ArcGIS nəticəsinə çox yaxındır.

2) Statistikanı hesablamadan əvvəl rasteri düzgün bir şəkildə "kəssəniz" də, PostGIS klipinin və ya rasterin kəsişməsinin yolu ArcGIS-dən fərqli ola bilər: kəsişən piksellər necə seçilir? Centroid ilə? Ümumi kəsişən əraziyə görə? Nodata dəyərləri necə işlənir? Əgər belədirsə, kəsmə / kəsişmədən əvvəl çoxbucaqlar rasterə necə çevrilir? Əgər vəziyyət belədirsə, rastrlar çoxbucaqlılara necə çevrilir? Bütün bunlar seçdiyiniz alətlərdən və bəzən əməliyyatları həyata keçirmə qaydasından asılıdır.

ArcGIS ilə fərqi anlamaq üçün məlumatlarınızın bir hissəsi ilə işləməyinizi və nəticələrini addım-addım görüntüləməyinizi / təhlil etməyinizi (məsələn OpenJump-da) təklif edirəm.

Bu təlimata da baxın.


Landsat rasterini verilənlər bazasına idxal etdiyiniz zaman, şəkli 300x300 piksel plitələrə bölən bir -t 300x300 parametri qurduğunuz vəziyyət ola bilər. Beləliklə, poliqonlarınızı verilənlər bazanıza idxal etdiyiniz bir şəkil şəklində olduğunu düşünün.

Çoxbucaqlı bir xüsusiyyəti əhatə etmək üçün birdən çox raster plitəyə ehtiyacınız ola bilər, buna görə də bu məsələnin qarşısını almaq üçün sorğunuza St_union funksiyasını əlavə etməlisiniz.

SEÇİN id, (St_SummaryStats (St_Union (ST_Clip (rast, 1, geom, true)))). * * Görüntüdən.L8_2015_09_12, analyz_results.zonal_stats_test HARADA st_intersects (rast, geom) GROUP BY id;

Ayrıca. * İstifadə edərək bonus olaraq cəmi, say, max, min, stddev əldə edirsiniz. Yalnız orta istəsəniz, ortalama olaraq dəyişdirin *, bu sətir belə nəticələnəcəkdir:

(St_SummaryStats (St_Union (ST_Clip (rast, 1, geom, true)))) deməkdir.

Videoya baxın: QGIS Tutorial: Connect to PostGIS database EN