CN104142624A - Time synchronization method and system based on waveform matching - Google Patents

Abstract

The invention provides a time synchronization method and system based on waveform matching. According to the method and system, based on the correlation principle, accurate time information is synchronized for a set of data not containing time scales or containing wrong time scales by utilizing a set of data containing reference time. The method comprises the steps that firstly, the data to be matched are segmented, correlation solution is performed on all segments and different sub-segments of reference data, when the correlation coefficient is maximum, the time difference between the segment of data and the reference data is calculated, and therefore time difference sequences of all the segments are obtained; secondly, the time difference sequences serve as observed values, time differences and time-scale drift factors serve as parameters to be estimated, and least square adjustment is performed to figure out the parameters to be estimated; finally, accurate time is synchronized for the data to be matched through the time differences and the time-scale drift factors. The time synchronization method and system can be widely applied to time synchronization during multi-sensor data fusion and can be further applied to correctness detection for time synchronization based on other methods.

Description

A kind of method for synchronizing time and system based on Waveform Matching

Technical field

The present invention relates to sensor application field, concrete is a kind of method for synchronizing time and system based on Waveform Matching.

Background technology

At present, sensor Data Fusion research is in the ascendant, and some proven technique has been applied in Practical Project, and has obtained good effect.The application of sensor Data Fusion not only can improve precision and the reliability of system, can also improve the measurement range of system, increases the confidence level of system, shortens the reaction time of system.But Fusion is a complicated uncertain information processing procedure, the precondition that can merge is that the information obtaining from each sensor is must be the description to the synchronization of same target.This comprises two aspects, first will guarantee that the information that each sensor obtains is the description that the same parameter of same target is carried out.The temporal information that secondly, guarantee the data that merge is synchronous.In dynamic duty environment, it is particularly outstanding that time synchronization problem shows.

In different engineering practices, there is the method for synchronizing time for particular problem.For example, utilize the method for curve to carry out the algorithm of time synchronized, while utilizing the method for sequence to solve, become the Time synchronization algorithm of observation, and adopt smothing filtering algorithm, by the methods such as time synchronized corresponding to the measurement data between each sensor.Before these algorithm application, each sensor is all stamped identical markers when starting to measure, and carries out timing afterwards according to sampling rate separately.Yet, for whole measuring process, all there is not identical markers in each sensor, or the processor that itself carries out timing the phenomenon of frequency marking drift occurs because of the impact of the factors such as temperature characterisitic, and above algorithm all can not correctly be realized the time synchronized between multisensor.

The more typical example of Fusion is exactly the integrated navigation system that GPS (Global Positioning System, GPS) and INS (Inertial Navigation System, inertial navigation system) form.The atomic clock that wherein GPS is equipped with degree of precision provides precise time for it, and INS can only carry out counter timing by the frequency of setting, and drift phenomenon often occurs INS frequency marking.GPS/INS integrated navigation system adopts Kalman filtering to carry out data fusion conventionally, and while only having GPS and INS subsystem data to put at one time, integrated navigation just has actual meaning.The most general, the effective practice is to utilize GPS pps pulse per second signal as the time synchronized benchmark of integrated navigation system at present, and integrated structure design realizes time synchronized.Yet when a IMU newly going out (Inertial Measurement Unit, Inertial Measurement Unit) need carry out integrated navigation test, though said method precision is high, complexity is high, cost is large, needs software and hardware support, inapplicable.In actual test process, the integrated navigation system that often carries one group of degree of precision is as with reference to system, and INS has wherein stamped high-precision GPS markers.

INS markers in integrated navigation system can be used as reference time system, only need the suitable algorithm of design, can realize the time synchronized between IMU to be tested and frame of reference, this algorithm is not only applicable to the time synchronized between inertial sensor, also can be applicable to the time synchronized between other sensors with reference time scale.

Summary of the invention

The present invention is directed to the problems referred to above, designed a kind of method for synchronizing time and system based on Waveform Matching; The method and system be to the whether identical not restriction of the original time system of different sensors data, and simple to operate, computing is quick, cost is low;

Technical scheme of the present invention is: a kind of method for synchronizing time based on Waveform Matching, comprises the steps:

Step 1: read reference data R ₀and data T to be matched ₀;

Step 2: unified reference data R ₀and data T to be matched ₀data processing form, extract the reference data R after conversion ₁and data T to be matched ₁;

Step 3: respectively to R ₁and T ₁down-sampled, to the reference data R after down-sampled ₂and data T to be matched ₂carry out simple crosscorrelation sequence and solve, corresponding interval while getting simple crosscorrelation sequence maximal value, obtains thick mistiming int_bias now;

Step 4: the thick mistiming int_bias according to trying to achieve, makes reference data R ₂time period cover data T to be matched ₂time period, and front and back respectively exceed data T to be matched ₂certain length, intercepting reference data R ₂and data T to be matched ₂; Enter step 5;

Step 5: unified R ₂and T ₂sampling rate, obtains the reference data R of identical sampling rate ₃with data T to be matched ₃;

Step 6: try to achieve R ₃and T ₃the poor bias of optimal time and markers drift parameter k.

Step 7: according to best bias and k combination, upgrade data T to be matched ₀time, and storage;

Step 8: finish.

Unified reference data R in described step 2 ₀and data T to be matched ₀the formula of data processing form as follows:

v＝Δs/Δt _s

Wherein, Δ s is gyro or the measured increment size of accelerometer between former and later two samplings, Δ t _sfor the time interval between former and later two samplings, v is the rate value that Δ s is corresponding;

Reference data R after described conversion ₁and data T to be matched ₁that the time row after changing and an axis data of carrying out Waveform Matching form;

In described step 3, comprise the steps:

Step 3.1: respectively to R ₁and T ₁down-sampled, obtain the reference data R after down-sampled ₂and data T to be matched ₂;

Step 3.2: to down-sampled rear R ₂and T ₂among shorter sequence trailing zero, until both are equal in length;

Step 3.3: calculate R ₂with T ₂simple crosscorrelation sequence; Simple crosscorrelation sequence computing formula is as follows:

${\hat{R}}_{\mathrm{xy}}\left(m\right)=\left\{\begin{array}{cc}\underset{K=0}{\overset{N-m-1}{\mathrm{\&Sigma;}}}{x}_{K+m}{y}_m& m\&GreaterEqual;0\\ {\hat{R}}_{\mathrm{xy}}(-m)& m<0\end{array}\right.$

Wherein, x and y are the sequence of carrying out correlativity processing, and x and y length are N, and K and m are respectively used to represent the ordinal number of the value in sequence x and y, m value is-(N-1) to+(N-1);

Step 3.4: get simple crosscorrelation sequence corresponding interval during maximal value m, obtains thick mistiming int_bias now, and the thick mistiming can utilize following formula to try to achieve:

int_bias＝R ₂(t ₀)-T ₂(t ₀)-T ₂_interval×m

Wherein, R ₂(t ₀) be R ₂start time, T ₂(t ₀) be T ₂start time, T ₂_ interval is the time span of the interval representative that staggers, is T ₂sampling interval, T ₂_ interval * m is the time that two column datas stagger;

In described step 6.15, comprise following sub-step:

Step 6.1: calculate sub-sampling T_inter interval time;

Step 6.2: take T_inter as factor, calculate the data to be matched δ s of translation distance to the right;

Step 6.3: judge whether δ s is the positive integer times of sampling interval; If so, forward step S6.15 to, if not, forward step S6.4 to;

Step 6.4: by data data T to be matched ₃translation δ s obtains T ₄;

Step 6.5: take Δ t as length of window, by T ₄be divided into limited n section, get present segment data T ₄(i), 1≤i≤n;

Step 6.6: judge whether data to be matched extract n section; If so, forward step S6.13 to; If not, forward step S6.7 to;

Step 6.7: set moving window W ₁, window is long is Δ t, starting point is the present segment zero hour;

Step 6.8: judge whether window sliding distance is greater than preseting length σ _s; As if so, forwarded step S6.11 to; Otherwise, with W ₁for window, get reference data R ₃(W ₁);

Step 6.9: calculate T ₄and R (i) ₃(W ₁) between Pearson's related coefficient; Pearson's Calculation of correlation factor as shown in the formula:

$r=\frac{\underset{i=1}{\overset{N}{\mathrm{\&Sigma;}}}(x_i-\stackrel{\&OverBar;}{x})(y_i-\stackrel{\&OverBar;}{y})}{\sqrt{\underset{i=1}{\overset{N}{\mathrm{\&Sigma;}}}{(x_i-\stackrel{\&OverBar;}{x})}^2\&CenterDot;\underset{i=1}{\overset{N}{\mathrm{\&Sigma;}}}{(y_i-\stackrel{\&OverBar;}{y})}^2}}$

Wherein for x serial mean, for y serial mean, N is the length of sequence x and y, and r is the Pearson's related coefficient between sequence x and y;

Step 6.10: judge that whether reference data reads last position, if so, forwards step S6.11 to; Otherwise, by window W ₁slide backward after the sampling interval of a reference data, forward step S6.8 to;

Step 6.11: calculate this section of mistiming Bias (i) and corresponding time time (i);

Step 6.12: extract next section of data to be matched.Take off one section of T ₄(i), forward step 6.6 to;

Step 6.13: with the poor bias of least square adjustment seeking time and markers drift parameter k; Using bias and k as solve for parameter, and Bias (i) is observed reading, and observation equation is as follows:

Bias(i)＝k×[time(i)-t ₀]+bias

Wherein the initial value Bias0 of bias is Bias (i) sequence intermediate value, and the corresponding time is t ₀, the initial value of k is 0;

Step 6.14: judge that whether δ s is greater than a sampling interval, if so, forwards step S6.15 to; Otherwise, by the δ s T_inter that doubles, forward step S6.3 to;

Step 6.15: select best bias and k combination; Utilize not bias and k on the same group, obtain the to be matched data of different time after synchronous;

In described step 6.15, if data are single shaft, that group bias while selecting covariance maximum and k are as best bias and k combination; If data are multiaxis, select the most believable axis data to carry out covariance and solve, select best bias and k combination.

In described step 7, utilize following formula to calculate and upgrade T ₀time:

t′＝(t-t ₀)×k+bias+t

Wherein, t is T ₀time row, t ₀for time corresponding to mistiming initial value of selecting in step 6.12, k is markers drift parameter, and bias is the mistiming, and t ' is the rear time for t upgrades.

A clock synchronization system based on Waveform Matching, is characterized in that, comprises as lower module:

Read module: for reading reference data R ₀and data T to be matched ₀;

Data processing format module: for unified reference data R ₀and data T to be matched ₀data processing form, extract the reference data R after conversion ₁and data T to be matched ₁;

Down-sampled module: for respectively to R ₁and T ₁down-sampled, to the reference data R after down-sampled ₂and data T to be matched ₂carry out simple crosscorrelation sequence and solve, corresponding interval while getting simple crosscorrelation sequence maximal value, obtains thick mistiming int_bias now;

Thick mistiming module: for according to the thick mistiming int_bias trying to achieve, make reference data R ₂time period cover data T to be matched ₂time period, and front and back respectively exceed data T to be matched ₂certain length, intercepting reference data R ₂and data T to be matched ₂;

Sampling rate computing module: for unified R ₂and T ₂sampling rate, obtains the reference data R of identical sampling rate ₃with data T to be matched ₃;

Computing module: for trying to achieve R ₃and T ₃the poor bias of optimal time and markers drift parameter k;

Update module: for according to best bias and k combination, upgrade data T to be matched ₀time, and storage.

Accompanying drawing explanation

Fig. 1 is schematic flow sheet of the present invention;

Fig. 2 is the schematic flow sheet of step 6 of the present invention;

Fig. 3-1st, the vertical deflection speed of EPSON gyro waveform schematic diagram;

Fig. 3-2nd, the vertical deflection speed of FSAS gyro waveform schematic diagram;

Fig. 4-1st, the maximum correlation coefficient schematic diagram that three axle segmentations of EPSON gyro are tried to achieve;

Fig. 4-2nd, the maximum correlation coefficient schematic diagram that three axle segmentations of EPSON accelerometer are tried to achieve;

Fig. 5-1st, the time scale difference sequence schematic diagram that three axle segmentations of EPSON gyro solve;

Fig. 5-2nd, the time scale difference sequence schematic diagram that three axle segmentations of EPSON accelerometer solve;

Fig. 6 is the vertical direction time match of gyro result schematic diagram;

Fig. 7 is the vertical direction time match of gyro result partial schematic diagram;

Fig. 8 is system architecture schematic diagram of the present invention.

Embodiment

Below in conjunction with drawings and Examples, the present invention is described in further detail.

Embodiment describes with the data of inertial sensor collection, and inertial sensor can gather multiaxis data, and the present embodiment take multiaxis data processing as example.

Data Source in embodiment: when vehicle-mounted test, carry low precision IMU EPSON and the GPS/INS integrated navigation system that contains high precision IMU FSAS simultaneously, in the data that wherein EPSON collects with counter timing, in integrated navigation system, FSAS process and GPS integrated design, stamped accurate gps time.

Data layout in embodiment: EPSON sampling rate is 125Hz, data are containing 7 row, wherein first classify time row as, 2nd～4 classify the three row angular speeds that gyroscope gathers as, be respectively forward direction, dextrad and vertical, 5th～7 classify the three row acceleration that accelerometer gathers as, are respectively forward direction, dextrad and vertical; FSAS sampling rate is 200Hz, and first classifies time row as, and 2nd～4 classify the three row angular velocity increments that gyroscope gathers as, are respectively forward direction, dextrad and vertical, and 5th～7 classify the three row acceleration increments that accelerometer gathers as, are respectively forward direction, dextrad and vertical.

The present invention can be without time mark or time and marks the synchronous precise time information of inaccurate data.As shown in Figure 1, Fig. 2 is the poor and markers drift parameter solution procedure of optimal time in Fig. 1 to implementation step.

S1: read reference data R ₀and data T to be matched ₀.Wherein FSAS data are as with reference to data R ₀, EPSON data are as data T to be matched ₀.

S2: unified reference data R ₀and data T to be matched ₀data processing form, extract the reference data R after conversion ₁and data T to be matched ₁(useful data row); According to existing prior imformation, judgement R ₀with T ₀whether represented parameter is identical.If different, according to the mutual relationship between parameter, change, finally with frequently-used data form, be as the criterion.In this example, EPSON is rate form, and FSAS is increment form, and conventional data processing form is rate format.It is as follows that increment is converted to the formula of speed:

v＝Δs/Δt _s

Wherein, Δ s is gyro or the measured increment size of accelerometer between former and later two samplings, Δ t _sfor the time interval between former and later two samplings, v is the rate value that Δ s is corresponding.

Utilize above formula to R ₀change the average that the corresponding time is former and later two sampling times.If carry out Waveform Matching with an axis data wherein, extract time row after conversion and form new reference data R when axis data ₁.In this example, the data of 6 axles are all carried out to Waveform Matching, final FSAS data R ₁for R ₀data after transforming, EPSON data T ₁for T ₀.

As Fig. 3-1 and Fig. 3-2 are depicted as R ₁with T ₁the vertical gyro data waveform of the 4th row, as can be seen from Figure, two waveforms have high consistency, and other axles also have same characteristic features, can carry out Waveform Matching according to the related coefficient of two data like this.

S3: to R ₁and T ₁down-sampled, to the reference data R after down-sampled ₂and data T to be matched ₂carry out simple crosscorrelation sequence and solve, corresponding interval while getting simple crosscorrelation sequence maximal value, obtains thick mistiming int_bias now;

First, treat matched data and reference data and carry out down-sampled processing, down-sampled rear frequency is the common factor of two groups of data original frequencies, meeting on the basis of thick mistiming precision, the smaller the better.In embodiment by R ₁and T ₁down-sampled is 5Hz, to down-sampled rear FSAS data R ₂and EPSON data T ₂carry out simple crosscorrelation sequence and solve, corresponding interval while getting simple crosscorrelation sequence maximal value, obtains thick mistiming int_bias now.Concrete steps are as follows:

S3.1: to R ₁and T ₁down-sampled, obtain the reference data R after down-sampled ₂and data T to be matched ₂;

S3.2: at R ₂and T ₂shorter sequence trailing zero is until both are equal in length;

S3.3: calculate R ₂with T ₂simple crosscorrelation sequence.It is as follows that simple crosscorrelation sequence is calculated General Principle formula:

${\hat{R}}_{\mathrm{xy}}\left(m\right)=\left\{\begin{array}{cc}\underset{K=0}{\overset{N-m-1}{\mathrm{\&Sigma;}}}{x}_{K+m}{y}_m& m\&GreaterEqual;0\\ {\hat{R}}_{\mathrm{xy}}(-m)& m<0\end{array}\right.$

Wherein x and y are the sequence of carrying out correlativity processing, and x and y length are N, and subscript K and m be for representing the ordinal number of value of sequence x and y, m value be-(N-1) arrive+(N-1). for simple crosscorrelation sequence, while being less than zero for m, utilize m to be greater than the simple crosscorrelation sequence that zero formula calculates.When m get-(N-1) to+in (N-1) during different value, all can utilize above formula to calculate and try to achieve a cross correlation value finally, obtain the cross correlation function sequential value that vector length is 2 * N-1.

Step 3.4: get simple crosscorrelation sequence corresponding interval during maximal value m, obtains thick mistiming int_bias now, in this example, and the R after end zero-padded length equates ₂with T ₂the 4th row be sequence x and the y that need carry out correlativity processing, utilize above formula to calculate and try to achieve simple crosscorrelation sequential value corresponding m when value is maximum, is two column data correlativitys when maximum, the interval of staggering.The thick mistiming can utilize following formula to try to achieve:

int_bias＝R ₂(t ₀)-T ₂(t ₀)-T ₂_interval×m

R ₂(t ₀) be R ₂start time, T ₂(t ₀) be T ₂start time, T ₂_ interval is the time span of the interval representative that staggers, is T in this example ₂sampling interval, T ₂_ interval * m is the time that two column datas stagger.

S4: the thick mistiming int_bias according to trying to achieve, makes reference data R ₂time period cover data T to be matched ₂time period, and front and back respectively exceed data T to be matched ₂certain length, intercepting reference data R ₂and data T to be matched ₂.The object of this step is according to the thick mistiming, finds data T to be matched ₂with reference data R ₂between the corresponding time period, make reference data R ₂comprise data T to be matched ₂, and front and back respectively exceed data T to be matched ₂the principle of certain length, intercepting reference data R ₂and data T to be matched ₂.In the present embodiment, thick match time of the precision of trying to achieve is 1 second level, therefore when getting FSAS data, make respectively to exceed FSAS data 0.5 second before and after it, enters S5.

S5: unified R ₂and T ₂sampling rate, obtains the reference data R of identical sampling rate ₃with data T to be matched ₃; FSAS sampling rate is 200Hz, and EPSON sampling rate is 125Hz, finally to T ₀carry out linear interpolation, making its sampling rate is 200Hz.

Through above-mentioned in steps after, FSAS data are R ₃, EPSON data are T ₃.

S6: try to achieve R ₃and T ₃the poor bias of optimal time and markers drift parameter k; Solving of mistiming and markers drift parameter is the core of this algorithm, and calculated amount is maximum, is also a most complicated step simultaneously.

S6.1: calculate sub-sampling T_inter interval time, with sampling interval 1/10～1/2 between be advisable.Sub-sampling interval is defined as the arbitrary value lower than sampling interval, in concrete enforcement, take user's accuracy requirement as criterion, suitably value.T_inter is too small, and precision is improved, but not raising of accuracy, and calculated amount increases; T_inter is excessive, does not reach and proposes high-precision requirement.For example, R in this example ₃with T ₃the sampling interval time be 0.005 second, sub-sampling interval T _ inter is defined as 1/5 second of sampling interval.

S6.2: take T_inter as factor, calculate the data to be matched δ s of translation distance to the right.At circulation time for the first time, δ s is 0 times of T_inter, i.e. δ s=0 * T_inter=0.

S6.3: judge whether δ s is the positive integer times of sampling interval.If so, forward step S6.15 to if not, forward step S6.4 to.

S6.4: by EPSON data translation δ s.First by T ₃time row carry out to right translation, translation distance is δ s, after translation, the time be listed as interpolated point, the numerical value T that other column datas linearity interpolations are made new advances ₄.For example there is original time t ₁, t ₂, corresponding numerical value is val ₁, val ₂, for new time t ₃(t ₁< t ₃< t ₂), corresponding numerical value val ₃by following formula linear interpolation out:

${\mathrm{val}}_3={\mathrm{val}}_1+(t_3-t_1)\&times;\frac{{\mathrm{val}}_2-{\mathrm{val}}_1}{t_2-t_1}$

S6.5: take Δ t as length of window, by T ₄be divided into limited n section, get present segment data T ₄(i), 1≤i≤n.In Δ t length of window, should be able to comprise abundant multidate information, to restrain the impact of random noise.In this example, Δ t is taken as 60 seconds, T ₄be divided into 61 sections, get first paragraph T ₄(1).

S6.6: judge whether data to be matched extract n section.If so, forward step S6.13 to.If not, forward step S6.7 to.

S6.7: set moving window W ₁, window is long is Δ t, starting point is the present segment zero hour.It is 60 seconds that window length is grown the same with the window in S6.5, W ₁at every turn along time shaft backward sampling interval of translation get reference data, window sliding distance is W ₁backward translation T.T., account form is: W ₁the number of translation sampling interval is backward multiplied by the sampling interval time.

S6.8: judge whether window sliding distance is greater than preseting length σ _s(σ _saccording to thick match time of precision and determining, be about twice thick match time of precision), as if so, forwarded step S6.11 to.Otherwise, with W ₁for window, get reference data R ₃(W ₁).In this example, σ _sit is 1 second.

S6.9: calculate T ₄and R (i) ₃(W ₁) between Pearson's related coefficient.Pearson's Calculation of correlation factor principle as shown in the formula:

$r=\frac{\underset{i=1}{\overset{N}{\mathrm{\&Sigma;}}}(x_i-\stackrel{\&OverBar;}{x})(y_i-\stackrel{\&OverBar;}{y})}{\sqrt{\underset{i=1}{\overset{N}{\mathrm{\&Sigma;}}}{(x_i-\stackrel{\&OverBar;}{x})}^2\&CenterDot;\underset{i=1}{\overset{N}{\mathrm{\&Sigma;}}}{(y_i-\stackrel{\&OverBar;}{y})}^2}}$

Wherein for x serial mean, for y serial mean, N is the length of sequence x and y, and r is the Pearson's related coefficient between sequence x and y.In this example, T ₄and R (i) ₃(W ₁) respectively containing six column datas, each respective column data are carried out separately related coefficient according to above formula and solved, each respective shaft all obtains a related coefficient, totally six.

S6.10: judge that whether reference data reads last position, if so, forwards step S6.11 to.Otherwise, by window W ₁slide backward after the sampling interval of a reference data, forward step S6.8 to.

S6.11: calculate this section of mistiming Bias (i) and corresponding time time (i).For this section of EPSON data, W ₁sliding process in, T ₄and R (i) ₃(W ₁) all can obtain corresponding facies relationship ordered series of numbers between each respective column data.For one group of respective column data wherein, when related coefficient is maximum, two column data waveform optimum matching, maximum correlation coefficient max_corr now, W ₁the sampling interval number sliding is to the right num.This axle can calculate a Bias (i) with following formula, and other axles in like manner.

Bias(i)＝R ₃(t ₀)-T ₄(t ₀)+R ₃_interval×num

R ₃(t ₀) be R ₃start time, T ₄(t ₀) be T ₄start time, R ₃_ interval is W ₁dextrad moves the time of a sampling interval representative, and this example is R ₃sampling interval, R ₃_ interval * num is W ₁total duration slides to the right.

Calculate after the Bias (i) of six axles, each axle max_corr of take is weight, and weighted mean obtains T ₄(i) unique Bias (i), the max_corr of this section (i) is the arithmetic mean of each axle max_corr, time (i) is T ₄(i) time in the moment in the middle of.

Mistiming between this section of data to be matched and reference data has calculated complete.

S6.12: extract next section of data to be matched.Take off one section of T ₄(i), forward step 6.6 to.

As shown in Fig. 4-1 and 4-2, be all sections of maximum correlation coefficients of trying to achieve of six axles of EPSON.Related coefficient is larger, and data waveform matching degree is higher.From figure, also can find out, after the 16th section, related coefficient, close to 1, show that EPSON and FSAS data waveform consistance are higher, and method identifies optimum matching section thus.

As shown in Fig. 5-1 and 5-2, the mistiming sequence being calculated by each section of maximum correlation coefficient, chooses the period of related coefficient when larger.In figure, can find out that the method can not only calculate the mistiming of EPSON, also can find that the frequency marking of EPSON exists drift.

S6.13: the poor bias of least square adjustment seeking time and markers drift parameter k.Using bias and k as solve for parameter, and Bias (i) is observed reading, and observation equation is as follows:

Bias(i)＝k×[time(i)-t ₀]+bias

Wherein the initial value Bias0 of bias is Bias (i) sequence intermediate value, and the corresponding time is t ₀, the initial value of k is 0.

The power of Bias (i) is max_corr (i).For removing the rough error in observed reading, set a correlation coefficient threshold σ _r, in this example, σ _rbe 0.9.Only have when max_corr (i) is greater than 0.9, observed reading is substitution observation equation, carries out least square adjustment and solves, and tries to achieve one group of bias and k.

S6.14: judge that whether δ s is greater than a sampling interval, if so, forwards step S6.15 to.Otherwise, by the δ s T_inter that doubles, forward step S6.3 to.

S6.15: select best bias and k combination.Utilize not bias and k on the same group, obtain the to be matched data of different time after synchronous.Calculate the covariance between data to be matched and reference data.If data are single shaft, that group bias while selecting covariance maximum and k are as best bias and k combination; If when data are multiaxis, be as the criterion with practical experience, select the most believable axis data to carry out covariance and solve, select best bias and k combination.

S7: according to best bias and k combination, upgrade data T to be matched ₀time, and storage.In this example, the most believable is the 4th column data, and after time synchronized, the 3rd group of covariance that data obtain is maximum.Utilize following formula to calculate and upgrade T ₀time:

t′＝(t-t ₀)×k+bias+t

Wherein, t is T ₀time row, t ₀for time corresponding to mistiming initial value of selecting in step 6.12, k is markers drift parameter, and bias is the mistiming, and t ' is the rear time for t upgrades.

By T ₀in t ' replacement for time row, and data are stored as and T again ₀the data file that data layout is consistent.

S8: finish.

As shown in Figure 6, for utilizing algorithm of the present invention to carry out after time synchronized EPSON and FSAS, vertical gyrobearing is marked lower oscillogram at one time.From then on figure can find out, this algorithm can correctly carry out time synchronized computing.

As shown in Figure 7, be the Waveform Matching partial, detailed view in certain period in Fig. 6.In figure, can find out, this algorithm can not only correctly carry out the time synchronized between data to be matched and reference data, and synchronization accuracy is higher.

A clock synchronization system based on Waveform Matching, is characterized in that, comprises as lower module:

Read module: for reading reference data R ₀and data T to be matched ₀;

Data processing format module: for unified reference data R ₀and data T to be matched ₀data processing form, extract the reference data R after conversion ₁and data T to be matched ₁;

Down-sampled module: for respectively to R ₁and T ₁down-sampled, to the reference data R after down-sampled ₂and data T to be matched ₂carry out simple crosscorrelation sequence and solve, corresponding interval while getting simple crosscorrelation sequence maximal value, obtains thick mistiming int_bias now;

Thick mistiming module: for according to the thick mistiming int_bias trying to achieve, make reference data R ₂time period cover data T to be matched ₂time period, and front and back respectively exceed data T to be matched ₂certain length, intercepting reference data R ₂and data T to be matched ₂;

Sampling rate computing module: for unified R ₂and T ₂sampling rate, obtains the reference data R of identical sampling rate ₃with data T to be matched ₃;

Computing module: for trying to achieve R ₃and T ₃the poor bias of optimal time and markers drift parameter k;

Update module: for according to best bias and k combination, upgrade data T to be matched ₀time, and storage.

Claims (7)

1. the method for synchronizing time based on Waveform Matching, is characterized in that, comprises the steps:

Step 1: read reference data R ₀and data T to be matched ₀;

Step 2: unified reference data R ₀and data T to be matched ₀data processing form, extract the reference data R after conversion ₁and data T to be matched ₁;

Step 3: respectively to R ₁and T ₁down-sampled, to the reference data R after down-sampled ₂and data T to be matched ₂carry out simple crosscorrelation sequence and solve, corresponding interval while getting simple crosscorrelation sequence maximal value, obtains thick mistiming int_bias now;

Step 4: the thick mistiming int_bias according to trying to achieve, makes reference data R ₂time period cover data T to be matched ₂time period, and front and back respectively exceed data T to be matched ₂certain length, intercepting reference data R ₂and data T to be matched ₂, enter step 5;

Step 5: unified R ₂and T ₂sampling rate, obtains the reference data R of identical sampling rate ₃with data T to be matched ₃;

Step 6: try to achieve R ₃and T ₃the poor bias of optimal time and markers drift parameter k;

Step 7: according to best bias and k combination, upgrade data T to be matched ₀time, and storage;

Step 8: finish.

2. a kind of method for synchronizing time based on Waveform Matching according to claim 1, is characterized in that: unified reference data R in described step 2 ₀and data T to be matched ₀the formula of data processing form as follows:

v＝Δs/Δt _s

Wherein, Δ s is gyro or the measured increment size of accelerometer between former and later two samplings, Δ t _sfor the time interval between former and later two samplings, v is the rate value that Δ s is corresponding;

Reference data R after described conversion ₁and data T to be matched ₁that the time row after changing and an axis data of carrying out Waveform Matching form.

3. a kind of method for synchronizing time based on Waveform Matching according to claim 1, is characterized in that: in described step 3, comprise the steps:

Step 3.1: respectively to R ₁and T ₁down-sampled, obtain the reference data R after down-sampled ₂and data T to be matched ₂;

Step 3.2: to down-sampled rear R ₂and T ₂among shorter sequence trailing zero, until both are equal in length;

Step 3.3: calculate R ₂with T ₂simple crosscorrelation sequence; Simple crosscorrelation sequence computing formula is as follows:

${\hat{R}}_{\mathrm{xy}}\left(m\right)=\left\{\begin{array}{cc}\underset{K=0}{\overset{N-m-1}{\mathrm{\&Sigma;}}}{x}_{K+m}{y}_m& m\&GreaterEqual;0\\ {\hat{R}}_{\mathrm{xy}}(-m)& m<0\end{array}\right.$

Wherein, x and y are the sequence of carrying out correlativity processing, and x and y length are N, and K and m are respectively used to represent the ordinal number of the value in sequence x and y, m value is-(N-1) to+(N-1);

Step 3.4: get simple crosscorrelation sequence corresponding interval during maximal value m, obtains thick mistiming int_bias now, and the thick mistiming can utilize following formula to try to achieve:

int_bias＝R ₂(t ₀)-T ₂(t ₀)-T ₂_interval×m

Wherein, R ₂(t ₀) be R ₂start time, T ₂(t ₀) be T ₂start time, T ₂_ interval is the time span of the interval representative that staggers, is T ₂sampling interval, T ₂_ interval * m is the time that two column datas stagger.

4. a kind of method for synchronizing time based on Waveform Matching according to claim 1, is characterized in that: in described step 6.15, comprise following sub-step:

Step 6.1: calculate sub-sampling T_inter interval time;

Step 6.2: take T_inter as factor, calculate the data to be matched δ s of translation distance to the right;

Step 6.3: judge whether δ s is the positive integer times of sampling interval; If so, forward step S6.15 to, if not, forward step S6.4 to;

Step 6.4: by data data T to be matched ₃translation δ s obtains T ₄;

Step 6.5: take Δ t as length of window, by T ₄be divided into limited n section, get present segment data T ₄(i), 1≤i≤n;

Step 6.6: judge whether data to be matched extract n section; If so, forward step S6.13 to; If not, forward step S6.7 to;

Step 6.7: set moving window W ₁, window is long is Δ t, starting point is the present segment zero hour;

Step 6.8: judge whether window sliding distance is greater than preseting length σ _s; As if so, forwarded step S6.11 to;

Otherwise, with W ₁for window, get reference data R ₃(W ₁);

Step 6.9: calculate T ₄and R (i) ₃(W ₁) between Pearson's related coefficient; Pearson's Calculation of correlation factor as shown in the formula:

$r=\frac{\underset{i=1}{\overset{N}{\mathrm{\&Sigma;}}}(x_i-\stackrel{\&OverBar;}{x})(y_i-\stackrel{\&OverBar;}{y})}{\sqrt{\underset{i=1}{\overset{N}{\mathrm{\&Sigma;}}}{(x_i-\stackrel{\&OverBar;}{x})}^2\&CenterDot;\underset{i=1}{\overset{N}{\mathrm{\&Sigma;}}}{(y_i-\stackrel{\&OverBar;}{y})}^2}}$

Wherein for x serial mean, for y serial mean, N is the length of sequence x and y, and r is the Pearson's related coefficient between sequence x and y;

Step 6.10: judge that whether reference data reads last position, if so, forwards step S6.11 to; Otherwise,

By window W ₁slide backward after the sampling interval of a reference data, forward step S6.8 to;

Step 6.11: calculate this section of mistiming Bias (i) and corresponding time time (i);

Step 6.12: extract next section of data to be matched.Take off one section of T ₄(i), forward step 6.6 to;

Step 6.13: with the poor bias of least square adjustment seeking time and markers drift parameter k; Using bias and k as solve for parameter, and Bias (i) is observed reading, and observation equation is as follows:

Bias(i)＝k×[time(i)-t ₀]+bias

Wherein the initial value Bias0 of bias is Bias (i) sequence intermediate value, and the corresponding time is t ₀, the initial value of k is 0;

Step 6.14: judge that whether δ s is greater than a sampling interval, if so, forwards step S6.15 to; Otherwise, by the δ s T_inter that doubles, forward step S6.3 to;

Step 6.15: select best bias and k combination; Utilize not bias and k on the same group, obtain the to be matched data of different time after synchronous.

5. a kind of method for synchronizing time based on Waveform Matching according to claim 4, is characterized in that: in described step 6.15, if data are single shaft, that group bias while selecting covariance maximum and k are as best bias and k combination; If data are multiaxis, select the most believable axis data to carry out covariance and solve, select best bias and k combination.

6. a kind of method for synchronizing time based on Waveform Matching according to claim 1, is characterized in that: in described step 7, utilize following formula to calculate and upgrade T ₀time:

t′＝(t-t ₀)×k+bias+t

Wherein, t is T ₀time row, t ₀for time corresponding to mistiming initial value of selecting in step 6.12, k is markers drift parameter, and bias is the mistiming, and t ' is the rear time for t upgrades.

7. the clock synchronization system based on Waveform Matching, is characterized in that, comprises as lower module:

Read module: for reading reference data R ₀and data T to be matched ₀;

Data processing format module: for unified reference data R ₀and data T to be matched ₀data processing form, extract the reference data R after conversion ₁and data T to be matched ₁;

Down-sampled module: for respectively to R ₁and T ₁down-sampled, to the reference data R after down-sampled ₂and data T to be matched ₂carry out simple crosscorrelation sequence and solve, corresponding interval while getting simple crosscorrelation sequence maximal value, obtains thick mistiming int_bias now;

Interception module: for according to the thick mistiming int_bias trying to achieve, make reference data R ₂time period cover data T to be matched ₂time period, and front and back respectively exceed data T to be matched ₂certain length, intercepting reference data R ₂and data T to be matched ₂;

Sampling rate computing module: for unified R ₂and T ₂sampling rate, obtains the reference data R of identical sampling rate ₃with data T to be matched ₃;

Computing module: for trying to achieve R ₃and T ₃the poor bias of optimal time and markers drift parameter k;

Update module: for according to best bias and k combination, upgrade data T to be matched ₀time, and storage.