JP5359738B2 - Telemetry data analysis device and telemetry data analysis program for structure-variable flying object - Google Patents

Description

  The present invention relates to a telemetry data analysis apparatus and a telemetry data analysis program for analyzing telemetry data transmitted from a structure-changeable flying object such as a rocket and detecting a motion of a flying object suspected of being abnormal.

The state of rockets in flight, artificial satellites in orbit, unmanned spacecraft and manned spacecraft is constantly observed at observation stations on the ground through telemetry data transmitted from these flying objects. Telemetry data is time-series data transmitted from various listed devices that indicate the state of the flying object. Mainly, the voltage and current of various mounted devices, the on / off state of operation, and various parts of the flying object Data such as temperature and pressure of the aircraft, transmission / reception power of the transmission / reception device, information on the external environment of the flying object. Each of these data is collected by an encoder of a telemeter transmitter mounted on the flying object and transmitted to a ground observation base in a certain format. As shown in FIG. 10, in the observation base B, after receiving various telemetry data S _{1 to} S _N transmitted from the flying object (rocket R) by the receiving device A (A _{1 to} A _N ) such as an antenna. The signal demodulator / converter D (D _{1 to} D _N ) demodulates and performs A / D conversion and various physical quantity conversions. Then, in the telemetry data analysis device Z, the converted time series data is subjected to statistical analysis and frequency analysis, or the data is visualized and displayed on the display device P, and the result is past data, design values, or Compare with predicted values. Through such time series data analysis, it is confirmed that the flying object was able to fly as planned, abnormal movement of the flying object was detected, and if any accident occurred in the flying object, the cause To improve the installed equipment (Patent Document 1).

  In such conventional telemetry data analysis (time-series data analysis), if the difference is not so large when comparing the telemetry data to be analyzed with past data, design values, or predicted values, the difference is Usually recognized as random noise.

  However, in recent years, it has become clear that in the field of mathematical physics, seemingly random time-series data may follow deterministic laws in chaotic dynamical systems. A characteristic of chaotic dynamical systems is that their dynamics are based on some kind of non-linear causal relationship, so that the time-series data shows complex behavior similar to random noise. In order to determine whether certain time-series data is time-series data based on chaotic dynamical systems, nonlinear statistics such as fractal dimensions and Lyapunov spectra, which are features of chaotic dynamical systems, are calculated, or chaotic dynamics It is necessary to perform data analysis such as visualizing the attractor trajectory (strange attractor) peculiar to the system in the phase space. Such time-series data analysis is particularly called chaos time-series data analysis, and is currently practically used for measurement data of various engineering systems such as various physiological signals in the living body and state signals of plants, etc. (Non-Patent Documents) 1).

US Patent Publication No. 20090105892

Kazuyuki Aihara, Tohru Ikeguchi, Yasushi Yamada, Motomasa Komuro, Basics and Applications of Chaos Time Series Analysis, Sangyo Tosho, 2000.

  Therefore, in the conventional data analysis for telemetry data of flying objects such as rockets, if the telemetry data is seemingly random at first glance, it is a deterministic factor based on the nonlinearity hidden in it. In other words, there was a possibility of overlooking the existence of some chaotic dynamical system in the movement of the flying object. This means that in the conventional telemetry data analysis of the flying object, the suspicious operation of the flying object based on the chaotic dynamical system is overlooked.

  The present invention has been made in view of the above circumstances, and an object thereof is to detect a suspicious operation of a flying object based on a chaotic dynamical system from telemetry data transmitted from a structure-variable flying object such as a rocket. And

  In order to solve the above-mentioned problem, the invention described in claim 1 is a telemetry data analysis device for a variable structure flying object, which is a time series to be analyzed from a telemetry data group transmitted from the variable structure flying object. A data cutout unit that cuts out a data group, a data preprocessing unit that calculates an embedding dimension of the phase space that enables embedding of the cutout time series data group, and a time series embedded in the phase space having the calculated dimension A data visualization processing unit that calculates visualization information from the data group, a first surrogate data group that stores primary statistics from the extracted time-series data group, and a second surrogate data group that stores secondary statistics The surrogate data creation unit that creates the non-linear statistics from the visualized time-series data group during real-time analysis, and is visualized during post-flight analysis A nonlinear statistic estimation unit for estimating a nonlinear statistic from the time series data group and the first and second surrogate data groups, and a nonlinear statistic estimated from the first and second surrogate data groups during post-flight analysis A data test unit that performs the significance test of the nonlinear statistics estimated from the visualized time-series data group and determines whether or not the extracted time-series data group follows a chaotic dynamical system, and chaos Calculate the difference between the nonlinear statistics and visualization information of the time series data group determined to follow the dynamic system and the nonlinear statistics and visualization information of the past time series data group, and the abnormality of the flying object based on the chaotic dynamic system An abnormal suspicion detection unit that detects a suspicious operation, and a control unit that calculates an embedded dimension of the extracted time-series data group each time the structure of the flying object changes with respect to the data preprocessing unit, Characterized by comprising.

  According to the second aspect of the present invention, there is provided a data cut-out means for cutting out a time-series data group to be analyzed from a telemetry data group transmitted from a structure-variable flying object, and a cut-out time-series data group. Data preprocessing means for calculating the embedding dimension of the phase space enabling embedding, data visualization processing means for calculating the visualization information from the time-series data group embedded in the phase space having the calculated dimension, and Surrogate data creation means for creating a first surrogate data group for storing primary statistics from a time series data group and a second surrogate data group for storing secondary statistics, and at the time of visualization during real-time analysis Nonlinear statistics are estimated from the series data group, and at the time of post-flight analysis, the visualized time series data group and the first and second surrogate data A non-linear statistic estimation means for estimating a non-linear statistic from the above, and a non-linear statistic estimated from the first and second surrogate data groups at the time of post-flight analysis as a population, estimated from a visualized time series data group Data test means to determine whether the extracted time-series data group follows a chaotic dynamical system, and nonlinear statistics of the time-series data group determined to follow a chaotic dynamical system Calculating the difference between the quantity and visualization information and the nonlinear statistics and visualization information of the past time series data group, and detecting the suspected abnormality detecting means of the flying object based on the chaotic dynamical system, and A variable structure, characterized in that the data preprocessing means functions as a control means for calculating an embedding dimension of the extracted time-series data group whenever the structure of the flying object changes. Characterized in that it is a telemetry data analysis program projectile.

  According to the present invention, it is possible to systematically find deterministic factors based on a chaotic dynamical system from telemetry data that appears to be random at first glance from a structure-variable flying object. It is possible to detect a motion of a flying object that may have been suspected of being abnormal. Furthermore, when this device is applied to post-flight analysis, it is possible to efficiently identify the cause of a flying object failure that has already occurred and to predict a failure that may occur in the future. Therefore, it is possible to improve the design of the flying object after the next unit.

It is a figure which shows the telemetry data from the multistage rocket used as the analysis object in the telemetry data analysis apparatus which concerns on one embodiment of this invention. It is a figure which shows the characteristic of the telemetry data sent to an observation base from the multistage rocket shown in FIG. 1 is a block diagram of a telemetry data analysis apparatus according to an embodiment of the present invention. It is a flowchart which shows the whole flow of the telemetry data analysis performed in the telemetry data analysis apparatus shown in FIG. 5 is a flowchart showing details of data preprocessing shown in FIG. 4. 5 is a flowchart showing details of data visualization processing shown in FIG. 4. 5 is a flowchart showing details of the nonlinear statistic estimation process shown in FIG. 4. 5 is a flowchart showing details of surrogate data creation processing shown in FIG. 4. 5 is a flowchart showing details of the nonlinear statistic estimation process shown in FIG. 4. It is a block diagram which shows the main components of the observation base which receives the telemetry data from a rocket.

  Hereinafter, a telemetry data analysis apparatus and a telemetry data analysis program for a variable structure flying object according to an embodiment of the present invention will be described with reference to the drawings. FIG. 1 is a diagram showing telemetry data from a multistage rocket R to be analyzed in the telemetry data analysis apparatus 1 (FIG. 3) according to an embodiment of the present invention.

[Telemetry data to be analyzed]
The object of analysis in the present invention is telemetry data transmitted from a structure-variable flying vehicle R such as a multistage rocket as shown in FIG. A two-stage rocket will be described as an example of a multistage rocket. This rocket is a two-stage liquid fuel rocket equipped with a rocket engine with liquid oxygen and liquid hydrogen at each stage R ₁ , R _{2. The} first stage R ₁ has a large solid rocket booster for thrust enhancement. 2 is equipped with, are equipped with satellite R ₃ covered in the payload fairing of the second stage R ₂ tip. In order to further increase the thrust, a plurality of solid auxiliary rockets may be provided in the first stage.

Each of these parts R _{1 to} R ₃ is equipped with a telemeter transmitter to check the state of the rocket R in flight (onboard equipment), and these telemeter transmitters are equipped with parts R _{1 to} R. ₃ equipped with various mounting telemetry data group SG ₁ to SG ₃ representing the state of the device, respectively sent to ground research station B in real time. Each telemetry data groups _SG 1 to SG _3, including _{N 1} pieces of telemetry data _S 11 _~S 1N1, _{N 2} pieces of telemetry data _S 21 _~S 2N2, _{N 3} pieces of telemetry data _S 31 _{to S 3N3.} The number of telemetry data in these telemetry data groups SG _{1 to} SG ₃ is generally different for each part. Although there are a great variety of on-board equipment, the on-board equipment related to the communication system, measurement system, control system, timer / ignition system, power source / control system, wiring system, etc. is essential for the movement of the rocket R.

The parts R ₁ and R ₂ operate in order from the bottom during the flight of the rocket R, and when the rocket R finishes producing thrust to reach the required speed, it is separated from the main body one after another. only R ₃ is turned on orbit the earth. According to such separation of the parts R ₁ and R ₂ from the rocket R, the telemetry data transmitted to the observation base B changes with time as shown in FIG. FIG. 2 is a diagram showing the characteristics of the telemetry data sent from the multistage rocket R shown in FIG. Rocket R was launched at time t _0, when the parts R ₁ at time t ₁ is disconnected from the article body, after time t _1, telemetry data group SG ₁ that has been transmitted from the part R ₁ in the observation base B is It will not be sent. Then, the part R ₂ at time t ₂ is once disconnected from the body, the time t ₂ later, telemetry data group SG _2, which has been transmitted from the part R ₂ observations base B will not be transmitted. Finally, after that (current time t ₃ ), only the telemetry data group SG ₃ from the artificial satellite R ₀ put into the orbit is continuously transmitted to the observation base B. That is, telemetry data transmitted from a structure-variable flying object such as the multistage rocket R is received at the observation base B every time the structure of the flying object changes according to a command transmitted from the observation base B. It is characterized in that the number of telemetry data to be changed changes.

[Configuration of analyzer]
Next, the configuration of the telemetry data analysis apparatus 1 according to one embodiment of the present invention will be described. FIG. 3 is a block diagram of the telemetry data analysis apparatus 1 according to the embodiment of the present invention. The telemetry data analysis apparatus 1 includes a data reception unit 10, a data extraction unit 20, a data preprocessing unit 30, a data visualization processing unit 40, a nonlinear statistic estimation unit 50, a surrogate data creation unit 60, a data The verification unit 70, the suspected abnormality detection unit 80, the data transmission unit 90, the input unit 100, the storage unit 110, and the control unit 120 are provided.

The data receiving unit 10 is a telemetry data group received at the observation base B (antennas A _{1 to} A _N ), demodulated by the signal demodulation conversion device D (D _{1 to} D _N ), and subjected to A / D conversion and the like. Receive.

  The data cutout unit 20 performs a process of cutting out a time-series data group in the time domain for performing telemetry data analysis from the telemetry data group received by the data receiving unit 10.

  When the time-series data group cut out by the data cut-out section 20 is a one-dimensional time-series data group, the data preprocessing section 30 follows the trajectory drawn by the time-series data group by the error neighborhood point method. When the process of calculating the dimension (embedding dimension) m of the phase space H that can be embedded is performed and the extracted time-series data group is multi-dimensional, the dimension n is set as the embedding dimension of the time-series data group. Process to set as.

  The data visualization processing unit 40 binarizes the time series data group based on a predetermined norm that defines the distance between two points in the high-dimensional phase space H having the embedded dimension calculated by the data preprocessing unit 30. A recurrence plot to be visualized, and a unidirectional recurrence plot to visualize the binary data of a time series data based on a norm for evaluating the degree of variation between direction vectors at two points in the high-dimensional phase space H Based on the results of these two plots, a directional neighborhood plot is performed, and based on the results of the directional neighborhood plot, processing is performed to determine whether the time-series data group conforms to a deterministic law.

  The nonlinear statistic estimation unit 50 estimates a correlation dimension Dc, which is one of the fractal dimensions, from the time series data group by a GP (Grassberger-Proccia) method as a nonlinear statistic of the extracted time series data group, Further, processing for estimating the Lyapunov spectrum and the Lyapunov exponent from the time series data group is performed.

  The surrogate data creation unit 60 generates RS surrogate data for storing primary statistics from the extracted time-series data group by the RS (Random Shuffle) method when the data analysis is not real-time analysis but post-flight analysis. At the same time, the FT (Fourier Transform) method is used to create FT surrogate data for storing secondary statistics.

  Note that the nonlinear statistic estimation unit 50 is one of the fractal dimensions by the GP method for the FT surrogate data and the RS surrogate data created by the surrogate data creation unit 60 as the nonlinear statistic of the surrogate data. The correlation dimension Dc is estimated, and the Lyapunov spectrum and Lyapunov exponent are estimated from these surrogate data.

  The data test unit 70 compares the nonlinear statistics estimated from the extracted time series data group and the nonlinear statistics estimated from each surrogate data by the Monte Carlo significance test, and the time series data group is chaotic. To determine whether or not to follow general determinism.

  The suspected abnormality detection unit 80 uses a non-linear statistic of a past time series data group read from the storage unit 110 (to be described later) as a time series data group determined by the data verification unit 70 to follow chaotic determinism ( Compared with fractal dimensions, Lyapunov spectrum, etc.) and visualization data (recurrence plot, attractor trajectory in phase space, etc.), a process for detecting a suspicious operation of the rocket R is performed by obtaining the difference.

  The data transmission unit 90 includes a data extraction unit 20, a data preprocessing unit 30, a data visualization processing unit 40, a nonlinear statistic estimation unit 50, a surrogate data creation unit 60, a data verification unit 70, and a suspected abnormality detection unit 80. The obtained data is transmitted to the storage unit 110 and, if necessary, these data are transmitted to the display device P of the observation base B.

  In addition, the data receiving unit 10, the data extracting unit 20, the data preprocessing unit 30, the data visualization processing unit 40, the nonlinear statistic estimation unit 50, the surrogate data creation unit 60, the data verification unit 70, and the abnormality suspect detection unit 80, The data transmitter 90 constitutes a telemetry data analyzer 1A that is a main part of the telemetry data analyzer 1 according to the present embodiment.

  In addition, the input unit 100 prompts the user to input parameters (details will be described later) necessary for performing the above-described processing in the telemetry data analysis unit 1A.

  The storage unit 110 stores a sub-memory for storing various data obtained by the above processing in the telemetry data analysis unit 1A, and a computer-readable program for executing each step (described later) of the above processing. And main memory for storing. The storage unit 110 includes a RAM (Random Access Memory), a ROM (Read Only Memory), and the like. Furthermore, the sub memory and the main memory of the storage unit 110 may be configured separately.

  The control unit 120 includes a CPU (Central Processing Unit) that controls each operation of the telemetry data analysis unit 1A according to a program read from the storage unit 110. In particular, the control unit 120 is characterized in that the data preprocessing unit 30 is controlled to calculate the embedding dimension of the extracted time-series data group every time the structure of the flying object R changes.

  In the present embodiment, the telemetry data analysis device 1 is configured by integrating the telemetry data analysis unit 1A, the input unit 100, the storage unit 110, and the control unit 120. However, the storage unit 110 is an independent storage device. The telemetry data analyzing unit 1A, the input unit 100, and the control unit 120 may be separated from each other. Alternatively, the telemetry data analysis unit 1B may have a configuration in which the data reception unit 10 and the data transmission unit 90, which are data transmission / reception interfaces, are removed from the telemetry data analysis unit 1A.

  Here, a computer refers to a device that processes structured input in accordance with a predetermined rule and outputs the processed result in a structured manner, for example, a general-purpose computer, a supercomputer, a mainframe, a workstation, a micro computer. Computers, servers, etc. are included. In addition, a configuration (for example, distributed) composed of two or more computers connected via a communication network (for example, a communication network including an intranet, a local area network (LAN), a wide area network (WAN), and a combination thereof). Computer system).

[Data analysis method]
Next, an overall flow of telemetry data analysis performed in the telemetry data analysis apparatus 1 will be described. FIG. 4 is a flowchart showing an overall flow of telemetry data analysis performed in the telemetry data analysis apparatus shown in FIG.

When the data receiving unit 10 receives the telemetry data group from the signal demodulation conversion device D (D _{1 to} D _N ) of the observation base B, the control unit 120 performs the telemetry data analysis unit 1A of the telemetry data analysis unit 1A according to the program stored in the storage unit 110. Each unit is controlled so that each unit performs the following processing.

  In step S10, the data extraction unit 20 reads the telemetry data group received by the data reception unit 10, and in step S20, extracts the time-series time series data group for performing data analysis from the read telemetry data group. I do. At that time, the control unit 120 prompts the user to input this time region from the input unit 100.

  In step S30, the control unit 120 determines whether or not the time series data group extracted by the data extraction unit 20 is a one-dimensional time series data group, and determines that the time series data group is a one-dimensional time series data group. In this case, the data preprocessing unit 30 is caused to perform a process of calculating the dimension (embedding dimension) m of the phase space that can embed the time series data group by the error neighborhood point method, and the extracted time series If the data group is multi-dimensional, a process of setting the dimension n as an embedding dimension of the time series data group is performed. Here, the dimension of the time-series data group of the extracted time-series data group depends on the change in the structure of the flying object R. Therefore, the control unit 120 calculates the embedded dimension of the extracted time-series data group every time the structure of the flying object R changes, such as when each stage is disconnected in the multistage rocket. The processing unit 30 is controlled.

FIG. 5 is a flowchart showing details of the setting process of the embedding dimension. When the time-series data group cut out by the data cut-out unit 20 is the one-dimensional time-series data group x (t), the control unit 120 sends a delay time τ from the input unit 100 to the user in step S41. And the data preprocessing unit 30 is configured to form an m-dimensional phase space H for the input delay time τ. At this time, a point on the trajectory of the extracted time-series data group is the m-dimensional phase space H configured,

It is expressed. Further, the control unit 120 separately displays the trajectory y _m (t) of the time series data group in the m-dimensional phase space H on the display device P of the observation base B.

At this time, unless an appropriate delay time τ is set, it is impossible to determine whether or not the trajectory of the time-series data group in the m-dimensional phase space H is asymptotic to the attractor. Specifically, if the delay time τ is too small, the correlation becomes too large, and the trajectory of the time series data group is collapsed on the diagonal line of each quadrant in the m-dimensional phase space, while the delay time τ is large. If too much, the correlation becomes too small, and the trajectory of the time series data group draws a random trajectory in the m-dimensional phase space. Therefore, the control unit 120 causes the user to check the image of the trajectory y _m (t) of the time-series data group in the m-dimensional phase space displayed on the display device P, and selects an appropriate value from the input unit 100. Prompt to enter the delay time τ.

Next, in step S <b> 42, the data preprocessing unit 30 obtains the lowest embedding dimension by the error neighborhood method and sets it as the dimension m of the phase space H. First, the data preprocessing unit 30 is the nearest point of the point y _m (t) on the locus of the time-series data group cut out in the m-dimensional phase space H represented by the equation (1).

Is calculated by heap sort. Here, n (t, m) is an index of the nearest point. The data preprocessing unit 30 then calculates the ratio between the point-to-point processing in the m-dimensional phase space H and the point-to-point distance in the m + 1-dimensional phase space.

Calculate Here, the distance is defined by the Euclidean norm for both the numerator and denominator. The data preprocessing unit 30 then calculates the average value of a for all data points N (t = 0, 1,..., N).

And the index of the rate of change when the phase space dimension is changed from m dimension to m + 1 dimension

Is calculated.

  Next, in step S43, the data preprocessing unit 30 increases the dimension m to calculate the value of Expression (5), and determines whether or not the value converges. When the extracted time-series data group conforms to some deterministic law, the value of Equation (5) converges, so that the data preprocessing unit 30 changes the value of Equation (5) to within the threshold value. The value of m when it disappears is set as the lowest dimension that can be embedded. This threshold value prompts the control unit 120 to set it from the input unit 100 to the user.

  In step S44, the data preprocessing unit 30 performs a process of embedding the time-series data group cut out in the phase space H of the set embedding dimension m.

  In step S30, if the control unit 120 determines that the time-series data group cut out by the data cut-out unit 20 is a multi-dimensional time-series data group, the control unit 120 performs the data pre-processing unit 30. The dimension n is set as the embedding dimension of the time series data group.

  If no convergence solution is found in step S43, the control unit 120 prompts the user to set an appropriate embedding dimension m from the input unit 100 in step S45, and then performs the process in step S44. Let it be done.

  Next, in step S50, the data visualization processing unit 40 performs a visualization process for enabling the display device P to display the extracted time-series data group.

FIG. 6 is a flowchart showing details of this visualization processing. First, in step S51, when the data visualization processing unit 40 is cut out in the high-dimensional phase space H having the embedding dimension calculated by the data preprocessing unit 30 (m or n: hereinafter, it is assumed to be m). Distance between two points y _m (i) and y _m (j) on the locus of the series data group

Are calculated for all pairs (i, j) of i and j. Here, the distance between the two points is defined by the Euclidean norm. Both i and j are time indexes in the extracted time-series data group.

Here, the control unit 120, the user prompts from the input unit 100, the cut-out time series data group to enter the threshold value Th ₁ for visualizing by binarization to recurrence plots. Next, the data visualization processing unit 40 binarizes D _ij calculated by Expression (6) based on the input threshold Th ₁ . At this time, a point (i, j) where D _ij is equal to or smaller than a threshold value is set to True, and a point (i, j) equal to or higher than the threshold value is set to False. In addition, the control unit 120 displays the time series data group binarized by the data visualization processing unit 40 on the display device P as a recurrence plot on D _ij (i, i, i, j = 0, 1,..., N−1) can be transmitted to the display device P via the data transmission unit 90.

Next, in step S <b> 52, the data visualization processing unit 40 determines the direction vectors at two points y _m (i) and y _m (j) on the trajectory of the extracted time-series data group in the m-dimensional phase space H.

In order to evaluate the variation between

Is calculated.

Here, the control unit 120 provides the user with a threshold Th ₂ for binarizing the time-series data group cut out from the input unit 100 and visualizing it in the same-direction recurrence plot, and Expression (7). , (8) is prompted to input the parameter (minute time difference) Δ. Next, the data visualization processing unit 40 binarizes based on the threshold Th ₂ input the DV _ij calculated by Equation (8). At this time, a point (i, j) where DV _ij is equal to or smaller than a threshold value is set to True, and a point (i, j) equal to or higher than the threshold value is set to False. Further, the control unit 120 displays DV _ij (2) calculated by Expression (8) in order to display the time series data group binarized by the data visualization processing unit 40 on the display device P as the same direction recurrence plot. i, j = 0, 1,..., N−1) can be transmitted to the display device P via the data transmission unit 90.

Next, in step S53, the data visualization processing unit 40 determines that D _ij = (True or False) calculated by Expression (6) and DV _ij = (True or False) calculated by Expression (8). AND of

Calculate This IDNP also takes a binary value of True or False.

Next, in step S60, the data visualization processing unit 40

If the calculated value of R is close to 1, it is determined that the extracted time-series data group follows a deterministic law, and the process proceeds to step S70. Otherwise, The extracted time-series data group is determined to be non-deterministic, and the data analysis ends. Here, the control unit 120 prompts the user to input from the input unit 100 a threshold value that is a criterion for determining whether or not the value of R is close to 1.

  Next, in step S70, the control unit 120 determines whether or not the data analysis is real-time analysis. If the control unit 120 determines that the data analysis is real-time analysis, in step S80, the control unit 120 The estimation unit 50 is caused to estimate nonlinear statistics of the extracted time series data group.

FIG. 7 is a flowchart showing details of the nonlinear statistic calculation process. First, in step S81, the nonlinear statistic estimation unit 50 estimates a correlation dimension Dc, which is one of the fractal dimensions, from the time-series data group cut out by the GP method. Specifically, the nonlinear quantity estimation unit 50

Is calculated. Where I is the snake side function

It is. C _m (r) is a quantity called a correlation integral, and qualitatively, m of radius r centering on a point y _m (i) on the locus of the cut time series data group in the m-dimensional phase space H. Counting the number of points y _m (j) entering the dimensional hypersphere and calculating the ratio to the total number of data points N is performed for all points y _m (i) on the trajectory, and the average is calculated. It was taken.

Next, the nonlinear statistic estimation unit 50 calculates the appropriate range of the radius r.

A dimension d (m) that satisfies the above is calculated for each embedding dimension m. Specifically, when the equation (13) is plotted on a logarithmic graph, the slope of the graph in the logr region where the slope is constant is calculated. Here, the constant inclination means a case where the inclination variation is within a certain threshold over one decade (a section that is 10 times on the logarithmic frequency axis) or more. Therefore, at this time, the control unit 120 prompts the user to input this threshold value from the input unit 100. Then, the nonlinear statistic estimation unit 50 calculates the convergence destination slope accompanying the change in the embedding dimension m as Dc, and this Dc is the condition.

It is determined whether or not the above is satisfied. When this condition is satisfied, the nonlinear statistic estimation unit 50 estimates Dc as the correlation dimension. On the other hand, when this condition is not satisfied, the control unit 120 notifies the user via the display device P that the number of data points is increased and this data analysis is executed again.

  Next, in step S82, the nonlinear statistic estimation unit 50 performs processing for estimating the Lyapunov spectrum and the Lyapunov dimension from the extracted time-series data group.

First, the nonlinear statistic estimation unit 50 obtains each point y _m (t) (t = 0, 1,...) Of N data points on the trajectory of the time-series data group cut out in the m-dimensional phase space H. , N−1), and consider an m-dimensional hypersphere including M neighboring points y _m (k _i ) (i = 1, 2,..., M).

But

Processing for calculating the m × m matrix J _t satisfying. Here, s represents a unit step time. Therefore, first, the control unit 120 prompts the user to input, from the input unit 100, the number M of neighboring points y _m (k _i ) on the trajectory included in the hypersphere. Note that for the neighboring points y _m (k _i ) (i = 1, 2,..., M), the control unit 120 sorts D _ij calculated in step S61 in ascending order using the heap sort algorithm, Performs processing to call the index.

Next, the nonlinear statistic estimation unit 50 performs the operation at time t = 0.

Is calculated. Here, ε _ik and ε _il are the k-th and l-th components of the minute displacement vector ε _i (i = 1, 2,..., M) of Expression (15), respectively. Also, w _kl and c _kl are respectively

The components of the m × m matrix satisfying the above, W and C are the variance matrix and the covariance matrix, respectively. Incidentally, equation (18) is

Minimal condition to minimize

Is a relational expression obtained from

Then, the nonlinear statistic estimation unit 50 solves the right side of the equation (18) by the LU decomposition method using the calculated equation (17), and the J _{t at} the point y _m (t) (t = 0) on the orbit. Is calculated as an approximate Jacobian matrix. At this time, if there is no degeneracy with M ≧ m, the Jacobian matrix J _t can be calculated uniquely. Then, the nonlinear statistic estimation unit 50 performs this process on all data point sequences (that is, repeats N times), so that N Jacobian matrices J _t (t = 0, 1,... N-1) is calculated.

Next, the nonlinear statistic estimation unit 50 generates an m-dimensional orthonormal basis prepared in advance.

The Jacobian matrix at time t = 0 calculated by the above processing

Multiplied by

Map to. here,

It is.

Then, the nonlinear statistic estimation unit 50 orthogonalizes each column vector of z ⁰ _{i in} the equation (23) again by the Gram-Schmidt orthogonalization method,

Is calculated. At this time, the norm ‖v ¹ _iの of the i-th orthogonal vector is equal to the i-th eigenvalue of the Jacobian matrix J, and the exponent part approximates the Lyapunov spectrum.

Next, the nonlinear statistic estimation unit 50 normalizes V ¹ again to calculate W ¹ and repeats the above processing N times to calculate V ¹ , V ² ,..., V ^N. , Lyapunov spectrum by averaging the norms of N column vectors v ^j _i (j = 1, 2,..., N)

Is estimated.

Subsequently, the nonlinear statistic estimation unit 50 calculates the sum of the Lyapunov spectrum λ _i estimated by the equation (26).

For the largest integer P satisfying

According, to calculate the Lyapunov dimension _{D L.}

Next, in step S90, the control unit 120, if at least the maximum Lyapunov exponent λ ₁ in the Lyapunov spectrum λ _i calculated by the nonlinear statistic estimation unit 50 is positive, is a time-series data group that has been cut out. In step S100, the data obtained in the above process is stored in the storage unit 110. On the other hand, if the maximum Lyapunov exponent λ1 is negative, this data analysis is terminated assuming that the extracted time-series data group follows a deterministic law that is not chaotic.

  Here, it returns to the process of step S70 and continues description. In step S70, when the control unit 120 determines that the data analysis is not a real-time analysis but a post-flight analysis, in step S110, the control unit 120 extracts the time-series data group extracted from the surrogate data creation unit 60. RS surrogate data is created by RS (Random Shuffle) method and FT surrogate data is created by FT (Fourier Transform) method.

  FIG. 8 is a flowchart showing details of the surrogate data creation process. First, in step S111, the control unit 120 determines whether or not the time series data group embedded in the data preprocessing unit 30 was originally a one-dimensional time series data group, and the one-dimensional time series data group. In step S112, the control unit 120 causes the surrogate data creation unit 60 to change the arrangement of the data point sequences in the original one-dimensional time-series data group at random. A predetermined number of sets of data are created. In the RS surrogate data created at this time, primary statistics such as average value, variance, and median are stored.

  Next, in step S113, the surrogate data creation unit 60 performs fast Fourier transform on the original one-dimensional time-series data group to convert the one-dimensional time-series data group from real space to frequency space. . Then, the surrogate data creation unit 60 randomizes only the phase while preserving the amplitude of each wave in the fast Fourier transformed frequency space, and then performs inverse fast Fourier transform to inverse transform from the frequency space to the real space. Create a predetermined number of surrogate data. In the FT surrogate data created at this time, secondary statistics such as power spectral density are stored.

  On the other hand, when the control unit 120 determines in step S111 that the time series data group embedded in the data preprocessing unit 30 is originally a multidimensional time series data group, the control unit 120 In S114, the surrogate data creation unit 60 creates a predetermined number of sets of RS and TF surrogate data for each (one-dimensional) time-series data group in steps S112 and S113.

Note that, in step S112, the control unit 120, an input to the user from the input unit 100, the number of sets _{S 2} of the FT surrogate data created in the number of sets _{S 1} and step S113 in RS surrogate data generated in step S112 Encourage you to. Hereinafter, it is assumed that the total number of surrogate data sets is S (= S ₁ + S ₂ ).

Next, in step S120, the nonlinear statistic estimation unit 50 extracts the extracted time-series data group and S sets of surrogate data (S ₁ set of RS surrogate data and S ₂₎ created by the surrogate data creation unit 60. for each set of FT surrogate data), as a nonlinear statistics, correlation dimension Dc, Rya-flop Nof spectrum _{λ i (i = 1,2, ···} , N), and estimates the Lyapunov dimension _{D L.}

  FIG. 9 is a flowchart showing details of this nonlinear statistic estimation process. First, in step S121, the control unit 120 determines that each of the extracted original time series data group and the S sets of surrogate data created by the surrogate data creation unit 60 was originally one-dimensional time series data. If it is determined, in step S122, data is embedded in the phase space H of the embedding dimension m calculated by the data preprocessing unit 30.

  Next, in step S123, the nonlinear statistic estimation unit 50 performs the same process as the process performed in step S81 on the original time-series data group and each surrogate data embedded in the m-dimensional phase space H. By applying, the correlation dimension Dc which is one of the fractal dimensions is estimated.

Next, in step S124, the nonlinear statistic estimation unit 50 performs the same process as the process performed in step S82 on the original time-series data group and each surrogate data embedded in the m-dimensional phase space H. by subjecting, Rya-flop Nof spectrum _{λ i (i = 1,2, ···} , N) ( norm and the maximum Lyapunov exponent), and estimates the Lyapunov dimension _{D L.}

  On the other hand, when the control unit 120 determines in step S121 that each of the original time-series data group and the S sets of surrogate data created by the surrogate data creation unit 60 was originally multidimensional time-series data. In step S125, data is embedded in the phase space H with the dimension itself as the embedding dimension, and processing is recommended to step S123.

Next, in step S130, the data verification unit 70 uses the non-linear statistics Q ₀ estimated from the original time-series data group by the non-linear statistics estimation unit 50 and the non-linear statistics Q _l (l = 1, 2, ..., S) with a Monte Carlo significance test. At that time, the control unit 120 prompts the user to input the significance level α (for example, α = 0.05) in the Monte Carlo significance test from the input unit 100. Here, _{Q 0} and _{Q l} is correlation dimension Dc, Rya-flop Nof spectrum _{λ i (i = 1,2, ···} , N) ( norm and the maximum Lyapunov exponent), and a Lyapunov dimension _{D L.}

Next, in step S140, the control unit 120, _{Q 0} is largest from the direction (Bs + 1) · α / 2 -th _{Q l,} or sandwiched between towards smallest (Bs + 1) · α / 2 th and it is judged whether or not the that, in the case where Q ₀ is not pinched during this time, the original time series data group is determined that follow the chaotic rule, the process proceeds to step S150.

On the other hand, in step S140, if the Q ₀ is determined to be sandwiched during this time, the original time series data group as not in accordance with chaotic law, and terminates the data processing.

Next, in step S150, the suspected abnormality detection unit 80 determines the nonlinear statistic Q ₀ (correlation dimension Dc, Lyapunov spectrum λ) of the original time-series data group determined to follow the chaotic determinism by the data verification unit 70. _i (i = 1, 2,..., N) and Lyapunov dimension D _L ) and visualization data (recurrence plot, directional recurrence plot, directional neighborhood plot, attractor trajectory in phase space), The non-linear statistic of the past time-series data group read from the storage unit 110 and the Lyapunov dimension D _L ) and visualization data are compared with each other to detect the difference. Then, the control unit 120 presents the difference detected by the suspected abnormality detection unit 80 to the user using the display device P or the like. Accordingly, the user can detect an operation with a suspected abnormality based on the chaotic dynamical system of the structure-variable flying object R that has transmitted the time-series data group to be analyzed.

  Note that the present invention is not limited to the above-described embodiment, and can be embodied in other configurations without departing from the gist thereof.

  For example, as another embodiment, a telemetry data analysis unit 1C having a learning mechanism such as a three-layer neural network can be applied. In the telemetry data analysis unit 1C, by using this learning mechanism, for example, based on the back-propagation learning rule, a model of the chaotic dynamical system is directly constructed from the extracted time-series data group, so that the chaos dynamics is derived from the time-series data group. The system is identified. Then, using the constructed model of the chaotic dynamical system, non-linear statistics such as the Lyapunov spectrum are estimated from the time series data group to be analyzed, and compared with the past data stored in the storage unit 110. Thus, the suspected abnormality detection unit 80 can detect the motion of the flying object suspected of being abnormal based on the model of the chaotic dynamical system.

DESCRIPTION OF SYMBOLS 1 Telemetry data analysis apparatus 1A Telemetry data analysis part 10 Data receiving part 20 Data extraction part 30 Data pre-processing part 40 Data visualization processing part 50 Nonlinear statistic estimation part 60 Surrogate data creation part 70 Data verification part 80 Abnormal doubt detection part 90 Data transmission unit 100 Input unit 110 Storage unit 120 Control unit

Claims (2)

A data cutout unit that cuts out a time-series data group to be analyzed from a telemetry data group transmitted from a variable-structure flying object;
A data pre-processing unit that calculates the embedding dimension of the phase space that enables embedding the extracted time-series data group;
A data visualization processing unit for calculating visualization information from a time-series data group embedded in a phase space having a calculated dimension;
A surrogate data creation unit for creating a first surrogate data group for storing primary statistics and a second surrogate data group for storing secondary statistics from the extracted time-series data group;
Non-linear statistics are estimated from the visualized time-series data group during real-time analysis, and non-linear statistics are estimated from the visualized time-series data group and the first and second surrogate data groups during post-flight analysis A statistic estimator;
At the time of post-flight analysis, the nonlinear statistics estimated from the first and second surrogate data groups were used as a population, and the significance test of the nonlinear statistics estimated from the visualized time series data groups was performed and extracted. A data test unit that determines whether or not the time series data group follows a chaotic dynamical system;
Calculate the difference between the nonlinear statistics and visualization information of the time series data group determined to follow the chaotic dynamic system and the nonlinear statistics and visualization information of the past time series data group. A suspected abnormality detection unit for detecting suspected abnormal behavior;
A control unit that calculates the embedding dimension of the cut-out time-series data group each time the structure of the flying object changes with respect to the data preprocessing unit,
A telemetry data analysis device for a variable structure flying object characterized by comprising: Computer
Data extraction means for extracting a time-series data group to be analyzed from a telemetry data group transmitted from a structure-variable flying object;
Data pre-processing means for calculating the embedding dimension of the phase space that enables embedding the extracted time-series data group;
Data visualization processing means for calculating visualization information from a time-series data group embedded in a phase space having a calculated dimension;
Surrogate data creation means for creating a first surrogate data group for storing primary statistics and a second surrogate data group for storing secondary statistics from the extracted time-series data group;
Non-linear statistics are estimated from the visualized time-series data group during real-time analysis, and non-linear statistics are estimated from the visualized time-series data group and the first and second surrogate data groups during post-flight analysis Statistic estimation means;
At the time of post-flight analysis, the nonlinear statistics estimated from the first and second surrogate data groups were used as a population, and the significance test of the nonlinear statistics estimated from the visualized time series data groups was performed and extracted. Data verification means for determining whether or not the time series data group follows a chaotic dynamical system;
Calculate the difference between the nonlinear statistics and visualization information of the time series data group determined to follow the chaotic dynamic system and the nonlinear statistics and visualization information of the past time series data group. A suspected abnormality detecting means for detecting a suspected anomaly;
Control means for calculating the embedding dimension of the cut out time-series data group each time the structure of the flying object is changed with respect to the data preprocessing means,
A telemetry data analysis program for a structure-variable flying object, which is characterized by being made to function.

Images