前言

本文主要是處理葵花8（Himawari 8） L1級netCDF格式數據
葵花8（Himawari 8）衛星數據請自行百度教程，官網數據下載地址：https://www.eorc.jaxa.jp/ptree/registration_top.html

數據簡介：https://www.eorc.jaxa.jp/ptree/userguide.html

所需要的庫：

netCDF4
gdal

代碼

import numpy as np
import netCDF4 as nc
from osgeo import gdal, osr, ogr
import os
data = r"...\dataset\H8\NC_H08_20160902_1050_R21_FLDK.06001_06001.nc"
nc_data = nc.Dataset(data) # 讀入數據
print(nc_data)

輸出結果如下，可以看到頭文件信息，包括行列號、數據時間、范圍和網格大小（空間分辨率）

<class 'netCDF4._netCDF4.Dataset'>
root group (NETCDF4 data model, file format HDF5):
title: Himawari-8 AHI equal latitude-longitude map data
id: NC_H08_20160902_1050_R21_FLDK.06001_06001.nc
date_created: 2016-09-02T11:05:56Z
pixel_number: 6001
line_number: 6001
upper_left_latitude: 60.0
upper_left_longitude: 80.0
grid_interval: 0.02
band_number: 16
dimensions(sizes): latitude(6001), longitude(6001), band(16), time(1), geometry(17)
variables(dimensions): float32 latitude(latitude), float32 longitude(longitude), int32 band_id(band), float64 start_time(time), float64 end_time(time), float64 geometry_parameters(geometry), int16 albedo_01(latitude, longitude), int16 albedo_02(latitude, longitude), int16 albedo_03(latitude, longitude), int16 sd_albedo_03(latitude, longitude), int16 albedo_04(latitude, longitude), int16 albedo_05(latitude, longitude), int16 albedo_06(latitude, longitude), int16 tbb_07(latitude, longitude), int16 tbb_08(latitude, longitude), int16 tbb_09(latitude, longitude), int16 tbb_10(latitude, longitude), int16 tbb_11(latitude, longitude), int16 tbb_12(latitude, longitude), int16 tbb_13(latitude, longitude), int16 tbb_14(latitude, longitude), int16 tbb_15(latitude, longitude), int16 tbb_16(latitude, longitude), int16 SAZ(latitude, longitude), int16 SAA(latitude, longitude), int16 SOZ(latitude, longitude), int16 SOA(latitude, longitude), int16 Hour(latitude, longitude)
groups:

list(nc_data.variables.keys())

這樣可以清晰的看到數據集中包含哪些數據，可以把nc數據看成是一個大的字典，可以按keys來提取數據。

['latitude',
'longitude',
'band_id',
'start_time',
'end_time',
'geometry_parameters',
'albedo_01',
'albedo_02',
'albedo_03',
'sd_albedo_03',
'albedo_04',
'albedo_05',
'albedo_06',
'tbb_07',
'tbb_08',
'tbb_09',
'tbb_10',
'tbb_11',
'tbb_12',
'tbb_13',
'tbb_14',
'tbb_15',
'tbb_16',
'SAZ',
'SAA',
'SOZ',
'SOA',
'Hour']

# 查看溫度波段的頭文件信息，可以看到名稱、比例因子、偏移量、單位等信息。
nc_data['tbb_13']

<class 'netCDF4._netCDF4.Variable'>
int16 tbb_13(latitude, longitude)
long_name: Brightness temperature of band 13
units: K
valid_min: -27315
valid_max: 32767
scale_factor: 0.01
add_offset: 273.15
missing_value: -32768
unlimited dimensions:
current shape = (6001, 6001)
filling on, default _FillValue of -32767 used

提取12波段的亮溫數據

tbb_12 = np.asarray(nc_data['tbb_12'][:]) # 結果為物理量

nc原數據是一個整數，這裡netCDF4直接進行了計算，將其轉為了亮溫的開爾文物理量。

lat = nc_data['latitude'][:]
lon = nc_data['longitude'][:]
latMin, latMax, lonMin, lonMax = lat.min(), lat.max(), lon.min(), lon.max()
# 分辨率
lat_Res = (latMax - latMin) / (lat.shape[0]-1)
lon_Res = (lonMax - lonMin) / (lon.shape[0]-1)
# 保存路徑和名稱，先占個坑
driver = gdal.GetDriverByName('GTiff')
out_tif_name = r'...\dataset\H8\H8_20160902_tbb_12_temp.tif'
out_tif = driver.Create(out_tif_name, lon.shape[0], lat.shape[0], 1, gdal.GDT_Float32)
# 輸出結果
geotransform = (lonMin,lon_Res, 0, latMax, 0, -lat_Res)
out_tif.SetGeoTransform(geotransform)
srs = osr.SpatialReference()
srs.SetWellKnownGeogCS('WGS84')
out_tif.SetProjection(srs.ExportToWkt())
out_tif.GetRasterBand(1).WriteArray(tbb_12)
out_tif.FlushCache()
out_tif = None