WO2014080833A1 - 超音波診断装置、画像処理装置及び画像処理方法 - Google Patents
超音波診断装置、画像処理装置及び画像処理方法 Download PDFInfo
- Publication number
- WO2014080833A1 WO2014080833A1 PCT/JP2013/080813 JP2013080813W WO2014080833A1 WO 2014080833 A1 WO2014080833 A1 WO 2014080833A1 JP 2013080813 W JP2013080813 W JP 2013080813W WO 2014080833 A1 WO2014080833 A1 WO 2014080833A1
- Authority
- WO
- WIPO (PCT)
- Prior art keywords
- position information
- image data
- region
- unit
- interest
- Prior art date
Links
Images
Classifications
-
- A—HUMAN NECESSITIES
- A61—MEDICAL OR VETERINARY SCIENCE; HYGIENE
- A61B—DIAGNOSIS; SURGERY; IDENTIFICATION
- A61B8/00—Diagnosis using ultrasonic, sonic or infrasonic waves
- A61B8/08—Detecting organic movements or changes, e.g. tumours, cysts, swellings
- A61B8/0883—Detecting organic movements or changes, e.g. tumours, cysts, swellings for diagnosis of the heart
-
- A—HUMAN NECESSITIES
- A61—MEDICAL OR VETERINARY SCIENCE; HYGIENE
- A61B—DIAGNOSIS; SURGERY; IDENTIFICATION
- A61B8/00—Diagnosis using ultrasonic, sonic or infrasonic waves
- A61B8/52—Devices using data or image processing specially adapted for diagnosis using ultrasonic, sonic or infrasonic waves
- A61B8/5207—Devices using data or image processing specially adapted for diagnosis using ultrasonic, sonic or infrasonic waves involving processing of raw data to produce diagnostic data, e.g. for generating an image
-
- G—PHYSICS
- G06—COMPUTING; CALCULATING OR COUNTING
- G06T—IMAGE DATA PROCESSING OR GENERATION, IN GENERAL
- G06T7/00—Image analysis
- G06T7/20—Analysis of motion
- G06T7/246—Analysis of motion using feature-based methods, e.g. the tracking of corners or segments
-
- G—PHYSICS
- G06—COMPUTING; CALCULATING OR COUNTING
- G06T—IMAGE DATA PROCESSING OR GENERATION, IN GENERAL
- G06T2207/00—Indexing scheme for image analysis or image enhancement
- G06T2207/10—Image acquisition modality
- G06T2207/10132—Ultrasound image
-
- G—PHYSICS
- G06—COMPUTING; CALCULATING OR COUNTING
- G06T—IMAGE DATA PROCESSING OR GENERATION, IN GENERAL
- G06T2207/00—Indexing scheme for image analysis or image enhancement
- G06T2207/30—Subject of image; Context of image processing
- G06T2207/30004—Biomedical image processing
- G06T2207/30048—Heart; Cardiac
Definitions
- Embodiments described herein relate generally to an ultrasonic diagnostic apparatus, an image processing apparatus, and an image processing method.
- the endocardial contour of the myocardium is given as the initial contour position at the end diastole (first R wave time phase) or the temporary phase at the end systole.
- the initial contour position is automatically tracked by using local pattern matching processing and motion information (movement vector) obtained by the optical flow method in the remaining time phases. Get the contour position in the phase. Therefore, it is indispensable to accurately estimate the movement vector in order to realize accurate contour position tracking.
- the problem to be solved by the present invention is to provide an ultrasonic diagnostic apparatus, an image processing apparatus, and an image processing method capable of accurately obtaining a tracking result of a region of interest.
- the ultrasonic diagnostic apparatus includes a setting unit, an estimation unit, and a tracking unit.
- the setting unit performs two-dimensional or at least one time-phase ultrasonic image data among the moving image data of the two-dimensional or three-dimensional ultrasonic image data collected by imaging an area including a moving tissue.
- a three-dimensional region of interest is set.
- the estimation unit is a first position of a region estimated as the region of interest based on motion information in each ultrasonic image data in a remaining time phase other than the at least one time phase in the section in which the moving image data is collected Information and second position information of the region estimated as the region of interest based on information other than the motion information are acquired.
- the tracking unit acquires, as position information of the region of interest, position information obtained by combining the first position information and the second position information based on a reliability index related to the reliability of the motion information. To track.
- FIG. 1 is a block diagram illustrating a configuration example of the ultrasonic diagnostic apparatus according to the first embodiment.
- FIG. 2 is a diagram illustrating an example of moving image data acquired by the acquisition unit.
- FIG. 3 is a diagram illustrating an example of an initial contour set by the operator.
- FIG. 4 is a block diagram illustrating a configuration example of the tracking unit illustrated in FIG.
- FIG. 5A is a diagram (1) for explaining the tracking unit according to the first embodiment.
- FIG. 5B is a diagram (2) for explaining the tracking unit according to the first embodiment.
- FIG. 6 is a diagram (3) for explaining the tracking unit according to the first embodiment.
- FIG. 7 is a diagram (4) for explaining the tracking unit according to the first embodiment.
- FIG. 8 is a diagram (5) for explaining the tracking unit according to the first embodiment.
- FIG. 5A is a diagram (1) for explaining the tracking unit according to the first embodiment.
- FIG. 5B is a diagram (2) for explaining the tracking unit according to the first embodiment.
- FIG. 6
- FIG. 9 is a diagram (6) for explaining the tracking unit according to the first embodiment.
- FIG. 10 is a diagram (7) for explaining the tracking unit according to the first embodiment.
- FIG. 11 is a diagram (8) for explaining the tracking unit according to the first embodiment.
- FIG. 12 is a diagram (9) for explaining the tracking unit according to the first embodiment.
- FIG. 13A is a diagram (10) for explaining the tracking unit according to the first embodiment.
- FIG. 13B is a diagram (11) for explaining the tracking unit according to the first embodiment.
- FIG. 13C is a diagram (12) illustrating the tracking unit according to the first embodiment.
- FIG. 14 is a diagram (13) for explaining the tracking unit according to the first embodiment.
- FIG. 15 is a diagram (14) for explaining the tracking unit according to the first embodiment.
- FIG. 16 is a diagram (15) for explaining the tracking unit according to the first embodiment.
- FIG. 17 is a flowchart for explaining an overview of processing performed by the ultrasonic diagnostic apparatus according to the first embodiment.
- FIG. 18 is a flowchart for explaining an example of processing performed by the tracking unit according to the first embodiment.
- FIG. 19 is a flowchart for explaining an example of processing performed by the tracking unit according to the second embodiment.
- FIG. 20 is a diagram for explaining a modification according to the second embodiment.
- FIG. 1 is a block diagram illustrating a configuration example of the ultrasonic diagnostic apparatus according to the first embodiment.
- the ultrasound diagnostic apparatus according to the first embodiment includes an ultrasound probe 1, a monitor 2, an input device 3, an electrocardiograph 4, and a device body 10.
- the ultrasonic probe 1 includes a plurality of piezoelectric vibrators, and the plurality of piezoelectric vibrators generate ultrasonic waves based on a drive signal supplied from a transmission / reception unit 11 included in the apparatus main body 10 described later.
- the ultrasonic probe 1 receives a reflected wave from the subject P and converts it into an electrical signal.
- the ultrasonic probe 1 includes a matching layer provided in the piezoelectric vibrator, a backing material that prevents propagation of ultrasonic waves from the piezoelectric vibrator to the rear, and the like.
- the ultrasonic probe 1 is detachably connected to the apparatus main body 10.
- the transmitted ultrasonic waves are transmitted from the ultrasonic probe 1 to the subject P
- the transmitted ultrasonic waves are reflected one after another at the discontinuous surface of the acoustic impedance in the body tissue of the subject P
- the ultrasonic probe is used as a reflected wave signal. 1 is received by a plurality of piezoelectric vibrators.
- the amplitude of the received reflected wave signal depends on the difference in acoustic impedance at the discontinuous surface where the ultrasonic wave is reflected.
- the reflected wave signal when the transmitted ultrasonic pulse is reflected by the moving blood flow or the surface of the heart wall depends on the velocity component of the moving object in the ultrasonic transmission direction due to the Doppler effect. And undergoes a frequency shift.
- a 1D array probe in which a plurality of piezoelectric vibrators are arranged in a row is connected to the apparatus main body 10 as an ultrasonic probe 1 for two-dimensional scanning of the subject P.
- the 1D array probe as the ultrasonic probe 1 is a sector probe that performs sector scanning, a convex probe that performs offset sector scanning, a linear probe that performs linear scanning, or the like.
- a mechanical 4D probe or a 2D array probe may be connected to the apparatus main body 10 as the ultrasonic probe 1 for three-dimensional scanning of the subject P.
- the mechanical 4D probe is capable of two-dimensional scanning using a plurality of piezoelectric vibrators arranged in a line like a 1D array probe, and swings the plurality of piezoelectric vibrators at a predetermined angle (swing angle). By doing so, three-dimensional scanning is possible.
- the 2D array probe can be three-dimensionally scanned by a plurality of piezoelectric vibrators arranged in a matrix and can be two-dimensionally scanned by focusing and transmitting ultrasonic waves.
- the input device 3 includes a mouse, a keyboard, a button, a panel switch, a touch command screen, a foot switch, a trackball, a joystick, etc., receives various setting requests from an operator of the ultrasonic diagnostic apparatus, The various setting requests received are transferred.
- the setting information received from the operator by the input device 3 according to the first embodiment will be described in detail later.
- the monitor 2 displays a GUI (Graphical User Interface) for the operator of the ultrasonic diagnostic apparatus to input various setting requests using the input device 3, and displays ultrasonic image data generated in the apparatus main body 10. Or display.
- GUI Graphic User Interface
- the electrocardiograph 4 acquires an electrocardiogram (ECG) of the subject P as a biological signal of the subject P to be ultrasonically scanned.
- ECG electrocardiogram
- the electrocardiograph 4 transmits the acquired electrocardiogram waveform to the apparatus main body 10.
- the apparatus main body 10 is an apparatus that generates ultrasonic image data based on a reflected wave signal received by the ultrasonic probe 1.
- the apparatus main body 10 shown in FIG. 1 is an apparatus that can generate two-dimensional ultrasonic image data based on two-dimensional reflected wave data received by the ultrasonic probe 1.
- the apparatus main body 10 shown in FIG. 1 is an apparatus that can generate three-dimensional ultrasonic image data based on three-dimensional reflected wave data received by the ultrasonic probe 1.
- the three-dimensional ultrasound image data may be referred to as “volume data”.
- the apparatus body 10 includes a transmission / reception unit 11, a B-mode processing unit 12, a Doppler processing unit 13, an image generation unit 14, an image memory 15, an internal storage unit 16, and an image processing unit. 17 and a control unit 18.
- the transmission / reception unit 11 includes a pulse generator, a transmission delay unit, a pulser, and the like, and supplies a drive signal to the ultrasonic probe 1.
- the pulse generator repeatedly generates rate pulses for forming transmission ultrasonic waves at a predetermined rate frequency.
- the transmission delay unit generates a delay time for each piezoelectric vibrator necessary for focusing the ultrasonic wave generated from the ultrasonic probe 1 into a beam and determining transmission directivity. Give for each rate pulse.
- the pulser applies a drive signal (drive pulse) to the ultrasonic probe 1 at a timing based on the rate pulse. That is, the transmission delay unit arbitrarily adjusts the transmission direction of the ultrasonic wave transmitted from the piezoelectric vibrator surface by changing the delay time given to each rate pulse.
- the transmission / reception unit 11 has a function capable of instantaneously changing a transmission frequency, a transmission drive voltage, and the like in order to execute a predetermined scan sequence based on an instruction from the control unit 18 described later.
- the change of the transmission drive voltage is realized by a linear amplifier type transmission circuit capable of instantaneously switching the value or a mechanism for electrically switching a plurality of power supply units.
- the transmission / reception unit 11 includes a preamplifier, an A / D (Analog / Digital) converter, a reception delay unit, an adder, and the like.
- the transmission / reception unit 11 performs various processing on the reflected wave signal received by the ultrasonic probe 1 and reflects it. Generate wave data.
- the preamplifier amplifies the reflected wave signal for each channel.
- the A / D converter A / D converts the amplified reflected wave signal.
- the reception delay unit gives a delay time necessary for determining the reception directivity.
- the adder performs an addition process on the reflected wave signal processed by the reception delay unit to generate reflected wave data. By the addition processing of the adder, the reflection component from the direction corresponding to the reception directivity of the reflected wave signal is emphasized, and a comprehensive beam for ultrasonic transmission / reception is formed by the reception directivity and the transmission directivity.
- the transmission / reception unit 11 transmits a two-dimensional ultrasonic beam from the ultrasonic probe 1 when the subject P is two-dimensionally scanned. Then, the transmission / reception unit 11 generates two-dimensional reflected wave data from the two-dimensional reflected wave signal received by the ultrasonic probe 1. In addition, when the subject P is three-dimensionally scanned, the transmission / reception unit 11 transmits a three-dimensional ultrasonic beam from the ultrasonic probe 1. Then, the transmission / reception unit 11 generates three-dimensional reflected wave data from the three-dimensional reflected wave signal received by the ultrasonic probe 1.
- the form of the output signal from the transmitter / receiver 11 can be selected from various forms such as a signal including phase information called an RF (Radio Frequency) signal or amplitude information after envelope detection processing. Is possible.
- RF Radio Frequency
- the B-mode processing unit 12 receives the reflected wave data from the transmission / reception unit 11, performs logarithmic amplification, envelope detection processing, and the like, and generates data (B-mode data) in which the signal intensity is expressed by brightness. .
- the Doppler processing unit 13 performs frequency analysis on velocity information from the reflected wave data received from the transmission / reception unit 11, extracts blood flow, tissue, and contrast agent echo components due to the Doppler effect, and obtains moving body information such as velocity, dispersion, and power. Data extracted for multiple points (Doppler data) is generated.
- the B-mode processing unit 12 and the Doppler processing unit 13 can process both two-dimensional reflected wave data and three-dimensional reflected wave data. That is, the B-mode processing unit 12 generates two-dimensional B-mode data from the two-dimensional reflected wave data, and generates three-dimensional B-mode data from the three-dimensional reflected wave data.
- the Doppler processing unit 13 generates two-dimensional Doppler data from the two-dimensional reflected wave data, and generates three-dimensional Doppler data from the three-dimensional reflected wave data.
- the image generation unit 14 generates ultrasonic image data from the data generated by the B mode processing unit 12 and the Doppler processing unit 13. That is, the image generation unit 14 generates two-dimensional B-mode image data in which the intensity of the reflected wave is represented by luminance from the two-dimensional B-mode data generated by the B-mode processing unit 12. Further, the image generation unit 14 generates two-dimensional Doppler image data representing moving body information from the two-dimensional Doppler data generated by the Doppler processing unit 13. The two-dimensional Doppler image data is velocity image data, distributed image data, power image data, or image data obtained by combining these. Further, the image generation unit 14 can also generate a Doppler waveform in which blood flow and tissue velocity information is plotted in time series from the Doppler data generated by the Doppler processing unit 13.
- the image generation unit 14 generally converts (scan converts) a scanning line signal sequence of ultrasonic scanning into a scanning line signal sequence of a video format represented by a television or the like, and displays ultrasonic waves for display. Generate image data. Specifically, the image generation unit 14 generates ultrasonic image data for display by performing coordinate conversion in accordance with the ultrasonic scanning mode of the ultrasonic probe 1. In addition to the scan conversion, the image generation unit 14 performs various image processing, such as image processing (smoothing processing) for regenerating an average luminance image using a plurality of image frames after scan conversion, Image processing (edge enhancement processing) using a differential filter is performed in the image. In addition, the image generation unit 14 synthesizes character information, scales, body marks, and the like of various parameters with the ultrasound image data.
- image processing smoothing processing
- Image processing edge enhancement processing
- the B-mode data and the Doppler data are ultrasonic image data before the scan conversion process, and the data generated by the image generation unit 14 is ultrasonic image data for display after the scan conversion process.
- the B-mode data and Doppler data are also called raw data.
- the image generation unit 14 generates three-dimensional B-mode image data by performing coordinate conversion on the three-dimensional B-mode data generated by the B-mode processing unit 12. Further, the image generation unit 14 generates three-dimensional Doppler image data by performing coordinate conversion on the three-dimensional Doppler data generated by the Doppler processing unit 13. In other words, the image generation unit 14 generates “three-dimensional B-mode image data or three-dimensional Doppler image data” as “three-dimensional ultrasound image data (volume data)”.
- the image generation unit 14 performs a rendering process on the volume data in order to generate various two-dimensional image data for displaying the volume data on the monitor 2.
- the rendering processing performed by the image generation unit 14 includes processing for generating MPR image data from volume data by performing a cross-section reconstruction method (MPR: Multi Planer Reconstruction).
- MPR Multi Planer Reconstruction
- the rendering process performed by the image generation unit 14 includes a process of performing “Curved MPR” on volume data and a process of performing “Maximum Intensity Projection” on volume data.
- the rendering process performed by the image generation unit 14 includes a volume rendering (VR) process.
- VR volume rendering
- the image memory 15 is a memory for storing image data for display generated by the image generation unit 14.
- the image memory 15 can also store data generated by the B-mode processing unit 12 and the Doppler processing unit 13.
- the B-mode data and Doppler data stored in the image memory 15 can be called by an operator after diagnosis, for example, and become ultrasonic image data for display via the image generation unit 14.
- the image generation unit 14 associates the ultrasonic image data and the time of the ultrasonic scanning performed for generating the ultrasonic image data with the electrocardiographic waveform transmitted from the electrocardiograph 4.
- the internal storage unit 16 stores a control program for performing ultrasonic transmission / reception, image processing and display processing, diagnostic information (for example, patient ID, doctor's findings, etc.), various data such as a diagnostic protocol and various body marks. To do.
- the internal storage unit 16 is also used for storing image data stored in the image memory 15 as necessary. Data stored in the internal storage unit 16 can be transferred to an external device via an interface (not shown).
- the external device is, for example, a high-performance workstation for image processing, a PC (Personal Computer) used by a doctor who performs image diagnosis, a storage medium such as a CD or a DVD, a printer, or the like.
- the image processing unit 17 is installed in the apparatus main body 10 in order to perform computer-aided diagnosis (CAD).
- CAD computer-aided diagnosis
- the image processing unit 17 acquires ultrasonic image data stored in the image memory 15 and performs image processing for diagnosis support. Then, the image processing unit 17 stores the image processing result in the image memory 15 or the internal storage unit 16.
- the image processing unit 17 according to the first embodiment is installed to provide exercise information of a moving tissue.
- the moving tissue is a tissue that moves periodically, such as the heart.
- the image processing unit 17 according to the first embodiment includes an acquisition unit 171, a setting unit 172, a tracking unit 173, and an exercise information calculation unit 174. The processing performed by the image processing unit 17 will be described in detail later.
- the control unit 18 controls the entire processing of the ultrasonic diagnostic apparatus. Specifically, the control unit 18 is based on various setting requests input from the operator via the input device 3 and various control programs and various data read from the internal storage unit 16. Controls the processing of the processing unit 12, Doppler processing unit 13, image generation unit 14, and image processing unit 17. Further, the control unit 18 controls the monitor 2 to display the ultrasonic image data for display stored in the image memory 15 and the internal storage unit 16. Further, the control unit 18 performs control so that the processing result of the image processing unit 17 is displayed on the monitor 2 or is output to an external device.
- the image processing unit 17 included in the ultrasonic diagnostic apparatus according to the first embodiment is time-sequential in order to provide motion information of a moving tissue (for example, a periodically moving tissue).
- the region of interest is tracked with the ultrasonic data group (moving image data of the ultrasonic image data) along.
- the image processing unit 17 according to the first embodiment uses the moving image data of the ultrasound image data to provide the motion information of the heart wall that periodically moves, follow the outline of the myocardium.
- the moving image data to be processed by the image processing unit 17 may be a two-dimensional ultrasonic image data group or a three-dimensional ultrasonic image data group.
- the acquisition unit 171 illustrated in FIG. 1 acquires moving image data of two-dimensional or three-dimensional ultrasonic image data collected using a region including a moving tissue as an imaging target.
- the moving image data is a two-dimensional or three-dimensional image over a period of at least one cycle collected from a periodically moving tissue that is a moving tissue. Sound image data. That is, the acquisition unit 171 acquires moving image data of two-dimensional or three-dimensional ultrasonic image data over a section of at least one cycle collected using a region including a tissue that moves periodically as an imaging target.
- the operator performs a two-dimensional scan or a three-dimensional scan of a region including the heart of the subject P with a sector probe, and a moving image of two-dimensional or three-dimensional ultrasound image data in which the myocardium is depicted.
- the moving image data is, for example, a group of ultrasonic image data collected in the B mode.
- the image generation unit 14 generates myocardial moving image data and stores the generated moving image data in the image memory 15.
- Ultrasound image data of the left ventricular short-axis section at the papillary muscle level includes a part of the short-axis surface of the right ventricle (RV: Right Ventricular) and the entire short-axis surface of the left ventricle (LV: Left Ventricular). It is drawn.
- the left ventricular short-axis ultrasonic image data at the papillary muscle level depict the anterior papillary muscle (APM: Anterolateral Papillary Muscle) and the posterior papillary muscle (PPM: Posteromedial Papillary Muscle).
- APM Anterolateral Papillary Muscle
- PPM Posteromedial Papillary Muscle
- the operator sets, for example, one heartbeat section from the first end diastole (ED) to the next end diastole as the section to be tracked.
- ED0 the time phase of the first end diastole
- ED1 the time phase of the next end diastole
- FIG. 2 is a diagram illustrating an example of moving image data acquired by the acquisition unit.
- this embodiment is applicable even when the acquisition unit 171 acquires moving image data of three-dimensional ultrasonic image data.
- the present embodiment can be applied even when the section to be tracked is set as a 2-beat section or a 3-beat section.
- the setting unit 172 illustrated in FIG. 1 sets a two-dimensional or three-dimensional region of interest using at least one time-phase ultrasonic image data among the moving image data acquired by the acquisition unit 171.
- the setting unit 172 sets a two-dimensional region of interest.
- the setting unit 172 sets a 2D or 3D region of interest.
- the region of interest is the contour of the left ventricular membrane or the contour of the myocardium.
- the setting unit 172 sets a region of interest based on information manually set by the operator.
- FIG. 3 is a diagram illustrating an example of an initial contour set by the operator.
- the monitor 2 displays the image data “I (t1)” under the control of the control unit 18.
- the operator traces the position of the contour of the left ventricular membrane with reference to the image data “I (t1)” shown in the upper diagram of FIG.
- the operator traces the position of the contour of the myocardium by referring to the image data “I (t1)” shown in the upper diagram of FIG.
- the setting unit 172 sets the position of the contour of the left ventricular membrane traced by the operator as the initial contour position “C (t1)” of the left ventricular membrane, as shown in the lower diagram of FIG.
- the setting unit 172 sets the position of the contour of the myocardial epicardium traced by the operator as the initial contour position “D (t1)” of the myocardium.
- the outline traced by the operator is indicated by a thick curve.
- the contour position of the myocardial epicardium may be automatically generated as a position away from the endocardial position by a predetermined thickness (predetermined distance).
- the position of the contour of the middle layer of the myocardium may be set as the initial contour position. Further, the present embodiment is not limited to the case where the initial contour position manually set by the operator is used.
- the setting unit 172 determines the initial contour position based on the luminance of the ultrasonic image data. It is also possible to set automatically.
- the position of the initial region of interest is tracked by speckle tracking.
- the position of the region of interest of each ultrasonic image data is acquired.
- the degree of pattern change between frames and between volumes is large in the early diastole phase with a large heart strain rate and the peak systole phase.
- the tracking unit 173 illustrated in FIG. 1 performs a tracking process described below in order to accurately obtain the tracking result of the region of interest. That is, the tracking unit 173 illustrated in FIG. 1 acquires the first position information of the region estimated as the region of interest based on the motion information, using the remaining temporal phase ultrasonic image data. In addition, the tracking unit 173 illustrated in FIG. 1 acquires the second position information of the region estimated as the region of interest based on information other than the motion information (for example, shape information) in each ultrasonic image data in the remaining time phase. To do. That is, the tracking unit 173 functions as an estimation unit that acquires the first position information and the second position information. Then, the tracking unit 173 illustrated in FIG.
- the tracking unit 173 illustrated in FIG. 1 tracks the region of interest with each ultrasonic image data of the remaining time phase.
- the “reliability index” is also defined as “motion quality” which is the quality of motion information.
- the “reliability index” may be described as “motion quality”.
- FIG. 4 is a block diagram showing a configuration example of the tracking unit shown in FIG.
- the tracking unit 173 includes a first estimation unit 173a, a definition unit 173b, a second estimation unit 173c, and a combination unit 173d, for example, as shown in FIG.
- the 1st estimation part 173a and the 2nd estimation part 173c demonstrated below may be the case where it installs as an estimation part in the image process part 17, separately from the tracking part 173.
- the first estimation unit 173a performs a region of interest tracking process based on motion information, that is, a speckle tracking process. That is, the first estimation unit 173a includes first image data that is ultrasonic image data of the first time phase and second image that is ultrasonic image data of the second time phase that is temporally adjacent to the first time phase. Motion information is estimated with respect to the data, the region of interest of the first image data is moved based on the estimated motion information, and the first position information of the second image data is estimated.
- the first estimation unit 173a calculates the movement vector of each of a plurality of points in the region including the region of interest. More specifically, the first estimation unit 173a calculates a movement vector for each of a plurality of points by template matching processing. Then, the first estimation unit 173a, for each tracking point constituting the region of interest, removes a statistically abnormal movement vector from the movement vector group obtained in a segment of a predetermined size including the tracking point. A movement vector of the tracking point is estimated using a group. The first estimation unit 173a estimates the first position information using the movement vector of each tracking point obtained in this way. Then, the first estimation unit 173a estimates the movement vector of each tracking point constituting the region of interest using the remaining movement vector group obtained by removing the statistically abnormal movement vector from the obtained movement vector group. First position information is estimated.
- the definition unit 173b defines a reliability index (motion quality) based on at least one variable obtained by the processing of the first estimation unit 173a.
- the second estimation unit 173c estimates the second position information of the second image data based on information (for example, shape information) other than the motion information.
- the combination unit 173d acquires position information obtained by combining the first position information and the second position information based on the reliability index (motion quality) as the position information of the region of interest in the second image data.
- FIG. 5A to FIG. 16 are diagrams for explaining the tracking unit according to the first embodiment.
- the initial contour position “D (t1)” of the outer myocardium can be similarly tracked by the processing described below. For this reason, the description of the tracking process of the contour position of the myocardium is omitted.
- the first estimation unit 173a sets a plurality of tracking points at the initial contour position “C (t1)”.
- the first estimation unit 173a sets a plurality of points in the vicinity of the initial contour position “C (t1)”.
- Each of the plurality of points illustrated in FIG. 5B is a calculation point for calculating motion information (movement vector).
- the first estimating unit 173a sets the first time phase to “t1” and the second time phase to “t2”. Then, the first estimation unit 173a estimates motion information between “I (t1)” that is the first image data and “I (t2)” that is the second image data, and based on the estimated motion information. The initial contour position “C (t1)” of I (t1) is moved. Thereby, the first estimation unit 173a estimates the first position information “Cvector (t2)” of “I (t2)” as illustrated in FIG.
- the first estimation unit 173a sets template data for each of the plurality of calculation points set in “I (t1)”.
- the template data is composed of a plurality of pixels centered on the calculation point.
- the first estimation unit 173a searches for an area that most closely matches the speckle pattern of the template data between two frames, that is, between “I (t1)” and “I (t2)”. Track where the template data has moved in the next frame.
- the first estimation unit 173a calculates the movement vector of each of the plurality of calculation points of “I (t1)”.
- a set of movement vectors of each of the plurality of calculation points of “I (t1)” is denoted by ⁇ Mv0 ⁇ .
- ⁇ Mv0 ⁇ is a group of movement vectors calculated at each of a plurality of calculation points.
- FIG. 7 shows an example of the movement vector group ⁇ Mv0 ⁇ calculated by the first estimation unit 173a.
- white squares indicate calculation points, and line segments extending from the white squares indicate movement vectors.
- the movement vector group illustrated in FIG. 7 schematically shows the state of the movement vector group calculated between the image data in the early dilation in which the heart lumen expands.
- the 1st estimation part 173a divides
- the first estimation unit 173a uses a segmented region recommended by the American Echocardiographic Society or the American Heart Association to convert the left ventricular myocardium having a short-axis cross section into “ANT-SEPT (anterior wall). It is divided into six segments: “septum”, ANT (front wall), LAT (side wall), POST (rear wall), INF (lower wall), and SEPT (septum). Accordingly, each tracking point shown in FIG. 5A is included in one of the six segments illustrated in FIG. In addition, each calculation point illustrated in FIG. 5B is included in one of the six segments illustrated in FIG.
- the 1st estimation part 173a calculates the statistical value of each segment by carrying out statistical processing of the movement vector group contained in each segment, and specifies the movement vector from which a statistical value becomes an abnormal value in each segment. Then, the first estimation unit 173a sets ⁇ Mv ′ ⁇ to the identified set of movement vectors, that is, the set of movement vectors whose statistical values are abnormal values. Then, the first estimation unit 173a sets a set of movement vectors ⁇ Mv ⁇ obtained by removing ⁇ Mv ' ⁇ from ⁇ Mv0 ⁇ as a remaining movement vector group.
- the first estimation unit 173a calculates the variance “Vari (s)” of the movement vector group in each segment.
- the first estimation unit 173a has six variance values “Vari (s1),. ) ”.
- the first estimation unit 173a specifies a movement vector exceeding “ ⁇ ⁇ Vari (s1)” as an abnormal movement vector.
- ⁇ is a predetermined coefficient
- a value of “2” is a suitable example as the value of “ ⁇ ”.
- the first estimation unit 173a similarly performs this process for each of the segments “s2, s3, s4, s5, s6”.
- the first estimation unit 173a identifies ⁇ Mv ′ ⁇ that is an abnormal movement vector group and removes ⁇ Mv ′ ⁇ from ⁇ Mv0 ⁇ , whereby the remaining movement vector group ⁇ Mv ⁇ is acquired.
- FIG. 9 shows an example of the abnormal movement vector group specified by the first estimation unit 173a.
- four abnormal movement vectors identified by the first estimation unit 173a are indicated by dotted lines.
- the abnormal movement vector group ⁇ Mv ′ ⁇ illustrated in FIG. 9 is the movement vector group ⁇ Mv0 ⁇ illustrated in FIG. 7, that is, a movement vector group calculated between image data in the early diastole in which the heart lumen expands. Is an abnormal movement vector group identified from For this reason, the remaining movement vector group ⁇ Mv ⁇ illustrated in FIG. 9 is a movement vector that generally goes in the outward direction.
- the first estimation unit 173a estimates the movement vectors of the plurality of tracking points constituting the initial contour position “C (t1)” using the remaining movement vector group ⁇ Mv ⁇ .
- the movement vector (residual movement vector) of each element constituting ⁇ Mv ⁇ will be generalized as “Mv”.
- an arbitrary tracking point among the plurality of tracking points constituting “C (t1)” is set to “p (t1)”.
- the first estimation unit 173a averages several movement vectors located in the vicinity of p (t1) from the remaining movement vector group ⁇ Mv ⁇ to obtain the movement vector “Mv_p (t1)” of p (t1). get.
- p (t1) is indicated by a circle hatched with a point.
- a plurality of calculation points located in the vicinity of p (t1) are indicated by white circles, and the remaining movement vectors at the respective calculation points are indicated by arrows extending from the white circles.
- the first estimation unit 173a sets a circle with a predetermined radius centered on p (t1) (see the dotted circle shown in FIG. 10).
- the length of the predetermined radius may be a value set in advance in the system or a value set by the operator.
- the first estimating unit 173a averages a plurality of remaining movement vectors (four remaining movement vectors in FIG. 10) included in the circle set by p (t1), whereby the movement vector “p (t1)” Mv_p (t1) "is estimated.
- the 1st estimation part 173a performs the process mentioned above in all the tracking points using ⁇ Mv ⁇ .
- the 1st estimation part 173a estimates the position of the point p (t2) which moved p (t1) by Mv_p (t1) as 1st position information, as shown in the upper figure of FIG.
- the first estimation unit 173a acquires the first position information “Cvector (t2)” of “I (t2)” as illustrated in the lower diagram of FIG.
- an arbitrary tracking point is assumed to be “p (t2)” among a plurality of tracking points constituting “Cvector (t2)”.
- the definition unit 173b defines the motion quality “Q_p (t2)” that is the reliability index of “p (t2)” using “at least one variable obtained by the processing of the first estimation unit 173a”.
- the first estimation unit 173a estimates the movement vector by template matching. That is, the first estimation unit 173a moves at each calculation point by searching for a position where the most similar pattern exists in the next time phase with respect to the density pattern of the image signal value in the template in the current time phase. Detect vectors. Therefore, the definition unit 173b determines the reliability “R (Mv) of the individual remaining movement vectors“ Mv ”constituting the remaining movement vector group ⁇ Mv ⁇ from the variables obtained by the template matching process performed by the first estimation unit 173a. Is defined. The variables used for defining the reliability include the three variables described below.
- the first variable used for defining the reliability is the average signal value of the reference template in the template matching process. That is, when the signal value of the template is extremely small, this template is not a tissue region effective for matching, but is considered to be an intracardiac region or a region of white noise. Therefore, when the average value “Tm” of the signals in the template is small, it is appropriate to reduce the value of “R (Mv)”.
- the definition unit 173b defines “R (Mv)” by the following equation (1).
- “Mth” shown in the following equation (1) is a lower limit threshold set for the signal average value in the reference template.
- Mth is stored in advance in the internal storage unit 16.
- “Mth” is set by an operator, for example.
- the definition unit 173b defines the reliability as “1” for the remaining movement vector calculated using the reference template whose “Tm” is equal to or greater than “Mth”, and the reference whose “Tm” is smaller than “Mth”. For the residual movement vector calculated using the template, the reliability is defined as “0”.
- the second variable used for defining the reliability is the signal variance value of the reference template in the template matching process. That is, when the light and shade pattern of the signal value of the reference template is too uniform, the light and shade pattern in the next time phase search range is considered to be nearly uniform. In such a case, since the most similar position is likely to be detected everywhere in the second image data, it is impossible to accurately estimate the movement vector. Therefore, when the variance value “T ⁇ ” of the signal in the reference template is small, it is appropriate to decrease the value of “R (Mv)”.
- the definition unit 173b defines “R (Mv)” by the following equation (2).
- “Sth” shown in the following equation (2) is a lower limit threshold set for the signal variance value in the reference template.
- Sth is stored in advance in the internal storage unit 16.
- “Sth” is set by an operator, for example.
- the definition unit 173b defines the reliability as “1” for the remaining movement vector calculated using the reference template having “T ⁇ ” equal to or greater than “Sth”, and the criterion “T ⁇ ” is smaller than “Sth”. For the residual movement vector calculated using the template, the reliability is defined as “0”.
- the third variable used for defining the reliability is the similarity between templates in the template matching process.
- the first estimation unit 173a searches the second image data for the position most similar to the reference template by comparing the reference template with the search region of the next time phase. Therefore, it is possible to define a similarity indicating how much the area of the finally searched position is similar to the reference template.
- a typical similarity there is a cross-correlation coefficient “Ck”, and when the input signal is a real number, “Ck” becomes “0 ⁇ Ck ⁇ 1”.
- the input signal is a real number, and “0 ⁇ Ck ⁇ 1”.
- Lowering the similarity means that the pattern change of the tracking target part is large over time from the current phase to the next time phase. Therefore, when the similarity “Ck” is small, it is appropriate to reduce the value of “R (Mv)”.
- the definition unit 173b defines “R (Mv)” by the following equation (3).
- “Cth” shown in the following equation (3) is a lower limit threshold set for the similarity in template matching.
- “Cth” is stored in advance in the internal storage unit 16.
- “Cth” is set by an operator, for example.
- the definition unit 173b defines the reliability as “Ck” for the remaining movement vector in which “Ck” is equal to or greater than “Cth”, and for the remaining movement vector in which “Ck” is smaller than “Cth”, The reliability is defined as “0”.
- R (Mv) is simply given by “1” or “0” depending on the threshold value.
- the present embodiment may be a case where “R (Mv)” is defined as a continuous monotonically increasing function using the above control variables.
- “R (Mv)” can also be defined by a combination using at least one of the three control variables. In the present embodiment, it is preferable to give the final “R (Mv)” using all the three control variables. In such a case, “R (Mv)” is defined by a function having three control variables “Tm, T ⁇ , Ck”.
- the defining unit 173b performs quality evaluation (reliability evaluation) of each remaining movement vector by defining “R (Mv)” of each remaining movement vector as described above. Furthermore, as will be described below, the definition unit 173b performs quality evaluation (reliability evaluation) of a motion vector in an area (segment) having a predetermined size using spatial continuity of a plurality of motion vectors. Do.
- the first estimation unit 173a calculates the variance “Vari (s)” of the motion vector group in each segment, and the motion vector having an abnormal value exceeding a predetermined threshold from ⁇ Mv0 ⁇ .
- the group ⁇ Mv ′ ⁇ is removed as an abnormal movement vector group to obtain a remaining movement vector group ⁇ Mv ⁇ .
- the defining unit 173b sets the variable obtained by the statistical processing performed by the first estimating unit 173a as a variable for defining the motion quality (reliability index).
- the definition unit 173b uses the above six segments “s1,..., S6” as the neighborhood region. In such a case, the definition unit 173b obtains the variance “Vari (s)” of the movement vector group in each segment from the first estimation unit 173a. That is, the definition unit 173b acquires “Vari (s1),..., Vari (s6)” of each of the six segments “s1,.
- the definition unit 173b uses the above six segments “s1,..., S6” as the neighborhood region. In such a case, the definition unit 173b acquires, for example, the number of ⁇ Mv0 ⁇ in each segment “s” and the number of ⁇ Mv ⁇ in each segment “s” from the first estimation unit 173a. The density “ ⁇ (s)” of the effective movement vector in each segment “s” is calculated. For example, the definition unit 173b defines “ ⁇ (s)” by the following equation (4).
- the definition unit 173b defines a value obtained by dividing the number of ⁇ Mv0 ⁇ in “s” by the number of ⁇ Mv ⁇ in the segment “s” as “ ⁇ (s)”.
- the defining unit 173b acquires “ ⁇ (s1),..., ⁇ (s6)” of each of the six segments “s1,.
- the definition unit 173b uses each of the calculation points calculated by the first estimation unit 173a using at least one of “R (Mv)”, “ ⁇ (s)”, and “Vari (s)” described above.
- the motion quality “Q (s, Mv)” of the movement vector is calculated.
- the defining unit 173b defines “Q (s, Mv)” by the following equation (5).
- “Lth” shown in the following equation (5) is a lower limit threshold set for the spatial dispersion value of the movement vector.
- “Lth” is stored in advance in the internal storage unit 16.
- “Lth” is set by an operator, for example.
- the definition unit 173b sets the motion quality to “0” for the movement vector of the calculation point included in the segment in which “Vari (s)” is smaller than “Lth”. For example, when “Vari (s2)” is smaller than “Lth”, the definition unit 173b sets the motion quality to “0” for the movement vector of the calculation point included in “s2”.
- the definition unit 173b calculates motion quality for each movement vector for the movement vector of the calculation point included in the segment whose “Vari (s)” is “Lth” or more. For example, when “Vari (s5)” is equal to or greater than “Lth”, the definition unit 173b sets “ ⁇ (Mv)” to the reliability “R (Mv)” of the movement vector “Mv” of each calculation point included in “s5”. By multiplying (s5) ”,“ Q (s5, Mv) ”is calculated.
- the definition unit 173b defines the motion quality (reliability index) of each remaining movement vector used for the estimation of “Cvector (t2)”.
- the definition unit 173b may define the motion quality of each remaining movement vector using “R (Mv)”.
- the definition unit 173b may define the motion quality of each residual motion vector using at least one of “Vari (s)” and “ ⁇ (s)”.
- the definition unit 173b averages the motion qualities of several residual motion vectors located in the vicinity of p (t1), and thereby the motion that is a reliability index of “p (t2) moved from p (t1)” Quality “Q_p (t2)” is defined. For example, as illustrated in FIG. 12, the definition unit 173b sets the average value of the motion quality of each of the four motion vectors (residual motion vectors) used in the calculation of “Mv_p (t1)” as “Q_p (t2)”. Calculate.
- the process of the 1st estimation part 173a and the process of the definition part 173b are performed in parallel for every time phase used as tracking object.
- Such parallel processing can be realized by the control of the control unit 18.
- the number and size of segments used in the above processing can be changed to an arbitrary number and size by an operator, for example.
- the segment used in the process of the first estimation unit 173a may be different from the segment used in the process of the definition unit 173b.
- the second estimation unit 173c uses the position information obtained in the first process described below, the position information obtained in the second process, and the position information obtained in the third process to use the first position information. 2 Estimate position information.
- the second estimation unit 173c collates the second image data of the time phase “t” with the shape dictionary information “dict (t)” of the simultaneous phase, so that the position of the region of interest in the second image data
- the information “Cdict (t)” is estimated.
- the internal storage unit 16 stores in advance, for example, a shape dictionary of the left ventricular myocardium in all cardiac phases. Then, since the second image data is “I (t2)”, the second estimation unit 173c stores the cardiac phase shape dictionary “dict (t2)” corresponding to the time phase “t2” in the internal storage unit 16. Get from. Then, as shown in FIG.
- the second estimation unit 173c collates “I (t2)” with “dict (t2)” by the function of the discriminator, and the region of interest (“I (t2)”) ( (Contour) position information “Cdict (t2)” is estimated.
- I (t2) region of interest
- Contour region of interest
- the second estimation unit 173c estimates the position information “Ceng (t)” of the region of interest in the second image data of the time phase “t” based on the principle of shape energy minimization.
- the second estimation unit 173c uses the SNAKE method that is known as a kind of an active contour model (ACM).
- ACM active contour model
- the second estimation unit 173c sets a contour line near the boundary of the heart lumen in the second image data to be processed.
- the contour line to be set may be set based on C (t1) or may be set based on luminance information.
- the second estimation unit 173c repeats the correction of the position and shape of the contour line until the value of the shape energy evaluation function is minimized. For example, the second estimation unit 173c determines that the contour line is located where the luminance changes rapidly, the luminance value in the contour line is uniform, the contour line is short, and the shape of the contour line is smooth. Repeat the position and shape correction. Thereby, as shown in FIG. 13B, the second estimation unit 173c estimates the position information “Ceng (t2)” of the region of interest (outline) of “I (t2)” that is the second image data.
- the second estimation unit 173c estimates the position information “Cfit (t)” of the region of interest in the second image data of the time phase “t” by a process including fitting based on the least square method.
- a plurality of control points “r (t)” including position information that is accurate to some extent are required as fitting inputs.
- the 2nd estimation part 173c performs a 3rd process only with respect to the site
- the motion quality “Q_p (t2)” is defined at each tracking point p (t2) constituting “Cvector (t2)” of “I (t2)” that is the second image data. Therefore, the second estimation unit 173c identifies p (t2) where “Q_p (t2)” is equal to or less than a predetermined value. For example, as illustrated in the left diagram of FIG. 13C, the second estimation unit 173c identifies p (t2) included in a dotted circle. Then, the second estimation unit 173c sets a set of p (t2) other than p (t2) included in the dotted circle as a control point group “r (t2)” for performing the fitting process. Perform a least-squares fitting using “r (t2)”. Thereby, for example, as illustrated in the left diagram of FIG. 13C, the second estimation unit 173c obtains the position information “Cfit (t2)” of the region of interest (contour) of “I (t2)” that is the second image data. presume.
- the third process in the case where the motion quality is degraded due to the influence of local noise on Cvector (t), only the contour position of the part where the motion quality is degraded is determined as Cvector (t) having a high motion quality in the vicinity.
- Cfit (t) By appropriately interpolating with the position of “Cfit (t)”, “Cfit (t)” is obtained.
- the third process is a process of estimating the contour shape of the part where the motion quality (reliability index) is lowered based on the contour shape of the peripheral motion quality (reliability index).
- the second estimation unit 173c performs the second process only when the first process is performed to obtain “Cdict (t)” only at a local position where the motion quality is reduced on Cvector (t). To obtain “Ceng (t)”.
- the second estimation unit 173c uses the at least one of “Cdict (t)”, “Ceng (t)”, and “Cfit (t)” to store the second position information “Cmodel (t) of the second image data. ) ”. Since the second image data is “I (t2)”, the second estimation unit 173c uses at least one of “Cdict (t2)”, “Ceng (t2)”, and “Cfit (t2)”. , “I (t2)” second position information “Cmodel (t2)” is defined.
- the second estimation unit 173c sets “Cdict (t2)” to “Cmodel (t2)”.
- the second estimation unit 173c sets “Ceng (t2)” to “Cmodel (t2)”.
- the second estimation unit 173c sets “Cfit (t2)” to “Cmodel (t2)”.
- the second estimation unit 173c can estimate “Cfit (t)” if “Cvector (t)” and “Q_p (t)” are obtained. Further, the second estimation unit 173c estimates “Cdict (t)” and “Ceng (t)” even when “Cvector (t)” and “Q_p (t)” are not obtained. It is possible. Furthermore, the second estimation unit 173c can perform the estimation of “Cdict (t)” and the estimation of “Ceng (t)” independently of each other.
- the second estimation unit 173c can estimate these three pieces of position information. Therefore, the second estimation unit 173c may give “Cmodel (t)” by linearly adding all three pieces of position information or two pieces of position information selected from these three pieces of position information. For example, the second estimation unit 173c gives “Cmodel (t)” by, for example, averaging a plurality of contour positions.
- the second estimation unit 173c sets the position information obtained by the first process as the second position information in the specific time phase having the shape dictionary information.
- the second estimation unit 173c sets the position information obtained by the first process as the second position information in a specific time phase in which the degree of local deformation of the shape of the tissue (myocardium) is large.
- the second estimation unit 173c combines the position information obtained by the second process or the position information obtained by the second process and the position information obtained by the third process in a time phase other than the specific time phase. It is preferable that the obtained position information is the second position information.
- the specific time phase is a typical time phase in which the rate of deformation of the heart is large, a large early diastole (e ′) time phase that is likely to be untracked, and a contraction peak (s ′ ) Time phase.
- the second estimation unit 173c changes “Cdict (t)” to “Cmodel (t) by the first process in the specific time phase (e ′, s ′). ) ”.
- the second estimation unit 173c sets “Ceng (t)” to “Cmodel (t)” in a time phase other than the specific time phase.
- the second estimation unit 173c performs linear addition of “Ceng (t)” and “Cfit (t)” in a time phase other than the specific time phase.
- the contour position thus obtained is defined as “Cmodel (t)”.
- the second estimation unit 173c can be executed in parallel with the processes of the first estimation unit 173a and the definition unit 173b.
- the second estimation unit 173c can be executed independently of the processes of the first estimation unit 173a and the definition unit 173b. That is, the second estimation unit 173c may perform the first process and the second process in all remaining time phases when the acquisition unit 171 acquires the moving image data.
- the second estimation unit 173c be performed in parallel with the processes of the first estimation unit 173a and the definition unit 173b in synchronization with each time phase to be tracked.
- the progress of processing of each unit can be controlled by the control unit 18.
- the first estimation unit 173a generates a plurality of contour positions “Cvector (t2)” based on the contour position “C (t1)” and the motion information ⁇ Mv ⁇ , as shown in the upper diagram of FIG. As a set of tracking points “p (t2)”.
- the second estimation unit 173c determines that the contour position “Cmodel (t2)” based on information other than motion information (for example, shape information) is “p (t2)”. Obtained as a set of contour control points “q (t2)” whose position and number are independent.
- “Cvector (t2)” is the first position information
- “Cmodel (t2)” is the second position information.
- a motion quality (reliability index) is defined for each of the plurality of tracking points “p (t2)”.
- the combination unit 173d combines the first position information and the second position information based on the reliability index (motion quality). Specifically, the combination unit 173d weights and adds the first position information and the second position information based on the reliability index (motion quality). More specifically, the combination unit 173d weights and adds the first position information and the second position information so that the weight of the first position information becomes high when the reliability index (motion quality) is high.
- the combination unit 173d acquires a plurality of points “cross_p (t2)” corresponding to the plurality of tracking points “p (t2)” in “Cmodel (t2)”. For example, as illustrated in the upper diagram of FIG.
- the combination unit 173d sets the normal vector “Vnormal_p (t2)” for each “p (t2)” in the contour defined by “Cvector (t2)”. Ask. Then, as illustrated in the upper diagram of FIG. 14, the combination unit 173d determines the intersection of the normal vector “Vnormal_p (t2)” and the contour defined by “Cmodel (t2)” as “cross_p (t2)”. "
- the combination unit 173d converts “Cmodel (t2)” into “Cmodel (t2) ′” including a plurality of points “cross_p (t2)”, as illustrated in the lower diagram of FIG.
- the combination unit 173d acquires a plurality of points “cross_p (t2)” on Cmodel (t2) corresponding to each of the plurality of tracking points “p (t2)”.
- “Q_p (t2)” of “p (t2)” becomes “Q_cross_p (t2)” which is the motion quality of “cross_p (t2)”.
- the first position information of the time phase “t” is generalized as “Cvector (t)” composed of a plurality of tracking points “p (t)”, and each tracking point “p (t)” Is assumed to be “Q_p (t)”.
- the second position information of the time phase “t” is “Cmodel (t)” composed of a plurality of control points “q (t)”, and a plurality of points “cross_p (t2)”.
- the second position information obtained by converting “Cmodel (t)” is assumed to be “Cmodel (t) ′”.
- the combination unit 173d weights and adds “Cvector (t)” and “Cmodel (t) ′” based on “Q_p (t)”, for example, by the following equation (6).
- the contour position “C (t)” of the time phase “t” is obtained. Note that “Q_p (t)” is normalized to “0 ⁇ Q_p (t) ⁇ 1” according to the above definition.
- Equation (6) a point determined by a position vector obtained by adding these two position vectors is estimated as a point constituting the contour position “C (t2)” in “I (t2)” of the time phase “t2”. Is done.
- FIG. 15 illustrates a contour position “C (t)” in which the combination unit 173d combines “Cvector (t)” and “Cmodel (t)” based on “Q_p (t)” in the short-axis image.
- FIG. FIG. 15 shows an example in which the left wall has high motion quality and the right wall has low motion quality. That is, in FIG. 15, the position of “C (t)” is reconfigured as a position close to “Cvector (t)” on the left and close to “Cmodel (t)” on the right.
- the combination unit 173d notifies the first estimation unit 173a of a plurality of points constituting the contour position “C (t2)” estimated by the equation (6) as a plurality of tracking points “p (t2)”. To do. Thereby, each unit constituting the tracking unit 173 performs processing in which “t1” is replaced with “t2” and “t2” is replaced with “t3” in the above description.
- the first estimation unit 173a performs speckle tracking between “I (t2)” that is the first image data and “I (t3)” that is the second image data.
- the defining unit 173b defines the motion quality “Q_ (pt3)” of “p (t3)”.
- the second estimation unit 173c estimates “Cmodel (t3)” of “I (t3)”. Then, as illustrated in FIG. 16, the combination unit 173d converts “Cmodel (t3)” to “Cmodel (t3) ′”, and “Cvector (t3)” and “Cmodel” based on “Q_p (t3)”. (T3) "is combined to estimate the contour position" C (t3) ".
- the motion information calculation unit 174 illustrated in FIG. 1 calculates the motion information of the tissue using the position information of the region of interest of each moving image data. For example, the motion information calculation unit 174 calculates various physical quantities such as strain and strain change rate (Strain Rate) as wall motion information.
- strain Rate strain change rate
- the control unit 18 displays the exercise information on the monitor 2.
- FIG. 17 is a flowchart for explaining an outline of processing performed by the ultrasonic diagnostic apparatus according to the first embodiment
- FIG. 18 illustrates an example of processing performed by the tracking unit according to the first embodiment. It is a flowchart for.
- the acquisition unit 171 included in the ultrasonic diagnostic apparatus determines whether moving image data and a tracking target section are specified (step S101).
- the acquisition unit 171 waits until it is specified.
- the acquisition unit 171 acquires the moving image data of the tracking target section (Step S102). Then, the setting unit 172 determines whether or not the initial contour position has been set by the operator (step S103). Here, when it is not set (No at Step S103), the setting unit 172 waits until the initial contour position is set.
- Step S104 the tracking unit 173 acquires the initial contour position and tracks the position of the remaining time phase contour.
- the motion information calculation unit 174 calculates heart wall motion information (step S105), and the monitor 2 displays the heart wall motion information under the control of the control unit 18. (Step S106), and the process ends.
- the 1st estimation part 173a calculates ⁇ Mv0 ⁇ by ST from I (tk) to I (tk + 1) (step S202). Then, the first estimation unit 173a obtains ⁇ Mv ⁇ that is the remaining movement vector group by statistical processing (step S203), calculates Mv_p (tk) at each tracking point using ⁇ Mv ⁇ , and Cvector (Tk + 1) is estimated (step S204).
- the definition unit 173b calculates a motion quality Q_p (tk + 1) that is a reliability index of each tracking point p (tk + 1) constituting Cvector (tk + 1) while the process of step S204 is performed (step S204).
- the second estimating unit 173c estimates Cmodel (tk + 1) of I (tk + 1) while the processes of Steps S202 to S204 and Step S205 are being performed (Step S206).
- the combination unit 173d converts Cmodel (tk + 1) to Cmodel (tk + 1) ′ using Cvector (tk + 1) (step S207). Then, the combination unit 173d combines Cvector (tk + 1) and Cmodel (tk + 1) ′ based on the motion quality Q_p (tk + 1) of each tracking point p (tk + 1), and the contour position C (tk + 1) of I (tk + 1) ) Is acquired (step S208).
- Step S209 the combination unit 173d notifies the motion information calculation unit 174 that the contour tracking has been completed in all time phases, and starts the process of Step S105.
- the position of the region of interest is estimated based on information other than the motion information as well as the position of the region of interest based on the motion information.
- the quality (reliability index) of the motion information is defined from various information, and the estimated position based on the motion information and the estimated position based on information other than the motion information are defined as the quality (reliability index).
- the local position where the motion information is insufficient, the local position where the reliability of the motion information is low, and the local where the spatial and temporal consistency of the motion information is low Specific positions are identified from the motion quality (reliability index), and for these positions, the final contour position is determined by using the contour positions obtained from information other than the motion.
- the position of the region of interest without extreme tracking error is automatically acquired. Can do. Therefore, in the first embodiment, the tracking result of the region of interest can be obtained accurately. In the first embodiment, since the position of the region of interest without extreme tracking error can be automatically acquired, it is possible to easily provide a robust tissue motion evaluation.
- the tracking unit 173 may perform tracking of the region of interest in the backward direction (backward).
- the first estimation unit 173a can calculate the forward movement vector “MvF” from the time phase “tk” to the time phase “tk + 1”, and can calculate the time phase “t + 1” from the time phase “t + 1”.
- the movement vector “MvB” in the reverse direction toward “tk” can be calculated. That is, in the first embodiment, the time vector “MvF_p (tk ⁇ tk + 1)” is estimated from “MvF” with respect to the tracking point p (tk) of the time phase “tk”. The position of the tracking point p (tk + 1) of “tk + 1” was estimated.
- the first estimating unit 173a can calculate the movement vector “MvB” in the reverse direction from the time phase “tk + 1” to the time phase “tk”. Then, the first estimation unit 173a can estimate the movement vector “MvB_p (tk + 1 ⁇ tk)” with respect to the tracking point p (tk + 1) of the time phase “tk + 1” from “MvB”. Thereby, the first estimation unit 173a can estimate the position of the tracking point p (tk) of the time phase “tk”.
- the definition unit 173b can provide information defining motion quality, which is a reliability index, by using “MvF” and “MvB” together. That is, the position of p (tk) obtained from “MvB_p (tk + 1 ⁇ tk)” does not necessarily match the position of p (tk) moved to p (tk + 1) by “MvF_p (tk ⁇ tk + 1)”. This is because “MvF” and “MvB” may not match due to the influence of speckle noise due to different image data to be searched in template matching. On the other hand, when the influence of noise is small, it is expected that the degree of coincidence of both movement vectors is high.
- the first estimation unit 173a performs a region of interest between the first image data and the second image data. Forward tracking vectors from the first time phase to the second time phase of the tracking points constituting the region, and backward movement vectors from the second time phase to the first time phase of the tracking points constituting the region of interest. presume. Then, the definition unit 173b uses the degree of coincidence between the forward direction movement vector and the backward direction movement vector at each tracking point constituting the region of interest as “at least one variable obtained by the processing of the first estimation unit 173a”.
- the first estimation unit 173a performs the calculation of “MvF”, the estimation of “MvF_p (tk ⁇ tk + 1)”, and the position estimation of p (tk + 1). The process described in the above embodiment is performed. Then, the first estimation unit 173a performs the process described in the above embodiment on p (tk + 1) determined by “MvF” in the reverse direction from the time phase “tk + 1” to the time phase “tk”. Execute by reversing the time phase. As a result, the first estimation unit 173a calculates “MvB” and estimates “MvB_p (tk + 1 ⁇ tk)”.
- both “MvF_p (tk ⁇ tk + 1)” and “MvB_p (tk + 1 ⁇ tk)” can be variables at the tracking point p (tk + 1) in the time phase “tk + 1”. That is, when generalized, these two variables can be expressed as two variables “MvF_p (t)” and “MvB_p (t)” for the tracking point p (t) in the time phase “t”.
- the definition unit 173b calculates the degree of coincidence “MM_p (t)” of both the forward and backward movement vectors of the tracking point p (t) by the following equation (7). Can be defined.
- “ ⁇ ” indicates an inner product of vectors.
- the definition unit 173b defines the motion quality “Q (s, Mv)” by the following equation (8) using “MM_p (t)”.
- the definition unit 173b sets each residual movement vector “Mv” in the segment “s” including the tracking point “p (t)” from which the value “MM_p (t)” is “0” or more is obtained.
- the motion quality “Q (s, Mv)” is defined as “MM_p (t)”.
- the defining unit 173b determines that each remaining movement vector “Mv” in the segment “s” including the tracking point “p (t)” from which the value “MM_p (t)” is “negative” is obtained.
- the motion quality “Q (s, Mv)” is defined as “0”. Note that the motion quality defined by Expression (8) satisfies the condition “0 ⁇ Q (s, Mv) ⁇ 1”. Note that the definition unit 173b may define the motion quality “Q (s, Mv)” using Equation (5) and Equation (8) together.
- the definition unit 173b defines “Q (s, Mv)” defined by Expression (8) or “Q (s, Mv)” defined by the combination of Expression (5) and Expression (8). Is used to define “Q_P (t)”.
- the initial contour position may be set in a plurality of time phases.
- the first initial contour position and the second initial contour position may be set in the time phase in the early diastole and the time phase in the end systole where speckle tracking is difficult.
- the tracking unit 173 performs the tracking process using the first initial contour position and the tracking process using the second initial contour position in parallel.
- the process executed by the second estimation unit 173c is not limited to the first process, the second process, and the third process described above.
- the second estimation unit 173c may execute an AQ (Acoustic Quantification) method or an ASM (Active Shape Model) method.
- the second embodiment a control method for reducing the calculation processing amount will be described.
- the first embodiment the case where the second position information “Cmodel (t)” is calculated for all the remaining time phases has been described.
- the amount of calculation increases when all the remaining time phases are calculated from the second position information.
- the first process of estimating “Cdict (t)” using the classifier is a process that requires a discrimination time by the classifier although the accuracy in the estimation of the contour position is increased.
- the ultrasonic diagnostic apparatus according to the second embodiment has the same configuration as the ultrasonic diagnostic apparatus according to the first embodiment illustrated in FIG.
- the tracking unit 173 according to the second embodiment has the same configuration as the tracking unit 173 according to the first embodiment illustrated in FIG.
- the tracking unit 173 according to the second embodiment performs processing described below under the control of the control unit 18. That is, the tracking unit 173 according to the second embodiment outputs the first position information as the position information of the region of interest for the time phase in which the average value of the reliability index (motion quality) in the region of interest is equal to or greater than the predetermined threshold. To do. Further, the tracking unit 173 according to the second embodiment, for a time phase in which the above average value is smaller than a predetermined threshold, uses position information obtained by combining the first position information and the second position information as position information on the region of interest. Output as.
- the processing of the acquisition unit 171 and the setting unit 172 is performed as in the first embodiment.
- the first estimation unit 173a estimates the first position information “Cvector (t)” of “I (t2)”.
- the definition unit 173b defines the motion quality “Q_p (t2)” of each point “p (t2)” constituting “Cvector (t)”.
- the defining unit 173b calculates “Q (t2)” obtained by averaging the motion quality “Q_p (t2)” at each point “p (t2)”.
- the control unit 18 acquires the value of “Q (t2)” from the definition unit 173b, and compares the quality threshold value “Qth” with “Q (t2)”. “Qth” is stored in advance in, for example, the internal storage unit 16. Alternatively, “Qth” is set by an operator, for example.
- the control unit 18 stops the processes of the second estimation unit 173c and the combination unit 173d, and changes “Cvector (t2)” estimated by the first estimation unit 173a to “C”. (T2) ".
- the first estimation unit 173a outputs “Cvector (t2)” as “C (t2)” to the motion information calculation unit 174.
- the control unit 18 starts the processing of the second estimation unit 173c and the combination unit 173d. Accordingly, the second estimation unit 173c estimates “Cmodel (t2)”, and the combination unit 173d obtains “Cmodel (t2) ′” from “Cvector (t2)” and “Cmodel (t2)”. Then, the combination unit 173d weights and adds “Cvector (t2)” and “Cmodel (t2) ′” based on “Q_p (t2)”, acquires “C (t2)”, and exercises The information is output to the information calculation unit 174.
- FIG. 19 is a flowchart for explaining an example of processing performed by the tracking unit according to the second embodiment.
- the 1st estimation part 173a calculates ⁇ Mv0 ⁇ by ST from I (tk) to I (tk + 1) (step S302). Then, the first estimation unit 173a obtains ⁇ Mv ⁇ that is the remaining movement vector group by statistical processing (Step S303), calculates Mv_p (tk) at each tracking point using ⁇ Mv ⁇ , and Cvector (Tk + 1) is estimated (step S304).
- the definition unit 173b calculates the motion quality Q_p (tk + 1) of each tracking point p (tk + 1) constituting Cvector (tk + 1) while the process of step S304 is being performed (step S305). Then, the definition unit 173b calculates an average value Q (tk + 1) of the motion quality Q_p (tk + 1) of each tracking point p (tk + 1) (step S306).
- control part 18 determines whether Q (tk + 1) is smaller than Qth (step S307). If Q (tk + 1) is equal to or greater than Qth (No at Step S307), the tracking unit 173 acquires Cvector (tk + 1) as the contour position C (tk + 1) of the image data I (tk + 1) (Step S308).
- the second estimation unit 173c estimates Cmodel (tk + 1) of I (tk + 1) (Step S309).
- the combination unit 173d converts Cmodel (tk + 1) to Cmodel (tk + 1) ′ using Cvector (tk + 1) (step S310). Then, the combination unit 173d combines Cvector (tk + 1) and Cmodel (tk + 1) ′ based on the motion quality Q_p (tk + 1) of each tracking point p (tk + 1), and the contour position C (tk + 1) of I (tk + 1) ) Is acquired (step S311).
- step S312 the combination unit 173d notifies the motion information calculation unit 174 that the contour tracking has been completed in all time phases, and starts the process of step S105.
- the processing of the second estimation unit 173c and the combination unit 173d is omitted, so that the tracking result of the region of interest is accurate. And it can obtain efficiently.
- FIG. 20 is a diagram for explaining a modification according to the second embodiment.
- the tracking unit 173 sets Cvector (t) as the contour position C (t) of the image data I (t).
- the second estimating unit 173c performs the third process described with reference to FIG. 13C in accordance with an instruction from the control unit 18. That is, as shown in the right diagram of FIG. 20, the second estimation unit 173 c interpolates the contour in the dotted circle by the least square fitting process using surrounding points with high motion quality, thereby obtaining “I “Cfit (t)” of (t) ”is estimated. Then, the second estimation unit 173c outputs “Cfit (t)” as “Cmodel (t)” to the combination unit 173d.
- the third process is a process with a small amount of calculation compared to the first process and the second process. For this reason, by performing the above-described modification, it is avoided that Cvector (t) whose reliability is high overall but locally low is “C (t)”. “C (t)” having high reliability can be obtained quickly.
- the initial contour is set based on information other than the motion information described above. It may be performed by a method for estimating the position.
- the setting unit 172 may set the initial contour by various processes executed by the second estimation unit 173c. In such a configuration, accurate tracking of the region of interest can be further automated, and the operator can more easily perform tissue motion analysis.
- the initial contour is automatically set, it is preferable that the operator can accept a procedure for correcting the initial contour position.
- the image processing method described in the first and second embodiments can be executed in substantially real time together with the collection of moving image data. That is, the ultrasonic diagnostic apparatus according to the first and second embodiments performs reliability by performing the image processing method described in the first and second embodiments sequentially while collecting moving image data.
- a stable tissue motion method can be provided in substantially real time.
- the image processing methods described in the first and second embodiments are applicable even when tracking processing is performed on a plurality of two-dimensional moving image data in a plurality of cross sections. Further, as described above, the image processing methods described in the first and second embodiments are applicable even when tracking processing is performed on three-dimensional moving image data.
- the organ to which the above image processing method is applied is not limited to the heart, but may be an arterial blood vessel such as a carotid artery that repeatedly expands and contracts in synchronization with the cardiac cycle.
- the above image processing method can be used for analyzing the motion index such as strain and displacement to grasp the hardness and dynamics of the tissue such as “soft tissue such as liver and thyroid gland” and “muscle”. It may be applied to a moving tissue.
- the image processing method described above is applicable to two-dimensional or three-dimensional moving image data of medical image data capable of tracking a region of interest by template matching processing, such as an X-ray diagnostic apparatus, an X-ray CT apparatus, or an MRI apparatus. May be applied. That is, the image processing method described in the first and second embodiments may be performed by a medical image diagnostic apparatus other than the ultrasonic diagnostic apparatus. Further, the above image processing method may be executed by an image processing apparatus installed independently of the medical image diagnostic apparatus.
- each component of each illustrated apparatus is functionally conceptual, and does not necessarily need to be physically configured as illustrated.
- the specific form of distribution / integration of each device is not limited to the one shown in the figure, and all or a part of the distribution / integration is functionally or physically distributed in arbitrary units according to various loads or usage conditions. Can be integrated and configured.
- all or a part of each processing function performed in each device may be realized by a CPU and a program analyzed and executed by the CPU, or may be realized as hardware by wired logic.
- image processing methods described in the above embodiments and modifications can be realized by executing an image processing program prepared in advance on a computer such as a personal computer or a workstation.
- This image processing program can be distributed via a network such as the Internet.
- the image processing program can also be executed by being recorded on a computer-readable recording medium such as a hard disk, a flexible disk (FD), a CD-ROM, an MO, and a DVD, and being read from the recording medium by the computer. .
Landscapes
- Health & Medical Sciences (AREA)
- Life Sciences & Earth Sciences (AREA)
- Engineering & Computer Science (AREA)
- Physics & Mathematics (AREA)
- Molecular Biology (AREA)
- General Health & Medical Sciences (AREA)
- Pathology (AREA)
- Radiology & Medical Imaging (AREA)
- Biophysics (AREA)
- Biomedical Technology (AREA)
- Heart & Thoracic Surgery (AREA)
- Medical Informatics (AREA)
- Veterinary Medicine (AREA)
- Surgery (AREA)
- Animal Behavior & Ethology (AREA)
- Nuclear Medicine, Radiotherapy & Molecular Imaging (AREA)
- Public Health (AREA)
- Computer Vision & Pattern Recognition (AREA)
- Cardiology (AREA)
- Multimedia (AREA)
- General Physics & Mathematics (AREA)
- Theoretical Computer Science (AREA)
- Ultra Sonic Daignosis Equipment (AREA)
- Magnetic Resonance Imaging Apparatus (AREA)
- Image Analysis (AREA)
- Apparatus For Radiation Diagnosis (AREA)
Abstract
実施形態の超音波診断装置は、設定部(172)と、追跡部(173)とを備える。設定部(172)は、運動する組織を含む領域を撮影対象として収集された2次元又は3次元の超音波画像データの動画像データの中で、少なくとも1つの時相の超音波画像データにて2次元又は3次元の関心領域を設定する。推定部としての追跡部(173)は、前記動画像データが収集された区間のうち前記少なくとも1つの時相以外の残余時相の各超音波画像データで、動き情報に基づいて前記関心領域として推定した領域の第1位置情報と、前記動き情報以外の情報に基づいて前記関心領域として推定した領域の第2位置情報とを取得する。追跡部(173)は、前記動き情報の信頼度に関わる信頼度指標に基づいて前記第1位置情報と前記第2位置情報とを合成した位置情報を前記関心領域の位置情報として取得して、前記関心領域を追跡する。
Description
本発明の実施形態は、超音波診断装置、画像処理装置及び画像処理方法に関する。
近年、超音波画像の動画像に対して、2次元又は3次元のスペックルトラッキング(ST:Speckle Tracking)技術を適用することで、心臓の壁運動解析を行なうことが実用化されている。一般的なSTでは、拡張末期(最初のR波時相)、或いは、収縮末期の一時相にて心筋の内外膜の輪郭を初期輪郭位置として与える。そして、STでは、残りの時相では局所的なパターンマッチング処理や、オプティカルフロー法で得られる動き情報(移動ベクトル)を用いて、初期輪郭位置を自動的に追跡することで、必要な全時相での輪郭位置を取得する。従って、正確な輪郭位置の追跡の実現には、移動ベクトルを正確に推定することが欠かせない。
しかし、心臓の変形率(strain rate)が大きい拡張早期(e’)の時相や、収縮ピーク(s’)の時相では、フレーム間やボリューム間でのパターン変化の程度が大きい。このため、STでは、フレームレートやボリュームレートが不足すると移動ベクトルの正確な推定が困難となる。特に、2次元画像の場合には、追跡する輪郭が走査断面を通過する動き(through-plane)の影響を受けるので、フレーム間でのパターン変化の程度は、更に大きくなり、移動ベクトルの正確な推定がより困難となる。
更に、ノイズやアーティファクトにより動画像の画質が悪い場合、不要な成分が混入している箇所では、移動ベクトルの正確な推定が困難となる。上記の様々な要因により、動き情報が正確に推定できない場合には、追跡外れが生じ、その結果、正確な壁運動解析が行えなくなる。
このため、正確な追跡を行なうための様々な方法が提案されているが、これらの方法は、上記の様々な要因の全てに対応した方法ではなかった。このように、現状のSTでは、輪郭の追跡結果を正確に得ることができない場合があった。また、上記の課題は、超音波画像以外の医用画像の動画像を用いて関心領域の追跡を行なう場合でも、同様に生じる。
本発明が解決しようとする課題は、関心領域の追跡結果を正確に得ることができる超音波診断装置、画像処理装置及び画像処理方法を提供することである。
実施形態の超音波診断装置は、設定部と、推定部と、追跡部とを備える。設定部は、運動する組織を含む領域を撮影対象として収集された2次元又は3次元の超音波画像データの動画像データの中で、少なくとも1つの時相の超音波画像データにて2次元又は3次元の関心領域を設定する。推定部は、前記動画像データが収集された区間のうち前記少なくとも1つの時相以外の残余時相の各超音波画像データで、動き情報に基づいて前記関心領域として推定した領域の第1位置情報と、前記動き情報以外の情報に基づいて前記関心領域として推定した領域の第2位置情報とを取得する。追跡部は、前記動き情報の信頼度に関わる信頼度指標に基づいて前記第1位置情報と前記第2位置情報とを合成した位置情報を前記関心領域の位置情報として取得して、前記関心領域を追跡する。
以下、添付図面を参照して、超音波診断装置の実施形態を詳細に説明する。
(第1の実施形態)
まず、第1の実施形態に係る超音波診断装置の構成について説明する。図1は、第1の実施形態に係る超音波診断装置の構成例を示すブロック図である。図1に例示するように、第1の実施形態に係る超音波診断装置は、超音波プローブ1と、モニタ2と、入力装置3と、心電計4と、装置本体10とを有する。
まず、第1の実施形態に係る超音波診断装置の構成について説明する。図1は、第1の実施形態に係る超音波診断装置の構成例を示すブロック図である。図1に例示するように、第1の実施形態に係る超音波診断装置は、超音波プローブ1と、モニタ2と、入力装置3と、心電計4と、装置本体10とを有する。
超音波プローブ1は、複数の圧電振動子を有し、これら複数の圧電振動子は、後述する装置本体10が有する送受信部11から供給される駆動信号に基づき超音波を発生する。また、超音波プローブ1は、被検体Pからの反射波を受信して電気信号に変換する。また、超音波プローブ1は、圧電振動子に設けられる整合層と、圧電振動子から後方への超音波の伝播を防止するバッキング材等を有する。なお、超音波プローブ1は、装置本体10と着脱自在に接続される。
超音波プローブ1から被検体Pに超音波が送信されると、送信された超音波は、被検体Pの体内組織における音響インピーダンスの不連続面で次々と反射され、反射波信号として超音波プローブ1が有する複数の圧電振動子にて受信される。受信される反射波信号の振幅は、超音波が反射される不連続面における音響インピーダンスの差に依存する。なお、送信された超音波パルスが、移動している血流や心臓壁等の表面で反射された場合の反射波信号は、ドプラ効果により、移動体の超音波送信方向に対する速度成分に依存して、周波数偏移を受ける。
例えば、本実施形態では、被検体Pの2次元走査用に、複数の圧電振動子が一列で配置された1Dアレイプローブが超音波プローブ1として装置本体10と接続される。例えば、超音波プローブ1としての1Dアレイプローブは、セクタ走査を行なうセクタプローブや、オフセットセクタ走査を行なうコンベックスプローブ、リニア走査を行なうリニアプローブ等である。或いは、例えば、本実施形態では、被検体Pの3次元走査用に、メカニカル4Dプローブや2Dアレイプローブが超音波プローブ1として装置本体10と接続されても良い。メカニカル4Dプローブは、1Dアレイプローブのように一列で配列された複数の圧電振動子を用いて2次元走査が可能であるとともに、複数の圧電振動子を所定の角度(揺動角度)で揺動させることで3次元走査が可能である。また、2Dアレイプローブは、マトリックス状に配置された複数の圧電振動子により3次元走査が可能であるとともに、超音波を集束して送信することで2次元走査が可能である。
入力装置3は、マウス、キーボード、ボタン、パネルスイッチ、タッチコマンドスクリーン、フットスイッチ、トラックボール、ジョイスティック等を有し、超音波診断装置の操作者からの各種設定要求を受け付け、装置本体10に対して受け付けた各種設定要求を転送する。なお、第1の実施形態に係る入力装置3が操作者から受け付ける設定情報については、後に詳述する。
モニタ2は、超音波診断装置の操作者が入力装置3を用いて各種設定要求を入力するためのGUI(Graphical User Interface)を表示したり、装置本体10において生成された超音波画像データ等を表示したりする。
心電計4は、超音波走査される被検体Pの生体信号として、被検体Pの心電波形(ECG: Electrocardiogram)を取得する。心電計4は、取得した心電波形を装置本体10に送信する。
装置本体10は、超音波プローブ1が受信した反射波信号に基づいて超音波画像データを生成する装置である。図1に示す装置本体10は、超音波プローブ1が受信した2次元の反射波データに基づいて2次元の超音波画像データを生成可能な装置である。また、図1に示す装置本体10は、超音波プローブ1が受信した3次元の反射波データに基づいて3次元の超音波画像データを生成可能な装置である。以下、3次元の超音波画像データを「ボリュームデータ」と記載する場合がある。
装置本体10は、図1に示すように、送受信部11と、Bモード処理部12と、ドプラ処理部13と、画像生成部14と、画像メモリ15と、内部記憶部16と、画像処理部17と、制御部18とを有する。
送受信部11は、パルス発生器、送信遅延部、パルサ等を有し、超音波プローブ1に駆動信号を供給する。パルス発生器は、所定のレート周波数で、送信超音波を形成するためのレートパルスを繰り返し発生する。また、送信遅延部は、超音波プローブ1から発生される超音波をビーム状に集束し、かつ送信指向性を決定するために必要な圧電振動子ごとの遅延時間を、パルス発生器が発生する各レートパルスに対し与える。また、パルサは、レートパルスに基づくタイミングで、超音波プローブ1に駆動信号(駆動パルス)を印加する。すなわち、送信遅延部は、各レートパルスに対し与える遅延時間を変化させることで、圧電振動子面から送信される超音波の送信方向を任意に調整する。
なお、送受信部11は、後述する制御部18の指示に基づいて、所定のスキャンシーケンスを実行するために、送信周波数、送信駆動電圧等を瞬時に変更可能な機能を有している。特に、送信駆動電圧の変更は、瞬間にその値を切り替え可能なリニアアンプ型の発信回路、又は、複数の電源ユニットを電気的に切り替える機構によって実現される。
また、送受信部11は、プリアンプ、A/D(Analog/Digital)変換器、受信遅延部、加算器等を有し、超音波プローブ1が受信した反射波信号に対して各種処理を行って反射波データを生成する。プリアンプは、反射波信号をチャネル毎に増幅する。A/D変換器は、増幅された反射波信号をA/D変換する。受信遅延部は、受信指向性を決定するために必要な遅延時間を与える。加算器は、受信遅延部によって処理された反射波信号の加算処理を行なって反射波データを生成する。加算器の加算処理により、反射波信号の受信指向性に応じた方向からの反射成分が強調され、受信指向性と送信指向性とにより超音波送受信の総合的なビームが形成される。
送受信部11は、被検体Pを2次元走査する場合、超音波プローブ1から2次元の超音波ビームを送信させる。そして、送受信部11は、超音波プローブ1が受信した2次元の反射波信号から2次元の反射波データを生成する。また、送受信部11は、被検体Pを3次元走査する場合、超音波プローブ1から3次元の超音波ビームを送信させる。そして、送受信部11は、超音波プローブ1が受信した3次元の反射波信号から3次元の反射波データを生成する。
なお、送受信部11からの出力信号の形態は、RF(Radio Frequency)信号と呼ばれる位相情報が含まれる信号である場合や、包絡線検波処理後の振幅情報である場合等、種々の形態が選択可能である。
Bモード処理部12は、送受信部11から反射波データを受信し、対数増幅、包絡線検波処理等を行なって、信号強度が輝度の明るさで表現されるデータ(Bモードデータ)を生成する。
ドプラ処理部13は、送受信部11から受信した反射波データから速度情報を周波数解析し、ドプラ効果による血流や組織、造影剤エコー成分を抽出し、速度、分散、パワー等の移動体情報を多点について抽出したデータ(ドプラデータ)を生成する。
なお、第1の実施形態に係るBモード処理部12及びドプラ処理部13は、2次元の反射波データ及び3次元の反射波データの両方について処理可能である。すなわち、Bモード処理部12は、2次元の反射波データから2次元のBモードデータを生成し、3次元の反射波データから3次元のBモードデータを生成する。また、ドプラ処理部13は、2次元の反射波データから2次元のドプラデータを生成し、3次元の反射波データから3次元のドプラデータを生成する。
画像生成部14は、Bモード処理部12及びドプラ処理部13が生成したデータから超音波画像データを生成する。すなわち、画像生成部14は、Bモード処理部12が生成した2次元のBモードデータから、反射波の強度を輝度で表した2次元Bモード画像データを生成する。また、画像生成部14は、ドプラ処理部13が生成した2次元のドプラデータから、移動体情報を表す2次元ドプラ画像データを生成する。2次元ドプラ画像データは、速度画像データ、分散画像データ、パワー画像データ、又は、これらを組み合わせた画像データである。また、画像生成部14は、ドプラ処理部13が生成したドプラデータから、血流や組織の速度情報を時系列に沿ってプロットしたドプラ波形を生成することも可能である。
ここで、画像生成部14は、一般的には、超音波走査の走査線信号列を、テレビ等に代表されるビデオフォーマットの走査線信号列に変換(スキャンコンバート)し、表示用の超音波画像データを生成する。具体的には、画像生成部14は、超音波プローブ1による超音波の走査形態に応じて座標変換を行なうことで、表示用の超音波画像データを生成する。また、画像生成部14は、スキャンコンバート以外に種々の画像処理として、例えば、スキャンコンバート後の複数の画像フレームを用いて、輝度の平均値画像を再生成する画像処理(平滑化処理)や、画像内で微分フィルタを用いる画像処理(エッジ強調処理)等を行なう。また、画像生成部14は、超音波画像データに、種々のパラメータの文字情報、目盛り、ボディーマーク等を合成する。
すなわち、Bモードデータ及びドプラデータは、スキャンコンバート処理前の超音波画像データであり、画像生成部14が生成するデータは、スキャンコンバート処理後の表示用の超音波画像データである。なお、Bモードデータ及びドプラデータは、生データ(Raw Data)とも呼ばれる。
更に、画像生成部14は、Bモード処理部12が生成した3次元のBモードデータに対して座標変換を行なうことで、3次元Bモード画像データを生成する。また、画像生成部14は、ドプラ処理部13が生成した3次元のドプラデータに対して座標変換を行なうことで、3次元ドプラ画像データを生成する。すなわち、画像生成部14は、「3次元のBモード画像データや3次元ドプラ画像データ」を「3次元超音波画像データ(ボリュームデータ)」として生成する。
更に、画像生成部14は、ボリュームデータをモニタ2にて表示するための各種の2次元画像データを生成するために、ボリュームデータに対してレンダリング処理を行なう。画像生成部14が行なうレンダリング処理としては、断面再構成法(MPR:Multi Planer Reconstruction)を行なってボリュームデータからMPR画像データを生成する処理がある。また、画像生成部14が行なうレンダリング処理としては、ボリュームデータに対して「Curved MPR」を行なう処理や、ボリュームデータに対して「Maximum Intensity Projection」を行なう処理がある。また、画像生成部14が行なうレンダリング処理としては、ボリュームレンダリング(VR:Volume Rendering)処理がある。
画像メモリ15は、画像生成部14が生成した表示用の画像データを記憶するメモリである。また、画像メモリ15は、Bモード処理部12やドプラ処理部13が生成したデータを記憶することも可能である。画像メモリ15が記憶するBモードデータやドプラデータは、例えば、診断の後に操作者が呼び出すことが可能となっており、画像生成部14を経由して表示用の超音波画像データとなる。
なお、画像生成部14は、超音波画像データと当該超音波画像データを生成するために行なわれた超音波走査の時間とを、心電計4から送信された心電波形に対応付けて画像メモリ15に格納する。後述する画像処理部17や制御部18は、画像メモリ15に格納されたデータを参照することで、超音波画像データを生成するために行なわれた超音波走査時の心時相を取得することができる。
内部記憶部16は、超音波送受信、画像処理及び表示処理を行なうための制御プログラムや、診断情報(例えば、患者ID、医師の所見等)や、診断プロトコルや各種ボディーマーク等の各種データを記憶する。また、内部記憶部16は、必要に応じて、画像メモリ15が記憶する画像データの保管等にも使用される。また、内部記憶部16が記憶するデータは、図示しないインターフェースを経由して、外部装置へ転送することができる。なお、外部装置は、例えば、画像処理用の高性能なワークステーションや、画像診断を行なう医師が使用するPC(Personal Computer)、CDやDVD等の記憶媒体、プリンター等である。
画像処理部17は、コンピュータ支援診断(Computer-Aided Diagnosis:CAD)を行なうために、装置本体10に設置される。画像処理部17は、画像メモリ15に格納された超音波画像データを取得して、診断支援のための画像処理を行なう。そして、画像処理部17は、画像処理結果を、画像メモリ15や内部記憶部16に格納する。
具体的には、第1の実施形態に係る画像処理部17は、運動する組織の運動情報を提供するために設置される。上記の運動する組織は、例えば心臓のように、周期的に運動する組織である。第1の実施形態に係る画像処理部17は、図1に例示するように、取得部171と、設定部172と、追跡部173と、運動情報演算部174とを有する。なお、画像処理部17が行なう処理については、後に詳述する。
制御部18は、超音波診断装置の処理全体を制御する。具体的には、制御部18は、入力装置3を介して操作者から入力された各種設定要求や、内部記憶部16から読込んだ各種制御プログラム及び各種データに基づき、送受信部11、Bモード処理部12、ドプラ処理部13、画像生成部14及び画像処理部17の処理を制御する。また、制御部18は、画像メモリ15や内部記憶部16が記憶する表示用の超音波画像データをモニタ2にて表示するように制御する。また、制御部18は、画像処理部17の処理結果をモニタ2に表示するように制御したり、外部装置に出力したりするように制御する。
以上、第1の実施形態に係る超音波診断装置の全体構成について説明した。かかる構成のもと、第1の実施形態に係る超音波診断装置が有する画像処理部17は、運動する組織(例えば、周期的に運動する組織)の運動情報を提供するために、時系列に沿った超音波データ群(超音波画像データの動画像データ)にて、関心領域の追跡を行なう。一例として、第1の実施形態に係る画像処理部17は、周期的に運動する心壁の運動情報を提供するために、超音波画像データの動画像データにて、心筋の内膜の輪郭や心筋の外膜の輪郭の追跡を行なう。
なお、画像処理部17の処理対象となる動画像データは、2次元超音波画像データ群であっても、3次元超音波画像データ群であっても良い。
まず、図1に例示する取得部171は、運動する組織を含む領域を撮影対象として収集された2次元又は3次元の超音波画像データの動画像データを取得する。具体的には、第1の実施形態では、動画像データは、運動する組織である周期的に運動する組織を撮影対象として収集された少なくとも1周期以上の区間に渡る2次元又は3次元の超音波画像データである。すなわち、取得部171は、周期的に運動する組織を含む領域を撮影対象として収集された少なくとも1周期以上の区間に渡る2次元又は3次元の超音波画像データの動画像データを取得する。例えば、操作者は、セクタプローブにより、被検体Pの心臓を含む領域を2次元走査、又は、3次元走査を行なって、心筋が描出された2次元又は3次元の超音波画像データの動画像データの撮影を行なう。動画像データは、例えば、Bモードで収集された超音波画像データ群である。これにより、画像生成部14は、心筋の動画像データを生成し、生成した動画像データを画像メモリ15に格納する。
以下では、一例として、操作者が、乳頭筋レベルの左室短軸断面の2次元走査を行なうことで、複数心拍に渡る2次元超音波画像データの動画像データが収集される場合について説明する。乳頭筋レベルの左室短軸断面の超音波画像データには、右心室(RV:Right Ventricular)の短軸面の一部と、左心室(LV:Left Ventricular)の短軸面の全体とが描出される。また、乳頭筋レベルの左室短軸断面の超音波画像データには、左心室内の前乳頭筋(APM : Anterolateral Papillary Muscle)と後乳頭筋(PPM : Posteromedial Papillary Muscle)とが描出される。
そして、操作者は、追跡対象の区間として、例えば、最初の拡張末期(ED:End Diastole)から次の拡張末期までの1心拍区間を設定する。以下、最初の拡張末期の時相を「ED0」と記載し、次の拡張末期の時相を「ED1」と記載する。「ED0」は、追跡開始時相となり、「ED1」は、追跡終了時相となる。
取得部171は、操作者が設定した1心拍区間(t=ED0~ED1)の動画像データを画像メモリ15から取得する。図2は、取得部が取得する動画像データの一例を示す図である。例えば、取得部171は、図2に示すように、1心拍区間(t=ED0~ED1)の2次元動画像データとして、「時相:t1=E0」の画像データ「I(t1)」~「時相:tn=E1」の画像データ「I(tn)」を画像メモリ15から取得する。なお、本実施形態は、取得部171が3次元超音波画像データの動画像データを取得する場合であっても適用可能である。また、本実施形態は、追跡対象の区間が2心拍区間や3心拍区間として設定される場合であっても適用可能である。
図1に例示する設定部172は、取得部171が取得した動画像データの中で、少なくとも1つの時相の超音波画像データにて2次元又は3次元の関心領域を設定する。動画像データが2次元超音波画像データ群である場合、設定部172は、2次元の関心領域を設定する。また、動画像データが3次元超音波画像データ群である場合、設定部172は、2次元又は3次元の関心領域を設定する。
また、周期的に運動する心壁の運動情報を操作者に提供する本実施形態では、関心領域は、左室内膜の輪郭や心筋外膜の輪郭となる。本実施形態では、設定部172は、操作者が用手的に設定した情報から、関心領域を設定する。図3は、操作者により設定される初期輪郭の一例を示す図である。
例えば、操作者は、入力装置3を用いて、追跡開始時相「t1=E0」の画像データ「I(t1)」の表示要求を行なう。モニタ2は、制御部18の制御により、画像データ「I(t1)」を表示する。操作者は、図3の上図に示す画像データ「I(t1)」を参照して、左室内膜の輪郭の位置をトレースする。また、操作者は、図3の上図に示す画像データ「I(t1)」を参照して、心筋外膜の輪郭の位置をトレースする。これにより、設定部172は、図3の下図に示すように、操作者がトレースした左室内膜の輪郭の位置を、左室内膜の初期輪郭位置「C(t1)」として設定する。また、設定部172は、図3の下図に示すように、操作者がトレースした心筋外膜の輪郭の位置を、心筋外膜の初期輪郭位置「D(t1)」として設定する。図3では、操作者がトレースした輪郭線を、太い曲線で示している。ここで、心筋外膜の輪郭位置は、心内膜位置から所定の厚さ(所定の距離)だけ離れた位置として自動的に発生させても良い。
なお、本実施形態は、心筋の中層の輪郭の位置が初期輪郭位置として設定される場合であっても良い。また、本実施形態は、操作者により用手的に設定された初期輪郭位置を用いる場合に限定されるものではなく、設定部172が超音波画像データの輝度等に基づいて、初期輪郭位置を自動的に設定する場合であっても良い。
そして、図1に例示する追跡部173は、動画像データが収集された区間のうち、初期関心領域が設定された時相以外の残余時相の各超音波画像データで関心領域の追跡を行なう。具体的には、追跡部173は、追跡対象の区間のうち初期関心領域が設定された時相以外の残余時相の各超音波画像データで、関心領域の追跡を行なう。例えば、追跡部173は、設定部172が追跡開始時相「t1=E0」で設定した初期輪郭位置を、追跡開始時相「t1=E0」以外の残余時相『「t2」~「tn=E1」』それぞれの超音波画像データ「I(t2)~I(tn)」にて追跡する。
従来では、初期関心領域の位置(初期輪郭位置)の追跡は、スペックルトラッキングにより行なわれていた。すなわち、従来では、局所的なパターンマッチング処理や、オプティカルフロー法で得られる動き情報(移動ベクトル)を用いて、初期輪郭位置を自動的に追跡することで、残余時相「t2~tn=E1」の各超音波画像データの関心領域の位置を取得していた。しかし、心臓の変形率(strain rate)が大きい拡張早期時相や、収縮ピーク時相では、フレーム間やボリューム間でのパターン変化の程度が大きい。
このため、スペックルトラッキングでは、フレームレートやボリュームレートが不足すると移動ベクトルの正確な推定が困難となる。特に、2次元画像データの場合には、追跡する輪郭が走査断面を通過する動き(through-plane)の影響を受けるので、フレーム間でのパターン変化の程度は、更に大きくなり、移動ベクトルの正確な推定がより困難となる。
更に、ノイズやアーティファクトにより動画像データの画質が悪い場合、不要な成分が混入している箇所では、移動ベクトルの正確な推定が困難となる。上記の様々な要因により、動き情報が正確に推定できない場合には、追跡外れが生じ、その結果、正確な壁運動解析が行えなくなる。
これに対して、追跡結果により得られた輪郭の位置と初期輪郭位置との間で、動きが矛盾している追跡点を検知し、検知した追跡点の位置を修正する手法が知られている。しかし、この方法は、検知した追跡点の周辺の追跡点の位置が正確であるという前提に基づくものである。しかし、パターン変化の程度が大きくなる上記の時相では、広範囲で動き情報の推定が破綻する。すなわち、拡張早期時相や、収縮ピーク時相では、この前提が満たされず、高精度な追跡を行なえない。
また、画像内でノイズ成分が混入している領域をユーザーが指定し、ユーザーが指定した領域の輪郭位置を、この領域以外で追跡された輪郭の位置情報から与える手法が知られている。しかし、この方法は、ユーザーが用手的に追跡外れ位置での輪郭を修正する必要があり、手間を要する。
また、画像データに基づく3次元輪郭形状とオプティカルフローに基づく3次元輪郭形状との線形和により最終的な3次元の輪郭を決定する方法が知られている。しかし、この方法は、3次元動画像データを用いた輪郭追跡に適用できるが、2次元動画像データを用いた輪郭追跡に適用する場合には、上述した「through-plane」の影響を排除できないため性能が制約される。また、線形和を如何なる条件に基づいて行なうかを明確に示す手法は、現時点で知られていない。
そこで、図1に例示する追跡部173は、関心領域の追跡結果を正確に得るために、以下に説明する追跡処理を行なう。すなわち、図1に例示する追跡部173は、残余時相の各超音波画像データで、動き情報に基づいて関心領域として推定した領域の第1位置情報を取得する。また、図1に例示する追跡部173は、残余時相の各超音波画像データで、動き情報以外の情報(例えば、形状情報)に基づいて関心領域として推定した領域の第2位置情報を取得する。すなわち、追跡部173は、第1位置情報及び第2位置情報を取得する推定部として機能する。そして、図1に例示する追跡部173は、動き情報の信頼度に関わる信頼度指標に基づいて、第1位置情報と、第2位置情報とを合成した位置情報を、関心領域の位置情報として取得する。これにより、図1に例示する追跡部173は、残余時相の各超音波画像データで、関心領域を追跡する。ここで、「信頼度指標」は、動き情報の品質である「動き品質」とも定義される。以下では、「信頼度指標」を「動き品質」と記載する場合がある。
図4は、図1に示す追跡部の構成例を示すブロック図である。上記の追跡処理を行なうため、追跡部173は、例えば、図4に示すように、第1推定部173aと、定義部173bと、第2推定部173cと、組み合わせ部173dとを有する。なお、以下に説明する第1推定部173aと第2推定部173cとは、追跡部173とは別に、画像処理部17において、推定部として設置される場合であっても良い。
第1推定部173aは、動き情報に基づく関心領域の追跡処理、すなわち、スペックルトラッキング処理を行なう。すなわち、第1推定部173aは、第1時相の超音波画像データである第1画像データと、第1時相と時間的に隣接する第2時相の超音波画像データである第2画像データとの間で動き情報を推定し、推定した動き情報に基づいて第1画像データの関心領域を移動して、第2画像データの第1位置情報を推定する。
具体的には、第1推定部173aは、関心領域を含む領域内の複数の点それぞれの移動ベクトルを演算する。より具体的には、第1推定部173aは、テンプレートマッチング処理により、複数の点それぞれの移動ベクトルを演算する。そして、第1推定部173aは、関心領域を構成する各追跡点について、追跡点を含む所定のサイズのセグメント内で得られた移動ベクトル群から統計的に異常な移動ベクトルを除去した残存移動ベクトル群を用いて該追跡点の移動ベクトルを推定する。これにより得られた各追跡点の移動ベクトルを用いて、第1推定部173aは、第1位置情報を推定する。そして、第1推定部173aは、得られた移動ベクトル群から統計的に異常な移動ベクトルを除去した残存移動ベクトル群を用いて、関心領域を構成する各追跡点の移動ベクトルを推定することで、第1位置情報を推定する。
定義部173bは、第1推定部173aの処理で得られる少なくとも1つの変量に基づいて、信頼度指標(動き品質)を定義する。第2推定部173cは、第2画像データの第2位置情報を動き情報以外の情報(例えば、形状情報)により推定する。組み合わせ部173dは、信頼度指標(動き品質)に基づいて第1位置情報と第2位置情報とを組み合わせた位置情報を、第2画像データにおける関心領域の位置情報として取得する。
以下、上述した各部が行なう具体的な処理の一例について、図5A~図16を用いて、詳細に説明する。図5A~図16は、第1の実施形態に係る追跡部を説明するための図である。以下では、図3の下図に示す左室内膜の初期輪郭位置「C(t1)」を残余時相「t2~tn=E1」の各超音波画像データで追跡する場合について説明する。以下に説明する処理により、心筋外膜の初期輪郭位置「D(t1)」についても、同様に、追跡可能である。このため、心筋外膜の輪郭位置の追跡処理については、説明を省略する。
まず、第1推定部173aは、図5Aに示すように、初期輪郭位置「C(t1)」に複数の追跡点を設定する。また、第1推定部173aは、図5Bに示すように、初期輪郭位置「C(t1)」の近傍に、複数の点を設定する。図5Bに例示する複数の点それぞれは、動き情報(移動ベクトル)を演算するための演算点である。
そして、第1推定部173aは、第1時相を「t1」とし、第2時相を「t2」とする。そして、第1推定部173aは、第1画像データである「I(t1)」と第2画像データである「I(t2)」との間で動き情報を推定し、推定した動き情報に基づいてI(t1)の初期輪郭位置「C(t1)」を移動する。これにより、第1推定部173aは、図6に示すように、「I(t2)」の第1位置情報「Cvector(t2)」を推定する。
具体的には、第1推定部173aは、「I(t1)」で設定した複数の演算点それぞれに対して、テンプレートデータを設定する。テンプレートデータは、演算点を中心とする複数のピクセルから構成される。そして、第1推定部173aは、2つのフレーム間、すなわち、「I(t1)」と「I(t2)」との間で、テンプレートデータのスペックルパターンと最も一致する領域を探索することで、テンプレートデータが次のフレームでどの位置に移動したかを追跡する。これにより、第1推定部173aは、「I(t1)」の複数の演算点それぞれの移動ベクトルを演算する。以下、「I(t1)」の複数の演算点それぞれの移動ベクトルの集合を{Mv0}とする。{Mv0}は、複数の演算点それぞれで演算された移動ベクトル群となる。
図7は、第1推定部173aが演算した移動ベクトル群{Mv0}の一例を示している。図7では、白抜きの正方形が演算点を示し、白抜きの正方形から延びる線分が移動ベクトルを示している。なお、図7に例示する移動ベクトル群は、心内腔が拡大する拡張早期の画像データ間で演算された移動ベクトル群の様子を模式的に示したものである。
そして、第1推定部173aは、図5Aに示す各追跡点を、所定のサイズのセグメント単位で分割する。例えば、第1推定部173aは、図8に示すように、アメリカ心エコー図学会やアメリカ心臓協会が推奨する分割領域を用いて、短軸断面の左心室心筋を、「ANT-SEPT(前壁中隔)、ANT(前壁)、LAT(側壁)、POST(後壁)、INF(下壁)、SEPT(中隔)」の6つのセグメントに分ける。これにより、図5Aに示す各追跡点は、図8に例示する6つのセグメントのいずれかに含まれることになる。また、図5Bに示す各演算点は、図8に例示する6つのセグメントのいずれかに含まれることになる。
そして、第1推定部173aは、各セグメントに含まれる移動ベクトル群を統計処理することで、各セグメントの統計値を演算し、各セグメントで統計値が異常値となる移動ベクトルを特定する。そして、第1推定部173aは、特定した移動ベクトルの集合、すなわち、統計値が異常値となる移動ベクトルの集合を{Mv’}とする。そして、第1推定部173aは、{Mv0}から{Mv’}を除去した移動ベクトルの集合{Mv}を残存移動ベクトル群とする。
例えば、第1推定部173aは、各セグメントでの移動ベクトル群の分散値「Vari(s)」を演算する。ここで、図8に例示する6つのセグメントを「s1、・・・・、s6」とすると、第1推定部173aは、6つの分散値「Vari(s1)、・・・・、Vari(s6)」を演算する。そして、第1推定部173aは、例えば、セグメント「s1」に含まれる移動ベクトル群において、「α×Vari(s1)」を超える移動ベクトルを異常移動ベクトルとして特定する。なお、「α」は、所定の係数であり、「α」の値としては、「2」程度が好適な例となる。第1推定部173aは、かかる処理を、セグメント「s2、s3、s4、s5、s6」それぞれで同様に行なう。これにより、第1推定部173aは、図9に示すように、異常移動ベクトル群である{Mv’}を特定し、{Mv0}から{Mv’}を除去することで、残存移動ベクトル群{Mv}を取得する。
図9は、第1推定部173aが特定して異常移動ベクトル群の一例を示している。図9では、第1推定部173aが特定した4つの異常移動ベクトルを点線で示している。なお、図9に例示する異常移動ベクトル群{Mv’}は、図7で例示した移動ベクトル群{Mv0}、すなわち、心内腔が拡大する拡張早期の画像データ間で演算された移動ベクトル群から特定された異常移動ベクトル群である。このため、図9で例示する残存移動ベクトル群{Mv}は、概して外向きの方向へ向かう移動ベクトルとなる。
そして、第1推定部173aは、残存移動ベクトル群{Mv}を用いて、初期輪郭位置「C(t1)」を構成する複数の追跡点それぞれの移動ベクトルを推定する。以下では、{Mv}を構成する各要素の移動ベクトル(残存移動ベクトル)を、「Mv」と一般化して説明する。また、「C(t1)」を構成する複数の追跡点の中で、任意の追跡点を「p(t1)」とする。第1推定部173aは、残存移動ベクトル群{Mv}の中からp(t1)近傍に位置する幾つかの移動ベクトルを平均することで、p(t1)の移動ベクトル「Mv_p(t1)」を取得する。
上記方法の一例について、図10を用いて説明する。図10では、p(t1)を、点でハッチングした丸で示している。また、図10では、p(t1)の近傍に位置する複数の演算点を白丸で示し、各演算点の残存移動ベクトルを白丸から延びる矢印で示している。
例えば、第1推定部173aは、p(t1)を中心とする所定半径の円を設定する(図10に示す点線の丸を参照)。所定半径の長さは、システム上、予め設定された値であっても、操作者により設定された値であっても良い。そして、第1推定部173aは、p(t1)で設定した円内に含まれる複数の残存移動ベクトル(図10では4つの残存移動ベクトル)を平均することで、p(t1)の移動ベクトル「Mv_p(t1)」を推定する。第1推定部173aは、上述した処理を、{Mv}を用いて、全ての追跡点において行なう。
そして、第1推定部173aは、図11の上図に示すように、p(t1)をMv_p(t1)により移動させた点p(t2)の位置を第1位置情報として推定する。これにより、第1推定部173aは、図11の下図に示すように、「I(t2)」の第1位置情報「Cvector(t2)」を取得する。以下では、「Cvector(t2)」を構成する複数の追跡点の中で、任意の追跡点を「p(t2)」とする。
そして、定義部173bは、「第1推定部173aの処理で得られる少なくとも1つの変量」を用いて、「p(t2)」の信頼度指標である動き品質「Q_p(t2)」を定義する。
第1推定部173aは、上述したように、テンプレートマッチングにより移動ベクトルを推定する。すなわち、第1推定部173aは、現時相でのテンプレートにおける画像信号値の濃淡パターンに対して、次の時相で最も似ているパターンが存在する位置を探索することで各演算点での移動ベクトルを検出する。そこで、定義部173bは、第1推定部173aが行なったテンプレートマッチング処理で得られる変量から、残存移動ベクトル群{Mv}を構成する個々の残存移動ベクトル「Mv」の信頼度「R(Mv)」を定義する。信頼度の定義に用いる変量は、以下に説明する3つの変量が挙げられる。
信頼度の定義に用いる第1の変量は、テンプレートマッチング処理における基準テンプレートの平均信号値である。すなわち、テンプレートの信号値が極端に小さい場合には、このテンプレートは、マッチングに有効な組織領域ではなく、心腔内やホワイトノイズの領域と考えられる。従って、テンプレート内の信号の平均値「Tm」が小さい時には「R(Mv)」の値を小さくするのが妥当となる。
そこで、定義部173bは、一例として、以下に示す式(1)により、「R(Mv)」を定義する。なお、以下の式(1)で示す「Mth」は、基準テンプレート内の信号平均値に対して設定された下限閾値である。「Mth」は、例えば、内部記憶部16に予め格納される。或いは、「Mth」は、例えば、操作者により設定される。
すなわち、定義部173bは、「Tm」が「Mth」以上の基準テンプレートを用いて演算された残存移動ベクトルについては、信頼度を「1」と定義し、「Tm」が「Mth」より小さい基準テンプレートを用いて演算された残存移動ベクトルについては、信頼度を「0」と定義する。
信頼度の定義に用いる第2の変量は、テンプレートマッチング処理における基準テンプレートの信号分散値である。すなわち、基準テンプレートの信号値の濃淡パターンが一様過ぎる場合には、次の時相の探索範囲における濃淡パターンも一様に近いと考えられる。かかる場合、最も似ている位置が、第2画像データにて至る所で検出される可能性が高いため、正確な移動ベクトルの推定ができなくなる。従って、基準テンプレート内の信号の分散値「Tσ」が小さい時には「R(Mv)」の値を小さくするのが妥当となる。
そこで、定義部173bは、一例として、以下に示す式(2)により、「R(Mv)」を定義する。なお、以下の式(2)で示す「Sth」は、基準テンプレート内の信号分散値に対して設定された下限閾値である。「Sth」は、例えば、内部記憶部16に予め格納される。或いは、「Sth」は、例えば、操作者により設定される。
すなわち、定義部173bは、「Tσ」が「Sth」以上の基準テンプレートを用いて演算された残存移動ベクトルについては、信頼度を「1」と定義し、「Tσ」が「Sth」より小さい基準テンプレートを用いて演算された残存移動ベクトルについては、信頼度を「0」と定義する。
信頼度の定義に用いる第3の変量は、テンプレートマッチング処理におけるテンプレート間の類似度である。第1推定部173aは、基準テンプレートと次の時相の探索領域との比較により、第2画像データで基準テンプレートと最も似ている位置を探索する。従って、最終的に探索した位置の領域が、基準テンプレートとどの程度似ているかを示す類似度を定義することができる。代表的な類似度としては、相互相関係数「Ck」があり、入力信号が実数の場合、「Ck」は、「0≦Ck≦1」となる。画像データ間の比較を行なう本実施形態では、入力信号が実数となり、「0≦Ck≦1」となる。類似度が低くなることは、現時相から次の時相までの時間経過において、追跡対象部位のパターン変化が大きいことを意味する。従って、類似度である「Ck」が小さい時には「R(Mv)」の値を小さくするのが妥当となる。
そこで、定義部173bは、一例として、以下に示す式(3)により、「R(Mv)」を定義する。なお、以下の式(3)で示す「Cth」は、テンプレートマッチングでの類似度に対して設定された下限閾値である。「Cth」は、例えば、内部記憶部16に予め格納される。或いは、「Cth」は、例えば、操作者により設定される。
すなわち、定義部173bは、「Ck」が「Cth」以上であった残存移動ベクトルについては、信頼度を「Ck」と定義し、「Ck」が「Cth」より小さかった残存移動ベクトルについては、信頼度を「0」と定義する。
なお、上記の式(1)及び式(2)では、単純に、閾値により「1」又は「0」で「R(Mv)」を与えている。しかし、本実施形態は、「R(Mv)」を、上記の制御変数を用いた連続的な単調増加の関数として定義する場合であっても良い。また、「R(Mv)」は、上記の3つの制御変数のうち、少なくとも1つ以上を用いた組合せで定義することも可能である。本実施形態では、上記の3つの制御変数を全て用いて、最終的な「R(Mv)」を与えることが好適である。かかる場合、「R(Mv)」は、「Tm、Tσ、Ck」の3つを制御変数とする関数により定義される。
定義部173bは、個々の残存移動ベクトルの「R(Mv)」を上記のように定義することで、個々の残存移動ベクトルの品質評価(信頼性の評価)を行なう。更に、定義部173bは、以下に説明するように、複数の移動ベクトルの空間的な連続性を用いて、所定のサイズを有する領域(セグメント)における移動ベクトルの品質評価(信頼性の評価)を行なう。
ここで、生体組織の動きでは、極端な空間的不連続は起こりにくい。一方、スペックルトラッキングでは、スペックルノイズの影響により不当な移動ベクトルが検出されることがあり得る。このため、第1推定部173aは、上述したように、各セグメントでの移動ベクトル群の分散値「Vari(s)」を演算し、{Mv0}から所定の閾値を越える異常値を有する移動ベクトル群{Mv’}を異常移動ベクトル群として取り除き、残存移動ベクトル群{Mv}を得ている。
そこで、定義部173bは、第1推定部173aが行なった統計処理で得られる変量を、動き品質(信頼度指標)を定義するための変量とする。
統計処理で得られる第1の変量としては、複数の追跡点それぞれの近傍領域内での移動ベクトルの分散値がある。ここで、近傍領域内における移動ベクトルの空間的分散の値が小さいほど、動きの品質は低いと考えられる。近傍領域として、定義部173bは、例えば、上記の6つのセグメント「s1、・・・・、s6」を用いる。かかる場合、定義部173bは、第1推定部173aから、各セグメントでの移動ベクトル群の分散値「Vari(s)」を得る。すなわち、定義部173bは、6つのセグメント「s1、・・・・、s6」それぞれの「Vari(s1)、・・・・、Vari(s6)」を取得する。
統計処理で得られる第2の変量としては、複数の追跡点それぞれの近傍領域において残存移動ベクトル群として採用された移動ベクトルの密度がある。ここで、近傍領域内における有効な移動ベクトル(残存移動ベクトル)の密度の値が小さいほど、動きの品質は低いと考えられる。近傍領域として、定義部173bは、例えば、上記の6つのセグメント「s1、・・・・、s6」を用いる。かかる場合、定義部173bは、例えば、第1推定部173aから、各セグメント「s」内での{Mv0}の数と、各セグメント「s」内での{Mv}の数とを取得し、各セグメント「s」における有効移動ベクトルの密度「ρ(s)」を演算する。定義部173bは、一例として、以下に示す式(4)により、「ρ(s)」を定義する。
すなわち、定義部173bは、「s」内での{Mv0}の数を、セグメント「s」内での{Mv}の数で除算した値を、「ρ(s)」として定義する。定義部173bは、式(4)を用いて、6つのセグメント「s1、・・・・、s6」それぞれの「ρ(s1)、・・・・、ρ(s6)」を取得する。
そして、定義部173bは、上述した「R(Mv)」、「ρ(s)」及び「Vari(s)」の少なくとも1つ以上を用いて、第1推定部173aが演算した各演算点の移動ベクトルの動き品質「Q(s,Mv)」を演算する。例えば、定義部173bは、「Q(s,Mv)」を以下の式(5)により定義する。なお、以下の式(5)で示す「Lth」は、移動ベクトルの空間的分散値に対して設定された下限閾値である。「Lth」は、例えば、内部記憶部16に予め格納される。或いは、「Lth」は、例えば、操作者により設定される。
すなわち、定義部173bは、「Vari(s)」が「Lth」より小さいセグメント内に含まれる演算点の移動ベクトルについては、動き品質を「0」とする。例えば、定義部173bは、「Vari(s2)」が「Lth」より小さい場合、「s2」に含まれる演算点の移動ベクトルについては、動き品質を「0」とする。
また、定義部173bは、「Vari(s)」が「Lth」以上であったセグメント内に含まれる演算点の移動ベクトルについては、個々の移動ベクトルについて、動き品質の演算を行なう。例えば、定義部173bは、「Vari(s5)」が「Lth」以上である場合、「s5」に含まれる各演算点の移動ベクトル「Mv」の信頼度「R(Mv)」に、「ρ(s5)」を乗算することで、「Q(s5,Mv)」を演算する。
これにより、定義部173bは、「Cvector(t2)」の推定に用いられた各残存移動ベクトルの動き品質(信頼度指標)を定義する。なお、定義部173bは、「R(Mv)」を用いて、各残存移動ベクトルの動き品質を定義しても良い。また、定義部173bは、「Vari(s)」及び「ρ(s)」の少なくとも1つを用いて、各残存移動ベクトルの動き品質を定義しても良い。
そして、定義部173bは、p(t1)近傍に位置する幾つかの残存移動ベクトルの動き品質を平均することで、「p(t1)から移動したp(t2)」の信頼度指標である動き品質「Q_p(t2)」を定義する。例えば、定義部173bは、図12に示すように、「Mv_p(t1)」の演算に用いた4つの移動ベクトル(残存移動ベクトル)それぞれの動き品質の平均値を、「Q_p(t2)」として演算する。
なお、第1推定部173aの処理と、定義部173bとの処理は、追跡対象となる時相ごとに同期させて、並列で行なわれることが好適である。かかる並列処理は、制御部18の制御により実現可能である。また、上記の処理で用いられたセグメントの数及びサイズは、例えば、操作者により任意の数及びサイズに変更可能である。また、第1推定部173aの処理で用いられるセグメントと、定義部173bの処理で用いられるセグメントとは、異なる場合であっても良い。
そして、第2推定部173cは、以下に説明する第1処理で得られる位置情報、第2処理で得られる位置情報及び第3処理で得られる位置情報の少なくとも1つの位置情報を用いて、第2位置情報を推定する。
第1処理では、第2推定部173cは、時相「t」の第2画像データと同時相の形状辞書情報「dict(t)」と照合することで、第2画像データにおける関心領域の位置情報「Cdict(t)」を推定する。第1処理を行なう場合、内部記憶部16は、例えば、全ての心時相における左室心筋の形状辞書を予め記憶する。そして、第2推定部173cは、第2画像データが「I(t2)」であることから、時相「t2」に対応する心時相の形状辞書「dict(t2)」を内部記憶部16から取得する。そして、第2推定部173cは、図13Aに示すように、「I(t2)」と「dict(t2)」とを識別器の機能により照合して、「I(t2)」の関心領域(輪郭)の位置情報「Cdict(t2)」を推定する。識別器としては、形状辞書の学習時や識別時の両方でランダム木を用いるアルゴリズムが高速な手法として知られており、この手法が本実施形態における好適な例である。
第2処理では、第2推定部173cは、形状エネルギー最小化原理に基づいて、時相「t」の第2画像データにおける関心領域の位置情報「Ceng(t)」を推定する。一例として、第2推定部173cは、動的輪郭モデル(ACM:Active Contour Model)の一種として知られているSNAKE法を用いる。第2推定部173cは、処理対象となる第2画像データ内に、心内腔の境界近傍の輪郭線を設定する。設定する輪郭線は、例えば、C(t1)に基づいて設定される場合であっても、輝度情報に基づいて設定される場合であっても良い。
そして、第2推定部173cは、形状エネルギーの評価関数の値が最小となるまで、輪郭線の位置及び形状の修正を繰り返す。例えば、第2推定部173cは、輝度が急激に変化する箇所に輪郭線が位置し、輪郭線内の輝度値が均一となり、輪郭線が短く、輪郭線の形状が滑らかとなるまで、輪郭線の位置及び形状の修正を繰り返す。これにより、第2推定部173cは、図13Bに示すように、第2画像データである「I(t2)」の関心領域(輪郭)の位置情報「Ceng(t2)」を推定する。
第3処理では、第2推定部173cは、最小2乗法に基づくフィッティングを含む処理により、時相「t」の第2画像データにおける関心領域の位置情報「Cfit(t)」を推定する。ただし、最小2乗法に基づくフィッティング処理で、信頼性の高い位置情報を得るためには、フィッティングの入力として、ある程度正確な位置情報を含む複数の制御点「r(t)」が必要となる。このため、第2推定部173cは、動き品質(信頼度指標)が低下する部位のみに対して第3処理を行なう。
上述したように、第2画像データである「I(t2)」の「Cvector(t2)」を構成する各追跡点p(t2)では、動き品質「Q_p(t2)」が定義されている。そこで、第2推定部173cは、「Q_p(t2)」が所定の値以下となるp(t2)を特定する。第2推定部173cは、例えば、図13Cの左図に示すように、点線の円内に含まれるp(t2)を特定する。そして、第2推定部173cは、点線の円内に含まれるp(t2)以外のp(t2)の集合を、フィッティング処理を行なうための制御点群「r(t2)」とし、制御点群「r(t2)」を用いて最小2乗フィッティングを実施する。これにより、第2推定部173cは、例えば、図13Cの左図に示すように、第2画像データである「I(t2)」の関心領域(輪郭)の位置情報「Cfit(t2)」を推定する。
第3処理では、Cvector(t)上で局所的なノイズの影響で動き品質が低下するような場合において、動き品質が低下する部位の輪郭位置のみを、周辺の動き品質の高いCvector(t)の位置により適切に補間することで、「Cfit(t)」を得る。換言すると、第3処理は、動き品質(信頼度指標)が低下する部位の輪郭形状を、周辺の動き品質(信頼度指標)の高い輪郭形状に基づいて、推定する処理となる。なお、第2推定部173cは、Cvector(t)上で動き品質が低下する局所的な位置においてのみ、第1処理を行なって「Cdict(t)」を得る場合であっても、第2処理を行なって「Ceng(t)」を得る場合であっても良い。
そして、第2推定部173cは、「Cdict(t)」、「Ceng(t)」及び「Cfit(t)」の少なくとも1つを用いて、第2画像データの第2位置情報「Cmodel(t)」を定義する。第2画像データが「I(t2)」であることから、第2推定部173cは、「Cdict(t2)」、「Ceng(t2)」及び「Cfit(t2)」の少なくとも1つを用いて、「I(t2)」の第2位置情報「Cmodel(t2)」を定義する。
第1処理のみを行なう場合、第2推定部173cは、「Cdict(t2)」を「Cmodel(t2)」とする。また、第2処理のみを行なう場合、第2推定部173cは、「Ceng(t2)」を「Cmodel(t2)」とする。また、第3処理のみを行なう場合、第2推定部173cは、「Cfit(t2)」を「Cmodel(t2)」とする。
ここで、第2推定部173cは、「Cvector(t)」及び「Q_p(t)」が得られていれば、「Cfit(t)」を推定可能である。また、第2推定部173cは、「Cvector(t)」及び「Q_p(t)」が得られていない状態でも、「Cdict(t)」の推定や、「Ceng(t)」の推定を行なうことが可能である。更に、第2推定部173cは、「Cdict(t)」の推定と、「Ceng(t)」の推定とは、各々独立に行なうことが可能である。
すなわち、第2推定部173cは、これら3つの位置情報を推定することが可能である。そこで、第2推定部173cは、これら3つの位置情報の全て、又は、これら3つの位置情報から選択した2つの位置情報を線形加算することで、「Cmodel(t)」を与えても良い。例えば、第2推定部173cは、例えば、複数の輪郭位置を平均することで、「Cmodel(t)」を与える。
ただし、「Cdict(t)」を得るためには、形状辞書を全ての時相で用意する必要がある。臨床での心臓の形状には様々なバリエーションがあるため、高精度な形状辞書を構成するためには一時相だけであっても一般的に数百例以上のデータが必要であり、心周期に渡る全時相での形状辞書を用いるのは実際的でない。
そこで、第2推定部173cは、形状辞書情報を有している特定時相では、第1処理で得られる位置情報を第2位置情報とする。一例として、第2推定部173cは、組織(心筋)の形状の局所変形の程度が大きい特定時相では、第1処理で得られる位置情報を第2位置情報とする。また、第2推定部173cは、上記の特定時相以外の時相では、第2処理で得られる位置情報、又は、第2処理で得られる位置情報及び第3処理で得られる位置情報を組み合わせた位置情報を第2位置情報とするのが好適である。
ここで、特定時相は、心臓の変形率が大きいとされる代表的な時相であり、追跡外れの起こりやすい時相である大きい拡張早期(e’)時相や、収縮ピーク(s’)時相である。特定位相に応じて第2位置情報を変更する場合、第2推定部173cは、特定時相(e’、s’)においては、第1処理により、「Cdict(t)」を「Cmodel(t)」とする。また、特定位相に応じて第2位置情報を変更する場合、第2推定部173cは、特定時相以外の時相においては、「Ceng(t)」を「Cmodel(t)」とする。或いは、特定位相に応じて第2位置情報を変更する場合、第2推定部173cは、特定時相以外の時相においては、「Ceng(t)」と「Cfit(t)」との線形加算した輪郭位置を「Cmodel(t)」とする。
ここで、第1処理及び第2処理については、第2推定部173cは、第1推定部173a及び定義部173bの処理と並列して実行することが可能である。或いは、第1処理及び第2処理については、第2推定部173cは、第1推定部173a及び定義部173bの処理と独立して実行することが可能である。すなわち、第2推定部173cは、取得部171が動画像データを取得した時点で、残余時相全てで、第1処理及び第2処理を行なっても良い。また、第3処理を行なう場合、第2推定部173cは、追跡対象となる時相ごとに同期させて、第1推定部173a及び定義部173bの処理と並列で行なわれることが好適である。各部の処理の進行については、制御部18により制御可能である。
上記の処理により、第1推定部173aは、図14の上図に示すように、輪郭位置「C(t1)」と動き情報{Mv}とに基づく輪郭位置「Cvector(t2)」を、複数の追跡点「p(t2)」の集合として取得する。一方、第2推定部173cは、図14の上図に示すように、動き情報以外の情報(例えば、形状情報)に基づく輪郭位置「Cmodel(t2)」を、「p(t2)」とは位置や個数が独立な輪郭制御点「q(t2)」の集合として取得する。「Cvector(t2)」は、第1位置情報であり、「Cmodel(t2)」は、第2位置情報となる。そして、複数の追跡点である「p(t2)」それぞれには、動き品質(信頼度指標)が定義されている。
そして、組み合わせ部173dは、信頼度指標(動き品質)に基づいて第1位置情報と第2位置情報とを組み合わせる。具体的には、組み合わせ部173dは、信頼度指標(動き品質)に基づいて第1位置情報と第2位置情報とを重み付け加算する。より具体的には、組み合わせ部173dは、信頼度指標(動き品質)が高い場合に第1位置情報の重みが高くなるよう、第1位置情報と第2位置情報とを重み付け加算する。かかる処理を行なうために、組み合わせ部173dは、「Cmodel(t2)」において、複数の追跡点「p(t2)」それぞれに対応する複数の点「cross_p(t2)」を取得する。例えば、組み合わせ部173dは、図14の上図に示すように、「Cvector(t2)」により定義される輪郭線にて、「p(t2)」それぞれで法線ベクトル「Vnormal_p(t2)」を求める。そして、組み合わせ部173dは、図14の上図に示すように、法線ベクトル「Vnormal_p(t2)」と、「Cmodel(t2)」により定義される輪郭線との交点を、「cross_p(t2)」とする。
そして、組み合わせ部173dは、図14の下図に示すように、「Cmodel(t2)」を、複数の点「cross_p(t2)」で構成された「Cmodel(t2)’」に変換する。
これにより、組み合わせ部173dは、複数の追跡点「p(t2)」それぞれに対応するCmodel(t2)上にある複数の点「cross_p(t2)」を取得する。「p(t2)」の「Q_p(t2)」は、「cross_p(t2)」の動き品質である「Q_cross_p(t2)」となる。
ここで、一般化して、時相「t」の第1位置情報を、複数の追跡点「p(t)」で構成された「Cvector(t)」とし、各追跡点「p(t)」の動き品質を「Q_p(t)」とする。また、一般化して、時相「t」の第2位置情報を、複数の制御点「q(t)」で構成された「Cmodel(t)」とし、複数の点「cross_p(t2)」で「Cmodel(t)」を変換した第2位置情報を、「Cmodel(t)’」とする。かかる場合、組み合わせ部173dは、例えば、以下の式(6)により、「Q_p(t)」に基づいて、「Cvector(t)」と「Cmodel(t)’」とを重み付け加算することで、時相「t」の輪郭位置「C(t)」を求める。なお、「Q_p(t)」は、上記の定義により、「0≦Q_p(t)≦1」に正規化されている。
すなわち、式(6)の第1項では、「Cvector(t2)」上の追跡点「p(t2)」の位置ベクトルに重み「Q_p(t2)」が乗算された位置ベクトルが求められる。また、式(6)の第2項では、「Cmodel(t2)’」上の点「cross_p(t2)」の位置ベクトルに重み「1-Q_p(t2)」が乗算された位置ベクトルが求められる。そして、式(6)では、これら2つの位置ベクトルを加算した位置ベクトルにより定まる点が、時相「t2」の「I(t2)」における輪郭位置「C(t2)」を構成する点として推定される。このようにすれば、各追跡点p(t)に対して、動きの品質が十分に高い場合(Q_p(t)=1)では、従来通り、動き情報に基づくCvector(t)により最終的な輪郭位置C(t)が決まる。一方、動きの品質が極端に低い位置や時相の場合(Q_p(t)=0)では、動き情報以外の情報に基づくCmodel(t)’によりC(t)が決まるようになる。更に、動き品質が中庸な場合(0<Q_p(t)<1)では、Q_p(t)の値に応じて双方の輪郭位置が合成されて、動き品質に応じて自動的に確からしい輪郭位置が与えられるようになる。また、C(t)を構成する追跡点p(t)の個数は一定なので、輪郭位置を決めると同時に追跡動作が達成されている。
図15は、組み合わせ部173dが、短軸像において、「Q_p(t)」に基づいて「Cvector(t)」と「Cmodel(t)」とを合成した輪郭位置「C(t)」を模擬した模式図である。図15では、左側の壁では動き品質が高く、右側の壁では動き品質が低い場合の例を示している。すなわち、図15では、「C(t)」の位置は、左方では「Cvector(t)」に近く、右方では「Cmodel(t)」に近い位置として再構成されている。
そして、組み合わせ部173dは、式(6)により推定された輪郭位置「C(t2)」を構成する複数の点を、複数の追跡点「p(t2)」として、第1推定部173aに通知する。これにより、追跡部173を構成する各部は、上記の説明において、「t1」を「t2」に置き換え、「t2」を「t3」に置き換えた処理を行なう。端的には、第1推定部173aは、図16に示すように、第1画像データである「I(t2)」と第2画像データである「I(t3)」との間でスペックルトラッキングを行なって、Cvector(t3)を推定し、定義部173bは、「p(t3)」の動き品質「Q_(pt3)」を定義する。また、第2推定部173cは、「I(t3)」の「Cmodel(t3)」を推定する。そして、組み合わせ部173dは、図16に示すように、「Cmodel(t3)」を「Cmodel(t3)’」に変換し、「Q_p(t3)」に基づいて「Cvector(t3)」と「Cmodel(t3)」とを合成することで、輪郭位置「C(t3)」を推定する。
追跡部173は、内膜の輪郭位置について、上記の処理を、時相「tn=E1」まで繰り返して行なう。これにより、追跡部173は、全ての時相において、内膜の輪郭位置「C(t1)、C(t2)、・・・、C(tn)」を得る。また、追跡部173は、外膜の輪郭位置についても、上記の処理を、時相「tn=E1」まで繰り返して行なう。これにより、追跡部173は、全ての時相において、外膜の輪郭位置「D(t1)、D(t2)、・・・、D(tn)」を得る。
そして、図1に例示する運動情報演算部174は、動画像データそれぞれの関心領域の位置情報を用いて組織の運動情報を演算する。例えば、運動情報演算部174は、壁運動情報として歪み(Strain)や、歪みの時間変化率(Strain Rate)等、種々の物理量を演算する。操作者は、心機能の解析を行なううえで、評価したい局所的な壁運動情報や全体的な壁運動情報を任意に選択して、運動情報演算部174に演算させることが可能である。そして、制御部18は、運動情報をモニタ2に表示させる。
次に、図17及び図18を用いて、第1の実施形態に係る超音波診断装置の処理について説明する。図17は、第1の実施形態に係る超音波診断装置が行なう処理の概要を説明するためのフローチャートであり、図18は、第1の実施形態に係る追跡部が行なう処理の一例を説明するためのフローチャートである。
図17に例示するように、第1の実施形態に係る超音波診断装置が有する取得部171は、動画像データ及び追跡対象区間が指定されたか否かを判定する(ステップS101)。ここで、指定されない場合(ステップS101否定)、取得部171は、指定されるまで待機する。
一方、動画像データ及び追跡対象区間が指定された場合(ステップS101肯定)、取得部171は、追跡対象区間の動画像データを取得する(ステップS102)。そして、設定部172は、初期輪郭位置が操作者により設定されたか否かを判定する(ステップS103)。ここで、設定されない場合(ステップS103否定)、設定部172は、初期輪郭位置が設定されるまで待機する。
一方、設定された場合(ステップS103肯定)、追跡部173は、初期輪郭位置を取得して、残余時相の輪郭の位置を追跡する(ステップS104)。そして、残余時相の輪郭追跡が終了すると、運動情報演算部174は、心壁の運動情報を演算し(ステップS105)、モニタ2は、制御部18の制御により、心壁の運動情報を表示して(ステップS106)、処理を終了する。
ここで、図17に示すステップS104の処理は、図18に例示するフローチャートの処理となる。すなわち、第1推定部173aは、画像データI(tk)の輪郭位置C(tk)を取得する(ステップS201)。最初の処理では、「tk=t1」となる。
そして、第1推定部173aは、I(tk)からI(tk+1)へのSTにより、{Mv0}を演算する(ステップS202)。そして、第1推定部173aは、統計処理により、残存移動ベクトル群である{Mv}を取得し(ステップS203)、{Mv}を用いて各追跡点でMv_p(tk)を算出して、Cvector(tk+1)を推定する(ステップS204)。
また、定義部173bは、ステップS204の処理が行なわれている間に、Cvector(tk+1)を構成する各追跡点p(tk+1)の信頼度指標である動き品質Q_p(tk+1)を演算する(ステップS205)。第2推定部173cは、ステップS202~ステップS204及びステップS205の処理が行なわれている間に、I(tk+1)のCmodel(tk+1)を推定する(ステップS206)。
そして、組み合わせ部173dは、Cvector(tk+1)を用いて、Cmodel(tk+1)をCmodel(tk+1)’に変換する(ステップS207)。そして、組み合わせ部173dは、各追跡点p(tk+1)の動き品質Q_p(tk+1)に基づいて、Cvector(tk+1)とCmodel(tk+1)’とを組み合わせて、I(tk+1)の輪郭位置C(tk+1)を取得する(ステップS208)。
そして、組み合わせ部173dは、「k+1=n」であるか否かを判定する(ステップS209)。ここで、「k+1<n」である場合(ステップS209否定)、組み合わせ部173dは、輪郭追跡を行なっていない時相があると判定し、「k=k+1」として設定し(ステップS210)、そして、C(tk+1)をC(tk)として、追跡部173は、ステップS202以降の処理を繰り返す。
一方、「k+1=n」である場合(ステップS209肯定)、組み合わせ部173dは、全ての時相で輪郭追跡が終了したことを、運動情報演算部174に通知して、ステップS105の処理を開始させる。
上述したように、第1の実施形態では、関心領域の追跡を行なう際に、動き情報に基づく関心領域の位置推定とともに、動き情報以外の情報に基づく関心領域の位置推定を行なう。そして、第1の実施形態では、動き情報の品質(信頼度指標)を様々な情報から定義し、動き情報に基づく推定位置と動き情報以外の情報に基づく推定位置とを、品質(信頼度指標)に応じて合成する。すなわち、第1の実施形態では、スペックルトラッキング処理において、動き情報が不足している局所的な位置、動き情報の信頼性が低い局所的な位置、動き情報の時空間的整合性が低い局所的な位置を、動き品質(信頼度指標)から同定し、これらの位置については、動き以外の情報で求めた輪郭位置を併用して最終的な輪郭位置を決める。
これにより、第1の実施形態では、画像ノイズの存在やフレームレート不足により動き情報の推定性能が低下する場合であっても、極端な追跡外れの無い関心領域の位置を自動的に取得することができる。従って、第1の実施形態では、関心領域の追跡結果を正確に得ることができる。また、第1の実施形態では、極端な追跡外れの無い関心領域の位置を自動的に取得することができるので、ロバストな組織運動の評価を簡便に提供することが可能となる。
なお、上記の第1の実施形態では、説明を簡単にするため、追跡部173が、関心領域の追跡を時間的に順方向(forward)で行なう場合について説明した。しかし、追跡部173が、関心領域の追跡を時間的に逆方向(backward)で行なう場合であっても良い。
換言すると、第1推定部173aは、時相「tk」から時相「tk+1」へ向かう順方向の移動ベクトル「MvF」を演算することが可能であるとともに、時相「t+1」から時相「tk」へ向かう逆方向の移動ベクトル「MvB」を演算することが可能である。すなわち、上記の第1の実施形態では、「MvF」から、時相「tk」の追跡点p(tk)に対して移動ベクトル「MvF_p(tk→tk+1)」を推定することで、時相「tk+1」の追跡点p(tk+1)の位置を推定していた。
一方、第1推定部173aは、時相「tk+1」から時相「tk」へ向かう逆方向の移動ベクトル「MvB」を演算することができる。そして、第1推定部173aは、「MvB」から、時相「tk+1」の追跡点p(tk+1)に対して移動ベクトル「MvB_p(tk+1→tk)」を推定することができる。これにより、第1推定部173aは、時相「tk」の追跡点p(tk)の位置を推定することができる。
この際、定義部173bは、「MvF」と「MvB」とを併用することで、信頼度指標である動き品質を定義する情報を与えることが可能である。すなわち、「MvB_p(tk+1→tk)」より得られるp(tk)の位置は、「MvF_p(tk→tk+1)」によりp(tk+1)に移動されたp(tk)の位置と必ずしも一致しない。これは、テンプレートマッチングの際に探索対象となる画像データが異なることで、スペックルノイズの影響により「MvF」と「MvB」とが一致しないことがあり得るためである。一方で、ノイズの影響が少ない場合には双方の移動ベクトルの一致度は高いことが期待される。
そこで、「MvF」と「MvB」とから信頼度指標である動き品質を定義する変形例を行なう場合、第1推定部173aは、第1画像データと第2画像データとの間で、関心領域を構成する各追跡点の第1時相から第2時相へ向かう順方向移動ベクトルと、関心領域を構成する各追跡点の第2時相から第1時相へ向かう逆方向移動ベクトルとを推定する。そして、定義部173bは、関心領域を構成する各追跡点における順方向移動ベクトルと逆方向移動ベクトルとの一致度を、「第1推定部173aの処理で得られる少なくとも1つの変量」として用いる。
第1推定部173aは、時相「tk」から時相「tk+1」へ向かう順方向において、「MvF」の演算、「MvF_p(tk→tk+1)」の推定及びp(tk+1)の位置推定を、上記の実施形態で説明した処理により行なう。そして、第1推定部173aは、時相「tk+1」から時相「tk」へ向かう逆方向においては、「MvF」で決まったp(tk+1)に対して、上記の実施形態で説明した処理の時相を逆転させて実行する。これにより、第1推定部173aは、「MvB」の演算、「MvB_p(tk+1→tk)」の推定を行なう。
ここで、「MvF_p(tk→tk+1)」及び「MvB_p(tk+1→tk)」は、いずれも、時相「tk+1」における追跡点p(tk+1)における変数とすることができる。すなわち、これら2つの変数は、一般化すると、時相「t」における追跡点p(t)に対する2つの変数「MvF_p(t)」及び「MvB_p(t)」と表現することができる。かかる表現を用いて、定義部173bは、一例として、以下に示す式(7)により、追跡点p(t)の順方向と逆方向の双方の移動ベクトルの一致度「MM_p(t)」を定義することができる。なお、式(7)における「・」は、ベクトルの内積を示す。
式(7)により「MM_p(t)」を定義すると、追跡点p(t)の順方向移動ベクトルの逆ベクトルが、追跡点p(t)の逆方向移動ベクトルと一致する場合には、「MM_p(t)=1」となる。また、式(7)の定義では、追跡点p(t)の順方向移動ベクトルと追跡点p(t)の逆方向移動ベクトルとが直交する場合には、「MM_p(t)=0」となる。また、「MM_p(t)」の値が「1より小さい正値」となる場合には、追跡点p(t)の順方向移動ベクトルと逆方向移動ベクトルとの成す角が0度より大きく90度より小さい範囲にあり、少なくとも互いに同じ向きのベクトルが検出されていることを意味する。また、「MM_p(t)」の値が「負」となる場合には、追跡点p(t)の順方向移動ベクトルと逆方向移動ベクトルとの成す角が90度を超えており、互いに逆向きのベクトルが検出されていることを意味する。従って、「MM_p(t)」の値が極性も加味して小さいほど、信頼度指標である動き品質は、低いと考えることができる。
そこで、定義部173bは、「MM_p(t)」を用いて、動き品質「Q(s,Mv)」を以下の式(8)により定義する。
すなわち、定義部173bは、「MM_p(t)」が「0」以上となる値が得られた追跡点「p(t)」が含まれるセグメント「s」内の各残存移動ベクトル「Mv」の動き品質「Q(s,Mv)」を、「MM_p(t)」と定義する。また、定義部173bは、「MM_p(t)」の値が「負」となる値が得られた追跡点「p(t)」が含まれるセグメント「s」内の各残存移動ベクトル「Mv」の動き品質「Q(s,Mv)」を「0」と定義する。なお、式(8)で定義される動き品質は、「0≦Q(s,Mv)≦1」の条件を満たしている。なお、定義部173bは、式(5)と式(8)とを併用して、動き品質「Q(s,Mv)」を定義しても良い。
そして、定義部173bは、式(8)で定義された「Q(s,Mv)」、又は、式(5)と式(8)との併用で定義された「Q(s,Mv)」を用いて、「Q_P(t)」を定義する。
なお、上記では、初期輪郭位置が、1つの時相で設定される場合について説明した。しかし、初期輪郭位置は、複数の時相で設定される場合であっても良い。例えば、本実施形態は、スペックルトラッキングが困難とされる拡張早期の時相と収縮末期の時相とで、第1初期輪郭位置及び第2初期輪郭位置を設定しても良い。かかる場合、追跡部173は、第1初期輪郭位置を用いた追跡処理と、第2初期輪郭位置を用いた追跡処理とを並行して行なう。
また、第2推定部173cが実行する処理は、上記で説明した第1処理、第2処理及び第3処理に限定されるものではない。例えば、第2推定部173cは、AQ(Acoustic Quantification)法や、ASM(Active Shape Model)法を実行しても良い。
(第2の実施形態)
第2の実施形態では、演算処理量を軽減するための制御方法について説明する。第1の実施形態では、残余時相の全時相について、第2位置情報「Cmodel(t)」を演算する場合について説明した。しかし、画質によらず安定して正確な関心領域(輪郭)の位置を自動検出するためには、第2位置情報を残余時相の全時相を演算すると演算量が増加する。特に、識別器を用いて「Cdict(t)」を推定する第1処理は、輪郭位置の推定における正確性は増すが、識別器による識別時間を要する処理となる。
第2の実施形態では、演算処理量を軽減するための制御方法について説明する。第1の実施形態では、残余時相の全時相について、第2位置情報「Cmodel(t)」を演算する場合について説明した。しかし、画質によらず安定して正確な関心領域(輪郭)の位置を自動検出するためには、第2位置情報を残余時相の全時相を演算すると演算量が増加する。特に、識別器を用いて「Cdict(t)」を推定する第1処理は、輪郭位置の推定における正確性は増すが、識別器による識別時間を要する処理となる。
一方で、スペックルトラッキング処理では、殆どの場合で、第1位置情報「Cvector(t)」による関心領域(輪郭)の決定が行われるが、信頼度指標である動きの品質が低下する場合にのみ「Cmodel(t)」が必要となる。そこで、第2の実施形態では、動きの品質が低下する時相のみで「Cmodel(t)」を演算して、演算時間を低減する場合について説明する。
第2の実施形態に係る超音波診断装置は、図1に例示した第1の実施形態に係る超音波診断装置と同様の構成となる。また、第2の実施形態に係る追跡部173は、図4に例示した第1の実施形態に係る追跡部173と同様の構成となる。
ただし、第2の実施形態に係る追跡部173は、制御部18の制御により、以下に説明する処理を行なう。すなわち、第2の実施形態に係る追跡部173は、信頼度指標(動き品質)の関心領域における平均値が所定の閾値以上の時相については、第1位置情報を関心領域の位置情報として出力する。また、第2の実施形態に係る追跡部173は、上記の平均値が所定の閾値より小さい時相については、第1位置情報と第2位置情報とを組み合わせた位置情報を関心領域の位置情報として出力する。
以下、第1画像データが「I(t1)」であり、第2画像データが「I(t2)」である場合で、第2の実施形態で行なわれる処理について説明する。第2の実施形態でも、第1の実施形態と同様に、取得部171及び設定部172の処理が行なわれる。そして、第1の実施形態で説明した処理により、第1推定部173aは、「I(t2)」の第1位置情報「Cvector(t)」を推定する。また、第1の実施形態で説明した処理により、定義部173bは、「Cvector(t)」を構成する各点「p(t2)」の動き品質「Q_p(t2)」を定義する。
そして、定義部173bは、各点「p(t2)」での動き品質「Q_p(t2)」を平均した「Q(t2)」を演算する。
制御部18は、「Q(t2)」の値を定義部173bから取得し、品質閾値「Qth」と「Q(t2)」との比較を行なう。「Qth」は、例えば、内部記憶部16に予め格納される。或いは、「Qth」は、例えば、操作者により設定される。
そして、「Q(t2)≧Qth」の場合、制御部18は、第2推定部173c及び組み合わせ部173dの処理を停止し、第1推定部173aが推定した「Cvector(t2)」を「C(t2)」とする。第1推定部173aは、「Cvector(t2)」を「C(t2)」として、運動情報演算部174に出力する。
一方、「Q(t2)<Qth」の場合、制御部18は、第2推定部173c及び組み合わせ部173dの処理を開始させる。これにより、第2推定部173cは、「Cmodel(t2)」を推定し、組み合わせ部173dは、「Cvector(t2)」及び「Cmodel(t2)」から「Cmodel(t2)’」を得る。そして、組み合わせ部173dは、「Cvector(t2)」と「Cmodel(t2)’」とを、「Q_p(t2)」に基づいて、重み付け加算して、「C(t2)」を取得し、運動情報演算部174に出力する。
なお、第2推定部173c及び組み合わせ部173dの処理を平均値に基づいて行なうか否かを判定する点以外は、第1の実施形態で説明した内容は、第2の実施形態においても適用可能である。
次に、図19を用いて、第2の実施形態に係る超音波診断装置の処理について説明する。図19は、第2の実施形態に係る追跡部が行なう処理の一例を説明するためのフローチャートである。
すなわち、図17に示すステップS104の処理は、第2の実施形態では、図19に例示するフローチャートの処理となる。すなわち、第1推定部173aは、画像データI(tk)の輪郭位置C(tk)を取得する(ステップS301)。最初の処理では、「tk=t1」となる。
そして、第1推定部173aは、I(tk)からI(tk+1)へのSTにより、{Mv0}を演算する(ステップS302)。そして、第1推定部173aは、統計処理により、残存移動ベクトル群である{Mv}を取得し(ステップS303)、{Mv}を用いて各追跡点でMv_p(tk)を算出して、Cvector(tk+1)を推定する(ステップS304)。
また、定義部173bは、ステップS304の処理が行なわれている間に、Cvector(tk+1)を構成する各追跡点p(tk+1)の動き品質Q_p(tk+1)を演算する(ステップS305)。そして、定義部173bは、各追跡点p(tk+1)の動き品質Q_p(tk+1)の平均値Q(tk+1)を演算する(ステップS306)。
そして、制御部18は、Q(tk+1)がQthより小さいか否かを判定する(ステップS307)。ここで、Q(tk+1)がQth以上である場合(ステップS307否定)、追跡部173は、Cvector(tk+1)を画像データI(tk+1)の輪郭位置C(tk+1)として取得する(ステップS308)。
一方、Q(tk+1)がQthより小さい場合(ステップS307肯定)、第2推定部173cは、I(tk+1)のCmodel(tk+1)を推定する(ステップS309)。
そして、組み合わせ部173dは、Cvector(tk+1)を用いて、Cmodel(tk+1)をCmodel(tk+1)’に変換する(ステップS310)。そして、組み合わせ部173dは、各追跡点p(tk+1)の動き品質Q_p(tk+1)に基づいて、Cvector(tk+1)とCmodel(tk+1)’とを組み合わせて、I(tk+1)の輪郭位置C(tk+1)を取得する(ステップS311)。
そして、ステップS308又はステップS311の処理の後、組み合わせ部173dは、「k+1=n」であるか否かを判定する(ステップS312)。ここで、「k+1<n」である場合(ステップS312否定)、組み合わせ部173dは、輪郭追跡を行なっていない時相があると判定し、「k=k+1」として設定し(ステップS313)、そして、C(tk+1)をC(tk)として、追跡部173は、ステップS302以降の処理を繰り返す。
一方、「k+1=n」である場合(ステップS312肯定)、組み合わせ部173dは、全ての時相で輪郭追跡が終了したことを、運動情報演算部174に通知して、ステップS105の処理を開始させる。
上述したように、第2の実施形態では、動き情報の品質(信頼度)が高い場合には、第2推定部173c及び組み合わせ部173dの処理を省略するので、関心領域の追跡結果を正確、かつ、効率的に得ることができる。
なお、第2の実施形態は、以下に説明する変形例により行なわれる場合であっても良い。図20は、第2の実施形態に係る変形例を説明するための図である。「Q(t)≧Qth」の場合、全ての「Q_p(t)」が「Qth」以上である第1の場合と、局所的に「Q_p(t)」が「Qth」より小さい第2の場合とに大別できる。そこで、第1の場合には、追跡部173は、Cvector(t)を画像データI(t)の輪郭位置C(t)とする。
一方、第2の場合、すなわち、「Q(t)≧Qth」であるが局所的に「Q_p(t)<Qth」となる点(図20の左図の点線の丸内を参照)が存在する場合には、制御部18の指示により、第2推定部173cは、図13Cを用いて説明した第3処理を行なう。すなわち、第2推定部173cは、図20の右図に示すように、点線の丸内の輪郭を、動き品質の高い周囲の点を用いた最小2乗フィッティング処理により補間することで、「I(t)」の「Cfit(t)」を推定する。そして、第2推定部173cは、「Cfit(t)」を「Cmodel(t)」として、組み合わせ部173dに出力する。
第3処理は、第1処理や第2処理と比較して、演算量が少ない処理である。このため、上記の変形例を行なうことで、全体的には信頼性が高いが、局所的に信頼性が低くなっているCvector(t)を「C(t)」とすることを回避し、迅速に信頼性が高い「C(t)」を得ることができる。
なお、上記の第1及び第2の実施形態では、初期輪郭の設定が、操作者により用手的に行なわれる場合について説明したが、初期輪郭の設定は、上述した動き情報以外の情報により輪郭位置を推定する方法により行なわれる場合であっても良い。換言すると、第1及び第2の実施形態では、設定部172は、第2推定部173cが実行する種々の処理により、初期輪郭を設定しても良い。かかる構成では、関心領域の正確な追跡を、更に、自動化することができ、操作者は、より簡便に組織の運動解析を行なうことが可能となる。なお、初期輪郭を自動設定する場合には、操作者により、初期輪郭位置の修正手続きを受け付け可能な構成とすることが好適である。
また、初期輪郭を自動設定する場合には、第1及び第2の実施形態で説明した画像処理方法を、動画像データの収集とともに略リアルタイムで実行することが可能となる。すなわち、第1及び第2の実施形態に係る超音波診断装置は、動画像データを収集しながら、逐次、第1及び第2の実施形態で説明した画像処理方法を行なうことで、信頼性が安定した組織運動方法を、略リアルタイムで提供することができる。
また、上記の第1及び第2の実施形態では、1つの断面における2次元動画像データに対する追跡処理を行なう場合について説明した。しかし、上記の第1及び第2の実施形態で説明した画像処理方法は、複数の断面における複数の2次元動画像データに対する追跡処理を行なう場合であっても適用可能である。また、上述したように、上記の第1及び第2の実施形態で説明した画像処理方法は、3次元動画像データに対する追跡処理を行なう場合であっても適用可能である。
また、上記の画像処理方法が適用される臓器は、心臓に限らず、心周期に同期して拡張と収縮を繰り返す頸動脈等の動脈血管であっても良い。更には、上記の画像処理方法は、strainや変位といった動き指標を解析することで組織の硬さや動態を把握する用途のために、「肝臓や甲状腺のような軟部組織」や「筋肉」等の運動する組織へ適用しても良い。
また、上記の画像処理方法は、X線診断装置や、X線CT装置、MRI装置等、テンプレートマッチング処理により関心領域の追跡が可能な医用画像データの2次元又は3次元の動画像データに対して適用される場合であっても良い。すなわち、上記の第1及び第2の実施形態で説明した画像処理方法は、超音波診断装置以外の医用画像診断装置で行なわれる場合であっても良い。また、上記の画像処理方法は、医用画像診断装置とは独立に設置された画像処理装置により実行される場合であっても良い。
また、上記の実施形態及び変形例において、図示した各装置の各構成要素は機能概念的なものであり、必ずしも物理的に図示の如く構成されていることを要しない。すなわち、各装置の分散・統合の具体的形態は図示のものに限られず、その全部又は一部を、各種の負荷や使用状況等に応じて、任意の単位で機能的又は物理的に分散・統合して構成することができる。更に、各装置にて行なわれる各処理機能は、その全部又は任意の一部が、CPU及び当該CPUにて解析実行されるプログラムにて実現され、或いは、ワイヤードロジックによるハードウェアとして実現され得る。
また、上記の実施形態及び変形例で説明した画像処理方法は、あらかじめ用意された画像処理プログラムをパーソナルコンピュータやワークステーション等のコンピュータで実行することによって実現することができる。この画像処理プログラムは、インターネット等のネットワークを介して配布することができる。また、この画像処理プログラムは、ハードディスク、フレキシブルディスク(FD)、CD-ROM、MO、DVD等のコンピュータで読み取り可能な記録媒体に記録され、コンピュータによって記録媒体から読み出されることによって実行することもできる。
以上、説明したとおり、第1及び第2の実施系形態、並びに、変形例によれば、関心領域の追跡結果を正確に得ることができる。
本発明のいくつかの実施形態を説明したが、これらの実施形態は、例として提示したものであり、発明の範囲を限定することは意図していない。これら実施形態は、その他の様々な形態で実施されることが可能であり、発明の要旨を逸脱しない範囲で、種々の省略、置き換え、変更を行うことができる。これら実施形態やその変形は、発明の範囲や要旨に含まれると同様に、請求の範囲に記載された発明とその均等の範囲に含まれるものである。
Claims (14)
- 運動する組織を含む領域を撮影対象として収集された2次元又は3次元の超音波画像データの動画像データの中で、少なくとも1つの時相の超音波画像データにて2次元又は3次元の関心領域を設定する設定部と、
前記動画像データが収集された区間のうち前記少なくとも1つの時相以外の残余時相の各超音波画像データで、動き情報に基づいて前記関心領域として推定した領域の第1位置情報と、前記動き情報以外の情報に基づいて前記関心領域として推定した領域の第2位置情報とを取得する推定部と、
前記動き情報の信頼度に関わる信頼度指標に基づいて前記第1位置情報と前記第2位置情報とを合成した位置情報を前記関心領域の位置情報として取得して、前記関心領域を追跡する追跡部と、
を備えたことを特徴とする超音波診断装置。 - 前記動画像データは、前記運動する組織である周期的に運動する組織を撮影対象として収集された少なくとも1周期以上の区間に渡る2次元又は3次元の超音波画像データである、請求項1に記載の超音波診断装置。
- 前記推定部は、
第1時相の超音波画像データである第1画像データと、前記第1時相と時間的に隣接する第2時相の超音波画像データである第2画像データとの間で前記動き情報を推定し、推定した動き情報に基づいて前記第1画像データの関心領域を移動して、前記第2画像データの前記第1位置情報を推定する第1推定部と、
前記第2画像データの前記第2位置情報を推定する第2推定部と、
を有し、
前記追跡部は、
前記第1推定部の処理で得られる少なくとも1つの変量に基づいて、前記信頼度指標を定義する定義部と、
前記信頼度指標に基づいて前記第1位置情報と前記第2位置情報とを組み合わせた位置情報を、前記第2画像データにおける関心領域の位置情報として取得する組み合わせ部と、
を有する、請求項1に記載の超音波診断装置。 - 前記第1推定部は、前記関心領域を含む領域内の複数の点それぞれの移動ベクトルを演算し、得られた移動ベクトル群から統計的に異常な移動ベクトルを除去した残存移動ベクトル群を用いて、前記関心領域を構成する各追跡点の移動ベクトルを推定することで、前記第1位置情報を推定し、
前記定義部は、複数の追跡点それぞれの近傍領域において前記残存移動ベクトル群として採用された移動ベクトルの密度、及び、前記複数の追跡点それぞれの近傍領域内での移動ベクトルの分散値の少なくとも1つを、前記少なくとも1つの変量として用いる、請求項3に記載の超音波診断装置。 - 前記第1推定部は、前記第1画像データと前記第2画像データとの間で、前記関心領域を構成する各追跡点の前記第1時相から前記第2時相へ向かう順方向移動ベクトルと、前記関心領域を構成する各追跡点の前記第2時相から前記第1時相へ向かう逆方向移動ベクトルとを推定し、
前記定義部は、前記関心領域を構成する各追跡点における順方向移動ベクトルと逆方向移動ベクトルとの一致度を、前記少なくとも1つの変量として用いる、請求項3に記載の超音波診断装置。 - 前記第1推定部は、テンプレートマッチング処理により、前記関心領域を含む領域内の複数の点それぞれの移動ベクトルを演算し、
前記定義部は、前記テンプレートマッチング処理におけるテンプレート間の類似度、前記テンプレートマッチング処理における基準テンプレートの信号分散値、前記基準テンプレートの平均信号値の少なくとも1つを、前記少なくとも1つの変量として用いる、請求項3に記載の超音波診断装置。 - 前記第2推定部は、前記第2画像データと同時相の形状辞書情報とを照合する第1処理で得られる位置情報、形状エネルギー最小化原理に基づく第2処理で得られる位置情報、及び、最小2乗法に基づくフィッティングを含む第3処理で得られる位置情報の少なくとも1つの位置情報を用いて、前記第2位置情報を推定する、請求項3に記載の超音波診断装置。
- 前記組み合わせ部は、前記信頼度指標が高い場合に前記第1位置情報の重みが高くなるよう、前記第1位置情報と第2位置情報とを重み付け加算する、請求項3に記載の超音波診断装置。
- 前記信頼度指標の前記関心領域における平均値が所定の閾値以上の時相については、前記第1位置情報を前記関心領域の位置情報として前記追跡部に出力させ、当該平均値が前記所定の閾値より小さい時相については、前記第1位置情報と前記第2位置情報とを組み合わせた位置情報を前記関心領域の位置情報として前記追跡部に出力させる制御部、
を更に備える、請求項1に記載の超音波診断装置。 - 前記第2推定部は、前記形状辞書情報を有している特定時相では、前記第1処理で得られる位置情報を前記第2位置情報とし、前記特定時相以外の時相では、前記第2処理で得られる位置情報、又は、前記第2処理で得られる位置情報及び前記第3処理で得られる位置情報を組み合わせた位置情報を前記第2位置情報とする、請求項7に記載の超音波診断装置。
- 前記第2推定部は、前記信頼度指標が低下する部位に対して前記第2位置情報を推定する、請求項3に記載の超音波診断装置。
- 前記動画像データそれぞれの関心領域の位置情報を用いて前記組織の運動情報を演算する運動情報演算部と、
前記運動情報を表示部に表示させる制御部と、
を更に備える、請求項1に記載の超音波診断装置。 - 運動する組織を含む領域を撮影対象として収集された2次元又は3次元の超音波画像データの動画像データの中で、少なくとも1つの時相の超音波画像データにて2次元又は3次元の関心領域を設定する設定部と、
前記動画像データが収集された区間のうち前記少なくとも1つの時相以外の残余時相の各超音波画像データで、動き情報に基づいて前記関心領域として推定した領域の第1位置情報と、前記動き情報以外の情報に基づいて前記関心領域として推定した領域の第2位置情報とを取得する推定部と、
前記動き情報の信頼度に関わる信頼度指標に基づいて前記第1位置情報と前記第2位置情報とを合成した位置情報を前記関心領域の位置情報として取得して、前記関心領域を追跡する追跡部と、
を備える、画像処理装置。 - 設定部が、運動する組織を含む領域を撮影対象として収集された2次元又は3次元の超音波画像データの動画像データの中で、少なくとも1つの時相の超音波画像データにて2次元又は3次元の関心領域を設定し、
推定部が、前記動画像データが収集された区間のうち前記少なくとも1つの時相以外の残余時相の各超音波画像データで、動き情報に基づいて前記関心領域として推定した領域の第1位置情報と、前記動き情報以外の情報に基づいて前記関心領域として推定した領域の第2位置情報とを取得し、
追跡部が、前記動き情報の信頼度に関わる信頼度指標に基づいて前記第1位置情報と前記第2位置情報とを合成した位置情報を前記関心領域の位置情報として取得して、前記関心領域を追跡する、
ことを含む、画像処理方法。
Priority Applications (2)
Application Number | Priority Date | Filing Date | Title |
---|---|---|---|
CN201380003169.6A CN104080408B (zh) | 2012-11-21 | 2013-11-14 | 超声波诊断装置、图像处理装置以及图像处理方法 |
US14/717,646 US10376236B2 (en) | 2012-11-21 | 2015-05-20 | Ultrasound diagnostic apparatus, image processing apparatus, and image processing method |
Applications Claiming Priority (4)
Application Number | Priority Date | Filing Date | Title |
---|---|---|---|
JP2012255579 | 2012-11-21 | ||
JP2012-255579 | 2012-11-21 | ||
JP2013235938A JP6173886B2 (ja) | 2012-11-21 | 2013-11-14 | 超音波診断装置、画像処理装置及び画像処理方法 |
JP2013-235938 | 2013-11-14 |
Related Child Applications (1)
Application Number | Title | Priority Date | Filing Date |
---|---|---|---|
US14/717,646 Continuation US10376236B2 (en) | 2012-11-21 | 2015-05-20 | Ultrasound diagnostic apparatus, image processing apparatus, and image processing method |
Publications (1)
Publication Number | Publication Date |
---|---|
WO2014080833A1 true WO2014080833A1 (ja) | 2014-05-30 |
Family
ID=50776016
Family Applications (1)
Application Number | Title | Priority Date | Filing Date |
---|---|---|---|
PCT/JP2013/080813 WO2014080833A1 (ja) | 2012-11-21 | 2013-11-14 | 超音波診断装置、画像処理装置及び画像処理方法 |
Country Status (4)
Country | Link |
---|---|
US (1) | US10376236B2 (ja) |
JP (1) | JP6173886B2 (ja) |
CN (1) | CN104080408B (ja) |
WO (1) | WO2014080833A1 (ja) |
Cited By (1)
Publication number | Priority date | Publication date | Assignee | Title |
---|---|---|---|---|
CN112512435A (zh) * | 2018-08-24 | 2021-03-16 | 深圳迈瑞生物医疗电子股份有限公司 | 超声图像处理设备及方法及计算机可读存储介质 |
Families Citing this family (18)
Publication number | Priority date | Publication date | Assignee | Title |
---|---|---|---|---|
CN106030657B (zh) * | 2014-02-19 | 2019-06-28 | 皇家飞利浦有限公司 | 医学4d成像中的运动自适应可视化 |
JP6098677B2 (ja) * | 2014-10-10 | 2017-03-22 | キヤノンマーケティングジャパン株式会社 | 医用画像処理装置、医用画像処理装置に搭載可能なプログラム、及び医用画像処理方法 |
JP6723173B2 (ja) * | 2017-02-10 | 2020-07-15 | 富士フイルム株式会社 | 画像処理装置、方法およびプログラム |
US10302751B2 (en) * | 2017-03-09 | 2019-05-28 | Russell H. Dewey | A-mode ultrasonic classifier |
JP6745237B2 (ja) * | 2017-03-21 | 2020-08-26 | 株式会社日立製作所 | 超音波診断装置 |
JP6815259B2 (ja) | 2017-03-31 | 2021-01-20 | キヤノンメディカルシステムズ株式会社 | 超音波診断装置、医用画像処理装置及び医用画像処理プログラム |
CN107527316B (zh) * | 2017-08-14 | 2019-10-18 | 马鞍山雨甜医疗科技有限公司 | 二维超声影像序列上的任意点构建点云数据的方法及*** |
JP6958202B2 (ja) * | 2017-10-03 | 2021-11-02 | コニカミノルタ株式会社 | 動態画像処理装置及びプログラム |
CN111542827B (zh) | 2017-12-29 | 2021-09-28 | 上海联影医疗科技股份有限公司 | 一种患者定位的***和方法 |
CN116777901A (zh) * | 2017-12-29 | 2023-09-19 | 上海联影医疗科技股份有限公司 | 一种用于医学成像中的生成对象的三维模型的*** |
JP6976869B2 (ja) | 2018-01-15 | 2021-12-08 | キヤノンメディカルシステムズ株式会社 | 超音波診断装置及びその制御プログラム |
US11145054B2 (en) * | 2018-07-09 | 2021-10-12 | Canon U.S.A., Inc. | Methods and systems for displaying intraluminal images |
JP7071898B2 (ja) | 2018-08-09 | 2022-05-19 | 富士フイルムヘルスケア株式会社 | 超音波診断装置、プログラムおよび超音波診断装置の動作方法 |
JP7155277B2 (ja) * | 2019-03-15 | 2022-10-18 | 株式会社東芝 | 構造物評価システム、構造物評価装置及び構造物評価方法 |
JP2019171102A (ja) * | 2019-05-30 | 2019-10-10 | コニカミノルタ株式会社 | 画像処理装置、画像処理方法、およびプログラム |
EP3754607B1 (en) * | 2019-06-19 | 2023-10-04 | Canon Medical Systems Corporation | Ultrasound diagnosis apparatus and ultrasound diagnosis apparatus controlling method |
CN111449684B (zh) * | 2020-04-09 | 2023-05-05 | 济南康硕生物技术有限公司 | 心脏超声标准扫查切面快速获取方法及*** |
JP6973570B1 (ja) * | 2020-06-19 | 2021-12-01 | 沖電気工業株式会社 | 画像処理装置、画像処理プログラム、及び画像処理方法 |
Citations (4)
Publication number | Priority date | Publication date | Assignee | Title |
---|---|---|---|---|
JP2003250804A (ja) * | 2002-03-05 | 2003-09-09 | Toshiba Corp | 画像処理装置及び超音波診断装置 |
JP2008073423A (ja) * | 2006-09-25 | 2008-04-03 | Toshiba Corp | 超音波診断装置、診断パラメータ計測装置及び診断パラメータ計測方法 |
WO2009013871A1 (ja) * | 2007-07-20 | 2009-01-29 | Panasonic Corporation | 超音波診断装置 |
WO2011052602A1 (ja) * | 2009-10-27 | 2011-05-05 | 株式会社 日立メディコ | 超音波イメージング装置、超音波イメージング方法および超音波イメージング用プログラム |
Family Cites Families (8)
Publication number | Priority date | Publication date | Assignee | Title |
---|---|---|---|---|
JP3725442B2 (ja) | 2001-04-18 | 2005-12-14 | 株式会社東芝 | 医用画像診断装置及びその方法 |
DE10132308A1 (de) * | 2001-07-06 | 2003-01-30 | Aventis Behring Gmbh | Kombinationspräparat zur Therapie von immunologischen Erkrankungen |
JP4079690B2 (ja) * | 2002-05-23 | 2008-04-23 | 株式会社東芝 | 物体追跡装置及び方法 |
JP5491836B2 (ja) | 2009-01-30 | 2014-05-14 | 株式会社東芝 | 超音波診断装置、超音波画像処理装置、医用画像診断装置及び医用画像処理装置 |
JP5295805B2 (ja) * | 2009-02-05 | 2013-09-18 | 株式会社東芝 | 超音波診断装置、超音波画像処理装置、及び超音波画像処理プログラム |
US8594398B2 (en) * | 2009-06-26 | 2013-11-26 | International Business Machines Corporation | Systems and methods for cardiac view recognition and disease recognition |
CN102724918A (zh) * | 2010-01-29 | 2012-10-10 | 株式会社日立医疗器械 | 超声波诊断装置及其测量点追踪方法 |
JP5645146B2 (ja) * | 2011-01-31 | 2014-12-24 | 日本電気株式会社 | 情報処理システム、情報処理方法、情報処理装置及びその制御方法とその制御プログラム |
-
2013
- 2013-11-14 WO PCT/JP2013/080813 patent/WO2014080833A1/ja active Application Filing
- 2013-11-14 JP JP2013235938A patent/JP6173886B2/ja active Active
- 2013-11-14 CN CN201380003169.6A patent/CN104080408B/zh active Active
-
2015
- 2015-05-20 US US14/717,646 patent/US10376236B2/en active Active
Patent Citations (4)
Publication number | Priority date | Publication date | Assignee | Title |
---|---|---|---|---|
JP2003250804A (ja) * | 2002-03-05 | 2003-09-09 | Toshiba Corp | 画像処理装置及び超音波診断装置 |
JP2008073423A (ja) * | 2006-09-25 | 2008-04-03 | Toshiba Corp | 超音波診断装置、診断パラメータ計測装置及び診断パラメータ計測方法 |
WO2009013871A1 (ja) * | 2007-07-20 | 2009-01-29 | Panasonic Corporation | 超音波診断装置 |
WO2011052602A1 (ja) * | 2009-10-27 | 2011-05-05 | 株式会社 日立メディコ | 超音波イメージング装置、超音波イメージング方法および超音波イメージング用プログラム |
Non-Patent Citations (1)
Title |
---|
KIYOTSUGU SEKIOKA ET AL.: "Assessment of Regional Intramyocardial Layer Function from Ultrasonic RF Echo Signals Using Hierarchical Correlation Method with Confidence", THE TRANSACTIONS OF THE INSTITUTE OF ELECTRONICS, INFORMATION AND COMMUNICATION ENGINEERS, vol. J87-D-II, no. 1, 1 January 2004 (2004-01-01), pages 98 - 108 * |
Cited By (1)
Publication number | Priority date | Publication date | Assignee | Title |
---|---|---|---|---|
CN112512435A (zh) * | 2018-08-24 | 2021-03-16 | 深圳迈瑞生物医疗电子股份有限公司 | 超声图像处理设备及方法及计算机可读存储介质 |
Also Published As
Publication number | Publication date |
---|---|
US20150257731A1 (en) | 2015-09-17 |
CN104080408A (zh) | 2014-10-01 |
CN104080408B (zh) | 2016-06-01 |
JP6173886B2 (ja) | 2017-08-02 |
JP2014121541A (ja) | 2014-07-03 |
US10376236B2 (en) | 2019-08-13 |
Similar Documents
Publication | Publication Date | Title |
---|---|---|
JP6173886B2 (ja) | 超音波診断装置、画像処理装置及び画像処理方法 | |
JP6640922B2 (ja) | 超音波診断装置及び画像処理装置 | |
JP6041350B2 (ja) | 超音波診断装置、画像処理装置及び画像処理方法 | |
JP6132614B2 (ja) | 超音波診断装置、画像処理装置及び画像処理方法 | |
JP6125281B2 (ja) | 医用画像診断装置、医用画像処理装置及び制御プログラム | |
JP5624258B2 (ja) | 超音波診断装置、超音波画像処理装置及び超音波画像処理プログラム | |
US9254115B2 (en) | Ultrasonic diagnostic apparatus for cardiac wall movement measurements by re-tracking the cardiac wall | |
US20080304730A1 (en) | Ultrasonic image processing apparatus and method for processing ultrasonic image | |
WO2005115249A1 (ja) | 超音波診断装置、超音波画像処理装置、及び超音波画像処理方法 | |
JP6381972B2 (ja) | 医用画像処理装置および医用画像診断装置 | |
JP7375140B2 (ja) | 超音波診断装置、医用画像診断装置、医用画像処理装置及び医用画像処理プログラム | |
US9877698B2 (en) | Ultrasonic diagnosis apparatus and ultrasonic image processing apparatus | |
JP6925824B2 (ja) | 超音波診断装置、画像処理装置、及び画像処理プログラム | |
JP5788230B2 (ja) | 超音波診断装置、超音波画像処理装置、超音波画像処理プログラム | |
JP6815259B2 (ja) | 超音波診断装置、医用画像処理装置及び医用画像処理プログラム | |
JP2007222533A (ja) | 超音波診断装置及び超音波画像処理方法 | |
JP6430558B2 (ja) | 超音波診断装置、画像処理装置及び画像処理方法 | |
JP4745455B2 (ja) | 超音波診断装置、超音波画像処理装置、及び超音波信号処理プログラム | |
JP7483519B2 (ja) | 超音波診断装置、医用画像処理装置、及び医用画像処理プログラム | |
JP7019287B2 (ja) | 超音波診断装置及び画像処理プログラム | |
JP2014176749A (ja) | 超音波診断装置及び超音波画像処理装置 | |
JP2020049212A (ja) | 装置、医用情報処理装置、及びプログラム |
Legal Events
Date | Code | Title | Description |
---|---|---|---|
121 | Ep: the epo has been informed by wipo that ep was designated in this application |
Ref document number: 13856695 Country of ref document: EP Kind code of ref document: A1 |
|
NENP | Non-entry into the national phase |
Ref country code: DE |
|
122 | Ep: pct application non-entry in european phase |
Ref document number: 13856695 Country of ref document: EP Kind code of ref document: A1 |