CN112438702B - Photoacoustic endoscopic imaging method and system for biological cavity - Google Patents

Photoacoustic endoscopic imaging method and system for biological cavity Download PDF

Info

Publication number
CN112438702B
CN112438702B CN201910806275.7A CN201910806275A CN112438702B CN 112438702 B CN112438702 B CN 112438702B CN 201910806275 A CN201910806275 A CN 201910806275A CN 112438702 B CN112438702 B CN 112438702B
Authority
CN
China
Prior art keywords
light absorption
absorption energy
sound
sound velocity
photoacoustic
Prior art date
Legal status (The legal status is an assumption and is not a legal conclusion. Google has not performed a legal analysis and makes no representation as to the accuracy of the status listed.)
Active
Application number
CN201910806275.7A
Other languages
Chinese (zh)
Other versions
CN112438702A (en
Inventor
孙正
段爽
Current Assignee (The listed assignees may be inaccurate. Google has not performed a legal analysis and makes no representation or warranty as to the accuracy of the list.)
North China Electric Power University
Original Assignee
North China Electric Power University
Priority date (The priority date is an assumption and is not a legal conclusion. Google has not performed a legal analysis and makes no representation as to the accuracy of the date listed.)
Filing date
Publication date
Application filed by North China Electric Power University filed Critical North China Electric Power University
Priority to CN201910806275.7A priority Critical patent/CN112438702B/en
Publication of CN112438702A publication Critical patent/CN112438702A/en
Application granted granted Critical
Publication of CN112438702B publication Critical patent/CN112438702B/en
Active legal-status Critical Current
Anticipated expiration legal-status Critical

Links

Images

Classifications

    • A—HUMAN NECESSITIES
    • A61—MEDICAL OR VETERINARY SCIENCE; HYGIENE
    • A61B—DIAGNOSIS; SURGERY; IDENTIFICATION
    • A61B5/00—Measuring for diagnostic purposes; Identification of persons
    • A61B5/0093—Detecting, measuring or recording by applying one single type of energy and measuring its conversion into another type of energy
    • A61B5/0095—Detecting, measuring or recording by applying one single type of energy and measuring its conversion into another type of energy by applying light and detecting acoustic waves, i.e. photoacoustic measurements
    • G—PHYSICS
    • G16—INFORMATION AND COMMUNICATION TECHNOLOGY [ICT] SPECIALLY ADAPTED FOR SPECIFIC APPLICATION FIELDS
    • G16H—HEALTHCARE INFORMATICS, i.e. INFORMATION AND COMMUNICATION TECHNOLOGY [ICT] SPECIALLY ADAPTED FOR THE HANDLING OR PROCESSING OF MEDICAL OR HEALTHCARE DATA
    • G16H30/00—ICT specially adapted for the handling or processing of medical images
    • G16H30/40—ICT specially adapted for the handling or processing of medical images for processing medical images, e.g. editing

Landscapes

  • Health & Medical Sciences (AREA)
  • Life Sciences & Earth Sciences (AREA)
  • Public Health (AREA)
  • Engineering & Computer Science (AREA)
  • Physics & Mathematics (AREA)
  • General Health & Medical Sciences (AREA)
  • Medical Informatics (AREA)
  • Acoustics & Sound (AREA)
  • Biomedical Technology (AREA)
  • Epidemiology (AREA)
  • Radiology & Medical Imaging (AREA)
  • Nuclear Medicine, Radiotherapy & Molecular Imaging (AREA)
  • Biophysics (AREA)
  • Pathology (AREA)
  • Primary Health Care (AREA)
  • Heart & Thoracic Surgery (AREA)
  • Molecular Biology (AREA)
  • Surgery (AREA)
  • Animal Behavior & Ethology (AREA)
  • Veterinary Medicine (AREA)
  • Ultra Sonic Daignosis Equipment (AREA)

Abstract

The invention discloses a photoacoustic endoscopic imaging method and a photoacoustic endoscopic imaging system for a biological cavity, which are used for acquiring an initial value A of light absorption energy at a position r of the biological cavity 0 (r) and initial value of sound velocity c s,0 (r) calculating light absorption energy A of kth time at position r k (r) and (k + 1) th light absorption energy A k+1 (r), speed of sound c of kth s,k (r) and sound velocity c of k +1 st time s,k+1 (r); by judging the absolute difference epsilon of light absorption energy A,k Sum sound velocity absolute difference
Figure DDA0002183766250000011
And whether a convergence condition is met or not is determined, and whether a light absorption energy distribution graph and an acoustic velocity distribution graph are continuously iterated or directly constructed. The method considers the complexity of biological tissues, comprehensively considers the difference of the speeds, combines light absorption energy and different sound velocities to construct a photoacoustic image, reduces errors caused by uneven sound velocity distribution in the tissues to be detected, and obtains a high-quality cavity cross section light absorption energy distribution map.

Description

Photoacoustic endoscopic imaging method and system for biological cavity
Technical Field
The invention relates to the technical field of medical imaging, in particular to a photoacoustic endoscopic imaging method and a photoacoustic endoscopic imaging system for a biological cavity.
Background
A biological Photoacoustic endoscopic (PAE) imaging method is a novel medical functional imaging method, a special imaging catheter is directly inserted into a target cavity, a catheter probe emits short pulse laser to irradiate surrounding tissues, the tissues absorb light energy and then are heated and expanded, ultrasonic waves are excited, and the ultrasonic waves are transmitted to the surfaces of the tissues, namely Photoacoustic signals. An ultrasound probe at the tip of the imaging catheter receives photoacoustic signals from various directions during the intraluminal scan of the surrounding tissue. After the photoacoustic pressure time sequences collected at different measurement positions are sent to a computer, a light absorption energy distribution image of the cross section of the cavity can be inverted and reconstructed by adopting a proper image reconstruction algorithm, and the morphological structure and the functional components of the tissue are reflected.
Given the complexity of biological tissue, there is often a large difference in the speed of ultrasound as it propagates through tissue having different compositions. In many cases, the sound velocity distribution in the tissue cannot be predicted before PAE imaging is performed, so that the assumption of a constant sound velocity in the process of reconstructing an image may cause serious problems such as acoustic distortion, artifacts, blurring, and object dislocation in the reconstructed image.
Disclosure of Invention
The invention aims to provide a photoacoustic endoscopic imaging method and a photoacoustic endoscopic imaging system for a biological cavity, which are used for reducing errors caused by uneven sound velocity distribution in tissues to be measured and obtaining a high-quality light absorption energy distribution image of the cross section of the cavity.
In order to achieve the purpose, the invention provides the following scheme:
a method of photoacoustic endoscopic imaging of a biological cavity, the method comprising:
obtaining the initial value A of the light absorption energy at the position r of the biological cavity 0 (r) and initial value of sound velocity c s,0 (r)；
Respectively calculating the light absorption energy A of the kth time at the position r of the biological cavity k (r) and (k + 1) th light absorption energy A k+1 (r); k is iteration times, and k is more than or equal to 1;
respectively calculating the sound velocity c of the kth time at the position r of the biological cavity s,k (r) and sound velocity c of k +1 st time s,k+1 (r)；
Calculating the light absorption energy A of the k time k (r) and the light absorption energy A of the k +1 th order k+1 Absolute difference epsilon of light absorption energy of (r) A,k And the speed of sound c of the kth time s,k (r) speed of sound c of the k +1 th time s,k+1 Absolute difference in sound velocity of (r)
Figure BDA0002183766230000021
Judging the absolute difference epsilon of the light absorption energy A,k And the absolute difference of sound velocity
Figure BDA0002183766230000022
Whether the convergence condition is met or not is judged, and a convergence judgment result is obtained;
if the convergence is judgedAs a result, the (k + 1) th light absorption energy A at the output position r k+1 (r) and speed of sound c s,k+1 (r) constructing a light absorption energy distribution map and an acoustic velocity distribution map, respectively;
if the convergence judgment result is negative, updating the iteration times, and returning to the step of respectively calculating the light absorption energy A of the kth time at the position r of the biological cavity k (r) and (k + 1) th light absorption energy A k+1 (r)。
Optionally, according to formula A k (r)＝A k-1 (r)-W -1 (A k-1 (r))G 1 '(r,c s,k-1 (r),A k-1 (r)) calculating the light absorption energy A of the k-th time k (r) and the light absorption energy A of the k +1 th order k+1 (r)；
Wherein,
Figure BDA0002183766230000023
argmin[·]is that g (A) k (r))+λΦ TV (A k (r)) A at the minimum value k (r)；
W(A k-1 (r)) is an approximate blackplug matrix of the first photoacoustic relationship function;
W -1 (A k-1 (r)) is W (A) k-1 (r)) an inverse matrix;
a is a light absorption energy distribution matrix; a is more than or equal to 0, which means the light absorption energy which is more than or equal to 0 in the light absorption energy distribution matrix;
g(A k (r))＝||p m (r)-H(c s,k-1 (r))·A k (r)|| 2 is a second photoacoustic relationship function;
p m (r) is a measurement of the photoacoustic signal at location r;
H(c s,k-1 (r)) is a first operator related to the speed of sound;
| | · | is a 2-norm;
Φ TV (A k (r)) is a TV regularization term,
Figure BDA0002183766230000024
λ is the TV regularization parameter;
L k is about A k (r) a sparse matrix of variance characteristics;
eta > 0 is a constant;
G 1 '(r,A k-1 (r),c s,k-1 (r)) is the gradient of the photoacoustic relationship function.
Optionally, according to a formula
Figure BDA0002183766230000031
Respectively calculating the sound velocity c of the k-th time s,k (r) and the sound speed c of the k +1 st time s,k+1 (r)；
Wherein,
Figure BDA0002183766230000032
▽f(q k ) Is f (q) k ) A gradient of (a);
f(q k )＝||p m (r)-H(q k )·A k-1 (r)|| 2 ；
Figure BDA0002183766230000033
γ k-1 is the weight parameter, γ, after the k-1 iteration k-2 Is the weight parameter after the k-2 iteration;
weight parameter gamma after the kth iteration k Is composed of
Figure BDA0002183766230000034
H(q k ) Is a second operator related to the speed of sound;
d is f (c) s,k-1 (r)) the prestz constant of the derivative;
f(c s,k-1 (r))＝||p m (r)-H(c s,k-1 (r))·A k-1 (r)|| 2 。
optionally, the convergence condition is epsilon A,k ＜ε A And is provided with
Figure BDA0002183766230000035
Wherein epsilon A For the convergence tolerance of the light absorption energy,
Figure BDA0002183766230000036
is the convergence tolerance of the speed of sound.
Optionally, the constructing the light absorption energy distribution image and the sound velocity distribution image specifically includes:
the light absorption energy distribution matrix A under the polar coordinates is normalized and subjected to gray scale processing by adopting the following formula:
Figure BDA0002183766230000037
b (i, j) is a value obtained by normalizing and graying a (i, j);
i is the abscissa of the biological cavity position r under the corresponding polar coordinate system, and j is the ordinate of the biological cavity position r under the corresponding polar coordinate system;
a (i, j) is the element value of the ith row and the jth column in the light absorption energy distribution matrix A;
min (A) is the minimum value of the elements in the light absorption energy distribution matrix A, and max (A) is the maximum value of the elements in the light absorption energy distribution matrix A;
adopting the following formula to align the sound velocity distribution matrix c under polar coordinates s Carrying out normalization and gray level processing;
Figure BDA0002183766230000041
e (i, j) is for c s (i, j) normalized and grayed values;
c s (i, j) is the sound velocity distribution matrix c s The value of the element in the ith row and the jth column;
min(c s ) Is the sound velocity distribution matrix c s Minimum value of (1); max (c) s ) Is the sound velocity distribution matrix c s Maximum value of (1);
converting said B (i, j) into B '(x, y) in planar rectangular coordinates according to the formula B' (x, y) = B (jcosi, jsini);
converting said E (i, j) into E '(x, y) in planar rectangular coordinates according to the formula E' (x, y) = E (jcosi, jsini);
wherein i belongs to [0,360], j belongs to [0,d ];
d is the maximum value of the polar diameter in the theta-l polar coordinate system;
b' (x, y) is a gray value of a point (x, y) in a rectangular coordinate system in the light absorption energy distribution map; e' (x, y) is the grayscale value of the point (x, y) in the rectangular coordinate system in the sound velocity map.
A bio-cavity photoacoustic endoscopic imaging system, the system comprising:
an initial value acquisition module for acquiring an initial value A of light absorption energy at a position r of the biological cavity 0 (r) and initial value of sound velocity c s,0 (r)；
A light absorption energy calculation module for calculating the light absorption energy A of the kth time at the position r of the biological cavity respectively k (r) and (k + 1) th light absorption energy A k+1 (r); k is iteration times, and k is more than or equal to 1;
a sound velocity calculation module for respectively calculating the kth sound velocity c at the position r of the biological cavity s,k (r) and sound velocity c of (k + 1) th order s,k+1 (r)；
An absolute difference calculation module for calculating the light absorption energy A of the kth time k (r) and the light absorption energy A of the k +1 th order k+1 Absolute difference of light absorption energy ε of (r) A,k And the speed of sound c of the kth time s,k (r) speed of sound c at the k +1 th time s,k+1 Absolute difference in sound velocity of (r)
Figure BDA0002183766230000051
A convergence judgment result module for judging the light absorption energyAbsolute difference epsilon A,k And the absolute difference of sound velocity
Figure BDA0002183766230000052
Whether the convergence condition is met or not is judged, and a convergence judgment result is obtained;
an image construction module connected with the convergence judgment result module and used for absorbing energy A according to the k +1 th time light at the position r when the convergence judgment result module is yes k+1 (r) and speed of sound c s,k+1 (r) respectively constructing a light absorption energy distribution map and an acoustic velocity distribution map;
and the iteration updating module is respectively connected with the convergence judgment result module and the light absorption energy calculating module, and is used for updating the iteration times and returning to the light absorption energy calculating module when the convergence judgment result module is negative.
Optionally, according to formula A k (r)＝A k-1 (r)-W -1 (A k-1 (r))G 1 '(r,c s,k-1 (r),A k-1 (r)) calculating the light absorption energy A of the k-th time k (r) and the light absorption energy A of the k +1 th order k+1 (r)；
Wherein,
Figure BDA0002183766230000053
argmin[·]is that g (A) k (r))+λΦ TV (A k (r)) A at the minimum k (r)；
W(A k-1 (r)) is an approximate blackplug matrix of the first photoacoustic relationship function;
W -1 (A k-1 (r)) is W (A) k-1 (r)) an inverse matrix;
a is a light absorption energy distribution matrix; a is more than or equal to 0, which means the light absorption energy which is more than or equal to 0 in the light absorption energy distribution matrix;
g(A k (r))＝||p m (r)-H(c s,k-1 (r))·A k (r)|| 2 is a second photoacoustic relationship function;
p m (r) is a measurement of the photoacoustic signal at location r;
H(c s,k-1 (r)) is a first operator related to the speed of sound;
| | · | is a 2-norm;
Φ TV (A k (r)) is a TV regularization term,
Figure BDA0002183766230000054
λ is a TV regularization parameter;
L k is about A k (r) a sparse matrix of variance characteristics;
η > 0 is a constant;
G 1 '(r,A k-1 (r),c s,k-1 (r)) is the gradient of the photoacoustic relationship function.
Optionally according to a formula
Figure BDA0002183766230000061
Respectively calculating the sound velocity c of the k-th time s,k (r) and the speed of sound c after the (k + 1) th iteration s,k+1 (r)；
Wherein,
Figure BDA0002183766230000062
▽f(q k ) Is f (q) k ) A gradient of (a);
f(q k )＝||p m (r)-H(q k )·A k-1 (r)|| 2 ；
Figure BDA0002183766230000063
γ k-1 is the weight parameter, γ, after the k-1 iteration k-2 Is the weight parameter after the k-2 iteration;
weight parameter gamma after the kth iteration k Is composed of
Figure BDA0002183766230000064
H(q k ) Is a second operator related to the speed of sound;
d is f (c) s,k-1 (r)) the petz constant of the derivative;
f(c s,k-1 (r))＝||p m (r)-H(c s,k-1 (r))·A k-1 (r)|| 2 。
optionally, the convergence condition is epsilon A,k ＜ε A And is
Figure BDA0002183766230000065
Wherein epsilon A For the convergence tolerance of the light absorption energy,
Figure BDA0002183766230000066
is the convergence tolerance of the speed of sound.
Optionally, the image construction module specifically includes:
the light absorption energy distribution matrix processing unit is used for carrying out normalization and gray scale processing on the light absorption energy distribution matrix A under the polar coordinates by adopting the following formula:
Figure BDA0002183766230000071
b (i, j) is a value obtained by normalizing and graying a (i, j);
i is the abscissa of the biological cavity position r under the corresponding polar coordinate system, and j is the ordinate of the biological cavity position r under the corresponding polar coordinate system;
a (i, j) is the element value of the ith row and the jth column in the light absorption energy distribution matrix A;
min (A) is the minimum value of the elements in the light absorption energy distribution matrix A, max (A) is the maximum value of the elements in the light absorption energy distribution matrix A;
a sound velocity distribution matrix processing unit for aligning the sound velocity distribution matrix c under polar coordinates by using the following formula s Carrying out normalization and gray level processing;
Figure BDA0002183766230000072
e (i, j) is for c s (i, j) normalized and grayed values;
c s (i, j) is the sound velocity distribution matrix c s The value of the element in the ith row and the jth column;
min(c s ) Is the sound velocity distribution matrix c s Minimum value of (1); max (c) s ) Is the sound velocity distribution matrix c s Maximum value of (1);
a light absorption energy distribution map gray value conversion unit for converting the B (i, j) into B '(x, y) in planar rectangular coordinates according to a formula B' (x, y) = B (jcosi, jsini);
a sound velocity profile gradation value conversion unit for converting the E (i, j) into E '(x, y) in planar rectangular coordinates according to a formula E' (x, y) = E (jcosi, jsini);
wherein, i belongs to [0,360], j belongs to [0,d ];
d is the maximum value of the polar diameter in the theta-l polar coordinate system;
b' (x, y) is a gray value of a point (x, y) in a rectangular coordinate system in the light absorption energy distribution map; e' (x, y) is the grayscale value of the point (x, y) in the rectangular coordinate system in the sound velocity map.
According to the specific embodiment provided by the invention, the invention discloses the following technical effects: the method obtains an initial value A of light absorption energy at a position r of a biological cavity 0 (r) and initial value of sound velocity c s,0 (r) and then calculates the light absorption energy A of the kth time at the position r, respectively k (r) and (k + 1) th light absorption energy A k+1 (r), speed of sound c of kth s,k (r) and sound velocity c of (k + 1) th order s,k+1 (r); by judging the absolute difference epsilon of light absorption energy A,k And absolute difference of sound velocity
Figure BDA0002183766230000073
Whether the convergence condition is satisfied, whether to continue iteration or to directly construct lightAn absorption energy distribution image and a sound speed distribution image. The invention considers the complexity of biological tissues, comprehensively considers the difference of the speeds, and reconstructs photoacoustic images by combining light absorption energy and different sound velocities, thereby reducing errors caused by uneven sound velocity distribution in tissues to be detected and obtaining high-quality light absorption energy distribution images of the cross section of the cavity.
Drawings
In order to more clearly illustrate the embodiments of the present invention or the technical solutions in the prior art, the drawings required in the embodiments will be briefly described below, it is obvious that the drawings in the following description are only some embodiments of the present invention, and it is obvious for those skilled in the art that other drawings can be obtained according to the drawings without creative efforts.
FIG. 1 is a schematic view of a photoacoustic endoscopic imaging method for a biological cavity according to the present invention;
FIG. 2 is a block diagram of a photoacoustic endoscopic imaging system for a biological cavity according to the present invention;
FIG. 3 is a schematic cross-sectional view of a biological cavity in a rectangular XOY coordinate system;
FIG. 4 is a schematic diagram of the multi-layer cavity wall structure in a polar θ -l coordinate system.
Detailed Description
The technical solutions in the embodiments of the present invention will be clearly and completely described below with reference to the drawings in the embodiments of the present invention, and it is obvious that the described embodiments are only a part of the embodiments of the present invention, and not all of the embodiments. All other embodiments, which can be derived by a person skilled in the art from the embodiments given herein without making any creative effort, shall fall within the protection scope of the present invention.
The invention aims to provide a photoacoustic endoscopic imaging method and a photoacoustic endoscopic imaging system for a biological cavity, which are more accurate in constructing sub-images by considering the complexity of biological tissues and combining the conditions of different sound velocities in the process of constructing a light absorption energy distribution image.
In order to make the aforementioned objects, features and advantages of the present invention comprehensible, embodiments accompanied with figures are described in further detail below.
A method of photoacoustic endoscopic imaging of a biological cavity, the method comprising:
obtaining the initial value A of the light absorption energy at the position r of the biological cavity 0 (r) and initial value of sound velocity c s,0 (r)。
In this embodiment, the initial value of the light absorption energy is A 0 (r) =0, initial value of sound velocity c s,0 (r)＝1600m/s。
Respectively calculating the light absorption energy A of the kth time at the position r of the biological cavity k (r) and (k + 1) th light absorption energy A k+1 (r); k is iteration times, and k is more than or equal to 1;
specifically, this embodiment is based on formula A k (r)＝A k-1 (r)-W -1 (A k-1 (r))G 1 '(r,c s,k-1 (r),A k-1 (r)) calculating the light absorption energy A of the k-th time k (r) and the light absorption energy A of the k +1 th order k+1 (r)；
Wherein,
Figure BDA0002183766230000091
argmin[·]is that g (A) k (r))+λΦ TV (A k (r)) A at the minimum k (r)；
W(A k-1 (r)) is an approximate blackplug matrix of the first photoacoustic relationship function;
W -1 (A k-1 (r)) is W (A) k-1 (r)) an inverse matrix;
a is a light absorption energy distribution matrix; a is more than or equal to 0, which means the light absorption energy which is more than or equal to 0 in the light absorption energy distribution matrix;
g(A k (r))＝||p m (r)-H(c s,k-1 (r))·A k (r)|| 2 is a second photoacoustic relationship function;
p m (r) is a measurement of the photoacoustic signal at location r;
H(c s,k-1 (r)) is related to the speed of soundA first operator;
Figure BDA0002183766230000092
p (r) is the theoretical value of the photoacoustic signal at location r;
| | · | is a 2-norm;
Φ TV (A k (r)) is a TV regularization term,
Figure BDA0002183766230000093
λ is the TV regularization parameter;
L k is about A k (r) a sparse matrix of variance characteristics;
eta > 0 is a constant;
G 1 '(r,A k-1 (r),c s,k-1 (r)) is the gradient of the photoacoustic relationship function.
In the specific treatment process, the light energy deposition in the cavity wall tissue under the short pulse laser irradiation is simulated by adopting a Monte Carlo simulation method. Then, solving a discrete photoacoustic wave equation under a polar coordinate system by adopting a finite difference time domain algorithm to obtain a theoretical value of a photoacoustic signal generated by a tissue:
Figure BDA0002183766230000101
wherein, (i, j) is the coordinate of a point r on the cross section of the cavity in a theta-l polar coordinate system; Δ θ and Δ l are unit lengths on the θ axis and the l axis, respectively; Δ t is a discrete time interval; n is a discrete time; p is a radical of (n+1) (i, j) is a theoretical value of a photoacoustic signal generated by a particle point with a position of (i, j) at the time n + 1;
Figure BDA0002183766230000102
and
Figure BDA0002183766230000103
the particle at position (i, j) is in the θ direction at time nAnd vibration speed in the l direction; c. C s (i, j) is the ultrasound at location (i, j); the propagation velocity of (c); β is the isobaric expansion coefficient of the tissue; c p Is the specific heat capacity of the tissue; ρ is a unit of a gradient 0 Is the density of the tissue; i is (n) Is the value of the laser pulse function at time n; a (i, j) is the light absorption energy at location (i, j).
Respectively calculating the sound velocity c of the kth time at the position r of the biological cavity s,k (r) and sound velocity c of (k + 1) th order s,k+1 (r)；
In the present embodiment, the first and second electrodes are,
Figure BDA0002183766230000104
after k iteration is carried out on the third photoacoustic relation function, the sound velocity at the position r is obtained
Figure BDA0002183766230000105
After the solution is carried out, the solution is obtained,
Figure BDA0002183766230000106
according to the formula
Figure BDA0002183766230000107
Respectively calculating the sound velocity c of the k-th time s,k (r) and the sound speed c of the k +1 st time s,k+1 (r)。
Wherein,
Figure BDA0002183766230000111
▽f(q k ) Is f (q) k ) A gradient of (a);
f(q k )＝||p m (r)-H(q k )·A k-1 (r)|| 2 ；
Figure BDA0002183766230000112
γ k-1 is the weight parameter, γ, after the k-1 iteration k-2 Is the weight parameter after the k-2 iteration;
weight parameter gamma after the kth iteration k Is composed of
Figure BDA0002183766230000113
H(q k ) Is a second operator related to the speed of sound;
d is f (c) s,k-1 (r)) the prestz constant of the derivative;
f(c s,k-1 (r))＝||p m (r)-H(c s,k-1 (r))·A k-1 (r)|| 2 ；
c s sound velocity distribution matrix, c s And > 0 is the sound velocity of 0 or more in the sound velocity distribution matrix.
Calculating the light absorption energy A of the k time k (r) and the light absorption energy A of the k +1 th order k+1 Absolute difference of light absorption energy ε of (r) A,k And the speed of sound c of the kth time s,k (r) speed of sound c at the k +1 th time s,k+1 Absolute difference in sound velocity of (r)
Figure BDA0002183766230000114
The convergence condition is epsilon A,k ＜ε A And is
Figure BDA0002183766230000115
Wherein epsilon A For the convergence tolerance of the light absorption energy,
Figure BDA0002183766230000116
is the convergence tolerance of the speed of sound.
In the present embodiment, the convergence tolerance ε of light absorption energy A Taking 0.01, the convergence tolerance of the speed of sound
Figure BDA0002183766230000117
0.01 is taken.
After iterative computationIn the process, if all the iteration values do not satisfy the convergence condition, the convergence tolerance epsilon of the light absorption energy can be changed A Convergence tolerance to speed of sound
Figure BDA0002183766230000118
And (4) adjusting the convergence condition to complete the iterative process.
Judging the absolute difference epsilon of the light absorption energy A,k And the absolute difference of sound velocity
Figure BDA0002183766230000119
Whether the convergence condition is met or not is judged, and a convergence judgment result is obtained;
if the convergence judgment result is yes, outputting the (k + 1) th light absorption energy A at the position r k+1 (r) and speed of sound c s,k+1 (r) respectively constructing a light absorption energy distribution map and an acoustic velocity distribution map;
if the convergence judgment result is negative, updating the iteration times, and returning to the step of respectively calculating the light absorption energy A of the kth time at the position r of the biological cavity k (r) and (k + 1) th light absorption energy A k+1 (r)。
The constructing of the light absorption energy distribution image and the sound velocity distribution image specifically includes:
the light absorption energy distribution matrix A under the polar coordinates is normalized and subjected to gray scale processing by adopting the following formula:
Figure BDA0002183766230000121
b (i, j) is a value obtained by normalizing and graying a (i, j);
i is the abscissa of the biological cavity position r under the polar coordinate system, and j is the ordinate of the biological cavity position r under the polar coordinate system;
a (i, j) is the element value of the ith row and the jth column in the light absorption energy distribution matrix A;
min (A) is the minimum value of the elements in the light absorption energy distribution matrix A, max (A) is the maximum value of the elements in the light absorption energy distribution matrix A;
adopting the following formula to align the sound velocity distribution matrix c under polar coordinates s Carrying out normalization and gray level processing;
Figure BDA0002183766230000122
e (i, j) is for c s (i, j) normalized and grayed values;
c s (i, j) is the sound velocity distribution matrix c s The element value of the ith row and the jth column;
min(c s ) Is the sound velocity distribution matrix c s Minimum value of (1); max (c) s ) Is the sound velocity distribution matrix c s Maximum value of (1);
converting said B (i, j) into B '(x, y) in planar rectangular coordinates according to the formula B' (x, y) = B (jcosi, jsini);
converting said E (i, j) into E '(x, y) in planar cartesian coordinates according to the formula E' (x, y) = E (jcosi, jsini);
wherein i belongs to [0,360], j belongs to [0,d ];
d is the maximum value of the polar diameter in the theta-l polar coordinate system;
b' (x, y) is a gray value of a point (x, y) in a rectangular coordinate system in the light absorption energy distribution map; e' (x, y) is the grayscale value of the point (x, y) in the rectangular coordinate system in the sound velocity map.
The invention also discloses a photoacoustic endoscopic imaging system for the biological cavity, which comprises:
an initial value acquisition module for acquiring an initial value A of light absorption energy at a position r of the biological cavity 0 (r) and initial value of sound velocity c s,0 (r)。
In this embodiment, the initial value of the light absorption energy is A 0 (r) =0, initial value of sound velocity c s,0 (r)＝1600m/s。
A light absorption energy calculation module for calculating the light absorption energy A of the kth time at the position r of the biological cavity respectively k (r) and (k + 1) th light absorption energy A k+1 (r); k is an iterationThe number of times; k is more than or equal to 1.
In this embodiment, the formula A is shown k (r)＝A k-1 (r)-W -1 (A k-1 (r))G 1 '(r,c s,k-1 (r),A k-1 (r)) calculating the light absorption energy A of the k-th time k (r) and the light absorption energy A of the k +1 th order k+1 (r)。
Wherein,
Figure BDA0002183766230000131
argmin[·]is that g (A) k (r))+λΦ TV (A k (r)) A at the minimum k (r)；
W(A k-1 (r)) is an approximate blackplug matrix of the first photoacoustic relationship function;
W -1 (A k-1 (r)) is W (A) k-1 (r)) an inverse matrix;
a is a light absorption energy distribution matrix; a is more than or equal to 0, which means the light absorption energy which is more than or equal to 0 in the light absorption energy distribution matrix;
g(A k (r))＝||p m (r)-H(c s,k-1 (r))·A k (r)|| 2 is a second photoacoustic relationship function;
p m (r) is a measurement of the photoacoustic signal at location r;
H(c s,k-1 (r)) is a first operator related to the speed of sound;
| | · | is a 2-norm;
Φ TV (A k (r)) is a TV regularization term,
Figure BDA0002183766230000132
λ is the TV regularization parameter;
L k is about A k (r) a sparse matrix of variance characteristics;
η > 0 is a constant;
G 1 '(r,A k-1 (r),c s,k-1 (r)) a ladder of said photoacoustic correlation functionAnd (4) degree.
In the specific treatment process, the light energy deposition in the cavity wall tissue under the short pulse laser irradiation is simulated by adopting a Monte Carlo simulation method. Then, solving a discrete photoacoustic wave equation under a polar coordinate system by adopting a finite difference time domain algorithm to obtain a theoretical value of a photoacoustic signal generated by a tissue:
namely, it is
Figure BDA0002183766230000141
Wherein, (i, j) is the coordinate of a point r on the cross section of the cavity in a theta-l polar coordinate system; Δ θ and Δ l are unit lengths on the θ axis and the l axis, respectively; Δ t is a discrete time interval; n is a discrete time; p is a radical of (n+1) (i, j) is the theoretical value of the photoacoustic signal generated by the particle at position (i, j) at time n + 1;
Figure BDA0002183766230000142
and
Figure BDA0002183766230000143
the vibration speeds of the mass point with the position (i, j) along the theta direction and the l direction at the time n respectively; c. C s (i, j) is the ultrasonic wave at position (i, j); the propagation velocity of (c); β is the isobaric expansion coefficient of the tissue; c p Is the specific heat capacity of the tissue; rho 0 Is the density of the tissue; i is (n) Is the value of the laser pulse function at time n; a (i, j) is the light absorption energy at location (i, j).
A sound velocity calculation module for respectively calculating the kth sound velocity c at the position r of the biological cavity s,k (r) and sound velocity c of (k + 1) th order s,k+1 (r)。
In the present embodiment of the present invention,
Figure BDA0002183766230000144
after k iteration is carried out on the third photoacoustic relation function, the sound velocity at the position r is obtained
Figure BDA0002183766230000145
After the solution is carried out, the solution is obtained,
Figure BDA0002183766230000146
according to the formula
Figure BDA0002183766230000147
Respectively calculating the sound velocity c of the k-th time s,k (r) and the sound speed c of the k +1 st time s,k+1 (r)。
Wherein,
Figure BDA0002183766230000151
▽f(q k ) Is f (q) k ) A gradient of (a);
f(c s,k-1 (r))＝||p m (r)-H(c s,k-1 (r))·A k-1 (r)|| 2 ；
Figure BDA0002183766230000152
γ k-1 is the weight parameter, γ, after the k-1 iteration k-2 Is the weight parameter after the k-2 iteration;
weight parameter gamma after the kth iteration k Is composed of
Figure BDA0002183766230000153
H(q k ) Is a second operator related to the speed of sound;
d is f (c) s,k-1 (r)) the prestz constant of the derivative;
f(c s,k-1 (r))＝||p m (r)-H(c s,k-1 (r))·A k-1 (r)|| 2 。
c s sound velocity distribution matrix, c s And > 0 is the sound velocity of 0 or more in the sound velocity distribution matrix.
An absolute difference calculation module for calculating the light absorption energy A of the kth time k (r) and the light absorption energy A of the k +1 th order k+1 Absolute difference epsilon of light absorption energy of (r) A,k And the speed of sound c of the kth time s,k (r) speed of sound c of the k +1 th time s,k+1 Absolute difference in sound velocity of (r)
Figure BDA0002183766230000154
The convergence condition is epsilon A,k ＜ε A And is
Figure BDA0002183766230000155
Wherein epsilon A For the convergence tolerance of the light absorption energy,
Figure BDA0002183766230000156
is the convergence tolerance of the speed of sound.
In particular, the convergence tolerance ε of the light absorption energy A Take 0.01, convergence tolerance of speed of sound
Figure BDA0002183766230000157
Take 0.01.
In the iterative calculation process, if all the iteration values do not meet the convergence condition, the convergence tolerance epsilon of the light absorption energy can be changed A Convergence tolerance to speed of sound
Figure BDA0002183766230000158
And (4) adjusting the convergence condition to complete the iterative process.
A convergence judgment result module for judging the absolute difference epsilon of the light absorption energy A,k And the absolute difference of sound velocity
Figure BDA0002183766230000159
Whether the convergence condition is met or not is judged, and a convergence judgment result is obtained;
the image construction module is connected with the convergence judgment result module and is used for judging whether the convergence judgment result module is yes or not according to the judgment resultThe (k + 1) th light absorption energy A at the position r k+1 (r) and speed of sound c s,k+1 (r) respectively constructing a light absorption energy distribution map and an acoustic velocity distribution map;
and the iteration updating module is respectively connected with the convergence judgment result module and the light absorption energy calculating module, and is used for updating the iteration times and returning to the light absorption energy calculating module when the convergence judgment result module is negative.
The image construction module specifically comprises:
the light absorption energy distribution matrix processing unit is used for normalizing and carrying out gray scale processing on the light absorption energy distribution matrix A under the polar coordinates by adopting the following formula:
Figure BDA0002183766230000161
b (i, j) is a value obtained by normalizing and graying A (i, j);
i is the abscissa of the biological cavity position r under the polar coordinate system, and j is the ordinate of the biological cavity position r under the polar coordinate system;
a (i, j) is the element value of the ith row and the jth column in the light absorption energy distribution matrix A;
min (A) is the minimum value of the elements in the light absorption energy distribution matrix A, max (A) is the maximum value of the elements in the light absorption energy distribution matrix A;
a sound velocity distribution matrix processing unit for aligning the sound velocity distribution matrix c under polar coordinates by using the following formula s Carrying out normalization and gray level processing;
Figure BDA0002183766230000162
e (i, j) is for c s (i, j) normalized and grayed values;
c s (i, j) is the sound velocity distribution matrix c s The value of the element in the ith row and the jth column;
min(c s ) Is the sound velocity distribution matrixc s Minimum value of (1); max (c) s ) Is the sound velocity distribution matrix c s Maximum value of (1);
a light absorption energy distribution map gray value conversion unit for converting the B (i, j) into B '(x, y) in planar rectangular coordinates according to a formula B' (x, y) = B (jcosi, jsini);
a sound velocity distribution image gradation value conversion unit for converting the E (i, j) into E '(x, y) in planar rectangular coordinates according to a formula E' (x, y) = E (jcosi, jsini);
wherein i belongs to [0,360], j belongs to [0,d ];
d is the maximum value of the polar diameter in the theta-l polar coordinate system;
b' (x, y) is a gray value of a point (x, y) in a rectangular coordinate system in the light absorption energy distribution map; e' (x, y) is the grayscale value of the point (x, y) in the rectangular coordinate system in the sound velocity map.
For the system disclosed by the embodiment, the description is relatively simple because the system corresponds to the method disclosed by the embodiment, and the relevant points can be referred to the method part for description.
The principles and embodiments of the present invention have been described herein using specific examples, which are provided only to help understand the method and the core concept of the present invention; meanwhile, for a person skilled in the art, according to the idea of the present invention, the specific embodiments and the application range may be changed. In view of the above, the present disclosure should not be construed as limiting the invention.

Claims (10)

1. A photoacoustic endoscopic imaging method for a biological cavity, the method comprising:
obtaining the initial value A of the light absorption energy at the position r of the biological cavity 0 (r) and initial value of sound velocity c s,0 (r)；
Respectively calculating the light absorption energy A of the kth time at the position r of the biological cavity k (r) and (k + 1) th light absorption energy A k+1 (r); k is iteration times, and k is more than or equal to 1;
respectively calculating the kth sound at the position r of the biological cavitySpeed c s,k (r) and sound velocity c of (k + 1) th order s,k+1 (r)；
Calculating the light absorption energy A of the k time k (r) and the light absorption energy A of the k +1 th order k+1 Absolute difference epsilon of light absorption energy of (r) A,k And the speed of sound c of the kth time s,k (r) speed of sound c at the k +1 th time s,k+1 Absolute difference in sound velocity of (r)
Figure FDA0003919830850000013
Judging the absolute difference epsilon of the light absorption energy A,k And the absolute difference of sound velocity
Figure FDA0003919830850000012
Whether the convergence condition is met or not is judged, and a convergence judgment result is obtained;
if the convergence judgment result is yes, outputting the (k + 1) th light absorption energy A at the position r k+1 (r) and speed of sound c s,k+1 (r) respectively constructing a light absorption energy distribution map and an acoustic velocity distribution map;
if the convergence judgment result is negative, updating the iteration times, and returning to the step of respectively calculating the light absorption energy A of the kth time at the position r of the biological cavity k (r) and (k + 1) th light absorption energy A k+1 (r)。
2. The photoacoustic endoscopic imaging method for the biological cavity of claim 1, wherein formula A is given k (r)＝A k-1 (r)-W -1 (A k-1 (r))G 1 '(r,A k-1 (r),c s,k-1 (r)) calculating the light absorption energy A of the k-th time k (r) and the light absorption energy A of the k +1 th order k+1 (r)；
Wherein,
Figure FDA0003919830850000011
is a first photoacoustic relationship function;
argmin[·]is that g (A) k (r))+λΦ TV (A k (r)) A at the minimum k (r)；
W(A k-1 (r)) is an approximate blackplug matrix of the first photoacoustic relationship function;
W -1 (A k-1 (r)) is W (A) k-1 (r)) an inverse matrix;
a is a light absorption energy distribution matrix; a is more than or equal to 0, which means the light absorption energy which is more than or equal to 0 in the light absorption energy distribution matrix;
g(A k (r))＝||p m (r)-H(c s,k-1 (r))·A k (r)|| 2 is a second photoacoustic relationship function;
p m (r) is a measurement of the photoacoustic signal at location r;
H(c s,k-1 (r)) is a first operator related to the speed of sound;
| | · | is a 2-norm;
Φ TV (A k (r)) is a TV regularization term,
Figure FDA0003919830850000021
λ is the TV regularization parameter;
L k is about A k (r) a sparse matrix of variance characteristics;
η > 0 is a constant;
G 1 '(r,A k-1 (r),c s,k-1 (r)) is the gradient of the first photoacoustic relationship function.
3. The photoacoustic endoscopic imaging method for the biological cavity according to claim 2, wherein the formula is given
Figure FDA0003919830850000022
Respectively calculating the sound velocity c of the k-th time s,k (r) and the sound speed c of the k +1 st time s,k+1 (r)；
Wherein,
Figure FDA0003919830850000023
Figure FDA0003919830850000024
is f (q) k ) A gradient of (a);
f(q k )＝||p m (r)-H(q k )·A k-1 (r)|| 2 ；
Figure FDA0003919830850000025
γ k-1 is the weight parameter, γ, after the k-1 iteration k-2 Is the weight parameter after the k-2 iteration;
weight parameter gamma after the kth iteration k Is composed of
Figure FDA0003919830850000026
H(q k ) Is a second operator related to the speed of sound;
d is f (c) s,k-1 (r)) the prestz constant of the derivative;
f(c s,k-1 (r))＝||p m (r)-H(c s,k-1 (r))·A k-1 (r)|| 2 。
4. the photoacoustic endoscopic imaging method according to claim 1, wherein the convergence condition is ε A,k ＜ε A And is provided with
Figure FDA0003919830850000031
Wherein epsilon A For the convergence tolerance of the light absorption energy,
Figure FDA0003919830850000032
is the convergence tolerance of the speed of sound.
5. The photoacoustic endoscopic imaging method according to claim 1, wherein the constructing the light absorption energy distribution image and the sound velocity distribution image specifically comprises:
the light absorption energy distribution matrix A under the polar coordinates is normalized and subjected to gray scale processing by adopting the following formula:
Figure FDA0003919830850000033
b (i, j) is a value obtained by normalizing and graying a (i, j);
i is the abscissa of the biological cavity position r under the corresponding polar coordinate system, and j is the ordinate of the biological cavity position r under the corresponding polar coordinate system;
a (i, j) is the element value of the ith row and the jth column in the light absorption energy distribution matrix A;
min (A) is the minimum value of the elements in the light absorption energy distribution matrix A, max (A) is the maximum value of the elements in the light absorption energy distribution matrix A;
adopting the following formula to align the sound velocity distribution matrix c under polar coordinates s Carrying out normalization and gray level processing;
Figure FDA0003919830850000034
e (i, j) is for c s (i, j) normalized and grayed values;
c s (i, j) is the sound velocity distribution matrix c s The element value of the ith row and the jth column;
min(c s ) Is the sound velocity distribution matrix c s Minimum value of (1); max (c) s ) Is the sound velocity distribution matrix c s Maximum value of (1);
converting said B (i, j) into B '(x, y) in planar rectangular coordinates according to the formula B' (x, y) = B (jcosi, jsini);
converting said E (i, j) into E '(x, y) in planar rectangular coordinates according to the formula E' (x, y) = E (jcosi, jsini);
wherein, i belongs to [0,360], j belongs to [0,d ];
d is the maximum value of the polar diameter in the theta-l polar coordinate system;
b' (x, y) is a gray value of a point (x, y) in a rectangular coordinate system in the light absorption energy distribution map; e' (x, y) is a gray value of a point (x, y) in a rectangular coordinate system in the sound velocity map.
6. A biological cavity photoacoustic endoscopic imaging system, the system comprising:
an initial value acquisition module for acquiring an initial value A of light absorption energy at a position r of the biological cavity 0 (r) and initial value of sound velocity c s,0 (r)；
A light absorption energy calculation module for respectively calculating the light absorption energy A of the kth time at the position r of the biological cavity k (r) and (k + 1) th light absorption energy A k+1 (r); k is iteration times, and k is more than or equal to 1;
a sound velocity calculation module for respectively calculating the kth sound velocity c at the position r of the biological cavity s,k (r) and sound velocity c of (k + 1) th order s,k+1 (r)；
An absolute difference calculation module for calculating the light absorption energy A of the kth time k (r) and the light absorption energy A of the (k + 1) th order k+1 Absolute difference of light absorption energy ε of (r) A,k And the speed of sound c of the kth time s,k (r) speed of sound c at the k +1 th time s,k+1 (r) absolute difference of sound velocity ε cs , k ；
A convergence judgment result module for judging the absolute difference epsilon of the light absorption energy A,k And the absolute difference of sound velocity ε cs,k Whether the convergence condition is met or not is judged, and a convergence judgment result is obtained;
an image construction module connected with the convergence judgment result module and used for absorbing energy A according to the (k + 1) th light at the position r when the convergence judgment result is yes k+1 (r) and speed of sound c s,k+1 (r) respectively constructing a light absorption energy distribution map and an acoustic velocity distribution map;
and the iteration updating module is respectively connected with the convergence judgment result module and the light absorption energy calculating module, and is used for updating the iteration times and returning the iteration times to the light absorption energy calculating module when the convergence judgment result is negative.
7. The bio-cavity photoacoustic endoscopic imaging system according to claim 6, wherein formula A is given as the basis k (r)＝A k-1 (r)-W -1 (A k-1 (r))G 1 '(r,A k-1 (r),c s,k-1 (r)) calculating the light absorption energy A of the kth time k (r) and the light absorption energy A of the k +1 th order k+1 (r)；
Wherein,
Figure FDA0003919830850000041
is a first photoacoustic relationship function;
argmin[·]is that g (A) k (r))+λΦ TV (A k (r)) A at the minimum k (r)；
W(A k-1 (r)) is an approximate blackplug matrix of the first photoacoustic relationship function;
W -1 (A k-1 (r)) is W (A) k-1 (r)) an inverse matrix;
a is a light absorption energy distribution matrix; a is more than or equal to 0, which means the light absorption energy which is more than or equal to 0 in the light absorption energy distribution matrix;
g(A k (r))＝||p m (r)-H(c s,k-1 (r))·A k (r)|| 2 is a second photoacoustic relationship function;
p m (r) is a measurement of the photoacoustic signal at location r;
H(c s,k-1 (r)) is a first operator related to the speed of sound;
| | · | is a 2-norm;
Φ TV (A k (r)) is a TV regularization term,
Figure FDA0003919830850000051
λ is the TV regularization parameter;
L k is offIn A k (r) a sparse matrix of variance characteristics;
η > 0 is a constant;
G 1 '(r,A k-1 (r),c s,k-1 (r)) is the gradient of the first photoacoustic relationship function.
8. The bio-cavity photoacoustic endoscopic imaging system according to claim 7, wherein the formula is given
Figure FDA0003919830850000052
Respectively calculating the sound velocity c of the k-th time s,k (r) and the speed of sound c after the (k + 1) th iteration s,k+1 (r)；
Wherein,
Figure FDA0003919830850000053
Figure FDA0003919830850000054
is f (q) k ) A gradient of (a);
f(q k )＝||p m (r)-H(q k )·A k-1 (r)|| 2 ；
Figure FDA0003919830850000055
γ k-1 is the weight parameter after the k-1 iteration, γ k-2 Is the weight parameter after the k-2 iteration;
weight parameter gamma after the kth iteration k Is composed of
Figure FDA0003919830850000061
H(q k ) Is a second operator related to the speed of sound;
d is f (c) s,k-1 (r)) the prestz constant of the derivative;
f(c s,k-1 (r))＝||p m (r)-H(c s,k-1 (r))·A k-1 (r)|| 2 。
9. the photoacoustic endoscopic imaging system of claim 6, wherein the convergence condition is ε A,k ＜ε A And is
Figure FDA0003919830850000062
Wherein epsilon A For the convergence tolerance of the light absorption energy,
Figure FDA0003919830850000063
is the convergence tolerance of the speed of sound.
10. The bio-cavity photoacoustic endoscopic imaging system according to claim 6, wherein the image construction module specifically comprises:
the light absorption energy distribution matrix processing unit is used for carrying out normalization and gray scale processing on the light absorption energy distribution matrix A under the polar coordinates by adopting the following formula:
Figure FDA0003919830850000064
b (i, j) is a value obtained by normalizing and graying A (i, j);
i is the abscissa of the biological cavity position r under the corresponding polar coordinate system, and j is the ordinate of the biological cavity position r under the corresponding polar coordinate system;
a (i, j) is the element value of the ith row and the jth column in the light absorption energy distribution matrix A;
min (A) is the minimum value of the elements in the light absorption energy distribution matrix A, max (A) is the maximum value of the elements in the light absorption energy distribution matrix A;
a sound velocity distribution matrix processing unit for aligning the sound velocity distribution matrix c under polar coordinates by using the following formula s Performing normalization and gray scale processingC, trimming;
Figure FDA0003919830850000065
e (i, j) is for c s (i, j) normalized and grayed values;
c s (i, j) is the sound velocity distribution matrix c s The value of the element in the ith row and the jth column;
min(c s ) Is the sound velocity distribution matrix c s Minimum value of (1); max (c) s ) Is the sound velocity distribution matrix c s Maximum value of (2);
a light absorption energy distribution map gray value conversion unit for converting the B (i, j) into B '(x, y) in a plane rectangular coordinate according to a formula B' (x, y) = B (jcosi, jsini);
a sound velocity profile gradation value conversion unit for converting the E (i, j) into E '(x, y) in planar rectangular coordinates according to a formula E' (x, y) = E (jcosi, jsini);
wherein, i belongs to [0,360], j belongs to [0,d ];
d is the maximum value of the polar diameter in a theta-l polar coordinate system;
b' (x, y) is a gray value of a point (x, y) in a rectangular coordinate system in the light absorption energy distribution map; e' (x, y) is a gray value of a point (x, y) in a rectangular coordinate system in the sound velocity map.
CN201910806275.7A 2019-08-29 2019-08-29 Photoacoustic endoscopic imaging method and system for biological cavity Active CN112438702B (en)

Priority Applications (1)

Application Number Priority Date Filing Date Title
CN201910806275.7A CN112438702B (en) 2019-08-29 2019-08-29 Photoacoustic endoscopic imaging method and system for biological cavity

Applications Claiming Priority (1)

Application Number Priority Date Filing Date Title
CN201910806275.7A CN112438702B (en) 2019-08-29 2019-08-29 Photoacoustic endoscopic imaging method and system for biological cavity

Publications (2)

Publication Number Publication Date
CN112438702A CN112438702A (en) 2021-03-05
CN112438702B true CN112438702B (en) 2022-12-27

Family

ID=74741903

Family Applications (1)

Application Number Title Priority Date Filing Date
CN201910806275.7A Active CN112438702B (en) 2019-08-29 2019-08-29 Photoacoustic endoscopic imaging method and system for biological cavity

Country Status (1)

Country Link
CN (1) CN112438702B (en)

Citations (5)

* Cited by examiner, † Cited by third party
Publication number Priority date Publication date Assignee Title
CN102946809A (en) * 2010-02-01 2013-02-27 松下电器产业株式会社 Ultrasonic probe and ultrasonic examination device using same
CN105139382A (en) * 2015-08-03 2015-12-09 华北电力大学（保定） Elastic registration method of intracoronary ultrasonic image sequence
CN105249993A (en) * 2015-11-16 2016-01-20 南京大学 Method for selecting optimum sound velocity group to optimize ultrasonic imaging through photoacoustic imaging
CN107049252A (en) * 2017-03-29 2017-08-18 华北电力大学（保定） A kind of biological magneto-optic sound joint endoscopic imaging method
CN109662720A (en) * 2018-11-27 2019-04-23 中国科学院深圳先进技术研究院 Multiband optoacoustic Woundless blood sugar concentration prediction system based on deep learning

Family Cites Families (4)

* Cited by examiner, † Cited by third party
Publication number Priority date Publication date Assignee Title
JP4739363B2 (en) * 2007-05-15 2011-08-03 キヤノン株式会社 Biological information imaging apparatus, biological information analysis method, and biological information imaging method
US9445786B2 (en) * 2011-11-02 2016-09-20 Seno Medical Instruments, Inc. Interframe energy normalization in an optoacoustic imaging system
JP6257190B2 (en) * 2013-07-09 2018-01-10 キヤノン株式会社 Subject information acquisition device and laser device
JP6685555B2 (en) * 2015-02-12 2020-04-22 国立大学法人京都工芸繊維大学 Ultrasonic particle size measuring device and ultrasonic measuring device

Patent Citations (5)

* Cited by examiner, † Cited by third party
Publication number Priority date Publication date Assignee Title
CN102946809A (en) * 2010-02-01 2013-02-27 松下电器产业株式会社 Ultrasonic probe and ultrasonic examination device using same
CN105139382A (en) * 2015-08-03 2015-12-09 华北电力大学（保定） Elastic registration method of intracoronary ultrasonic image sequence
CN105249993A (en) * 2015-11-16 2016-01-20 南京大学 Method for selecting optimum sound velocity group to optimize ultrasonic imaging through photoacoustic imaging
CN107049252A (en) * 2017-03-29 2017-08-18 华北电力大学（保定） A kind of biological magneto-optic sound joint endoscopic imaging method
CN109662720A (en) * 2018-11-27 2019-04-23 中国科学院深圳先进技术研究院 Multiband optoacoustic Woundless blood sugar concentration prediction system based on deep learning

Non-Patent Citations (2)

* Cited by examiner, † Cited by third party
Title
A computer-based simulator for intravascular photoacoustic images;Sun Zheng, Yuan Yuan, Han Duoduo;《Computers in Biology and Medicine》;20171231;全文 *
生物声学成像中声速不均匀性解决方法的研究进展;孙正，贾艺璇;《声学技术》;20181031;第37卷(第5期);全文 *

Also Published As

Publication number Publication date
CN112438702A (en) 2021-03-05

Similar Documents

Publication Publication Date Title
Karamalis et al. Fast ultrasound image simulation using the westervelt equation
US8197412B2 (en) Ultrasonic diagnostic apparatus
KR101610874B1 (en) Module for Processing Ultrasonic Signal Based on Spatial Coherence and Method for Processing Ultrasonic Signal
CN105249993A (en) Method for selecting optimum sound velocity group to optimize ultrasonic imaging through photoacoustic imaging
CN110840484B (en) Ultrasonic imaging method and device for adaptively matching optimal sound velocity and ultrasonic equipment
CN107970042B (en) Ultrasonic nonlinear quantitative real-time imaging method and system
WO2014055410A1 (en) Shear wave attenuation from k-space analysis system
CN111956180A (en) Method for reconstructing photoacoustic endoscopic tomography image
CN106940883B (en) Ultrasonic imaging method based on point spread function simulation and compressed sensing of ultrasonic system
CN106373103A (en) Ultrasonic data compounding method and apparatus
CN110536647A (en) For handling system, method and computer-readable medium with composite ultraphonic image in the presence of motion
CN112438702B (en) Photoacoustic endoscopic imaging method and system for biological cavity
Treeby et al. Equivalent-source acoustic holography for projecting measured ultrasound fields through complex media
CN107204021B (en) Ultrasonic imaging method based on Gaussian function probe response model and compressed sensing
CN101313858B (en) Inhibition method for complete type-B ultrasonic instrument supersonic wave aberration based on transmission iterative method
CN112927145B (en) Correction method for planar ultrasonic craniocerebral imaging
Pan et al. Improving axial resolution based on the deconvolution recovery method combined with adaptive weighting techniques for ultrasound imaging
CN109671130B (en) Method and system for reconstructing endoscopic photoacoustic tomography image by using sparse measurement data
Horeh et al. Regularized tracking of shear-wave in ultrasound elastography
CN106997045B (en) Ultrasonic imaging method based on ultrasonic system point spread function measurement and compressed sensing
CN111820868A (en) Biological photoacoustic endoscopic image reconstruction method and system
CN111481168B (en) Photoacoustic endoscopic imaging image reconstruction method and system
CN103126720B (en) The ultrasonic imaging method of many spread speeds and device
Jaeger et al. Learnt correction for regularization-related biases in pulse-echo speed-of-sound imaging
CN116482035B (en) Photoacoustic tomography method and device based on flexible ultrasonic probe

Legal Events

Date Code Title Description
PB01 Publication
PB01 Publication
SE01 Entry into force of request for substantive examination
SE01 Entry into force of request for substantive examination
GR01 Patent grant
GR01 Patent grant