GEE（python）将图像转为numpy，并通过gdal重新输出为图像

为了实现这个功能，真的找了很多很多的资料。这里感谢吴秋生老师在github和YouTube上的分享。
接下来分享我的代码（ps：使用的编译器是jupyter lab，使用的数据是上一篇博文GEE（Python）逐像元线性拟合的结果）：

1.将image类型的图像转为numpy

这里用到了geemap的ee_to_numpy的函数

#Edited by Xinglu Cheng 2021.12.8
import ee
import os
from ipygee import Map 
from IPython.display import Image
import numpy as np
import matplotlib.pyplot as plt
from osgeo import gdal

#初始化gee，并设置底图
ee.Initialize()
import geemap
import geemap.colormaps as cm
gee_Map = geemap.Map(center=[32.5, 113.5], zoom=6)
gee_Map.add_basemap('Gaode.Normal')
#这里可以把gee_Map先显示一下，方便后面画矩形

#首先在地图上绘制研究区矩形
#将画的矩形转为feature
feature=gee_Map.draw_last_feature
roi=feature.geometry()

#将图像转为numpy，并输出大小
linearScale=linear.select('scale')#这里的linear是一个有两个波段的图像，选择其中一个波段
#print(linearScale.projection().nominalScale())#获取原图像像元大小
linear_repro=linearScale.reproject(crs='EPSG:4326',scale=4638)
im_data = geemap.ee_to_numpy(linear_repro,region=roi,default_value=-999).squeeze()
print(im_data.shape)

#显示numpy的图像
im_data=np.ma.masked_equal(im_data,-999)#掩膜掉-999的无效值
plt.imshow(im_data,cmap='gray')#通过plt显示灰度值
plt.show()

贴上效果图（这里输出了矩阵的大小，并用plt进行矩阵的显示）

2.通过gdal重新输出为tif图像

（有关gdal的tif输入输出函数，可以参考我之前的博文python（GDAL）读取、输出tif数据）

#输出图像
def write_tif(newpath,im_data,im_Geotrans,im_proj,datatype):
    if len(im_data.shape)==3:
        im_bands, im_height, im_width = im_data.shape
    else:
        im_bands, (im_height, im_width) = 1, im_data.shape
    diver = gdal.GetDriverByName('GTiff')
    new_dataset = diver.Create(newpath, im_width, im_height, im_bands, datatype)
    new_dataset.SetGeoTransform(im_Geotrans)
    new_dataset.SetProjection(im_proj)

    if im_bands == 1:
        new_dataset.GetRasterBand(1).WriteArray(im_data)
    else:
        for i in range(im_bands):
            new_dataset.GetRasterBand(i+1).WriteArray(im_data[i])
    del new_dataset

#获取影像的经纬度
def getPosition(image):
  image = image.addBands(ee.Image.pixelLonLat())
  data = image.reduceRegion(
      reducer = ee.Reducer.toList(),
      geometry = study_area,
      scale = 4638,#这里的像元大小和区域需要自行修改
      maxPixels = 1e13
  )
  lat = np.array(ee.Array(data.get("latitude")).getInfo())
  lon = np.array(ee.Array(data.get("longitude")).getInfo())
  return lat, lon
lat,lon=getPosition(linearScale)#调用获取经纬度的函数
lat_min=lat.max()
lon_min=lon.min()#获得图像左上角的坐标
print(lat_min,lon_min)

#设置投影和仿射变换
im_Geotrans=(lon_min, 0.0416638628774634, 0.0, lat_min, 0.0, -0.0416638628774634)
im_proj='GEOGCS["WGS 84",DATUM["WGS_1984",SPHEROID["WGS 84",6378137,298.257223563,AUTHORITY["EPSG","7030"]],AUTHORITY["EPSG","6326"]],PRIMEM["Greenwich",0,AUTHORITY["EPSG","8901"]],UNIT["degree",0.0174532925199433,AUTHORITY["EPSG","9122"]],AXIS["Latitude",NORTH],AXIS["Longitude",EAST],AUTHORITY["EPSG","4326"]]'
#调用输出函数
newpath=os.path.expanduser('~/gee_data')#这里是指输出到本地用户文件夹下的gee_data文件夹
newpath = os.path.join(newpath, 'try.tif')#输出名为try
write_tif(newpath,im_data,im_Geotrans,im_proj,gdal.GDT_Float32)#调用输出函数

接下来贴上结果：
这个是输出的图像左上角经纬度坐标

图像输出了到用户文件夹下

接下来可以通过将图像输出到谷歌云端，再将图像输入到gee中；或者使用吴秋生老师，geemap的add_raster函数将本地栅格直接加入到gee中（后者我本人实践后，发现图层并未显示，后面还会多加尝试）
以上！OWO