nc转tif文件，并计算单个栅格的平均值

1、将nc文件转为tif文件

若是一个nc文件中有多年数据，则用MATLAB：

(8条消息) NC文件读取及批量转为TIFF-史上最详细讲解-含代码（ArcGIS/MATLAB）_nc转tif_xiatiandexia123的博客-CSDN博客

一个文件夹里多个nc文件（nc.m）

注：首先用 arcgis-多维工具-将nc文件输出为一个栅格图层 作为参考

clc
clear
ncdisp('F:\DATA\IRR\IWU.ens.Y2011.nc')
sfc=ncread('F:\DATA\IRR\IWU.ens.Y2011.nc','IWUens');
%查看影像的缩略图
% imshow(sfc(:,:,1));
%如果数据需要旋转，使用rot函数进行相应修改
sfc=rot90(sfc);
%加载任意一个tif文件，tif文件需要与转换出来的tif空间边界一致（很重要）
[lizi,r]=geotiffread('F:\DATA\IRR\tif\201101.tif');
%设置输出文件夹
outpath='F:\DATA\IRR\tif\nctif\';
for year = 2011:2017
%逐一读取nc数据
sfc=ncread(['F:\DATA\IRR\IWU.ens.Y',int2str(year),'.nc'],'IWUens');
%循环单个nc文件中的文件数
for t=1:12
%获取第t层数据
sfc0=sfc(:,:,t);
%设置tif文件名称
filename=['IWUens',int2str(year),'_',int2str(t)];
%存储文件
geotiffwrite([outpath,filename,'.tif'],sfc0,r);

end

end

2、tif批量进行裁剪

利用arcgis中模型构建器进行批量裁剪

Arcgis批量掩膜、栅格计算—模型构建器_哔哩哔哩_bilibili

3、计算栅格平均值

arcpy进行计算

(8条消息) ArcPy批量计算栅格数据平均值_arcgis批量计算栅格平均值_bestYou_的博客-CSDN博客ArcPy批量计算栅格数据平均值、最大和最小值_哔哩哔哩_bilibili

# -*- coding: UTF-8 -*-
import os
import glob
import arcpy
from arcpy.sa import *

'''
功能：
    计算输入文件夹inws内，所有栅格数据的平均值
'''

arcpy.CheckOutExtension("ImageAnalyst")  # 检查许可
arcpy.CheckOutExtension("spatial")

# 输入路径  应该注意，中文路径，会导致读不出文件；路径尽量不要有空格，写文件时会报错
inws = r"F:\DATA\IRR\tif\2clip"

OutputFile = open('IRRaverage1.csv', 'w')

# 利用glob包，将inws下的所有tif文件读存放到rasters中
rasters = glob.glob(os.path.join(inws, "*.tif"))

whereClause = "VALUE = 0"  # 去除异常值

# 循环rasters中的所有影像，进行“求平均值”操作
for ras in rasters:
    outSetNull = SetNull(ras, ras, whereClause)  # 去除异常值

    meanValueInfo = arcpy.GetRasterProperties_management(outSetNull, 'MEAN')
    meanValue = meanValueInfo.getOutput(0)
    maxValueInfo = arcpy.GetRasterProperties_management(outSetNull, 'MAXIMUM')
    maxValue = maxValueInfo.getOutput(0)
    minValueInfo = arcpy.GetRasterProperties_management(outSetNull, 'MINIMUM')
    minValue = minValueInfo.getOutput(0)
    print os.path.basename(ras).split('_')[0] + ', mean:' + str(meanValue)+' max: '+str(maxValue)+' min: '+str(minValue) + '\n'

    OutputFile.write(os.path.basename(ras).split('_')[0] + ',' + str(meanValue)+ ',' + str(maxValue)+',' + str(minValue) + '\n')

OutputFile.close()
print("All project is OK！")

arcpy基本可行，但是有时候会遇见栅格统计不可计算的问题，则可用Spyder（python）

Spyder进行计算

(8条消息) 【Python】基于Python统计文件目录内每景栅格图像均值并输出成表_栅格统计值成表_陈修一的博客-CSDN博客

# -*- coding: utf-8 -*-
"""
Created on Wed Jun 28 10:52:55 2023

@author: Administrator
"""

from osgeo import gdal
import glob
import numpy as np
import csv

# 计算栅格图像均值
def calc_grid_mean(tiffile, csvpath):
    open_tif = gdal.Open(tiffile) # 打开tif文件
    # 打开栅格图像的波段
    band = open_tif.GetRasterBand(1).ReadAsArray()
    count_value = [] # 存放像元值的数组
    # 遍历图像矩阵，剔除异常值，读取有效值
    for i in range(len(band)):
        for j in range(len(band[i])):
            if band[i][j] >= -1000:
                count_value.append(band[i][j]) # 将有效值放入新数组
    # 把文件名和参数写入csv
    mycsv = open(str(csvpath), 'a', newline='') # 新建并打开csv
    csv_write = csv.writer(mycsv, dialect='excel')
    tif_info = (str(tif[-18:]), float(np.mean(count_value))) # 把文件名称+均值整合为一行
    csv_write.writerow(tif_info) # 把文件名称+均值写入csv
    print(f'{str(tif[-18:])}, {np.mean(count_value)}') # 在终端中打印出文件名称+均值


csv_path = r'F:\DATA\\IRR\tif\2clip\clipfile\csv.csv' # csv输出路径
tif_dir = glob.glob(r'F:\DATA\IRR\tif\2clip\clipfile\*.tif') # tif目录

for tif in tif_dir: # 遍历文件夹
    calc_grid_mean(tif, csv_path) # 计算均值