Bu atölyede SRTM raster verileri işlenerek bazı temel analizler gerçekleştirilecektir.
Kullanılan kütüphaneler: rasterio, geopandas, rasterstats
- rasterio: mekansal raster veri işlemleri için
- geopandas: mekansal veri operasyonları için
- rasterstats: vektör geometriler ile raster verilerden istatistik çıkarmak için
- contextily: internetten altık haritaları çekmek için
-
https://search.earthdata.nasa.gov/search/granules?p=C1546314043-LPDAAC_ECS adresinden NASADEM Merged DEM Global 1 arc second V001 veri seti araması ile sağ taraftan indireceğimiz bölgenin alanını kare içine alarak arama yapıyoruz ve ilgili bölgenin .zip uzantılı .hgt yükseklik veri setini indiriyoruz.
-
İndirdiğimiz bölgeye ait ilçe sınırlarını vektör veri formatında https://gadm.org/ sitesinden indiriyoruz. Level-2 olan veri ilçeleri içermektedir. GeoJSON formatında Level-2 veri bu çalışmada yeterli olacaktır.
-
Level-2 Türkiye GeoJSON Verisi: https://geodata.ucdavis.edu/gadm/gadm4.1/json/gadm41_TUR_2.json.zip
-
İndirdiğimiz ilçe verilerini okuyup görselleştirelim.
# read turkey level-2 geojson
ilceler = gpd.read_file('gadm41_TUR_2.json')
ilceler.plot()
- İstanbul, Bursa ve Kocaeli illerini seçelim
secileniller = ilceler[(ilceler['NAME_1'] == 'Istanbul') | (ilceler['NAME_1'] == 'Bursa') | (ilceler['NAME_1'] == 'Kocaeli')]
- Bursa'nın, Karamürsel, Orhangazi ve İznik ilçelerini seçelim.
bursailceler = ilceler[(ilceler['NAME_2'] == 'Karamürsel') | (ilceler['NAME_2'] == 'Orhangazi') | (ilceler['NAME_2'] == 'İznik')]
- Seçilen ilçeleri çizdirelim. Bu ilçelere göre sayısal yükseklik modelimizi keseceğiz.
bursailceler.plot()
- Şimdi sayısal yükseklik modeli verimizi öncelikle rasterio ile okuyup, rioshow ile gösterip, contextily'dan alacağımız arkaplan harita ile sunalım. Üzerine'de seçtiğimiz illerin katmanı nı ekleyelim.
r = rasterio.open("n40e029.hgt")
f, ax = plt.subplots(1, figsize=(9, 9))
ax = rioshow(r, alpha=0.5, zorder=2, ax=ax)
cx.add_basemap(ax, crs=r.crs)
secileniller.plot(facecolor="none", edgecolor="#F93822", zorder=1, ax=ax)
- Seçtiğimiz Bursa ilçeleri için görüntüyü kesmeden önce bir kere daha görselleştirelim.
r = rasterio.open("n40e029.hgt")
f, ax = plt.subplots(1, figsize=(9, 9))
ax = rioshow(r, alpha=0.5, zorder=2, ax=ax)
cx.add_basemap(ax, crs=r.crs)
bursailceler.plot(facecolor="none", edgecolor="#F93822", zorder=1, ax=ax)
- Şimdi bölgesel istatistik hesabı için rasterstats kütüphanesinden zonal_stats fonksiyonu ile ilçelerimiz için min, max, median gibi değerleri hesaplayalım.
from rasterstats import zonal_stats
elevations2 = zonal_stats(bursailceler, 'n40e029.hgt', stats="count min mean max median")
elevations2 = pandas.DataFrame(elevations2)
-
Son aşamada sayısal yükseklik verisini bu ilçe sınırlarına göre keserek bir bir .tif uzantılı dosya üreteceğiz.
-
Bunun için öncelikle ilçeler verimizi geojson olarak kaydedelim.
bursailceler.to_file("bursa.geojson", driver="GeoJSON")
- rio mask komutu ile görüntümüzü bu vektör dosya sınırına göre keserek yeni bir dosya yazdıralım.
! rio mask n40e029.hgt \
bursa.tif \
--crop \
--geojson-mask\
bursa.geojson
- Sonrasında görselleştirelim.
f, ax = plt.subplots(1, figsize=(9, 9))
r = rasterio.open("bursa.tif")
rioshow(r, zorder=1, alpha=0.75, ax=ax)
cx.add_basemap(ax, alpha=0.5, crs=r.crs)