JP6738876B2 - Correction parameter creating device, correction parameter creating method, and correction parameter creating program - Google Patents

Description

本開示は、測定位置を地殻変動に基づいて補正するための補正パラメータを作成する補正パラメータ作成装置、方法およびプログラムに関する。 The present disclosure relates to a correction parameter creating device, method, and program for creating a correction parameter for correcting a measurement position based on crustal movement.

特許文献１および特許文献２には、全地球航法衛星システムを用いて位置情報を取得し、取得した位置情報が示す位置を地図上に表示する技術が知られている。
衛星測位分野では、大きく分けて、基準点の不要な単独測位と、基準点が必要な相対測位とが行われおり、一般的には後者の方が高精度である。実用化されている高精度な測位方法としては、相対測位であるＲＴＫ法があり、測量分野をはじめとして様々な分野で利用されている。ＲＴＫは、Real-Time Kinematicの略である。 Patent Documents 1 and 2 disclose a technique of acquiring position information using a global navigation satellite system and displaying the position indicated by the acquired position information on a map.
In the field of satellite positioning, there are broadly divided into single positioning that does not require a reference point and relative positioning that requires a reference point, and the latter is generally more accurate. As a highly accurate positioning method that has been put into practical use, there is an RTK method that is a relative positioning, which is used in various fields including the field of surveying. RTK is an abbreviation for Real-Time Kinematic.

しかし近年では、ＰＰＰ、ＰＰＰ−ＲＴＫ、ＰＰＰ−ＡＲなどの高精度な単独測位が発展してきている。ＰＰＰは、Precise Point Positioningの略である。ＡＲは、Ambiguity Resolutionの略である。 However, in recent years, highly accurate independent positioning such as PPP, PPP-RTK, and PPP-AR has been developed. PPP is an abbreviation for Precise Point Positioning. AR is an abbreviation for Ambiguity Resolution.

日本においては、２０１８年度より実運用が開始される準天頂衛星の補強信号として採用されている、センチメータ級測位補強サービス（ＣＬＡＳ）と、ＭＡＤＯＣＡとがその一例である。ＣＬＡＳは、Centimeter-Level Augmentation Serviceの略である。ＭＡＤＯＣＡは、Multi-GNSS Advanced Demonstration tool for Orbit and Clock Analysisの略である。高精度単独測位は、基準局を必要としないため、自動走行などの先端分野から期待が高まっている。なお、基準局が必要な場合には、基準局までの距離（すなわち、基線長）、および、基準局の切り替えなどの問題がある。 In Japan, an example is the Centimeter-class Positioning Augmentation Service (CLAS), which is adopted as an augmentation signal for the Quasi-Zenith Satellite, which will be put into practical use in FY2018, and MADOCA. CLAS is an abbreviation for Centimeter-Level Augmentation Service. MADOCA is an abbreviation for Multi-GNSS Advanced Demonstration tool for Orbit and Clock Analysis. High-accuracy independent positioning does not require a reference station, so expectations are rising from advanced fields such as autonomous driving. When a reference station is required, there are problems such as the distance to the reference station (that is, the base line length) and switching of the reference station.

地図の座標値は、測量成果に基づいており、時間情報である測量成果の元期に紐づけられている。このような座標を元期座標という。一方、単独測位によって得られる座標値は、観測日に紐づけられている。このような座標を今期座標という。一般的に、元期座標と今期座標とで、紐づけられた時間は異なる。つまり、時間経過に伴う影響、特に地殻変動の影響を大きく受ける。よって、単独測位によって得られた座標値と地図の座標値とは、時間経過とともに乖離するという問題が生じる。 The map coordinate values are based on the survey results, and are associated with the epoch of the survey results, which is time information. Such coordinates are called primal coordinates. On the other hand, the coordinate values obtained by independent positioning are associated with the observation date. Such coordinates are called current term coordinates. In general, the associated time differs between the original period coordinate and the current period coordinate. In other words, it is greatly affected by the passage of time, especially the effect of crustal deformation. Therefore, there arises a problem that the coordinate values obtained by the independent positioning and the coordinate values of the map deviate with the passage of time.

上記の問題の解決方法として、元期座標から今期座標の間に生じる地殻変動による移動量を補正する方法が考えられる。実例として、日本においては、国土地理院が公開している地殻変動補正パラメータ（以下、ＳＤパラメータ）を使用して地殻変動量を補正している。この補正方法をセミ・ダイナミック補正といい、地殻変動による位置の乖離を解消し、地図と整合性のとれた座標値を得ることができる。 As a solution to the above problem, a method of correcting the movement amount due to the crustal movement that occurs between the original period coordinate and the current period coordinate can be considered. As an example, in Japan, the amount of crustal movement is corrected using the crustal movement correction parameter (hereinafter referred to as SD parameter) published by the Geographical Survey Institute. This correction method is called semi-dynamic correction, and it is possible to eliminate the positional deviation due to crustal movement and obtain coordinate values that are consistent with the map.

特開２００１−０５９７２６号公報JP 2001-059726 A 特開２００４−１４４７１１号公報JP, 2004-144711, A

しかし、ＳＤパラメータの更新頻度は１年である。そして、地殻変動が活発な場所では年間約１０ｃｍも移動するため、高精度測位の精度よりも、１年間での地殻変動量の方が大きくなってしまう。また、ＳＤパラメータは毎年１月１日の位置情報を利用しているため、地殻変動補正の効果が時間経過とともに単調に減少する。
また、作成時点の基準点の位置から計算したパラメータは提供時には既に劣化し、使用時に最適に補正されない。また、基準点座標を連続観測していると、なかには地殻変動に依らないノイズ（システム雑音や観測雑音）を含む座標値が発生するので、それがパラメータに影響することがある。 However, the update frequency of the SD parameter is one year. Further, in a place where crustal movements are active, the amount of crustal movements in one year is larger than the accuracy of high-precision positioning because the distance moves about 10 cm per year. Further, since the SD parameter uses the position information on January 1st every year, the effect of the crustal deformation correction decreases monotonically with the passage of time.
Also, the parameters calculated from the position of the reference point at the time of creation are already deteriorated at the time of provision and are not optimally corrected at the time of use. Further, when the coordinates of the reference point are continuously observed, a coordinate value including noise (system noise or observation noise) that does not depend on the crustal movement is generated, which may affect the parameter.

本開示は、地殻変動補正の効果を向上させることを目的とする。 The present disclosure aims to improve the effect of crustal deformation correction.

本開示の一態様は、単独測位により測定された位置を地殻変動に基づいて補正するための補正パラメータを作成する補正パラメータ作成装置である。
そして、本開示の補正パラメータ作成装置は、基準日設定部と、有効期間設定部と、使用日設定部と、基準日位置算出部と、使用日位置算出部と、変動量算出部と、推定量算出部とを備える。ここで、基準日とは、地図等の位置の基準になる日を意味する。パラメータ有効期間とは、補正パラメータが有効な期間を意味する。使用日とは、パラメータ有効期間内の任意の日を意味する。 One aspect of the present disclosure is a correction parameter creation device that creates a correction parameter for correcting a position measured by independent positioning based on crustal movement.
Then, the correction parameter creation device of the present disclosure includes a reference date setting unit, a valid period setting unit, a use date setting unit, a reference date position calculation unit, a use date position calculation unit, a variation amount calculation unit, and an estimation unit. And an amount calculation unit. Here, the reference date means a date that serves as a reference for a position on a map or the like. The parameter valid period means a period during which the correction parameter is valid. The date of use means any day within the valid period of the parameter.

基準日設定部は、基準日を設定するように構成される。有効期間設定部は、補正パラメータが有効な期間であるパラメータ有効期間を設定するように構成される。使用日設定部は、パラメータ有効期間内に、使用日を設定するように構成される。 The reference date setting unit is configured to set the reference date. The valid period setting unit is configured to set a parameter valid period in which the correction parameter is valid. The usage date setting unit is configured to set the usage date within the parameter valid period.

基準日位置算出部は、補正パラメータを作成する対象となる地域である対象地域内に存在する複数の基準点のそれぞれについて、基準日を含む一定期間内の基準点の座標値に基づき、基準日における基準点の位置を基準日位置として算出するように構成される。 The reference date position calculation unit, based on the coordinate value of the reference point within a certain period including the reference date, for each of the plurality of reference points existing in the target area that is the area for which the correction parameter is to be created, The position of the reference point in is calculated as the reference date position.

使用日位置算出部は、対象地域内に存在する複数の基準点のそれぞれについて、算出時点に得られる一定期間内の基準点の座標値に基づき、使用日における基準点の位置を使用日位置として算出するように構成される。 The date-of-use position calculation unit determines the position of the reference point on the date of use as the date-of-use position based on the coordinate values of the reference point within a certain period obtained at the time of calculation for each of the plurality of reference points existing in the target area. Is configured to calculate.

変動量算出部は、対象地域内に存在する複数の基準点のそれぞれについて、基準日位置算出部により算出された基準日位置と、使用日位置算出部により算出された使用日位置との差を、地殻変動量として算出するように構成される。 The fluctuation amount calculation unit calculates the difference between the reference date position calculated by the reference date position calculation unit and the use date position calculated by the use date position calculation unit for each of the plurality of reference points existing in the target area. , Crustal movement amount is calculated.

推定量算出部は、対象地域を格子状に区切って形成された複数のメッシュのそれぞれの格子点にあたる地点について、空間補間法を用いて、その地点を囲む複数の基準点の地殻変動量に基づき、地殻変動推定量を算出し、算出した地殻変動推定量を補正パラメータとするように構成される。 The estimation amount calculation unit uses a spatial interpolation method for the points corresponding to the respective grid points of the plurality of meshes formed by dividing the target area into a grid pattern, and based on the crustal movements of the plurality of reference points surrounding the points. The crustal movement estimated amount is calculated, and the calculated crustal movement estimated amount is used as a correction parameter.

このように構成された本開示の補正パラメータ作成装置は、パラメータ有効期間内の使用日の地殻変動量に基づいて、補正パラメータを作成する。このため、本開示の補正パラメータ作成装置は、補正パラメータによる地殻変動補正の効果が、パラメータ有効期間内で最大となるようにすることができる。すなわち、本開示の補正パラメータ作成装置は、パラメータ有効期間の始期から使用日までは、地殻変動補正の効果が単調に増加し、使用日からパラメータ有効期間の終期までは、地殻変動補正の効果が単調に減少するようにすることが期待できる。これにより、本開示の補正パラメータ作成装置は、パラメータ有効期間全体として地殻変動補正の効果を向上させることができる。 The correction parameter creation device of the present disclosure configured as described above creates the correction parameter based on the amount of crustal movement on the day of use within the parameter valid period. Therefore, the correction parameter creation device according to the present disclosure can maximize the effect of the crustal movement correction by the correction parameter within the parameter valid period. That is, the correction parameter creation device of the present disclosure, the effect of crustal movement correction monotonically increases from the beginning of the parameter valid period to the use date, and the effect of crustal movement correction from the use date to the end of the parameter valid period. It can be expected to decrease monotonically. Accordingly, the correction parameter creation device of the present disclosure can improve the effect of crustal movement correction during the entire parameter valid period.

さらに、本開示の補正パラメータ作成装置は、一定期間内の基準点の座標値に基づき、基準日および使用日における基準点の位置（すなわち、基準日位置および使用日位置）を算出して、補正パラメータを作成する。これにより、本開示の補正パラメータ作成装置は、補正パラメータに対して、地殻変動によらないノイズによる基準点位置の変動の影響が及ぶのを抑制することができる。 Further, the correction parameter creation device of the present disclosure calculates the position of the reference point on the reference date and the use date (that is, the reference date position and the use date position) based on the coordinate value of the reference point within a certain period, and corrects the position. Create parameters. Accordingly, the correction parameter creation device of the present disclosure can suppress the influence of the fluctuation of the reference point position due to noise that is not due to crustal deformation on the correction parameter.

本開示の一態様では、使用日設定部は、パラメータ有効期間の中間に使用日を設定するようにしてもよい。これにより、本開示の補正パラメータ作成装置は、パラメータ有効期間内で最適に補正パラメータによる補正をさせることが期待できる。 In one aspect of the present disclosure, the usage date setting unit may set the usage date in the middle of the parameter valid period. As a result, the correction parameter creation device of the present disclosure can be expected to optimally perform correction using the correction parameter within the parameter valid period.

本開示の別の態様は、単独測位により測定された位置を地殻変動に基づいて補正するための補正パラメータを作成する補正パラメータ作成方法である。 Another aspect of the present disclosure is a correction parameter creation method that creates a correction parameter for correcting a position measured by independent positioning based on crustal movement.

そして、本開示の補正パラメータ作成方法は、基準日設定手順と、有効期間設定手順と、使用日設定手順と、基準日位置算出手順と、使用日位置算出手順と、変動量算出手順と、推定量算出手順とを備える。 Then, the correction parameter creation method of the present disclosure includes a reference date setting procedure, an effective period setting procedure, a use date setting procedure, a reference date position calculating procedure, a use date position calculating procedure, a variation amount calculating procedure, and an estimation. And an amount calculation procedure.

基準日設定手順は、基準日を設定する。有効期間設定手順は、補正パラメータが有効な期間であるパラメータ有効期間を設定する。使用日設定手順は、パラメータ有効期間内に、使用日を設定する。 The reference date is set in the reference date setting procedure. The valid period setting procedure sets a parameter valid period during which the correction parameter is valid. In the usage date setting procedure, the usage date is set within the parameter valid period.

基準日位置算出手順は、補正パラメータを作成する対象となる地域である対象地域内に存在する複数の基準点のそれぞれについて、基準日を含む一定期間内の基準点の座標値に基づき、基準日における基準点の位置を基準日位置として算出する。 The reference date position calculation procedure is based on the coordinate value of the reference point within a certain period including the reference date for each of the plurality of reference points existing in the target area, which is the area where the correction parameter is created. The position of the reference point in is calculated as the reference date position.

使用日位置算出手順は、対象地域内に存在する複数の基準点のそれぞれについて、算出時点に得られる一定期間内の基準点の座標値に基づき、使用日における基準点の位置を使用日位置として算出する。 The date-of-use position calculation procedure uses the position of the reference point on the day of use as the date-of-use position based on the coordinate values of the reference point within a certain period obtained at the time of calculation for each of the plurality of reference points existing in the target area. calculate.

変動量算出手順は、対象地域内に存在する複数の基準点のそれぞれについて、基準日位置算出手順により算出された基準日位置と、使用日位置算出手順により算出された使用日位置との差を、地殻変動量として算出する。 The variation calculation procedure calculates the difference between the reference date position calculated by the reference date position calculation procedure and the use date position calculated by the use date position calculation procedure for each of the plurality of reference points existing in the target area. , Calculated as the amount of crustal movement.

推定量算出手順は、対象地域を格子状に区切って形成された複数のメッシュのそれぞれの格子点にあたる地点について、空間補間法を用いて、その地点を囲む複数の基準点の地殻変動量に基づき、地殻変動推定量を算出し、算出した地殻変動推定量を補正パラメータとする。 The estimated amount calculation procedure is based on the amount of crustal movements of multiple reference points surrounding the points, using spatial interpolation method, for the points corresponding to the respective lattice points of the multiple meshes formed by dividing the target area into a grid pattern. A crustal deformation estimated amount is calculated, and the calculated crustal deformation estimated amount is used as a correction parameter.

本開示の補正パラメータ作成方法は、本開示の補正パラメータ作成装置にて実行される方法であり、当該方法を実行することで、本開示の補正パラメータ作成装置と同様の効果を得ることができる。 The correction parameter creation method of the present disclosure is a method executed by the correction parameter creation device of the present disclosure, and by executing the method, the same effect as that of the correction parameter creation device of the present disclosure can be obtained.

本開示の更に別の態様は、コンピュータを、基準日設定部、有効期間設定部、使用日設定部、基準日位置算出部、使用日位置算出部、変動量算出部、及び、推定量算出部として機能させるための補正パラメータ作成プログラムである。 Yet another aspect of the present disclosure is to use a computer as a reference date setting unit, a valid period setting unit, a use date setting unit, a reference date position calculation unit, a use date position calculation unit, a variation amount calculation unit, and an estimated amount calculation unit. It is a correction parameter creating program for functioning as.

本開示の補正パラメータ作成プログラムによって制御されるコンピュータは、本開示の補正パラメータ作成装置の一部を構成することができ、本開示の補正パラメータ作成装置と同様の効果を得ることができる。 The computer controlled by the correction parameter creation program of the present disclosure can configure a part of the correction parameter creation device of the present disclosure, and can obtain the same effect as the correction parameter creation device of the present disclosure.

補正パラメータ作成装置の構成を示すブロック図である。It is a block diagram which shows the structure of a correction parameter production apparatus. 第１実施形態の補正パラメータ作成処理を示すフローチャートである。It is a flow chart which shows correction parameter creation processing of a 1st embodiment. 有効期間中間日の算出方法を説明する図である。It is a figure explaining the calculation method of an effective period middle day. 第１実施形態においてパラメータ基準日での基準点の座標を算出する方法を説明する図である。It is a figure explaining the method of calculating the coordinates of the reference point on the parameter reference day in the first embodiment. 第１実施形態において有効期間中間日での基準点の座標を算出する方法を説明する図である。It is a figure explaining the method of calculating the coordinate of the reference point on the effective period middle day in 1st Embodiment. 地殻変動量ベクトルを算出する方法を説明する図である。It is a figure explaining the method of calculating a crustal movement amount vector. 経験バリオグラムを示す図である。It is a figure which shows an experience variogram. 第２実施形態の補正パラメータ作成処理を示すフローチャートである。It is a flow chart which shows correction parameter creation processing of a 2nd embodiment. ＳＤパラメータの線形補間方法を説明する図である。It is a figure explaining the linear interpolation method of SD parameter. 第２実施形態においてパラメータ基準日での基準点の座標を算出する方法を説明する図である。It is a figure explaining the method of calculating the coordinates of the reference point on the parameter reference day in the second embodiment. 第２実施形態において有効期間中間日での基準点の座標を算出する方法を説明する図である。It is a figure explaining the method of calculating the coordinate of the reference point on the effective period middle day in 2nd Embodiment. 有効期間中間日における地殻変動量を算出する方法を説明する図である。It is a figure explaining the method of calculating the amount of crustal movement on the middle day of the effective period.

（第１実施形態）
以下に本開示の第１実施形態を図面とともに説明する。
本実施形態の補正パラメータ作成装置１は、図１に示すように、表示部１１と、操作入力部１２と、データ記憶部１３と、データ入出力部１４と、制御部１５とを備える。 (First embodiment)
Hereinafter, a first embodiment of the present disclosure will be described with reference to the drawings.
As shown in FIG. 1, the correction parameter creating device 1 of the present embodiment includes a display unit 11, an operation input unit 12, a data storage unit 13, a data input/output unit 14, and a control unit 15.

表示部１１は、図示しない表示装置を備え、表示装置の表示画面に各種画像を表示する。
操作入力部１２は、図示しないキーボードおよびマウスを介して使用者が行った入力操作を特定するための入力操作情報を出力する。 The display unit 11 includes a display device (not shown) and displays various images on the display screen of the display device.
The operation input unit 12 outputs input operation information for specifying an input operation performed by the user via a keyboard and a mouse (not shown).

データ記憶部１３は、各種データを記憶するための記憶装置である。
データ入出力部１４は、有線または無線で接続された外部機器との間でデータの入出力を行う。 The data storage unit 13 is a storage device for storing various data.
The data input/output unit 14 inputs/outputs data to/from an external device connected by wire or wirelessly.

制御部１５は、ＣＰＵ、ＲＯＭおよびＲＡＭ等を備えたマイクロコンピュータを中心に構成される。マイクロコンピュータの各種機能は、ＣＰＵが非遷移的実体的記録媒体に格納されたプログラムを実行することにより実現される。この例では、ＲＯＭが、プログラムを格納した非遷移的実体的記録媒体に該当する。また、このプログラムの実行により、プログラムに対応する方法が実行される。なお、ＣＰＵが実行する機能の一部または全部を、一つあるいは複数のＩＣ等によりハードウェア的に構成してもよい。また、制御部１５を構成するマイクロコンピュータの数は１つでも複数でもよい。 The control unit 15 is mainly composed of a microcomputer including a CPU, a ROM, a RAM and the like. Various functions of the microcomputer are realized by the CPU executing the programs stored in the non-transitional substantive recording medium. In this example, the ROM corresponds to the non-transitional substantive recording medium storing the program. Further, by executing this program, the method corresponding to the program is executed. It should be noted that some or all of the functions executed by the CPU may be configured as hardware by one or a plurality of ICs. Further, the number of microcomputers forming the control unit 15 may be one or plural.

次に、制御部１５が実行する補正パラメータ作成処理の手順を説明する。補正パラメータ作成処理は、補正パラメータ作成処理を実行するために制御部１５に記憶された補正パラメータ作成プログラム２０を使用者の入力操作により起動することで実行される。なお、補正パラメータ作成プログラム２０は、補正パラメータ作成装置１に予めインストールされていてもよいし、記録媒体またはネットワークを介してインストールされるようにしてもよい。記録媒体としては、例えば光ディスク、磁気ディスクおよび半導体メモリなどが挙げられる。 Next, the procedure of the correction parameter creation process executed by the control unit 15 will be described. The correction parameter creating process is executed by activating the correction parameter creating program 20 stored in the control unit 15 to execute the correction parameter creating process by an input operation of the user. The correction parameter creating program 20 may be installed in the correction parameter creating apparatus 1 in advance, or may be installed via a recording medium or a network. Examples of the recording medium include an optical disc, a magnetic disc, and a semiconductor memory.

補正パラメータ作成処理が実行されると、制御部１５は、図２に示すように、まず、Ｓ１０にて、日付計算処理を実行する。日付計算処理では、制御部１５は、まず、有効期間開始日と有効期間日数を示す有効期間情報を入力するための画像を表示部１１の表示画面に表示する。そして、制御部１５は、有効期間情報が操作入力部１２から入力されるまで待機する。そして、有効期間情報が入力されると、制御部１５は、図３に示すように、式（１）により、有効期間中間日を算出する。式（１）におけるＤＭは有効期間中間日、ＤＢは有効期間開始日、ＤＵＲは有効期間日数である。また、式（１）における［＊］は、ガウス記号である。なお、ガウス記号内の＊は数字である。 When the correction parameter creation process is executed, the control unit 15 first executes a date calculation process in S10 as shown in FIG. In the date calculation process, the control unit 15 first displays on the display screen of the display unit 11 an image for inputting effective period start date and effective period information indicating the number of effective period days. Then, the control unit 15 waits until the valid period information is input from the operation input unit 12. Then, when the valid period information is input, the control unit 15 calculates the valid period intermediate date by the formula (1) as shown in FIG. In the formula (1), DM is the effective period intermediate date, DB is the effective period start date, and DUR is the effective period number of days. [*] in the formula (1) is a Gauss symbol. Note that * in the Gauss symbol is a number.

Ｓ１０の処理が終了すると、図２に示すように、Ｓ２０にて、制御部１５は、メッシュ取込処理を実行する。メッシュ取込処理では、制御部１５は、まず、補正パラメータを作成する対象となる地域（以下、対象地域）を示す対象地域情報を入力するための画像を表示部１１の表示画面に表示する。そして制御部１５は、対象地域情報が操作入力部１２から入力されるまで待機する。そして、対象地域情報が入力されると、制御部１５は、対象地域情報のメッシュデータをデータ記憶部１３から取得する。メッシュは、対象地域を緯度および経度に基づいて格子状に区切って形成される。そしてメッシュデータは、格子状に形成されるメッシュの格子点（以下、メッシュ頂点）における緯度および経度を示すデータである。なお、本実施形態では、メッシュデータは、補正パラメータ作成処理を開始する前に予めデータ記憶部１３に記憶される。

When the process of S10 is completed, as shown in FIG. 2, in S20, the control unit 15 executes a mesh loading process. In the mesh loading process, the control unit 15 first displays an image for inputting target area information indicating a target area (hereinafter, target area) for which correction parameters are created on the display screen of the display unit 11. Then, the control unit 15 waits until the target area information is input from the operation input unit 12. Then, when the target area information is input, the control unit 15 acquires the mesh data of the target area information from the data storage unit 13. The mesh is formed by dividing the target area into a grid shape based on latitude and longitude. The mesh data is data indicating latitude and longitude at grid points (hereinafter, mesh vertices) of the mesh formed in a grid pattern. In the present embodiment, the mesh data is stored in the data storage unit 13 in advance before starting the correction parameter creating process.

またＳ３０にて、制御部１５は、位置情報取得処理を実行する。位置情報取得処理では、制御部１５は、対象地域内に存在するＣＯＲＳの位置を示す基準点位置情報をデータ記憶部１３から取得する。ＣＯＲＳは、測量における基準点であり、測位衛星からの信号を受信することにより自身の位置が計測される。ＣＯＲＳは、Continuously Operating Reference Stationの略である。基準点位置情報は、少なくとも、ＣＯＲＳの緯度、経度、楕円体高、年月日および時刻を示す情報である。 Further, in S30, control unit 15 executes a position information acquisition process. In the position information acquisition process, the control unit 15 acquires the reference point position information indicating the position of the CORS existing in the target area from the data storage unit 13. CORS is a reference point in surveying, and its position is measured by receiving a signal from a positioning satellite. CORS is an abbreviation for Continuously Operating Reference Station. The reference point position information is information that indicates at least the latitude, longitude, ellipsoidal height, date, and time of CORS.

位置情報取得処理では、具体的に、制御部１５は、まず、後述するパラメータ基準日と、基準点位置情報の数（以下、位置取得数）とを示す取得情報を入力するための画像を表示部１１の表示画面に表示する。そして制御部１５は、取得情報が操作入力部１２から入力されるまで待機する。そして、取得情報が入力されると、制御部１５は、複数のＣＯＲＳ毎に、パラメータ基準日を中心とした位置取得数分の基準点位置情報を取得する。さらに制御部１５は、複数のＣＯＲＳ毎に、年月日および時刻が新しい順に位置取得数分の基準点位置情報を取得する。位置取得数は、例えば、１００〜２００である。 In the position information acquisition process, specifically, the control unit 15 first displays an image for inputting acquisition information indicating a parameter reference date described later and the number of reference point position information (hereinafter, the number of position acquisitions). It is displayed on the display screen of the section 11. Then, the control unit 15 waits until the acquired information is input from the operation input unit 12. Then, when the acquisition information is input, the control unit 15 acquires, for each of the plurality of CORS, reference point position information corresponding to the number of position acquisitions centered on the parameter reference date. Further, the control unit 15 acquires the reference point position information for each of the plurality of CORSs in the order of the latest date and time and the number of acquired positions. The number of acquired positions is, for example, 100 to 200.

次にＳ４０にて、制御部１５は、位置決定処理を実行する。位置決定処理では、制御部１５は、Ｓ３０で基準点位置情報を取得した複数のＣＯＲＳ毎に、パラメータ基準日におけるＣＯＲＳの位置と、有効期間中間日におけるＣＯＲＳの位置とを算出する。 Next, in S40, the control unit 15 executes a position determination process. In the position determination process, the control unit 15 calculates the position of the CORS on the parameter reference date and the position of the CORS on the intermediate date of the valid period for each of the plurality of CORSs for which the reference point position information is acquired in S30.

位置決定処理では、具体的には、制御部１５は、まず、複数のＣＯＲＳのうち、パラメータ基準日および有効期間中間日における位置が算出されていない１つのＣＯＲＳを選択する。そして制御部１５は、選択したＣＯＲＳにおいて、パラメータ基準日を中心とした位置取得数分の基準点位置情報の中から、最新の年月日および時刻の基準点位置情報（以下、第１最新基準点位置情報）と、最古の年月日および時刻の基準点位置情報（以下、第１最古基準点位置情報）とを抽出する。さらに制御部１５は、第１最古基準点位置情報における年月日および時刻を第１算出期間開始時ＤＵＰＢ１とし、第１最新基準点位置情報における年月日および時刻を第１算出期間終了時ＤＵＰＥ１とする。なお、第１算出期間開始時ＤＵＰＢ１から第１算出期間終了時ＤＵＰＥ１までの期間を第１算出期間ＤＵＰ１という。 In the position determination process, specifically, the control unit 15 first selects one of the plurality of CORSs for which the position on the parameter reference date and the intermediate date of the effective period has not been calculated. Then, in the selected CORS, the control unit 15 selects the reference point position information of the latest date and time from the reference point position information of the number of position acquisitions centered on the parameter reference date (hereinafter referred to as the first latest reference). (Point position information) and the reference point position information of the oldest date and time (hereinafter, the first oldest reference point position information) are extracted. Further, the control unit 15 sets the date and time in the first oldest reference point position information to DUPB1 at the start of the first calculation period, and sets the date and time in the first latest reference point position information at the end of the first calculation period. DUPE1. The period from the start of the first calculation period DUPB1 to the end of the first calculation period DUPE1 is called the first calculation period DUP1.

そして制御部１５は、第１算出期間ＤＵＰ１内の基準点位置情報を用い、緯度、経度および楕円体高のそれぞれについて線形回帰直線を算出する。さらに制御部１５は、緯度、経度および楕円体高のそれぞれで算出された線形回帰直線について、パラメータ基準日における緯度、経度および楕円体高を算出する。図４に示すように、点Ｐ１，Ｐ２，Ｐ３，Ｐ４，Ｐ５，Ｐ６，Ｐ７，Ｐ８は、第１算出期間ＤＵＰ１内の基準点位置情報における緯度、経度または楕円体高を示す。直線Ｌ１は、線形回帰直線を示す。点Ｐ９は、線形回帰直線に基づいて算出されたパラメータ基準日ＤＲにおける緯度、経度または楕円体高である。なお、線形回帰直線の代わりに、状態空間推定やディープラーニングを利用したリカレントネットワーク等の時系列推定を利用してもよい。 Then, the control unit 15 uses the reference point position information within the first calculation period DUP1 to calculate a linear regression line for each of latitude, longitude, and ellipsoidal height. Further, the control unit 15 calculates the latitude, the longitude and the ellipsoidal height on the parameter reference date for the linear regression lines calculated for the latitude, the longitude and the ellipsoidal height. As shown in FIG. 4, points P1, P2, P3, P4, P5, P6, P7, and P8 indicate the latitude, longitude, or ellipsoidal height in the reference point position information within the first calculation period DUP1. The straight line L1 shows a linear regression line. The point P9 is the latitude, longitude, or ellipsoidal height on the parameter reference day DR calculated based on the linear regression line. Instead of the linear regression line, state space estimation or time series estimation such as recurrent network using deep learning may be used.

次に制御部１５は、選択したＣＯＲＳにおいて、年月日および時刻が新しい順に取得した位置取得数分の基準点位置情報の中から、最新の年月日および時刻の基準点位置情報（以下、第２最新基準点位置情報）と、最古の年月日および時刻の基準点位置情報（以下、第２最古基準点位置情報）とを抽出する。さらに制御部１５は、第２最古基準点位置情報における年月日および時刻を第２算出期間開始時ＤＵＰＢ２とし、第２最新基準点位置情報における年月日および時刻を第２算出期間終了時ＤＵＰＥ２とする。なお、第２算出期間開始時ＤＵＰＢ２から第２算出期間終了時ＤＵＰＥ２までの期間を第２算出期間ＤＵＰ２という。 Next, in the selected CORS, the control unit 15 selects the reference point position information of the latest date and time from the reference point position information corresponding to the number of position acquisitions acquired in the order of latest date and time (hereinafter, reference point position information (hereinafter, referred to as “reference point position information”). The second latest reference point position information) and the reference point position information of the oldest date and time (hereinafter, second oldest reference point position information) are extracted. Further, the control unit 15 sets the date and time in the second oldest reference point position information as DUPB2 at the start of the second calculation period, and sets the date and time in the second latest reference point position information at the end of the second calculation period. DUPE2. The period from the start of the second calculation period DUPB2 to the end of the second calculation period DUPE2 is referred to as the second calculation period DUP2.

そして制御部１５は、第２算出期間ＤＵＰ２内の基準点位置情報を用い、緯度、経度および楕円体高のそれぞれについて線形回帰直線を算出する。そして、緯度、経度および楕円体高のそれぞれで算出された線形回帰直線について、有効期間中間日における緯度、経度および楕円体高を算出する。図５に示すように、点Ｐ１１，Ｐ１２，Ｐ１３，Ｐ１４，Ｐ１５は、第２算出期間ＤＵＰ２内の基準点位置情報における緯度、経度または楕円体高を示す。直線Ｌ２は、線形回帰直線を示す。点Ｐ１６は、線形回帰直線に基づいて算出された有効期間中間日ＤＭにおける緯度、経度または楕円体高である。 Then, the control unit 15 uses the reference point position information within the second calculation period DUP2 to calculate a linear regression line for each of latitude, longitude, and ellipsoidal height. Then, the latitude, longitude, and ellipsoid height on the intermediate day of the effective period are calculated for the linear regression lines calculated for each of the latitude, longitude, and ellipsoid height. As shown in FIG. 5, points P11, P12, P13, P14, and P15 indicate the latitude, longitude, or ellipsoidal height in the reference point position information within the second calculation period DUP2. The straight line L2 shows a linear regression line. The point P16 is the latitude, longitude, or ellipsoidal height on the effective period middle day DM calculated based on the linear regression line.

そして制御部１５は、パラメータ基準日および有効期間中間日における位置が算出されていないＣＯＲＳが存在しているか否かを判断する。ここで、位置が算出されていないＣＯＲＳが存在している場合には、ＣＯＲＳを新たに選択して、上述の処理を繰り返す。一方、位置が算出されていないＣＯＲＳが存在しない場合には、位置決定処理を終了する。 Then, the control unit 15 determines whether or not there is a CORS whose position has not been calculated on the parameter reference date and the effective period intermediate date. Here, when there is a CORS whose position has not been calculated, the CORS is newly selected and the above processing is repeated. On the other hand, if there is no CORS whose position has not been calculated, the position determination process ends.

Ｓ４０の処理が終了すると、図２に示すように、制御部１５は、Ｓ５０にて、変動量算出処理を実行する。具体的には、制御部１５は、図６に示すように、Ｓ３０で基準点位置情報を取得した複数のＣＯＲＳ毎に、式（２）により、パラメータ基準日におけるＣＯＲＳの座標値ＣＲ（すなわち、緯度、経度および楕円体高）を始点とし、有効期間中間日におけるＣＯＲＳの座標値ＣＭを終点とする地殻変動量ベクトルＶＭを算出する。 When the processing of S40 ends, as shown in FIG. 2, the control unit 15 executes the fluctuation amount calculation processing in S50. Specifically, as shown in FIG. 6, the control unit 15 uses the equation (2) for each of the plurality of CORSs for which the reference point position information is acquired in S30, and then uses the coordinate value CR (that is, Latitude, longitude and ellipsoidal height) are the starting points, and the crustal movement amount vector VM having the coordinate value CM of the CORS on the middle day of the effective period as the ending point is calculated.

Ｓ５０の処理が終了すると、図２に示すように、制御部１５は、Ｓ６０にて、経験バリオグラム作成処理を実行する。

When the process of S50 ends, as shown in FIG. 2, the control unit 15 executes the experience variogram creation process in S60.

経験バリオグラム作成処理では、具体的に、制御部１５は、まず、ビン幅を示すビン幅情報を入力するための画像を表示部１１の表示画面に表示する。そして制御部１５は、ビン幅情報が操作入力部１２から入力されるまで待機する。そして、ビン幅情報が入力されると、制御部１５は、Ｓ３０で基準点位置情報を取得した複数のＣＯＲＳの中から、後述する非類似度が算出されていない２つのＣＯＲＳの組み合わせ（以下、ＣＯＲＳ対）を選択する。そして制御部１５は、選択したＣＯＲＳ対において、二点間の距離（以下、点間距離）を算出する。２つのＣＯＲＳ間の距離として、平面投影した座標値や地心座標系における直線距離、測地線長を利用してもよい。 In the experience variogram creation process, specifically, the control unit 15 first displays an image for inputting bin width information indicating the bin width on the display screen of the display unit 11. Then, the control unit 15 waits until the bin width information is input from the operation input unit 12. Then, when the bin width information is input, the control unit 15 selects a combination of two CORSs whose dissimilarities, which will be described later, have not been calculated from the plurality of CORSs for which the reference point position information has been acquired in S30 (hereinafter, CORS pair). Then, the control unit 15 calculates the distance between two points (hereinafter, the point distance) in the selected CORS pair. As the distance between the two CORS, a coordinate value projected on a plane, a straight line distance in a geocentric coordinate system, or a geodesic line length may be used.

そして制御部１５は、選択したＣＯＲＳ対における各成分（すなわち、緯度、経度および楕円体高）の非類似度ＤＩＳを、式（３）により算出する。式（３）におけるＳＭ１は、選択したＣＯＲＳ対の一方における地殻変動量ベクトルＶＭの各成分である。式（３）におけるＳＭ２は、選択したＣＯＲＳ対の他方における地殻変動量ベクトルＶＭの各成分である。 Then, the control unit 15 calculates the dissimilarity DIS of each component (that is, the latitude, the longitude, and the ellipsoidal height) in the selected CORS pair by the formula (3). SM1 in the equation (3) is each component of the crustal movement amount vector VM in one of the selected CORS pairs. SM2 in Formula (3) is each component of the crustal movement amount vector VM in the other of the selected CORS pairs.

そして制御部１５は、図７に示すように、地殻変動量ベクトルＶＭの各成分について、点間距離を横軸とし、非類似度ＤＩＳを縦軸とした二次元直交座標系に、算出した点間距離と非類似度ＤＩＳとを示す座標点を配置することにより、経験バリオグラムを作成する。図７は、点間距離を示す横軸をビン幅で分割したＢＮ１〜ＢＮ７内に座標点Ｐ２１〜Ｐ３２が配置された経験バリオグラムを示している。

Then, as shown in FIG. 7, the control unit 15 calculates, for each component of the crustal movement amount vector VM, a calculated point in a two-dimensional orthogonal coordinate system with the distance between points as the horizontal axis and the dissimilarity DIS as the vertical axis. An empirical variogram is created by arranging coordinate points indicating the inter-distance and the dissimilarity DIS. FIG. 7 shows an empirical variogram in which coordinate points P21 to P32 are arranged in BN1 to BN7 in which the horizontal axis indicating the distance between points is divided by the bin width.

次に制御部１５は、非類似度ＤＩＳが算出されていないＣＯＲＳ対が存在しているか否かを判断する。ここで、非類似度ＤＩＳが算出されていないＣＯＲＳ対が存在している場合には、ＣＯＲＳ対を新たに選択して、上述の処理を繰り返す。一方、非類似度ＤＩＳが算出されていないＣＯＲＳ対が存在しない場合には、経験バリオグラム作成処理を終了する。 Next, the control unit 15 determines whether or not there is a CORS pair for which the dissimilarity DIS has not been calculated. Here, if there is a CORS pair for which the dissimilarity DIS has not been calculated, a CORS pair is newly selected and the above processing is repeated. On the other hand, if there is no CORS pair for which the dissimilarity DIS has not been calculated, the empirical variogram creation processing ends.

Ｓ６０の処理が終了すると、図２に示すように、制御部１５は、Ｓ７０にて、モデル推定処理を実行する。モデル推定処理では、制御部１５は、経験バリオグラムから、式（４）で表されるバリオグラム関数γ（ｒ）を推定する。式（４）におけるＣ_３は、式（５）で表される。 When the process of S60 ends, as shown in FIG. 2, the control unit 15 executes the model estimation process in S70. In the model estimation process, the control unit 15 estimates the variogram function γ(r) represented by the equation (4) from the empirical variogram. _{C 3} in formula (4) is represented by the formula (5).

モデル推定処理では、具体的に、制御部１５は、まず、Ｓ６０で作成された経験バリオグラムを用いて、推定用データ列（ｘ_ｉ，ｙ_ｉ）を作成する。ｘ_ｉは、ビンｉの階級の中心である。ｙ_ｉは、ビンｉに含まれるＮｉ個の非類似度ＤＩＳの平均値である。 In the model estimation process, specifically, the control unit 15 first creates an estimation data string (x _i , y _i ) using the empirical variogram created in S60. x _i is the center of the bin of bin i. y _i is the average value of the Ni dissimilarities DIS contained in the bin i.

例えば、図７において、ビンＢＮ１，ＢＮ２，ＢＮ３，ＢＮ４，ＢＮ５，ＢＮ６，ＢＮ７はそれぞれ、ビン１，２，３，４，５，６，７に相当する。そして、ビン１に相当するビンＢＮ１内には、座標点Ｐ２１に対応する１個の非類似度ＤＩＳが含まれる。また、ビン５に相当するビンＢＮ５内には、座標点Ｐ２７，Ｐ２８に対応する２個の非類似度ＤＩＳが含まれる。 For example, in FIG. 7, bins BN1, BN2, BN3, BN4, BN5, BN6, BN7 correspond to bins 1, 2, 3, 4, 5, 6, 7, respectively. Then, the bin BN1 corresponding to the bin 1 includes one dissimilarity DIS corresponding to the coordinate point P21. Further, the bin BN5 corresponding to the bin 5 includes two dissimilarities DIS corresponding to the coordinate points P27 and P28.

ｘ_ｉ，ｙ_ｉはそれぞれ、式（６），（７）で表される。式（６）におけるＢＷは、ビン幅である。 x _i and y _i are represented by equations (6) and (7), respectively. BW in equation (6) is the bin width.

そして制御部１５は、最小二乗法によって、作成された複数の推定用データ列（ｘ_ｉ，ｙ_ｉ）をバリオグラム関数γ（ｒ）で近似することにより、バリオグラム関数γ（ｒ）のパラメータａ，ｂを算出する。

Then, the control unit 15 approximates the created plurality of estimation data strings (x _i , y _i ) by the variogram function γ(r) by the least squares method, thereby obtaining the parameters a, Calculate b.

制御部１５は、上記のパラメータａ，ｂの算出を、地殻変動量ベクトルＶＭの各成分について実行する。
Ｓ７０の処理が終了すると、図２に示すように、制御部１５は、Ｓ８０にて、クリギング法計算処理を実行する。クリギング法計算処理では、具体的に、制御部１５は、まず、通常型クリギング法の係数行列Ａを算出する。係数行列Ａは、式（８）で表される。式（８）におけるｎは、Ｓ３０で基準点位置情報を取得したＣＯＲＳの数である。式（８）におけるｒ_ｉ，ｊは、ｉで特定されるＣＯＲＳとｊで特定されるＣＯＲＳとの点間距離である。ｉ，ｊは、１〜ｎの整数である。クリギング法は、空間補間法の一種であり、国土地理院が提供しているＳＤパラメータの算出にも使用されている。 The control unit 15 executes the above calculation of the parameters a and b for each component of the crustal movement amount vector VM.
When the process of S70 ends, as shown in FIG. 2, the control unit 15 executes the Kriging method calculation process in S80. In the Kriging method calculation process, specifically, the control unit 15 first calculates the coefficient matrix A of the normal Kriging method. The coefficient matrix A is expressed by equation (8). N in Formula (8) is the number of CORS which acquired the reference point position information in S30. R _i,j in the equation (8) is a point-to-point distance between the CORS specified by i and the CORS specified by j. i and j are integers of 1 to n. The Kriging method is a kind of spatial interpolation method, and is also used for calculating SD parameters provided by GSI.

そして制御部１５は、係数行列ＡのＬＵ分解を行う。すなわち、制御部１５は、係数行列Ａを、下三角行列Ｌと、上三角行列Ｕとの積に分解する。なお、ＬＵ分解のアルゴリズムとして、Ｃｒｏｕｔ法およびＤｏｏｌｉｔｔｌｅ法が知られている。

Then, the control unit 15 performs LU decomposition of the coefficient matrix A. That is, the control unit 15 decomposes the coefficient matrix A into the product of the lower triangular matrix L and the upper triangular matrix U. The Crout method and the Doolittle method are known as LU decomposition algorithms.

制御部１５は、上記の係数行列Ａの算出と係数行列ＡのＬＵ分解とを、地殻変動量ベクトルＶＭの各成分について実行する。
Ｓ８０の処理が終了すると、制御部１５は、Ｓ９０にて、地殻変動量推定処理を実行する。 The control unit 15 executes the above calculation of the coefficient matrix A and LU decomposition of the coefficient matrix A for each component of the crustal movement amount vector VM.
When the process of S80 ends, the control unit 15 executes the crustal deformation amount estimation process in S90.

地殻変動量推定処理では、具体的に、制御部１５は、まず、Ｓ２０にてメッシュデータが取得された複数のメッシュ頂点の中から、後述する地殻変動推定量が算出されていないメッシュ頂点を一つ選択する。そして制御部１５は、選択したメッシュ頂点のメッシュデータ（すなわち、メッシュ頂点における緯度および経度）を取得する。次に制御部１５は、係数ベクトルｖを算出する。係数ベクトルｖは、ｖ_１，ｖ_２，ｖ_３，・・・・，ｖ_ｎ，ｖ_ｎ＋１を成分とするベクトルである。ｎは、Ｓ３０で基準点位置情報を取得したＣＯＲＳの数である。ｉ＝１，２，３，・・・，ｎである場合において、ｖ_ｉ＝γ（ｒ_ｉ）である。ｒ_ｉは、選択したメッシュ頂点と、ｉで特定されるＣＯＲＳとの間の距離である。また、ｖ_ｎ＋１＝１である。 In the crustal deformation amount estimation process, specifically, the control unit 15 first selects one of the mesh vertices for which the crustal deformation estimated amount, which will be described later, has not been calculated from the plurality of mesh vertices for which the mesh data was acquired in S20. Choose one. Then, the control unit 15 acquires the mesh data of the selected mesh vertex (that is, the latitude and longitude at the mesh vertex). Next, the control unit 15 calculates the coefficient vector v. The coefficient vector v is a vector whose components are v ₁ , v ₂ , v ₃ ,..., V _n , v _n+1 . n is the number of CORS whose reference point position information was acquired in S30. In the case where i=1, 2, 3,..., N, v _i =γ(r _i ). r _i is the distance between the selected mesh vertex and the CORS identified by i. Also, v _n+1 =1.

そして制御部１５は、変数ベクトルｗを用いて、式（９）に示す方程式を解く。変数ベクトルｗは、式（９）に示すように、ｗ_１，ｗ_２，ｗ_３，・・・・，ｗ_ｎ，λを成分とするベクトルである。λは、ラグランジュの未定乗数である。すなわち、制御部１５は、ラグランジュの未定乗数法により、式（９）を解く。なお、制御部１５は、式（９）を解く際に、係数行列ＡのＬＵ分解を行う。したがって、制御部１５は、Ａｗ＝（ＬＵ）ｗ＝ｖを解くことにより、ｗ_１，ｗ_２，ｗ_３，・・・・，ｗ_ｎ，λを算出する。そして制御部１５は、上記のｗ_１，ｗ_２，ｗ_３，・・・・，ｗ_ｎ，λの算出を各成分（すなわち、緯度、経度および楕円体高）について実行する。 Then, the control unit 15 uses the variable vector w to solve the equation shown in Expression (9). The variable vector w is a vector whose components are w ₁ , w ₂ , w ₃ ,..., W _n , λ, as shown in Expression (9). λ is a Lagrange's undetermined multiplier. That is, the control unit 15 solves the equation (9) by the Lagrange's undetermined multiplier method. The control unit 15 performs LU decomposition of the coefficient matrix A when solving the equation (9). Therefore, the control unit 15 calculates w ₁ , w ₂ , w ₃ ,..., W _n , λ by solving Aw=(LU)w=v. Then, the control unit 15 executes the above-described calculation of w ₁ , w ₂ , w ₃ ,..., W _n , λ for each component (that is, latitude, longitude, and ellipsoid height).

ｗ_１，ｗ_２，ｗ_３，・・・・，ｗ_ｎを算出した後に、制御部１５は、式（１０）により、選択したメッシュ頂点の地殻変動推定量ｚを算出する。式（１０）におけるｚ_ｉは、ｉで特定されるＣＯＲＳにおける各成分（すなわち、緯度、経度および楕円体高）の地殻変動量である。

After calculating w ₁ , w ₂ , w ₃ ,..., W _n , the control unit 15 calculates the crustal movement estimated amount z of the selected mesh vertex by the equation (10). Z _i in Expression (10) is the amount of crustal movement of each component (that is, latitude, longitude, and ellipsoidal height) in CORS specified by i.

さらに制御部１５は、式（１１）により、メッシュ頂点における地殻変動推定量zの予測誤差σを算出する。なお、式（１１）における「・」はベクトルの内積を表す。

Further, the control unit 15 calculates the prediction error σ of the crustal movement estimated amount z at the mesh vertices by the formula (11). In addition, “·” in Expression (11) represents an inner product of vectors.

そして制御部１５は、上記の地殻変動推定量ｚおよび予測誤差σの算出を、各成分（すなわち、緯度、経度および楕円体高）について実行する。

Then, the control unit 15 executes the calculation of the crustal movement estimated amount z and the prediction error σ for each component (that is, latitude, longitude, and ellipsoid height).

各成分の地殻変動推定量ｚおよび予測誤差σの算出が終了すると、制御部１５は、地殻変動推定量ｚおよび予測誤差σが算出されていないメッシュ頂点が存在しているか否かを判断する。ここで、地殻変動推定量ｚおよび予測誤差σが算出されていないメッシュ頂点が存在している場合には、メッシュ頂点を新たに選択して、上述の処理を繰り返す。一方、地殻変動推定量ｚおよび予測誤差σが算出されていないメッシュ頂点が存在しない場合には、地殻変動量推定処理を終了する。 When the calculation of the crustal movement estimated amount z and the prediction error σ of each component is completed, the control unit 15 determines whether there is a mesh vertex for which the crustal movement estimated amount z and the prediction error σ are not calculated. Here, when there is a mesh vertex for which the crustal movement estimated amount z and the prediction error σ have not been calculated, a mesh vertex is newly selected and the above process is repeated. On the other hand, if there are no mesh vertices for which the estimated crustal movement amount z and the prediction error σ have not been calculated, the crustal movement amount estimation processing ends.

Ｓ９０の処理が終了すると、制御部１５は、Ｓ１００にて、地殻変動推定量出力処理を実行する。地殻変動推定量出力処理では、具体的に、制御部１５は、Ｓ９０で算出した地殻変動推定量ｚを示す地殻変動推定量情報を、補正パラメータ作成装置１の外部へ出力する。そして、Ｓ１００の処理が終了すると、制御部１５は、Ｓ１１０にて、予測誤差出力処理を実行する。予測誤差出力処理では、具体的に、制御部１５は、Ｓ９０で算出した予測誤差σを示す予測誤差情報を、補正パラメータ作成装置１の外部へ出力する。 When the process of S90 ends, the control unit 15 executes the crustal deformation estimated amount output process in S100. In the crustal movement estimated amount output processing, specifically, the control unit 15 outputs the crustal movement estimated amount information indicating the crustal movement estimated amount z calculated in S90 to the outside of the correction parameter creation device 1. Then, when the process of S100 ends, the control unit 15 executes the prediction error output process in S110. In the prediction error output process, specifically, the control unit 15 outputs the prediction error information indicating the prediction error σ calculated in S90 to the outside of the correction parameter creation device 1.

そして、Ｓ１１０の処理が終了すると、制御部１５は、補正パラメータ作成処理を終了する。
このように構成された補正パラメータ作成装置１は、単独測位により測定された位置を地殻変動に基づいて補正するための補正パラメータを作成する。 Then, when the process of S110 ends, the control unit 15 ends the correction parameter creation process.
The correction parameter creating device 1 configured as described above creates a correction parameter for correcting the position measured by the single positioning based on the crustal movement.

そして補正パラメータ作成装置１は、パラメータ基準日を設定する。また補正パラメータ作成装置１は、有効期間を設定する。また補正パラメータ作成装置１は、有効期間中間日を設定する。ここで、パラメータ基準日とは、地図等の位置の基準になる日を意味する。有効期間とは、補正パラメータが有効な期間を意味する。 Then, the correction parameter creation device 1 sets the parameter reference date. Further, the correction parameter creating device 1 sets the valid period. Further, the correction parameter creation device 1 sets the effective period intermediate date. Here, the parameter reference date means a date that serves as a reference for a position on a map or the like. The valid period means a period during which the correction parameter is valid.

そして補正パラメータ作成装置１は、補正パラメータを作成する対象となる対象地域内に存在する複数のＣＯＲＳのそれぞれについて、基準日を含む一定期間内のＣＯＲＳの座標値に基づき、パラメータ基準日におけるＣＯＲＳの位置（以下、基準日位置）を算出する。 Then, the correction parameter creation device 1 determines the CORS of the parameter reference date based on the coordinate values of the CORS within a certain period including the reference date for each of the plurality of CORS existing in the target area for which the correction parameter is to be created. The position (hereinafter referred to as the reference date position) is calculated.

また補正パラメータ作成装置１は、対象地域内に存在する複数のＣＯＲＳのそれぞれについて、算出時点に得られる一定期間内のＣＯＲＳの座標値に基づき、有効期間中間日におけるＣＯＲＳの位置（以下、使用日位置）を算出する。 Further, the correction parameter creation device 1 determines the position of the CORS on the intermediate day of the effective period (hereinafter, used date) for each of the plurality of CORS existing in the target area based on the coordinate values of the CORS within a certain period obtained at the time of calculation. Position).

さらに補正パラメータ作成装置１は、対象地域内に存在する複数のＣＯＲＳのそれぞれについて、算出された基準日位置と、算出された使用日位置との差を、地殻変動量ベクトルＶＭとして算出する。 Further, the correction parameter creation device 1 calculates the difference between the calculated reference date position and the calculated use date position for each of the plurality of CORS existing in the target area as the crustal movement amount vector VM.

そして補正パラメータ作成装置１は、対象地域を格子状に区切って形成された複数のメッシュのそれぞれの格子点にあたる地点について、クリギング法を用いて、その地点を囲む複数のＣＯＲＳの地殻変動量ベクトルＶＭに基づき、地殻変動推定量ｚを算出し、算出した地殻変動推定量ｚを補正パラメータとする。 Then, the correction parameter creation device 1 uses the kriging method for the points corresponding to the respective grid points of the plurality of meshes formed by partitioning the target area into a grid pattern, and uses the kriging method for the plurality of CORS crustal movement vector VMs surrounding the points. The crustal movement estimated amount z is calculated based on the above, and the calculated crustal movement estimated amount z is used as a correction parameter.

このように補正パラメータ作成装置１は、有効期間内の使用日の地殻変動量に基づいて、補正パラメータを作成する。このため、補正パラメータ作成装置１は、補正パラメータによる地殻変動補正の効果が、有効期間内で最大となるようにすることができる。すなわち、補正パラメータ作成装置１は、有効期間の始期から有効期間中間日までは、地殻変動補正の効果が単調に増加し、有効期間中間日から有効期間の終期までは、地殻変動補正の効果が単調に減少するようにすることが期待できる。これにより、補正パラメータ作成装置１は、有効期間全体として地殻変動補正の効果を向上させることができる。 In this way, the correction parameter creation device 1 creates the correction parameter based on the amount of crustal movement on the day of use within the effective period. Therefore, the correction parameter creating device 1 can maximize the effect of the crustal deformation correction by the correction parameter within the effective period. That is, in the correction parameter creating device 1, the effect of the crustal movement correction monotonically increases from the beginning of the effective period to the intermediate day of the effective period, and the effect of the crustal deformation correction increases from the intermediate day of the effective period to the end of the effective period. It can be expected to decrease monotonically. As a result, the correction parameter creation device 1 can improve the effect of crustal deformation correction over the entire effective period.

さらに補正パラメータ作成装置１は、一定期間内のＣＯＲＳの座標値に基づき、パラメータ基準日および有効期間中間日におけるＣＯＲＳの位置（すなわち、基準日位置および使用日位置）を算出して、補正パラメータを作成する。これにより、補正パラメータ作成装置１は、補正パラメータに対して、地殻変動によらないノイズによるＣＯＲＳ位置の変動の影響が及ぶのを抑制することができる。 Further, the correction parameter creation device 1 calculates the position of the CORS (that is, the reference date position and the use date position) on the parameter reference date and the intermediate date of the effective period based on the coordinate values of the CORS within a certain period, and sets the correction parameter. create. As a result, the correction parameter creation device 1 can suppress the influence of the fluctuation of the CORS position due to the noise that does not depend on the crustal movement on the correction parameter.

また補正パラメータ作成装置１は、有効期間の中間に有効期間中間日を設定する。これにより、補正パラメータ作成装置１は、有効期間内で最適に補正パラメータによる補正をさせることが期待できる。 Further, the correction parameter creation device 1 sets the valid period intermediate date to the middle of the valid period. As a result, the correction parameter creation device 1 can be expected to perform optimum correction with the correction parameter within the effective period.

以上説明した実施形態において、Ｓ３０は基準日設定部および基準日設定手順としての処理に相当し、Ｓ１０は有効期間設定部および有効期間設定手順としての処理に相当し、Ｓ１０は使用日設定部および使用日設定手順としての処理に相当する。 In the embodiment described above, S30 corresponds to the process as the reference date setting unit and the reference date setting procedure, S10 corresponds to the process as the valid period setting unit and the valid period setting procedure, and S10 is the use date setting unit and This corresponds to the process as the use date setting procedure.

また、Ｓ４０は基準日位置算出部および基準日位置算出手順としての処理に相当し、Ｓ４０は使用日位置算出部および使用日位置算出手順としての処理に相当し、Ｓ５０は変動量算出部および変動量算出手順としての処理に相当し、Ｓ６０〜Ｓ９０は推定量算出部および推定量算出手順としての処理に相当する。 Further, S40 corresponds to a process as a reference date position calculation unit and a reference date position calculation procedure, S40 corresponds to a process as a use date position calculation unit and a use date position calculation procedure, and S50 denotes a variation amount calculation unit and variation. This corresponds to the process as the amount calculation procedure, and S60 to S90 correspond to the process as the estimated amount calculation unit and the estimated amount calculation procedure.

また、パラメータ基準日は基準日に相当し、有効期間はパラメータ有効期間に相当し、有効期間中間日は使用日に相当し、ＣＯＲＳは基準点に相当し、地殻変動量ベクトルＶＭは地殻変動量に相当する。 Further, the parameter reference date corresponds to the reference date, the effective period corresponds to the parameter effective period, the effective period intermediate day corresponds to the use day, CORS corresponds to the reference point, and the crustal movement amount vector VM indicates the crustal movement amount. Equivalent to.

（第２実施形態）
以下に本開示の第２実施形態を図面とともに説明する。なお第２実施形態では、第１実施形態と異なる部分を説明する。共通する構成については同一の符号を付す。 (Second embodiment)
The second embodiment of the present disclosure will be described below with reference to the drawings. In the second embodiment, parts different from the first embodiment will be described. The same components are designated by the same reference numerals.

第２実施形態の補正パラメータ作成装置１は、補正パラメータ作成処理が変更された点が第１実施形態と異なる。
第２実施形態の補正パラメータ作成処理が実行されると、制御部１５は、図８に示すように、まず、Ｓ２１０にて、日付計算処理を実行する。日付計算処理では、制御部１５は、まず、有効期間開始日と有効期間日数を示す有効期間情報を入力するための画像を表示部１１の表示画面に表示する。そして、制御部１５は、有効期間情報が操作入力部１２から入力されるまで待機する。そして、有効期間情報が入力されると、制御部１５は、Ｓ１０と同様にして、式（１）により、有効期間中間日を算出する。さらに制御部１５は、有効期間中間日の年における１月１日をパラメータ基準日とする。但し、制御部１５は、有効期間中間日が１月１日から３月３１日の場合は、有効期間中間日の前年における１月１日をパラメータ基準日とする。例えば、有効期間中間日が２０１８年４月１３日である場合には、パラメータ基準日は２０１８年１月１日である。 The correction parameter creation device 1 of the second embodiment is different from the first embodiment in that the correction parameter creation process is changed.
When the correction parameter creation process of the second embodiment is executed, the control unit 15 first executes a date calculation process in S210, as shown in FIG. In the date calculation process, the control unit 15 first displays on the display screen of the display unit 11 an image for inputting effective period start date and effective period information indicating the number of effective period days. Then, the control unit 15 waits until the valid period information is input from the operation input unit 12. Then, when the valid period information is input, the control unit 15 calculates the valid period intermediate date by the formula (1), as in S10. Further, the control unit 15 sets January 1 of the year on the intermediate date of the effective period as the parameter reference date. However, when the effective period intermediate date is January 1 to March 31, the control unit 15 sets January 1 of the previous year of the effective period intermediate date to the parameter reference date. For example, when the effective period middle date is April 13, 2018, the parameter reference date is January 1, 2018.

Ｓ２１０の処理が終了すると、制御部１５は、Ｓ２２０にて、ＳＤパラメータ取得処理を実行する。
ＳＤパラメータ取得処理では、制御部１５は、まず、補正パラメータを作成する対象となる地域（以下、対象地域）を示す対象地域情報を入力するための画像を表示部１１の表示画面に表示する。そして制御部１５は、対象地域情報が操作入力部１２から入力されるまで待機する。そして、対象地域情報が入力されると、制御部１５は、対象地域情報のＳＤパラメータをデータ記憶部１３から取得する。 When the process of S210 ends, the control unit 15 executes the SD parameter acquisition process in S220.
In the SD parameter acquisition process, the control unit 15 first displays on the display screen of the display unit 11 an image for inputting target area information indicating a target area (hereinafter, target area) for which correction parameters are to be created. Then, the control unit 15 waits until the target area information is input from the operation input unit 12. Then, when the target area information is input, the control unit 15 acquires the SD parameter of the target area information from the data storage unit 13.

ここで、ＳＤパラメータについて説明する。セミ・ダイナミック補正では、相対測位による測位結果に対して、地殻変動による歪みの影響を補正する。ＳＤパラメータは、今期座標と対象地図の測地系の元期座標との差であり、ここで、今期座標は、１年度（４月１日から次年の３月３１日まで）を通して固定のものが使用される。 Here, the SD parameter will be described. In semi-dynamic correction, the effect of distortion due to crustal movement is corrected on the positioning result of relative positioning. The SD parameter is the difference between the current period coordinates and the original period coordinates of the geodetic system of the target map, where the current period coordinates are fixed throughout the first year (April 1 to March 31 of the following year). Is used.

電子基準点は、国土地理院により全国約１３００箇所に設置されている。そして国土地理院は、電子基準点の観測データに基づいて、電子基準点の位置を高精度に計測している。また国土地理院は、電子基準点の日々の座標値を国土地理院のホームページで公開している。電子基準点の日々の座標値のうち、約２週間後に提供される、最終的に解析された座標値はＦ３解と呼ばれる。 Electronic control points have been set up by the Geographical Survey Institute at about 1,300 locations nationwide. The Geographical Survey Institute accurately measures the position of the electronic reference point based on the observation data of the electronic reference point. The Geographical Survey Institute publishes the daily coordinate values of electronic control points on the Geographical Survey Institute website. Of the daily coordinate values of the electronic reference point, the finally analyzed coordinate value provided after about two weeks is called the F3 solution.

また、ＳＤパラメータは、地域メッシュコード毎に設定されている。地域メッシュは、日本国内の地域を緯度および経度に基づいて格子状に区切って形成される。そして地域メッシュコードは、格子状に形成される地域メッシュの南西の格子点（以下、地域メッシュ頂点）における緯度および経度を示すデータである。なお、本実施形態では、ＳＤパラメータおよび地域メッシュコードは、補正パラメータ作成処理を開始する前に予めデータ記憶部１３に記憶される。 The SD parameter is set for each area mesh code. The regional mesh is formed by dividing an area in Japan into a grid pattern based on latitude and longitude. The regional mesh code is data indicating the latitude and longitude at the southwest grid point (hereinafter, regional mesh apex) of the regional mesh formed in a grid pattern. In the present embodiment, the SD parameter and the area mesh code are stored in the data storage unit 13 in advance before starting the correction parameter creating process.

Ｓ２２０の処理が終了すると、制御部１５は、Ｓ２３０にて、補間計算処理を実行する。補間計算処理では、制御部１５は、まず、対象地域内の地域メッシュ頂点のうち、ＳＤパラメータが設定されていない地域メッシュ頂点を選択し、選択した地域メッシュ頂点の地域メッシュコードを取得する。 When the processing of S220 ends, the control unit 15 executes the interpolation calculation processing in S230. In the interpolation calculation process, the control unit 15 first selects a regional mesh vertex whose SD parameter is not set among the regional mesh vertices in the target region, and acquires the regional mesh code of the selected regional mesh vertex.

そして、図９に示すように、制御部１５は、選択した地域メッシュ頂点ＭＲ１の東西に配置されている地域メッシュ頂点ＭＲ２，ＭＲ３を選択し、選択した地域メッシュ頂点ＭＲ２，ＭＲ３に設定されているＳＤパラメータを取得する。ここで、地域メッシュ頂点ＭＲ２，ＭＲ３のＳＤパラメータを取得することができた場合には、制御部１５は、地域メッシュ頂点ＭＲ２，ＭＲ３のＳＤパラメータの線形補間によって、選択した地域メッシュ頂点ＭＲ１のＳＤパラメータを算出する。 Then, as shown in FIG. 9, the control unit 15 selects the regional mesh vertices MR2 and MR3 arranged east and west of the selected regional mesh vertex MR1 and sets the selected regional mesh vertices MR2 and MR3. Get SD parameters. Here, when the SD parameters of the regional mesh vertices MR2 and MR3 can be acquired, the control unit 15 linearly interpolates the SD parameters of the regional mesh vertices MR2 and MR3, and thereby the SD of the selected regional mesh vertices MR1. Calculate the parameters.

一方、地域メッシュ頂点ＭＲ２，ＭＲ３の少なくとも一方のＳＤパラメータを取得することができなかった場合には、制御部１５は、選択した地域メッシュ頂点ＭＲ１の南北に配置されている地域メッシュ頂点ＭＲ４，ＭＲ５を選択し、選択した地域メッシュ頂点ＭＲ４，ＭＲ５に設定されているＳＤパラメータを取得する。ここで、地域メッシュ頂点ＭＲ４，ＭＲ５のＳＤパラメータを取得することができた場合には制御部１５は、地域メッシュ頂点ＭＲ４，ＭＲ５のＳＤパラメータの線形補間によって、選択した地域メッシュ頂点ＭＲ１のＳＤパラメータを算出する。一方、地域メッシュ頂点ＭＲ４，ＭＲ５の少なくとも一方のＳＤパラメータを取得することができなかった場合には、制御部１５は、地域メッシュ頂点ＭＲ１のＳＤパラメータの算出を行わない。 On the other hand, when the SD parameter of at least one of the regional mesh vertices MR2 and MR3 cannot be acquired, the control unit 15 controls the regional mesh vertices MR4 and MR5 located in the north and south of the selected regional mesh vertex MR1. Is selected, and the SD parameters set in the selected regional mesh vertices MR4 and MR5 are acquired. Here, when the SD parameters of the regional mesh vertices MR4 and MR5 can be acquired, the control unit 15 performs the linear interpolation of the SD parameters of the regional mesh vertices MR4 and MR5 to the SD parameters of the selected regional mesh vertices MR1. To calculate. On the other hand, when the SD parameter of at least one of the regional mesh vertices MR4 and MR5 cannot be acquired, the control unit 15 does not calculate the SD parameter of the regional mesh vertex MR1.

次に制御部１５は、ＳＤパラメータが設定されていない地域メッシュ頂点の中から、これまでに選択されていない地域メッシュ頂点が存在しているか否かを判断する。ここで、選択されていない地域メッシュ頂点が存在している場合には、地域メッシュ頂点を新たに選択して、上述の処理を繰り返す。一方、選択されていない地域メッシュ頂点が存在しない場合には、補間計算処理を終了する。 Next, the control unit 15 determines whether or not there is a regional mesh vertex that has not been selected so far among the regional mesh vertices for which the SD parameter is not set. If there is an unselected area mesh vertex, the area mesh vertex is newly selected and the above-described processing is repeated. On the other hand, when there are no unselected regional mesh vertices, the interpolation calculation process ends.

Ｓ２３０の処理が終了すると、図８に示すように、制御部１５は、Ｓ２４０にて、Ｆ３解取得処理を実行する。
Ｆ３解取得処理では、制御部１５は、対象地域内に存在する電子基準点のＦ３解をデータ記憶部１３から取得する。 When the process of S230 ends, as shown in FIG. 8, the control unit 15 executes the F3 solution acquisition process in S240.
In the F3 solution acquisition process, the control unit 15 acquires the F3 solution of the electronic reference point existing in the target area from the data storage unit 13.

Ｆ３解取得処理では、具体的に、制御部１５は、まず、Ｆ３解の数（以下、Ｆ３解取得数）を示す取得情報を入力するための画像を表示部１１の表示画面に表示する。そして制御部１５は、取得情報が操作入力部１２から入力されるまで待機する。そして、取得情報が入力されると、制御部１５は、複数の電子基準点毎に、Ｓ２１０で算出されたパラメータ基準日を中心としたＦ３解取得数分のＦ３解を取得する。さらに制御部１５は、複数の電子基準点毎に、年月日および時刻が新しい順にＦ３解取得数分のＦ３解を取得する。Ｆ３解取得数は、例えば、１００〜２００である。 In the F3 solution acquisition processing, specifically, the control unit 15 first displays an image for inputting acquisition information indicating the number of F3 solutions (hereinafter, the F3 solution acquisition number) on the display screen of the display unit 11. Then, the control unit 15 waits until the acquired information is input from the operation input unit 12. Then, when the acquisition information is input, the control unit 15 acquires, for each of the plurality of electronic reference points, F3 solutions for the number of F3 solutions acquired centered on the parameter reference date calculated in S210. Further, the control unit 15 acquires, for each of the plurality of electronic reference points, F3 solutions for the number of F3 solutions acquired in the order of latest date and time. The number of F3 solution acquisitions is, for example, 100 to 200.

次にＳ２５０にて、制御部１５は、Ｆ３解位置決定処理を実行する。Ｆ３解位置決定処理では、制御部１５は、Ｓ２４０でＦ３解を取得した複数の電子基準点毎に、パラメータ基準日における電子基準点の位置と、有効期間中間日における電子基準点の位置とを算出する。 Next, in S250, the control unit 15 executes F3 solution position determination processing. In the F3 solution position determination process, the control unit 15 determines the position of the electronic reference point on the parameter reference date and the position of the electronic reference point on the intermediate day of the valid period for each of the plurality of electronic reference points for which the F3 solution was acquired in S240. calculate.

Ｆ３解位置決定処理では、具体的には、制御部１５は、まず、複数の電子基準点のうち、パラメータ基準日および有効期間中間日における位置が算出されていない１つの電子基準点を選択する。そして制御部１５は、選択した電子基準点において、式（１２），（１３），（１４）により、第１算出期間開始日ＤＵＦＢ１および第１算出期間終了日ＤＵＦＥ１を算出する。式（１２）〜（１４）におけるＤＲは、パラメータ基準日である。式（１２）〜（１４）におけるＮＵＦ３は、Ｆ３解取得数である。また、式（１３）は、Ｆ３解取得数が偶数である場合に用いられ、式（１４）は、Ｆ３解取得数が奇数である場合に用いられる。なお、第１算出期間開始日ＤＵＦＢ１から第１算出期間終了日ＤＵＦＥ１までの期間を第１算出期間ＤＵＦ１という。 In the F3 solution position determination process, specifically, the control unit 15 first selects, from among the plurality of electronic reference points, one electronic reference point for which the position on the parameter reference date and the intermediate date of the effective period is not calculated. .. Then, the control unit 15 calculates the first calculation period start date DUFB1 and the first calculation period end date DUFE1 by the equations (12), (13), and (14) at the selected electronic reference point. DR in Expressions (12) to (14) is a parameter reference date. NUF3 in Expressions (12) to (14) is the F3 solution acquisition number. The equation (13) is used when the F3 solution acquisition number is an even number, and the equation (14) is used when the F3 solution acquisition number is an odd number. The period from the first calculation period start date DUFB1 to the first calculation period end date DUFE1 is referred to as the first calculation period DUF1.

そして制御部１５は、第１算出期間ＤＵＦ１内のＦ３解を用い、緯度、経度および楕円体高のそれぞれについて線形回帰直線を算出する。さらに制御部１５は、緯度、経度および楕円体高のそれぞれで算出された線形回帰直線について、パラメータ基準日における緯度、経度および楕円体高を算出する。図１０に示すように、点Ｐ４１，Ｐ４２，Ｐ４３，Ｐ４４，Ｐ４５，Ｐ４６，Ｐ４７，Ｐ４８は、第１算出期間ＤＵＦ１内のＦ３解における緯度、経度または楕円体高を示す。直線Ｌ３は、線形回帰直線を示す。点Ｐ４９は、線形回帰直線に基づいて算出されたパラメータ基準日ＤＲにおける緯度、経度または楕円体高である。

Then, the control unit 15 uses the F3 solution in the first calculation period DUF1 to calculate a linear regression line for each of the latitude, longitude, and ellipsoidal height. Further, the control unit 15 calculates the latitude, the longitude and the ellipsoidal height on the parameter reference day for the linear regression line calculated for each of the latitude, the longitude and the ellipsoidal height. As shown in FIG. 10, points P41, P42, P43, P44, P45, P46, P47, P48 indicate the latitude, longitude, or ellipsoidal height in the F3 solution within the first calculation period DUF1. The straight line L3 shows a linear regression line. The point P49 is the latitude, longitude, or ellipsoidal height on the parameter reference day DR calculated based on the linear regression line.

次に制御部１５は、選択した電子基準点において、式（１５），（１６）により、第２算出期間開始日ＤＵＦＢ２および第２算出期間終了日ＤＵＦＥ２を算出する。式（１５）におけるＤＬは、Ｓ２４０で取得した最新のＦ３解の日付である。式（１６）におけるＮＵＦ３は、Ｆ３解取得数である。なお、第２算出期間開始日ＤＵＦＢ２から第２算出期間終了日ＤＵＦＥ２までの期間を第２算出期間ＤＵＦ２という。 Next, the control unit 15 calculates the second calculation period start date DUFB2 and the second calculation period end date DUFE2 using the equations (15) and (16) at the selected electronic reference point. DL in Expression (15) is the latest F3 solution date acquired in S240. NUF3 in Expression (16) is the F3 solution acquisition number. The period from the second calculation period start date DUFB2 to the second calculation period end date DUFE2 is referred to as the second calculation period DUF2.

そして制御部１５は、第２算出期間ＤＵＦ２内のＦ３解を用い、緯度、経度および楕円体高のそれぞれについて線形回帰直線を算出する。さらに制御部１５は、緯度、経度および楕円体高のそれぞれで算出された線形回帰直線について、有効期間中間日における緯度、経度および楕円体高を算出する。図１１に示すように、点Ｐ５１，Ｐ５２，Ｐ５３，Ｐ５４，Ｐ５５は、第２算出期間ＤＵＦ２内のＦ３解における緯度、経度または楕円体高を示す。直線Ｌ４は、線形回帰直線を示す。点Ｐ５６は、線形回帰直線に基づいて算出された有効期間中間日ＤＭにおける緯度、経度または楕円体高である。

Then, the control unit 15 uses the F3 solution in the second calculation period DUF2 to calculate a linear regression line for each of the latitude, longitude, and ellipsoidal height. Further, the control unit 15 calculates the latitude, the longitude and the ellipsoidal height on the intermediate day of the effective period for the linear regression lines calculated for the latitude, the longitude and the ellipsoidal height. As shown in FIG. 11, points P51, P52, P53, P54, and P55 indicate the latitude, longitude, or ellipsoidal height in the F3 solution within the second calculation period DUF2. The straight line L4 shows a linear regression line. The point P56 is the latitude, longitude, or ellipsoid height on the effective period middle day DM calculated based on the linear regression line.

そして制御部１５は、パラメータ基準日および有効期間中間日における位置が算出されていない電子基準点が存在しているか否かを判断する。ここで、位置が算出されていない電子基準点が存在している場合には、電子基準点を新たに選択して、上述の処理を繰り返す。一方、位置が算出されていない電子基準点が存在しない場合には、Ｆ３解位置決定処理を終了する。 Then, the control unit 15 determines whether or not there is an electronic reference point whose position has not been calculated on the parameter reference date and the valid period intermediate date. Here, if there is an electronic reference point whose position has not been calculated, the electronic reference point is newly selected and the above-described processing is repeated. On the other hand, if there is no electronic reference point whose position has not been calculated, the F3 solution position determination processing ends.

Ｓ２５０の処理が終了すると、図８に示すように、制御部１５は、Ｓ２６０にて、Ｆ３解変動量算出処理を実行する。具体的には、制御部１５は、Ｓ５０と同様にして、Ｓ２４０でＦ３解を取得した複数の電子基準点毎に、パラメータ基準日における電子基準点の座標値（すなわち、緯度、経度および楕円体高）を始点とし、有効期間中間日における電子基準点の座標値を終点とする地殻変動量ベクトルＶＭを算出する。 When the process of S250 ends, as shown in FIG. 8, the control unit 15 executes the F3 solution variation amount calculation process in S260. Specifically, similarly to S50, the control unit 15 sets, for each of the plurality of electronic reference points for which the F3 solution is acquired in S240, the coordinate values of the electronic reference points on the parameter reference date (that is, the latitude, the longitude, and the ellipsoid height). ) As the starting point, and the coordinate value of the electronic reference point on the intermediate day of the effective period as the ending point is used to calculate the crustal movement vector VM.

次に制御部１５は、Ｓ２７０にて、経験バリオグラム作成処理を実行する。具体的には、制御部１５は、ＣＯＲＳの代わりに電子基準点を用いる点以外はＳ６０と同様にして、経験バリオグラムを作成する。 Next, the control unit 15 executes the experience variogram creation process in S270. Specifically, the control unit 15 creates an empirical variogram in the same manner as S60 except that an electronic reference point is used instead of CORS.

そして制御部１５は、Ｓ２８０にて、Ｓ７０と同様にして、モデル推定処理を実行する。さらに制御部１５は、Ｓ２９０にて、ＣＯＲＳの代わりに電子基準点を用いる点以外はＳ８０と同様にして、クリギング法計算処理を実行する。 Then, in S280, the control unit 15 executes the model estimation process as in S70. Further, in S290, the control unit 15 executes the Kriging method calculation process in the same manner as in S80 except that the electronic reference point is used instead of CORS.

次に制御部１５は、Ｓ３００にて、地殻変動量推定処理を実行する。
地殻変動量推定処理では、具体的に、制御部１５は、まず、対象地域内の地域メッシュ頂点の中から、後述する地殻変動推定量が算出されていない地域メッシュ頂点を一つ選択する。そして制御部１５は、選択した地域メッシュ頂点の地域メッシュコード（すなわち、メッシュ頂点における緯度および経度）を取得する。次に制御部１５は、係数ベクトルｖを算出する。係数ベクトルｖは、ｖ_１，ｖ_２，ｖ_３，・・・・，ｖ_ｎ，ｖ_ｎ＋１を成分とするベクトルである。ｎは、Ｓ２４０でＦ３解を取得した電子基準点の数である。ｉ＝１，２，３，・・・，ｎである場合において、ｖ_ｉ＝γ（ｒ_ｉ）である。ｒ_ｉは、選択した地域メッシュ頂点と、ｉで特定される電子基準点との間の距離である。また、ｖ_ｎ＋１＝１である。 Next, the control unit 15 executes a crustal deformation amount estimation process in S300.
In the crustal deformation amount estimation process, specifically, the control unit 15 first selects one regional mesh vertex from which a crustal deformation estimated amount to be described later has not been calculated, from the regional mesh vertices in the target area. Then, the control unit 15 acquires the regional mesh code of the selected regional mesh vertex (that is, the latitude and longitude at the mesh vertex). Next, the control unit 15 calculates the coefficient vector v. The coefficient vector v is a vector whose components are v ₁ , v ₂ , v ₃ ,..., V _n , v _n+1 . n is the number of electronic reference points for which the F3 solution was acquired in S240. In the case where i=1, 2, 3,..., N, v _i =γ(r _i ). r _i is the distance between the selected regional mesh vertex and the electronic reference point identified by i. Also, v _n+1 =1.

そして制御部１５は、変数ベクトルｗを用いて、式（９）に示す方程式を解く。変数ベクトルｗは、式（９）に示すように、ｗ_１，ｗ_２，ｗ_３，・・・・，ｗ_ｎ，λを成分とするベクトルである。λは、ラグランジュの未定乗数である。すなわち、制御部１５は、ラグランジュの未定乗数法により、式（９）を解く。なお、制御部１５は、式（９）を解く際に、係数行列ＡのＬＵ分解を行う。したがって、制御部１５は、Ａｗ＝（ＬＵ）ｗ＝ｖを解くことにより、ｗ_１，ｗ_２，ｗ_３，・・・・，ｗ_ｎ，λを算出する。そして制御部１５は、上記のｗ_１，ｗ_２，ｗ_３，・・・・，ｗ_ｎ，λの算出を各成分（すなわち、緯度、経度および楕円体高）について実行する。 Then, the control unit 15 uses the variable vector w to solve the equation shown in Expression (9). The variable vector w is a vector whose components are w ₁ , w ₂ , w ₃ ,..., W _n , λ, as shown in Expression (9). λ is a Lagrange's undetermined multiplier. That is, the control unit 15 solves the equation (9) by the Lagrange's undetermined multiplier method. The control unit 15 performs LU decomposition of the coefficient matrix A when solving the equation (9). Therefore, the control unit 15 calculates w ₁ , w ₂ , w ₃ ,..., W _n , λ by solving Aw=(LU)w=v. Then, the control unit 15 executes the above-described calculation of w ₁ , w ₂ , w ₃ ,..., W _n , λ for each component (that is, latitude, longitude, and ellipsoid height).

ｗ_１，ｗ_２，ｗ_３，・・・・，ｗ_ｎを算出した後に、制御部１５は、式（１０）により、選択した地域メッシュ頂点の地殻変動推定量ｚを算出する。式（１０）におけるｚ_ｉは、ｉで特定される電子基準点における各成分（すなわち、緯度、経度および楕円体高）の地殻変動量である。 After calculating w ₁ , w ₂ , w ₃ ,..., W _n , the control unit 15 calculates the crustal movement estimator z of the selected regional mesh vertices by the equation (10). Z _i in Expression (10) is the amount of crustal movement of each component (that is, latitude, longitude, and ellipsoidal height) at the electronic reference point specified by i.

さらに制御部１５は、式（１１）により、選択した地域メッシュ頂点における地殻変動推定量の予測誤差σを算出する。
そして制御部１５は、上記の地殻変動推定量ｚおよび予測誤差σの算出を各成分（すなわち、緯度、経度および楕円体高）について実行する。 Further, the control unit 15 calculates the prediction error σ of the crustal movement estimation amount at the selected regional mesh vertices by the equation (11).
Then, the control unit 15 executes the calculation of the crustal movement estimated amount z and the prediction error σ for each component (that is, latitude, longitude, and ellipsoid height).

各成分の地殻変動推定量ｚおよび予測誤差σの算出が終了すると、制御部１５は、地殻変動推定量ｚおよび予測誤差σが算出されていない地域メッシュ頂点が存在しているか否かを判断する。ここで、地殻変動推定量ｚおよび予測誤差σが算出されていない地域メッシュ頂点が存在している場合には、地域メッシュ頂点を新たに選択して、上述の処理を繰り返す。一方、地殻変動推定量ｚおよび予測誤差σが算出されていない地域メッシュ頂点が存在しない場合には、地殻変動量推定処理を終了する。 When the calculation of the crustal movement estimated amount z and the prediction error σ of each component is completed, the control unit 15 determines whether or not there is a regional mesh vertex for which the crustal movement estimated amount z and the prediction error σ are not calculated. .. Here, when there is a regional mesh vertex for which the crustal movement estimator z and the prediction error σ have not been calculated, a new regional mesh vertex is selected and the above process is repeated. On the other hand, when there is no regional mesh vertex for which the crustal movement estimated amount z and the prediction error σ have not been calculated, the crustal movement estimation process ends.

Ｓ３００の処理が終了すると、制御部１５は、Ｓ３１０にて、ＳＤ／Ｒパラメータ算出処理を実行する。ＳＤ／Ｒパラメータ算出処理では、制御部１５は、対象地域内の複数の地域メッシュ毎に、地殻変動推定量とＳＤパラメータとを合成して、ＳＤ／Ｒパラメータを算出する。 When the process of S300 ends, the control unit 15 executes the SD/R parameter calculation process in S310. In the SD/R parameter calculation process, the control unit 15 calculates the SD/R parameter by synthesizing the crustal movement estimated amount and the SD parameter for each of a plurality of regional meshes in the target region.

ＳＤ／Ｒパラメータ算出処理では、具体的に、制御部１５は、まず、対象地域内の地域メッシュ頂点の中から、ＳＤ／Ｒパラメータが算出されていない地域メッシュ頂点を一つ選択する。そして制御部１５は、選択した地域メッシュ頂点の地域メッシュコードを取得する。次に制御部１５は、取得した地域メッシュコードに対応するＳＤパラメータを取得する。さらに制御部１５は、取得した地域メッシュコードに対応する地殻変動推定量ｚを取得する。そして制御部１５は、式（１７）により、ＳＤ／Ｒパラメータを算出する。式（１７）におけるＶＭＳＤＲは、ＳＤ／Ｒパラメータである。式（１７）におけるＶＭＳＤは、ＳＤパラメータである。式（１７）におけるＶＭＥは、地殻変動推定量ｚ（緯度、経度、楕円体高）である。 In the SD/R parameter calculation process, specifically, the control unit 15 first selects one area mesh vertex whose SD/R parameter is not calculated from the area mesh vertices in the target area. Then, the control unit 15 acquires the regional mesh code of the selected regional mesh vertex. Next, the control unit 15 acquires the SD parameter corresponding to the acquired regional mesh code. Further, the control unit 15 acquires the crustal movement estimated amount z corresponding to the acquired regional mesh code. Then, the control unit 15 calculates the SD/R parameter by the equation (17). VMSDR in Expression (17) is an SD/R parameter. VMSD in Expression (17) is an SD parameter. VME in the equation (17) is the crustal movement estimated amount z (latitude, longitude, ellipsoid height).

図１２に示すように、有効期間中間日ＤＭにおける地殻変動量は、パラメータ基準日ＤＲにおける地殻変動量（すなわち、ＳＤパラメータ）と、地殻変動推定量とを加算することにより算出される。

As shown in FIG. 12, the crustal movement amount on the effective period intermediate day DM is calculated by adding the crustal movement amount (that is, the SD parameter) on the parameter reference date DR and the crustal movement estimated amount.

次に制御部１５は、ＳＤ／Ｒパラメータが算出されていない地域メッシュ頂点が存在しているか否かを判断する。ここで、ＳＤ／Ｒパラメータが算出されていない地域メッシュ頂点が存在している場合には、地域メッシュ頂点を新たに選択して、上述の処理を繰り返す。一方、ＳＤ／Ｒパラメータが算出されていない地域メッシュ頂点が存在しない場合には、ＳＤ／Ｒパラメータ算出処理を終了する。 Next, the control unit 15 determines whether or not there are regional mesh vertices for which SD/R parameters have not been calculated. If there is a regional mesh vertex for which the SD/R parameter has not been calculated, a new regional mesh vertex is selected and the above process is repeated. On the other hand, if there is no regional mesh vertex for which the SD/R parameter has not been calculated, the SD/R parameter calculation process ends.

Ｓ３１０の処理が終了すると、制御部１５は、図８に示すように、Ｓ３２０にて、ＳＤ／Ｒパラメータ出力処理を実行する。ＳＤ／Ｒパラメータ出力処理では、制御部１５は、Ｓ３１０で算出したＳＤ／Ｒパラメータを示すＳＤ／Ｒパラメータ情報を、補正パラメータ作成装置１の外部へ出力する。 When the process of S310 is completed, the control unit 15 executes the SD/R parameter output process in S320, as shown in FIG. In the SD/R parameter output process, the control unit 15 outputs the SD/R parameter information indicating the SD/R parameter calculated in S310 to the outside of the correction parameter creation device 1.

Ｓ３２０の処理が終了すると、制御部１５は、Ｓ３３０にて、予測誤差出力処理を実行する。予測誤差出力処理では、制御部１５は、Ｓ３００で算出した予測誤差σを示す予測誤差情報を、補正パラメータ作成装置１の外部へ出力する。 When the process of S320 ends, the control unit 15 executes the prediction error output process in S330. In the prediction error output process, the control unit 15 outputs the prediction error information indicating the prediction error σ calculated in S300 to the outside of the correction parameter creation device 1.

そして、Ｓ３３０の処理が終了すると、制御部１５は、補正パラメータ作成処理を終了する。
このように構成された補正パラメータ作成装置１は、単独測位により測定された位置を地殻変動に基づいて補正するための補正パラメータを作成する。 Then, when the process of S330 ends, the control unit 15 ends the correction parameter creation process.
The correction parameter creating device 1 configured as described above creates a correction parameter for correcting the position measured by the single positioning based on the crustal movement.

そして補正パラメータ作成装置１は、パラメータ基準日を設定する。また補正パラメータ作成装置１は、有効期間を設定する。また補正パラメータ作成装置１は、有効期間中間日を設定する。 Then, the correction parameter creation device 1 sets the parameter reference date. Further, the correction parameter creating device 1 sets the valid period. Further, the correction parameter creation device 1 sets the effective period intermediate date.

そして補正パラメータ作成装置１は、補正パラメータを作成する対象となる対象地域内に存在する複数の電子基準点のそれぞれについて、基準日を含む一定期間内の電子基準点の座標値に基づき、パラメータ基準日における電子基準点の位置（以下、基準日位置）を算出する。 Then, the correction parameter creation device 1 determines the parameter reference based on the coordinate values of the electronic reference points within a certain period including the reference date for each of the plurality of electronic reference points existing in the target area where the correction parameter is created. The position of the electronic reference point on the day (hereinafter referred to as the reference date position) is calculated.

また補正パラメータ作成装置１は、対象地域内に存在する複数の電子基準点のそれぞれについて、算出時点に得られる一定期間内の電子基準点の座標値に基づき、有効期間中間日における電子基準点の位置（以下、使用日位置）を算出する。 Further, the correction parameter creating device 1 determines the electronic reference points on the intermediate day of the effective period based on the coordinate values of the electronic reference points within a certain period obtained at the time of calculation for each of the plurality of electronic reference points existing in the target area. The position (hereinafter referred to as the use date position) is calculated.

さらに補正パラメータ作成装置１は、対象地域内に存在する複数の電子基準点のそれぞれについて、算出された基準日位置と、算出された使用日位置との差を、地殻変動量ベクトルＶＭとして算出する。 Further, the correction parameter creation device 1 calculates the difference between the calculated reference date position and the calculated use date position for each of the plurality of electronic reference points existing in the target area as the crustal movement amount vector VM. ..

そして補正パラメータ作成装置１は、対象地域を格子状に区切って形成された複数の地域メッシュのそれぞれについて、クリギング法を用いて、複数の電子基準点の地殻変動量ベクトルＶＭに基づき、地域メッシュ頂点における地殻変動推定量ｚを算出し、算出した地殻変動推定量ｚを補正パラメータとする。 Then, the correction parameter creation device 1 uses the kriging method for each of a plurality of area meshes formed by dividing the target area into a grid pattern, and based on the crustal movement amount vector VM of the plurality of electronic reference points, the area mesh vertices. The crustal deformation estimated amount z at is calculated, and the calculated crustal deformation estimated amount z is used as a correction parameter.

このように第２実施形態の補正パラメータ作成装置１は、有効期間内の時点での地殻変動量に基づいて、補正パラメータを作成する。このため、第２実施形態の補正パラメータ作成装置１は、第１実施形態の補正パラメータ作成装置１と同様の効果を得ることができる。 As described above, the correction parameter creation device 1 of the second embodiment creates the correction parameter based on the crustal movement amount at the time point within the effective period. Therefore, the correction parameter creation device 1 of the second embodiment can obtain the same effect as the correction parameter creation device 1 of the first embodiment.

以上説明した実施形態において、Ｓ２１０は基準日設定部および基準日設定手順としての処理に相当し、Ｓ２１０は有効期間設定部および有効期間設定手順としての処理に相当し、Ｓ２１０は使用日設定部および使用日設定手順としての処理に相当する。 In the embodiment described above, S210 corresponds to a process as a reference date setting unit and a reference date setting procedure, S210 corresponds to a process as an effective period setting unit and an effective period setting procedure, and S210 denotes a use date setting unit and This corresponds to the process as the use date setting procedure.

また、Ｓ２５０は基準日位置算出部および基準日位置算出手順としての処理に相当し、Ｓ２５０は使用日位置算出部および使用日位置算出手順としての処理に相当し、Ｓ２６０は変動量算出部および変動量算出手順としての処理に相当し、Ｓ２７０〜Ｓ３００は推定量算出部および推定量算出手順としての処理に相当する。また、電子基準点は基準点に相当する。 Further, S250 corresponds to a process as a reference date position calculation unit and a reference date position calculation procedure, S250 corresponds to a process as a use date position calculation unit and a use date position calculation procedure, and S260 denotes a variation amount calculation unit and variation. This corresponds to the process as the amount calculation procedure, and S270 to S300 correspond to the process as the estimated amount calculation unit and the estimated amount calculation procedure. The electronic reference point corresponds to the reference point.

以上、本開示の一実施形態について説明したが、本開示は上記実施形態に限定されるものではなく、種々変形して実施することができる。
例えば上記実施形態では、使用日が有効期間の中間に設定される形態を示したが、有効期間の中間に限定されるものではなく、有効期間内の時点であればよい。 Although one embodiment of the present disclosure has been described above, the present disclosure is not limited to the above embodiment, and various modifications can be implemented.
For example, in the above-described embodiment, the usage date is set to the middle of the valid period, but the present invention is not limited to the middle of the valid period, and may be any time within the valid period.

上記実施形態では、クリギング法を用いて地殻変動推定量を算出する形態を示した。しかし、クリギング法以外の空間補間法を用いて地殻変動推定量を算出するようにしてもよい。 In the above-described embodiment, the mode in which the crustal deformation estimator is calculated using the Kriging method has been shown. However, the crustal movement estimator may be calculated using a spatial interpolation method other than the Kriging method.

上記実施形態では、日付を示すパラメータ基準日および有効期間中間日を設定する形態を示した。しかし、パラメータ基準日および有効期間中間日の代わりに、日付だけではなく時刻も含んだパラメータ基準日時および有効期間中間日時を設定するようにしてもよい。 In the above embodiment, the parameter reference date indicating the date and the valid period intermediate date are set. However, instead of the parameter reference date and the valid period intermediate date, the parameter reference date and the valid period intermediate date and time including not only the date but also the time may be set.

また、上記実施形態における１つの構成要素が有する機能を複数の構成要素に分担させたり、複数の構成要素が有する機能を１つの構成要素に発揮させたりしてもよい。また、上記実施形態の構成の一部を省略してもよい。また、上記実施形態の構成の少なくとも一部を、他の上記実施形態の構成に対して付加、置換等してもよい。 Further, the function of one constituent element in the above-described embodiment may be shared by a plurality of constituent elements, or the function of a plurality of constituent elements may be exerted by one constituent element. Moreover, you may omit a part of structure of the said embodiment. Further, at least a part of the configuration of the above embodiment may be added to or replaced with the configuration of the other above embodiment.

上述した補正パラメータ作成装置１の他、当該補正パラメータ作成装置１を構成要素とするシステム、当該補正パラメータ作成装置１としてコンピュータを機能させるためのプログラム、このプログラムを記録した媒体、補正パラメータ作成方法など、種々の形態で本開示を実現することもできる。 In addition to the correction parameter creation device 1 described above, a system including the correction parameter creation device 1 as a component, a program for causing a computer to function as the correction parameter creation device 1, a medium in which this program is recorded, a correction parameter creation method, etc. The present disclosure can be realized in various forms.

１…補正パラメータ作成装置、１５…制御部、２０…補正パラメータ作成プログラム DESCRIPTION OF SYMBOLS 1... Correction parameter creation device, 15... Control part, 20... Correction parameter creation program

Claims (4)

単独測位により測定された位置を地殻変動に基づいて補正するための補正パラメータを作成する補正パラメータ作成装置であって、
前記補正パラメータを作成する基準となる日である基準日を設定するように構成された基準日設定部と、
前記補正パラメータが有効な期間であるパラメータ有効期間を設定するように構成された有効期間設定部と、
前記パラメータ有効期間内に、前記パラメータ有効期間内の任意の日である使用日を設定するように構成された使用日設定部と、
前記補正パラメータを作成する対象となる地域である対象地域内に存在する複数の基準点のそれぞれについて、前記基準日を含む一定期間内の前記基準点の座標値に基づき、前記基準日における前記基準点の位置を基準日位置として算出するように構成された基準日位置算出部と、
前記対象地域内に存在する複数の前記基準点のそれぞれについて、算出時点に得られる一定期間内の前記基準点の座標値に基づき、前記使用日における前記基準点の位置を使用日位置として算出するように構成された使用日位置算出部と、
前記対象地域内に存在する複数の前記基準点のそれぞれについて、前記基準日位置算出部により算出された前記基準日位置と、前記使用日位置算出部により算出された前記使用日位置との差を、地殻変動量として算出するように構成された変動量算出部と、
前記対象地域を格子状に区切って形成された複数のメッシュのそれぞれの格子点にあたる地点について、空間補間法を用いて、前記地点を囲む複数の前記基準点の前記地殻変動量に基づき、地殻変動推定量を算出し、算出した前記地殻変動推定量を前記補正パラメータとするように構成された推定量算出部と
を備える補正パラメータ作成装置。 A correction parameter creating device for creating a correction parameter for correcting a position measured by independent positioning based on crustal movement,
A reference date setting unit configured to set a reference date which is a reference date for creating the correction parameter ,
An effective period setting unit configured to set a parameter effective period, wherein the correction parameter is an effective period,
Within the parameter validity period, a usage date setting unit configured to set a usage date that is any day within the parameter validity period ,
For each of the plurality of reference points existing in the target area that is the area for which the correction parameter is created, based on the coordinate value of the reference point within a certain period including the reference date, the reference on the reference day A reference date position calculation unit configured to calculate the position of the point as the reference date position,
For each of the plurality of reference points existing in the target area, the position of the reference point on the use day is calculated as a use day position based on the coordinate value of the reference point within a certain period obtained at the time of calculation. A date-of-use position calculator configured as described above,
For each of the plurality of reference points existing in the target area, the difference between the reference date position calculated by the reference date position calculation unit and the use date position calculated by the use date position calculation unit, , A fluctuation amount calculation unit configured to calculate as a crustal deformation amount,
For the points corresponding to the respective grid points of a plurality of meshes formed by dividing the target area into a grid pattern, using spatial interpolation, based on the crustal movements of the plurality of reference points surrounding the points, crustal movements An estimated amount calculation unit configured to calculate an estimated amount and use the calculated crustal movement estimated amount as the correction parameter. 請求項１に記載の補正パラメータ作成装置であって、
前記使用日設定部は、前記パラメータ有効期間の中間に前記使用日を設定する補正パラメータ作成装置。 The correction parameter creating device according to claim 1,
The usage date setting unit is a correction parameter creation device that sets the usage date in the middle of the parameter valid period. 単独測位により測定された位置を地殻変動に基づいて補正するための補正パラメータを作成する補正パラメータ作成方法であって、
前記補正パラメータを作成する基準となる日である基準日を設定する基準日設定手順と、
前記補正パラメータが有効な期間であるパラメータ有効期間を設定する有効期間設定手順と、
前記パラメータ有効期間内に、前記パラメータ有効期間内の任意の日である使用日を設定する使用日設定手順と、
前記補正パラメータを作成する対象となる地域である対象地域内に存在する複数の基準点のそれぞれについて、前記基準日を含む一定期間内の前記基準点の座標値に基づき、前記基準日における前記基準点の位置を基準日位置として算出する基準日位置算出手順と、
前記対象地域内に存在する複数の前記基準点のそれぞれについて、算出時点に得られる一定期間内の前記基準点の座標値に基づき、前記使用日における前記基準点の位置を使用日位置として算出する使用日位置算出手順と、
前記対象地域内に存在する複数の前記基準点のそれぞれについて、前記基準日位置算出手順により算出された前記基準日位置と、前記使用日位置算出手順により算出された前記使用日位置との差を、地殻変動量として算出する変動量算出手順と、
前記対象地域を格子状に区切って形成された複数のメッシュのそれぞれの格子点にあたる地点について、空間補間法を用いて、前記地点を囲む複数の前記基準点の前記地殻変動量に基づき、地殻変動推定量を算出し、算出した前記地殻変動推定量を前記補正パラメータとする推定量算出手順と
を備える補正パラメータ作成方法。 A correction parameter creating method for creating a correction parameter for correcting a position measured by independent positioning based on crustal movement,
A reference date setting procedure for setting a reference date which is a reference date for creating the correction parameter ,
A valid period setting procedure for setting a parameter valid period in which the correction parameter is valid,
Within the parameter validity period, a usage date setting procedure for setting a usage date that is any day within the parameter validity period ,
For each of a plurality of reference points existing in the target area that is the area for which the correction parameter is created, based on the coordinate value of the reference point within a certain period including the reference date, the reference on the reference day A reference date position calculation procedure that calculates the position of the point as the reference date position,
For each of the plurality of reference points existing in the target area, the position of the reference point on the use day is calculated as a use day position based on the coordinate value of the reference point within a certain period obtained at the time of calculation. Use date position calculation procedure,
For each of the plurality of reference points existing in the target area, the difference between the reference date position calculated by the reference date position calculation procedure and the use date position calculated by the use date position calculation procedure is calculated. , A variation calculation procedure for calculating the amount of crustal movement,
For the points corresponding to the respective grid points of the plurality of meshes formed by dividing the target area in a grid pattern, using spatial interpolation, based on the crustal movements of the plurality of reference points surrounding the points, crustal movements An estimated amount calculation procedure in which an estimated amount is calculated and the calculated crustal movement estimated amount is used as the correction parameter. 単独測位により測定された位置を地殻変動に基づいて補正するための補正パラメータを作成するために、コンピュータを、
前記補正パラメータを作成する基準となる日である基準日を設定するように構成された基準日設定部、
前記補正パラメータが有効な期間であるパラメータ有効期間を設定するように構成された有効期間設定部、
前記パラメータ有効期間内に、前記パラメータ有効期間内の任意の日である使用日を設定するように構成された使用日設定部、
前記補正パラメータを作成する対象となる地域である対象地域内に存在する複数の基準点のそれぞれについて、前記基準日を含む一定期間内の前記基準点の座標値に基づき、前記基準日における前記基準点の位置を基準日位置として算出するように構成された基準日位置算出部、
前記対象地域内に存在する複数の前記基準点のそれぞれについて、算出時点に得られる一定期間内の前記基準点の座標値に基づき、前記使用日における前記基準点の位置を使用日位置として算出するように構成された使用日位置算出部、
前記対象地域内に存在する複数の前記基準点のそれぞれについて、前記基準日位置算出部により算出された前記基準日位置と、前記使用日位置算出部により算出された前記使用日位置との差を、地殻変動量として算出するように構成された変動量算出部、及び、
前記対象地域を格子状に区切って形成された複数のメッシュのそれぞれの格子点にあたる地点について、空間補間法を用いて、前記地点を囲む複数の前記基準点の前記地殻変動量に基づき、地殻変動推定量を算出し、算出した前記地殻変動推定量を前記補正パラメータとするように構成された推定量算出部
として機能させるための補正パラメータ作成プログラム。 A computer is used to create a correction parameter for correcting the position measured by independent positioning based on the crustal movement.
A reference date setting unit configured to set a reference date which is a reference date for creating the correction parameter ,
An effective period setting unit configured to set a parameter effective period in which the correction parameter is an effective period,
Within the parameter validity period, a usage date setting unit configured to set a usage date that is any day within the parameter validity period ,
For each of the plurality of reference points existing in the target area that is the area for which the correction parameter is created, based on the coordinate value of the reference point within a certain period including the reference date, the reference on the reference day A reference date position calculation unit configured to calculate the position of the point as the reference date position,
For each of the plurality of reference points existing in the target area, the position of the reference point on the use day is calculated as a use day position based on the coordinate value of the reference point within a certain period obtained at the time of calculation. A date-of-use position calculator configured as
For each of the plurality of reference points existing in the target area, the difference between the reference date position calculated by the reference date position calculation unit and the use date position calculated by the use date position calculation unit, A fluctuation amount calculation unit configured to calculate a crustal deformation amount, and
For the points corresponding to the respective grid points of a plurality of meshes formed by dividing the target area into a grid pattern, using spatial interpolation, based on the crustal movements of the plurality of reference points surrounding the points, crustal movements A correction parameter creation program for causing an estimated amount to be calculated, and for causing the calculated crustal deformation estimated amount to serve as the correction parameter.

Images