当前位置：首页 > news >正文

SMOS数据处理，投影变换，‘EPSG:6933‘转为‘EPSG:4326‘

news 2025/8/9 19:37:14

在处理SMOS数据时，遇到了读取nc数据并存为tif后，影像投影无法改变，因此全球数据无法重叠。源数据的投影为'EPSG:6933'，希望转为'EPSG:4326'。

解决代码。

```python
import os
import netCDF4 as nc
import numpy as np
from osgeo import gdal, gdalconst
# *************************导入netCDF4
file_path = r'I:\\'
input_name = file_path + 'SM_SCIE_MIR_L4AGBB_20110101T000000_20211231T235959_100_001_8.nc'
out_path = r'I:\out\\'
res = 0.25 #目标数据的分辨率
top = 57 #目标数据的top坐标
geo_transform = [-179, res, 0.00, top, 0.00, -res] #投影转换os.makedirs(out_path, exist_ok=True)with nc.Dataset(input_name) as DS:data = DS.variables['AGB'] #nc文件数据data_arr1 = np.asarray(data) #转为数组src_proj = DS.proj4text #nc原始文件投影信息for j in range(0, 11):year = str(2011 + j)data_arr3 = np.where(data_arr1[j] < 0, 0, data_arr1[j])data_arr4=np.flipud(data_arr3) #影像上下变换# 获取数组的行数和列数rows, cols = data_arr4.shape# 创建一个新的Dataset对象driver = gdal.GetDriverByName('MEM')dataset = driver.Create('', cols, rows, 1, gdal.GDT_Float32)# 将数组的数据写入Dataset对象band = dataset.GetRasterBand(1)band.WriteArray(data_arr4)# 设置地理变换信息dataset.SetGeoTransform(geo_transform)  # 这里的参数需要设置一下，尽管这个参数在warp后被抹去。# 设置投影信息dataset.SetProjection(src_proj)  # 源文件的投影。option = gdal.WarpOptions(format='GTiff', srcSRS=src_proj, dstSRS='EPSG:4326', outputType=gdal.GDT_Float32, resampleAlg=gdalconst.GRA_NearestNeighbour)out_name = out_path + '**' + year + '.tif'result_data=gdal.Warp(out_name, dataset, options=option,  xRes=res, yRes=res, creationOptions=["COMPRESS=LZW"], dstNodata=0)result_data.SetGeoTransform(geo_transform) #warp后geo_transform信息被抹去，这里需要再次设定。result_data.FlushCache()del result_dataprint(out_name)
'''