CN110084185A - A kind of bullet train slightly crawls the rapid extracting method of operation characteristic - Google Patents
A kind of bullet train slightly crawls the rapid extracting method of operation characteristic Download PDFInfo
- Publication number
- CN110084185A CN110084185A CN201910336833.8A CN201910336833A CN110084185A CN 110084185 A CN110084185 A CN 110084185A CN 201910336833 A CN201910336833 A CN 201910336833A CN 110084185 A CN110084185 A CN 110084185A
- Authority
- CN
- China
- Prior art keywords
- matrix
- bullet train
- slightly
- signal
- formula
- Prior art date
- Legal status (The legal status is an assumption and is not a legal conclusion. Google has not performed a legal analysis and makes no representation as to the accuracy of the status listed.)
- Granted
Links
- 238000000034 method Methods 0.000 title claims abstract description 66
- 239000011159 matrix material Substances 0.000 claims abstract description 135
- 230000004927 fusion Effects 0.000 claims abstract description 12
- 238000000354 decomposition reaction Methods 0.000 claims abstract description 9
- 241000270295 Serpentes Species 0.000 claims description 23
- 238000000605 extraction Methods 0.000 claims description 18
- 230000010355 oscillation Effects 0.000 claims description 18
- 239000010977 jade Substances 0.000 claims description 13
- 230000009467 reduction Effects 0.000 claims description 11
- 238000005070 sampling Methods 0.000 claims description 9
- 230000007704 transition Effects 0.000 claims description 9
- 238000012545 processing Methods 0.000 claims description 6
- 230000017105 transposition Effects 0.000 claims description 4
- 238000012935 Averaging Methods 0.000 claims description 3
- 238000010276 construction Methods 0.000 claims description 3
- 238000001914 filtration Methods 0.000 claims description 3
- 238000010606 normalization Methods 0.000 claims description 3
- 230000009466 transformation Effects 0.000 claims description 3
- 238000012952 Resampling Methods 0.000 claims description 2
- 102000002274 Matrix Metalloproteinases Human genes 0.000 claims 1
- 108010000684 Matrix Metalloproteinases Proteins 0.000 claims 1
- 238000012706 support-vector machine Methods 0.000 abstract description 4
- 230000006870 function Effects 0.000 description 12
- 230000001133 acceleration Effects 0.000 description 10
- 239000000284 extract Substances 0.000 description 7
- 238000012360 testing method Methods 0.000 description 7
- 238000012544 monitoring process Methods 0.000 description 5
- 230000008569 process Effects 0.000 description 5
- 230000009193 crawling Effects 0.000 description 4
- 238000011160 research Methods 0.000 description 4
- 238000010586 diagram Methods 0.000 description 3
- 238000004458 analytical method Methods 0.000 description 2
- 230000008901 benefit Effects 0.000 description 2
- 238000011161 development Methods 0.000 description 2
- 238000002474 experimental method Methods 0.000 description 2
- 238000010801 machine learning Methods 0.000 description 2
- 238000005457 optimization Methods 0.000 description 2
- 238000012549 training Methods 0.000 description 2
- 241001269238 Data Species 0.000 description 1
- 241000271559 Dromaiidae Species 0.000 description 1
- 102000008297 Nuclear Matrix-Associated Proteins Human genes 0.000 description 1
- 108010035916 Nuclear Matrix-Associated Proteins Proteins 0.000 description 1
- 230000002159 abnormal effect Effects 0.000 description 1
- 230000002776 aggregation Effects 0.000 description 1
- 238000004220 aggregation Methods 0.000 description 1
- 238000004364 calculation method Methods 0.000 description 1
- 230000008859 change Effects 0.000 description 1
- 238000013480 data collection Methods 0.000 description 1
- 230000004069 differentiation Effects 0.000 description 1
- 230000000694 effects Effects 0.000 description 1
- 230000008030 elimination Effects 0.000 description 1
- 238000003379 elimination reaction Methods 0.000 description 1
- 238000005516 engineering process Methods 0.000 description 1
- 239000012634 fragment Substances 0.000 description 1
- 238000005259 measurement Methods 0.000 description 1
- 239000000203 mixture Substances 0.000 description 1
- 210000000299 nuclear matrix Anatomy 0.000 description 1
- 238000005096 rolling process Methods 0.000 description 1
Classifications
-
- G—PHYSICS
- G06—COMPUTING; CALCULATING OR COUNTING
- G06F—ELECTRIC DIGITAL DATA PROCESSING
- G06F2218/00—Aspects of pattern recognition specially adapted for signal processing
- G06F2218/02—Preprocessing
-
- G—PHYSICS
- G06—COMPUTING; CALCULATING OR COUNTING
- G06F—ELECTRIC DIGITAL DATA PROCESSING
- G06F2218/00—Aspects of pattern recognition specially adapted for signal processing
- G06F2218/08—Feature extraction
-
- G—PHYSICS
- G06—COMPUTING; CALCULATING OR COUNTING
- G06F—ELECTRIC DIGITAL DATA PROCESSING
- G06F2218/00—Aspects of pattern recognition specially adapted for signal processing
- G06F2218/12—Classification; Matching
Landscapes
- Train Traffic Observation, Control, And Security (AREA)
- Measurement Of Mechanical Vibrations Or Ultrasonic Waves (AREA)
Abstract
It slightly crawls the rapid extracting method of operation characteristic the invention discloses a kind of bullet train.Aiming at the problem that existing method can not quickly identify the small size hunting state of bullet train (operate normally, slightly restrain, slightly dissipate and substantially crawl), the present invention decomposes pretreated signal using average empirical mode decomposition (EEMD) method, its result is converted into energy matrix, then the joint approximate diagonalization (JAD) under nonstationary condition is carried out to energy matrix to handle, the energy matrix for merging multiple sensors obtains fused bullet train and slightly crawls operating status feature.It is trained and identifies by the way that fusion feature is put into least square method support vector machines, demonstrate this method can quickly and accurately by bullet train operate normally, slightly convergence, slightly diverging and substantially crawl four kinds of operating statuses separate, to ensure train operating safety.
Description
Technical field
The invention belongs to train operation state Feature Extraction Technology fields, and in particular to a kind of bullet train slightly crawls fortune
The rapid extracting method of row feature.
Background technique
Bullet train is rapid in Chinese development, ends 2018, and China is high-speed rail mileage open to traffic and in Jianli (CV 11) journey in the world
Most countries, high-speed rail brings convenience to people's trip, but with the raising of high-speed rail speed, operational safety stability becomes
Focus of attention.
During bullet train operation, the study found that after the travel speed of vehicle is more than a certain critical value, high speed
Train can generate unstable hunting.Snakes model can generate biggish lateral wheel force, and it is de- to will lead to train when serious
Rail.Therefore, the significant obstacle of Snakes model bullet train safe operation.Currently, to vehicle Snakes model, there is no unite both at home and abroad
One judgment criteria, most of Overseas Correlative Standards be by wheel shaft cross force in train traveling process, wheel-rail lateral force and
Vehicle lateral acceleration etc. is to determine whether there are Snakes models.The railroad coach traffic safety monitoring standard in China also uses structure
This index of frame transverse acceleration evaluates the lateral stability of rolling stock, " high-speed EMUs complete vehicle test specification " rule
Determine train bogie transverse acceleration after 10Hz is filtered, meet or exceed limiting value continuous 6 times, then can be determined that column
Vehicle loses hunting stability.The existing monitoring of peak method in China is all based on above-mentioned specification and is monitored to bullet train, with true
Determine in train travelling process with the presence or absence of Snakes model phenomenon.
The discovery when being analyzed and processed to experimental data, train also easily occur small size hunting in high-speed cruising,
Slightly crawl.Slightly snake refers to that during train operation the amplitude of framework transverse acceleration signal is not up to safe pole
The part of limit value.Slightly snake has many hazards, for example influences the comfort that passenger takes, causes wheel-rail wear and tear to aggravate, reduce
The active time of train, and violent harm can be brought to the safe operation of train, it is violent Snakes model that slightly snake is abnormal
Sign.Present invention discovery when slightly the state of crawling is studied to bullet train, according to the development of the small size hunting of train,
Small size snake can be divided into two classes: slightly convergence and slightly diverging.The small size convergence state of bullet train can slowly converge to just
Normal operating status, and the oscillation crosswise amplitude of the small size divergent state of bullet train then can constantly expand, and lose until developing into snake
Surely, to be caused huge harm to bullet train operation, therefore, the two stages slightly crawled are differentiated in time, train is transported
Row stability is most important.
Such as Fig. 1, Fig. 1 (a) is that bullet train slightly crawls convergence graph, and Fig. 1 (b) is that bullet train slightly crawls scatter graph,
As can be seen that bogie lateral amplitude of vibration is gradually reduced to normal operating condition in Fig. 1 (a), and bogie transverse-vibration in Fig. 1 (b)
Width is then gradually increased to Snakes model state.Therefore, it is necessary that slightly snake two states are differentiated in time.
It slightly crawls for bullet train differentiation, researcher has carried out relevant research to it.Ning J proposes one
Kind of the bullet train based on multiple dimensioned arrangement entropy and local tangent space alignment slightly crawls Characteristics of Evolution extracting method, and the method is quasi-
The operation characteristic of multiple operating statuses really is extracted, but the operation time for needing to consume is longer.The speed of service of bullet train
It is exceedingly fast, too long operation time makes feature extraction result that can not rapidly be reflected to train, and train is caused to miss best tune
Whole window.To solve the above-mentioned problems, it slightly crawls the rapidly extracting side of operation characteristic the invention proposes a kind of bullet train
Method.
At the same time, in test when carrying out data collection to running state of high speed strain, what is generally used is pass more
The mode of sensor network system, and divide to bullet train Snakes model phenomenon and bullet train phenomenon of slightly crawling
When analysis processing, certain one piece of data information that many researchers have selected in a certain position sensor is researched and analysed, and this
When show whether research conclusion can reflect that whole fault messages of train need to be discussed.For this purpose, being directed to High Speed Train in China
Urgent problem to be solved in hunting research, the present invention is combined with the theoretical analysis of and engineering practice, with the evolution Feature that slightly crawls
Extracting method is research core, combines the vibration data of the multiple position sensors of train, steering when running according to bullet train
The measurement data of framework frame acceleration signal extracts the small size snake Characteristics of Evolution of bullet train, is High Speed Train in China state
Monitoring and safe early warning provide reference frame.
Summary of the invention
A kind of slightly the crawl rapid extracting method of operation characteristic of bullet train provided by the invention solves the prior art
The feature extraction of high speed train operation state not exclusively, can not rapidly extract train slightly snakelike Characteristics of Evolution the problem of.
In order to achieve the above object of the invention, the technical solution adopted by the present invention are as follows: a kind of bullet train slightly crawls operation
The rapid extracting method of feature:
Original vibration signal when S1, extraction bullet train operation, and it is pre-processed;
S2, pretreated vibration signal is decomposed into several IMF components by EEMD;
S3, according to IMF component, construct the corresponding energy matrix of vibration signal;
S4, energy matrix is handled to obtain fusion feature matrix by the JADE method of unstable condition, i.e., at a high speed
The small size snake operation characteristic of train.
Further, the original vibration signal in the step S1 includes 1 axle box oscillation crosswise signal, 1 framework cross
To vibration signal and 4 framework oscillation crosswise signals;
Further, the step S2 specifically:
S21, the Gaussian sequence ω (t) for pretreated vibration signal a (t) being added identical amplitude, obtain total
Body signal a'(t);
Wherein, overall signal a'(t) are as follows:
A'(t)=a (t)+ω (t)
S22, the overall signal a'(t after each addition white noise) is decomposed according to EMD method, determines each rank point
Measure ci:
Wherein, overall signal a'(t) decomposition formula are as follows:
Wherein: ciFor each rank IMF component;
I indicates i-th of IMF component;
R is residual error item;
N is the IMF number decomposited.
S23, the different white noise sequence ω that identical amplitude is added every timej(t), step S21-S22 is repeated, is obtained:
In formula: aj' (t) be overall signal after the Gaussian sequence that jth time is added;
N is the number that Gaussian sequence is added;
ωj(t) Gaussian sequence being added for jth time;
cijI-th of the IMF component decomposited when Gaussian sequence is added for jth time;
rjThe residual values decomposited when Gaussian sequence is added for jth time.
S24, the zero-mean principle according to white Gaussian noise frequency eliminate influence of the white Gaussian noise to IMF component, obtain
The corresponding IMF component of pretreated original vibration signal are as follows:
In formula: ciIt (t) is the corresponding IMF component of i-th of original vibration signal;
cij(t) i-th of IMF component to be decomposed to original vibration signal progress EEMD.
Further, the step S3 specifically:
S31, IMF energy square is constituted according to IMF component are as follows:
In formula: P is total sampling number;
M is sampled point;
ciFor the corresponding IMF component of vibration signal of selection;
Δ t is the sampling period.
S32, according to IMF energy square, the feature vector M of the IMF energy square after construction normalization are as follows:
In formula: q=1 ..., Q is the sample size of the original vibration signal of a position;
N is the IMF number decomposited;
S33, according to step S31-S32, by the corresponding feature of the IMF component of all vibration signals of vibration signal of the same race
Vector M forms energy matrix E0;
S34, according to step S31-S33, respectively obtain 1 axle box oscillation crosswise signal, 1 framework oscillation crosswise signal and
The corresponding energy matrix E of 4 upper framework oscillation crosswise signals1, energy matrix E2With energy matrix E3;
S35, by energy matrix E1, energy matrix E2, energy matrix E3It joins together, obtains the corresponding energy of vibration signal
Matrix E.
Further, the step S4 specifically:
S41, time series T is resolved into K sections of time interval TK, and then generate corresponding K covariance matrix
Wherein, time series T is corresponding timing node when extracting original vibration signal;
S42, the corresponding diagonalizable matrix of K covariance matrix is determined
Wherein, diagonalizable matrixAre as follows:
In formula, k=2 ..., K;
For first covariance matrixDiagonalizable matrix, andV is first association side
Poor matrixEigenmatrix, first covariance matrix of ΛFeature vector, subscript H be conjugate transposition operator;
S43, according to diagonalizable matrixDetermine that a unitary matrice U makes following formula value maximum;
Wherein: function diag indicates diagonal element;
P is the dimension after dimensionality reduction;
The b row and d of b and d representing matrix arrange.
S44, given threshold ε, according to spin matrix to unitary matrice U and covariance matrixIt is iterated update, Zhi Daogeng
Unitary matrice U and covariance matrix after newIn the value of all off diagonal elements when being respectively less than threshold epsilon, obtain final
Unitary matriceIt completes to K covariance matrixJoint diagonalization;
Wherein, k=2 ... K;
Spin matrix G (i, j, θ) are as follows:
In formula: i and j is respectively the ith row and jth column of spin matrix;
θ is intermediate computations parameter;
S45, according to final unitary matriceDetermine transition matrix A;
Wherein, transition matrix A are as follows:
In formula: subscript # is pseudo-inverse transformation operator;
S46, the fusion feature matrix Z according to transition matrix A, after being decomposed;
Z=AE.
Further, which is characterized in that in the step S44, to the iterative formula of unitary matrice U are as follows:
U←UG(1,2,θ)
In formula: G (1,2, θ) is the element value of spin matrix G (i, j, θ) the 1st row the 2nd column;
To covariance matrixIterative formula are as follows:
Further, the small size snake operation characteristic of the bullet train in the step S4 includes normal operating condition spy
Sign, small size convergence state feature, small size divergent state feature and substantially snakelike state feature.
The rapid extracting method of operation characteristic the invention has the benefit that bullet train provided by the invention slightly crawls
Bullet train run signal is analyzed, is successfully operated normally bullet train, slightly convergence and substantially snakelike four kinds of operations
State separates, and quickly and accurately detects the operating status of bullet train, and using the LSSVM using machine learning to bullet train
Operating status feature extracts, to predict whether bullet train will appear the instability status substantially to crawl, helps to improve
The accuracy and timeliness of train monitoring.
Detailed description of the invention
Fig. 1 is background of invention high speed train slightly snakelike operation characteristic schematic diagram.
Fig. 2 is that bullet train in the present invention slightly crawls the rapid extracting method flow chart of operation characteristic.
Fig. 3 is the time-domain diagram that signal x (t), y (t), z (t) are emulated in embodiment provided by the invention.
Fig. 4 is emulation signal JADE clustering effect picture.
Fig. 5 is train speed curve, bogie frame transverse acceleration curve graph in embodiment provided by the invention.
Fig. 6 is train normal operation signal EEMD decomposition result figure in embodiment provided by the invention.
Fig. 7 is that train slightly restrains signal EEMD decomposition result figure in embodiment provided by the invention.
Fig. 8 is the small size diverged signal EEMD decomposition result figure of train in embodiment provided by the invention.
Fig. 9 is train substantially snakelike signal EEMD decomposition result figure in embodiment provided by the invention.
Figure 10 is the JADE feature extraction figure of position 1 in embodiment provided by the invention.
Figure 11 is the JADE feature extraction figure of position 2 in embodiment provided by the invention.
Figure 12 is the JADE feature extraction figure of position 3 in embodiment provided by the invention.
Figure 13 is that multiposition merges JADE feature extraction figure in embodiment provided by the invention.
Specific embodiment
A specific embodiment of the invention is described below, in order to facilitate understanding by those skilled in the art this hair
It is bright, it should be apparent that the present invention is not limited to the ranges of specific embodiment, for those skilled in the art,
As long as various change is in the spirit and scope of the present invention that the attached claims limit and determine, these variations are aobvious and easy
See, all are using the innovation and creation of present inventive concept in the column of protection.
The rapid extracting method of operation characteristic as shown in Fig. 2, a kind of bullet train slightly crawls, comprising the following steps:
Original vibration signal when S1, extraction bullet train operation, and it is pre-processed;
Wherein, original vibration signal includes 1 axle box oscillation crosswise signal, 1 framework oscillation crosswise signal and 4 upper structures
Frame oscillation crosswise signal;
Pretreated method is carried out to original vibration signal specifically: successively adopted again to the original vibration signal of extraction
The processing such as sample processing, bandpass filtering treatment, zero averaging processing and removal trend term.
S2, pretreated vibration signal is decomposed into several IMF components by EEMD;
S3, according to IMF component, construct the corresponding energy matrix of vibration signal;
S4, energy matrix is handled to obtain fusion feature data matrix by the JADE method of unstable condition, i.e.,
The small size snake operation characteristic of bullet train.
The small size snake operation characteristic of bullet train therein includes normal operating condition feature, small size convergence state spy
Sign, small size divergent state feature and substantially snakelike state feature.
Wherein, step S2 specifically:
S21, the Gaussian sequence ω (t) for pretreated vibration signal a (t) being added identical amplitude, obtain total
Body signal a'(t);
Wherein, overall signal a'(t) are as follows:
A'(t)=a (t)+ω (t)
S22, the overall signal a'(t after each addition white noise) is decomposed according to EMD method, determines each rank point
Measure ci:
Wherein, overall signal a'(t) decomposition formula are as follows:
In formula: ciFor each rank IMF component;
I indicates i-th of IMF component;
R is residual error item;
N is the IMF number decomposited.
S23, the different white noise sequence ω that identical amplitude is added every timej(t), step S21-S22 is repeated, is obtained:
In formula: a'j(t) overall signal after the Gaussian sequence being added for jth time;
N is the number that Gaussian sequence is added;
ωj(t) Gaussian sequence being added for jth time;
cijI-th of the IMF component decomposited when Gaussian sequence is added for jth time;
rjThe residual values decomposited when Gaussian sequence is added for jth time.
S24, the zero-mean principle according to white Gaussian noise frequency eliminate influence of the white Gaussian noise to IMF component, obtain
The corresponding IMF component of pretreated original vibration signal are as follows:
In formula: ciIt (t) is the corresponding IMF component of i-th of original vibration signal;
cij(t) i-th of IMF component to be decomposed to original vibration signal progress EEMD.
When bullet train operation, signal is generally non-stationary signal, and the method based on traditional IMF energy is embodying letter
Result is unsatisfactory when the feature hidden in number, and present invention uses IMF energy Moment Methods to extract bullet train energy spy
Sign, therefore step S3 specifically:
S31, IMF energy square is constituted according to IMF component are as follows:
In formula: P is total sampling number;
M is sampled point;
ciFor the corresponding IMF component of vibration signal of selection;
Δ t is the sampling period.
S32, according to IMF energy square, the feature vector M of the IMF energy square after construction normalization are as follows:
In formula: q=1 ... Q is the sample size of the original vibration signal of a position;
N is the number of the IMF component decomposed;
S33, according to step S31-S32, by the corresponding feature of the IMF component of all vibration signals of vibration signal of the same race
Vector M forms energy matrix Ε0;
S34, according to step S31-S33, respectively obtain 1 axle box oscillation crosswise signal, 1 framework oscillation crosswise signal and
The corresponding energy matrix E of 4 upper framework oscillation crosswise signals1, energy matrix E2With energy matrix E3;
S35, by energy matrix E1, energy matrix E2With energy matrix E3It joins together, obtains the corresponding energy of vibration signal
Matrix E.
The feature vector of the original vibration signal sample of three positions is [M1;M2;...;MQ], it can composition energy square
Battle array E1、E2And E3, the dimension of energy matrix is higher, and feature is unobvious.It is non-stationary when in view of bullet train operation, in order to
Bullet train vibration signal characteristics are preferably extracted, herein using the JADE algorithm used under unstable condition to three positions
Energy matrix carry out fusion feature extraction.
It is non-stationary in view of bullet train signal, in order to use JADE method, this hair under non-stationary environment
It is bright that entire time series T is resolved into K sections of time interval T1,...,TK, and then K covariance matrix can be generated
The method of the present invention will carry out Joint diagonalization to this K covariance matrix, can be by K covariance matrix diagonalization to find one
Unitary matrice U, and then to energy matrix E dimensionality reduction and extract dimensionality reduction characteristic;Therefore step S4 specifically:
S41, time series T is resolved into K sections of time interval TK, and then generate corresponding K covariance matrix
Wherein, time series T is corresponding timing node when extracting original vibration signal;
For time interval T1,...,TK, covariance matrix are as follows:
Wherein: E () is desired operator;
EtFor time t ∈ TkInterior energy matrix E;
S42, K covariance matrix is determinedCorresponding diagonalizable matrix
Want the multiple matrixes of simultaneous diagonalizationFirstly the need of first matrix of diagonalization again to remaining K-1
A matrix carries out diagonalization, and matrix W is covariance matrixDiagonalizable matrix;
The wherein diagonalizable matrix of first covariance matrix are as follows:
Wherein, V is first covariance matrixEigenmatrix, first covariance matrix of ΛFeature vector,
Subscript H is conjugate transposition operator;
Therefore, diagonalizable matrixAre as follows:
In formula, k=2 ..., K;
For first covariance matrixDiagonalizable matrix;
S43, according to diagonalizable matrixDetermine that a unitary matrice U makes following formula value maximum;
It can be converted to the problem of joint diagonalization matrix in this step and find a unitary matrice U and make following formula minimum:
In formula: the quadratic sum of function off expression off diagonal element;
B, d respectively indicates the b row and d column of data matrix.
Since the quadratic sum of above formula remains unchanged when multiplied by orthogonal matrix, we can equally be maximized diagonally
The quadratic sum of element, i.e. problem, which are transformed into, to be found a unitary matrice U and makes following formula maximum, therefore is had:
Wherein: function diag indicates diagonal element;
P is the dimension after dimensionality reduction;
B and d is the b row and d column of matrix.
S44, given threshold ε, according to spin matrix to unitary matrice U and covariance matrixIt is iterated update, Zhi Daogeng
Unitary matrice U and covariance matrix after newIn the value of all off diagonal elements when being respectively less than threshold epsilon, obtain final
Unitary matriceIt completes to K covariance matrixJoint diagonalization;
The method of the present invention is rotated using Jacobi is changed into diagonal matrix, the base of Jacobi method for multiple covariance matrixes
This thought be by a sub-orthogonal transforation, by the off-diagonal element of a pair of of non-zero in objective matrix be melted into zero and make it is non-right
The quadratic sum of angle element reduces.The above process is repeated, the quadratic sum of the off-diagonal element of transformed matrix is made to go to zero,
To make the approximate matrix diagonal matrix, All Eigenvalues and feature vector are obtained.
The above method specifically: substituted into unit matrix I as the initial value of matrix U, the value of θ is then by matrix's
ElementIt is calculated,Using as follows
Matrix group is transformed into diagonal form by spin matrix G (i, j, θ);
In formula: i and j is respectively the ith row and jth column of spin matrix;
θ is intermediate computations parameter.
Wherein, to the iterative formula of unitary matrice U are as follows:
U←UG(1,2,θ)
In formula, G (1,2, θ) is the element value of spin matrix G (i, j, θ) the 1st row the 2nd column;
To covariance matrixIterative formula are as follows:
S45, according to final unitary matriceDetermine transition matrix A;
Wherein, transition matrix A are as follows:
In formula, subscript # is pseudo-inverse transformation operator;
S46, the fusion feature matrix Z according to transition matrix A, after being decomposed;
Z=AE.
The fusion feature matrix of extraction is verified using least square method support vector machines (LSSVM) in the present invention,
It is input to (LSSVM) to be trained and test, the accuracy of test feature data, reinforces the accuracy of prediction;The side LSSVM
The core of method replaces quadratic programming problem with solution system of linear equations, to avoid insensitive loss function, this method is reduced
Computational complexity, increases operation efficiency.
In this process, taking dimensionality reduction characteristic matrix Z is sample, because the first dimension data in Z includes most letter
Breath, therefore the first n dimensional vector n z of dimensionality reduction matrix Z is taken, have: { zq,yq, q=1 ..., Q, y ∈ [- 1,1], z are dimensionality reduction characteristic,
Y is tag along sort, and Q is sample size, and the function of the hyperplane of least square method support vector cassification is as follows:
Wherein:For adjustable weighted vector;
G is biasing;
φ (z) is nonlinear function, can handle sample linearly inseparable problem, and basic thought is the input space
Interior linear inseparable data are mapped to the feature vector in high-dimensional feature space, so that data become to divide.
The optimization problem of least square method support vector machines may be expressed as:
In formula: z is input vector, i.e. dimensionality reduction characteristic;
ξl> 0 is slack variable, for measuring departure degree;
γ is penalty factor;
Subscript H is conjugate transposition.
Above formula can be obtained using Lagrange multiplier:
Wherein: β is Lagrange factor;
And then it is right respectivelyAsk local derviation that can obtain following formula:
It is obtained by above formula:
Wherein, I is unit matrix;
Ωql=yqylφ(zl)Hφ(zq)=yqylD(zq,zl), q, l=1 ..., η is nuclear matrix;
For the kernel function for meeting Mercer theorem, q, o indicate the element of sample z, σ table
Show the parameter of control kernel function width.
By above formula, LSSVM optimization problem can be exchanged into the Solve problems of system of linear equations, therefore LSSVM classifier
Hyperplane function is writeable are as follows:
In formula: f (z) is least square method support vector machines objective function;
Sgn () is jump function;
βlFor Lagrange factor, l=1 ..., η are βlNumber;
ylFor tag along sort, y ∈ [- 1,1];
D(z,zl) it is the kernel function for meeting Mercer theorem;
G is biasing.
The different motion state of bullet train is extracted in prediction in time by the above process, prevents bullet train from snakelike mistake occurs
Surely, the driving safety of train has been ensured.
In the first embodiment of the present invention, in order to verify the superiority of proposition method, herein to this method emulation testing.
Selection analog signal x (t), y (t), z (t) analyze this feature extracting method, are shown below.Emulation signal x (t),
Y (t), z (t) are to be added by the original signal amplitude modulation that three frequencies are 10Hz, 20Hz, 30Hz, and on this basis, give
The white Gaussian noise that three emulation signal addition signal-to-noise ratio are 10.Sample frequency is 500Hz, the time-domain diagram of three obtained signal
It is illustrated in fig. 3 shown below.
Emulation signal is handled, and dimension-reduction treatment, available low-dimensional feature according to said frame, and Projection Character is existed
Two-dimensional space is illustrated in fig. 4 shown below;
The respective feature height aggregation for emulating signal x (t), y (t), z (t), can only see a point in two-dimensional space,
Also, the feature permutation position difference between difference emulation signal is larger, can tell the feature of unlike signal well, know
Other excellent, and then confirm the correctness of the method.
In the second embodiment of the present invention, provides and slightly crawl operating status with the method for the present invention to bullet train
The example of feature extraction:
(1) obtain experimental data: the Mr. Yu's type EMU scientific experiment of experimental data source, signal sampling frequencies are
2500Hz, train speed information are provided by onboard wireless GPS, and acceleration information is by 2 vehicle 1 axle box lateral acceleration sensor, 2
Multiple sensors such as 1 framework lateral acceleration sensor of vehicle, 2 vehicle, 4 axle box upper framework transverse accelerations provide.Route uses
The rich format non-fragment orbit technology of Germany, rail are the domestic rail of fixed length 100m.Wherein train time-rate curve and bogie
When framework shown in m- accelerating curve such as Fig. 5 (a), time overall length is 2491 seconds.It is observed that the experiment train speed 330~
When 350km/h, occurs Snakes model repeatedly.
The initial data obtained from bogie frame sensor is relatively rough, needs to do it certain pretreatment.Due to
The frequency of bullet train hunting is 2~12.07Hz, carries out the bandpass filtering of 2~12.07Hz, to initial data with effective
A large amount of interference signal in original signal is eliminated on ground, improves the accuracy and efficiency of feature identification, and according to Shannon's sampling theorem,
The resampling that sample frequency is 250Hz is carried out to initial data, later to data zero averaging, elimination trend term and smooth place
Reason, to obtain smooth accurate initial data.
(2) EEMD is decomposed: the initial data of a position chooses the bullet train speed of service in 330~350km/h,
It operates normally, small size convergence, slightly dissipate, substantially crawl each 10 samples of four kinds of states, total 40 samples, three position originals
Beginning data choose altogether 120 samples.Sample data is too long, can increase calculation amount and operation time;Sample data is short, then not
The operating status of train can be completely shown, taking into account the above, in conjunction with reality, selecting sample length herein is 500, the sampling time
It is 2 seconds.
The sample that each initial data obtains separately is handled, bullet train is operated normally, slightly convergence, small size hair
Scattered, four kinds of states of substantially crawling signals do EEMD decomposition, and Fig. 6 to Fig. 9 is that bullet train is operated normally, slightly restrained, slightly
The IMF component that diverging, four kinds of operating statuses of substantially crawling are decomposed through EEMD, it can be seen that original signal energy focuses primarily upon
Preceding several IMF components.
(3) IMF after decomposing EEMD extracts the special medical treatment of IMF energy, is converted to IMF energy square, and it is special to extract IMF energy
Sign, the energy matrix E of available three positions1、E2、E3, and join together to obtain energy matrix E, due to energy feature dimension
It is 9, dimension is too high to cause the bad identification of energy feature.
(4) matrix signal, joint pair are vibrated using the bullet train of the JADE method fusion different location of unstable condition
The multiple energy matrixs of angling, fusion feature extract to obtain the feature of bullet train difference operating status.Data are dropped into three-dimensional, and
By Projection Character to Three Dimensional Interface, as shown in figure 13.As a comparison, Figure 10-12 is respectively axle box position (position 1), 1 framework
Position (position 2), 4 framework positions (position 3) JADE feature extraction figure.Also, set forth herein the fast of method in order to protrude
Victory, table 1 are that the runing time of multiple methods compares, and PC platform CPU is Intel Core i5-4460,12GB memory, video card
For NVIDIA GeForce GT720;
Table 1: the runing time comparison of multiple methods
(5) in order to which the accuracy for verifying feature extraction verifies feature extraction using dimensionality reduction feature as the input of LSSVM
Accuracy.Split data into training, two classes of test, every each 20 groups of data of class (5 groups of normal operations, 5 groups slightly convergence, 5 groups it is small size
Diverging, 5 groups substantially crawl).As a comparison, the feature that other methods generate also is put into together inside LLSVM herein, half is made
For training, the other half is as test.Since multiple position bullet train vibration datas being utilized herein, and utilize unstable condition
Under the JADE method fusion dimensionality reduction data of multiple positions, the discrimination of method proposed in this paper be it is highest, reach
100%, this also demonstrates the accuracy and superiority of feature extracting method provided by the invention.
The rapid extracting method of operation characteristic the invention has the benefit that bullet train provided by the invention slightly crawls
Bullet train run signal is analyzed, is successfully operated normally bullet train, slightly convergence and substantially snakelike four kinds of operations
State separates, and detects the operating status of bullet train in time and accurately, and using the LSSVM using machine learning to bullet train
Slightly snake operation characteristic extracts, to predict whether bullet train will appear the instability status substantially to crawl, facilitates
Improve the accuracy and timeliness of train monitoring.
Claims (7)
- The rapid extracting method of operation characteristic 1. a kind of bullet train slightly crawls, which comprises the following steps:Original vibration signal when S1, extraction bullet train operation, and it is pre-processed;S2, pretreated vibration signal is decomposed into several IMF components by EEMD;S3, according to IMF component, construct the corresponding energy matrix of vibration signal;S4, energy matrix is handled to obtain fusion feature matrix, i.e. bullet train by the JADE method of unstable condition Small size snake operation characteristic.
- 2. the operation characteristic rapid extracting method that bullet train according to claim 1 slightly crawls, which is characterized in that institute Stating the original vibration signal in step S1 includes 1 axle box oscillation crosswise signal, 1 framework oscillation crosswise signal and 4 upper structures Frame oscillation crosswise signal;Pretreated method is carried out to original vibration signal specifically: successively carry out at resampling to the original vibration signal of extraction Reason, bandpass filtering treatment, zero averaging processing and removal trend term processing.
- The rapid extracting method of operation characteristic 3. bullet train according to claim 2 slightly crawls, which is characterized in that institute State step S2 specifically:S21, the Gaussian sequence ω (t) for pretreated vibration signal a (t) being added identical amplitude, obtain overall letter Number a'(t);Wherein, overall signal a'(t) are as follows:A'(t)=a (t)+ω (t)S22, the overall signal a'(t after each addition white noise) is decomposed according to EMD method, determines each order component ci:Wherein, overall signal a'(t) decomposition formula are as follows:In formula: ciFor each rank IMF component;I indicates i-th of IMF component;R is residual error item;N is the IMF number decomposited;S23, the different white noise sequence ω that identical amplitude is added every timej(t), step S21-S22 is repeated, is obtained:In formula: a'j(t) overall signal after the Gaussian sequence being added for jth time;N is the number that Gaussian sequence is added;ωj(t) Gaussian sequence being added for jth time;cijI-th of the IMF component decomposited when Gaussian sequence is added for jth time;rjThe residual values decomposited when Gaussian sequence is added for jth time;S24, the zero-mean principle according to white Gaussian noise frequency eliminate influence of the white Gaussian noise to IMF component, obtain pre- place The corresponding IMF component of original vibration signal after reason are as follows:In formula;ciIt (t) is the corresponding IMF component of i-th of original vibration signal;cij(t) i-th of IMF component to be decomposed to original vibration signal progress EEMD.
- The rapid extracting method of operation characteristic 4. bullet train according to claim 3 slightly crawls, which is characterized in that institute State step S3 specifically:S31, IMF energy square is constituted according to IMF component are as follows:In formula: P is total sampling number;M is sampled point;ciFor the corresponding IMF component of vibration signal of selection;Δ t is the sampling period;S32, according to IMF energy square, the feature vector M of the IMF energy square after construction normalization are as follows:In formula: q=1 ..., Q is the sample size of the original vibration signal of a position;N is the IMF number decomposited;S33, according to step S31-S32, by the corresponding feature vector of the IMF component of all vibration signals of vibration signal of the same race M forms energy matrix E0;S34, according to step S31-S33, respectively obtain 1 axle box oscillation crosswise signal, 1 framework oscillation crosswise signal and 4 The corresponding energy matrix E of upper framework oscillation crosswise signal1, energy matrix E2, energy matrix E3;S35, by energy matrix E1, energy matrix E2With energy matrix E3It joins together, obtains the corresponding energy matrix of vibration signal E。
- The rapid extracting method of operation characteristic 5. bullet train according to claim 4 slightly crawls, which is characterized in that institute State step S4 specifically:S41, time series T is resolved into K sections of time interval TK, and then generate corresponding K covariance matrixWherein, time series T is corresponding timing node when extracting original vibration signal;S42, the corresponding diagonalizable matrix of K covariance matrix is determinedWherein, diagonalizable matrixAre as follows:In formula: k=2 ..., K;For first covariance matrixDiagonalizable matrix, andV is first covariance MatrixEigenmatrix, first covariance matrix of ΛFeature vector, subscript H be conjugate transposition operator;S43, according to diagonalizable matrixDetermine that a unitary matrice U makes following formula value maximum;Wherein: function diag indicates diagonal element;P --- the dimension after dimensionality reduction;B, the b row and d of d --- matrix arrange;S44, given threshold ε, according to spin matrix to unitary matrice U and covariance matrixIt is iterated update, after updating Unitary matrice U and covariance matrixIn the value of all off diagonal elements when being respectively less than threshold epsilon, obtain final tenth of the twelve Earthly Branches square Battle arrayIt completes to K covariance matrixJoint diagonalization;Wherein, k=2 ... K;Spin matrix G (i, j, θ) are as follows:In formula: i and j is respectively the ith row and jth column of spin matrix;θ is intermediate computations parameter;S45, according to final unitary matriceDetermine transition matrix A;Wherein, transition matrix A are as follows:In formula: subscript # is pseudo-inverse transformation operator;S46, the fusion feature matrix Z according to transition matrix A, after being decomposed;Z=AE.
- The rapid extracting method of operation characteristic 6. bullet train according to claim 5 slightly crawls, which is characterized in that institute It states in step S44, to the iterative formula of unitary matrice U are as follows:U←UG(1,2,θ)In formula: G (1,2, θ) is the element value of spin matrix G (i, j, θ) the 1st row the 2nd column;To covariance matrixIterative formula are as follows:
- The rapid extracting method of operation characteristic 7. bullet train according to claim 1 slightly crawls, which is characterized in that institute State the bullet train in step S4 small size snake operation characteristic include normal operating condition feature, small size convergence state feature, Small size divergent state feature and substantially snakelike state feature.
Priority Applications (1)
Application Number | Priority Date | Filing Date | Title |
---|---|---|---|
CN201910336833.8A CN110084185B (en) | 2019-04-25 | 2019-04-25 | Method for rapidly extracting small-amplitude snaking operation characteristics of high-speed train |
Applications Claiming Priority (1)
Application Number | Priority Date | Filing Date | Title |
---|---|---|---|
CN201910336833.8A CN110084185B (en) | 2019-04-25 | 2019-04-25 | Method for rapidly extracting small-amplitude snaking operation characteristics of high-speed train |
Publications (2)
Publication Number | Publication Date |
---|---|
CN110084185A true CN110084185A (en) | 2019-08-02 |
CN110084185B CN110084185B (en) | 2021-03-16 |
Family
ID=67416623
Family Applications (1)
Application Number | Title | Priority Date | Filing Date |
---|---|---|---|
CN201910336833.8A Expired - Fee Related CN110084185B (en) | 2019-04-25 | 2019-04-25 | Method for rapidly extracting small-amplitude snaking operation characteristics of high-speed train |
Country Status (1)
Country | Link |
---|---|
CN (1) | CN110084185B (en) |
Cited By (4)
Publication number | Priority date | Publication date | Assignee | Title |
---|---|---|---|---|
CN112948981A (en) * | 2021-04-08 | 2021-06-11 | 西南交通大学 | Interval prediction method for small-amplitude snaking evolution trend of high-speed train |
CN113624521A (en) * | 2021-08-20 | 2021-11-09 | 厦门物之联智能科技有限公司 | Method and system for monitoring snake-shaped instability of train based on axle box vibration |
CN114997252A (en) * | 2022-08-05 | 2022-09-02 | 西南交通大学 | Vehicle-mounted detection method for wheel polygon based on inertia principle |
CN115017965A (en) * | 2022-08-09 | 2022-09-06 | 西南交通大学 | Snake classification method based on HHT energy and maximum Lyapunov exponent |
Citations (12)
Publication number | Priority date | Publication date | Assignee | Title |
---|---|---|---|---|
US20060182455A1 (en) * | 2005-02-17 | 2006-08-17 | Nobuhito Matsushiro | Image processing method, image processing apparatus, and image forming apparatus |
US20100057260A1 (en) * | 2008-08-29 | 2010-03-04 | Behrooz Fallahi | Self-tuning vibration absorber |
CN101900708A (en) * | 2010-08-18 | 2010-12-01 | 哈尔滨工业大学 | Vibration and audio signal-based high-speed train track defect detecting method |
CN102175768A (en) * | 2011-02-22 | 2011-09-07 | 哈尔滨工业大学 | Method and device for detecting defects and failures of high-speed rail based on vibration signals |
CN103163505A (en) * | 2013-01-31 | 2013-06-19 | 西安电子科技大学 | Time-varying narrow-band interference suppression method based on joint approximate diagonalization of eigen-matrices (JADE) |
CN104246850A (en) * | 2012-04-25 | 2014-12-24 | 丰田自动车株式会社 | Drift-assessment device |
CN104502126A (en) * | 2014-12-28 | 2015-04-08 | 华东交通大学 | Modal intervals-based high-speed train bogie fault diagnosis method |
CN105699080A (en) * | 2015-12-18 | 2016-06-22 | 华北电力大学(保定) | Wind turbine generator set bearing fault feature extraction method based on vibration data |
CN105890744A (en) * | 2016-03-30 | 2016-08-24 | 中车青岛四方机车车辆股份有限公司 | Method for determining snakelike instability of rail vehicle |
CN107192553A (en) * | 2017-06-28 | 2017-09-22 | 石家庄铁道大学 | Gear-box combined failure diagnostic method based on blind source separating |
CN107563403A (en) * | 2017-07-17 | 2018-01-09 | 西南交通大学 | A kind of recognition methods of bullet train operating condition |
CN108254179A (en) * | 2017-08-08 | 2018-07-06 | 常州路航轨道交通科技有限公司 | A kind of bullet train wheel set bearing method for diagnosing faults based on MEEMD arrangement entropys |
-
2019
- 2019-04-25 CN CN201910336833.8A patent/CN110084185B/en not_active Expired - Fee Related
Patent Citations (12)
Publication number | Priority date | Publication date | Assignee | Title |
---|---|---|---|---|
US20060182455A1 (en) * | 2005-02-17 | 2006-08-17 | Nobuhito Matsushiro | Image processing method, image processing apparatus, and image forming apparatus |
US20100057260A1 (en) * | 2008-08-29 | 2010-03-04 | Behrooz Fallahi | Self-tuning vibration absorber |
CN101900708A (en) * | 2010-08-18 | 2010-12-01 | 哈尔滨工业大学 | Vibration and audio signal-based high-speed train track defect detecting method |
CN102175768A (en) * | 2011-02-22 | 2011-09-07 | 哈尔滨工业大学 | Method and device for detecting defects and failures of high-speed rail based on vibration signals |
CN104246850A (en) * | 2012-04-25 | 2014-12-24 | 丰田自动车株式会社 | Drift-assessment device |
CN103163505A (en) * | 2013-01-31 | 2013-06-19 | 西安电子科技大学 | Time-varying narrow-band interference suppression method based on joint approximate diagonalization of eigen-matrices (JADE) |
CN104502126A (en) * | 2014-12-28 | 2015-04-08 | 华东交通大学 | Modal intervals-based high-speed train bogie fault diagnosis method |
CN105699080A (en) * | 2015-12-18 | 2016-06-22 | 华北电力大学(保定) | Wind turbine generator set bearing fault feature extraction method based on vibration data |
CN105890744A (en) * | 2016-03-30 | 2016-08-24 | 中车青岛四方机车车辆股份有限公司 | Method for determining snakelike instability of rail vehicle |
CN107192553A (en) * | 2017-06-28 | 2017-09-22 | 石家庄铁道大学 | Gear-box combined failure diagnostic method based on blind source separating |
CN107563403A (en) * | 2017-07-17 | 2018-01-09 | 西南交通大学 | A kind of recognition methods of bullet train operating condition |
CN108254179A (en) * | 2017-08-08 | 2018-07-06 | 常州路航轨道交通科技有限公司 | A kind of bullet train wheel set bearing method for diagnosing faults based on MEEMD arrangement entropys |
Non-Patent Citations (7)
Title |
---|
JARI MIETTINEN 等: "Blind Source Separation Based on Joint Diagonalization in R: The Packages JADE and BSSasymp", 《JOURNAL OF STATISTICAL SOFTWARE》 * |
JING NING 等: "A multi-sensor fusion framework for detecting small amplitude hunting of high-speed trains", 《JOURNAL OF VIBRATION AND CONTROL》 * |
YONGBIN LIU 等: "Feature fusion using kernel joint approximate diagonalization of eigen-matrices for rolling bearing fault identification", 《JOURNAL OF SOUND AND VIBRATION》 * |
冉伟 等: "基于EEMD-SVD-LTSA的高速列车蛇行演变特征提取框架", 《电子测量技术》 * |
叶运广 等: "高速列车转向架蛇行失稳的MEEMD-LSSVM预测模型", 《铁道学报》 * |
崔万里 等: "EMD-ISOMAP高速列车小幅蛇行异常特征提取", 《中国测试》 * |
徐智超 等: "JADE 盲源分离算法在变压器振动信号监测中的应用", 《传感器与微***》 * |
Cited By (7)
Publication number | Priority date | Publication date | Assignee | Title |
---|---|---|---|---|
CN112948981A (en) * | 2021-04-08 | 2021-06-11 | 西南交通大学 | Interval prediction method for small-amplitude snaking evolution trend of high-speed train |
CN112948981B (en) * | 2021-04-08 | 2022-06-21 | 西南交通大学 | Interval prediction method for small-amplitude snaking evolution trend of high-speed train |
CN113624521A (en) * | 2021-08-20 | 2021-11-09 | 厦门物之联智能科技有限公司 | Method and system for monitoring snake-shaped instability of train based on axle box vibration |
CN113624521B (en) * | 2021-08-20 | 2024-04-02 | 厦门物之联智能科技有限公司 | Train serpentine instability monitoring method and system based on axle box vibration |
CN114997252A (en) * | 2022-08-05 | 2022-09-02 | 西南交通大学 | Vehicle-mounted detection method for wheel polygon based on inertia principle |
CN115017965A (en) * | 2022-08-09 | 2022-09-06 | 西南交通大学 | Snake classification method based on HHT energy and maximum Lyapunov exponent |
CN115017965B (en) * | 2022-08-09 | 2022-11-01 | 西南交通大学 | Snake classification method based on HHT energy and maximum Lyapunov exponent |
Also Published As
Publication number | Publication date |
---|---|
CN110084185B (en) | 2021-03-16 |
Similar Documents
Publication | Publication Date | Title |
---|---|---|
CN110084185A (en) | A kind of bullet train slightly crawls the rapid extracting method of operation characteristic | |
Ye et al. | Fault diagnosis of high-speed train suspension systems using multiscale permutation entropy and linear local tangent space alignment | |
Wang et al. | An optimal vibration control strategy for a vehicle's active suspension based on improved cultural algorithm | |
Cao et al. | The fault diagnosis of a switch machine based on deep random forest fusion | |
Shuanbao et al. | Optimization design for aerodynamic elements of high speed trains | |
CN106250613A (en) | A kind of wheel service state security domain is estimated and method for diagnosing faults | |
CN100485342C (en) | Integrated supporting vector machine mixed intelligent diagnosing method for mechanical fault | |
CN106021789A (en) | Fuzzy-intelligence-based rail car suspension system fault classification method and system | |
CN106596149A (en) | Method for monitoring and diagnosing flat wheel state of rail vehicle based on VMD | |
Ito et al. | Predicting traffic congestion using driver behavior | |
Liu et al. | Modeling crossing behavior of drivers at unsignalized intersections with consideration of risk perception | |
CN104688252A (en) | Method for detecting fatigue status of driver through steering wheel rotation angle information | |
CN108664673A (en) | A kind of wheel rail force load identification feature-based data model method for building up and device | |
Shi et al. | Applying a WNN-HMM based driver model in human driver simulation: Method and test | |
Liu et al. | Vehicle state estimation based on unscented Kalman filtering and a genetic-particle swarm algorithm | |
Wang et al. | A novel optimal fuzzy integrated control method of active suspension system | |
Liang et al. | 1d convolutional neural networks for fault diagnosis of high-speed train bogie | |
Xi et al. | A hierarchical ensemble learning framework for energy-efficient automatic train driving | |
Zhang et al. | Stacking‐based ensemble learning method for the recognition of the preceding vehicle lane‐changing manoeuvre: A naturalistic driving study on the highway | |
Yu et al. | Driving risk of road vehicle shielded by bridge tower under strong crosswind | |
Gong et al. | Modal vibration decomposition method and its application on multi-mode vibration control of high-speed railway car bodies | |
CN116467570B (en) | Heavy-duty locomotive coupler swing angle quantitative identification method based on data driving | |
CN116946183A (en) | Commercial vehicle driving behavior prediction method considering driving capability and vehicle equipment | |
Qin et al. | Bogie fault identification based on EEMD information entropy and manifold learning | |
CN112308824B (en) | Curve radius classification identification method and device based on track geometric detection data |
Legal Events
Date | Code | Title | Description |
---|---|---|---|
PB01 | Publication | ||
PB01 | Publication | ||
SE01 | Entry into force of request for substantive examination | ||
SE01 | Entry into force of request for substantive examination | ||
GR01 | Patent grant | ||
GR01 | Patent grant | ||
CF01 | Termination of patent right due to non-payment of annual fee |
Granted publication date: 20210316 |
|
CF01 | Termination of patent right due to non-payment of annual fee |