METHOD FOR SEPARATING AND ESTIMATING MULTIPLE MOTION PARAMETERS IN XRAY ANGIOGRAM IMAGE

0Associated
Cases 
0Associated
Defendants 
0Accused
Products 
3Forward
Citations 
0
Petitions 
1
Assignment
First Claim
1. A method for separating and estimating multiple motion parameters in an Xray angiogram image, the method comprising:
 (1) tracing structure feature points of vessels in a Xray angiogram image sequence whereby obtaining tracing curves of said feature points s_{i}(n), i=1, . . . , I, n=1, . . . , N;
wherein I represents the number of the feature points, and N represents the number of image frames in said Xray angiogram image sequence;
(2) obtaining a simulated translation curve d_{α}(n) according to a variation frame sequence {N_{t}, t=1, . . . , T} and a slope angle sequence {α
_{t}, t=1, . . . , T−
1} of translational motion, where T represents the number of variations in motion directions;
(3) determining a cardiac motion signal cycle N_{c }according to said Xray angiogram image sequence, obtaining a multimotion synthetic motion curve ŝ
_{iα}^{1}(n)=s_{i}(n)−
d_{α}(n) without translational motion according to said tracing curve s_{i}(n) and said simulated translation curve d_{α}(n), and processing said synthetic motion curve ŝ
_{iα}^{1}(n) via Fourier frequencydomain filtering according to said cardiac motion signal cycle N_{c }thereby obtaining a cardiac motion curve ĉ
_{iα}(n);
(4) obtaining a residual motion curve ŝ
_{1α}^{2}(n)=ŝ
_{1α}^{1 }(n)−
ĉ
_{iα}(n) with no translational motion signal or cardiac motion signal according to said synthetic motion curve ŝ
_{iα}^{1}(n) and said cardiac motion curve ĉ
_{iα}(n), processing said residual motion curve ŝ
_{iα}(n) via Fourier frequencydomain filtering according to each respiratory motion signal cycle in a cycle range of [3N_{c}, 10N_{c}], thereby obtaining a corresponding respiratory motion curve, obtaining an optimum respiratory motion curve {circumflex over (r)}_{iα}(n) and an optimum respiratory motion signal cycle N_{αr }with respect to a current simulated translation curve using a fitting curve {circumflex over (r)}′
_{iα}(n) closest to said residual motion curve ŝ
_{iα}^{2}(n) as an optimum criteria;
(5) detecting whether an amplitude of a curve ĥ
′
_{iα}(n)=ŝ
_{iα}^{2}(n)−
{circumflex over (r)}′
_{iα}(n) obtained according to said residual motion curve ŝ
_{iα}^{2}(n) and said fitting curve {circumflex over (r)}′
_{iα}(n) is less than three pixels, determining there is no highfrequency component if yes, otherwise processing said residual motion curve ŝ
_{iα}^{2}(n) via Fourier frequencydomain filtering according to each highfrequency motion signal cycle in a cycle range of [1/7N_{c}, 5/7N_{c}], thereby obtaining a corresponding highfrequency motion curve, obtaining an optimum highfrequency motion curve ĥ
_{iα}(n) and an optimum highfrequency motion signal cycle {circumflex over (N)}_{αh }with respect to a current simulated translation curve using a fitting curve ĥ
′
_{iα}(n) closest to said highfrequency motion curve as an optimum criteria;
(6) obtaining a synthetic motion estimation curve ŝ
_{iα}(n)=d_{α}(n)+{circumflex over (r)}_{iα}(n)+ĉ
_{iα}(n)+ĥ
_{iα}(n) according to said simulated translation curve d_{α}(n), said respiratory motion curve {circumflex over (r)}_{iα}(n), said cardiac motion curve ĉ
_{iα}(n), and said highfrequency motion curve ĥ
_{iα}(n), and obtaining an optimum translational motion curve, a cardiac motion curve, a respiratory motion curve, and a highfrequency motion curve using a tracing curve s_{i}(n) closest to said synthetic motion estimation curve ŝ
_{iα}(n) as an optimum criteria.
1 Assignment
0 Petitions
Accused Products
Abstract
A method for separating and estimating multiple motion parameters in an Xray angiogram image. The method includes: determining a cardiac motion signal cycle and a variation frame sequence of translational motion according to an angiogram image sequence, tracing structure feature points of vessels in the angiogram image sequence whereby obtaining a motion sequence, processing the motion sequence via multivariable optimization and Fourier frequencydomain filtering, separating an optimum translational motion curve, a cardiac motion curve, a respiratory motion curve and a highfrequency motion curve according to the variation frame sequence of translational motion, a cycle of the cardiac motion signal, a range of a respiratory motion signal cycle, and a range of a highfrequency motion signal cycle.
9 Citations
View as Search Results
INFANT MONITORING SYSTEM  
Patent #
US 20180035082A1
Filed 03/27/2017

Current Assignee
Chigru Innovations OPC Private Limited

Sponsoring Entity
Chigru Innovations OPC Private Limited

Infant monitoring system  
Patent #
US 10,447,972 B2
Filed 03/27/2017

Current Assignee
Chigru Innovations OPC Private Limited

Sponsoring Entity
Chigru Innovations OPC Private Limited

Movementdependent reconstruction of magnetic resonance images  
Patent #
US 10,726,586 B2
Filed 07/06/2018

Current Assignee
Siemens Healthcare GMBH

Sponsoring Entity
Siemens Healthcare GMBH

METHOD OF COMPENSATION OF RESPIRATORY MOTION IN CARDIAC IMAGING  
Patent #
US 20100310140A1
Filed 04/27/2010

Current Assignee
Siemens Healthcare GMBH

Sponsoring Entity
Siemens AG

METHOD OF DETECTION AND COMPENSATION FOR RESPIRATORY MOTION IN RADIOGRAPHY CARDIAC IMAGES SYNCHRONIZED WITH AN ELECTROCARDIOGRAM SIGNAL  
Patent #
US 20080240536A1
Filed 03/25/2008

Current Assignee
General Electric Company

Sponsoring Entity
General Electric Company

Method for measuring myocardial motion and the like  
Patent #
US 6,171,241 B1
Filed 06/11/1998

Current Assignee
NcneilPpc Incorporated

Sponsoring Entity
NcneilPpc Incorporated

Method for computed tomography motion estimation and compensation  
Patent #
US 8,224,056 B2
Filed 12/15/2009

Current Assignee
General Electric Company

Sponsoring Entity
General Electronics Company

Automatic Registration of Image Series with Varying Contrast Based on Synthetic Images Derived from Intensity Behavior Model  
Patent #
US 20120189183A1
Filed 09/20/2011

Current Assignee
Siemens Healthcare GMBH

Sponsoring Entity
Siemens AG

PHASE SENSITIVE T1 MAPPING IN MAGNETIC RESONANCE IMAGING  
Patent #
US 20130272591A1
Filed 04/17/2013

Current Assignee
Siemens Healthcare GMBH

Sponsoring Entity
Siemens AG

1 Claim
 1. A method for separating and estimating multiple motion parameters in an Xray angiogram image, the method comprising:
(1) tracing structure feature points of vessels in a Xray angiogram image sequence whereby obtaining tracing curves of said feature points s_{i}(n), i=1, . . . , I, n=1, . . . , N;
wherein I represents the number of the feature points, and N represents the number of image frames in said Xray angiogram image sequence;(2) obtaining a simulated translation curve d_{α}(n) according to a variation frame sequence {N_{t}, t=1, . . . , T} and a slope angle sequence {α
_{t}, t=1, . . . , T−
1} of translational motion, where T represents the number of variations in motion directions;(3) determining a cardiac motion signal cycle N_{c }according to said Xray angiogram image sequence, obtaining a multimotion synthetic motion curve ŝ
_{iα}^{1}(n)=s_{i}(n)−
d_{α}(n) without translational motion according to said tracing curve s_{i}(n) and said simulated translation curve d_{α}(n), and processing said synthetic motion curve ŝ
_{iα}^{1}(n) via Fourier frequencydomain filtering according to said cardiac motion signal cycle N_{c }thereby obtaining a cardiac motion curve ĉ
_{iα}(n);(4) obtaining a residual motion curve ŝ
_{1α}^{2}(n)=ŝ
_{1α}^{1 }(n)−
ĉ
_{iα}(n) with no translational motion signal or cardiac motion signal according to said synthetic motion curve ŝ
_{iα}^{1}(n) and said cardiac motion curve ĉ
_{iα}(n), processing said residual motion curve ŝ
_{iα}(n) via Fourier frequencydomain filtering according to each respiratory motion signal cycle in a cycle range of [3N_{c}, 10N_{c}], thereby obtaining a corresponding respiratory motion curve, obtaining an optimum respiratory motion curve {circumflex over (r)}_{iα}(n) and an optimum respiratory motion signal cycle N_{αr }with respect to a current simulated translation curve using a fitting curve {circumflex over (r)}′
_{iα}(n) closest to said residual motion curve ŝ
_{iα}^{2}(n) as an optimum criteria;(5) detecting whether an amplitude of a curve ĥ
′
_{iα}(n)=ŝ
_{iα}^{2}(n)−
{circumflex over (r)}′
_{iα}(n) obtained according to said residual motion curve ŝ
_{iα}^{2}(n) and said fitting curve {circumflex over (r)}′
_{iα}(n) is less than three pixels, determining there is no highfrequency component if yes, otherwise processing said residual motion curve ŝ
_{iα}^{2}(n) via Fourier frequencydomain filtering according to each highfrequency motion signal cycle in a cycle range of [1/7N_{c}, 5/7N_{c}], thereby obtaining a corresponding highfrequency motion curve, obtaining an optimum highfrequency motion curve ĥ
_{iα}(n) and an optimum highfrequency motion signal cycle {circumflex over (N)}_{αh }with respect to a current simulated translation curve using a fitting curve ĥ
′
_{iα}(n) closest to said highfrequency motion curve as an optimum criteria;(6) obtaining a synthetic motion estimation curve ŝ
_{iα}(n)=d_{α}(n)+{circumflex over (r)}_{iα}(n)+ĉ
_{iα}(n)+ĥ
_{iα}(n) according to said simulated translation curve d_{α}(n), said respiratory motion curve {circumflex over (r)}_{iα}(n), said cardiac motion curve ĉ
_{iα}(n), and said highfrequency motion curve ĥ
_{iα}(n), and obtaining an optimum translational motion curve, a cardiac motion curve, a respiratory motion curve, and a highfrequency motion curve using a tracing curve s_{i}(n) closest to said synthetic motion estimation curve ŝ
_{iα}(n) as an optimum criteria.
1 Specification
This application is a continuationinpart of International Patent Application No. PCT/CN2014/085724 with an international filing date of Sep. 2, 2014, designating the United States, now pending, and further claims priority benefits to Chinese Patent Application No. 201310752751.4 filed Dec. 31, 2013. The contents of all of the aforementioned applications, including any intervening amendments thereto, are incorporated herein by reference. Inquiries from the public to applicants or assignees concerning this document or the related applications should be directed to: Matthias Scholl P. C., Attn.: Dr. Matthias Scholl Esq., 245 First Street, 18th Floor, Cambridge, Mass. 02142.
1. Field of the Invention
The invention relates to a method for separating and estimating multiple motion parameters in an Xray angiogram image.
2. Description of the Related Art
Respiratory motion includes inspiration motion and expiration motion formed by rhythmic inflation and deflation of a thorax. The respiratory motion may cause translational motion of the heart in a threedimensional space. In an Xray angiography system, twodimensional motion of coronary vessels may occur on angiography plane. Moreover, an angiography sequence may cover the entire coronary vessels during angiography, causing twodimensional translational motion between different frames of the imaging sequence.
Therefore, a coronary angiogram image operates to record a projection of cardiac motion on a twodimensional plane, and superposition of twodimensional translational motion of coronary on the angiography plane caused by respiratory motion of a human body, and twodimensional translational motion of the patient. To obtain an accurate twodimensional angiogram, it is beneficial to subtract cardiac motion, the respiratory motion and translational motion of the patient.
A conventional method for separating the translational motion is to conduct image registration between frames. Medical image registration methods reported in some literatures normally comprise an external registration method and an internal registration method. The external registration method sets some obvious reference points prior to angiography, and traces them during the angiography. However, this kind of labeling method is normally invasive. The internal registration method is divided into a labeling method, a segmentation method and so on after angiography is completed. The labeling method selects some anatomical structure points from an angiogram for registration, but not all angiograms have such anatomical structure points. The segmentation method can enable registration via anatomical structure lines that are segmented, and can be applicable for both rigid models and deformation models. However, it can only extract motion of the anatomical structure line, and cannot extract rigid translational motion alone. However, motion of vessel extracted from the angiogram images contains translational motion of a patient, along with physiological motion and respiratory motion of the vessel itself. Even if some angiogram images contain skeletons with no physiological motion (such as spines) as internal labeling points and it is possible to enable registration via the internal labeling method, not all the angiogram images have such a feature for registration. As no such feature exists in an angiogram image, it is difficult to extract the translational motion. Moreover, extraction of the translational motion is seldom reported in the medical image field.
It has been suggested to set labeling points before extracting human'"'"'s respiratory motion, and tracing those points in sequence. One feature of the respiratory motion is that as one person breaths, some organs in his body, such as a heart, a diaphragm and so on may move along with a lung in a threedimensional space. Therefore, a motion curve obtained by tracing structure feature points that are not to move along with the heart can be equivalent to respiratory motion. Two commonused methods in Xray angiography comprise manually tracing noncardiac structure feature points in the angiogram image, and simultaneously recording motion of these points during angiography. Obviously, both these two methods have defects. The former one features poor applicability: since it is impossible to ensure that a labeling point meeting the abovementioned requirement exists in each frame of the angiogram image, and to find the labeling point requires good understanding of human'"'"'s anatomical structure, it is hard to extract the respiratory motion as there is no feature point in the angiogram image. Implementation of the latter one requires large amount of experiment for controlling, which is inappropriate for normal clinical application. In addition, there is still another method that is implemented under a dualarm Xray angiogram, and a principle thereof for separating the cardiac motion from the respiratory motion is to conduct threedimensional reconstruction on coronary vessels of two angiogram images at different projection angles at the same time thereby obtaining threedimensional distribution of the vessel at that time. After all angiogram images in one respiratory cycle are matched and reconstructed, a group of threedimensional structure sequences can be obtained, and spatial displacement vectors therebetween constitute the respiratory motion. This method can obtain comparatively reliable respiratory motion, but cannot be widely applied to medical practice due to limitation of conditions of the dualarm Xray angiogram.
A typical method for extracting respiratory motion parameters from xray angiogram image discloses Fourier frequencydomain separation. However, a residual motion curve can be regarded as respiratory motion only after the cardiac motion is extracted. In addition, it is possible to extract translational motion and highfrequency motion, and as a result, respiratory motion is inaccurate, and robustness of the method is insufficient.
In view of the abovementioned problems, it is an objective of the invention to provide a method for separating and estimating multiple motion parameters in an Xray angiogram image that is based on frequencydomain filtering and multivariable optimization, and capable of accurately separating various motion parameters in a simple, fast and undamaged manner.
To achieve the above objective, in accordance with one embodiment of the invention, there is provided a method for separating and estimating multiple motion parameters in an Xray angiogram image, comprising steps of:
(1) tracing structure feature points of vessels in a Xray angiogram image sequence thereby obtaining tracing curves of the feature points s_{i}(n), i=1, . . . , I, n=1, . . . , N; wherein I represents the number of the feature points, and N represents the number of image frames in the Xray angiogram image sequence;
(2) obtaining a simulated translation curve d_{α}(n) according to a variation frame sequence {N_{t}, t=1, . . . , T} and a slope angle sequence {α_{t}, t=1, . . . , T−1} of translational motion, where T represents the number of variations in motion directions;
(3) determining a cardiac motion signal cycle N_{c }according to the Xray angiogram image sequence, obtaining a multimotion synthetic motion curve ŝ_{iα}^{1}(n)=s_{i}(n)−d_{α}(n) without translational motion according to the tracing curve s_{i}(n) and the simulated translation curve d_{α}(n), and processing the synthetic motion curve ŝ_{iα}^{1}(n) via Fourier frequencydomain filtering according to the cardiac motion signal cycle N_{c }thereby obtaining a cardiac motion curve ĉ_{iα}(n);
(4) obtaining a residual motion curve ŝ_{1α}^{2}(n)=ŝ_{1α}^{1 }(n)−ĉ_{iα}(n) with no translational motion signal or cardiac motion signal according to the synthetic motion curve ŝ_{1α}^{1}(n) and the cardiac motion curve ĉ_{iα}(n), processing the residual motion curve ŝ_{1α}^{2 }(n) via Fourier frequencydomain filtering according to each respiratory motion signal cycle in a cycle range of [3N_{c}, 10N_{c}], thereby obtaining a corresponding respiratory motion curve, obtaining an optimum respiratory motion curve {circumflex over (r)}_{iα}(n) and an optimum respiratory motion signal cycle {circumflex over (N)}_{αr}, with respect to a current simulated translation curve using a fitting curve {circumflex over (r)}′_{iα}(n) closest to the residual motion curve ŝ_{iα}^{2 }(n) as an optimum criteria;
(5) detecting whether an amplitude of a curve ĥ′_{iα}(n)=ŝ_{iα}^{2}(n)−{circumflex over (r)}′_{iα}(n) obtained according to the residual motion curve ŝ_{iα}^{2}(n) and the fitting curve {circumflex over (r)}′_{iα}(n) is less than three pixels, determining there is no highfrequency component if yes, otherwise processing the residual motion curve ŝ_{iα}^{2}(n) via Fourier frequencydomain filtering according to each highfrequency motion signal cycle in a cycle range of [1/7N_{c}, 5/7N_{c}], thereby obtaining a corresponding highfrequency motion curve, obtaining an optimum highfrequency motion curve ĥ_{iα}(n) and an optimum highfrequency motion signal cycle N_{αh }with respect to a current simulated translation curve using a fitting curve ĥ′_{iα}(n) closest to the highfrequency motion curve as an optimum criteria;
(6) obtaining a synthetic motion estimation curve ŝ_{iα}(n)=d_{α}(n)+{circumflex over (r)}_{iα}(n)+ĉ_{iα}(n)+ĥ_{iα}(n) according to the simulated translation curve d_{α}(n), the respiratory motion curve {circumflex over (r)}_{iα}(n), the cardiac motion curve ĉ_{iα}(n), and the highfrequency motion curve ĥ_{iα}(n), and obtaining an optimum translational motion curve, a cardiac motion curve, a respiratory motion curve, and a highfrequency motion curve using a tracing curve s_{i}(n) closest to the synthetic motion estimation curve ŝ_{iα}(n) as an optimum criteria.
To summarize, the invention provides separated estimation of multiple motion parameters and takes an Xray angiogram image that reflects different kinds of motion signals such as a translational motion signal, a cardiac motion signal, a respiratory motion signal, a highfrequency motion signal and so on into account, all motions that are separated are more accurate over the prior art, which provides accurate and effective help for medical diagnosis.
In
For clear understanding of the objectives, features and advantages of the invention, detailed description of the invention will be given below in conjunction with accompanying drawings and specific embodiments. It should be noted that the embodiments are only meant to explain the invention, and not to limit the scope of the invention.
Specifically, the invention relates to a data separation method that integrates selection of structure feature points of vessels and Fourier series expansion, and is capable of solving a problem of separating motion signals of human bodies, such as a translational motion signal, a cardiac motion signal, a respiratory motion signal and so on.
An objective of the invention is to provide a method for separating and estimating multiple motion parameters in an Xray angiogram image that is based on frequencydomain filtering and multivariable optimization and capable of accurately separating all motion parameters in a simple, fast and undamaged manner, and can be widely applied to singlearm angiogram and dualarm angiogram.
An Xray angiogram image records a projection of cardiac motion on a coronary, and superposition of twodimensional motion of coronary on a angiography plane caused by respiratory motion of a human body, and twodimensional translational motion of the patient.
The invention integrates the Fourier series expansion and multivariable optimization for separating a translational signal, a cardiac signal, and a respiratory motion signal, and comprises the following steps:
1. Tracing structure feature points of vessels in an Xray angiogram image sequence thereby obtaining a tracing curve s_{i}(n), i=1, . . . , I, n=1, . . . , N, where I represents the number of the feature points, and N represents the number of image frames in the Xray angiogram image sequence. It is required that the structure feature points of vessels can really and comprehensively reflect motion information of the vessels, and thus based on teaching of physiology and anatomy, feature points that are selected are bifurcation points of vessels. To obtain the structure feature points of vessels, the Xray angiogram image needs to be preprocessed, segmented and thinned thereby generating a vascular skeleton, and then the structure feature points of vessels are selected via the vascular skeleton. A tracing method of the invention obtains a coordinate of each feature point at each frame by establishing a relationship for each bifurcation point of different frames in a skeleton sequence of the Xray angiogram image.
2. Assuming a variation frame sequence of translational motion {N_{t}, t=1, . . . , T}, and a slope angle sequence {α_{t}, t=1, . . . . , T=1} of translational motion, where T represents the number of variations in motion directions; then a simulated translation curve d_{α}(n) can be obtained according to the following formula:
It is possible to determine an Xaxis variation frame sequence of translational motion {N_{t}, t_{x}=1, . . . , T_{s}} via the angiogram image, where T_{x }represents the number of variations in an Xaxis motion direction, and an Xaxis simulated translation curve d_{α}_{x}(n) at a range of (−π/2, π/2) of the slope angle sequence {α_{t}, t=1, . . . , T−1}. Also, it is possible to determine a Yaxis variation frame sequence of translational motion {N_{t}, t_{y}=1, . . . , T_{y}} via the angiogram image, where T_{y }represents the number of variations in a Yaxis motion direction, and a Yaxis simulated translation curve d_{α}_{y}(n) at a range of (−π/2, π/2) of the slope angle sequence {α_{t}, t_{y}=1, . . . . , T_{x}−1}. d_{α}(n)=(d_{α}_{x}(n), d_{α}_{y}(n)) represents translational motion of a patient.
3. Determining a cardiac motion signal cycle N_{c }according to the Xray angiogram image sequence, obtaining a multimotion synthetic motion curve ŝ_{iα}^{1}(n)=s_{i}(n)−d_{α}(n) without translational motion, and processing the synthetic motion curve ŝ_{iα}^{1}(n) via Fourier frequencydomain filtering according to the cardiac motion signal cycle N_{c }thereby obtaining a cardiac motion curve ĉ_{iα}(n);
4. Obtaining a residual component ŝ_{iα}^{2}(n) with no translational motion signal or cardiac motion signal, processing each respiratory motion signal ŝ_{iα}^{2}(n) via Fourier frequencydomain filtering according to each respiratory motion signal cycle in a cycle range of N_{αr}=m_{c}N_{c}, M_{c}ε[3, 10], obtaining an optimum respiratory motion signal curve {circumflex over (r)}_{iα}(n) and an optimum respiratory motion signal cycle {circumflex over (N)}_{αr}, with respect to a current simulated translation curve using a fitting curve {circumflex over (r)}_{iα}(n) closest to the residual motion curve ŝ_{iα}^{2}(n) as an optimum criteria;
5. Detecting whether an amplitude of a curve ĥ′_{iα}(n)=ŝ_{iα}^{2}(n)−{circumflex over (r)}′_{iα}(n) obtained according to the residual component ŝ_{iα}^{2}(n) and the fitting curve {circumflex over (r)}′_{iα}(n) is less than three pixels, determining there is no highfrequency component if yes, otherwise processing each highfrequency component motion cycle ŝ_{iα}(n) via Fourier frequencydomain filtering according to each highfrequency motion signal cycle in a cycle range of N_{αh}=M_{h}N_{c}, m_{h}ε[1/7, 5/7], thereby obtaining an optimum highfrequency motion signal curve ĥ_{iα}(n) and an optimum highfrequency motion signal cycle {circumflex over (N)}_{αh }with respect to a current translation curve using a fitting curve ĥ′_{iα}(n) closest to the highfrequency motion curve as an optimum criteria;
6. Obtaining a synthetic motion estimation curve ŝ_{iα}(n)=d_{α}(n)+{circumflex over (r)}_{iα}(n)+ĉ_{iα}(n)+(n)+ĥ_{iα}(n) iε[1, I], and obtaining an optimum translational motion curve, a cardiac motion curve, a respiratory motion curve, and a highfrequency motion curve using a tracing curve s_{i}(n) closest to the synthetic motion estimation curve ŝ_{iα}(n) as an optimum criteria.
The method of the invention takes into account that an appropriate feature point (such as an intersection point between ribs) might not be found in every angiogram image under a condition of Xray angiogram, and integrates selection of structure feature points of vessels, Fourier frequencydomain filtering, and multivariable optimization for extracting multiple motion parameters, and thus the invention features wider applicability and flexibility over a traditional method that uses manual tracing to obtain respiratory motion or translational motion, and can be applied to almost all angiogram image sequences. Meanwhile, in comparison with a conventional method that sets labeling points in the vicinity of a heart and then traces the points via an imaging method, which may cause damage to a patient'"'"'s health since labeling objects in a body thereof is invasive, and may bring inevitable trouble and errors in practice due to a fact that adding, imaging, and removal of the objects and extraction of respiratory motion are very complex, the method of the invention features higher safety and feasibility.
The invention will be described hereinafter in conjunction with accompanying drawings and examples.
The invention provides a method for directly extracting multiple motion parameters such as a translational motion signal, a respiratory motion signal, and a cardiac motion signal from an Xray angiogram image sequence, the methods selects some feature points on coronary vessels and traces the points thereby obtaining a motion curve thereof with time, then obtains signals, like a translational signal, a cardiac signal, a respiratory motion signal and so on, by separating the motion curve via Fourier series expansion, frequencydomain filtering, and multivariable optimization. The invention features good clinical applicability.
Normal cardiac motion measured in clinical practice is 60 to 100 times per minute, and a cycle thereof is 0.6 to 1.0 s. Periodic motion having a frequency higher than 140 times per second and a cycle thereof less than 0.42 s is diagnosed as pathologically abnormal motion. If it is able to detect a frequency of highfrequency motion is higher than that of normal cardiac motion by analyzing the coronary angiogram image, and an amplitude of a highfrequency signal is greater than 2 to 3 pixels, it can be reasonably deduced that the highfrequency signal is somewhat a useful signal (or a disease diagnosis signal). The reason for this is, some errors (normally not greater than three pixels) may occur during procession of the coronary angiogram image due to inaccuracy of segmentation or skeleton extraction. Therefore, as it is detected that the highfrequency motion exists and amplitude thereof is greater than three pixels, a new and useful disease diagnosis signal can be provided for a doctor. A frequency of the respiratory motion signal is far less than that of each of those two kinds of motion referred above, normally 3 to 6 s, and may be much higher in a quiet condition.
The invention separates motion components (comprising the translational signal, the cardiac motion signal, the respiratory motion signal, highfrequency motion signal and so on) in the angiogram image by integrating Fourier series expansion with frequencydomain filtering based on features thereof as follows.
Next discrete Fourier series transformation will be introduced hereinafter.
1) Basis of Discrete Fourier Series Transformation
According to principle of Fourier series expansion, continuous temporal periodic signals can be expressed by Fourier series. Correspondingly, a periodic sequence can also be expressed by discrete Fourier series. If a cycle of a periodic sequence x_{p}(n) is N_{p}, then
x_{p}(n)=x_{p}(n+rN_{p}) (1)
where r can be any integer,
Then Fourier series transformation is conducted on x_{p}(n) according to the following formulae (2) and (3):
In formula (3),
represents a fundamental frequency component of the periodic sequence,
represents a k^{th }harmonic component, and X_{p}(k) represents a harmonic coefficient. Since
the number of independent harmonic components is N_{p}, the number of terms in sum of the Fourier series is k=0, . . . , N_{p}−1. Formula (2) is a summation formula that determines the coefficient X_{p}(k) by x_{p}(n).
2) Separation Algorithm
Assuming in an angiogram image sequence, a motion curve of an i^{th }feature point P_{i}(x, y) on a coronary vessel along an Xaxis is x_{i}(n) (where n represents the number of angiography frames that are discrete in time, and N represents a length of the sequence), a motion curve along a Yaxis is y_{i}(n), and s_{i}(n)=(x_{i}(n), y_{i}(n)), then s_{i}(n) can be expressed as follows:
s_{i}(n)=d_{i}(n)+r_{i}(n)+c_{i}(n)+h_{i}(n)+t_{i}(n)iε[1,I] (4)
where I represents the number of feature points, d(n)=(x_{d}(n), y_{d}(n)) represents rigid translational motion of a patient that varies with time, r_{i}(n)=(x_{ri}(n), y_{ri}(n)) represents motion caused by respiratory motion, c_{i}(n)=(x_{ci}(n), y_{ci}(n)) represents motion of each point on the vessel caused by cardiac motion, h_{i}(n)=(x_{hi}(n), y_{hi}(n)) represents motion caused by possible highfrequency motion of the patient, and t_{i}(n)=(x_{ti}(n), y_{ti}(n)) represents other noise components caused by vessel segmentation or thinning. For the purpose of easy explanation, s_{i}(n) represents a coordinate variation curve of extracted feature points on a vessel along both the Xaxis and the Yaxis, d(n), r_{i}(n), c_{i}(n), h_{i}(n) and t_{i}(n) respectively represents a coordinate variation curve of a translational motion signal, a respiratory motion signal, a cardiac motion signal, a highfrequency motion signal, and a noise component along the Xaxis and the Yaxis in the content hereinafter.
Assuming a cardiac motion signal cycle is N_{c}, a respiratory motion signal cycle is N_{r}, a highfrequency motion signal cycle is N_{h}, then according to formula (3), c(n) r(n) and h(n) are:
where
represents a harmonic component of the cardiac motion signal,
represents a harmonic component of the respiratory motion signal,
represents a harmonic component of the highfrequency motion signal, C(k), R(k) and H(k) respectively represents a harmonic coefficient of c(n), r(n) and h(n). To completely separate c(n), r(n), and h(n), harmonic components thereof cannot be overlapped to each other (assuming the translational motion d(n) is already eliminated), namely:
where n_{1}, n_{2 }and n_{3 }are timedomain independent variables, and k_{1}, k_{2 }and k_{3 }are frequencydomain independent variables.
To satisfy formula (6), it is required that N_{c}, N_{r }and N_{h }are coprime. In addition, to separate c(n), r(n) and h(n), a length of s(n) should be greater than or equal to a cycle of s(n), namely being greater than N_{c}·N_{r}·N_{h}.
However, the above two requirements are very difficult to be met: even if N_{c}, N_{r }and N_{h }are coprime, it is almost impossible to meet the second requirement. The reason for this is that time contrast agent staying in a human body is not long enough (not greater than 10 seconds, and normally 6 seconds for an angiogram image sequence), assuming a frame rate is 80 ms/frame (or even less), then N_{c}≈10, N_{r}>30, N_{c}·N_{r}·N_{h}>N_{c}·N_{r}>300>10/0.08, and it is difficult to obtain a complete cyclic sequence. Moreover, in terms of a translational motion curve, a requirement for determining the respiratory motion signal cycle is comparatively high.
Although it is impossible to perfectly separate the cardiac motion curve, the respiratory motion curve, the highfrequency motion curve and the translational motion curve, it is feasible to approximately separate them via Fourier series expansion, frequencydomain filtering, and multivariable optimization.
The separation algorithm is based on the following assumption:
(1) assuming a cycle is N_{a}, as a length of a sequence is an integer multiple of the cycle and is not an integer multiple of other motion signal cycles, a component of motion only exists at the frequency
which is proven by formula (3), and other components of motion are smooth at the frequency
This assumption is feasible since the respiratory motion signal, the cardiac motion signal and the highfrequency motion signal are smooth, and it is almost impossible for any one of them to become a fundamental frequency component of the other.
(2) since the translational motion signal is nonperiodic rigid motion varies with time, it is estimated and assumed that the translational motion signal is in polyline motion expressed as follows:
where {N_{t}, t=1, . . . , T} represents a variation frame of the translational motion signal, {α_{t}=1, . . . , T−1} represents a slope angle of the translational motion signal, T represents the number of variations in motion directions.
Step 1: tracing structure feature points of vessels in an Xray angiogram image sequence thereby obtaining a tracing curve s_{i}(n), i=1, . . . , I, n=1, . . . , N, where I represents the number of the feature points, and N represents the number of image frames in the Xray angiogram image sequence. It is required that the structure feature points of vessels can really and comprehensively reflect motion information of the vessels, and thus based on teaching of physiology and anatomy, feature points that are selected are bifurcation points of vessels. As shown in
Ideally, s_{i}(n) can be expressed as:
s_{i}(n)≈d(n)+r_{i}(n)+c_{i}(n)+h_{i}(n)iε[1,I] (8)
Step 2: Estimation of a Translational Motion Curve:
If translational motion of a patient exists, it is possible to determine an Xaxis variation frame sequence of translational motion {N_{t}_{x}, t_{x}=1, . . . , T_{x}} via the angiogram image, where T_{x }represents the number of variations in an Xaxis motion direction, and a Yaxis variation frame sequence of translational motion {N_{t}_{y}, t_{y}=1, . . . , T_{y}} via the angiogram image, where T_{y }represents the number of variations in a Yaxis motion direction. Moreover, it is possible to determine an Xaxis simulated translation curve d_{α}_{x}(n) at a range of (−π/2, π/2) of the slope angle sequence, and a Yaxis simulated translation curve d_{α}_{y}(n) at a range of (−π/2, π/2) of the slope angle sequence {α_{t}, t_{y}=1, . . . , T−1}. d_{α}(n)=(d_{α}_{x}(n), d_{α}_{y}(n)) represents translational motion of a patient.
A multimotion synthetic motion curve without translational motion is represented by ŝ_{iα}^{1}(n)=s_{i}(n)−d_{α}(n).
Step 3: Estimation of the Cardiac Motion Curve:
A cardiac motion signal cycle N_{c }is determined according to the angiogram image sequence. Normally N_{c}≈0.8/f, where N_{c }represents the number of angiogram image frames in a cardiac cycle, 0.8 represents time of the cardiac cycle (in second), and f represents a frame rate, namely the number of angiogram images that are taken per second. By choosing a sequence ŝ_{iα}^{11 }(n) having a length of n_{sc }from ŝ_{iα}^{1}(n), where n_{sc }is an integer multiple of N_{c}, ŝ_{iα}^{11}(n) is divided into a motion curve along the Xaxis x_{ic}(n), and a motion curve along the Yaxis y_{ic}(n), then discrete Fourier expansion is conducted thereon thereby obtaining harmonic coefficients X_{ic}(k) and Y_{ic}(k), k=0, . . . , n_{sc}−1;
A frequencydomain response of a residual motion curve without the translational motion signal and the cardiac motion signal is:
A frequencydomain response of the cardiac motion is:
C_{iα}^{x}(k)=X_{i}(k)−Ŝ_{iα}^{x}(k)(k=1, . . . ,n_{sc}−1)
Discrete Fourier inverse transformation is conducted on C_{iα}^{x}(k) thereby obtaining an Xaxis cardiac motion curve ĉ_{iα}^{x}(n), then ĉ_{iα}^{y}(n), and finally a cardiac motion curve ĉ_{iα}(n). Then discrete Fourier inverse transformation is conducted on the residual motion curve without the translational motion signal and the cardiac motion signal, thereby obtaining ŝ_{iα}^{2}(n).
Step 4: Estimation of a Respiratory Motion Curve
For the respiratory motion signal cycle N_{αr}=m_{c}N_{c}, m_{c}ε[3, 10], as a length of an original motion signal is greater than the respiratory motion signal cycle N_{αr}, Fourier frequencydomain filtering is conducted on ŝ_{iα}^{2}, (n) via the method of separating the cardiac motion signal thereby obtaining a respiratory motion curve r_{iαN}_{αr}(n); as the respiratory motion signal cycle is greater than length of the original motion signal, it is determined there it does not contain a respiratory motion signal cycle, and the following steps are taken to separate the respiratory motion signal:
a. selecting a sequence ŝ_{iN}_{αr}_{j}^{21}(n), j=1, . . . , n_{sc}−N_{αr}/2+1 with a length of N_{αr}/2 from the j^{th }point of the sequence ŝ_{iα}^{2}(n), conducting Fourier series expansion thereon thereby obtaining harmonic coefficients Ŝ_{iN}_{αr}_{j}^{21}(k), k=0, . . . , N_{αr}−1.
b. processing a curve having a length half than that of the respiratory periodic signal (the curve comprising a motion curve with a cycle of N_{αr}/2, and only has components at a frequency of
via Fourier frequencydomain filtering thereby obtaining a motion curve r_{iαN}_{αr}_{j}(n), j=1, . . . , n_{sc}−N_{αr}/2+1;
c. averaging the r_{iαN}_{αr}_{j}(n), j=1, . . . , n_{sc}−N_{αr}/2+1 thereby obtaining an approximate respiratory motion curve r_{iαN}_{αr}(n).
Then ŝ_{iα}^{2}(n) is processed via curve fitting thereby obtaining an estimation curve {circumflex over (r)}′_{iα}(n) of the respiratory motion signal (since a cycle of the respiratory motion signal is the largest and far greater than that of the highfrequency motion signal, the curve is an approximate shape of the respiratory motion signal). The following formula (10) is used to obtain an optimum respiratory motion signal cycle {circumflex over (N)}_{αr }and an optimum respiratory motion curve {circumflex over (r)}_{iα}(n) for an optimum criteria.
Step 5: Estimation of a HighFrequency Motion Curve
For a curve ĥ′_{iα}(n) obtained by substracting the curve {circumflex over (r)}′_{iα}(n) with the residual component ŝ_{iα}^{2 }(n), namely ĥ′_{iα}(n)=ŝ_{iα}^{2}(n)−{circumflex over (r)}′iα(n), if an amplitude of ĥ′_{iα}(n) is not greater than 3 pixels, it is determined there is no highfrequency component (this is because that some errors that normally are not greater than three pixels may occur during procession of the coronary angiogram image due to inaccuracy of segmentation or skeleton extraction), otherwise ŝ_{iα}^{2}(n) is processed via Fourier frequencydomain filtering (which is identical to the cardiac motion signal) according to each highfrequency motion signal cycle in a cycle range of N_{αh}=m_{h}N_{c}, m_{h}ε[1/7, 5/7], thereby obtaining a corresponding highfrequency motion curve h_{iαN}_{αr}(n).
The following formula (11) is used to obtain an optimum highfrequency motion signal cycle {circumflex over (N)}_{αh }and an optimum highfrequency motion curve ĥ_{iα}(n) according to the optimum criteria:
Step 6: Obtaining a Synthetic Motion Estimation Curve:
ŝ_{iα}(n)=d_{α}(n)+{circumflex over (r)}_{iα}(n)+ĉ_{iα}(n)+ĥ_{iα}(n) iε[1, I], and using the following formula (12) to obtain an optimum translational motion curve {circumflex over (d)}(n), a cardiac motion curve ĉ_{i}(n), a respiratory motion curve {circumflex over (r)}_{i}(n), and a highfrequency motion curve ĥ_{i}(n) according to the optimum criteria:
It is determined that the cardiac motion signal cycle is 11 frames after observing the angiogram image sequence, and
It can be seen from
After separation, the cardiac motion curve features good periodicity. As respiratory motion curve for all labeling points are compared with each other, it can be found that there is small difference between respiratory motion curves for feature points in the same angiography plane (as shown in
While particular embodiments of the invention have been shown and described, it will be obvious to those skilled in the art that changes and modifications may be made without departing from the invention in its broader aspects, and therefore, the aim in the appended claims is to cover all such changes and modifications as fall within the true spirit and scope of the invention.