311 lines
16 KiB
Python
311 lines
16 KiB
Python
"""
|
||
@Project ：microproduct
|
||
@File ：BackScatteringXmlInfo.PY
|
||
@Function ：主函数
|
||
@Author ：LMM
|
||
@Date ：2021/10/19 14:39
|
||
@Version ：1.0.0
|
||
"""
|
||
import os
|
||
from xml.etree.ElementTree import ElementTree, Element
|
||
import xml.dom.minidom
|
||
from lxml import etree
|
||
import shutil
|
||
from tool.algorithm.image.ImageHandle import ImageHandler
|
||
from tool.algorithm.algtools.PreProcess import PreProcess as pp
|
||
from osgeo import gdal
|
||
import numpy as np
|
||
import datetime
|
||
from PIL import Image
|
||
|
||
|
||
class CreateDict:
|
||
"""根据影像/DEM的属性信息添加到字典中"""
|
||
def __init__(self):
|
||
self.ImageHandler = ImageHandler()
|
||
pass
|
||
|
||
def calu_nature(self, image_path, image_pair, out_path1, out_path2):
|
||
"""
|
||
将productinfo节点需要填写的信息存入字典中
|
||
image_path:影像路径
|
||
image_pair:输入的压缩包中的极化对 例：hh,hv,vh,vv=【1，1，1，1】
|
||
out_path1：地理转平面的输出路径
|
||
out_path2：平面转地理的输出路径
|
||
"""
|
||
|
||
para_dict = {}
|
||
imageinfo_width = self.ImageHandler.get_img_width(image_path)
|
||
para_dict.update({"imageinfo_width":imageinfo_width})
|
||
imageinfo_height = self.ImageHandler.get_img_height(image_path)
|
||
para_dict.update({"imageinfo_height":imageinfo_height})
|
||
para_dict.update({"imageinfo_EarthModel": "WGS84"})
|
||
para_dict.update({"imageinfo_ProjectModel": "UTM"})
|
||
proj = self.ImageHandler.get_projection(image_path) # 输出的影像若是投影坐标系则先转成地理坐标系
|
||
keyword = proj.split("[", 2)[0] # 若是地理坐标系则pass
|
||
if keyword == "GEOGCS":
|
||
pass
|
||
elif keyword == "PROJCS":
|
||
pp.trans_projcs2geogcs(out_path2, image_path)
|
||
image_path = out_path2
|
||
elif len(keyword) == 0 or keyword.strip() == "" or keyword.isspace() is True:
|
||
raise Exception('image projection is missing!')
|
||
|
||
pp.trans_geogcs2projcs(out_path1, image_path) # 坐标投影, 地理转平面投影坐标
|
||
imageinfo_widthspace = self.ImageHandler.get_geotransform(out_path1)[1] # 投影后的分辨率
|
||
imageinfo_heightspace = -self.ImageHandler.get_geotransform(out_path1)[5] # 投影后的分辨率
|
||
para_dict.update({"imageinfo_widthspace":imageinfo_widthspace})
|
||
para_dict.update({"imageinfo_heightspace":imageinfo_heightspace})
|
||
para_dict.update({"NominalResolution":imageinfo_widthspace})
|
||
|
||
WidthInMeters = imageinfo_width*imageinfo_widthspace # 投影后的分辨率×宽度
|
||
para_dict.update({"WidthInMeters":WidthInMeters})
|
||
|
||
image_array = self.ImageHandler.get_band_array(image_path)
|
||
a2 = np.where(np.isnan(image_array), 999999, image_array)
|
||
MinValue = np.min(a2)
|
||
a3 = np.where(np.isnan(image_array), -999999, image_array)
|
||
MaxValue = np.max(a3)
|
||
|
||
para_dict.update({"MaxValue":MaxValue})
|
||
para_dict.update({"MinValue":MinValue})
|
||
|
||
get_scope = self.ImageHandler.get_scope(image_path)
|
||
point_upleft, point_upright, point_downleft, point_downright=get_scope[0], get_scope[1], get_scope[2], get_scope[3]
|
||
para_dict.update({"imageinfo_corner_topLeft_latitude": point_upleft[1]})
|
||
para_dict.update({"imageinfo_corner_topLeft_longitude": point_upleft[0]})
|
||
para_dict.update({"imageinfo_corner_topRight_latitude": point_upright[1]})
|
||
para_dict.update({"imageinfo_corner_topRight_longitude": point_upright[0]})
|
||
para_dict.update({"imageinfo_corner_bottomLeft_latitude": point_downleft[1]})
|
||
para_dict.update({"imageinfo_corner_bottomLeft_longitude": point_downleft[0]})
|
||
para_dict.update({"imageinfo_corner_bottomRight_latitude": point_downright[1]})
|
||
para_dict.update({"imageinfo_corner_bottomRight_longitude": point_downright[0]})
|
||
|
||
longitude_max=np.array([point_upleft[0], point_upright[0], point_downleft[0], point_downright[0]]).max()
|
||
longitude_min=np.array([point_upleft[0], point_upright[0], point_downleft[0], point_downright[0]]).min()
|
||
latitude_max=np.array([point_upleft[1], point_upright[1], point_downleft[1], point_downright[1]]).max()
|
||
latitude_min=np.array([point_upleft[1], point_upright[1], point_downleft[1], point_downright[1]]).min()
|
||
imageinfo_center_latitude=(latitude_max+latitude_min)/2
|
||
imageinfo_center_longitude=(longitude_max+longitude_min)/2
|
||
para_dict.update({"imageinfo_center_latitude": imageinfo_center_latitude})
|
||
para_dict.update({"imageinfo_center_longitude": imageinfo_center_longitude})
|
||
|
||
# self.para_dict.update({"productType": "GTC"}) # 设置产品类型
|
||
para_dict.update({"productFormat": "TIF"})
|
||
productGentime = datetime.datetime.now()
|
||
para_dict.update({"productGentime": productGentime})
|
||
para_dict.update({"unit": "none"}) # 设置单位
|
||
para_dict.update({"NoDataValue": "nan"})
|
||
para_dict.update({"productLevel": "4"}) # 设置图像位深度
|
||
|
||
image_array = self.ImageHandler.get_band_array(image_path)
|
||
try: # 设置图像位深度
|
||
gdal_dtypes = {
|
||
'int8': gdal.GDT_Byte,
|
||
'unit16': gdal.GDT_UInt16,
|
||
'int16': gdal.GDT_Int16,
|
||
'unit32': gdal.GDT_UInt32,
|
||
'int32': gdal.GDT_Int32,
|
||
'float32': gdal.GDT_Float32,
|
||
'float64': gdal.GDT_Float64,
|
||
}
|
||
bit_dtypes = {
|
||
'int8': 8,
|
||
'unit16': 16,
|
||
'int16': 16,
|
||
'unit32': 32,
|
||
'int32': 32,
|
||
'float32': 32,
|
||
'float64': 64,
|
||
}
|
||
if not gdal_dtypes.get(image_array.dtype.name, None) is None:
|
||
bit_num = str(bit_dtypes[image_array.dtype.name])
|
||
datatype=bit_num+"bit"
|
||
else:
|
||
datatype = str(32) + "bit"
|
||
# datatype = str(gdal.GDT_Float32)+"bit"
|
||
para_dict.update({"imagebit": datatype})
|
||
except Exception:
|
||
para_dict.update({"imagebit": "None"})
|
||
|
||
HH, HV, VH ,VV= image_pair[0], image_pair[1], image_pair[2], image_pair[3]
|
||
if HH == 0:
|
||
HH = "delete"
|
||
else:
|
||
HH = "NULL"
|
||
para_dict.update({"imageinfo_QualifyValue_HH": HH})
|
||
if HV==0:
|
||
HV = "delete"
|
||
else:
|
||
HV = "NULL"
|
||
para_dict.update({"imageinfo_QualifyValue_HV": HV})
|
||
if VH==0:
|
||
VH = "delete"
|
||
else:
|
||
VH = "NULL"
|
||
para_dict.update({"imageinfo_QualifyValue_VH": VH})
|
||
if VV==0:
|
||
VV = "delete"
|
||
else:
|
||
VV = "NULL"
|
||
para_dict.update({"imageinfo_QualifyValue_VV": VV})
|
||
|
||
return para_dict
|
||
|
||
def calu_dem_nature(self, dem_path, out_dem_path1, out_dem_path2, sampling_f, para_A_arr):
|
||
"""
|
||
正射需要单独加上dem影像的信息
|
||
dem_path：mergedDEM.tif路径
|
||
out_dem_path1：将mergedDEM.tif由地理坐标转化为平面坐标的保存路径
|
||
out_dem_path2：将mergedDEM.tif由平面坐标转化为地理坐标的保存路径
|
||
sampling_f： 采样率
|
||
para_A_arr：四次多项式模型的参数数组
|
||
"""
|
||
para_dict2 = {}
|
||
proj = self.ImageHandler.get_projection(dem_path) # 输出的影像若是投影坐标系则先转成地理坐标系
|
||
keyword = proj.split("[", 2)[0] # 若是地理坐标系则pass
|
||
if keyword == "GEOGCS":
|
||
pass
|
||
elif keyword == "PROJCS":
|
||
pp.trans_projcs2geogcs(out_dem_path2, dem_path)
|
||
dem_path = out_dem_path2
|
||
elif len(keyword) == 0 or keyword.strip() == "" or keyword.isspace() is True:
|
||
raise Exception('image projection is missing!')
|
||
pp.trans_geogcs2projcs(out_dem_path1, dem_path) # 坐标投影, 地理转平面投影坐标
|
||
DEM_widthspace = self.ImageHandler.get_geotransform(out_dem_path1)[1] # 投影后的分辨率
|
||
DEM_heightspace = -self.ImageHandler.get_geotransform(out_dem_path1)[5] # 投影后的分辨率
|
||
para_dict2.update({"DEM_widthspace":DEM_widthspace})
|
||
para_dict2.update({"DEM_heightspace":DEM_heightspace})
|
||
# tree = ElementTree() # 获取DEMProduct：dem产品来源、DEMDate：dem对应的日期
|
||
# tree.parse(dem_meta) # 影像头文件
|
||
# root = tree.getroot()
|
||
# productinfo = root.find("metadata")
|
||
# DEMProduct = list(productinfo)[0].tag
|
||
# para_dict2.update({"DEM_DEMProduct":DEMProduct})
|
||
#
|
||
# DEMDate = root.find("metadata").find(DEMProduct).text
|
||
# para_dict2.update({"DEM_DEMDate": DEMDate})
|
||
get_scope = self.ImageHandler.get_scope(dem_path) # 获取mergedDEM.tif数据的四个角的经纬度信息
|
||
point_upleft, point_upright, point_downleft, point_downright=get_scope[0], get_scope[1], get_scope[2], get_scope[3]
|
||
para_dict2.update({"DEM_corner_topLeft_latitude": point_upleft[1]})
|
||
para_dict2.update({"DEM_corner_topLeft_longitude": point_upleft[0]})
|
||
para_dict2.update({"DEM_corner_topRight_latitude": point_upright[1]})
|
||
para_dict2.update({"DEM_corner_topRight_longitude": point_upright[0]})
|
||
para_dict2.update({"DEM_corner_bottomLeft_latitude": point_downleft[1]})
|
||
para_dict2.update({"DEM_corner_bottomLeft_longitude": point_downleft[0]})
|
||
para_dict2.update({"DEM_corner_bottomRight_latitude": point_downright[1]})
|
||
para_dict2.update({"DEM_corner_bottomRight_longitude": point_downright[0]})
|
||
#para_dict2.update({"orthoModel_samplingrate": sampling_f})
|
||
|
||
para_dict2.update({"satalliteOrbitModel_parameter_X_a0": para_A_arr[0, 0]}) # 获取四次多项式模型6个参数的数值
|
||
para_dict2.update({"satalliteOrbitModel_parameter_X_a1": para_A_arr[1, 0]})
|
||
para_dict2.update({"satalliteOrbitModel_parameter_X_a2": para_A_arr[2, 0]})
|
||
para_dict2.update({"satalliteOrbitModel_parameter_X_a3": para_A_arr[3, 0]})
|
||
para_dict2.update({"satalliteOrbitModel_parameter_X_a4": para_A_arr[4, 0]})
|
||
|
||
para_dict2.update({"satalliteOrbitModel_parameter_Y_b0": para_A_arr[0, 1]})
|
||
para_dict2.update({"satalliteOrbitModel_parameter_Y_b1": para_A_arr[1, 1]})
|
||
para_dict2.update({"satalliteOrbitModel_parameter_Y_b2": para_A_arr[2, 1]})
|
||
para_dict2.update({"satalliteOrbitModel_parameter_Y_b3": para_A_arr[3, 1]})
|
||
para_dict2.update({"satalliteOrbitModel_parameter_Y_b4": para_A_arr[4, 1]})
|
||
|
||
para_dict2.update({"satalliteOrbitModel_parameter_Z_c0": para_A_arr[0, 2]})
|
||
para_dict2.update({"satalliteOrbitModel_parameter_Z_c1": para_A_arr[1, 2]})
|
||
para_dict2.update({"satalliteOrbitModel_parameter_Z_c2": para_A_arr[2, 2]})
|
||
para_dict2.update({"satalliteOrbitModel_parameter_Z_c3": para_A_arr[3, 2]})
|
||
para_dict2.update({"satalliteOrbitModel_parameter_Z_c4": para_A_arr[4, 2]})
|
||
|
||
para_dict2.update({"satalliteOrbitModel_parameter_Vx_d0": para_A_arr[0, 3]})
|
||
para_dict2.update({"satalliteOrbitModel_parameter_Vx_d1": para_A_arr[1, 3]})
|
||
para_dict2.update({"satalliteOrbitModel_parameter_Vx_d2": para_A_arr[2, 3]})
|
||
para_dict2.update({"satalliteOrbitModel_parameter_Vx_d3": para_A_arr[3, 3]})
|
||
para_dict2.update({"satalliteOrbitModel_parameter_Vx_d4": para_A_arr[4, 3]})
|
||
|
||
para_dict2.update({"satalliteOrbitModel_parameter_Vy_e0": para_A_arr[0, 4]})
|
||
para_dict2.update({"satalliteOrbitModel_parameter_Vy_e1": para_A_arr[1, 4]})
|
||
para_dict2.update({"satalliteOrbitModel_parameter_Vy_e2": para_A_arr[2, 4]})
|
||
para_dict2.update({"satalliteOrbitModel_parameter_Vy_e3": para_A_arr[3, 4]})
|
||
para_dict2.update({"satalliteOrbitModel_parameter_Vy_e4": para_A_arr[4, 4]})
|
||
|
||
para_dict2.update({"satalliteOrbitModel_parameter_Vz_f0": para_A_arr[0, 5]})
|
||
para_dict2.update({"satalliteOrbitModel_parameter_Vz_f1": para_A_arr[1, 5]})
|
||
para_dict2.update({"satalliteOrbitModel_parameter_Vz_f2": para_A_arr[2, 5]})
|
||
para_dict2.update({"satalliteOrbitModel_parameter_Vz_f3": para_A_arr[3, 5]})
|
||
para_dict2.update({"satalliteOrbitModel_parameter_Vz_f4": para_A_arr[4, 5]})
|
||
return para_dict2
|
||
|
||
|
||
class CreateStadardXmlFile:
|
||
"""读取字典中的属性值，生成一个标准的xml文件"""
|
||
def __init__(self, xml_path, para_xml_path, par_dict, par_dict2, path):
|
||
"""
|
||
xml_path:模板路径
|
||
para_xml_path:算法配置文件的路径
|
||
par_dict:字典
|
||
path:xml模板输出路径
|
||
"""
|
||
self.par_dict = par_dict
|
||
self.par_dict2 = par_dict2
|
||
self.path = path
|
||
shutil.copy(xml_path, path)
|
||
pass
|
||
|
||
def create_standard_xml(self):
|
||
"""将字典中的信息写入到copy的xml文件中"""
|
||
tree = ElementTree()
|
||
tree.parse(self.path) # 影像头文件
|
||
root = tree.getroot()
|
||
|
||
productinfo = root.find("productinfo")
|
||
for key, value in self.par_dict.items():
|
||
if key.split("_")[0] != "imageinfo":
|
||
productinfo.find(key).text = str(value)
|
||
elif key.split("_")[0] == "imageinfo":
|
||
imageinfo = productinfo.find("imageinfo")
|
||
if key.split("_")[1] in ["EarthModel", "ProjectModel", "width", "height", "widthspace", "heightspace"]:
|
||
imageinfo.find(key.split("_")[1]).text = str(value)
|
||
elif key.split("_")[1] == "center":
|
||
center = imageinfo.find("center")
|
||
center.find(key.split("_")[2]).text = str(value)
|
||
elif key.split("_")[1] == "corner":
|
||
corner = imageinfo.find("corner")
|
||
corner.find(key.split("_")[2]).find(key.split("_")[3]).text = str(value)
|
||
elif key.split("_")[1] == "QualifyValue":
|
||
QualifyValue = imageinfo.find("QualifyValue")
|
||
if value =="delete":
|
||
element_QualifyValue = list(QualifyValue)
|
||
for i in element_QualifyValue:
|
||
if i.tag == key.split("_")[2]:
|
||
QualifyValue.remove(i)
|
||
else:
|
||
QualifyValue.find(key.split("_")[2]).text = str(value)
|
||
pass
|
||
orthoModel = root.find("processinfo").find("orthoModel") # 写入四次多项式模型
|
||
for key, value in self.par_dict2.items():
|
||
if key.split("_")[0] == "satalliteOrbitModel":
|
||
satalliteOrbitModel = orthoModel.find("satalliteOrbitModel")
|
||
satalliteOrbitModel.find(key.split("_")[1]).find(key.split("_")[2]).find(key.split("_")[3]).text = str(value)
|
||
elif key.split("_")[0] == "DEM": # 写入dem四个角坐标
|
||
DEM= orthoModel.find("DEM")
|
||
if key.split("_")[1] == "corner":
|
||
corner = DEM.find("corner")
|
||
corner.find(key.split("_")[2]).find(key.split("_")[3]).text = str(value)
|
||
elif key.split("_")[1] == "widthspace" or key.split("_")[1] == "heightspace":
|
||
DEM.find(key.split("_")[1]).text = str(value)
|
||
elif key.split("_")[1] == "samplingrate":
|
||
orthoModel.find(key.split("_")[1]).text = str(value)
|
||
tree.write(self.path, encoding="utf-8", xml_declaration=True)
|
||
|
||
#
|
||
# if __name__ == '__main__':
|
||
#
|
||
# xml_path = "./model_meta.xml"
|
||
# tem_folder= r"E:\microproduct\测试"
|
||
# image_path = r"E:\microproduct\测试\GF3_MYN_QPSI_011437_E98_HH_AtmosphericDelay.tif" # 输入影像
|
||
# out_path = os.path.join(tem_folder, "trans_geo_projcs.tif")
|
||
# image_pair=[1, 1, 1, 0]
|
||
# par_dict = CreateDict(image_path, image_pair, out_path).calu_nature() # 计算属性字典
|
||
#
|
||
# out_xml_path = os.path.join(tem_folder, "creat_standard.meta.xml") # 输出xml路径
|
||
# CreateStadardXmlFile(xml_path, par_dict, out_xml_path).create_standard_xml()
|