Raster を複数使用した安定計算。
解像度の違いでエラーを吐きます。
Arc で ポリゴンからラスターへ変換していたのですが、解像度が完全には一致しないようでした。仕方がないので、リサンプリングすることに。
以前、GDAL を使って Raster を読み込んでいたので、そこに追記。
# Open Raster
def open_src(file,band):
src = gdal.Open(file, gdal.GA_ReadOnly)
if src is None:
print('file not found')
else:
print(src.GetDescription() )
xy=src.GetGeoTransform()#
print(xy)
print(src.GetProjection())
print(src.GetMetadata() )
arr_b1=src.GetRasterBand(band).ReadAsArray(0)
print(arr_b1.shape)
nodata=src.GetRasterBand(band).GetNoDataValue())
print(nodata)
return src,xy,arr_b1,nodata
#図化
def plot_raster(arr_b1, nodata):
arr_b1 = np.where(arr_b1==nodata, np.nan, arr_b1)
print('max:',np.nanmax(arr_b1))
print('min:',np.nanmin(arr_b1))
print(arr_b1.shape)
plt.imshow(arr_b1)
plt.colorbar()
return arr_b1
#測地系変換
#EPSG:4612 から 4326に変換して保存
ds = gdal.Warp(file_out, file_in, srcSRS="EPSG:4612", dstSRS="EPSG:4326")
ds=None
#2つのラスターの領域、解像度を(ほぼ)揃える
#切り出し領域を決める
minx=136
miny=35
maxx=137
maxy=36
#領域 をforium で確認
fmap = folium.Map(location=[(miny+maxy)/2,(minx+maxx)/2], zoom_start=8)
line = [[miny,minx],[miny,maxx],[maxy,maxx],[maxy,minx],[miny,minx]]
folium.PolyLine(locations=line).add_to(fmap)
fmap
#ピクセル数を決める
width=18000
height=11100
#resampling
cut_file0="./output/cut0.tif"
cut_file1="./output/cut1.tif"
ds=gdal.Translate(cut_file0, file_in0, projWin=[minx,maxy,maxx,miny], width=width, height=height, noData=nodata)
ds=gdal.Translate(cut_file1, file_in1, projWin=[minx,maxy,maxx,miny], width=width, height=height, noData=nodata)
ds=None
#保存
def save_dst(rasterfile, arr_new, nodata, src):
driver = gdal.GetDriverByName("GTiff")
dst = driver.Create(rasterfile,
xsize=arr_new.shape[1],
ysize=arr_new.shape[0],
bands=1,
eType = gdal.GDT_Float32)
dst.SetGeoTransform(src.GetGeoTransform())
dst.SetProjection(src.GetProjection())
band1=dst_ds.GetRasterBand(1)
band1.WriteArray(arr_new)
band1.SetNoDataValue(nodata)
dst.FlushCache()
Arc でラスターを2枚作成・書き出していましたので、マルチバンドにしてから1枚出力すれば解像度が完全に一致したのですが。
https://pro.arcgis.com/ja/pro-app/2.8/tool-reference/data-management/composite-bands.htm
ま、今後使えるかもしれませんので残しておきましょう。