102 lines
3.9 KiB
Python
102 lines
3.9 KiB
Python
|
'''
|
|||
|
|
https://blog.csdn.net/vonuo/article/details/74783291
|
|||
|
|
本方法是使用asd的方法来进行corning的辐亮度定标（resonon的方法感觉和此方法差不多）
|
|||
|
|
'''
|
|||
|
|
from osgeo import gdal
|
|||
|
|
import os, math
|
|||
|
|
import sys
|
|||
|
|
import easygui
|
|||
|
|
import numpy as np
|
|||
|
|
import pandas as pd
|
|||
|
|
import xlwt
|
|||
|
|
|
|||
|
|
# 读写影像类
|
|||
|
|
class Grid(object):
|
|||
|
|
|
|||
|
|
#读图像文件
|
|||
|
|
@classmethod
|
|||
|
|
def read_img(cls, filename, xoff=0, yoff=0, im_width=None, im_height=None):
|
|||
|
|
try:
|
|||
|
|
dataset = gdal.Open(filename) # 打开文件
|
|||
|
|
if im_width == None:
|
|||
|
|
im_width = dataset.RasterXSize # 栅格矩阵的列数
|
|||
|
|
if im_height == None:
|
|||
|
|
im_height = dataset.RasterYSize # 栅格矩阵的行数
|
|||
|
|
num_bands = dataset.RasterCount # 栅格矩阵的波段数
|
|||
|
|
im_geotrans = dataset.GetGeoTransform() # 仿射矩阵
|
|||
|
|
im_proj = dataset.GetProjection() # 地图投影信息
|
|||
|
|
im_data = dataset.ReadAsArray(xoff, yoff, im_width, im_height) # 将数据写成数组，对应栅格矩阵
|
|||
|
|
del dataset
|
|||
|
|
return im_proj, im_geotrans, im_data
|
|||
|
|
except:
|
|||
|
|
sys.exit()
|
|||
|
|
|
|||
|
|
#写文件，以写成tif为例
|
|||
|
|
@classmethod
|
|||
|
|
def write_img(cls, dst_filename, data):
|
|||
|
|
format = "ENVI"
|
|||
|
|
driver = gdal.GetDriverByName(format)
|
|||
|
|
RasterXSize = data.shape[2] # 遥感影像的sample（列数）
|
|||
|
|
RasterYSize = data.shape[1] # 遥感影像的line（行数）
|
|||
|
|
band = data.shape[0]
|
|||
|
|
# driver.Create()函数中RasterXSize代表影像的sample（列数），RasterYSize代表影像的line（行数）
|
|||
|
|
dst_ds = driver.Create(dst_filename, RasterXSize, RasterYSize, band, gdal.GDT_Float32)
|
|||
|
|
for i in range(band):
|
|||
|
|
dst_ds.GetRasterBand(i + 1).WriteArray(data[i, :, :]) # gdal的band从1开始，所以dst_ds.GetRasterBand(i+1)
|
|||
|
|
dst_ds = None
|
|||
|
|
|
|||
|
|
# 是否转辐亮度：0→不转，1→转
|
|||
|
|
rad_switch = 0
|
|||
|
|
|
|||
|
|
img = r'D:\py_program\corning410\record_system_v24\baiban_record'
|
|||
|
|
img_baiban = r'D:\py_program\corning410\record_system_v24\baiban'
|
|||
|
|
img_gain = r'D:\py_program\corning410\corning410_radiance_calibration\jfq_dn_gain'
|
|||
|
|
img_offset = r'D:\py_program\corning410\record_system_v24\dark'
|
|||
|
|
|
|||
|
|
|
|||
|
|
dirpath = os.path.splitext(img)
|
|||
|
|
|
|||
|
|
# 读取影像
|
|||
|
|
proj, geotrans, data = Grid.read_img(img)
|
|||
|
|
proj_baiban, geotrans_baiban, data_baiban = Grid.read_img(img_baiban)
|
|||
|
|
proj_gain, geotrans_gain, data_gain = Grid.read_img(img_gain)
|
|||
|
|
proj_offset, geotrans_offset, data_offset = Grid.read_img(img_offset)
|
|||
|
|
|
|||
|
|
data_baiban = np.mean(data_baiban, axis=1)
|
|||
|
|
data_offset = np.mean(data_offset, axis=1)
|
|||
|
|
|
|||
|
|
if rad_switch == 1:
|
|||
|
|
# 计算辐射定标参数
|
|||
|
|
cal_it = 6059
|
|||
|
|
target_it = 200004
|
|||
|
|
gain_scale = cal_it / target_it
|
|||
|
|
data_gain_adjust = data_gain * gain_scale
|
|||
|
|
|
|||
|
|
# 影像和白板：1、扣除暗电流；2、转换成辐亮度；
|
|||
|
|
data_baiban = data_baiban - data_offset # 白板扣除暗电流
|
|||
|
|
data_baiban_rad = data_baiban * data_gain_adjust[:, 0, :] # 白板转辐亮度
|
|||
|
|
|
|||
|
|
data_rad = np.empty(data.shape)
|
|||
|
|
for i in range(data.shape[1]):
|
|||
|
|
data_rad[:, i, :] = (data[:, i, :] - data_offset) * data_gain_adjust[:, 0, :] # 转辐亮度
|
|||
|
|
Grid.write_img(dirpath[0] + '_rad', data_rad)
|
|||
|
|
|
|||
|
|
# 转换成反射率
|
|||
|
|
data_rad_ref = np.empty(data.shape)
|
|||
|
|
for i in range(data.shape[1]):
|
|||
|
|
data_rad_ref[:, i, :] = data_rad[:, i, :] / data_baiban_rad # 转反射率
|
|||
|
|
|
|||
|
|
Grid.write_img(dirpath[0] + '_rad_ref', data_rad_ref)
|
|||
|
|
elif rad_switch == 0:
|
|||
|
|
# 影像和白板扣除暗电流
|
|||
|
|
data_baiban = data_baiban - data_offset
|
|||
|
|
data_rmdark = np.empty(data.shape)
|
|||
|
|
for i in range(data.shape[1]):
|
|||
|
|
data_rmdark[:, i, :] = (data[:, i, :] - data_offset)
|
|||
|
|
|
|||
|
|
# 转反射率
|
|||
|
|
data_ref = np.empty(data.shape)
|
|||
|
|
for i in range(data.shape[1]):
|
|||
|
|
data_ref[:, i, :] = data_rmdark[:, i, :] / data_baiban # 转反射率
|
|||
|
|
Grid.write_img(dirpath[0] + '_ref', data_ref)
|