添加四极化后向散射系数到特征图中

.gitignore

@@ -35,3 +35,7 @@
 
 /soilSalinity/tool/
 /landcover_c_sar/tool/
+/atmosphericDelay-C-SAR/tool/
+/deformation-C-SAR/tool/
+/dem-C-SAR/tool/
+/soilMoistureTop/tool/

tool/algorithm/polsarpro/AHVToPolsarpro.py

@@ -12,6 +12,8 @@ import os
 import numpy as np
 import glob
 import struct
+
+from tool.algorithm.algtools.MetaDataHandler import MetaDataHandler
 from tool.algorithm.image.ImageHandle import ImageHandler
 
 
@@ -461,6 +463,65 @@ class AHVToPolsarpro:
             self._hh_hv_vh_vv_path_list = []
         return out_dir
 
+    @staticmethod
+    def sar_backscattering_sigma(in_sar_tif, meta_file_path, out_sar_tif, replece_VV=False, is_DB=True):
+
+        # 读取原始SAR影像
+        proj, geotrans, in_data = ImageHandler.read_img(in_sar_tif)
+
+        # 计算强度信息
+        I = np.array(in_data[0], dtype="float32")
+        Q = np.array(in_data[1], dtype="float32")
+
+        where_9999_0 = np.where(I == -9999)
+        where_9999_1 = np.where(Q == -9999)
+        I[where_9999_0] = 1.0
+        Q[where_9999_1] = 1.0
+
+        I2 = np.square(I)
+        Q2 = np.square(Q)
+        intensity_arr = I2 + Q2
+
+        # 获取极化类型
+        if 'HH' in os.path.basename(in_sar_tif):
+            polarization = 'HH'
+        elif 'HV' in os.path.basename(in_sar_tif):
+            polarization = 'HV'
+        elif 'VH' in os.path.basename(in_sar_tif):
+            polarization = 'VH'
+        elif 'VV' in os.path.basename(in_sar_tif):
+            polarization = 'VV'
+            if replece_VV:
+                polarization = 'HV'  # 土壤水分算法中可能会用HV替换VV
+        elif 'DH' in os.path.basename(in_sar_tif):
+            polarization = 'HH'
+        else:
+            raise Exception('there are not HH、HV、VH、VV in path:', in_sar_tif)
+
+        # 获取参数
+        QualifyValue = MetaDataHandler.get_QualifyValue(meta_file_path, polarization)
+        Kdb = MetaDataHandler.get_Kdb(meta_file_path, polarization)
+
+        # 计算后向散射系数
+        # 对数形式
+        coef_arr = 10 * (np.log10(intensity_arr * ((QualifyValue / 32767) ** 2))) - Kdb
+        coef_arr[np.isnan(coef_arr)] = -9999
+        coef_arr[np.isinf(coef_arr)] = -9999
+        coef_arr[where_9999_0] = -9999
+        coef_arr[where_9999_1] = -9999
+        ## 输出的SAR后向散射系数产品
+        # ImageHandler.write_img(out_sar_tif, proj, geotrans, coef_arr, 0)
+
+        tif_array = np.power(10.0, coef_arr / 10.0)  # dB --> 线性值  后向散射系数
+
+        tif_array[np.isnan(tif_array)] = 0
+        tif_array[np.isinf(tif_array)] = 0
+        tif_array[where_9999_0] = 0
+        tif_array[where_9999_1] = 0
+
+        ImageHandler.write_img(out_sar_tif, proj, geotrans, tif_array, 0)
+        return True
+
 
 
 if __name__ == '__main__':

vegetationPhenology/VegetationPhenology.xml

@@ -37,7 +37,7 @@
         <ParaType>File</ParaType>
         <DataType>tar.gz</DataType>
 		<ParaSource>Man</ParaSource>
-		<ParaValue>E:\VP-lk\GF3C_MYC_QPSI_006270_E100.4_N27.0_20230615_L1A_AHV_L10000158764-ortho.tar.gz</ParaValue>
+		<ParaValue>G:\VegetationPhenology-likun\rusuoces\GF3C_MYC_QPSI_006270_E100.4_N27.0_20230615_L1A_AHV_L10000158764-ortho.tar.gz</ParaValue>
         <MaxValue>DEFAULT</MaxValue>
         <MinValue>DEFAULT</MinValue>
         <OptionValue>DEFAULT</OptionValue>
@@ -51,7 +51,7 @@
         <ParaType>File</ParaType>
         <DataType>csv</DataType>
 		<ParaSource>Man</ParaSource>
-        <ParaValue>E:\VP-lk\VegetationPhe_samples.csv</ParaValue>
+        <ParaValue>G:\VegetationPhenology-likun\rusuoces\VegetationPhe_samples.csv</ParaValue>
          <MaxValue>DEFAULT</MaxValue>
         <MinValue>DEFAULT</MinValue>
         <OptionValue>DEFAULT</OptionValue>
@@ -65,7 +65,7 @@
         <ParaType>File</ParaType>
         <DataType>tif</DataType>
 		<ParaSource>Man</ParaSource>
-		<ParaValue>E:\VP-lk\landcover.tif</ParaValue>
+		<ParaValue>G:\VegetationPhenology-likun\rusuoces\landcover.tif</ParaValue>
         <EnModification>True</EnModification>
         <EnMultipleChoice>False</EnMultipleChoice>
         <Control>File</Control>

vegetationPhenology/VegetationPhenologyMain.py

@@ -26,7 +26,7 @@ from tool.algorithm.xml.CreatMetafile import CreateMetafile
 from tool.algorithm.algtools.ROIAlg import ROIAlg as alg
 from VegetationPhenologyXmlInfo import CreateDict, CreateStadardXmlFile
 from VegetationPhenologyAuxData import PhenoloyMeasCsv_geo
-
+from tool.algorithm.polsarpro.AHVToPolsarpro import AHVToPolsarpro
 from tool.algorithm.xml.CreateMetaDict import CreateMetaDict, CreateProductXml
 from tool.file.fileHandle import fileHandle
 from tool.algorithm.algtools.ROIAlg import ROIAlg as roi
@@ -95,7 +95,7 @@ class PhenologyMain:
 
         self.__workspace_path = self.__alg_xml_handler.get_workspace_path()
         self.__create_work_space()
-        self.__processing_paras = InitPara.init_processing_paras(self.__input_paras)
+        self.__processing_paras = InitPara.init_processing_paras(self.__input_paras, self.__workspace_preprocessing_path)
         self.__processing_paras.update(InitPara(DEBUG).get_mult_tar_gz_infs(self.__processing_paras, self.__workspace_preprocessing_path))
         SrcImagePath = self.__input_paras["AHVS"]['ParaValue']
         paths = SrcImagePath.split(';')
@@ -295,23 +295,38 @@ class PhenologyMain:
         # train_data, test_data, type_map = csvh.trans_VegePhenology_measdata_dic(csvh.readcsv(self.__processing_paras['MeasuredData']), ori_sim_path)
         logger.info("read phenology Measure.csv success!")
 
+        # 添加四极化后向散射系数到特征图中
+        feature_tif_dir = os.path.join(self.__workspace_processing_path, name, 'features')
+        if not os.path.exists(feature_tif_dir):
+            os.makedirs(feature_tif_dir)
+        origin_xml = self.__processing_paras[name + "_Origin_META"]
+
         # 特征分解
         hh_hv_vh_vv_list = ["","","",""]
         for key in key_list:
             if "_HH" in key :
                 hh_hv_vh_vv_list[0] = self.__preprocessed_paras[key]
+                outPath = os.path.join(feature_tif_dir, 'cal_HH.tif')
+                AHVToPolsarpro.sar_backscattering_sigma(hh_hv_vh_vv_list[0], origin_xml, outPath)
             if "_HV" in key:
                 hh_hv_vh_vv_list[1] = self.__preprocessed_paras[key]
+                outPath = os.path.join(feature_tif_dir, 'cal_HV.tif')
+                AHVToPolsarpro.sar_backscattering_sigma(hh_hv_vh_vv_list[1], origin_xml, outPath)
             if "_VH" in key:
                 hh_hv_vh_vv_list[2] = self.__preprocessed_paras[key]
+                outPath = os.path.join(feature_tif_dir, 'cal_VH.tif')
+                AHVToPolsarpro.sar_backscattering_sigma(hh_hv_vh_vv_list[2], origin_xml, outPath)
             if "_VV" in key:
                 hh_hv_vh_vv_list[3] = self.__preprocessed_paras[key]
+                outPath = os.path.join(feature_tif_dir, 'cal_VV.tif')
+                AHVToPolsarpro.sar_backscattering_sigma(hh_hv_vh_vv_list[3], origin_xml, outPath)
 
         cols = self.imageHandler.get_img_width(hh_hv_vh_vv_list[0])
         rows = self.imageHandler.get_img_height(hh_hv_vh_vv_list[0])
 
         featureInput = self.__getInputFeatures()
         feature_dir = CreateFeature.decompose_single_tar(hh_hv_vh_vv_list, self.__workspace_processing_path, self.__workspace_preprocessing_path, name, self._env_str, rows, cols, FILTER_SIZE=3, debug=DEBUG, FeatureInput=featureInput)
+
         feature_geo_dir = self.features_geo(feature_dir, paramter, sim_ori_path, name)
         # # 获取训练集提取特征的信息
         # ids = []