2

GDAL と Python を使用して ENVI ファイルを配列として読み込もうとしています

画像情報は次のとおりです。

Driver: ENVI/ENVI .hdr Labelled
Files: IMG-VV-ALPSRP248213250-P1.1__D_pwr_geo_sigma
     IMG-VV-ALPSRP248213250-P1.1__D_pwr_geo_sigma.hdr
Size is 1659, 2775
Coordinate System is:
PROJCS["UTM Zone 16, Northern Hemisphere",
    GEOGCS["WGS 84",
        DATUM["WGS_1984",
        SPHEROID["WGS 84",6378137,298.257223563,
            AUTHORITY["EPSG","7030"]],
        TOWGS84[0,0,0,0,0,0,0],
        AUTHORITY["EPSG","6326"]],
    PRIMEM["Greenwich",0,
        AUTHORITY["EPSG","8901"]],
    UNIT["degree",0.0174532925199433,
        AUTHORITY["EPSG","9108"]],
    AUTHORITY["EPSG","4326"]],
PROJECTION["Transverse_Mercator"],
PARAMETER["latitude_of_origin",0],
PARAMETER["central_meridian",-87],
PARAMETER["scale_factor",0.9996],
PARAMETER["false_easting",500000],
PARAMETER["false_northing",0],
UNIT["Meter",1]]
Origin = (332125.000000000000000,2017650.000000000000000)
Pixel Size = (25.000000000000000,-25.000000000000000)
Image Structure Metadata:
  INTERLEAVE=PIXEL
Corner Coordinates:
Upper Left  (  332125.000, 2017650.000) ( 88d35'16.06"W, 18d14'29.96"N)
Lower Left  (  332125.000, 1948275.000) ( 88d34'55.98"W, 17d36'53.41"N)
Upper Right (  373600.000, 2017650.000) ( 88d11'44.12"W, 18d14'40.22"N)
Lower Right (  373600.000, 1948275.000) ( 88d11'29.00"W, 17d37' 3.30"N)
Center      (  352862.500, 1982962.500) ( 88d23'21.23"W, 17d55'47.10"N)
Band 1 Block=1659x1 Type=Float32, ColorInterp=Undefined

私のコードは次のとおりです。

driver = gdal.GetDriverByName('ENVI')
driver.Register()

#Mind the suffix (It is an ENVI file)    
file = 'C:/img.1__A_pwr_geo_sigma

raster = gdal.Open(file,gdal.GA_ReadOnly)

raster_array = raster.ReadAsArray()

print raster_array

出力：

>>>array([[ nan, nan, nan, ..., nan, nan, nan], [ nan, nan, nan, ..., nan, nan, nan], [ nan, nan, nan, ..., nan, nan, nan], ..., [ nan, nan, nan, ..., nan, nan, nan], [ nan, nan, nan, ..., nan, nan, nan], [ nan, nan, nan, ..., nan, nan, nan]], dtype=float32)

イメージ内に浮動小数点 32 ビット値があることはわかっていますが、配列内のすべての値は NaN です (ENVI ソフトウェアで確認)。

ここで何が間違っていますか？それともサフィックスに問題がありますか？

また、gdal_translate を使用して ENVI 形式を Geotiff に変換しようとしましたが、GeoTiff も同じ配列を生成します。

4

1 に答える 1

2

朗報です。すべての値が nan というわけではありません。GDAL ドライバー (v.1.9.1) はうまく機能し、有効なデータは GDAL に依存するソフトウェア (QGIS など) で表示されます。

値はnanNoData を表すために使用されます。(通常、ほとんどの人は、NoData に有限の「ダミー」値、たとえば 1e+20 を使用します)。

import numpy as np
from osgeo import gdal
fname = r'C:\path\to\envi_data\IMG-HH-ALPSRP136260330-H1.1__A_pwr_geo_sigma'
ds = gdal.Open(fname, gdal.GA_ReadOnly)
band = ds.GetRasterBand(1)
print band.GetNoDataValue()
# None ## normally this would have a finite value, e.g. 1e+20
ar = band.ReadAsArray()
print np.isnan(ar).all()
# False
print '%.1f%% masked' % (np.isnan(ar).sum() * 100.0 / ar.size)
# 43.0% masked

Python/NumPy で欠損値のある配列を表す最良の方法は、マスクされた配列を使用することです。

mar = np.ma.masked_array(ar, np.isnan(ar))
print mar.min(), np.median(mar), mar.mean(), mar.std(), mar.max()
# 0.000715672 0.148093774915 0.0921848740388 0.0700106260235 5.0
于 2013-02-02T23:02:08.703 に答える