SRTM GEBCO mozaik

From GMT Türkiye Wiki
(Sürümler arası farklar)
Jump to: navigation, search
 
(3 kullanıcı tarafından yapılan 40 ara revizyon gösterilmiyor)
1. satır: 1. satır:
 +
[[Resim:marmozaik.jpg]]
 +
 
<bash>
 
<bash>
########################################################################
 
## SRTM veri indirme, GMT grid'e dönüştürme, birleştirme ve plot etme ##
 
########################################################################
 
 
# SRTM verileri aşağıdaki ftp sitesinden indirilebilir
 
 
ftp 152.61.4.9
 
# kullanıcı anonymous
 
# şifre email adresi
 
 
# /srtm/Documentation klasörünün altında Continent_def.gif isimli dosyayı indirip ilgilenen bölge
 
# hangi kategoriye giriyor öğren.
 
 
cd  /srtm/Documentation
 
get Continent_def.gif
 
 
# Farzedelim Eurasia katagorisi olsun. Bu durumda 90 m'lik yani 3 ark saniyelik SRTM verinizin
 
# bulunduğu yer şudur
 
 
cd /srt/version2/SRTM3/Eurasia
 
 
# burada 1x1 derecelik çerçevelere bölünmüş durumda  srtm verileri. Herbir çerçevenin nereyi
 
# kapsadığı isminden anlaşılabilir. Çerçeve isminlerinde o çerçevenin sol alt köşe kordinatı
 
# verilmektedir. Örneğin 40-41 kuzey enlemleri ve 29-30 boylamları arasında kalan bölgenin
 
# veririni içeren dosyanın adı N40E029.hgt.zip'dir. Bu çeçeve Bursa-Gemlik civarını kapsamaktadır.
 
 
# Şimdi bunu ve yandaki çerçeveyi indirelim.
 
 
get N40E029.hgt.zip
 
get N40E028.hgt.zip
 
 
# zipleri açalım
 
 
unzip N40E029.hgt.zip
 
unzip N40E028.hgt.zip
 
 
 
# SRTM verileri 2 byte binary formatındadır. Yani her bir yüksekik verisi için ayrılan yer 2 byte'dır.
 
# Bu durumda her bir çerçevenin ne kadarlık bir yer kapladığını hesap edebiliriz.
 
 
# 1 derece = 3600 saniye'dir
 
 
# indirdiğimiz SRTM 3 arc saniye ( 1 arc saniye verisi sadece kuzey amerika için var sadece)
 
 
# Her bir piksel 3 saniyelik bir uzunluğa sahip olduğuna göre bir derecelik mesafe için gerekli
 
# olan piksel sayısı 3600/3 = 1200
 
# Dolayısıyla her bir çerçevede 1x1 derecelik bir alan olduğuna göre piksel sayımız 1200x1200
 
# olması gerekiyor.  Ancak grid değilde piksel kaydı yapıldığ için birer piksel fazla oluyor
 
# yani 1201x1201. 
 
# peki niye SRTM 90 m diyoruz???
 
# Çünkü 1 derece = 111 km  ve 1200 piksel var olduğuna göre 111000/1200 = ~ 90 m
 
 
 
# Herbir piksel 2 byte lik bir yer tutuğuna göre çerçevenin boyutu 1201x1201x2=2884802 byte
 
# olması gerekir. Bunu doğrulamak için
 
 
ls -la
 
 
# yapınca aşağıdaki gibi bir çıktı görürüsünüz
 
 
# -rw-r--r--  1 ziyadin  staff  2884802 May 27  2005 N40E029.hgt
 
 
# ve buradanda boyutun gerçektende  2884802 olduğu anlaşılır. Aksi türlü bir yanlışlık var demektir.
 
 
# Şimdi bu dosyaları xyz2grd komutu ile gridleyip birleştirelim.
 
 
# -R ile dem'in sınırları verilir: -R/lonMin/lonMax/latMin/latMax
 
# -I ile derece cinsinden piksel büyüklüğü: 1/1200 = 0.000833333
 
# -Zh ile verinin 2 byte binary olduğu belirtilir
 
# -G ile output grd dosya ismi verilir
 
 
xyz2grd N40E029.hgt -R29/30/40/41 -I0.00083333 -Zh  -GN40E029.grd
 
xyz2grd N40E028.hgt -R28/29/40/41 -I0.00083333 -Zh  -GN40E028.grd
 
 
# Şimdi grdinfo programı  ile yükseklik değerlerine bir bakalım mantıklımı
 
 
grdinfo N40E029.grd
 
 
#  deyince şu satırlar gözükür
 
#  N40E029.grd: Title: N40E029.grd
 
#  N40E029.grd: Command: xyz2grd -R29/30/40/41 -I0.00083333 N40E029.hgt -Zh -GN40E029.grd
 
#  N40E029.grd: Remark:
 
#  N40E029.grd: Gridline node registration used
 
#  N40E029.grd: Grid file format: nf (# 18)
 
#  N40E029.grd: x_min: 29 x_max: 30 x_inc: 0.000833333 name: x nx: 1201
 
#  N40E029.grd: y_min: 40 y_max: 41 y_inc: 0.000833333 name: y ny: 1201
 
#  N40E029.grd: z_min: -32768 z_max: 32521 name: z
 
#  N40E029.grd: scale_factor: 1 add_offset: 0
 
 
# Görüldğü gibi bende z_min ve z_max değerleri -32768 ve 32521
 
# Burada bir hata var!!!
 
 
# Sizdede böyle çıkıyorsa o zaman bunun düzeltilmesi gerekiyor.
 
# Bunun sebebi binary dosyada byte-order'ın  (big veya little endiean) yanlış olmasıdır.
 
# Swap etmek için -Z  opsiyonuna w eklenmesi gerekir.
 
 
xyz2grd N40E029.hgt -R29/30/40/41 -I0.00083333 -Zhw  -GN40E029.grd
 
xyz2grd N40E028.hgt -R28/29/40/41 -I0.00083333 -Zhw  -GN40E028.grd
 
 
# Yeniden bakalım düzelmişmi
 
 
grdinfo N40E029.grd
 
 
#dersek şunları görürüz
 
 
#  N40E029.grd: z_min: -32768 z_max: 2531 name: z
 
 
 
# x_max düzelmiş ama z_min halen kötü. Bu değeri sıfırlamamız gerekiyor.
 
# Sıfırlamak için grdclip komutunu kullanabiliriz.
 
 
# -Sb0/NaN opsiyonu ile 0 den küçük tüm değerler NaN (Not A Number) olarak atanır ve
 
# yeni grid -G opsiyonunda verilir. Deniz seviyesinin olduğ yerler var malum bu durumda
 
# 0 yerine -1000 felan verilebilir.
 
 
grdclip N40E029.grd -Sb0/NaN -GN40E029C.grd
 
grdclip N40E028.grd -Sb0/NaN -GN40E028C.grd
 
 
# Yeniden bakalım
 
 
grdinfo N40E029C.grd
 
 
#dersek
 
#
 
#N40E029C.grd: z_min: 0 z_max: 2531 name: z
 
#
 
 
grdinfo N40E028C.grd
 
#.
 
#N40E028C.grd: z_min: 0 z_max: 1251 name: z
 
#
 
 
# bunlar yanyana olduğu için yine GMT grdpaste ile birleştirilip -G opsiyonu ile birleştirilen
 
# dosyaya mozaiksrtm.grd ismi verilebilir.
 
 
#########################  birleştirme #################################
 
grdpaste N40E029C.grd N40E028C.grd -Gmozaiksrtm.grd
 
 
grdinfo mozaiksrtm.grd
 
 
#  mozaiksrtm.grd: Title: mozaiksrtm.grd
 
#  mozaiksrtm.grd: Command: grdpaste N40E029C.grd N40E028C.grd -Gmozaiksrtm.grd
 
#  mozaiksrtm.grd: Remark:
 
#  mozaiksrtm.grd: Gridline node registration used
 
#  mozaiksrtm.grd: Grid file format: nf (# 18)
 
#  mozaiksrtm.grd: x_min: 28 x_max: 30 x_inc: 0.000833333 name: x nx: 2401
 
#  mozaiksrtm.grd: y_min: 40 y_max: 41 y_inc: 0.000833333 name: y ny: 1201
 
#  mozaiksrtm.grd: z_min: 0 z_max: 2531 name: z
 
#  mozaiksrtm.grd: scale_factor: 1 add_offset: 0
 
 
# Görüldüğü gibi doğu batı yönünde piksel adeti ikiye katlandı.
 
 
 
 
#################### Şimdi bunu plot edip görelim #####################
 
# bir grd değişkeni tanımlayıp tektek her seferinde yazmaktan kurtulalım
 
 
set grd = mozaiksrtm.grd
 
 
# Renkli palet yapalım
 
grd2cpt -Crainbow -Z $grd > ! color.cpt
 
 
# Gölge dosyası oluşturalım
 
grdgradient  $grd -A315 -Nt -Gshade.grd
 
 
# bazı gmt ayarları
 
gmtset PLOT_DEGREE_FORMAT +D MEASURE_UNIT inch
 
 
# dosya ismi verelim
 
set name = marmara.ps
 
 
# alanı ve projeksiyonu (JM = markatör) setleyelim
 
set R = "-R28/30/40/41 -JM6"
 
 
psbasemap $R -B0.2 -K > ! $name
 
grdimage $R -Ccolor.cpt  $grd -Ishade.grd -O -K  >> $name
 
pscoast -Df -O $R -W -Ia >> $name
 
xv $name
 
 
 
#gri ton olsun ama yükseltiye göre değişmesin tonlama, gölgeye göre değişşin.
 
# Palet oluşturalım
 
echo -10000 180 180 180 10000 180 180 180 > ! gray.cpt
 
 
psbasemap $R -B0.2 -K -P> ! $name
 
grdimage $R -Cgray.cpt $grd -Ishade.grd -O -K>> $name
 
pscoast -Df -S120 -W -O -R -JM -Ia/2/blue -K >> $name
 
psbasemap -Lf29.8/40.1/29.8/40.1/20 -R -JM -O -B.2  >> $name
 
xv $name &
 
 
  
 
############################################
 
############################################
## GEBCO veri indirme srtm ile mozaikleme ##
+
## GEBCO veri indirme SRTM ile mozaikleme ##
 
############################################
 
############################################
  
# Web arayüyü kullanarak istenilen diktörtgen bir alan seçilip veya kordinatları verilip indirilebilir.  
+
# Web arayüzü kullanarak istenilen dikdörtgen bir alan seçilip veya koordinatları verilip indirilebilir.  
 
# Default olarak  indirilen dosya aslında hazır GMT grd'si ancak uzantısı nc olarak veriliyor her  
 
# Default olarak  indirilen dosya aslında hazır GMT grd'si ancak uzantısı nc olarak veriliyor her  
 
# nedense. Bu nedenle kafa karışmaması için soyadını grd'ye çevirmeniz yeterli.
 
# nedense. Bu nedenle kafa karışmaması için soyadını grd'ye çevirmeniz yeterli.
  
# Tüm grid de indirilebilir vey 1 GB'ın üstünde.  Bu dosya grdinfo yaparsak şunu görebiliriz
+
# 1 GB'ın üstünde büyüklüğü sahip tüm grid indirilebilir ve indirilen bu dosyaya grdinfo yaparsak
# (bendeki dosya halen nc soyadı taşıyor :(
+
# şunu görebiliriz (bendeki dosya halen nc soyadı taşıyor :(
  
 
grdinfo GEBCO_08.nc  
 
grdinfo GEBCO_08.nc  
 +
  
 
# GEBCO_08.nc: Title: GEBCO_08 Grid
 
# GEBCO_08.nc: Title: GEBCO_08 Grid
214. satır: 27. satır:
 
# GEBCO_08.nc: scale_factor: 1 add_offset: 0
 
# GEBCO_08.nc: scale_factor: 1 add_offset: 0
  
# Görüldüğü gibi bu grid tüm dünyayı kapsıyor.
+
# Görüldüğü üzere bu grid tüm dünyayı kapsıyor.
  
  
 
# Tüm grid indirilmişşe kesilmesi gerekiyor ve bunun için de grdcut komutu kullanılabilir.  
 
# Tüm grid indirilmişşe kesilmesi gerekiyor ve bunun için de grdcut komutu kullanılabilir.  
# Gebco'yu yukarıda mozaiklediğimiz SRTM gridinin boyutunda keselim aşağıdaki komutla
+
# Gebco'yu aşağıdaki komutla yukarıda mozaiklediğimiz SRTM gridinin boyutunda keselim
  
 
grdcut -R28/30/40/41 GEBCO_08.nc -GGebco_cut.grd -fg
 
grdcut -R28/30/40/41 GEBCO_08.nc -GGebco_cut.grd -fg
  
# Gebco verisinde piksel aralığı 30 saniye yani SRTM in 3 katı. Dolayısıyla bunu SRTM in  
+
# Gebco verisinde piksel aralığı 30 saniye yani SRTM'in 3 katı. Dolayısıyla bunu SRTM in  
# piksel boyutuna düşürelim grdsample komtu ile.
+
# piksel boyutuna düşürelim: grdsample komutu ile.
 
+
  
 
grdsample -I0.000833333 Gebco_cut.grd -Ggebco_cutS.grd -T
 
grdsample -I0.000833333 Gebco_cut.grd -Ggebco_cutS.grd -T
 +
</bash>
  
# grdinfo gebco_cutS.grd yaparsak aşağıda görüldüğü gibi piksel büyüklüğü (0.000833333),  
+
grdinfo gebco_cutS.grd yaparsak aşağıda görüldüğü gibi piksel büyüklüğü (0.000833333),  
# alanı (28/30/40/41) ve dolayısıyla piksel sayısı (nx ve ny) aynı olur SRTM mozaiksrtm.grd ile.
+
alanı (28/30/40/41) ve dolayısıyla piksel sayısı (nx ve ny) aynı olur SRTM mozaiksrtm.grd ile mozaiksrtm.grd nasıl elde ediliyor [[SRTM çizimi]] sayfasına bakınız.
  
 +
<bash>
 
# gebco_cutS.grd: Title: gebco_cutS.grd
 
# gebco_cutS.grd: Title: gebco_cutS.grd
 
# gebco_cutS.grd: Command: grdsample -I0.000833333 Gebco_cut.grd -Ggebco_cutS.grd -T
 
# gebco_cutS.grd: Command: grdsample -I0.000833333 Gebco_cut.grd -Ggebco_cutS.grd -T
254. satır: 68. satır:
 
# alanlarını NaN yapıp grdmath komutunu kullanabiliriz.
 
# alanlarını NaN yapıp grdmath komutunu kullanabiliriz.
  
# karaların NaN denizlerin 0 olduğu bir mask dosyası yapalım -N0/NaN ile
+
# karaların NaN, denizlerin 0 olduğu bir mask dosyası yapalım -N0/NaN ile
 
grdlandmask -Df gebco_cutS.grd -N0/NaN -Gmask.grd -R28/30/40/41 -I0.000833333
 
grdlandmask -Df gebco_cutS.grd -N0/NaN -Gmask.grd -R28/30/40/41 -I0.000833333
  
#bu mask uygulayalım gebco verisine
+
#bu maskgebco verisine uygulayalım
 
grdmath mask.grd gebco_cutS.grd ADD = gebmasked.grd
 
grdmath mask.grd gebco_cutS.grd ADD = gebmasked.grd
  
# srtm ile maskelenmiş gebco verisini mozaikleyelim
+
# SRTM ile maskelenmiş Gebco verisini mozaikleyelim
 
grdmath gebmasked.grd mozaiksrtm.grd AND = mozaik.grd
 
grdmath gebmasked.grd mozaiksrtm.grd AND = mozaik.grd
  
278. satır: 92. satır:
 
xv $name
 
xv $name
 
</bash>
 
</bash>
 +
--[[Kullanıcı:Ziyadin|Ziyadin]] 22:33, 9 Nisan 2009 (CEST)

22:33, 9 Nisan 2009 itibarı ile sayfanın şu anki hâli

Marmozaik.jpg

 
 
############################################
## GEBCO veri indirme SRTM ile mozaikleme ##
############################################
 
# Web arayüzü kullanarak istenilen dikdörtgen bir alan seçilip veya koordinatları verilip indirilebilir. 
# Default olarak  indirilen dosya aslında hazır GMT grd'si ancak uzantısı nc olarak veriliyor her 
# nedense. Bu nedenle kafa karışmaması için soyadını grd'ye çevirmeniz yeterli.
 
# 1 GB'ın üstünde büyüklüğü sahip tüm grid indirilebilir ve indirilen bu dosyaya grdinfo yaparsak
# şunu görebiliriz (bendeki dosya halen nc soyadı taşıyor :(
 
grdinfo GEBCO_08.nc 
 
 
# GEBCO_08.nc: Title: GEBCO_08 Grid
# GEBCO_08.nc: Command: 20090202
# GEBCO_08.nc: Remark: 
# GEBCO_08.nc: Pixel node registration used
# GEBCO_08.nc: Grid file format: cs (# 8)
# GEBCO_08.nc: x_min: -180 x_max: 180 x_inc: 0.00833333 name: user_x_unit nx: 43200
# GEBCO_08.nc: y_min: -90 y_max: 90 y_inc: 0.00833333 name: user_y_unit ny: 21600
# GEBCO_08.nc: z_min: -10977 z_max: 8685 name: user_z_unit
# GEBCO_08.nc: scale_factor: 1 add_offset: 0
 
# Görüldüğü üzere bu grid tüm dünyayı kapsıyor.
 
 
# Tüm grid indirilmişşe kesilmesi gerekiyor ve bunun için de grdcut komutu kullanılabilir. 
# Gebco'yu aşağıdaki komutla yukarıda mozaiklediğimiz SRTM gridinin boyutunda keselim  
 
grdcut -R28/30/40/41 GEBCO_08.nc -GGebco_cut.grd -fg
 
# Gebco verisinde piksel aralığı 30 saniye yani SRTM'in 3 katı. Dolayısıyla bunu SRTM in 
# piksel boyutuna düşürelim: grdsample komutu ile.
 
grdsample -I0.000833333 Gebco_cut.grd -Ggebco_cutS.grd -T
 

grdinfo gebco_cutS.grd yaparsak aşağıda görüldüğü gibi piksel büyüklüğü (0.000833333), alanı (28/30/40/41) ve dolayısıyla piksel sayısı (nx ve ny) aynı olur SRTM mozaiksrtm.grd ile mozaiksrtm.grd nasıl elde ediliyor SRTM çizimi sayfasına bakınız.

 
# gebco_cutS.grd: Title: gebco_cutS.grd
# gebco_cutS.grd: Command: grdsample -I0.000833333 Gebco_cut.grd -Ggebco_cutS.grd -T
# gebco_cutS.grd: Remark: 
# gebco_cutS.grd: Gridline node registration used
# gebco_cutS.grd: Grid file format: nf (# 18)
# gebco_cutS.grd: x_min: 28 x_max: 30 x_inc: 0.000833333 name: longitude [degrees_east] nx: 2401
# gebco_cutS.grd: y_min: 40 y_max: 41 y_inc: 0.000833333 name: latitude [degrees_north] ny: 1201
# gebco_cutS.grd: z_min: -1250.12 z_max: 2442.18 name: z
# gebco_cutS.grd: scale_factor: 1 add_offset: 0
 
# mozaiksrtm.grd: Title: mozaiksrtm.grd
# mozaiksrtm.grd: Command: grdpaste N40E029C.grd N40E028C.grd -Gmozaiksrtm.grd
# mozaiksrtm.grd: Remark: 
# mozaiksrtm.grd: Gridline node registration used
# mozaiksrtm.grd: Grid file format: nf (# 18)
# mozaiksrtm.grd: x_min: 28 x_max: 30 x_inc: 0.000833333 name: x nx: 2401
# mozaiksrtm.grd: y_min: 40 y_max: 41 y_inc: 0.000833333 name: y ny: 1201
# mozaiksrtm.grd: z_min: 0 z_max: 2531 name: z
# mozaiksrtm.grd: scale_factor: 1 add_offset: 0
 
# Şimdi iki veri seti mozaiklenebilir çünkü tıpa tıp aynı iki grd. Bunun için Gebco gridinin kara 
# alanlarını NaN yapıp grdmath komutunu kullanabiliriz.
 
# karaların NaN, denizlerin 0 olduğu bir mask dosyası yapalım -N0/NaN ile
grdlandmask -Df gebco_cutS.grd -N0/NaN -Gmask.grd -R28/30/40/41 -I0.000833333
 
#bu mask'ı gebco verisine uygulayalım 
grdmath mask.grd gebco_cutS.grd ADD = gebmasked.grd
 
# SRTM ile maskelenmiş Gebco verisini mozaikleyelim
grdmath gebmasked.grd mozaiksrtm.grd AND = mozaik.grd
 
#Geçmiş olsun mozaik hazırdır !!!!!
 
 
#Plot et
gmtset PLOT_DEGREE_FORMAT +D MEASURE_UNIT inch
set grd = mozaik.grd 
grd2cpt -Crainbow -Z $grd > ! color.cpt
grdgradient  $grd -A315 -Nt -Gshade.grd
set name = marmozaik.ps
set R = "-R28/30/40/41 -JM6"
psbasemap $R -B0.2 -K > ! $name
grdimage $R -Ccolor.cpt  $grd -Ishade.grd -O -K  >> $name
pscoast -Df -O $R -W -Ia >> $name 
xv $name
 

--Ziyadin 22:33, 9 Nisan 2009 (CEST)

Personal tools