JP5567805B2 - Flying object detection method, system, and program - Google Patents

Description

  The present invention relates to a method and system for detecting flying rockets, missiles, and the like. More specifically, the present invention detects an object flying in the sky of several kilometers to several tens of kilometers (hereinafter, referred to as a “flying object”). The present invention relates to a method and system for estimating a direction, a destination, and the like.

  In recent years, if a flying object such as a rocket or missile has caused an accident such as an accidental crash or explosion, not only the human damage and property damage caused by the crash but also the crash point. And damage to the impact point is expected to become extremely large, and may cause serious damage.

  For this reason, a conventional system for monitoring a vehicle or a ship equipped with a radar system on the ground or the sea to monitor whether the flying object is flying normally or where the flying object is flying from where. Construction has been done.

  However, since the range of action of vehicles deployed on the ground or ships deployed on the sea is mainly limited within Japan, the monitoring area is limited from the global scale. Therefore, it was very difficult to detect a projectile that was launched from a distance that crossed countries and regions. In addition, when a vehicle equipped with a radar system is deployed on the ground, there is a problem that an obstacle such as a mountain narrows the monitoring area depending on the terrain.

  Therefore, an attempt has been made to construct a system using artificial satellites to complement the flying object monitoring system using vehicles and ships. For example, research and development of an early warning satellite system such as Non-Patent Document 1 has been conducted. ing.

  The technology disclosed in Non-Patent Document 1 uses two artificial satellites having a long elliptical orbit and equipped with an observation device. In general, since an artificial satellite having a long elliptical orbit flies at a high altitude (for example, an altitude of 40,000 km), it is assumed that many areas on the earth can be under surveillance.

  FIG. 19 shows a region that can be monitored when two artificial satellites are orbited (Satellite 1 orbit, Satellite 2 orbit). In FIG. 19, Satlite 1 and Satlite 2 are orbiting back and forth between the northern and southern hemispheres alternately. Specifically, Satelite1 flies in the northern hemisphere (monitors the northern hemisphere) and Satelite2 flies in the southern hemisphere (monitors the southern hemisphere) for 12 hours a day, while Satelite1 flies in the southern hemisphere for the next 12 hours. (Monitoring the southern hemisphere) and Satellite 2 flies in the northern hemisphere (monitors the northern hemisphere).

http://www.lockheedmartin.com/products/SpaceBasedInfraredSystemHigh/ (Lockheed Martin website) http://en.wikipedia.org/wiki/Space-Based_Infrared_System ("Space-Based Infrared System" section from Wikipedia, USA)

  However, in this method, since one artificial satellite orbits the northern hemisphere and the southern hemisphere, the area where two aircraft can simultaneously monitor (stereoscopic view) is only the hatched portion of FIG. In other words, although the presence of a flying object can be confirmed by monitoring one artificial satellite, the area in which the three-dimensional position of the flying object can be accurately identified by stereo vision using two artificial satellites and the flying direction are accurately determined. The area that can be predicted is limited.

  In addition, there is an attempt to secure a wide monitoring area using a plurality of geostationary satellites (Satelite 3 and Satelite 4) as shown in FIG. In this case, the hatched area shown in FIG.

  However, since the position of the geostationary satellite is easily specified in advance, it is easy for a third party to launch the flying object outside the monitoring area of the geostationary satellite. In addition, since geostationary satellites usually fly on the equator, the location of high latitudes on the earth becomes an image through clouds even within the area that can be monitored. There is. In addition, although the area that can be positioned by stereo vision is wide as described above, since the geostationary satellite is located on the equator, there is a problem that the position specifying accuracy becomes worse as the latitude becomes higher.

  Alternatively, a method of monitoring with a plurality of low-altitude satellites is also conceivable. However, because the altitude is low, the monitoring area is extremely narrow and remains at most 500 km square, and it can detect a flying object from a wide monitoring area. To do it quickly, dozens of low-altitude satellites are required. Furthermore, in order to accurately estimate the three-dimensional position of the flying object by observing with these low-altitude satellites in stereo, the number of low-altitude satellites that are double the accuracy is required. It is not realistic in thinking about construction.

  Accordingly, an object of the present invention is to solve various problems in the prior art. That is, the present invention is a flying object detection system composed of a plurality of artificial satellites, and the three-dimensional flying object is obtained by viewing the data (for example, captured image data) output from the observer in stereo (stereoscopic). It is an object of the present invention to provide a flying object detection method and system for accurately estimating or specifying a position.

The present invention includes a step of photographing a flying object flying on the ground by at least two artificial satellites and transmitting it as photographed image data, and position information, velocity information of each artificial satellite by the at least two artificial satellites, Transmitting at least two pieces of image data transmitted from the at least two artificial satellites by a base station equipped with a receiving antenna, position information, velocity information, and Receiving the attitude information, and calculating the vertex angle and azimuth angle of each artificial satellite by analyzing at least two pieces of image data received by the antenna by the position analysis device in the base station; A flying object including the step of analyzing the position of the flying object using the apex angle and the azimuth angle, and the position information, velocity information, and attitude information of each artificial satellite. A method of generating an observation vector from an apex angle and an azimuth of each artificial satellite, and generating a coordinate transformation matrix based on position information, velocity information, and attitude information of each artificial satellite A step of generating a coefficient from the initial position value of the flying object, position information of the artificial satellite, and the coordinate transformation matrix when the initial position value of the flying object is X _ref , the coefficient, the observation vector, And a step of specifying a flying object position by executing a parameter estimation process defined by:

In addition, the present invention provides at least two apparatuses including an observer capable of photographing a flying object flying on the ground, and a transmitter that transmits the photographed image data and each position information, speed information, and posture information. An artificial satellite, an antenna for receiving at least two pieces of image data transmitted from the at least two artificial satellites, position information, velocity information, and attitude information of each artificial satellite; and reception by the antenna A flying object detection comprising: a base station comprising a position analysis device that analyzes the position of the flying object based on the at least two pieces of image data and the position information, velocity information, and attitude information of each of the image data In the system, the position analysis device analyzes the at least two pieces of image data to calculate a vertex angle and an azimuth angle of each artificial satellite, and each of the artificial satellites. An observation vector generation unit that generates an observation vector from an angle and an azimuth, a coordinate conversion matrix generation unit that generates a coordinate conversion matrix based on the position information, velocity information, and attitude information of each artificial satellite, and a flying object position When the initial value is X _ref , a parameter defined by the coefficient generation unit that generates a coefficient from the initial position value of the flying object, the position information of the artificial satellite, and the coordinate transformation matrix, the coefficient, and the observation vector And a parameter estimation processing unit that executes the estimation process.

  According to the flying object detection method and system and the like according to the present invention, the flying object shown in a plurality of imaging data transmitted from a plurality of artificial satellites is analyzed stereoscopically, compared with the prior art. The three-dimensional position of the flying object can be specified accurately. In addition, by continuously receiving and analyzing the imaging data at a predetermined time interval (Δt), it is possible to estimate not only the three-dimensional position of the flying object but also the velocity vector, so that the flying object may be accidentally dropped. In the event of an accident such as an explosion or explosion, it is possible to accurately predict the crash point and landing point of the flying object, and prevent or minimize the damage at the crash point and landing point.

It is explanatory drawing explaining the structural example of the flying object detection system concerning one Embodiment of this invention. It is explanatory drawing explaining the structural example of the position analysis apparatus in the flying object detection system concerning one Embodiment of this invention. It is explanatory drawing explaining the example of the monitoring area | region of the flying body detection method and system concerning one Embodiment of this invention. It is explanatory drawing explaining a mode that the earth was looked down from the artificial satellite in the flying body detection method and system concerning one Embodiment of this invention. It is explanatory drawing explaining the satellite body coordinate system in the flying body detection method and system concerning one Embodiment of this invention. It is explanatory drawing explaining the satellite orbit coordinate system (O) and the geocentric direction orbit coordinate system (EDO) in the flying object detection method and system concerning one Embodiment of this invention. It is explanatory drawing explaining the relationship between the imaging screen imaged in the flying body detection method and system concerning one Embodiment of this invention, and a geocentric direction orbit coordinate system (EDO). It is explanatory drawing explaining the picked-up image imaged from two artificial satellites in the flying object detection method and system concerning one Embodiment of this invention. It is explanatory drawing explaining the relationship between the satellite body coordinate system (S) and a flying object vector in the flying object detection method and system concerning one Embodiment of this invention. It is explanatory drawing explaining the relationship between the earth center fixed coordinate system (ECEF) and the NED (North-East-Down) coordinate system in the flying object detection method and system concerning one Embodiment of this invention. It is explanatory drawing explaining the relationship between the captured image imaged in the flying object detection method and system concerning one Embodiment of this invention, and a NED (North-East-Down) coordinate system. It is explanatory drawing explaining the picked-up image imaged from three artificial satellites in the flying body detection method and system concerning other embodiment of this invention. It is explanatory drawing explaining the coordinate system conversion theory in the position-analysis apparatus of the flying object detection system concerning one Embodiment (when 2 satellites and satellite attitude information exist) of this invention. It is explanatory drawing explaining the coordinate system conversion theory in the position analysis apparatus of the flying object detection system concerning other embodiment (in the case of 2 satellites, satellite attitude information absence) of this invention. It is explanatory drawing explaining the coordinate system conversion theory in the position analyzer of the flying object detection system concerning other embodiment (in the case of 3 satellites) of this invention. It is explanatory drawing explaining the arithmetic processing flow in the position-analysis apparatus of the flying object detection system concerning one Embodiment (when 2 satellites and satellite attitude information exist) of this invention. It is explanatory drawing explaining the arithmetic processing flow in the position analysis apparatus of the flying object detection system concerning other embodiment (in the case of 2 satellites, satellite attitude information absence) of this invention. It is explanatory drawing explaining the arithmetic processing flow in the position-analysis apparatus of the flying object detection system concerning other embodiment (in the case of 3 satellites) of this invention. It is explanatory drawing explaining the example of the monitoring area | region in the conventional monitoring system which uses two artificial satellites. It is explanatory drawing explaining the example of the monitoring area | region in the conventional monitoring system using two geostationary satellites.

  Hereinafter, embodiments for carrying out the flying object detection method and system according to the present invention will be described in detail.

[System configuration example]
FIG. 1 shows a configuration example of a flying object detection system according to an embodiment of the present invention. The system 100 captures images of target flying objects transmitted from two (or three) artificial satellites 101a to 101b (101c) and artificial satellites 101a to 101b (101c) that orbit the earth. And a base station 102 provided with an antenna 103 for receiving various data such as the image data and the position and orientation of each artificial satellite itself. The base station 102 is provided with a position analysis device 105 for analyzing the image data received by the antenna 103 and analyzing the position of the target flying object. The antenna 103 and the position analysis device 105 are They are connected by a wired signal line 104. Note that part or all of the wired signal line 104 may be a wireless interface (including a wireless LAN, for example). Further, although the base station 102 has been described as a building, it may be a moving body such as an aircraft, a ship, or a vehicle. Alternatively, a portable position analysis device 105 equipped with a small antenna 103 may be carried (in this case, the position analysis device 105 provided with the antenna 103 itself becomes a base station).

  The artificial satellite 101a to the artificial satellite 101c include an observer (not shown) for observing the surface of the earth and photographing the image data, and the captured image data and the position information and attitude information of the artificial satellite itself on the earth. At least a transmitter (not shown) for transmitting to the base station is provided.

  FIG. 2 shows a configuration example of the position analysis apparatus in the flying object detection system according to the embodiment of the present invention of the position analysis apparatus 105. The position analysis device 105 receives image data of a target flying object transmitted from two or three artificial satellites orbiting the earth and various data such as the position and attitude of the artificial satellite itself. The antenna 103 and the receiving unit 201 are connected via a wired or wireless interface 104.

  In addition, the position analysis device 105 includes a communication interface 1051 for transmitting signals received by the antenna 103 and the reception unit 201 to the device main body as data, and a central processing for performing various numerical calculations, logical operations, device control, and the like by a program A device (CPU) 1052, a memory 1053 configured with a RAM and / or a ROM, and a storage disk 1054 configured with a magnetic or optical medium or the like are included. The software and / or program necessary for implementing the present invention is normally installed or stored in the storage disk 1055, and when executing the software and / or program, all or a part thereof is stored in the RAM or the like in the memory 1053 as necessary. These software modules are read out and executed by the CPU 1052.

  Furthermore, if necessary, a display unit 1056 configured with a digital signal processor (DSP) 1054 designed to process a specific matrix operation or the like at high speed, a display device for displaying processing results, etc., a keyboard, You may comprise so that the input part 1057 which consists of readers, such as a mouse | mouth, a pen, OCR, optics, or a magnet, may be included.

  Although not shown in FIG. 2, each unit in the position analysis device 105 is appropriately connected by a signal line such as a bus. Further, the position analysis device 105 can be installed with a program capable of executing calculation processing or logical operation processing according to the present invention, which will be described later, and has an interface that can be appropriately connected to the antenna 103 and the receiving unit 201. For example, a general-purpose computer or a portable information terminal (including a device having a communication function such as a cellular phone) can be used as the position analysis device 105.

  Further, in the present specification, the flying object detection method and system according to the present invention will be described based on some embodiments. However, each of calculations and operations including those exemplified by mathematical formulas are mainly positions. This is based on arithmetic processing executed by the analysis apparatus 105 and the computer program in cooperation. In the following, detailed description will be given of mounting processing based on the theoretical technical background (for example, various orbital systems for satellite orbits) on which the flying object detection method and system according to the present invention is premised.

  FIG. 3 shows a configuration of an artificial satellite (and an example of a circular orbit) implemented in a flying object detection method and system according to an embodiment of the present invention. Artificial satellites for monitoring the earth's surface are composed of three aircraft (Satelite 1 to Satelite 3), which is a flying object detection system mainly using the northern hemisphere as a monitoring area. Although three artificial satellites are used in the figure, only two artificial satellites (for example, Satelite 1 and Satelite 2) can be used. Alternatively, two of the three artificial satellites shown in FIG. 3 can be appropriately switched and used.

  The orbits of each of the three satellites from Satelite 1 to Satelite 3 are equally arranged with an eccentricity of 0.075, an orbital inclination angle of 45 degrees, and an ascending intersection red longitude of 120 degrees. The near point separation angle is 270 degrees.

  In order to set the northern hemisphere as a monitoring region, it is desirable that the trajectory fly at a relatively slow speed (relative to the southern hemisphere) over the northern hemisphere. Based on this, in one embodiment, the eccentricity is set to 0.75 and the near point separation angle is 270 degrees. Here, if the orbital inclination angle is set to 45 degrees, the artificial satellite also flies over Japan with high latitudes, etc., and flies over Japan as shown in FIG. In the case of a satellite), it is possible to monitor the area beyond the North Pole.

  Furthermore, in order to monitor the northern hemisphere using three artificial satellites, it is necessary to adjust the relationship between each other's orbits so that the three aircraft sequentially fly over the northern hemisphere. In this case, in one embodiment, it is desirable that the ascending intersection red diameters are equally arranged by 120 degrees.

  According to the configuration provided with the above conditions, the range that can be viewed in stereo by two of the three artificial satellites is as shown in the hatched portion of FIG. 3 and is wide.

[Earth center fixed coordinate system definition]
Next, the earth fixed coordinate system (ECEF) will be described. The Earth Center Fixed Coordinate System (ECEF) has the Earth's center of gravity as the origin, the z-axis passes through a conventional earth pole fixed to the earth's crust, the x-axis passes through the intersection of the equatorial plane and the average Greenwich meridian, and the y-axis Is a coordinate system that forms a right-handed system with these x and z axes.

[Satellite body coordinate system definition]
Next, the coordinate system (S) of the satellite aircraft will be described with reference to FIG. Now, the position vector of the satellite as seen from the origin of the earth fixed coordinate system (the origin of the xyz coordinates shown in the figure)
And the position vector of the flying object as seen from the center of the earth
And the vector from the flying object to the satellite
And At this time, the relationship of r ^S , the satellite position vector R ^S, and the flying object position vector R ⁰ is expressed by the following equation.

・ ・ ・ ・ ・ ・ (Formula 1)

The unit vector in the direction of the solar battery panel of the satellite
The unit vector in the sensor axis direction of the satellite
The unit vector that forms the right-handed system with YS and PS is
The flying object apex angle is
It is represented by

[Satellite orbit coordinate system definition]
Next, the satellite orbit coordinate system (O) will be described with reference to FIG. If there is a satellite S that orbits the elliptical orbit at the position shown in FIG. 6A, the orbital coordinate system (O) in the satellite S is a unit vector in the velocity vector direction in the orbital plane.
Is the unit vector orthogonal to the x-axis, the velocity vector direction and the orbital out-of-plane direction vector orthogonal thereto
Is a coordinate system with the y-axis and the x- and y-axes and the axis that forms the right hand system as the z-axis.

[Derivation of observation equations]
Next, FIG. 8 shows the state of a captured image when observation devices mounted on two artificial satellites simultaneously observe. FIG. 8A is an image of the earth taken by one artificial satellite, and FIG. 8B is an image of the earth taken by the other artificial satellite. 8 (a) and 8 (b)
Indicates a flying object. In principle, a stereo view (stereoscopic view) with respect to the earth is realized by analyzing the two photographed image data. Here, the observation amount directly obtained from the imaging data shown in FIGS. 8A and 8B is an angle (vertical angle) from the image center.
And the image rotation angle (azimuth angle) based on the vertical and horizontal directions of the image.
Are four parameters. At this time, as shown in FIG. 9, the position of the flying object viewed from the satellite body coordinate system (S) is expressed by the following equation using these two angles and the vector r ^S from the flying object to the satellite. The

・ ・ ・ ・ ・ ・ (Formula 2)

Here, | _S attached to the vector r ^S means that the basis of the vector is the satellite body coordinate system (S).

The vector to be finally obtained is the three-dimensional position vector of the flying object in the Earth Center Fixed Coordinate System (ECEF) with the center of the Earth as the origin.
However, R ⁰ is an image of each satellite based on the relationship between the coordinate transformation matrix C _S2ECEF from the satellite body coordinate system (S) to the earth center fixed coordinate system (ECEF) and (Equation 1) and (Equation 2). Observables in images
Is obtained by solving the observation equation shown by the following equation.

・ ・ ・ ・ ・ ・ (Formula 3)

<Coordinate system conversion with 2 satellites and satellite attitude information>
Here, the coordinate transformation matrix C _S2ECEF from the satellite body coordinate system (S) to the earth center fixed coordinate system (ECEF) is performed via the satellite orbit coordinate system (O) as shown in FIG.

Known roll attitude angle for each satellite
Pitch attitude angle
Yaw attitude angle
Transformation matrix from satellite body coordinate system (S) to orbital coordinate system (O) based on
And coordinates from the orbital coordinate system (O) to the earth center fixed coordinate system (ECEF) based on the known three-dimensional vector R ^S (s = 1, 2) and velocity vector v ^S (s = 1, 2) of each satellite. Transformation matrix
Thus, a coordinate transformation matrix C _S2ECEF from the satellite body coordinate system (S) to the earth center fixed coordinate system (ECEF) is determined as follows.

・ ・ ・ ・ ・ ・ (Formula 4)

・ ・ ・ ・ ・ ・ (Formula 5)

Here, C ₁ (θ _R ), C ₂ (θ _P ), and C ₃ (θ _Y ) are respectively

・ ・ ・ ・ ・ ・ (Formula 6a)
・ ・ ・ ・ ・ ・ (Formula 6b)
... (Formula 6c)
It is represented by

And C _O2ECEF is

・ ・ ・ ・ ・ ・ (Formula 7)
In ^{_{expressed, e s xO, e s yO}} , for e ^s _zO, respectively,

・ ・ ・ ・ ・ ・ (Formula 8a)
... (Formula 8b)
... (Formula 8c)

It is represented by The operator x means a vector product.

<Calculation process flow with 2 satellites and satellite attitude information>
Further, an example of arithmetic processing based on the above theory and the characteristics of the computer will be described in detail. In other words, since the observation equation (Equation 3) is a nonlinear equation, in order to obtain an operational solution, the first-order Taylor series around x _ref (which is the initial value of the flying object position) for the variable x It is realistic to obtain a solution by performing linearization by approximation and repeatedly performing arithmetic processing.

[Solution of general nonlinear equations]
Now, when a general nonlinear equation as shown in (Equation 9a) is given, when approximation is performed by first-order Taylor series expansion, (Equation 9b) is obtained, and when (Equation 9b) is solved for the variable x, (Equation 9c) ) The specific calculation of the variable x is obtained by giving an initial value x _ref as shown in (Equation 9d) and performing iterative calculation until it converges to a predetermined range.

... (Formula 9a)
... (Formula 9b)
... (Formula 9c)

... (Formula 9d)
here,

[Aircraft position calculation algorithm]
Therefore, when the solution of the nonlinear observation equation represented by (Equation 3) is applied to the flying object detection method or system according to the present invention, the processing flow shown in FIG. 16 is obtained. Hereinafter, a description will be given based on FIG.

  First, processing is started in S1601. Specifically, the flying object apex angle and azimuth angle are obtained from the captured images obtained from two of the satellites 101a to 101c by observation amount calculation processing. In other words, the captured image data transmitted from the artificial satellite is received by the antenna 103 and the receiving unit 201 and recorded in the memory 1053 or the disk 1055 via the communication interface 1051 under the control of the CPU 1052 in the position analysis device 105. The captured image is analyzed and the flying object apex angle and azimuth angle are obtained (S1602). At this time, a program necessary for the processing is read from the disk 1055 to the memory 1053 as appropriate and executed. In addition, attitude information, velocity information, and three-dimensional position information of each satellite are transmitted from each satellite, received by the antenna 103 and the receiving unit 201 as appropriate, and taken into the position analysis device 105 (S1605, S1606, S1607). .

  Next, in S1603, the observation vector y is calculated and generated from the flying object apex angle and azimuth by the observation amount calculation process represented by the following equation (S1604).

・ ・ ・ ・ ・ ・ (Formula 10)

  Next, in S1608, the coordinate transformation matrix shown in Equations 4, 5, and 7 is calculated from the position, velocity, and attitude angle of each satellite.

Next, an initial value of the flying object position (see x _ref related to the variable x described above) is determined (S1609), and the coordinate transformation matrix (Equations 4, 5, and 7), the position of the satellite, and the initial value of the flying object are determined. Each coefficient represented by the following equation is calculated (S1610).

... (Formula 11a)

... (Formula 11b)

... (Formula 11c)

  Next, the position vector of the flying object finally obtained is calculated by repeatedly calculating the following equation using the obtained coefficient.

・ ・ ・ ・ ・ ・ (Formula 12)

  That is, when the parameter estimation process shown in (Equation 12) is executed in S1611, the flying object position is calculated each time (S1612), and it is determined whether this position is within a predetermined range (S1613). . If it is not within the predetermined range (No in S1613), the process returns to S1610 to perform coefficient generation processing. If it is determined in S1613 that it is within the predetermined range (Yes), the process proceeds to S1614 and the process is terminated.

  The coefficient shown in (Expression 11c) is a partial differentiation of the function f (•) by x, and is derived from the following expression.

... (Formula 13)
Where L and I are

... (Formula 14a)
... (Formula 14b)
It is represented by

[Calculation of flying direction (speed) of flying object]
The flying direction (velocity) of the flying object is the position vector of the flying object at time t obtained from (Equation 12), where Δt is the observation time interval of the observation device mounted on the artificial satellite.
And the position vector of the flying object at time t-Δt
And is calculated from the following equation.

・ ・ ・ ・ ・ ・ (Formula 15)

  As a result of the processing as described above, the position and velocity direction (vector) of the flying object can be specified from the satellite image (satellite photograph) obtained by shooting the flying object.

  In the first embodiment, it is assumed that the observation times of two artificial satellites are “simultaneous”. However, when an estimation algorithm incorporating a dynamic state transition such as a Kalman filter is used, the following is shown. In addition, the observation times do not necessarily have to be “simultaneous”. This will be described in detail below.

In general, the Kalman filter has two processing configurations of observation update and time propagation, and each is represented by the following equation.
<Observation update processing>

・ ・ ・ ・ ・ ・ (Formula 16)

  Here, f (•) is an observation function, x is a state variable, y is an observation amount, P is a covariance matrix of the state variable, R is a variance of the observation amount, and H is a partial differentiation of the observation function by the state variable. The subscript i indicates that it is at the i-th time, and the superscript indicates that “−” is before update and “+” is after update.

<Time propagation processing>
・ ・ ・ ・ ・ ・ (Formula 17)

Here, g (•) is a time propagation function, Φ _i is a state transition matrix of the state variable x, and Q is a process noise of the state variable x.

  When the Kalman filter described above is applied to the flying object detection method or system according to the present invention, the elements such as the observation function f (•), the state variable x, and the observation amount y may be as follows.

... (Formula 18a)
... (Formula 18b)
... (Formula 18c)
... (Formula 18d)
・ ・ ・ ・ ・ ・ (Formula 18e)
... (Formula 18f)
・ ・ ・ ・ ・ ・ (Formula 18g)
・ ・ ・ ・ ・ ・ (Formula 18h)

  In the first and second embodiments, assuming that the attitude angle of the artificial satellite is known, coordinate transformation from the satellite body coordinate system (S) to the earth center fixed coordinate system (ECEF) is performed in the orbit coordinate system (O However, the present invention is not limited to this. In other words, even if the posture angle is not known, the present invention can be implemented, and for example, arithmetic processing may be performed via an alternative coordinate system (such as a NED coordinate system described later).

  In addition, even if the attitude angle of the satellite is known or unknown, the attitude of the satellite does not need to be precisely controlled, and the flying object is not captured in the image taken from the satellite. It is enough if it is included. For example, the center of the image in the captured image does not necessarily need to be centered on the screen, and even if the center of the earth protrudes greatly from the captured image, the center position can be estimated and calculated by image processing. . This means that the roll angle and pitch angle of the artificial satellite can be estimated from the image.

For example, the captured image may be rotated as long as a part of the image of the earth is captured. In this case, the north direction and the east direction of the earth can be estimated by a pattern matching processing technique between the observed image data and the terrain data on the map. This means that the angle corresponding to the yaw angle of the artificial satellite can be estimated from the image. More specifically, the center of the earth and the north and east directions of the earth are specified by image processing and correspond to the roll angle, pitch angle, and yaw angle of each satellite.
Can be estimated, a coordinate transformation from the satellite body coordinate system (S) to the earth center fixed coordinate system (ECEF) may be calculated via a NED (North-east-down) coordinate system shown in FIG.

[Definition of NED coordinate system and relationship with fixed center of earth]
In FIG. 10, the position of the artificial satellite S in the xyz coordinate system (ECEF coordinate system) with the earth center O as the origin (azimuth angle λ ^s rotated from the x axis around the z axis, from the xy plane) If it is at an elevation angle θ ^s ), the NED coordinate system is represented by three coordinate axes N, E, and D with the satellite S as the origin.

[Relationship between captured image and NED coordinate system]
Next, the relationship between the NED coordinate system on the imaging screen and the satellite roll angle, pitch angle, and yaw angle will be described with reference to FIG. In the drawing, the center O _{1 of the} imaging screen and the center O _{2 of the} earth are slightly shifted. And the pitch axis in the satellite body coordinate system
Matches the abscissa (X axis) with O ₂ as the origin, and the roll axis in the satellite body coordinate system
Is coincident with the ordinate (Y axis) with O ₂ as the origin.

Also, if the east axis in the NED coordinate system is E and the north axis in the coordinate system is N,
The pitch angle of the satellite
The roll angle of the satellite
The yaw angle of the satellite
It becomes. A clockwise arrow shown at the lower right in the figure indicates a positive direction when the angle is defined.

  In FIG. 11, it is assumed that the whole earth is contained in the image, but as described above, the present invention is not limited to this. It is only necessary that a part of the image of the earth is accommodated, and the center, north direction, and east direction of the earth can be estimated even if the image is a partial image of the earth.

<2 satellites, coordinate system conversion without satellite attitude information>
The coordinate transformation matrix C _S2ECEF from the satellite body coordinate system (S) to the earth center fixed coordinate system (ECEF) is performed via the NED coordinate system as shown in FIG.

Coordinate transformation matrix from satellite body coordinate system (S) to NED coordinate system
And a known 3D position vector for each satellite
・ ・ ・ ・ ・ ・ (Formula 19)
Based on the coordinate transformation matrix C _NED2ECEF (R ^S ) from the NED coordinate system to the earth center fixed coordinate system (ECEF), the coordinate transformation matrix C from the satellite body coordinate system (S) to the earth center fixed coordinate system (ECEF) _S2ECEF is determined as follows.

... (Formula 20a)
.... (Formula 20b)
here,

... (Formula 21a)
... (Formula 21b)
... (Formula 21c)
... (Formula 21d)
It is. In addition,

... (Formula 22a)
... (Formula 22b)
It becomes.

<Calculation process flow for 2 satellites and no satellite attitude information>
A specific calculation processing example in this embodiment will be described.

[Aircraft position calculation algorithm]
Next, FIG. 17 shows a processing flow based on the realistic calculation processing of the nonlinear observation equation represented by (Equation 3). Hereinafter, a description will be given with reference to FIG.

  First, processing is started in S1701. Specifically, the flying object apex angle and azimuth angle are obtained from the captured images obtained from two of the satellites 101a to 101c by observation amount calculation processing. In other words, the transmitted image is received by the antenna 103 and the receiving unit 201 and recorded in the memory 1053 or the disk 1055 via the communication interface 1051 under the control of the CPU 1052 in the position analysis device 105, and the captured image is analyzed. Thus, the flying object apex angle and azimuth angle are obtained (S1702). At this time, a program necessary for the processing is read from the disk 1055 to the memory 1053 as appropriate and executed.

  Next, in S1703, an observation vector y is calculated and generated from the flying object apex angle and azimuth by an observation amount calculation process represented by the following equation (S1704).

... (Formula 23)

Next, in S1705, the terrain shape of the earth surface in the captured image is compared with the terrain shape in other map data by the pattern matching processing technique, and the attitude difference (roll angle, pitch angle, yaw angle) of the artificial satellite itself is calculated. calculate. As a result, the roll angle, pitch angle, and yaw angle of the two satellites
Is obtained and stored in the memory 1053 (S1706).

  In addition, the velocity information and the three-dimensional position information of each satellite are transmitted from each satellite, and are appropriately received by the antenna 103 and the receiving unit 201 and taken into the position analysis device 105 (S1707, S1708).

  Next, in S1709, the coordinate transformation matrix shown in (Expressions 20a, 20b, 21d) is calculated from the position, velocity, and attitude angle of each satellite.

Next, an initial value of the flying object position (see x _ref related to the variable x described above) is determined (S1710), and the coordinate transformation matrix (Equations 20a, 20b, 21d), the position of the satellite, and the initial value of the flying object are determined. Each coefficient indicated by the following equation is calculated (S1711).

... (Formula 24a)
... (Formula 24b)
... (Formula 24c)

  Next, the position vector of the flying object finally obtained is calculated by repeatedly calculating the following equation using the obtained coefficient.

・ ・ ・ ・ ・ ・ (Formula 25)

  That is, when the parameter estimation process shown in (Equation 25) is executed in S1712, the flying object position is calculated each time (S1713), and it is determined whether this position is within a predetermined range (S1714). . If it is not within the predetermined range (No in S1714), the process returns to S1711 to perform coefficient generation processing. If it is determined in S1714 that it is within the predetermined range (Yes), the process proceeds to S1715 and the process is terminated.

[Calculation of flying direction (speed) of flying object]
The flying direction (velocity) of the flying object is the position vector of the flying object at time t obtained from (Equation 25), where Δt is the observation time interval of the observation device mounted on the artificial satellite.
And the position vector of the flying object at time t-Δt
And is calculated from the following equation.

・ ・ ・ ・ ・ ・ (Formula 31)

  As a result of the processing as described above, the position and velocity direction (vector) of the flying object can be specified from the satellite image (satellite photograph) obtained by shooting the flying object.

  In the first to third embodiments, the three-dimensional position and flight direction of the flying object using the stereo vision (stereoscopic) algorithm as described above for two pieces of image data observed by two artificial satellites. An example of calculating (speed) has been described. In the fourth embodiment, a method for calculating a three-dimensional position and a flight direction (speed) of a flying object using three pieces of image data observed by three artificial satellites will be described.

One of the characteristics when using three satellites is the yaw angle attitude around the yaw axis of the satellite.
Can also be estimated at the same time. Therefore, even if the roll angle, pitch angle, and yaw angle attitude of the satellite are unknown, or even if the center of the earth cannot be accurately grasped without accurately performing image pattern matching, the center of the earth is photographed. If it can be estimated from the image, the three-dimensional position of the target flying object can be specified. FIG. 12 shows image data observed from the observation equipment mounted on these three artificial satellites.

In this case, the directly obtained observation amount is an angle (vertical angle) from the center of the captured image.
And the image rotation angle (azimuth angle) based on the vertical and horizontal directions of the image.
Are six parameters.

At this time, the three-dimensional position of the target flying object
... (Formula 32)
As described in the case of two artificial satellites, the coordinate transformation matrix C _S2ECEF from the satellite body coordinate system (S) to the earth center fixed coordinate system (ECEF) and the above (formula 1) and (formula 2) Therefore, the amount of observation on the imaging screen of each satellite
The following relationship can be derived:

・ ・ ・ ・ ・ ・ (Formula 33)

  When three artificial satellites are used and the yaw angle attitude around the yaw axis of the artificial satellite is estimated, the coordinate conversion from the satellite body coordinate system (S) to the earth center fixed coordinate system (ECEF) This is done via the directional orbit coordinate system (EDO).

[Definition of satellite orbital coordinate system (EDO)]
First, the geocentric orbit coordinate system (EDO) will be described based on the relationship between the satellite orbit coordinate system (O) and the satellite's geocentric orbit coordinate system (EDO) based on FIG. Generally, since the eccentricity of the satellite orbit is not zero (that is, the satellite orbit is not a perfect circle), the satellite orbit coordinate system (O) and the satellite's geocentric orbit coordinate system (EDO) are different.

Now, assuming that there is a satellite S orbiting the elliptical orbit at the position shown in FIG. 6B, the satellite S in the orbital coordinate system (O) is a unit vector in the velocity vector direction in the orbital plane.
And the unit vector orthogonal to the velocity vector direction and the out-of-plane direction vector orthogonal to the direction
It is expressed as On the other hand, the satellite S in the geocentric orbit coordinate system (EDO) is a unit vector in the geocentric direction.
And a unit vector perpendicular to the geocentric direction and in the orbital plane
It is expressed as

[Relationship between captured image and geocentric orbit coordinate system (EDO)]
Next, FIG. 7 shows the relationship between the captured image and the geocentric orbit coordinate system (EDO). Now, if the rectangular outer frame in the figure is the outer frame of the imaging screen, the center O _{1 of the} imaging screen is slightly shifted from the center O _{2 of the} earth. And the pitch axis in the satellite body coordinate system
Matches the abscissa (X axis) with O ₂ as the origin, and the roll axis in the satellite body coordinate system
Is the same as the ordinate (Y axis) with O ₂ as the origin, and the unit vector projected onto the RS-PS plane in the vector direction perpendicular to the geocentric direction and in the orbital plane is
And the unit vector projected onto the RS-PS plane in the unit vector direction forming the right-handed system with xEDO and zEDO
given that,
The pitch angle of the satellite
The roll angle of the satellite
The yaw angle of the satellite
It becomes. In addition, the clockwise arrow shown at the lower right in the figure indicates a positive direction when the angle is defined.

<Coordinate system conversion for 3 satellites>
The coordinate transformation matrix C _S2ECEF from the satellite body coordinate system (S) to the earth center fixed coordinate system (ECEF) is performed via the geocentric orbit coordinate system as shown in FIG.

Coordinate transformation matrix C _S2EDO from satellite body coordinate system (S) to geocentric orbit coordinate system (EDO) and coordinate transformation matrix C from geocentric orbit coordinate system (EDO) to earth center fixed coordinate system (ECEF) It is expressed by the following formula from _EDO2ECEF .

・ ・ ・ ・ ・ ・ (Formula 34a)
... (Formula 34b)
here,
... (Formula 34c)
... (Formula 34d)
... (Formula 34e)
It is.

The coordinate transformation matrix C _S2EDO from the satellite body coordinate system (S) to the geocentric orbit coordinate system (EDO) is the roll angle and pitch angle of each satellite obtained by estimating the earth center from the image.
And the yaw angle treated as an estimated parameter
Based on the above, it is expressed by the following equation using a rotation matrix around each axis.

・ ・ ・ ・ ・ ・ (Formula 35)
here,

・ ・ ・ ・ ・ ・ (Formula 36a)
... (Formula 36b)
... (Formula 36c)
It is.

[Derivation of observation equations]
From the above, when estimating the yaw axis of each satellite by stereo viewing (stereoscopic viewing) of the three satellites, the observation equation shown below is solved.
... (Formula 37)
3D position of the flying object
... (Formula 38)
And the yaw axis of each satellite
Can be calculated.

<Calculation processing flow for 3 satellites>
A specific calculation processing example in this embodiment will be described.

[Aircraft position calculation algorithm]
Next, the processing flow based on the realistic processing of the nonlinear observation equation represented by (Equation 37) is as shown in FIG. Hereinafter, a description will be given based on FIG.

  First, the process is started in S1801. Specifically, the flying object apex angle and azimuth angle are obtained by the observation amount calculation process from the captured images obtained from the three satellites 101a to 101c. In other words, the transmitted image is received by the antenna 103 and the receiving unit 201 and recorded in the memory 1053 or the disk 1055 via the communication interface 1051 under the control of the CPU 1052 in the position analysis device 105, and the captured image is analyzed. Thus, the flying object apex angle and azimuth angle are obtained (S1802). At this time, a program necessary for the processing is read from the disk 1055 to the memory 1053 as appropriate and executed.

Next, in S1803, an observation vector y is calculated and generated from the flying object apex angle and azimuth by an observation amount calculation process represented by the following equation (S1804).
・ ・ ・ ・ ・ ・ (Formula 39)

Next, in S1805, the terrain shape of the earth surface in the captured image is compared with the terrain shape in other map data using a pattern matching processing technique, and the attitude difference (roll angle, pitch angle) of the artificial satellite itself is calculated. As a result, the roll angle and pitch angle of the three satellites
Is obtained and stored in the memory 1053 (S1806).

  In addition, speed information and three-dimensional position information of each satellite are transmitted from each satellite, and are appropriately received by the antenna 103 and the receiving unit 201 and taken into the position analysis apparatus 105 (S1807, S1808).

  Next, in S1809, the coordinate transformation matrix shown in (Expressions 34b, 36b, 36c) is calculated from the position, velocity, and attitude angle of each satellite.

Next, an initial value of the flying object position (see x _ref related to the variable x described above) is determined (S1810), and the coordinate transformation matrix (Equations 34b, 36b, 36c), the position of the satellite, and the initial value of the flying object are determined. Each coefficient indicated by the following equation is calculated (S1811).
... (Formula 40a)
.. (Formula 40b)
・ ・ ・ ・ ・ ・ (Formula 40c)
... (Formula 40d)
.. (Formula 40e)

  The position vector of the flying object finally obtained is calculated by repeatedly calculating the following equation using the obtained coefficient.

... (Formula 41)

  That is, when the parameter estimation process shown in (Equation 41) is executed in S1812, the flying object position and satellite attitude are calculated each time (S1813), and it is determined whether this position is within a predetermined range. (S1814). If it is not within the predetermined range (No in S1814), the process returns to S1811 to perform coefficient generation processing. If it is determined in S1814 that it is within the predetermined range (Yes), the process proceeds to S1815 and the process is terminated.

The coefficients shown in (Equations 40c, 40d, 40e) are the three-dimensional position of the flying object which is the estimation parameter of the function f (•).
... (Formula 42)
By partial differential and yaw axis of each satellite
And is derived from the following equation.

・ ・ ・ ・ ・ ・ (Formula 43)
here,

... (Formula 44)
And

... (Formula 45a)
... (Formula 45b)
... (Formula 45c)
It is. Also,

... (Formula 46)
It is.

[Calculation of flying direction (speed) of flying object]
The flying direction (velocity) of the flying object is the position vector of the flying object at time t obtained from (Equation 41), where Δt is the observation time interval of the observation device mounted on the artificial satellite.
And the position vector of the flying object at time t-Δt
And is calculated from the following equation.

・ ・ ・ ・ ・ ・ (Formula 47)

  As described above, the position and velocity vector of the flying object can be specified from images taken from two artificial satellites (Examples 1 to 3), and images taken from three artificial satellites are used. Then, it becomes possible to further estimate the yaw angle posture (Example 4).

  In addition, the present invention is not necessarily completed only when the entire system shown in FIG. 1 is configured. Images taken by two or three artificial satellites (and, if necessary, positional information of the artificial satellites taken). Etc.) It is also possible to use a program that uses only data and specifies the three-dimensional position and velocity vector of the target flying object, or a computer that executes such a program.

  It should be noted that all technical elements and methods or processing steps described in the claims, the description, the abstract, the drawings, and the method or processing step are at least of these elements and / or steps. Arbitrary combinations other than combinations that are mutually exclusive can be components of the method, system, and program according to the present invention.

  In addition, this invention is not restrict | limited to any individual specific detailed description of embodiment mentioned above. Further, the technical scope of the present invention is not limited only to the above description, but the extension thereof is specified by the description of the scope of claims, and substitutions or modifications equivalent to the scope of the claims are also included in the scope of the present invention. It is within the technical scope.

DESCRIPTION OF SYMBOLS 100 Flying object detection system 101a-c Artificial satellite 102 Base station 103 Antenna 104 Signal line 105 Position analysis apparatus 201 Reception part 1051 Communication interface 1052 CPU
1053 Memory 1054 DSP
1055 Disc 1056 Display unit 1057 Input unit

Claims (8)

Photographing at least two artificial satellites flying on the ground and transmitting them as photographed image data;
Transmitting the position information, velocity information, and attitude information of each satellite by the at least two satellites;
Receiving at least two pieces of image data transmitted from the at least two artificial satellites, and positional information, velocity information, and attitude information of each of the artificial satellites by a base station equipped with a receiving antenna;
Analyzing the at least two pieces of image data received by the antenna by the position analysis device in the base station and calculating the apex angle and azimuth angle of each artificial satellite;
A flying object detection method including the step of analyzing the position of the flying object using the apex angle and the azimuth angle, and the position information, velocity information, and attitude information of each of the artificial satellites,
The apex angle and azimuth angle of each artificial satellite
And generate an observation vector y defined by the following equation:
On the other hand, a coordinate transformation matrix is generated based on the following expression from the position information, velocity information, and attitude information of each artificial satellite,
On the other hand, when the initial position of the flying object is X _ref ,
A coefficient defined by the following equation is generated from the initial position value of the flying object, the position information of the artificial satellite, and the coordinate transformation matrix (Equation 4), (Equation 5), and (Equation 7).
Further, a parameter estimation process defined by the following equation from the coefficient defined by the (expression 11a) to (expression 11c) and the observation vector y:
A method characterized by specifying a flying object position by executing At least two artificial satellites comprising an observation device capable of photographing a flying object flying on the ground, and a transmitter that transmits the photographed image data and position information, velocity information, and attitude information of each of the photographed image data;
An antenna for receiving at least two pieces of image data transmitted from the at least two satellites and position information, velocity information, and attitude information of each of the satellites; and at least the two received by the antennas A flying object detection system comprising a base station including a piece of image data and a position analysis device that analyzes the position of the flying object based on each position information, velocity information, and attitude information. ,
The position analysis device
A calculation unit for analyzing the at least two pieces of image data and calculating an apex angle and an azimuth angle of each artificial satellite;
The apex angle and azimuth angle of each artificial satellite
The observation vector defined by
An observation vector generation unit for generating
A coordinate transformation matrix based on the following equation from the position information, velocity information, and attitude information of each artificial satellite
A coordinate transformation matrix generation unit for generating
When the initial position of the flying object is X _ref
A coefficient defined by the following equation from the initial position value of the flying object, position information of the artificial satellite, and the coordinate transformation matrix (Equation 4), (Equation 5), and (Equation 7)
A coefficient generation unit for generating
A parameter estimation process defined by the following equation from the coefficient defined by the (expression 11a) to (expression 11c) and the observation vector y:
And a parameter estimation processing unit for executing the above. At least two artificial satellites equipped with an observer capable of capturing a flying object flying on the ground, and a transmitter for transmitting the captured image data and each position information, velocity information, and attitude information;
An antenna for receiving at least two pieces of image data transmitted from the at least two artificial satellites and position information, velocity information, and attitude information of each of the artificial satellites, and at least two pieces received by the antennas Position analysis apparatus in a flying object detection system comprising a base station provided with a position analysis apparatus that analyzes the position of the flying object based on the image data of the image and each position information, velocity information, and attitude information Because
A calculation unit for analyzing the at least two pieces of image data and calculating an apex angle and an azimuth angle of each artificial satellite;
The apex angle and azimuth angle of each artificial satellite
The observation vector defined by
An observation vector generation unit for generating
A coordinate transformation matrix based on the following equation from the position information, velocity information, and attitude information of each artificial satellite
A coordinate transformation matrix generation unit for generating
When the initial position of the flying object is X _ref
A coefficient defined by the following equation from the initial position value of the flying object, the position information of the artificial satellite, and the coordinate transformation matrix (Equation 4), (Equation 5), (7)
A coefficient generation unit for generating
A parameter estimation process defined by the following equation from the coefficient defined by the (expression 11a) to (expression 11c) and the observation vector y:
And a parameter estimation processing unit for executing the above. At least two artificial satellites comprising an observation device capable of photographing a flying object flying on the ground, and a transmitter that transmits the photographed image data and position information, velocity information, and attitude information of each of the photographed image data;
An antenna for receiving at least two pieces of image data transmitted from the at least two satellites and position information, velocity information, and attitude information of each of the satellites; and at least the two received by the antennas Position analysis of a flying object detection system comprising a base station including a piece of image data and a position analysis device that analyzes the position of the flying object based on the position information, velocity information, and attitude information of each A computer program executed in an apparatus,
A calculation step of analyzing the at least two pieces of image data and calculating an apex angle and an azimuth angle of each artificial satellite;
The apex angle and azimuth angle of each artificial satellite
The observation vector defined by
An observation vector generation step for generating
A coordinate transformation matrix based on the following equation from the position information, velocity information, and attitude information of each artificial satellite
A coordinate transformation matrix generation step for generating
When the initial position of the flying object is X _ref
A coefficient defined by the following equation from the initial position value of the flying object, position information of the artificial satellite, and the coordinate transformation matrix (Equation 4), (Equation 5), and (Equation 7)
A coefficient generation step for generating
A parameter estimation process defined by the following equation from the coefficient defined by the (expression 11a) to (expression 11c) and the observation vector y:
A computer program comprising: a parameter estimation processing step for executing When Δt is the observation time interval of the observation device mounted on the artificial satellite, the position vector of the flying object at time t obtained from (Equation 12)
And the position vector of the flying object at time t-Δt
And use
The method of claim 1, further comprising the step of calculating a velocity vector of the flying object by. When Δt is the observation time interval of the observation device mounted on the artificial satellite, the position vector of the flying object at time t obtained from (Equation 12)
And the position vector of the flying object at time t-Δt
And use
The system according to claim 2 , further comprising a velocity vector calculation unit that calculates a velocity vector of the flying object according to the above. When Δt is the observation time interval of the observation device mounted on the artificial satellite, the position vector of the flying object at time t obtained from (Equation 12)
And the position vector of the flying object at time t-Δt
And use
The apparatus according to claim 3 , further comprising: a velocity vector calculation unit that calculates a velocity vector of the flying object. When Δt is the observation time interval of the observation device mounted on the artificial satellite, the position vector of the flying object at time t obtained from (Equation 12)
And the position vector of the flying object at time t-Δt
And use
The computer program according to claim 4 , further comprising a velocity vector calculation step of calculating a velocity vector of the flying object by the method.