GRASS GIS / Poznámky pro družicová data Landsat
Tato stránka obsahuje poznámky pro stažení dat LandSat Thematic Mapper a následné zpracování těchto dat v prostřední GRASS GIS.
- Časová řada
- 80. léta (co nejstarší)
- přelom 80. a 90. let
- kolem roku 2000
- nejnovější
Postup pro manuální stažení dat
- Data jsou dostupná na adrese: http://eros.usgs.gov (bližší informace zde)
- Data mohou být stažena také z adres: http://glovis.usgs.gov nebo http://edcsns17.cr.usgs.gov/EarthExplorer
Návod pro http://glovis.usgs.gov:
Import dat
Data je nejprve nutno dekomprimovat
tar xvzf LT51910252009213MOR00.tar.gz
Archiv obsahuje několik souborů ve formátu geotiff (jeden soubor na kanál), seznam identických bodů (GCP, ground control points) a metadata v textové podobě (MTL).
ls -w1 L5191025_02520090801_B10.TIF L5191025_02520090801_B20.TIF L5191025_02520090801_B30.TIF L5191025_02520090801_B40.TIF L5191025_02520090801_B50.TIF L5191025_02520090801_B60.TIF L5191025_02520090801_B70.TIF L5191025_02520090801_GCP.txt L5191025_02520090801_MTL.txt README.GTF
Nástroj knihovny GDAL gdalinfo poslouží pro kontrolu metadat, např.
gdalinfo L5191025_02520090801_B10.TIF
Driver: GTiff/GeoTIFF Files: L5191025_02520090801_B10.TIF Size is 8441, 7441 Coordinate System is: PROJCS["unnamed", GEOGCS["WGS 84", DATUM["WGS_1984", SPHEROID["WGS 84",6378137,298.257223563, AUTHORITY["EPSG","7030"]], AUTHORITY["EPSG","6326"]], PRIMEM["Greenwich",0], UNIT["degree",0.0174532925199433], AUTHORITY["EPSG","4326"]], PROJECTION["Transverse_Mercator"], PARAMETER["latitude_of_origin",0], PARAMETER["central_meridian",15], PARAMETER["scale_factor",0.9996], PARAMETER["false_easting",500000], PARAMETER["false_northing",0], UNIT["unknown",1], AUTHORITY["EPSG","32633"]] Origin = (398699.999999999941792,5679000.000000000000000) Pixel Size = (30.000000000000000,-30.000000000000000) Metadata: AREA_OR_POINT=Point Image Structure Metadata: INTERLEAVE=BAND Corner Coordinates: Upper Left ( 398700.000, 5679000.000) ( 13d32'54.25"E, 51d15'12.06"N) Lower Left ( 398700.000, 5455770.000) ( 13d36'29.29"E, 49d14'46.73"N) Upper Right ( 651930.000, 5679000.000) ( 17d10'35.91"E, 51d14'31.60"N) Lower Right ( 651930.000, 5455770.000) ( 17d 5'13.62"E, 49d14'9.03"N) Center ( 525315.000, 5567385.000) ( 15d21'18.49"E, 50d15'29.05"N) Band 1 Block=8441x1 Type=Byte, ColorInterp=Gray
Vytvoření lokace
Hromadný import dat
Data lze naimportovat pomocí modulu r.in.gdal (File -> Import raster data -> Import raster data), ukázka z příkazové řadky (přepínač -e nastaví po importu region podle vstupní dat).
r.in.gdal -e input=L5191025_02520090801_B10.TIF output=b10 title="kanal 1"
Základní metadata vypíše r.info.
r.info map=b10
+----------------------------------------------------------------------------+ | Layer: b10 Date: Thu Sep 16 17:48:23 2010 | | Mapset: PERMANENT Login of Creator: martin | | Location: tm | | DataBase: /home/martin/grassdata | | Title: kanal 1 ( b10 ) | | Timestamp: none | |----------------------------------------------------------------------------| | | | Type of Map: raster Number of Categories: 0 | | Data Type: CELL | | Rows: 7441 | | Columns: 8441 | | Total Cells: 62809481 | | Projection: UTM (zone 33) | | N: 5679000 S: 5455770 Res: 30 | | E: 651930 W: 398700 Res: 30 | | Range of data: min = 0 max = 255 | | | | Data Description: | | generated by r.in.gdal | | | | Comments: | | r.in.gdal -e input="L5191025_02520090801_B10.TIF" output="b10" | | | +----------------------------------------------------------------------------+
Data lze hromadně naimportovat pomocí jednoduchého skriptu. Následuje ukázka pro Bash a Python.
Bash
Úkázka pro Bash.
#!/bin/bash
idx=1
for file in *.TIF; do
oname=`echo $file | cut -d'_' -f3 | cut -d'.' -f1`
g.message "Importuji $file -> $oname..."
r.in.gdal -e input=$file output=$oname title="kanal $idx"
idx=$(($idx+1))
done
Python
Úkázka pro Python.
#!/usr/bin/python
import os
import glob
import grass.script as grass
idx = 1
for file in glob.glob1('.', '*.TIF'):
oname = os.path.splitext(file)[0].split('_')[-1]
grass.message('Importuji %s -> %s...' % (file, oname))
grass.run_command('r.in.gdal',
flags = 'e',
input = file,
output = oname,
title = "kanal %d" % idx)
idx += 1
wxGUI
Dále lze data naimportovat pomocí wxGUI.
File -> Import raster data -> Bulk import of raster data
Vizualizace dat
Příklad vizualizace kompozitního snímku (kanály 453).
Seznam měst
Brno Břeclav Česká Lípa Česká Třebová České Budějovice Český Těšín Děčín Frýdek-Místek Cheb Chomutov Chrudim Jablonec Jablonec nad Nisou Jihlava Jirkov Kadaň Karviná Klatovy Kutná Hora Liberec Litvínov Náchod Ostrava Plzeň Přerov Příbram Rožnov pod Radhoštěm Sokolov Teplice Třebíč Ústí nad Labem Ústí nad Orlicí Valašské Meziříčí Vsetín Zlín Znojmo Žďár nad Sázavou