本发明属于陆地和海洋参数测量领域，具体涉及一种机会信号反射遥感的全极化单次反射仿真方法及系统。
背景技术：
根区土壤水分处于土壤-植物-大气连续体的核心位置，对干旱半干旱地区作物生长、植被恢复和水土流失过程均具有重要影响。根区土壤水分的获取对于农业以及粮食生产尤为重要，同时土壤水分剖面信息的获取能与当前的陆面/水文模式相契合，更好地支撑数据同化研究、天气与气候预报、水资源管理等应用需求。
微波遥感技术是获取土壤水分的有效监测手段。目前合成孔径雷达和微波辐射计是两种较为常见和成熟的微波遥感方式。土壤水分和海洋盐度观测卫星(smos)和土壤水分探测卫星(smap)都是工作在l波段的微波遥感卫星。但在l波段，只能穿透浅层地表，因此适用于土壤表层(约5cm)的土壤水分获取，而且由于穿透强度的限制，l波段只能穿透中等厚度的植被冠层。
对于根区土壤水分的获取，则需要使传感器工作在更低的频率范围，如频率范围是240-380mhz的p波段。其中，p波段可以穿透深度约40cm,同时可以很好的规避茂密植被的影响。2012年美国在开展了机载p波段实验，即airmoss(airbornemicrowaveobservatoryofsub-canopyandsubsurface)试验。但是p波段的星载观测则面临着天线尺寸过大、rfi(radio-frequencyinterference)影响和频段分配竞争等挑战。
机会信号反射遥感为p波段根区土壤水分获取提供了新的机遇。机会信号反射遥感是一种现在存在的遥感技术手段，它利用现有的通讯卫星或者导航卫星作为信号源，采用专门的反射信号接收机来接收机会卫星信号经地表反射的信号，对地物参数进行遥感监测。由于使用数字通信卫星作为信号源，不需要研制专门的发射机，遥感过程中，只需要研制反射信号接收机接收地表的反射信号，因此，这种遥感手段相比于不采用机会信号反射遥感的p波段的雷达，可以克服天线尺寸过大、rfi和频段分配竞争等问题，存在着成本低、功耗小、价格便宜、时空分辨率高等诸多优点。
机会信号反射遥感作为一种新型的根区土壤水分遥感手段，其前期的机理仿真分析模型为分析解释卫星观测数据、卫星数据的仿真模拟、卫星数据同化、发展地表参数的定量反演等提供关键的机理工具。然而，目前，现有的仿真分析系统多是针对辐射计或者单站雷达或者gnss-r遥感的，对于机会反射信号遥感，针对根区土壤水分，尚未发现有相应的仿真分析系统。
技术实现要素：
本发明的目的在于提供一种机会信号反射遥感的全极化单次反射仿真方法及系统，从而通过机会反射信号遥感技术在根区土壤水分发生变化时预测出相应的接收机波形的变化。
为了实现上述目的，本发明提供了一种机会信号反射遥感的全极化单次反射仿真方法，其工作于p波段，包括：
s1：将一待测土壤剖面划分为多个分层，并根据待测土壤剖面的输入参数建立待测土壤剖面的介电常数模型，得到待测土壤剖面的各个分层的介电常数；
s2：利用所述步骤s1中的待测土壤剖面的各个分层的介电常数来获取每个分层的透过率、反射率和每两个相邻分层之间的衰减；
s3：根据各个分层的粗糙程度，分别构建光滑表面反射率矩阵和粗糙表面反射率矩阵；
s4：根据所述步骤s3中的光滑表面反射率矩阵和粗糙表面反射率矩阵，利用极化合成的方法计算各个分层在各个极化状态上的反射特性，该反射特性为各个分层的双站雷达散射截面；
s5：根据所述步骤s2中的每个分层的透过率、反射率和每两个相邻分层之间的衰减以及所述步骤s4中的反射特性来获得单次反射的总接收机波形。
所述全极化单次反射仿真系统的输入参数包括：土壤质地、土壤温度、土壤块密度、地表粗糙度和土壤含水量。
所述待测土壤剖面的介电常数模型是基于混合介质介电常数模型修正得到的，所述介电常数的实部采用线性校正来校正，所述介电常数的虚部ξ"通过修正有效电导率来修正。
所述介电常数的实部的校正公式为：
ξ'＝1.15ξ"-0.68
式中，ξ'和ξ"分别是介电常数的实部和虚部；
有效电导率σeff的修正公式为：
σeff＝0.0467 0.2204ρb-0.4111*s 0.6614*c，
式中，σeff是有效电导率，ρb是土壤块密度，s和c分别是沙土含量和粘土含量。
所述每个分层的透过率、反射率和每两个相邻分层之间的衰减为：
其中，下标i和j分别是彼此相邻的分层中的第i层和第j层，tv和th分别为第i层的垂直和水平极化的透过率，rv和rh分别是第i层的垂直和水平极化的反射率；θ为入射角，εi为第i层的介电常数；λ为波长，α为第i层和第j层之间的衰减系数。
在所述步骤s3中，对于分层为光滑表面的情况，则构建光滑表面反射率矩阵，所述光滑表面反射率矩阵是通过假定每层只发生镜像反射，基于菲涅尔公式来构建的；针对分层为粗糙表面的情况，则构建粗糙表面反射率矩阵，所述粗糙表面反射率矩阵是通过建立高级积分方程模型计算粗糙表面的反射率而得到的。
在所述步骤s3中，极化合成的计算公式为：
其中，σqp为各个分层的双站雷达散射截面，下标q和p分别为接受和发射的极化状态；ψ和τ分别为椭倾角和椭率角；yq、yp为修改的stokes矢量；q为旋转矩阵；m为光滑表面反射率矩阵或者粗糙表面反射率矩阵。
在所述步骤s5中，所述总接收机波形为：
式中，τ为时间延迟，上标dir表示接收机能量中的直射部分，i表示分层中的第i层，n为层数，ρi为第i层的延迟，lca为gps的c/a码片长度，j为复数单位，
其中，第i层的振幅ui、相位φi和延迟ρi根据所述每个分层的透过率、反射率、每两个相邻分层之间的衰减以及双站雷达散射截面计算得到。
另一方面，本发明提供一种机会信号反射遥感的全极化单次反射仿真系统，其工作于p波段，包括：
土壤剖面介电常数计算模块，设置为将一待测土壤剖面划分为多个分层，并根据待测土壤剖面的输入参数建立待测土壤剖面的介电常数模型，得到待测土壤剖面的各个分层的介电常数；
层内反射率计算模块、层内透过率计算模块和层间衰减率计算模块，分别设置为所述各个分层的介电常数来获取每个分层的透过率、反射率和每两个相邻分层之间的衰减；
全极化反射率计算模块，根据所述层内反射率计算模块得到的反射率得到全极化的反射特性，其包括：镜像反射模块、非镜像反射模块以及与这两者均相连的极化合成模块，镜像反射模块和非镜像反射模块分别设置为根据分层的粗糙程度构建光滑表面反射率矩阵和粗糙表面反射率矩阵；极化合成模块设置为根据所述光滑表面反射率矩阵和粗糙表面反射率矩阵，利用极化合成的方法计算各个分层在各个极化状态上的反射特性，该反射特性为各个分层的双站雷达散射截面；
单次反射能量计算模块设置为根据所述每个分层的透过率、反射率和每两个相邻分层之间的衰减以及所述各个分层在各个极化状态上的反射特性来获得总接收机波形。
所述全极化单次反射仿真系统的输入参数包括：土壤质地、土壤温度、土壤块密度、地表粗糙度和土壤含水量。
本发明的机会信号反射遥感的全极化单次反射仿真方法通过利用混合介质介电常数模型计算土壤的各个分层的介电常数，对根区土壤水分建立相应的仿真，从而在当根区土壤水分发生变化时，预测出相应的接收机波形的变化，为分析解释卫星观测数据、卫星数据的仿真模拟、卫星数据同化、发展地表参数的定量反演等提供关键的机理工具；此外，本发明的全极化单次反射仿真方法采用极化合成的方法，使其可以应用在机会信号反射遥感中，并且可以计算各种极化的反射率和波形信息，即机会信号卫星发射右旋圆极化信号，接收到圆极化(右旋圆极化和左旋圆极化)和线极化(水平极化和垂直极化)的总接收机波形，由于电磁波的极化对目标的介电常数、物理特性、几何形状和取向等比较敏感，因而极化测量可以大大提高对目标各种信息的获取能力；再者，现有的模型只可以计算每个小分层内的镜像反射，而本发明通过计算每个分层内的入射过来的机会卫星信号的镜像散射和非镜像的漫散射，从而避免了对各个方向的漫散射的忽略，使计算结果更准确。
附图说明
图1为具有多个分层的土壤剖面的单次反射信号的追踪示意图。
图2为根据本发明的一个实施例的p波段的机会信号反射遥感的全极化单次反射仿真系统的结构示意图。
具体实施方式
以下结合具体实施例，对本发明做进一步说明。应理解，以下实施例仅用于说明本发明而非用于限制本发明的范围。
如图1所示，根据本发明的一个实施例的机会信号反射遥感的全极化单次反射仿真方法工作于p波段，其包括以下步骤：
步骤s1：采用一土壤剖面介电常数计算模块1将一待测土壤剖面划分为多个分层，并根据待测土壤剖面的输入参数建立待测土壤剖面的介电常数模型，得到待测土壤剖面的各个分层的介电常数；
由此，模型计算的时候将土壤剖面划分为不同的多个分层，每层内发生单次散射，具有多个分层的土壤剖面的具体示意图如图1所示。所述划分方式为包括等间隔或者任意间隔的划分。
所述全极化单次反射仿真系统的输入参数包括：土壤质地(沙土含量、粘土含量)、土壤温度、土壤块密度、地表粗糙度(相关长度、均方根高度)、土壤含水量。这些参数都可以通过实际测量获得，或者利用现有文献的输入参数，或者使用者根据经验得到的数据。
待测土壤剖面的介电常数模型是基于dobson等人的混合介质介电常数模型修正得到的。0.3-1.3ghz情况下，需要针对现有的混合介质介电常数模型进行修改，以建立待测土壤剖面的介电常数模型。
所述介电常数的实部通过采用线性校正来校正。
所述介电常数的实部的校正公式为：
ξ'＝1.15ξ"-0.68(1)
式中，ξ'和ξ"分别是介电常数的实部和虚部。
所述介电常数的虚部ξ"，则通过修正有效电导率来修正。
其中，有效电导率σeff的修正采用以下公式：
σeff＝0.0467 0.2204ρb-0.4111*s 0.6614*c(2)
式中，σeff是有效电导率，ρb是土壤块密度，s和c分别是沙土含量和粘土含量。其中，土壤块密度ρb、沙土含量s和粘土含量c可以是使用者根据相应的实验测得，或者根据经验值给定。
介电常数详细的计算方法参照参考文献【hallikainenmt.microwavedielectricbehaviorofwetsoil-parti:empiricalmodelsandexperimentalobservations[j].ieeetrans.geosci.remotesensing,1985,23.】。
步骤s2：分别采用层内反射率计算模块2、层内透过率计算模块3和层间衰减率计算模块4，利用所述步骤s1中的待测土壤剖面的各个分层的介电常数来获取每个分层的透过率、反射率和每两个相邻分层之间的衰减。
其中，待测土壤的每个分层的透过率、反射率和每两个相邻分层之间的衰减为：
其中，下标i和j分别是彼此相邻的分层中的第i层和第j层，tv和th分别为第i层的垂直和水平极化的透过率，rv和rh分别是第i层的垂直和水平极化的反射率；θ为入射角，εi为第i层的介电常数；λ为波长，α为第i层和第j层之间的衰减系数，其受介电常数的虚部影响。
所述步骤s2求得的待测土壤的每个分层的透过率、反射率和每两个相邻分层之间的衰减可以用来计算下文的步骤s5中的总接收机波形的振幅、延迟。
步骤s3：根据各个分层的粗糙程度，采用镜像反射模块51和非镜像反射模块52分别构建光滑表面反射率矩阵和粗糙表面反射率矩阵。由此，可以为后续步骤s4中的极化合成，也就是使得反射率可以计算各种极化的反射特性做准备。
其中，对于分层为光滑表面的情况，则构建光滑表面反射率矩阵，该分层所对应的光滑表面反射率矩阵是通过假定每层只发生镜像反射，基于菲涅尔公式来构建的。
所述光滑表面反射率矩阵为：
其中，下标v和h分别为垂直极化和水平极化，rv和rh分别是各个分层的垂直和水平极化的反射率，其可以通过上文中的公式(5)和公式(6)得到；re是取实部，im是虚部，*表示共轭。
其中，针对分层为粗糙表面的情况，则构建粗糙表面反射率矩阵，该分层所对应的粗糙表面反射率矩阵是通过建立高级积分方程模型计算粗糙地表的反射率而得到的。
其中，粗糙表面的反射率为：
其中，
基尔霍夫项
所述粗糙表面的反射率在各个极化状态上的计算结果分别为所述粗糙表面反射率矩阵的元素。
步骤s4：极化合成。采用极化合成模块53，根据所述步骤s3中的光滑表面反射率矩阵和粗糙表面反射率矩阵，利用极化合成的方法计算各个分层在各个极化状态上的反射特性，该反射特性为各个分层的双站雷达散射截面(brcs)
其中，极化合成的计算公式如下：
其中，σqp为各个分层的双站雷达散射截面，下标q和p分别为接受和发射的极化状态，极化状态具体包括圆极化(右旋圆极化和左旋圆极化)和线极化(水平极化和垂直极化)；ψ和τ分别为椭倾角和椭率角，不同的极化对应不同的椭倾角和椭率角；yq、yp为修改的stokes矢量；q为旋转矩阵；m为散射矩阵，即为上文所述的光滑表面反射率矩阵或者粗糙表面反射率矩阵，当该分层为光滑表面时，m为光滑表面反射率矩阵，该分层为粗糙表面时，m为粗糙表面反射率矩阵。
由此，将光滑表面反射率矩阵和粗糙表面反射率矩阵作为公式中的m输入，既可以得到光滑表面或者粗糙表面的双站雷达散射截面
极化合成的具体计算方法可以参考相关文献【ulabyft,elachic.radarpolarimetryforgeoscienceapplications[j].geocartointernational,1990,5(3):38-38.】。
由此，本发明通过构造4*4的实数矩阵，通过极化合成的方法，使得模型可以计算各种极化的反射率，并通过计算每个分层内的入射过来的机会卫星信号的镜像散射和非镜像的漫散射，从而避免了对各个方向的漫散射的忽略，使计算结果更准确。
步骤s5：根据所述步骤s2中的每个分层的透过率、反射率和每两个相邻分层之间的衰减以及所述步骤s4中的双站雷达散射截面来获得单次反射的总接收机波形。
其中，所述总接收机波形为：
式中，τ为时间延迟，上标dir表示接收机能量中的直射部分，i表示分层中的第i层，n为层数，ρi为第i层的延迟，lca为gps的c/a码片长度，j为复数单位，
其中，第i层的振幅ui、相位φi和延迟ρi根据所述每个分层的透过率、反射率、每两个相邻分层之间的衰减以及双站雷达散射截面计算得到，具体的表达式参见文献【boniface,k.,etal.,comparisonofsnowdataassimilationsystemwithgpsreflectometrysnowdepthinthewesternunitedstates.hydrologicalprocesses,2015.29(10):p.2425-2437.】以及【cardellach,e.,etal.,characterizationofdry-snowsub-structureusinggnssreflectedsignals.remotesensingofenvironment,2012.124(0):p.122-134.】
所得到的总接收机波形包括各种极化(rr、rl、rv、rh)波形。
由此，利用上述仿真分析方法，可以给出当根区土壤水分发生变化时，预测出相应的接收机波形的变化，进而可以得到针对根区土壤水分的全极化的单次反射仿真模型，从而其前期的机理仿真分析模型为分析解释卫星观测数据、卫星数据的仿真模拟、卫星数据同化、发展地表参数的定量反演等提供关键的机理工具。
如图2所示，本发明提供一种机会信号反射遥感的全极化单次反射仿真系统，其包括土壤剖面介电常数计算模块1、层内反射率计算模块2、层内透过率计算模块3、层间衰减率计算模块4、全极化反射率计算模块5和单次反射能量计算模块6。
其中，土壤剖面介电常数计算模块1设置为将一待测土壤剖面划分为多个分层，并根据待测土壤剖面的输入参数建立待测土壤剖面的介电常数模型，得到待测土壤剖面的各个分层的介电常数；所述全极化单次反射仿真系统的输入参数包括：土壤质地(沙土含量、粘土含量)、土壤温度、土壤块密度、地表粗糙度(相关长度、均方根高度)、土壤含水量。这些参数都可以通过实际测量获得，或者利用现有文献的输入参数，或者使用者根据经验得到的数据。
层内反射率计算模块2、层内透过率计算模块3和层间衰减率计算模块4分别设置为根据土壤剖面介电常数计算模块1得到的各个分层的介电常数来获取每个分层的透过率、反射率和每两个相邻分层之间的衰减。
全极化反射率计算模块5设置为根据所述层内反射率计算模块2得到的反射率得到全极化的反射特性，其包括镜像反射模块51、非镜像反射模块52以及与所述镜像反射模块51和非镜像反射模块52均相连的极化合成模块53，镜像反射模块51和非镜像反射模块52分别设置为根据分层的粗糙程度构建光滑表面反射率矩阵和粗糙表面反射率矩阵，具体地，当分层为光滑表面时，镜像反射模块51设置为构建光滑表面反射率矩阵，当分层为粗糙表面时，非镜像反射模块52设置为构建粗糙表面反射率矩阵；极化合成模块53设置为根据所述光滑表面反射率矩阵和粗糙表面反射率矩阵，利用极化合成的方法计算各个分层在各个极化状态上的反射特性(即双站雷达散射截面brcs)。
单次反射能量计算模块6设置为根据所述层内反射率计算模块2、层内透过率计算模块3和层间衰减率计算模块4得到的每个分层的透过率、反射率和每两个相邻分层之间的衰减以及所述全极化反射率计算模块5得到的各个分层在各个极化状态上的反射特性来获得总接收机波形；
由此，本发明针对根区土壤水分首先构造了一个土壤剖面，利用混合介质介电常数模型计算这个分层土壤剖面的各个分层的介电常数；然后利用菲涅尔反射系数和高级积分方程模型来分别计算每个分层的镜像散射和非镜像的漫散射；利用极化合成的方法使得模型可以计算各种极化状态上的反射系数；最后将各种极化状态上的反射系数作为模型的一部分放在总接收机波形的公式中进行计算，可以给出当根区土壤水分发生变化时，预测出相应的接收机波形的变化，进而可以得到针对根区土壤水分的全极化的单次反射仿真模型。
以上所述的，仅为本发明的较佳实施例，并非用以限定本发明的范围，本发明的上述实施例还可以做出各种变化。凡是依据本发明申请的权利要求书及说明书内容所作的简单、等效变化与修饰，皆落入本发明专利的权利要求保护范围。本发明未详尽描述的均为常规技术内容。
1.一种机会信号反射遥感的全极化单次反射仿真方法，其特征在于，其工作于p波段，包括：
步骤s1：将一待测土壤剖面划分为多个分层，并根据待测土壤剖面的输入参数建立待测土壤剖面的介电常数模型，得到待测土壤剖面的各个分层的介电常数；
步骤s2：利用所述步骤s1中的待测土壤剖面的各个分层的介电常数来获取每个分层的透过率、反射率和每两个相邻分层之间的衰减；
步骤s3：根据各个分层的粗糙程度，分别构建光滑表面反射率矩阵和粗糙表面反射率矩阵；
步骤s4：根据所述步骤s3中的光滑表面反射率矩阵和粗糙表面反射率矩阵，利用极化合成的方法计算各个分层在各个极化状态上的反射特性，该反射特性为各个分层的双站雷达散射截面；
步骤s5：根据所述步骤s2中的每个分层的透过率、反射率和每两个相邻分层之间的衰减以及所述步骤s4中的反射特性来获得单次反射的总接收机波形。
2.根据权利要求1所述的机会信号反射遥感的全极化单次反射仿真方法，其特征在于，所述全极化单次反射仿真系统的输入参数包括：土壤质地、土壤温度、土壤块密度、地表粗糙度和土壤含水量。
3.根据权利要求1所述的机会信号反射遥感的全极化单次反射仿真方法，其特征在于，所述待测土壤剖面的介电常数模型是基于混合介质介电常数模型修正得到的，所述介电常数的实部采用线性校正来校正，所述介电常数的虚部ξ"通过修正有效电导率来修正。
4.根据权利要求3所述的机会信号反射遥感的全极化单次反射仿真方法，其特征在于，所述介电常数的实部的校正公式为：
ξ'＝1.15ξ"-0.68
式中，ξ'和ξ"分别是介电常数的实部和虚部；
有效电导率σeff的修正公式为：
σeff＝0.0467 0.2204ρb-0.4111*s 0.6614*c，
式中，σeff是有效电导率，ρb是土壤块密度，s和c分别是沙土含量和粘土含量。
5.根据权利要求1所述的机会信号反射遥感的全极化单次反射仿真方法，其特征在于，所述每个分层的透过率、反射率和每两个相邻分层之间的衰减为：
其中，下标i和j分别是彼此相邻的分层中的第i层和第j层，tv和th分别为第i层的垂直和水平极化的透过率，rv和rh分别是第i层的垂直和水平极化的反射率；θ为入射角，εi为第i层的介电常数；λ为波长，α为第i层和第j层之间的衰减系数。
6.根据权利要求1所述的机会信号反射遥感的全极化单次反射仿真方法，其特征在于，在所述步骤s3中，对于分层为光滑表面的情况，则构建光滑表面反射率矩阵，所述光滑表面反射率矩阵是通过假定每层只发生镜像反射，基于菲涅尔公式来构建的；针对分层为粗糙表面的情况，则构建粗糙表面反射率矩阵，所述粗糙表面反射率矩阵是通过建立高级积分方程模型计算粗糙表面的反射率而得到的。
7.根据权利要求1所述的机会信号反射遥感的全极化单次反射仿真方法，其特征在于，在所述步骤s3中，极化合成的计算公式为：
其中，σqp为各个分层的双站雷达散射截面，下标q和p分别为接受和发射的极化状态；ψ和τ分别为椭倾角和椭率角；yq、yp为修改的stokes矢量；q为旋转矩阵；m为光滑表面反射率矩阵或者粗糙表面反射率矩阵。
8.根据权利要求1所述的机会信号反射遥感的全极化单次反射仿真方法，其特征在于，在所述步骤s5中，所述总接收机波形为：
式中，τ为时间延迟，上标dir表示接收机能量中的直射部分，i表示分层中的第i层，n为层数，ρi为第i层的延迟，lc/a为gps的c/a码片长度，j为复数单位，
其中，第i层的振幅ui、相位φi和延迟ρi根据所述每个分层的透过率、反射率、每两个相邻分层之间的衰减以及双站雷达散射截面计算得到。
9.一种机会信号反射遥感的全极化单次反射仿真系统，其特征在于，其工作于p波段，包括：
土壤剖面介电常数计算模块(1)，设置为将一待测土壤剖面划分为多个分层，并根据待测土壤剖面的输入参数建立待测土壤剖面的介电常数模型，得到待测土壤剖面的各个分层的介电常数；
层内反射率计算模块(2)、层内透过率计算模块(3)和层间衰减率计算模块(4)，分别设置为所述各个分层的介电常数来获取每个分层的透过率、反射率和每两个相邻分层之间的衰减；
全极化反射率计算模块(5)，根据所述层内反射率计算模块(2)得到的反射率得到全极化的反射特性，其包括：镜像反射模块(51)、非镜像反射模块(52)以及与这两者均相连的极化合成模块(53)，镜像反射模块(51)和非镜像反射模块(52)分别设置为根据分层的粗糙程度构建光滑表面反射率矩阵和粗糙表面反射率矩阵；极化合成模块(53)设置为根据所述光滑表面反射率矩阵和粗糙表面反射率矩阵，利用极化合成的方法计算各个分层在各个极化状态上的反射特性，该反射特性为各个分层的双站雷达散射截面；
单次反射能量计算模块(6)设置为根据所述每个分层的透过率、反射率和每两个相邻分层之间的衰减以及所述各个分层在各个极化状态上的反射特性来获得总接收机波形。
10.根据权利要求9所述的机会信号反射遥感的全极化单次反射仿真方法，其特征在于，所述全极化单次反射仿真系统的输入参数包括：土壤质地、土壤温度、土壤块密度、地表粗糙度和土壤含水量。
技术总结