1.本发明涉及定向钻井技术领域,尤其涉及单向量航姿信息提取方法。
背景技术:2.在石油、天然气勘探开发工程中,使用定向钻井技术可以根据地质条件、油气藏位置等具体情况采用合理的钻井轨迹以使社会经济效益最大化。在定向钻井过程中,为监测与控制实际钻井轨迹,需要获取钻具及不同井深处井眼轨迹的航姿信息。获取航姿信息的一种常用方法是使用三轴加速度计及三轴磁传感器测量重力加速度向量及地磁向量在航姿测量单元坐标系中的坐标,并通过处理这些坐标值获取一组欧拉角,即航向角、俯仰角及横滚角,这些角度即为航姿信息参数,由航姿参数可确定定向钻井中方位角、井斜角及工具面角等工程参数。
3.在实际应用中,井下复杂的工作环境会对传感器测量信号造成严重干扰,或造成传感器故障,这些不利影响会导致在航姿信息参数计算过程中,部分测量分量信号不可用。针对这一问题,文献1(范光第,蒲文学,赵国山.磁力随钻测斜仪轴向磁干扰校正方法[j]. 石油钻探技术,2017,45(4):121-126.)利用向量模值不变条件估计不可用分量。文献2(基于非完整测量向量的航姿解算方法,发明专利zl202110525236.7)进一步引入点积约束,提供了一种在缺失一个或两个测量向量分量情况下的航姿求解方法。文献3(一种测量向量分量缺失情况下欧拉角优化方法,发明专利zl202111230235.6)使用遍历搜索寻优的方法确定分量缺失情况下的欧拉角最优解。上述方法均需要来自两个向量的测量值分量,当仅有单向量测量值时无法提取航姿信息。
[0004]
本发明提供一种仅使用单向量测量值提取航姿信息的方法。应当指出的是,现有方法可以只使用重力加速度向量确定俯仰角与横滚角,这是由于重力加速度垂直于水平面。但其它向量,如地磁向量,并不具备这一特殊性,本发明提供的方法可以使用一般的单向量测量值缩小部分欧拉角的取值范围,提供了一种利用向量测量值获取航姿信息的新技术方法,可有效提高向量航姿信息的利用率。同时,该法使用解析计算而非遍历搜索方式确定欧拉角范围,具有计算量小且计算精度高的优点。
技术实现要素:[0005]
本发明的目的是为了解决现有技术中问题,而提出单向量航姿信息提取方法。
[0006]
为了实现上述目的,以地磁向量为可用向量的情况为例,详细说明本发明的实施方案,需要强调的是,针对地磁向量的处理方法具有一般性,可推广至其它单向量情况,具体来说,单向量航姿信息提取方法包括以下步骤:s1、计算初始变换矩阵,设重力加速度与地磁向量在参考坐标系中的坐标分别为与gr与mr,两向量的测量值分别为gb与mb。向量的分量使用脚标x、y、z表示,如mb的分量形式可写为[m
bx
,m
by
,m
bz
]
t
,t表示矩阵转置。归一化向量使用脚标n表示,如mb的归一化向量为m
bn
。使用式(1)可计算得到
一个满足mr与mb坐标变换关系的变换矩阵,称之为初始变换矩阵cs。
[0007]
上式中:注意,由于gb不可用,在以上计算中g
bn
可取任意单位向量,其计算结果均可满足地磁向量坐标变换关系。
[0008]
s2、获取一般变换矩阵表达式,以初始变换矩阵表示的航姿状态为初始状态,以地磁向量为轴,使航姿测量单元做定轴转动。在这一转动过程中,地磁向量在三轴磁传感器坐标系中的坐标保持不变,即转动过程中的任一航姿状态均满足地磁向量坐标变换关系。因此,一般变换矩阵cb的表达式可写为:上式中φ为定轴旋转角度,其取值范围为0至360度,[m
bn
×
]为m
bn
的反对称矩阵,其表达式为:注意,一般变换矩阵表达式已将所有满足地磁向量坐标变换关系的变换矩阵表示为旋转角度φ的函数。
[0009]
s3、确定俯仰角的取值范围,使用θ表示俯仰角,使用脚标i、j表示矩阵元素,如c
bij
表示变换矩阵cb的第i行第j列元素。俯仰角的计算式为:上式中在俯仰角完整取值范围内反正弦函数是单调函数,因此俯仰角的极值点与c
b23
的极值点相同。设俯仰角的极值点为φ
θ
,将式(7)对φ求导,并令导数等于零可得极值点计算式:
反正切函数的值域为-90度至90度,因此在φ的取值范围内可得到两个极值点φ
θ1
与 φ
θ2
,其计算式为:将极值点值代入式(6)可求得俯仰角的各个极值。设俯仰角最小值与最大值分别为θ
min
与θ
max
,则俯仰角的取值范围为[θ
min
,θ
max
],其计算式为:s4、确定航向角的取值范围,使用ψ表示航向角,航向角的计算式为:上式中:反正切函数是单调函数,因此航向角与-c
b21
/c
b22
具有相同的极值点。将-c
b21
/c
b22
对φ求导,并令导数等于零可得极值点φ
ψ1
与φ
ψ2
的计算式:上式中:
由于反正弦函数的值域为-90度至90度,因此在φ的取值范围内还可获得另外两个极值点φ
ψ3
与φ
ψ4
:将极值点值代入式(11)可求得航向角的各个极值。设航向角最小值与最大值分别为ψ
min
与ψ
max
,则航向角的取值范围为[ψ
min
,ψ
max
],其计算式为:当a
2 + c
2 ‑ꢀ
d2《 0时,式(14)无法求得实数极值点,说明航向角是单调变化的,考虑到当φ取0度与360度时,航向角的计算值是相同的,因此航向角将在完整范围内取值,即s5、确定横滚角的取值范围,使用γ表示横滚角,横滚角的计算式为:上式中:
反正切函数是单调函数,因此横滚角与-c
b13
/c
b33
具有相同的极值点。将-c
b13
/c
b33
对φ求导,并令导数等于零可得极值点φ
γ1
与φ
γ2
的计算式:上式中:由于反正弦函数的值域为-90度至90度,因此在φ的取值范围内还可获得另外两个极值点φ
γ3
与φ
γ4
:将极值点值代入式(21)可求得横滚角的各个极值。设横滚角最小值与最大值分别为γ
min
与γ
max
,则横滚角的取值范围为[γ
min
,γ
max
],其计算式为:
当e
2 + f
2 ‑ꢀ
h2《 0时,式(24)无法求得实数极值点,说明横滚角是单调变化的,考虑到当φ取0度与360度时,横滚角的计算值是相同的,因此横滚角将在完整范围内取值,即与现有技术相比,本发明的有益效果是:1、传统航姿参数估计方法需要至少两个向量的完整测量值,但干扰或传感器故障会造成向量测量值的部分分量不可用,从而导致传统方法不可用。针对这一问题,一些方法使用非完整测量值实现了航姿计算,但这些方法均要求测量值需包含两个向量的分量。本发明仅使用单一向量测量值提取航姿信息,提供了一种利用向量测量值获取航姿信息的新技术方法,可有效提高向量航姿信息的利用率。
[0010]
2、传统遍历搜索式方法按步长递增在转动角取值范围内确定最值,计算量与计算精度受步长影响,步长越小,计算量越大但计算精度越高;步长越大计算量越小但计算精度越低。本发明依据欧拉角计算函数给出了确定全部欧拉角极值点的解析计算式,可以在确保计算准确的同时明显降低计算量,具有计算成本低且计算精度高的优势。
附图说明
[0011]
图1为航姿测量示意图;图2为单向量航姿信息提取方法示意图;图3为航向角变化曲线;图4为俯仰角变化曲线;图5为横滚角变换曲线。
具体实施方式
[0012]
下面将结合本发明实施例中的技术方案进行清楚、完整地描述,显然,所描述的实施例仅仅是本发明一部分实施例,而不是全部的实施例。
[0013]
目前,在航姿测量中广泛使用的双向量方式如附图1所示,g为重力加速度向量,m为地磁向量,r-xyz坐标系为参考坐标系,在航姿测量单元中,b
1-xyz与b
2-xyz分别为三轴加速度计传感器坐标系及三轴磁传感器坐标系,且均与航姿测量单元坐标系平行,因此,可认为上述两种传感器的测量结果分别为重力加速度向量与地磁向量在航姿测量单元中的坐标。
[0014]
使用上述双向量测量值可以计算航姿参数,但在实际应用中外界干扰或传感器故障会导致某个向量测量值不可用,即仅有单个向量可用,现有方法无法在该条件下使用。需要指出的是,使用单个向量无法唯一确定航姿参数值,但单个向量测量值中仍含有部分航姿信息,本发明提供的方法可提取单个向量测量值中包含的航姿信息,信息提取结果的表现形式为缩小部分欧拉角的取值范围。下文将以地磁向量为可用向量的情况为例,详细说明该方法的实施方案,需要强调的是,下文具体针对地磁向量的处理方法具有一般性,可推广至其它单向量情况,单向量航姿信息提取方法如附图2所示:
s1、计算初始变换矩阵,设重力加速度与地磁向量在参考坐标系中的坐标分别为与gr与mr,两向量的测量值分别为gb与mb。向量的分量使用脚标x、y、z表示,如mb的分量形式可写为[m
bx
,m
by
,m
bz
]
t
,t表示矩阵转置。归一化向量使用脚标n表示,如mb的归一化向量为m
bn
。使用式(1)可计算得到一个满足mr与mb坐标变换关系的变换矩阵,称之为初始变换矩阵cs:上式中:注意,由于gb不可用,在以上计算中g
bn
可取任意单位向量,其计算结果均可满足地磁向量坐标变换关系。
[0015]
s2、获取一般变换矩阵表达式,以初始变换矩阵表示的航姿状态为初始状态,以地磁向量为轴,使航姿测量单元做定轴转动。在这一转动过程中,地磁向量在三轴磁传感器坐标系中的坐标保持不变,即转动过程中的任一航姿状态均满足地磁向量坐标变换关系。因此,一般变换矩阵cb的表达式可写为:上式中φ为定轴旋转角度,其取值范围为0至360度,[m
bn
×
]为m
bn
的反对称矩阵,其表达式为:注意,一般变换矩阵表达式已将所有满足地磁向量坐标变换关系的变换矩阵表示为旋转角度φ的函数。
[0016]
s3、确定俯仰角的取值范围,使用θ表示俯仰角,使用脚标i、j表示矩阵元素,如c
bij
表示变换矩阵cb的第i行第j列元素。俯仰角的计算式为:上式中
在俯仰角完整取值范围内反正弦函数是单调函数,因此俯仰角的极值点与c
b23
的极值点相同。设俯仰角的极值点为φ
θ
,将式(7)对φ求导,并令导数等于零可得极值点计算式:反正切函数的值域为-90度至90度,因此在φ的取值范围内可得到两个极值点φ
θ1
与 φ
θ2
,其计算式为:将极值点值代入式(6)可求得俯仰角的各个极值。设俯仰角最小值与最大值分别为θ
min
与θ
max
,则俯仰角的取值范围为[θ
min
,θ
max
],其计算式为:s4、确定航向角的取值范围,使用ψ表示航向角,航向角的计算式为:上式中:反正切函数是单调函数,因此航向角与-c
b21
/c
b22
具有相同的极值点。将-c
b21
/c
b22
对φ求导,并令导数等于零可得极值点φ
ψ1
与φ
ψ2
的计算式:上式中:由于反正弦函数的值域为-90度至90度,因此在φ的取值范围内还可获得另外两个极值点φ
ψ3
与φ
ψ4
:将极值点值代入式(11)可求得航向角的各个极值。设航向角最小值与最大值分别为ψ
min
与ψ
max
,则航向角的取值范围为[ψ
min
,ψ
max
],其计算式为:当a
2 + c
2 ‑ꢀ
d2《 0时,式(14)无法求得实数极值点,说明航向角是单调变化的,考虑到当φ取0度与360度时,航向角的计算值是相同的,因此航向角将在完整范围内取值,即s5、确定横滚角的取值范围,使用γ表示横滚角,横滚角的计算式为:上式中:
反正切函数是单调函数,因此横滚角与-c
b13
/c
b33
具有相同的极值点。将-c
b13
/c
b33
对φ求导,并令导数等于零可得极值点φ
γ1
与φ
γ2
的计算式:上式中:由于反正弦函数的值域为-90度至90度,因此在φ的取值范围内还可获得另外两个极值点φ
γ3
与φ
γ4
:将极值点值代入式(21)可求得横滚角的各个极值。设横滚角最小值与最大值分别为γ
min
与γ
max
,则横滚角的取值范围为[γ
min
,γ
max
],其计算式为:
ꢀ
当e
2 + f
2 ‑ꢀ
h2《 0时,式(24)无法求得实数极值点,说明横滚角是单调变化的,考虑到当φ取0度与360度时,横滚角的计算值是相同的,因此横滚角将在完整范围内取值,即最后需要指出的是,s3至s5分别提取了三个欧拉角的取值范围信息,这三步之间没有依赖关系,即可以按任意顺序计算或只计算其中部分欧拉角的取值范围。
[0017]
为验证上述方法的效果,给出1个计算实例,该实例中首先使用重力加速度向量与地磁向量测量值确定欧拉角的参考值,然后使用本发明方法提取地磁向量包含的航姿信息,最后使用遍历搜索法确定欧拉角最小值与最大值点,以进一步验证本发明方法的准确性。具体过程如下:实施例1在东北天参考坐标系下gr=[0, 0, 9.8]
t
,mr=[-3627, 29344, 44065]
t
,重力加速度测量值gb=[-5.39,-8.01, 0.89]
t
,地磁向量测量值mb=[913, 53204,
ꢀ‑
406]
t
,其中,重力加速度向量单位为米/秒2,地磁向量单位为纳特。使用以上双向量测量值及传统triad方法可得航向角参考值为:6.55度,俯仰角参考值为:-55.41度,横滚角参考值为:80.65度,航向角的完整取值范围为[-180,180]度,俯仰角的完整取值范围为[-90,90]度,横滚角的完整取值范围为[-180,180]度,以下使用本发明方法提取地磁向量航姿信息。
[0018]
第一步:设g
bn
=[0.58,0.58, 0.58]
t
,使用式(1)可计算初始变换矩阵:第二步:将归一化的地磁向量测量值与第一步计算结果代入式(4),可得一般变换矩阵表达式:第三步:将归一化的地磁向量测量值与第一步计算结果代入式(8),可得φ
θ
=69.69度,进而由式(9)可得极值点φ
θ1
=69.69度,φ
θ2
=249.69度。将极值点代入式(10),可得θ
min
=-57.21度,θ
max
=-55.06度,即俯仰角的取值范围为[-57.21,-55.06]度。
[0019]
第四步:将归一化的地磁向量测量值与第一步计算结果代入式(14),可得极值点φ
ψ1
=-21.91度,φ
ψ2
=-18.71度,进而由式(18)可得另外两个极值点φ
ψ3
=158.09度,φ
ψ4
=198.71度。将极值点代入式(19),可得ψ
min
=5.12度,ψ
max
=8.98度,即航向角的取值范围为[5.12,8.98]度。
[0020]
第五步:将归一化的地磁向量测量值与第一步计算结果代入式(24),可得e2=9.089
×
10-6
,f2=6.636
×
10-5
,h2=0.0963。显然e
2 + f
2 ‑ꢀ
h2《 0,由式(30)可得γ
min
=-180度,
γ
max
=180度,即横滚角的取值范围为[-180,180]度。
[0021]
由以上结果可知,在实例中本发明方法仅使用地磁向量测量值显著缩小了航向角与俯仰角的取值范围,且其确定的取值范围与相关角度的参考值相符。为进一步说明本方法确定欧拉角取值范围的准确性,使用遍历搜索法将φ由0度至360度按0.01度步长递增,计算相应的变换矩阵cb,并将结果转换为欧拉角,可得欧拉角变化曲线,结果如附图3、附图4、附图5所示。附图3中最小值点的坐标为(158.09,5.12),最大值点的坐标为(341.29,8.98)。附图4中最小值点的坐标为(69.69,-57.21),最大值点的坐标为(249.69,-55.06),附图5中最小值点的坐标为(360,-180),最大值点的坐标为(0,180)。对比各欧拉角最值点的纵坐标与本发明方法的计算结果可知,本方法确定的欧拉角取值范围是准确的。
[0022]
本发明提供了一种利用向量测量值获取航姿信息的新技术方案,仅使用单向量测量值构建航向角、俯仰角及横滚角的计算式,当已知该单向量测量值时,各欧拉角均表示为定轴转动角度值的一元函数。同时,利用反三角函数在特定区间上的单调性将求欧拉角极值问题转变为求一般变换矩阵元素函数的极值问题。
[0023]
使用解析计算而非遍历搜索方式获取欧拉角极值点,本发明给出了全部欧拉角极值点的解析计算式,并利用三角函数的周期性确定在定轴转动角度完整取值范围内的极值点。此外,根据极值点是否存在实数解判断欧拉角变化的单调性,在单调变化情况下直接给出欧拉角取值范围。
[0024]
以上所述,仅为本发明较佳的具体实施方式,但本发明的保护范围并不局限于此,任何熟悉本技术领域的技术人员在本发明揭露的技术范围内,根据本发明的技术方案及其发明构思加以等同替换或改变,都应涵盖在本发明的保护范围之内。