A multi-dimensional feature fusion photoacoustic super-resolution imaging method, device and equipment

Abstract

The application provides a multi-dimensional feature fusion photoacoustic super-resolution imaging method, device and equipment, which comprises the following steps: irradiating a target biological tissue based on a multi-waveband laser beam, exciting to generate an original photoacoustic signal, and collecting multi-dimensional photoacoustic signals based on a multi-modal detector array; performing time-frequency analysis on the multi-dimensional photoacoustic signals, extracting energy distribution characteristics of the multi-dimensional photoacoustic signals in different frequency bands, and obtaining a time-frequency characteristic map; calculating spatial position coordinates of a photoacoustic signal source corresponding to the target biological tissue based on a propagation time difference of the original photoacoustic signal received by the detector and a position distribution of each detector in the detector array, constructing a spatial distribution model of the photoacoustic signal source, and obtaining spatial characteristics; fusing the time-frequency characteristic map and the spatial characteristics to obtain a multi-dimensional feature vector; and performing super-resolution image reconstruction based on the multi-dimensional feature vector to obtain a super-resolution photoacoustic image corresponding to the target biological tissue. The application improves the imaging resolution.

Description

TECHNICAL FIELD

[0001] The present application relates to the technical field of computer, and particularly relates to a multi-dimensional feature fusion photoacoustic super-resolution imaging method, device and equipment. BACKGROUND

[0002] In the field of photoacoustic imaging, an existing photoacoustic method is widely used in multiple fields such as biomedicine, and can obtain the structure and function information inside the biological tissue through photoacoustic effect. However, the method has a major drawback, that is, the imaging resolution is limited. The traditional photoacoustic imaging is limited by multiple factors such as light scattering, sound wave propagation characteristics and detector performance, and when imaging the fine structure inside the biological tissue, it is difficult to clearly distinguish the structure characteristics with small size, such as biological structure smaller than a certain scale (such as 200 nm), which leads to that the imaging result cannot meet the demand of high-precision imaging, and affects the observation and analysis of the fine structure inside the biological tissue, and limits the application in application scenarios with extremely high resolution requirements such as early disease diagnosis and micro-biological process research. SUMMARY

[0003] The present application provides a multi-dimensional feature fusion photoacoustic super-resolution imaging method, device and equipment, which aims to effectively improve the imaging resolution, so that the imaging result can clearly show the fine structure inside the biological tissue, and meet the demand of high-precision imaging.

[0004] In a first aspect, the present application provides a multi-dimensional feature fusion photoacoustic super-resolution imaging method, comprising:

[0005] Irradiating a target biological tissue based on a multi-band laser beam to excite original photoacoustic signals, and collecting the original photoacoustic signals with different frequency ranges and different propagation directions based on a multi-modal detector array to obtain multi-dimensional photoacoustic signals;

[0006] Performing time-frequency analysis on the multi-dimensional photoacoustic signals to convert the time-domain photoacoustic signals to the frequency domain, extracting the energy distribution characteristics of the multi-dimensional photoacoustic signals in different frequency bands to obtain a time-frequency feature map;

[0007] Based on the propagation time difference of the original photoacoustic signals received by the detector and the position distribution of each detector in the detector array, calculating the spatial position coordinates of the photoacoustic signal source corresponding to the target biological tissue to construct a spatial distribution model of the photoacoustic signal source, and obtaining a spatial feature;

[0008] Fusing the time-frequency feature map and the spatial feature to obtain a multi-dimensional feature vector;

[0009] Based on the multi-dimensional feature vector, a super-resolution image is reconstructed to obtain a super-resolution photoacoustic image corresponding to the target biological tissue.

[0010] In a second aspect, the present application also provides a multi-dimensional feature fusion photoacoustic super-resolution imaging device, which is applied to the multi-dimensional feature fusion photoacoustic super-resolution imaging method as described in the first aspect; the multi-dimensional feature fusion photoacoustic super-resolution imaging device comprises:

[0011] a signal acquisition module, configured to irradiate a target biological tissue based on a multi-band laser beam, excite to generate an original photoacoustic signal, and acquire the original photoacoustic signal in different frequency ranges and different propagation directions based on a multi-modal detector array, to obtain a multi-dimensional photoacoustic signal;

[0012] a frequency domain analysis module, configured to perform time-frequency analysis on the multi-dimensional photoacoustic signal, to convert the time-domain photoacoustic signal to a frequency domain, extract energy distribution characteristics of the multi-dimensional photoacoustic signal in different frequency bands, and obtain a time-frequency feature map;

[0013] a spatial feature extraction module, configured to calculate spatial position coordinates of a photoacoustic signal source corresponding to the target biological tissue based on a propagation time difference of the original photoacoustic signal received by the detector and a position distribution of each detector in the detector array, to construct a spatial distribution model of the photoacoustic signal source, and obtain spatial features;

[0014] a feature fusion module, configured to fuse the time-frequency feature map and the spatial features, to obtain a multi-dimensional feature vector;

[0015] an image reconstruction module, configured to perform super-resolution image reconstruction based on the multi-dimensional feature vector, to obtain a super-resolution photoacoustic image corresponding to the target biological tissue.

[0016] In a third aspect, the present application also provides an electronic device, comprising: a memory, configured to store a computer software program; and a processor, configured to read and execute the computer software program, to implement the multi-dimensional feature fusion photoacoustic super-resolution imaging method as described in any one of the above aspects.

[0017] In a fourth aspect, the present application also provides a non-transitory computer readable storage medium, wherein the storage medium stores a computer software program, and the computer software program is executed by a processor to implement the multi-dimensional feature fusion photoacoustic super-resolution imaging method as described in any one of the above aspects.

[0018] In a fifth aspect, the present application also provides a computer program product, comprising a computer program, and the computer program is executed by a processor to implement the multi-dimensional feature fusion photoacoustic super-resolution imaging method as described in any one of the above aspects.

[0019] The multi-dimension feature fusion photoacoustic super-resolution imaging method provided by the embodiment of the present application collects photoacoustic signals from multiple dimensions through a multi-modal detector array, makes up for the deficiencies of a single detector in the frequency response range and the spatial detection angle, and obtains more comprehensive signal data. The time-frequency analysis extracts the energy distribution characteristics of the signal in different frequency bands, which can reflect the dynamic characteristics of different components and structures in the tissue, and the spatial feature extraction clearly defines the position information of the signal source. After the fusion of the two, the multi-dimension feature vector set contains more abundant and comprehensive tissue information. Therefore, in the image reconstruction, based on the fusion features, the correlation between the features can be used to more accurately enhance and repair the details in the image. For example, through the correlation between the spatial position and the frequency characteristics, the position and direction of the fine blood vessels in the tumor tissue can be accurately identified and highlighted, thereby effectively improving the imaging resolution and solving the problem of limited imaging resolution of the existing photoacoustic method. The imaging result can clearly show the fine structure inside the biological tissue, and meet the high-precision imaging requirement. BRIEF DESCRIPTION OF DRAWINGS

[0020] Figure 1 is a flowchart of the multi-dimension feature fusion photoacoustic super-resolution imaging method provided by the embodiment of the present application;

[0021] Figure 2 is a structural schematic diagram of the multi-dimension feature fusion photoacoustic super-resolution imaging device provided by the embodiment of the present application;

[0022] Figure 3 is an embodiment diagram of the electronic device provided by the embodiment of the present application;

[0023] Figure 4 is an embodiment diagram of the computer readable storage medium provided by the embodiment of the present application. DETAILED DESCRIPTION

[0024] The technical solutions in the embodiments of the present application will be clearly and completely described below with reference to the drawings in the embodiments of the present application. Obviously, the described embodiments are only a part of the embodiments of the present application, rather than all the embodiments. Based on the embodiments in the present application, all other embodiments obtained by those skilled in the art without creative labor fall within the scope of protection of the present application.

[0025] In the description of the present application, the terms "first" and "second" are used only for the purpose of description, and cannot be understood as indicating or implying relative importance or implicitly indicating the number of the indicated technical features. Therefore, the features defined with "first" and "second" can explicitly or implicitly include one or more of the features. In the description of the present application, the meaning of "multiple" is two or more, unless otherwise specifically limited.

[0026] In the description of the present application, the term "for example" is used to mean "serving as an example, instance, or illustration." Any embodiment described as "for example" in the present application is not necessarily to be construed as preferred or advantageous over other embodiments. The following description is presented to enable any person skilled in the art to make and use the application. In the following description, for the purposes of explanation, numerous specific details are set forth in order to provide a thorough understanding of the present application. It will be apparent, however, to one skilled in the art that the present application can be practiced without using these specific details. In other instances, well-known structures and processes are not elaborated upon in order to avoid unnecessary detail and so as to not obscure the present application. Thus, the present application is not intended to be limited by the embodiments shown, but is to be accorded the widest scope consistent with the principles and features disclosed herein.

[0027] Referring to Figure 1 shown, Figure 1 is a flowchart of a multi-dimensional feature fusion photoacoustic super-resolution imaging method provided by the present application. The execution subject of the multi-dimensional feature fusion photoacoustic super-resolution imaging method in the present application is an imaging device. Therefore, the multi-dimensional feature fusion photoacoustic super-resolution imaging method includes:

[0028] Step 10, irradiating the target biological tissue based on the multi-band laser beam, exciting to generate original photoacoustic signals, and collecting the original photoacoustic signals of different frequency ranges and different propagation directions based on the multi-modal detector array, to obtain multi-dimensional photoacoustic signals.

[0029] Optionally, the imaging system in the present application includes a multi-band laser and a multi-modal detector array. Therefore, the multi-band laser beam irradiates the target biological tissue, and the tissue generates thermal-elastic expansion after absorbing light energy, and then excites original photoacoustic signals. The multi-modal detector array is composed of multiple types of detectors, which can collect original photoacoustic signals of different frequency ranges and propagation directions. The detectors convert the pressure changes of the photoacoustic signals into electrical signals according to the piezoelectric effect principle, thereby obtaining multi-dimensional photoacoustic signals.

[0030] In an embodiment, the imaging system uses a multi-band laser capable of emitting 532 nm, 808 nm and 1064 nm wavelength lasers as a light source. The laser is focused by an optical system and irradiates a mouse liver tissue sample placed on an imaging platform. The detector array is composed of piezoelectric ceramic detectors of different bandwidths, and the detectors are distributed in a semi-spherical shape around the tissue sample for omnidirectional collection of photoacoustic signals.

[0031] For example, a group of detectors with a center frequency of 1MHz is responsible for collecting photoacoustic signals in the low frequency range, and a group of detectors with a center frequency of 10MHz collects photoacoustic signals in the high frequency range. When the laser irradiates the liver tissue, the tissue absorbs the laser energy to generate photoacoustic signals, and the detectors at different positions and frequency responses convert these signals into electrical signals and transmit them to the signal acquisition system, and finally obtain multi-dimensional photoacoustic signals containing information of different frequency ranges and propagation directions.

[0032] Step 20, time-frequency analysis is performed on the multi-dimensional photoacoustic signal to convert the time-domain photoacoustic signal to the frequency domain, extract the energy distribution characteristics of the multi-dimensional photoacoustic signal in different frequency bands, and obtain a time-frequency feature map.

[0033] Further, the time-frequency analysis aims to convert the time-domain photoacoustic signal to the frequency domain to extract the energy distribution characteristics of different frequency bands, wherein the time-frequency analysis method includes short-time Fourier transform (STFT), wavelet transform, etc. Therefore, the imaging system performs time-frequency analysis on the multi-dimensional photoacoustic signal to convert the time-domain photoacoustic signal to the frequency domain, so that the energy distribution of the photoacoustic signal with time and frequency can be obtained, and a time-frequency feature map is generated, as described in steps 201 to 204.

[0034] Step 30, based on the propagation time difference of the original photoacoustic signal received by the detector and the position distribution of each detector in the detector array, the spatial position coordinates of the photoacoustic signal source corresponding to the target biological tissue are calculated to construct a spatial distribution model of the photoacoustic signal source, and the spatial feature is obtained.

[0035] Further, the imaging system obtains the propagation time difference of the original photoacoustic signal received by the detector and the position distribution of each detector in the detector array, and calculates the spatial position coordinates of the photoacoustic signal source corresponding to the target biological tissue by using a positioning algorithm according to the propagation time difference of the original photoacoustic signal received by the detector and the position distribution of each detector in the detector array. Further, the imaging system can obtain the spatial position information of the photoacoustic signal source in the target biological tissue, i.e. the spatial feature, by constructing a spatial distribution model of the photoacoustic signal source, as described in steps 301 to 304.

[0036] Step 40, fusion based on the time-frequency feature map and the spatial feature to obtain a multi-dimensional feature vector.

[0037] Further, the imaging system fuses the frequency energy information contained in the time-frequency feature map obtained in step 20 with the spatial feature (spatial position coordinates of the photoacoustic signal source) obtained in step 30, wherein the embodiment of the present application can use feature splicing, weighted fusion and other methods to combine different types of features into a multi-dimensional feature vector, so as to be used for subsequent image reconstruction and other tasks. The feature vector can more comprehensively describe the photoacoustic characteristics of the target biological tissue.

[0038] In an embodiment, for mouse liver tissue, the energy mean value of each frequency segment in the time-frequency feature map is taken as a set of characteristic values, for example, obtaining the energy mean value characteristic values of 10 frequency segments. At the same time, the spatial coordinates (x, y, z) of the photoacoustic signal source are taken as another set of characteristic values. The two sets of characteristic values are spliced together in order to form a multi-dimensional feature vector. For example, if the time-frequency energy mean value characteristic values are [0.1, 0.3, 0.5, 0.7, 0.9, 1.1, 1.3, 1.5, 1.7, 1.9], and the spatial coordinates of a photoacoustic signal source are (2, 3, 4), then the fused multi-dimensional feature vector is [0.1, 0.3, 0.5, 0.7, 0.9, 1.1, 1.3, 1.5, 1.7, 1.9, 2, 3, 4].

[0039] Step 50, based on the multi-dimensional feature vector, super-resolution image reconstruction is performed to obtain the super-resolution photoacoustic image corresponding to the target biological tissue.

[0040] Further, the imaging system generates the super-resolution photoacoustic image corresponding to the target biological tissue using the super-resolution image reconstruction algorithm according to the multi-dimensional feature vector, improves the resolution of the photoacoustic image, and more clearly displays the internal structure and details of the target biological tissue, as described in steps 501 to 504.

[0041] The present application implements the multi-modal detector array to collect photoacoustic signals from multiple dimensions, which makes up for the deficiencies of a single detector in the frequency response range and the spatial detection angle, and obtains more comprehensive signal data. The time-frequency analysis extracts the energy distribution characteristics of the signal at different frequency segments, which can reflect the dynamic characteristics of different components and structures in the tissue, and the spatial feature extraction clearly defines the position information of the signal source. After the fusion of the two, the multi-dimensional feature vector set contains more comprehensive and rich tissue information. Therefore, in the image reconstruction, based on the fused features, the correlation between the features can be used to more accurately enhance and repair the details in the image. For example, through the correlation between the spatial position and the frequency characteristics, the position and direction of the fine blood vessels in the tumor tissue can be accurately identified and highlighted, thereby effectively improving the imaging resolution, making the imaging result clearly display the fine structure inside the biological tissue, and meeting the high-precision imaging requirement.

[0042] In an embodiment, steps 201 to 204 are described as follows:

[0043] Step 201, block the multi-dimensional photoacoustic signal in time sequence, and periodically expand each block signal to obtain the expanded signal.

[0044] Optionally, the multi-dimensional photoacoustic signal is a continuous signal varying with time, therefore, the imaging system divides the signal into blocks in time sequence, i.e. cutting the continuous signal into multiple shorter segments, each segment is called a block signal. Further, the imaging system periodically extends each block signal, i.e. repeating the end part of each block signal to the beginning to form a periodic signal.

[0045] In an embodiment, the multi-dimensional photoacoustic signal of mouse liver tissue is collected, with a time length of 1 second and a sampling frequency of 10 MHz. The signal is divided into 10000 blocks with each block having 1000 sampling points. For one block signal, for example, the length of the signal is 1000 sampling points, and the value is [x1, x2, …, x1000]. When performing periodic extension, the last 200 sampling points [x801, x802, …, x1000] of the signal are copied and spliced to the beginning of the signal to obtain the extended signal [x801, x802, …, x1000, x1, x2, …, x1000].

[0046] Step 202, determining the time-frequency window width based on the local variation characteristics of each extended signal, and constructing an adaptive time-frequency window based on the time-frequency window width and the length and instantaneous frequency estimation value at the preset frequency.

[0047] Further, the time-frequency window width in the embodiment of the present application is determined according to the local variation characteristics of the extended signal, wherein the time-frequency window width is set smaller in the region where the signal changes rapidly to better capture high-frequency details, and the time-frequency window width is set larger in the region where the signal changes slowly. The length and instantaneous frequency estimation value at the preset frequency are used to further adjust the shape and position of the time-frequency window. The instantaneous frequency estimation can use methods such as Hilbert transform. For example, the signal is s(t), its Hilbert transform is H{s(t)}, and the instantaneous frequency ω(t) = d / dt[arctan(H{s(t)} / s(t))]. Therefore, the imaging system constructs an adaptive time-frequency window according to the time-frequency window width and the length and instantaneous frequency estimation value at the preset frequency, and the shape and width of the adaptive time-frequency window will dynamically change with the local characteristics of the signal.

[0048] In an embodiment, for the above-mentioned periodic extension of the mouse liver photoacoustic signal block, first calculate its local variation characteristics. By calculating the sum of the absolute values of the differences between adjacent sampling points, if the sum of the absolute values of the differences is large in a certain signal, it indicates that the signal changes sharply. For example, the signal changes sharply in the first half and gently in the second half. The instantaneous frequency estimate ω(t) is calculated by Hilbert transform. For the preset frequency f0, the length of the time-frequency window at this frequency is determined according to the empirical formula L=k / f0 (k is a constant, k=10 here). According to the local variation of the signal and the instantaneous frequency, the time-frequency window width is set to 50 sampling points in the first half and 200 sampling points in the second half, thereby constructing an adaptive time-frequency window. The time-frequency window is narrower in the place where the signal changes sharply and wider in the place where the signal changes gently, and its position is dynamically adjusted according to the instantaneous frequency.

[0049] Step 203, modulate each extended signal based on the adaptive time-frequency window of each extended signal to obtain a modulated signal corresponding to each extended signal.

[0050] Further, the imaging system modulates each extended signal using the constructed adaptive time-frequency window, wherein the modulation process is to multiply the time-frequency window function with the extended signal, so that the part of the signal within the time-frequency window range is highlighted, and other parts are suppressed. In this way, the features of the signal within the local time and frequency range covered by the time-frequency window can be extracted. For example, the extended signal is x(t), the time-frequency window function is w(t), and the modulated signal y(t)=x(t)*w(t).

[0051] Taking the above mouse liver photoacoustic signal block with the constructed adaptive time-frequency window as an example, the time-frequency window function w(t) has a width of 50 sampling points in the first half and a width of 200 sampling points in the second half. Multiply the time-frequency window function with the extended signal. For a time t1 in the first half of the signal, the modulated signal y(t1)=x(t1)*w(t1), at this time, since the time-frequency window is narrow, only the signal features in the narrow frequency range near the time are highlighted. At a time t2 in the second half, the multiplication operation is also performed, and since the time-frequency window is wide, the signal features in the wide frequency range near the time t2 are highlighted. Through such modulation, the modulated signal corresponding to the extended signal is obtained, which contains the time-frequency features of the signal in different local regions.

[0052] Step 204, determine the time-frequency feature map based on the energy distribution characteristics of each modulated signal in different frequency bands.

[0053] Further, the imaging system determines the time-frequency feature map according to the energy distribution characteristics of each modulated signal in different frequency bands, as described in steps 2041 to 204.

[0054] The embodiment of the present application can obtain a time-frequency feature map accurately reflecting the time-frequency energy distribution of a multi-dimensional photoacoustic signal, provide rich and accurate time-frequency domain information for subsequent analysis of the structure and function of a target biological tissue, help to more deeply understand the response characteristics of the biological tissue to multi-band laser, for example, clearly distinguish the photoacoustic signal energy difference of different structures (such as blood vessels, parenchymal cells, etc.) in the liver tissue at different frequencies and times, thereby effectively improving the imaging resolution, enabling the imaging result to clearly show the fine structure inside the biological tissue, and meeting the high-precision imaging requirement.

[0055] In an embodiment, steps 2041 to 204 are described as follows:

[0056] Step 2041: performing generalized time-frequency transformation on each modulated signal to obtain a generalized time-frequency transformation result, and determining the energy density of different frequency bands on the time-frequency plane based on each generalized time-frequency transformation result.

[0057] Optionally, the generalized time-frequency transformation is a method that can more flexibly convert a time-domain signal to a time-frequency domain, and can better adapt to the complex characteristics of the signal compared to the traditional time-frequency transformation (such as short-time Fourier transformation). The generalized time-frequency transformation is performed on each modulated signal to obtain the distribution of the signal on the time-frequency plane. Taking the Wigner-Ville distribution as an example, the modulated signal is y(t), and the WVD is defined as:

[0058]

[0059] wherein y * (t) is the conjugate of y(t), and the generalized time-frequency transformation result W y (t, f) is obtained through the transformation. Then, the energy density of different frequency bands on the time-frequency plane is determined. For a given frequency band [f1, f2], the energy density E d (t) can be obtained by integrating the generalized time-frequency transformation result in the frequency band:

[0060]

[0061] In an embodiment, the imaging system obtains a modulated photoacoustic signal y(t) of a mouse liver tissue, and uses the built-in algorithm module to perform Wigner-Ville distribution transformation on y(t). For example, the sampling time range is 0 to T seconds, and the sampling interval is Δt. In the calculation of W y(t,f) is calculated according to the WVD formula above. For numerical calculation convenience, the approximation calculation is performed in a discrete form. For example, the frequency range is divided into frequency segments of 0-1 MHz, and the |W y (t,f)| 2 The integration in the frequency dimension is performed to obtain the energy density E d (t) of the frequency segment at different time points t.

[0062] Step 2042, enhancing the energy density of each frequency segment based on a nonlinear logarithm to obtain an enhanced energy density.

[0063] Further, in order to highlight the difference between the energy densities of different frequency segments, the energy density of each frequency segment is enhanced by a nonlinear logarithm. The logarithmic transformation can amplify the smaller energy density values and relatively compress the larger energy density values, so that the originally weak energy changes can be more obviously exhibited in the subsequent atlas. For example, the original energy density E d (t) is enhanced to an enhanced energy density E wherein a is a constant for adjusting the enhancement degree, which can be selected according to the actual signal condition. Continuing the above example, for the energy density E d (t) of the 0-1 MHz frequency segment at different time points t calculated above, the imaging system selects a = 100. For the E d (t) value of each time point, the calculation is performed according to the logarithmic enhancement formula above. For example, at a time point t0, the original energy density E d (t0) = 0.01, and after the logarithmic enhancement, E Through such enhancement operation, the part with smaller energy density in the frequency segment is amplified.

[0064] Step 2043, taking time as the horizontal coordinate and taking the frequency segment index as the vertical coordinate, generating an initial feature atlas by taking the enhanced energy density of each frequency segment as the grayscale value or color value of the atlas.

[0065] Further, a map is constructed with time as the horizontal coordinate and frequency band index as the vertical coordinate. The frequency band index represents different frequency bands, for example, the frequency range is divided into N frequency bands, and the frequency band index ranges from 1 to N. The enhanced energy density on each frequency band is taken as the gray value or color value of the corresponding position of the map. If gray values are used to represent, the larger the gray value, the stronger the energy of the frequency band at the corresponding time; if color values are used to represent, a color mapping table can be defined in advance, different energy density values correspond to different colors, and the stronger the energy, the brighter the color or the more biased to a high-energy representative color.

[0066] In an embodiment, the imaging system divides the frequency range of the mouse liver photoacoustic signal into 10 frequency bands, and the frequency band index ranges from 1 to 10. The time range of the horizontal coordinate is 0 to 1 second, which is also discretized, for example, the time interval is Δt = 0.01 second. For the frequency band index of 3 (for example, corresponding to the 2-3 MHz frequency band), at the time point t = 0.2 second, the enhanced energy density is 0.8. If gray values are used to represent, according to the pre-set gray value range (for example, 0-1 corresponds to black to white), 0.8 is mapped to the corresponding gray value, and the gray value is set at the corresponding position (horizontal coordinate 0.2 second, vertical coordinate frequency band index 3) in the map. If color representation is used, according to the color mapping table, for example, energy density in 0-0.5 corresponds to blue, 0.5-1 corresponds to red, and the energy density value of 0.8 will make the position show a color biased to red, thereby generating an initial feature map, which preliminarily shows the distribution of energy density of different frequency bands over time.

[0067] Step 2044, for each feature point in the initial feature map, the gray value or color value of each feature point is updated based on the average gray value or average color value in the local neighborhood centered on each feature point, to generate a time-frequency feature map.

[0068] Further, for each feature point in the initial feature map, the gray value or color value of the feature point is updated considering the average gray value or average color value in its local neighborhood. The local neighborhood is a small area centered on the feature point, for example, a 3*3 or 5*5 neighborhood. In this way, noise can be smoothed, and the continuity and stability of the feature can be enhanced. Let the gray value of the feature point (i, j) in the initial feature map be G ij , and its local neighborhood be N ij , the updated gray value G ij new is:

[0069] where |N ij | is the local neighborhood N ijThe number of feature points contained in the local neighborhood. If color values are used, the method of calculating the average color value in the local neighborhood is similar, which can be calculated by averaging the color values in the RGB color space in the neighborhood.

[0070] In an embodiment, the initial feature map generated by the imaging system is represented by gray values, and for a feature point with a position of (10, 5) in the map (for example, the horizontal coordinate i = 10 corresponds to a time t = 0.1 seconds, and the vertical coordinate j = 5 corresponds to a frequency band index of 5), its local neighborhood is set to 3*3. The local neighborhood contains 9 feature points (9, 4), (9, 5), (9, 6), (10, 4), (10, 5), (10, 6), (11, 4), (11, 5), (11, 6), and their gray values are G 9,4 = 0.6, G 9,5 = 0.7, G 9,6 = 0.5, G 10,4 = 0.8, G 10,5 = 0.9, G 10,6 = 0.7, G 11,4 = 0.6, G 11,5 = 0.8, G 11,6 = 0.7. According to the above formula, the updated gray value G 10,5 new = (0.6+0.7+0.5+0.8+0.9+0.7+0.6+0.8+0.7) / 9≈0.7. The update operation is performed on all feature points in the initial feature map, and finally the time-frequency feature map is generated.

[0071] The embodiment of the present application can generate a high-quality time-frequency feature map, so that the time-frequency feature map can clearly and accurately reflect the energy distribution of the multi-dimensional photoacoustic signal of the mouse liver tissue at different frequency bands and times, providing more reliable and more recognizable time-frequency domain information for subsequent analysis of the characteristics of biological tissues, thereby effectively improving the imaging resolution, so that the imaging result can clearly display the fine structure inside the biological tissue, and meet the high-precision imaging requirement.

[0072] In an embodiment, steps 301 to 304 are described as follows:

[0073] Step 301, based on the propagation time difference of the original photoacoustic signal received by the detector and the propagation speed of the photoacoustic signal in the biological tissue, the corresponding propagation distance difference when the signal is received by different detectors is determined.

[0074] Optionally, the acoustic signal propagates in the biological tissue at a certain speed. When the original photoacoustic signal is received by the detectors, due to the different distances from the signal source to different detectors, a propagation time difference is generated. Given the propagation speed of the photoacoustic signal in the biological tissue as v, and the time difference of the signal received by detector A and detector B as Δt, the corresponding propagation distance difference Δd can be calculated by the formula Δd = v * Δt, wherein the propagation distance difference reflects the distance difference information from the photoacoustic signal source to different detectors.

[0075] Step 302, determine the detector position vector of each detector in space based on the position distribution of each detector in the detector array, and construct an initial distance difference equation about the position coordinates of the photoacoustic signal source based on the propagation distance difference and the detector position vector.

[0076] Further, each detector in the detector array has a certain position in space. The position coordinates of detector i are (x i ,y i ,z i ), and its position vector is represented as The position coordinates of the photoacoustic signal source are (x, y, z), and the distance from the photoacoustic signal source to detector i is obtained by the distance formula between two points.

[0077] Further, combining the propagation distance difference and the detector position vector, for two detectors A and detector B, an initial distance difference equation can be constructed:

[0078] Step 303, optimize the initial distance difference equation based on the triangular geometric relationship between the detector and the photoacoustic signal source to obtain the target distance difference equation.

[0079] Further, considering the triangular geometric relationship between the detector and the photoacoustic signal source, the initial distance difference equation is optimized by using some geometric properties and mathematical transformations. For example, the square operation is used to remove the square root, and then the terms are moved and simplified appropriately. Taking two detectors A and detector B as an example, the initial distance difference equation is operated as follows:

[0080]

[0081] Square both sides to get: (x-x A ) 2 +(y-y A ) 2 +(z-z A ) 2 =(Δd) 2 +2Δd

[0082]

[0083] Then the terms are moved, the same terms are combined, and the terms containing the square root are squared again, and after a series of complex algebraic operations, the optimized target distance difference equation is obtained. The equation is relatively more concise in form, and is more convenient for subsequent solving of the photoacoustic signal source position coordinates.

[0084] In step 304, the spatial position coordinates of the photoacoustic signal source corresponding to the target biological tissue are calculated based on the target distance difference equation, to construct a spatial distribution model of the photoacoustic signal source, and to obtain the spatial characteristics.

[0085] Further, the spatial position coordinates of the photoacoustic signal source corresponding to the target biological tissue are calculated according to the target distance difference equation, to construct a spatial distribution model of the photoacoustic signal source, and to obtain the spatial characteristics, which are specifically described in steps 3041 to 3043.

[0086] The embodiment of the present application can accurately construct a spatial distribution model of the photoacoustic signal source in the target biological tissue, and obtain accurate spatial characteristics, so as to intuitively present the spatial position of the photoacoustic signal source in the biological tissue, help researchers clearly understand the position distribution (position and direction) of the photoacoustic signal generated by different structures (such as blood vessels, lesion areas, etc.) in the biological tissue, thereby effectively improving the imaging resolution, making the imaging result clearly show the fine structure inside the biological tissue, and meeting the high-precision imaging requirement.

[0087] In an embodiment, steps 3041 to 3043 are described as follows:

[0088] In step 3041, a plurality of target distance difference equations are combined to obtain an equation group, and the offset of the photoacoustic signal source position coordinates relative to the reference detector is obtained by solving the equation group.

[0089] Optionally, in step 303, a plurality of target distance difference equations are obtained, which describe the distance relationship of the photoacoustic signal source to different detectors. These equations are combined to form an equation group. Since the unknowns in the equation group are the position coordinates of the photoacoustic signal source, by solving the equation group, the offset of the photoacoustic signal source position coordinates relative to the reference detector can be obtained. For example, the position coordinates of the reference detector are (x ref ,y ref ,z ref ), the position coordinates of the photoacoustic signal source are (x, y, z), then the offset is Δx=x-x ref , Δy=y-y ref , and Δz=z-z ref . The iterative algorithm such as Newton-Raphson iteration method can be used to solve the equation group. For an equation group F(Δx, Δy, Δz)=0 containing n target distance difference equations, where F=[f1, f2,..., fn], and the initial value of the offset is (Δx0, Δy0, Δz0). The iterative algorithm is as follows:n ] T The Newton-Raphson iteration formula is:

[0090] where k is the iteration number, J F is the Jacobian matrix, whose elements

[0091] Step 3042, based on the position offset, the absolute spatial position coordinates of the photoacoustic signal source in space are calculated in combination with the actual position coordinates of the reference detector.

[0092] Further, after obtaining the offset of the photoacoustic signal source position coordinates relative to the reference detector, in combination with the actual position coordinates of the reference detector, the absolute spatial position coordinates of the photoacoustic signal source in space can be calculated. Given that the position coordinates of the reference detector are (x ref ,y ref ,z ref ), the offset is (Δx, Δy, Δz), and the absolute spatial position coordinates (x, y, z) of the photoacoustic signal source are:

[0093] (x = x ref + Δx, y = y ref + Δy, z = z ref + Δz).

[0094] Continuing with the above example, the position coordinates of the reference detector D1 are (10 mm, 0, 0), and the offset (Δx, Δy, Δz) = (-3.8, 4.5, 2.2) is obtained through step 3041. Then the absolute spatial position coordinates of the photoacoustic signal source are: x = 10 + (-3.8) = 6.2 mm; y = 0 + 4.5 = 4.5 mm; z = 0 + 2.2 = 2.2 mm.

[0095] Step 3043, based on the absolute spatial position coordinates of the multiple photoacoustic signal sources, a spatial distribution model of the photoacoustic signal sources in the target biological tissue is constructed, and the spatial characteristics are obtained.

[0096] Further, for multiple photoacoustic signal sources in the target biological tissue, steps 3041 and 3042 are repeated to obtain the absolute spatial position coordinates of each photoacoustic signal source. These coordinates are integrated to construct a spatial distribution model of the photoacoustic signal sources in the target biological tissue. The spatial distribution model can be visually displayed by a three-dimensional coordinate system, in which the coordinates of each photoacoustic signal source correspond to a point in the coordinate system, and the spatial distribution of the photoacoustic signal sources in the biological tissue is presented by the distribution of the points, which is the obtained spatial feature. For example, the points can be plotted using a software tool (such as the three-dimensional plotting function of Matlab) to form a visual spatial distribution model. In an embodiment, the imaging system detects multiple photoacoustic signal sources in the liver tissue of a mouse, for example, a total of 10 photoacoustic signal sources. For each photoacoustic signal source, the absolute spatial position coordinates are calculated according to the method of steps 3041 and 3042. For example, in addition to the calculated coordinate (6.2, 4.5, 2.2) of one photoacoustic signal source described above, another photoacoustic signal source has a calculated coordinate (8.5, 3.1, 5.6). The coordinates of the 10 photoacoustic signal sources are input into the Matlab software, and the scatter3 function is used for plotting to display the distribution of these photoacoustic signal sources in the liver tissue of the mouse in the form of three-dimensional coordinates. From the plotted graph, the spatial position distribution of the photoacoustic signal sources in the liver tissue can be clearly seen, for example, some regions have a higher density of photoacoustic signal sources, which may correspond to blood vessels or lesion regions in the liver, thereby obtaining the spatial feature.

[0097] The embodiment of the present application can accurately construct a spatial distribution model of photoacoustic signal sources in a target biological tissue and obtain accurate spatial features, thereby visually and clearly displaying the spatial positions of the photoacoustic signal sources in the target biological tissue, effectively improving the imaging resolution, and enabling the imaging result to clearly show the fine structures inside the biological tissue, thereby meeting the high-precision imaging requirement.

[0098] In an embodiment, steps 501 to 504 are described as follows:

[0099] Step 501: performing dimension decomposition on the multi-dimensional feature vector according to the association of the multi-dimensional feature vector with different characteristics of the photoacoustic signal to obtain multiple feature sub-vector groups.

[0100] Optionally, the multi-dimensional feature vector contains rich information about the photoacoustic signal, which is related to different characteristics of the photoacoustic signal, such as the frequency, intensity, spatial position, etc. of the signal. According to these characteristics, the multi-dimensional feature vector is dimensionally decomposed and divided into multiple feature sub-vector groups. For example, the multi-dimensional feature vector is For example, the dimensions [f1, f2,..., f kare divided into a group of frequency feature sub-vectors; the dimensions related to spatial position [f j1 j2 jl are divided into a group of spatial feature sub-vectors.

[0101] Step 502, establish the mapping relationship between each feature sub-vector group and the pixel points in the super-resolution image. The mapping relationship of each feature sub-vector group represents the influence range and manner of each feature sub-vector group on the reconstruction of the pixel points.

[0102] Further, for each feature sub-vector group, the mapping relationship between it and the pixel points in the super-resolution image needs to be established. This mapping relationship defines how the feature sub-vector group affects the reconstruction of the pixel points, including the range of influence (for example, whether it has an influence on the local area pixel points or on the pixel points of the entire image) and the manner of influence (for example, changing the value of the pixel points through a certain mathematical transformation). Taking the spatial feature sub-vector group as an example, assuming that the pixel point coordinates of the super-resolution image are (x, y), the mapping relationship can be defined by a function which determines the manner of the influence of the feature sub-vector group on the reconstruction of the pixel points according to the value of the spatial feature sub-vector group and the coordinates of the pixel points. The mapping relationship can be determined by learning a large amount of known photoacoustic signals and corresponding super-resolution image data, for example, using a regression model in machine learning to train the parameters of the mapping function.

[0103] Continuing the above example, for the frequency feature sub-vector group The imaging system establishes the mapping relationship between it and the pixel points in the super-resolution image by training a multi-layer perceptron (MLP) model. The training data set contains a large number of multi-dimensional feature vectors of photoacoustic signals of mouse liver tissues and corresponding known super-resolution images. The frequency feature sub-vector group is taken as the input of the MLP, and the color value (assuming it is a gray value) of each pixel point in the super-resolution image is taken as the output for training. After training, the mapping function For the spatial feature sub-vector group The imaging system adopts a convolutional neural network (CNN) to establish the mapping relationship. After a series of convolution, pooling and other operations on the spatial feature sub-vector group, the mapping function For example, for the pixel point with coordinates (10, 20) in the image, the influence value of the frequency feature sub-vector group on the reconstruction of the pixel point can be calculated by and respectively.

[0104] ​​Step 503, based on the mapping relationship of each feature sub-vector group, the initial predicted pixel value of each pixel point in the super-resolution image is predicted according to the value of the feature sub-vector group.

[0105] Further, the initial predicted pixel value of each pixel point in the super-resolution image is predicted based on the mapping relationship and the value of each feature sub-vector group. For each pixel point (x, y), the corresponding feature sub-vector group is substituted into its mapping function respectively to obtain the influence value of each feature sub-vector group on the pixel point, and the influence values are combined to obtain the initial predicted pixel value P init (x, y). For example, for the frequency feature sub-vector group the spatial feature sub-vector group and the intensity feature sub-vector group The initial predicted pixel value is represented as

[0106] Where C[] is a combination function, the specific form of which is determined by the definition of the mapping relationship, for example, it may be a complex nonlinear transformation to combine the three influence values.

[0107] Continuing the above example, for the pixel point with coordinates (30, 40) in the super-resolution image of the mouse liver tissue, the imaging system substitutes the frequency feature sub-vector group into the mapping function to obtain the frequency influence value v freq ; the spatial feature sub-vector group is substituted into the mapping function to obtain the spatial influence value v space ; the intensity feature sub-vector group is substituted into the mapping function to obtain the intensity influence value v inten . For example, the combination function C[] is a complex nonlinear function, and the initial predicted pixel value P init (30, 40) = C[v freq , v space , v inten ] = 0.6 (assuming the gray value range is 0-1). All pixel points in the super-resolution image are calculated in this way to obtain the initial predicted pixel value distribution of the entire image.

[0108] Step 504, based on the propagation model of the photoacoustic signal in the biological tissue, the initial predicted pixel value of each pixel point is reconstructed into a super-resolution image to obtain the super-resolution photoacoustic image corresponding to the target biological tissue.

[0109] Further, a propagation model of the photoacoustic signal in the biological tissue describes the process from the generation of the photoacoustic signal to the reception of the photoacoustic signal by the detector, including the attenuation, scattering and other characteristics of the signal. The initial predicted pixel value obtained in step 503 is subjected to super-resolution image reconstruction by using the propagation model, to obtain a super-resolution photoacoustic image corresponding to the target biological tissue, as described in steps 5041 to 5044.

[0110] The embodiment of the present application can effectively reconstruct a super-resolution photoacoustic image corresponding to the target biological tissue from the multi-dimensional feature vector, so as to clearly display the fine structure in the target biological tissue, such as the branch of a blood vessel, the boundary of a lesion area, etc., thereby effectively improving the imaging resolution, so that the imaging result can clearly display the fine structure in the biological tissue, and meet the high-precision imaging requirement.

[0111] In an embodiment, steps 5041 to 5044 are described as follows:

[0112] In step 5041, the initial predicted pixel value of each pixel point is subjected to constraint adjustment based on the propagation model of the photoacoustic signal in the biological tissue, to obtain a pixel value after constraint adjustment.

[0113] Optionally, the propagation model of the photoacoustic signal in the biological tissue includes a plurality of physical characteristics, such as the attenuation, scattering and absorption of the signal by the biological tissue, etc. The initial predicted pixel value of each pixel point is subjected to constraint adjustment based on these characteristics. For example, considering the attenuation characteristic of the photoacoustic signal in the propagation process, the pixel point farther away from the photoacoustic signal source should theoretically have a lower signal intensity. Assuming that the function describing the attenuation in the propagation model is A(d)=exp(-0.1d), where d is the distance of the pixel point to the photoacoustic signal source, and the initial predicted pixel value is P init (x,y), the pixel value P constrained (x,y) after constraint adjustment can be calculated by the following formula: P constrained (x,y)=P init (x,y)*A(d(x,y)), where d(x,y) is the distance of the pixel point (x,y) to the photoacoustic signal source, which can be calculated according to the spatial position coordinates of the photoacoustic signal source and the pixel point coordinates by using the Euclidean distance formula.

[0114] In step 5042, the pixel value after constraint adjustment of each pixel point and the pixel value of the neighborhood pixel point of each pixel point are optimized based on the correlation between the feature sub-vectors, to obtain an optimized pixel value.

[0115] Further, there is a certain correlation between different feature sub-vector groups, which reflects the internal relationship between different characteristics of photoacoustic signals. The correlation is used to optimize the adjusted constraint pixel value of each pixel point and the pixel values of its neighborhood pixel points. First, the correlation coefficient matrix R between the feature sub-vector groups is obtained by analyzing a large amount of data, and the element R ij represents the correlation between the i-th and j-th feature sub-vector groups. For pixel point (x, y), its neighborhood pixel point set is denoted as N(x, y). Let the optimization function be O, then the optimized pixel value P optimized (x, y) can be calculated by the following formula:

[0116] P optimized (x, y) = O(P constrained (x, y), {P constrained (m, n)} (m,n)∈N(x,y) , R).

[0117] Where the optimization function O can be a complex function based on deep learning, such as a convolutional neural network (CNN), which takes the adjusted constraint pixel value and its neighborhood pixel value and the correlation coefficient matrix as input, and outputs the optimized pixel value.

[0118] In an embodiment, the imaging system obtains the correlation coefficient matrix R between the feature sub-vector groups by analyzing a large amount of photoacoustic signal data of mouse liver tissue. For pixel point (15, 25), its adjusted constraint pixel value P constrained (15, 25) = 0.6, and the neighborhood pixel point set N(15, 25) includes (14, 24), (14, 25), (14, 26), (15, 24), (15, 26), (16, 24), (16, 25), (16, 26), etc. The adjusted constraint pixel values of these neighborhood pixel points are P constrained (14 / 24) = 0.5, P constrained (14, 25) = 0.55, etc. The imaging system uses a pre-trained convolutional neural network (CNN) as the optimization function O, inputs P constrained (15, 25), the adjusted constraint pixel values of the neighborhood pixel points, and the correlation coefficient matrix R into the CNN, and after a series of convolution, pooling, etc. operations of the network, outputs the optimized pixel value P optimized (15, 25) = 0.65. All pixel points in the image are optimized based on the correlation between the feature sub-vector groups.

[0119] Step 5043, based on the influence of the feature sub-vector groups at different scales on the pixel points, the optimized pixel value of each pixel point is fused and adjusted to obtain the fused pixel value.

[0120] Further, the feature sub-vector groups of different scales have different influences on the pixel points, for example, the feature sub-vector groups of large scales can reflect the overall structure information of the image, and the feature sub-vector groups of small scales can reflect the detail information of the image. Based on the influence difference, the fused pixel value of each pixel point is adjusted. The feature sub-vector groups are divided into different sets according to the scales, and the set of feature sub-vector groups of large scales is denoted as S large , and the set of feature sub-vector groups of small scales is denoted as S small . For each pixel point (x, y), the weight function of the influence of the feature sub-vector groups of large scales on the pixel point is denoted as W large (x, y), the weight function of the influence of the feature sub-vector groups of small scales on the pixel point is denoted as W small (x, y), and W large (x, y) + W small (x, y) = 1. The fused pixel value P fused (x, y) can be calculated by the following formula:

[0121]

[0122] wherein, and are the optimized pixel values obtained based on the feature sub-vector groups of large scales and small scales respectively. The weight functions W large (x, y) and W small (x, y) can be determined by calculating the correlation degree of the feature sub-vector groups of different scales and the pixel point, for example, by using the mutual information method.

[0123] In an embodiment, the imaging system divides the feature sub-vector groups into large-scale and small-scale sets. For the pixel point (20, 30), the correlation degree of the feature sub-vector groups of large scales and the pixel point is high, which is calculated by the mutual information method, and W large (20, 30) = 0.7, and W small (20, 30) = 0.3. The optimized pixel value P obtained based on the feature sub-vector groups of large scales is 0.7, and the optimized pixel value P obtained based on the feature sub-vector groups of small scales is 0.8. According to the fusion formula, the fused pixel value P fused (x, y) = 0.7*0.7 + 0.3*0.8 = 0.73. The fusion adjustment is performed on all the pixel points in the image, and the influences of the feature sub-vector groups of different scales on the pixel points are comprehensively considered.

[0124] In step 5044, the fused pixel values of the entire image are globally corrected based on the spatial distribution model of the photoacoustic signal source, so as to ensure that the relationship between the pixel points in the image matches the spatial distribution of the photoacoustic signal source, and an ultrahigh-resolution photoacoustic image corresponding to the target biological tissue is obtained.

[0125] Further, the spatial distribution model of the photoacoustic signal source describes the distribution of the photoacoustic signal source in the target biological tissue. Based on this model, the global correction is made to the fused pixel value of the entire image, ensuring that the relationship between each pixel point in the image matches the spatial distribution of the photoacoustic signal source. For example, in the area where the photoacoustic signal source is dense, the pixel value in the image should exhibit corresponding characteristics, such as higher contrast or more obvious structural features. Let the global correction function be G, then the corrected pixel value P final (x, y) can be calculated by the following formula:

[0126] P final (x, y) = G(P fused (x, y), spatial distribution model of the photoacoustic signal source).

[0127] Where the global correction function G can be a complex function based on the combination of physical model and deep learning, which adjusts the pixel value according to the spatial position of the photoacoustic signal source and the fused pixel value, so that the reconstructed image is consistent with the spatial distribution of the photoacoustic signal source as a whole.

[0128] In an embodiment, the imaging system constructs a spatial distribution model of the photoacoustic signal source of the mouse liver tissue, and it is known that the photoacoustic signal source is relatively dense in a certain area of the liver. For the pixel point (25, 35) in this area, the fused pixel value P fused (25, 35) = 0.75. The imaging system uses a pre-trained global correction function G based on the combination of physical model and deep learning, inputs P fused (25, 35) and the spatial distribution model of the photoacoustic signal source into the function, and calculates the corrected pixel value P final (25, 35) = 0.8, so that the pixel value in this area can better reflect the characteristics of the dense photoacoustic signal source. All pixel points in the entire super-resolution image are globally corrected in this way, and finally the super-resolution photoacoustic image corresponding to the target biological tissue is obtained.

[0129] The embodiment of the present application can comprehensively optimize the super-resolution image based on the initial predicted pixel value from multiple dimensions, and obtain a high-quality super-resolution photoacoustic image corresponding to the target biological tissue, so that the fine structure and the distribution relationship of the photoacoustic signal source in the target biological tissue can be more clearly and accurately displayed, thereby effectively improving the imaging resolution, so that the imaging result can clearly display the fine structure inside the biological tissue, and meet the high-precision imaging requirement.

[0130] Further, the multi-dimensional feature fusion photoacoustic super-resolution imaging system provided by the present application is described below, and the multi-dimensional feature fusion photoacoustic super-resolution imaging system described below can be correspondingly referred to the multi-dimensional feature fusion photoacoustic super-resolution imaging method described above. Optionally, referring to Figure 2 , Figure 2 is a structural schematic diagram of the multi-dimensional feature fusion photoacoustic super-resolution imaging device provided by the present application, and the multi-dimensional feature fusion photoacoustic super-resolution imaging device comprises:

[0131] The signal acquisition module 210 is configured to irradiate the target biological tissue based on the multi-band laser beam, excite to generate an original photoacoustic signal, and acquire the original photoacoustic signal of different frequency ranges and different propagation directions based on the multi-modal detector array to obtain a multi-dimensional photoacoustic signal.

[0132] The frequency domain analysis module 220 is configured to perform time-frequency analysis on the multi-dimensional photoacoustic signal to convert the time-domain photoacoustic signal to the frequency domain, extract the energy distribution characteristics of the multi-dimensional photoacoustic signal in different frequency bands, and obtain a time-frequency feature map.

[0133] The spatial feature extraction module 230 is configured to calculate the spatial position coordinates of the photoacoustic signal source corresponding to the target biological tissue based on the propagation time difference of the original photoacoustic signal received by the detector and the position distribution of each detector in the detector array, construct a spatial distribution model of the photoacoustic signal source, and obtain a spatial feature.

[0134] The feature fusion module 240 is configured to fuse the time-frequency feature map and the spatial feature to obtain a multi-dimensional feature vector.

[0135] The image reconstruction module 250 is configured to perform super-resolution image reconstruction based on the multi-dimensional feature vector to obtain a super-resolution photoacoustic image corresponding to the target biological tissue.

[0136] The embodiment of the present application acquires the photoacoustic signal from multiple dimensions through the multi-modal detector array, makes up for the deficiencies of a single detector in the frequency response range and the spatial detection angle, and obtains more comprehensive signal data. The time-frequency analysis extracts the energy distribution characteristics of the signal in different frequency bands, which can reflect the dynamic characteristics of different components and structures in the tissue. The spatial feature extraction determines the position information of the signal source. After the fusion of the two, the multi-dimensional feature vector set contains richer and more comprehensive tissue information. Therefore, during image reconstruction, the fusion features can be used to enhance and repair the details in the image more accurately based on the correlation between the features. For example, through the correlation between the spatial position and the frequency characteristics, the position and direction of the fine blood vessels in the tumor tissue can be accurately identified and highlighted, thereby effectively improving the imaging resolution and enabling the imaging result to clearly display the fine structure inside the biological tissue, thereby meeting the high-precision imaging requirement.

[0137] Embodiment of the electronic device provided by the embodiment of the present application is shown. As shown in Figure 3 , Figure 3 Embodiment of the electronic device provided by the embodiment of the present application is shown. As shown in Figure 3 The embodiment of the present application provides an electronic device 300, which comprises a memory 310, a processor 320 and a computer program 311 stored in the memory 310 and capable of running on the processor 320, and the processor 320 realizes the following steps when executing the computer program 311:

[0138] The target biological tissue is irradiated based on a multi-band laser beam, original photoacoustic signals are excited and generated, and the original photoacoustic signals of different frequency ranges and different propagation directions are collected based on a multi-modal detector array to obtain multi-dimensional photoacoustic signals;

[0139] The multi-dimensional photoacoustic signals are subjected to time-frequency analysis to convert the time-domain photoacoustic signals to the frequency domain, the energy distribution characteristics of the multi-dimensional photoacoustic signals in different frequency bands are extracted, and a time-frequency feature map is obtained;

[0140] Based on the propagation time difference of the original photoacoustic signals received by the detector and the position distribution of each detector in the detector array, the spatial position coordinates of the photoacoustic signal source corresponding to the target biological tissue are calculated to construct a spatial distribution model of the photoacoustic signal source, and a spatial feature is obtained;

[0141] The time-frequency feature map and the spatial feature are fused to obtain a multi-dimensional feature vector;

[0142] Based on the multi-dimensional feature vector, a super-resolution image is reconstructed to obtain a super-resolution photoacoustic image corresponding to the target biological tissue.

[0143] Embodiment of the computer readable storage medium provided by the embodiment of the present application is shown. As shown in Figure 4 , Figure 4 Embodiment of the computer readable storage medium provided by the embodiment of the present application is shown. As shown in Figure 4 The embodiment provides a computer readable storage medium 400, which stores a computer program 311, and the computer program 311 is executed by a processor to realize the following steps:

[0144] The target biological tissue is irradiated based on a multi-band laser beam, original photoacoustic signals are excited and generated, and the original photoacoustic signals of different frequency ranges and different propagation directions are collected based on a multi-modal detector array to obtain multi-dimensional photoacoustic signals;

[0145] The multi-dimensional photoacoustic signals are subjected to time-frequency analysis to convert the time-domain photoacoustic signals to the frequency domain, the energy distribution characteristics of the multi-dimensional photoacoustic signals in different frequency bands are extracted, and a time-frequency feature map is obtained;

[0146] ​​​Based on the propagation time difference of the original photoacoustic signal received by the detector and the position distribution of each detector in the detector array, the spatial position coordinates of the photoacoustic signal source corresponding to the target biological tissue are calculated to construct a spatial distribution model of the photoacoustic signal source, so as to obtain the spatial feature;

[0147] Based on the time-frequency feature map and the spatial feature, a multi-dimensional feature vector is obtained.

[0148] Based on the multi-dimensional feature vector, a super-resolution image is reconstructed to obtain a super-resolution photoacoustic image corresponding to the target biological tissue.

[0149] On the other hand, the present application also provides a computer program product, which comprises a computer program, the computer program can be stored on a non-transitory computer readable storage medium, and the computer program can be executed by a processor to enable the computer to perform the multi-dimensional feature fusion photoacoustic super-resolution imaging method provided by the above method, which comprises:

[0150] Based on the multi-waveband laser beam, the target biological tissue is irradiated to generate an original photoacoustic signal, and based on the multi-modal detector array, the original photoacoustic signals of different frequency ranges and different propagation directions are collected to obtain multi-dimensional photoacoustic signals.

[0151] The multi-dimensional photoacoustic signals are subjected to time-frequency analysis to convert the time-domain photoacoustic signals to the frequency domain, and the energy distribution characteristics of the multi-dimensional photoacoustic signals in different frequency bands are extracted to obtain a time-frequency feature map.

[0152] Based on the propagation time difference of the original photoacoustic signal received by the detector and the position distribution of each detector in the detector array, the spatial position coordinates of the photoacoustic signal source corresponding to the target biological tissue are calculated to construct a spatial distribution model of the photoacoustic signal source, so as to obtain the spatial feature;

[0153] Based on the time-frequency feature map and the spatial feature, a multi-dimensional feature vector is obtained.

[0154] Based on the multi-dimensional feature vector, a super-resolution image is reconstructed to obtain a super-resolution photoacoustic image corresponding to the target biological tissue.

[0155] The system embodiments described above are only illustrative, wherein the units illustrated as separate components can or can not be physically separated, and the components illustrated as units can or can not be physical units, i.e., they can be located in one place or distributed on multiple network units. Part or all of the modules can be selected to achieve the purpose of the embodiment scheme according to actual needs. Those skilled in the art can understand and implement without creative labor.

[0156] Those skilled in the art can clearly understand the implementation of the embodiments by means of software and necessary general hardware platforms through the description of the above embodiments, and of course, the embodiments can also be implemented by hardware. Based on such understanding, the above technical solutions can be embodied in the form of a software product, and the computer software product can be stored in a computer readable storage medium, such as a ROM / RAM, a magnetic disk, an optical disk, etc., and includes a plurality of instructions for causing a computer device (which can be a personal computer, a server, or a network device, etc.) to execute the method described in each embodiment or some parts of the embodiment.

[0157] Finally, it should be noted that: the above examples are only used to illustrate the technical solutions of the present application, and not to limit them; although the present application has been described in detail with reference to the foregoing examples, those skilled in the art should understand that: it can still modify the technical solutions recorded in the foregoing examples, or make equivalent replacement for some technical features; and these modifications or replacements do not make the essence of the corresponding technical solutions deviate from the spirit and scope of the technical solutions of the embodiments of the present application.

Claims

1. A photoacoustic super-resolution imaging method with multi-dimensional feature fusion, characterized in that, The method comprises the following steps: Irradiating a target biological tissue based on a multi-wavelength laser beam to excite original photoacoustic signals, and collecting original photoacoustic signals of different frequency ranges and different propagation directions based on a multi-modal detector array to obtain multi-dimensional photoacoustic signals; Performing time-frequency analysis on the multi-dimensional photoacoustic signals to convert the time-domain photoacoustic signals to the frequency domain, extract the energy distribution characteristics of the multi-dimensional photoacoustic signals in different frequency bands, and obtain a time-frequency feature map; Based on the propagation time difference of the original photoacoustic signals received by the detector and the position distribution of each detector in the detector array, the spatial position coordinates of the photoacoustic signal source corresponding to the target biological tissue are calculated to construct a spatial distribution model of the photoacoustic signal source and obtain spatial features; Fusing the time-frequency feature map and the spatial features to obtain a multi-dimensional feature vector; Based on the multi-dimensional feature vector, a super-resolution image is reconstructed to obtain a super-resolution photoacoustic image corresponding to the target biological tissue; Wherein, the super-resolution image reconstruction based on the multi-dimensional feature vector to obtain the super-resolution photoacoustic image corresponding to the target biological tissue comprises: The multi-dimensional feature vector is dimensionally decomposed according to its association with different characteristics of photoacoustic signals to obtain a plurality of feature sub-vector groups; A mapping relationship between each feature sub-vector group and a pixel point in the super-resolution image is established; the mapping relationship of each feature sub-vector group represents the influence range and manner of each feature sub-vector group on the reconstruction of the pixel point; Based on the mapping relationship of each feature sub-vector group, the initial predicted pixel value of each pixel point in the super-resolution image is predicted according to the value of the feature sub-vector group; Based on the propagation model of photoacoustic signals in biological tissues, the initial predicted pixel value of each pixel point is subjected to super-resolution image reconstruction to obtain the super-resolution photoacoustic image corresponding to the target biological tissue; The time-frequency analysis of the multi-dimensional photoacoustic signals to convert the time-domain photoacoustic signals to the frequency domain, extract the energy distribution characteristics of the multi-dimensional photoacoustic signals in different frequency bands, and obtain a time-frequency feature map comprises: The multi-dimensional photoacoustic signals are blocked in time sequence, and each block signal is periodically expanded to obtain an expanded signal; The time-frequency window width is determined based on the local variation characteristics of each expanded signal, and an adaptive time-frequency window is constructed based on the time-frequency window width and the length and instantaneous frequency estimation value at the preset frequency; Each expanded signal is modulated based on the adaptive time-frequency window of each expanded signal to obtain a modulated signal corresponding to each expanded signal; The time-frequency feature map is determined based on the energy distribution characteristics of each modulated signal in different frequency bands.

2. The multi-dimensional feature-fusion photoacoustic super-resolution imaging method of claim 1, wherein, The super-resolution image reconstruction based on the propagation model of photoacoustic signals in biological tissues to each pixel point initial predicted pixel value, to obtain the super-resolution photoacoustic image corresponding to the target biological tissue, comprises: Based on the propagation model of photoacoustic signals in biological tissues, the initial predicted pixel value of each pixel point is subjected to constraint adjustment to obtain a constraint adjusted pixel value; Optimize the constraint adjusted pixel value of each pixel point and the pixel value of the neighborhood pixel point of each pixel point based on the correlation between the feature sub-vector groups, to obtain an optimized pixel value; Fuse and adjust the optimized pixel value of each pixel point based on the influence of the feature sub-vector groups under different scales, to obtain a fused pixel value; Based on the spatial distribution model of the photoacoustic signal source, globally correct the fused pixel value of the entire image, ensure that the relationship between the pixel points in the image matches the spatial distribution of the photoacoustic signal source, and obtain the super-resolution photoacoustic image corresponding to the target biological tissue.

3. The multi-dimensional feature-fusion photoacoustic super-resolution imaging method of claim 1, wherein, Based on the propagation time difference of the original photoacoustic signal received by the detector and the position distribution of each detector in the detector array, the spatial position coordinates of the photoacoustic signal source corresponding to the target biological tissue are calculated to construct the spatial distribution model of the photoacoustic signal source, and the spatial features are obtained, including: Based on the propagation time difference of the original photoacoustic signal received by the detector and the propagation speed of the photoacoustic signal in the biological tissue, the corresponding propagation distance difference when the signal is received by different detectors is determined; Based on the position distribution of each detector in the detector array, the detector position vector of each detector in space is determined, and based on the propagation distance difference and the detector position vector, an initial distance difference equation about the position coordinates of the photoacoustic signal source is constructed; Based on the triangular geometric relationship between the detector and the photoacoustic signal source, the initial distance difference equation is optimized to obtain a target distance difference equation; Based on the target distance difference equation, the spatial position coordinates of the photoacoustic signal source corresponding to the target biological tissue are calculated to construct the spatial distribution model of the photoacoustic signal source, and the spatial features are obtained.

4. The multi-dimensional feature-fusion photoacoustic super-resolution imaging method of claim 3, wherein, Based on the target distance difference equation, the spatial position coordinates of the photoacoustic signal source corresponding to the target biological tissue are calculated to construct the spatial distribution model of the photoacoustic signal source, and the spatial features are obtained, including: The equation group is obtained by combining a plurality of target distance difference equations, and the offset of the photoacoustic signal source position coordinates relative to the reference detector is obtained by solving the equation group; Based on the position offset and the actual position coordinates of the reference detector, the absolute spatial position coordinates of the photoacoustic signal source in space are calculated; Based on the absolute spatial position coordinates of a plurality of photoacoustic signal sources, the spatial distribution model of the photoacoustic signal source in the target biological tissue is constructed, and the spatial features are obtained.

5. The multi-dimensional feature-fusion photoacoustic super-resolution imaging method of claim 1, wherein, Based on the energy distribution characteristics of each modulated signal in different frequency bands, the time-frequency feature map is determined, including: Perform generalized time-frequency transformation on each modulated signal to obtain a generalized time-frequency transformation result, and determine the energy density of different frequency bands on the time-frequency plane based on each generalized time-frequency transformation result; Enhance the energy density of each frequency band based on a nonlinear logarithm to obtain an enhanced energy density; Take time as the horizontal coordinate and frequency band index as the vertical coordinate, take the enhanced energy density of each frequency band as the gray value or color value of the map, and generate an initial feature map; For each feature point in the initial feature map, the gray value or color value of each feature point is updated based on the average gray value or average color value in the local neighborhood centered on each feature point, to generate the time-frequency feature map.

6. A multi-dimensional feature fusion photoacoustic super-resolution imaging device, characterized in that, The application is applied to the multi-dimensional feature fusion photoacoustic super-resolution imaging method of any one of claims 1 to 5; the multi-dimensional feature fusion photoacoustic super-resolution imaging device comprises: A signal acquisition module is configured to irradiate a target biological tissue based on a multi-waveband laser beam, excite original photoacoustic signals, and acquire original photoacoustic signals of different frequency ranges and different propagation directions based on a multi-modal detector array to obtain multi-dimensional photoacoustic signals. A frequency domain analysis module is configured to perform time-frequency analysis on the multi-dimensional photoacoustic signals to convert time-domain photoacoustic signals to a frequency domain, extract energy distribution characteristics of the multi-dimensional photoacoustic signals in different frequency bands, and obtain a time-frequency feature map. A spatial feature extraction module is configured to calculate spatial position coordinates of a photoacoustic signal source corresponding to the target biological tissue based on a propagation time difference of the original photoacoustic signals received by the detectors and a position distribution of the detectors in the detector array, construct a spatial distribution model of the photoacoustic signal source, and obtain spatial features. A feature fusion module is configured to fuse the time-frequency feature map and the spatial features to obtain a multi-dimensional feature vector. An image reconstruction module is configured to perform super-resolution image reconstruction based on the multi-dimensional feature vector to obtain a super-resolution photoacoustic image corresponding to the target biological tissue. The super-resolution image reconstruction based on the multi-dimensional feature vector to obtain the super-resolution photoacoustic image corresponding to the target biological tissue comprises: dimensionally decomposing the multi-dimensional feature vector according to its association with different characteristics of photoacoustic signals to obtain a plurality of feature sub-vector groups; establishing a mapping relationship between each feature sub-vector group and a pixel point in the super-resolution image, wherein the mapping relationship of each feature sub-vector group represents an influence range and a manner of each feature sub-vector group on the reconstruction of the pixel point; predicting an initial predicted pixel value of each pixel point in the super-resolution image based on the mapping relationship of each feature sub-vector group and the value of the feature sub-vector group; performing super-resolution image reconstruction on the initial predicted pixel value of each pixel point based on a propagation model of photoacoustic signals in biological tissues to obtain the super-resolution photoacoustic image corresponding to the target biological tissue. The time-frequency analysis on the multi-dimensional photoacoustic signals to convert time-domain photoacoustic signals to a frequency domain, extract energy distribution characteristics of the multi-dimensional photoacoustic signals in different frequency bands, and obtain a time-frequency feature map comprises: blocking the multi-dimensional photoacoustic signals in a time sequence and periodically expanding each block of signals to obtain expanded signals; determining a time-frequency window width based on the local variation characteristics of each expanded signal, and constructing an adaptive time-frequency window based on the time-frequency window width and the length and instantaneous frequency estimation value at the preset frequency; modulating each expanded signal based on the adaptive time-frequency window of each expanded signal to obtain a modulated signal corresponding to each expanded signal. The time-frequency feature atlas is determined based on energy distribution characteristics of each modulated signal in different frequency bands.

7. An electronic device, comprising: A memory for storing a computer software program; A processor for reading and executing the computer software program, wherein the processor, when executing the computer software program, implements the multi-dimensional feature fusion photoacoustic super-resolution imaging method according to any one of claims 1 to 5.

8. A non-transitory computer-readable storage medium having stored therein a computer software program, characterized in that, The computer software program, when executed by the processor, implements the multi-dimensional feature fusion photoacoustic super-resolution imaging method according to any one of claims 1 to 5.

Images