CN113393472A - Zoning composite phase unwrapping method based on digital holographic microscopic imaging - Google Patents
Zoning composite phase unwrapping method based on digital holographic microscopic imaging Download PDFInfo
- Publication number
- CN113393472A CN113393472A CN202110595634.6A CN202110595634A CN113393472A CN 113393472 A CN113393472 A CN 113393472A CN 202110595634 A CN202110595634 A CN 202110595634A CN 113393472 A CN113393472 A CN 113393472A
- Authority
- CN
- China
- Prior art keywords
- phase
- quality
- holographic
- image
- map
- Prior art date
- Legal status (The legal status is an assumption and is not a legal conclusion. Google has not performed a legal analysis and makes no representation as to the accuracy of the status listed.)
- Granted
Links
- 238000000034 method Methods 0.000 title claims abstract description 73
- 238000003384 imaging method Methods 0.000 title claims abstract description 31
- 239000002131 composite material Substances 0.000 title claims abstract description 24
- 238000013316 zoning Methods 0.000 title claims abstract description 16
- 238000010587 phase diagram Methods 0.000 claims abstract description 51
- 238000004422 calculation algorithm Methods 0.000 claims abstract description 9
- 230000005540 biological transmission Effects 0.000 claims abstract description 8
- 238000012847 principal component analysis method Methods 0.000 claims abstract description 5
- 239000011159 matrix material Substances 0.000 claims description 54
- 238000010586 diagram Methods 0.000 claims description 26
- 238000001228 spectrum Methods 0.000 claims description 24
- 238000009647 digital holographic microscopy Methods 0.000 claims description 13
- 238000000513 principal component analysis Methods 0.000 claims description 10
- 230000009466 transformation Effects 0.000 claims description 8
- 230000003287 optical effect Effects 0.000 claims description 7
- 230000001131 transforming effect Effects 0.000 claims description 3
- 239000013598 vector Substances 0.000 claims description 3
- 230000008569 process Effects 0.000 description 7
- 230000001427 coherent effect Effects 0.000 description 5
- 230000000694 effects Effects 0.000 description 5
- 230000011218 segmentation Effects 0.000 description 4
- 238000000638 solvent extraction Methods 0.000 description 3
- 230000009286 beneficial effect Effects 0.000 description 2
- 238000004364 calculation method Methods 0.000 description 2
- 238000005259 measurement Methods 0.000 description 2
- 238000012545 processing Methods 0.000 description 2
- 238000011084 recovery Methods 0.000 description 2
- 241000381592 Senegalia polyacantha Species 0.000 description 1
- 230000004075 alteration Effects 0.000 description 1
- 238000004458 analytical method Methods 0.000 description 1
- 238000003491 array Methods 0.000 description 1
- 238000006243 chemical reaction Methods 0.000 description 1
- 150000001875 compounds Chemical class 0.000 description 1
- 230000007547 defect Effects 0.000 description 1
- 238000005516 engineering process Methods 0.000 description 1
- 230000002349 favourable effect Effects 0.000 description 1
- 230000004927 fusion Effects 0.000 description 1
- 230000010354 integration Effects 0.000 description 1
- 238000005192 partition Methods 0.000 description 1
- 230000010287 polarization Effects 0.000 description 1
- 230000009467 reduction Effects 0.000 description 1
- 238000011160 research Methods 0.000 description 1
- 230000000717 retained effect Effects 0.000 description 1
- 238000004088 simulation Methods 0.000 description 1
- 238000013519 translation Methods 0.000 description 1
Images
Classifications
-
- G—PHYSICS
- G06—COMPUTING; CALCULATING OR COUNTING
- G06T—IMAGE DATA PROCESSING OR GENERATION, IN GENERAL
- G06T7/00—Image analysis
- G06T7/10—Segmentation; Edge detection
- G06T7/11—Region-based segmentation
-
- G—PHYSICS
- G06—COMPUTING; CALCULATING OR COUNTING
- G06F—ELECTRIC DIGITAL DATA PROCESSING
- G06F17/00—Digital computing or data processing equipment or methods, specially adapted for specific functions
- G06F17/10—Complex mathematical operations
- G06F17/14—Fourier, Walsh or analogous domain transformations, e.g. Laplace, Hilbert, Karhunen-Loeve, transforms
- G06F17/141—Discrete Fourier transforms
-
- G—PHYSICS
- G06—COMPUTING; CALCULATING OR COUNTING
- G06F—ELECTRIC DIGITAL DATA PROCESSING
- G06F17/00—Digital computing or data processing equipment or methods, specially adapted for specific functions
- G06F17/10—Complex mathematical operations
- G06F17/16—Matrix or vector computation, e.g. matrix-matrix or matrix-vector multiplication, matrix factorization
-
- G—PHYSICS
- G06—COMPUTING; CALCULATING OR COUNTING
- G06T—IMAGE DATA PROCESSING OR GENERATION, IN GENERAL
- G06T5/00—Image enhancement or restoration
- G06T5/50—Image enhancement or restoration using two or more images, e.g. averaging or subtraction
-
- G—PHYSICS
- G06—COMPUTING; CALCULATING OR COUNTING
- G06T—IMAGE DATA PROCESSING OR GENERATION, IN GENERAL
- G06T7/00—Image analysis
- G06T7/10—Segmentation; Edge detection
- G06T7/136—Segmentation; Edge detection involving thresholding
-
- G—PHYSICS
- G06—COMPUTING; CALCULATING OR COUNTING
- G06T—IMAGE DATA PROCESSING OR GENERATION, IN GENERAL
- G06T2207/00—Indexing scheme for image analysis or image enhancement
- G06T2207/20—Special algorithmic details
- G06T2207/20048—Transform domain processing
- G06T2207/20056—Discrete and fast Fourier transform, [DFT, FFT]
-
- G—PHYSICS
- G06—COMPUTING; CALCULATING OR COUNTING
- G06T—IMAGE DATA PROCESSING OR GENERATION, IN GENERAL
- G06T2207/00—Indexing scheme for image analysis or image enhancement
- G06T2207/20—Special algorithmic details
- G06T2207/20212—Image combination
- G06T2207/20221—Image fusion; Image merging
-
- G—PHYSICS
- G06—COMPUTING; CALCULATING OR COUNTING
- G06T—IMAGE DATA PROCESSING OR GENERATION, IN GENERAL
- G06T2207/00—Indexing scheme for image analysis or image enhancement
- G06T2207/30—Subject of image; Context of image processing
- G06T2207/30168—Image quality inspection
Landscapes
- Engineering & Computer Science (AREA)
- Physics & Mathematics (AREA)
- General Physics & Mathematics (AREA)
- Mathematical Physics (AREA)
- Theoretical Computer Science (AREA)
- Pure & Applied Mathematics (AREA)
- Mathematical Optimization (AREA)
- Mathematical Analysis (AREA)
- Data Mining & Analysis (AREA)
- Computational Mathematics (AREA)
- Algebra (AREA)
- Databases & Information Systems (AREA)
- Software Systems (AREA)
- General Engineering & Computer Science (AREA)
- Computer Vision & Pattern Recognition (AREA)
- Computing Systems (AREA)
- Discrete Mathematics (AREA)
- Holo Graphy (AREA)
Abstract
The invention provides a zoning composite phase unwrapping method based on digital holographic microscopic imaging, which comprises the following steps of: collecting an original holographic image of a tested sample; extracting a holographic phase diagram; carrying out phase dephasing on the holographic phase diagram by adopting a PCA principal component analysis method to obtain the holographic phase diagram of the sample after dephasing; calculating phase derivative variance based on path transmission unwrapping to obtain a quality map; dividing the distribution of the quality map according to the quality factor of the quality map, unwrapping different quality areas by adopting different phase unwrapping algorithms, and fusing and splicing the quality areas to obtain unwrapped phases. The invention maintains the accuracy of the unwrapped phase obtained by resolving the high-quality area and strengthens the accuracy of the phase obtained by resolving the low-quality area.
Description
Technical Field
The invention belongs to the optical microscopic imaging technology, and particularly relates to a zoning composite phase unwrapping method based on digital holographic microscopic imaging.
Background
The traditional phase unwrapping method based on digital holographic microscopic imaging mainly comprises two types of methods based on path tracking and minimum norm. Commonly used methods based on path tracking are the "cutch" method and the quality-oriented graph method proposed by Goldstein et al. The branch cutting method is mainly characterized in that a branch cutting line is established at a pixel point corresponding to two error phase values of the phase diagram to distinguish the phase diagram in an error area, so that the end points of the branch cutting line can be bypassed to perform unwrapping, and the accuracy is improved. Although the 'branch-and-cut' method can effectively avoid a 'pole' with a large error, a new 'pole' must be judged and found out first so as to establish an error segmentation line, so that the process is complex and time-consuming, and is not beneficial to unpacking in an image with overlarge pixels; the Quality-oriented unpacking method mainly utilizes a similar Quality Map (Quality Map) to accurately judge the Quality and the physical complexity of each point, starts from a high-Quality point and finishes at a low-Quality starting point, and unpacks the parcel based on path integration. The minimum unwrapping method based on norm is the least square unwrapping method which is a global unwrapping mode, error points cannot be eliminated in the unwrapping process, but the 'dead points' are bypassed, so that the dead points still exist in the unwrapped phase diagram, but tend to be smooth in numerical calculation, and the method has the main defect that a high-precision phase diagram cannot be unwrapped in a low-quality phase diagram.
Disclosure of Invention
The invention aims to provide a partitioned composite phase unwrapping method based on digital holographic microscopy imaging, which aims to solve the influence of 'dead spots' in a digital holographic unwrapping algorithm on unwrapping, maintain high-quality areas to obtain wrapping accuracy, simultaneously strengthen the wrapping accuracy in low-quality areas and effectively inhibit the propagation of errors.
The technical scheme for realizing the purpose of the invention is as follows: a zoning composite phase unwrapping method based on digital holographic microscopic imaging comprises the following specific steps:
and 5, dividing the distribution of the quality map according to the quality factor of the quality map, unwrapping different quality areas by adopting different phase unwrapping algorithms, and fusing and splicing the quality areas to obtain unwrapped phases.
Preferably, the original holographic image is fourier-transformed into a spectrogram, and the specific method for extracting the +1 level spectrum is as follows:
the light intensity of the original holographic image is:
in the formula IH(x, y) is the intensity information of the original image, IO(x, y) is the object light intensity information, IR(x, y) is the reference light intensity information, phiO(x, y) is object light phase information, phiRAnd (x, y) is reference light phase information.
And carrying out fast Fourier transform on the light intensity of the holographic image to obtain:
in the formula F [ ■]Denotes a fast Fourier transform, α (x, y) ═ IO(x,y)+IR(x, y) represents object light IO(x, y) and reference light IR(x, y) a direct current component after the correlation, a ═ F [ α (x, y)]The direct current component of the two phases of dry light is subjected to fast Fourier transform to obtain a result;
in order to be a real image,indicating the intensity of the real image, exp [ i φ ]O(x,y)-iφR(x,y)]Is the argument distribution of the real image, i denotes the imaginary unit, F+1Is the +1 order spectrum corresponding to the real image;
the +1 order spectrum is extracted separately from the fourier spectrogram directly by clipping the image.
Preferably, the extracted holographic phase map is:
wherein Im (■) represents f+1The imaginary part of Re (■) represents f+1Real part of (f)+1Is a spatial light intensity distribution.
Preferably, the specific method for performing phase dephasing on the holographic phase diagram by using the PCA principal component analysis method comprises the following steps:
(1) holographic phase pattern phiH(X, y) is denoted as matrix X;
(2) carrying out zero average value transformation on each row of the matrix X again;
(4) computing eigenvectors e of the covariance matrix C1,e2……etAnd a characteristic value λ, t ═ 1,2 … … N × M;
(5) sorting the eigenvalues from big to small, and selecting the largest k of the eigenvalues; respectively taking the corresponding k eigenvectors as row vectors to form an eigenvector matrix P; y is PX which is the data matrix after the dimension of the matrix X is reduced to k dimension;
(6) and performing least square fitting on the obtained reduced-dimension data matrix Y to remove the phase difference.
Preferably, the phase derivative variance is calculated as:
wherein (m, n) is a pixel gradient range of the subscript value (i, j) of each pixel point in a square having a size represented by k × k pixels in the central range,andis the gradient of the phase difference between the pixel wrap points of the phase along two different directions x and y respectively,andrespectively representing two averages of two phase-deflector function within the window. According to the formula, a quality map Q (M, N) corresponding to the holographic phase map is obtained, wherein the quality map Q (M, N) is a two-dimensional data matrix map with M rows and N columns, each element of the matrix is a quality factor, and the quality factor is a value of [0, 1]]The numerical value of (c).
Preferably, the specific method for partitioning the distribution of the quality map according to the quality factor of the quality map is as follows:
setting a quality threshold value, and determining the quality ratio of the quality factor in the quality map which is greater than or equal to the set threshold value;
in the quality map Q (M, N) of M multiplied by N, small rectangles of a multiplied by b are searched from four top corners of the quality map to the center, and a is more than or equal to 1 and less than or equal to M; b is more than or equal to 1 and less than or equal to N, when the mass ratio r 'of the small rectangle is more than or equal to x r, x is less than or equal to a random index, the searching is stopped, four small-mass diagram position coordinates determined after four vertexes are respectively searched are obtained, the areas of the four small-mass diagrams are compared, and one block with the smallest area is selected as a final zoning position to obtain a zoning mass diagram Q' (m, N); dividing the quality map Q (m, n) into four rectangular areas according to the position coordinates of the partitioned quality map Q' (m, n), reserving the overlapped part of L rows and L columns of each area, and determining the coordinate positions of the four areas;
respectively calculating the mass ratio of the four regions, and judging the region as a low-mass region when the mass ratio of the regions is more than or equal to the mass ratio x r of the mass diagram; and when the area quality ratio is smaller than the quality area quality ratio, judging the area as a high quality area.
Preferably, the different quality regions are unwrapped by using a different phase unwrapping algorithm, and the specific method comprises the following steps: the quality diagram Q (m, n) is divided into position coordinates of four areas and is brought into the phase diagram, the phase diagram is divided into four small phase diagrams with different qualities according to the dividing coordinates, the small phase diagrams are unwrapped in the low-quality areas by adopting a path tracking method, and the unwrapped phase diagrams in the high-quality areas are unwrapped by adopting a minimum norm method to obtain unwrapped phase diagrams of four different quality areas.
Preferably, the quality threshold is set to 0.5.
Preferably, the mass ratio is calculated by the formula:
compared with the prior art, the invention has the following remarkable advantages: (1) the invention maintains the accuracy of the unwrapped phase obtained by resolving the high-quality area and strengthens the accuracy of the phase obtained by resolving the low-quality area; (2) the invention effectively inhibits the propagation of errors and greatly optimizes the effect of the finally obtained unwrapped phase diagram.
The present invention is described in further detail below with reference to the attached drawing figures.
Drawings
FIG. 1 is a flow chart of a partitioned composite phase unwrapping method based on digital holographic microscopy imaging.
FIG. 2 is a system diagram of a partitioned composite phase unwrapping method based on digital holographic microscopy imaging.
Fig. 3 is a process diagram of phase information extracted from the digital hologram of HELA cell according to the present invention.
FIG. 4 is a comparison of HELA cell digital holograms before and after phase difference analysis by PCA principal component analysis.
FIG. 5 is a schematic sectional view of the HELA cell digital hologram sectional composite unwrapping method of the present invention.
FIG. 6 is a comparison of phase information obtained by different phase unwrapping methods.
FIG. 7 is a numerical processing diagram of cell segmentation using the phase map obtained by the present invention.
Fig. 8 is a diagram of final partitioning of the quality map.
Detailed Description
The invention mainly focuses on the research of the unwrapping method, provides zoning composite unwrapping and can obtain more accurate complete phase information.
As shown in fig. 1, a zonal compound phase unwrapping method based on digital holographic microscopy imaging is based on hologram acquisition by a microscope system built by a transmission type holographic optical path. Taking the light path based on transmission type holographic imaging as an example, the holographic microscope imaging system comprises a laser light source, a beam splitter, a reflector, a phase objective, a sample and a CCD camera. The sample selects living cells (HELA cells), the hologram of the living cells is collected, the phase information of the hologram is extracted according to the digital holographic imaging principle, the phase distortion is further corrected, and the hologram phase image after phase difference removal is used as the information of phase unwrapping for subsequent work.
The method comprises the following specific steps:
step 1: collecting a digital holographic microscopic image: a microscopic imaging system built by a transmission type holographic optical path is adopted, and a CCD camera is used for collecting a holographic image of a tested sample under a microscope. Two beams of interference light are introduced in the imaging process and are respectively marked as object light O and reference light R, and an image under the microscope lens is recorded by a CCD camera. As shown in fig. 2(a), after the microscopic imaging system is constructed, the following adjustment operations are performed on the optical path:
1. the CCD camera, the beam splitter prism and a 20cm collimating lens positioned in front of the CCD camera are connected with each other, and the camera is subjected to focusing adjustment.
2. Placing the sample at the position of the focal plane of the microscope objective;
3. adjusting a light path, wherein a laser emits a light beam, the light beam is divided into object light and reference light by a beam splitter, and the object light emits parallel light from a beam splitter prism through a reflector, a microscope objective and a sample;
4. the reference light emits parallel light from the beam splitter prism through the reflector;
5. after the assembly, a live cell digital holographic microscope is formed, as shown in fig. 2(b), and a digital holographic image is collected on the platform and marked as IH(x, y) as shown in FIG. 2 (c).
Step 2: extracting the wrapping phase of the holographic image: fourier transforming the digital holographic image into a spectrogram, extracting +1 level spectrum, moving to the center, and performing inverse Fourier transformation to obtain a holographic image, which is marked as IH"(x, y). And the holographic image is focused by using the angular spectrum propagation principle, and a holographic phase image is extracted and recorded as phiH。
wherein A isO(x,y)、AR(x, y) are the light intensity information of the object light and the reference light respectively,According to the principle of light wave superposition, the complex amplitude distribution of the phase information of the object light and the reference light respectively is the interference superposition of two coherent lights:
H(x,y)=UO(x,y)+UR(x,y)
the intensity of the two coherent light beams can be expressed by the following formula:
the first two terms in the formula represent the phase intensity of two coherent lights, and the third term represents a hidden interference information element, which includes phase intensity information and amplitude distribution information of the optical wavefront with the same frequency.
Considering that the recording process of the digital hologram is a process of interference fringe digitization, the intensity distribution of the digital hologram is IH(x, y) which converts this two-dimensional continuous light intensity distribution into an M x N matrix I by means of a digital imaging deviceH' (x ', y ') (x ' ═ 0, 1, … …, M-1; y ' ═ 0, 1, … …, N-1), where the pixel size of the digital detector can be expressed as Δ x Δ y. Therefore, the recording of the digital hologram has a spatial discretization process. To further analyze the individual information recorded by the digital hologram, its light intensity is marked by the following formula:
IH(x,y)=|O|2+|R|2+RO*+R*O
the first two terms of the above equation are the direct addition of the intensity of the object light and the reference light, referred to as the zero order image. Third part RO*Is a real image, or +1 level image, which represents information of the object light conjugate image. Fourth term R*O is a virtual image, or-1 order image, representing the information of the object light itself, and is most conveniently viewed directly. The twin image is composed of a real image and a virtual image, both of which can resolve information of object light, but the combination of the two can easily cause information loss, especially phase information. Therefore, in order to extract the phase information better, an off-axis light path recording method is adopted, and the method is common and can separate the required information well.
The light intensity distribution of the digital hologram is further expressed as:
IH(x,y)=|O|2+|R|2+|R|Oexp(-jksinθx)+|R|O*exp(jksinθx)
wherein | O |2Representing the luminous intensity information of the reference object, | R |2Representing reference light intensity information, | R | Oexp (-jksin θ x) is a real image RO*Expanded representation, | R | O*exp (jksin θ x) is the virtual image R*The expansion of O indicates that exp (-jksin theta x) indicates the phase factor and theta is the angle of the reference light relative to the object light.
It is easy to see that the effective information to be acquired is hidden in the real image (+1 level image) or the virtual image (-1 level image), and the phase factor is just to perform frequency domain mirror image translation on the fourier spectrum of the real image or the virtual image with respect to the zero-level image, when the θ angle is large enough, the zero-level image is located at the origin of the frequency domain, and the real image and the virtual image (also called as two interference images) present mirror image distribution with respect to the origin, so that one of the real image or the virtual image is screened out and then reconstructed from the frequency domain back to the spatial domain, and then the complete spatial light intensity distribution I capable of extracting the amplitude information and the phase information can be obtainedH”(x,y)。
The light intensity recorded by the digital hologram is expressed as:
in the formula IH(x, y) is the intensity information of the original image, IO(x, y) is the object light intensity information, IR(x, y) is the reference light intensity information, phiO(x, y) is object light phase information, phiRAnd (x, y) is reference light phase information.
The fast Fourier transform is carried out on the light intensity, so that the following can be obtained:
in the formula F [ ■]Denotes a fast Fourier transform, α (x, y) ═ IO(x,y)+IR(x, y) represents object light IO(x, y) and reference light IR(x, y) a direct current component after the correlation, a ═ F [ α (x, y)]The direct current component of the two phases of dry light is subjected to fast Fourier transform to obtain a result;
in order to be a real image,indicating the intensity of the real image, exp [ i φ ]O(x,y)-iφR(x,y)]Is the argument distribution of the real image, i denotes the imaginary unit, F+1Is the +1 order spectrum corresponding to the real image;
is a virtual image, F-1The result is obtained by carrying out fast Fourier transform on the virtual image, and is a-1 level spectrum corresponding to the virtual image;
the +1 level spectrum (or-1 level spectrum) can be extracted from the Fourier spectrogram directly by intercepting the image, and the real image or the virtual image of the digital hologram is separated. In the simulation of MATLAB, in order to separate a real image and a virtual image, a certain piece of spectrum information can be captured in a screenshot mode, but the information obtained in screenshot is a user-defined pixel size, is different from the original image and is not beneficial to the next calculation, so that the image obtained in screenshot needs to be subjected to pixel recovery through numerical value conversion. The most intense direct current component of a spectrogram is usually located at a dot (center) of an image, so that the + 1-level spectrogram after pixel recovery is moved to the dot, and a pixel point corresponding to the maximum value of the + 1-level spectrogram needs to be obtained.
Then, the spatial light intensity distribution I of the holographic image for extracting the amplitude information and the phase information can be obtained by carrying out inverse Fourier transform on the spatial light intensity distribution IH″(x,y):
IH″(x,y)=Φ-1[F+1]=Φ-1[F(f+1(x,y))]=f+1(x,y)
Wherein phi-1[■]Representing the inverse Fourier transform, and finally f+1And (x, y) performing positive and negative tangent function transformation to extract phase information of the digital hologram:
wherein Im (■) represents f+1The imaginary part of Re (■) represents f+1The real part of (a); the spectrogram obtained by fourier-transforming the phase information extracted from the HELA digital hologram of fig. 2(c) is shown in fig. 3(a), the spectrogram obtained by extracting the +1 level spectrum from the spectrogram at the center position is shown in fig. 3(b), and the spectrogram obtained by inverse fourier-transforming the +1 level spectrum to extract the phase information is shown in fig. 3 (c).
And step 3: generating a wrapped phase map: performing phase dephasing on the holographic phase diagram by adopting a PCA (principal component analysis) method to obtain a complete wrapped phase diagram, eliminating the influence of phase distortion, and obtaining the holographic phase diagram recorded as phi 'of the sample after the phase distortion is eliminated'H。
The dimensionality reduction operation is carried out by a Principal Component Analysis (PCA) method, and mainly the covariance matrix form is symmetrically triangulated. The assumed original data matrix is X, the corresponding covariance data matrix is C, the data matrix P is defined as a matrix formed by combining a group of original bases according to different step sequences, and the assumed matrix Y is the original matrix data after the matrix X is subjected to basis transformation. And the covariance matrix relationship of Y is represented by matrix C and matrix D, respectively, then this covariance matrix relationship can be represented as follows:
the PCA principal component analysis method is used for carrying out specific work of phase aberration removal:
(1) holographic phase pattern phiH(X, y) is denoted as matrix X, which is the phase diagram φHThe two-dimensional arrays of (X, y) are arranged into a matrix X with m rows and n columns;
(2) carrying out zero average value transformation on each row of the matrix X again;
(4) finding eigenvectors e of the covariance matrix C1,e2……etAnd a characteristic value λ, t ═ 1,2 … … N × M;
(5) sorting the eigenvalues from large to small, and selecting the largest k of the eigenvalues. Then, respectively taking the k eigenvectors corresponding to the k eigenvectors as row vectors to form an eigenvector matrix P; y is PX which is the data matrix after the dimension of the matrix X is reduced to k dimension;
(6) and performing least square fitting on the obtained new dimension-reduced data matrix Y to remove the phase difference.
The phase difference results obtained are shown in fig. 4(a), and the HELA cell phase compensation results obtained by correcting the phase difference are shown in fig. 4 (b).
And 4, step 4: calculating wrapped phase quality map phase: unwrapping based on path propagation, calculating phase derivative variance Z according to pixel size of wrapped phase hologrammnAnd obtaining a quality diagram.
The phase derivative variance may be defined by:
this expression indicates that (m, n) is the pixel gradient range of the subscript value (i, j) of each pixel point in a square of the size represented by k × k pixels in the central range,andis the gradient of the phase difference between the pixel wrap points of the phase along two different directions x and y respectively,andwhich are used to represent two averages of two phase-deflector-rate functions within the window, respectively. According to the formula, a mass diagram Q (M, N) corresponding to the holographic phase diagram can be obtained, wherein the mass diagram Q (M, N) is a two-dimensional data matrix with M rows and N columnsIn the figure, each element of the matrix is a quality factor, and the quality factor is a value of [0, 1]The numerical value of (c) is shown in FIG. 5 (b).
And 5: partitioning and compositely unwrapping according to the quality diagram: and dividing the quality map into a low-quality area and a high-quality area according to the quality factor of the quality map, wherein different phase unwrapping algorithms are adopted for different quality areas. Obtaining an unwrapped phase diagram in a low-quality region by adopting an unwrapping method based on path tracking; and obtaining unwrapped phase diagrams of different quality areas by adopting a phase unwrapping method based on the minimum norm in the high quality area. And fusing and splicing the quality areas to obtain the unwrapped phase.
And judging whether the quality is good or bad according to the quality factors in the obtained quality map Q (m, n). A quality map Q (M, N) solved by the phase derivative variance is a matrix map with a quality factor value of M × N between [0, 1], when the quality factor of the quality map Q (M, N) is greater than 0.5, it indicates that the point tends to have poor quality, and when a region in a quality map has a point with a dense Q (M, N) element greater than 0.5, the region is a low-quality region, which is not favorable for the traditional phase unwrapping method based on the minimum norm, the traditional phase unwrapping method based on the minimum norm will polish a local region to be the same as the whole region during unwrapping, actively drive the unwrapped phase to be "smooth" in the region, and greatly lose phase accuracy; therefore, in low quality areas, it is more recommended to adopt a unwrapping method based on path tracking, while in high quality areas, a phase unwrapping method based on minimum norm can be considered.
The zoning composite unwrapping principle based on the quality map is as follows:
step 5.1: setting the quality threshold value to be 0.5, and determining the quality proportion of the quality factor in the quality map to be greater than or equal to 0.5 as follows:
step 5.2: determining a high-quality region and a low-quality region in a quality map by taking the mass ratio r as a basic reference, wherein the specific method comprises the following steps: firstly, in an M multiplied by N quality map Q (M, N), small rectangles of a multiplied by b are searched from four top corners of the quality map to the center, and a is more than or equal to 1 and less than or equal to M; b is more than or equal to 1 and less than or equal to N, when the mass ratio r 'of the small rectangle is more than or equal to x r (x is more than or equal to 0.5 and less than or equal to 1 and is a random index), the searching is stopped, the position coordinates of four small-mass graphs determined after the four vertexes are respectively searched are obtained, the areas of the four small-mass graphs are compared, and one block with the smallest area is selected as the final zoning position to obtain a zoning mass graph Q' (m, N). Then, based on the position coordinates of the partitioned quality map Q' (m, n), the quality map Q (m, n) is divided into four rectangular regions, and each region is reserved with an overlapping portion of L rows and L columns, which can be used for subsequent splicing and fusion, so as to determine the coordinate positions of the four regions, as shown in fig. 8.
Step 5.3: calculating a mass ratio r' of each region, and judging the region as a low-quality region when the mass ratio of the region is more than or equal to the mass ratio x r of the mass diagram (x is a random index); and when the area quality ratio is smaller than the quality area quality ratio, judging the area as a high quality area. Finally, the position coordinates of the four regions divided into the quality diagram Q (m, n) are brought into the phase diagram, the phase diagram is divided into four small phase diagrams with different qualities according to the dividing coordinates, namely, a high-quality region and a low-quality region are divided according to the mode shown in the figure 5(a), and are divided according to the dividing mode in the phase diagram, as shown in the figure 5(b), and are divided according to the positions, as shown in the figure 5(c), and finally, the wrapping phases of the different quality regions are respectively unwrapped, as shown in the figure 5 (d);
step 5.4: after unwrapping, the four unwrapped images are merged and fused to form the final phase image, which is shown in fig. 6 (d). The invention maintains the accuracy of the unwrapped phase obtained by solving the high-quality area, enhances the accuracy of the phase obtained by solving the low-quality area, effectively inhibits the propagation of errors and greatly optimizes the effect of the finally obtained unwrapped phase diagram.
In order to compare the unwrapping effect of the partitioned composite unwrapping method based on the digital holographic microscopy imaging system, fig. 6 shows the effect of the phase diagram obtained by partitioned composite unwrapping and the phase diagram obtained by the traditional unwrapping method. As can be seen from the schematic diagram of mass division of the partition in fig. 5, the conventional mass diagram has an obvious distribution of high-quality and low-quality regions, wherein it can be obviously seen from the mass division position shown in fig. 5(b), the bright regions are low-quality regions (upper left, lower left, and lower right), and the dark regions are high-quality regions (upper right), and then the wrapped phase diagram (for the accuracy of the phase, which needs to be phase-deblurred before) is divided according to the schematic diagram of mass division to obtain the wrapped phase diagrams of the four sub-regions, as shown in fig. 6 (d). Since the unwrapping method based on the minimum norm has the characteristic of locally extending parts, if the method is adopted, the most accurate phase information cannot be obtained in the low-quality region, as shown in fig. 6 (a). The unwrapping method based on path tracking can better cope with low-quality areas, but if the phase is too large, the error propagation is more obvious, as shown in fig. 6(b), the phase diagram obtained by the branch-cut method retains more phase information but the central area has a more obvious wire-pulling phenomenon, as shown in fig. 6(c), the phase diagram obtained by the mass-guided graph method has an obvious dark and bright area, and the error cannot be well suppressed. In order to inhibit the propagation of errors, the method carries out partitioned composite unwrapping on the wrapped phase diagram corresponding to the low-quality area, carries out branch-cut unwrapping on the wrapped phase diagram corresponding to the low-quality area, carries out weighted least square unwrapping on the phase wrapped diagram corresponding to the high-quality area, and then splices and fuses the phase wrapped phase diagram, as shown in fig. 6(d), the wrapped phase accuracy is greatly optimized and unwrapped.
According to the method, only the areas are divided according to the quality factors of the quality images, the high-quality areas and the low-quality areas are respectively unwrapped, the accuracy of the unwrapped phase is obviously improved, the influence of 'dead spots' on an unwrapping algorithm can be reduced, the holographic phase information is correct, and meanwhile the integrity of the holographic phase information is greatly improved.
In order to verify the application effect of the partitioned composite unwrapping method based on the digital holographic microscopic imaging system, the solved Hela cell phase information graph is further subjected to numerical processing, and the cell segmentation is performed by using a flood back-off algorithm, as shown in FIG. 7, the solved phase information of the digital holographic cell graph is accurate, the cell segmentation can be performed well, and the cell morphology is retained to the greatest extent.
Claims (9)
1. A zoning composite phase unwrapping method based on digital holographic microscopic imaging is characterized by comprising the following specific steps:
step 1, collecting an original holographic image of a tested sample by adopting a microscopic imaging system built by a transmission type holographic optical path;
step 2, Fourier transforming the original holographic image into a spectrogram, extracting a +1 level spectrum, moving the spectrum to the center for inverse Fourier transformation to obtain a holographic image, focusing the holographic image by using an angular spectrum propagation principle, and extracting a holographic phase diagram;
step 3, performing phase dephasing on the holographic phase diagram by adopting a PCA principal component analysis method to obtain the holographic phase diagram of the sample after dephasing;
step 4, unwrapping based on path transmission, and calculating phase derivative variance according to the pixel size of the holographic phase diagram after phase distortion is eliminated to obtain a quality diagram;
and 5, dividing the distribution of the quality map according to the quality factor of the quality map, unwrapping different quality areas by adopting different phase unwrapping algorithms, and fusing and splicing the quality areas to obtain unwrapped phases.
2. The zoned composite phase unwrapping method based on digital holographic microscopy imaging as claimed in claim 1, wherein fourier transform of the original holographic image into a spectrogram is performed, and the specific method for extracting +1 level spectrum is as follows:
the light intensity of the original holographic image is:
in the formula IH(x, y) is the intensity information of the original image, IO(x, y) is the object light intensity information, IR(x, y) is the reference light intensity information, phiO(x, y) is object light phase information, phiRAnd (x, y) is reference light phase information.
And carrying out fast Fourier transform on the light intensity of the holographic image to obtain:
in the formula F [ ■]Denotes a fast Fourier transform, α (x, y) ═ IO(x,y)+IR(x, y) represents object light IO(x, y) and reference light IR(x, y) a direct current component after the correlation, a ═ F [ α (x, y)]The direct current component of the two phases of dry light is subjected to fast Fourier transform to obtain a result;
in order to be a real image,indicating the intensity of the real image, exp [ i φ ]O(x,y)-iφR(x,y)]Is the argument distribution of the real image, i denotes the imaginary unit, F+1Is the +1 order spectrum corresponding to the real image;
the +1 order spectrum is extracted separately from the fourier spectrogram directly by clipping the image.
4. The partitioned composite phase unwrapping method based on digital holographic microscopy imaging as claimed in claim 1, wherein the specific method for performing phase dephasing on the holographic phase map by PCA principal component analysis comprises:
(1) holographic phase pattern phiH(X, y) is denoted as matrix X;
(2) carrying out zero average value transformation on each row of the matrix X again;
(4) computing eigenvectors e of the covariance matrix C1,e2……etAnd a characteristic value λ, t ═ 1,2 … … N × M;
(5) sorting the eigenvalues from big to small, and selecting the largest k of the eigenvalues; respectively taking the corresponding k eigenvectors as row vectors to form an eigenvector matrix P; y is PX which is the data matrix after the dimension of the matrix X is reduced to k dimension;
(6) and performing least square fitting on the obtained reduced-dimension data matrix Y to remove the phase difference.
5. The zoned composite phase unwrapping method based on digital holographic microscopy imaging as claimed in claim 1, wherein the phase derivative variance is calculated by the formula:
wherein (m, n) is a pixel gradient range of the subscript value (i, j) of each pixel point in a square having a size represented by k × k pixels in the central range,andis the gradient of the phase difference between the pixel wrap points of the phase along two different directions x and y respectively,andrespectively representing two averages of two phase-deflector function within the window. According to the formula, a quality map Q (M, N) corresponding to the holographic phase map is obtained, wherein the quality map Q (M, N) is a two-dimensional data matrix map with M rows and N columns, each element of the matrix is a quality factor, and the quality factor is a value of [0, 1]]The numerical value of (c).
6. The zoning composite phase unwrapping method based on digital holographic microscopy imaging as claimed in claim 1, wherein the specific method for zoning the distribution of the quality map according to the quality factor of the quality map is as follows:
setting a quality threshold value, and determining the quality ratio of the quality factor in the quality map which is greater than or equal to the set threshold value;
in the quality map Q (M, N) of M multiplied by N, small rectangles of a multiplied by b are searched from four top corners of the quality map to the center, and a is more than or equal to 1 and less than or equal to M; b is more than or equal to 1 and less than or equal to N, when the mass ratio r 'of the small rectangle is more than or equal to x r, x is less than or equal to a random index, the searching is stopped, four small-mass diagram position coordinates determined after four vertexes are respectively searched are obtained, the areas of the four small-mass diagrams are compared, and one block with the smallest area is selected as a final zoning position to obtain a zoning mass diagram Q' (m, N); dividing the quality map Q (m, n) into four rectangular areas according to the position coordinates of the partitioned quality map Q' (m, n), reserving the overlapped part of L rows and L columns of each area, and determining the coordinate positions of the four areas;
respectively calculating the mass ratio of the four regions, and judging the region as a low-mass region when the mass ratio of the regions is more than or equal to the mass ratio x r of the mass diagram; and when the area quality ratio is smaller than the quality area quality ratio, judging the area as a high quality area.
7. The partitioned composite phase unwrapping method based on digital holographic microscopy imaging as claimed in claim 6, wherein different quality regions are unwrapped by using different phase unwrapping algorithms, and the specific method is as follows: the quality diagram Q (m, n) is divided into position coordinates of four areas and is brought into the phase diagram, the phase diagram is divided into four small phase diagrams with different qualities according to the dividing coordinates, the small phase diagrams are unwrapped in the low-quality areas by adopting a path tracking method, and the unwrapped phase diagrams in the high-quality areas are unwrapped by adopting a minimum norm method to obtain unwrapped phase diagrams of four different quality areas.
8. The zoned composite phase unwrapping method based on digital holographic microscopy imaging as recited in claim 6, wherein a quality threshold is set to 0.5.
Priority Applications (1)
Application Number | Priority Date | Filing Date | Title |
---|---|---|---|
CN202110595634.6A CN113393472B (en) | 2021-05-29 | 2021-05-29 | Zoning composite phase unwrapping method based on digital holographic microscopic imaging |
Applications Claiming Priority (1)
Application Number | Priority Date | Filing Date | Title |
---|---|---|---|
CN202110595634.6A CN113393472B (en) | 2021-05-29 | 2021-05-29 | Zoning composite phase unwrapping method based on digital holographic microscopic imaging |
Publications (2)
Publication Number | Publication Date |
---|---|
CN113393472A true CN113393472A (en) | 2021-09-14 |
CN113393472B CN113393472B (en) | 2022-08-16 |
Family
ID=77619521
Family Applications (1)
Application Number | Title | Priority Date | Filing Date |
---|---|---|---|
CN202110595634.6A Active CN113393472B (en) | 2021-05-29 | 2021-05-29 | Zoning composite phase unwrapping method based on digital holographic microscopic imaging |
Country Status (1)
Country | Link |
---|---|
CN (1) | CN113393472B (en) |
Cited By (2)
Publication number | Priority date | Publication date | Assignee | Title |
---|---|---|---|---|
CN115078413A (en) * | 2022-06-30 | 2022-09-20 | 广西大学 | Novel phase shift electronic holography based on transmission electron microscope |
WO2024055602A1 (en) * | 2022-09-13 | 2024-03-21 | 南京理工大学 | Lens-free single-frame phase recovery method based on partially coherent light-emitting diode illumination |
Citations (3)
Publication number | Priority date | Publication date | Assignee | Title |
---|---|---|---|---|
CN103063155A (en) * | 2012-12-12 | 2013-04-24 | 浙江师范大学 | Rapid package removal method of digital microscopic holographic phase diagram |
EP3163251A1 (en) * | 2015-10-30 | 2017-05-03 | KLA-Tencor Corporation | Method and system for regional phase unwrapping with pattern-assisted correction |
CN111561877A (en) * | 2020-04-24 | 2020-08-21 | 西安交通大学 | Variable resolution phase unwrapping method based on point diffraction interferometer |
-
2021
- 2021-05-29 CN CN202110595634.6A patent/CN113393472B/en active Active
Patent Citations (3)
Publication number | Priority date | Publication date | Assignee | Title |
---|---|---|---|---|
CN103063155A (en) * | 2012-12-12 | 2013-04-24 | 浙江师范大学 | Rapid package removal method of digital microscopic holographic phase diagram |
EP3163251A1 (en) * | 2015-10-30 | 2017-05-03 | KLA-Tencor Corporation | Method and system for regional phase unwrapping with pattern-assisted correction |
CN111561877A (en) * | 2020-04-24 | 2020-08-21 | 西安交通大学 | Variable resolution phase unwrapping method based on point diffraction interferometer |
Cited By (2)
Publication number | Priority date | Publication date | Assignee | Title |
---|---|---|---|---|
CN115078413A (en) * | 2022-06-30 | 2022-09-20 | 广西大学 | Novel phase shift electronic holography based on transmission electron microscope |
WO2024055602A1 (en) * | 2022-09-13 | 2024-03-21 | 南京理工大学 | Lens-free single-frame phase recovery method based on partially coherent light-emitting diode illumination |
Also Published As
Publication number | Publication date |
---|---|
CN113393472B (en) | 2022-08-16 |
Similar Documents
Publication | Publication Date | Title |
---|---|---|
CN113393472B (en) | Zoning composite phase unwrapping method based on digital holographic microscopic imaging | |
EP3065001B1 (en) | Holographic microscope and data processing method for high-resolution hologram image | |
Palacios et al. | 3D image reconstruction of transparent microscopic objects using digital holography | |
CN109916522B (en) | Aberration compensation method based on hologram continuation and implementation device thereof | |
CN105184295B (en) | A kind of holoscan space length extracting method based on wavelet transformation and connected domain | |
CN113418469B (en) | Spectrum confocal scanning common-path digital holographic measurement system and measurement method | |
CN106054570A (en) | Method for realizing large-phase reconstruction of single digital hologram by adopting intensity transmission equation | |
CN114660060A (en) | Macroscopic Fourier laminated super-resolution imaging method based on matrix scanning | |
He et al. | Robust phase aberration compensation in digital holographic microscopy by self-extension of holograms | |
Wu et al. | Autofocusing algorithm for pixel-super-resolved lensfree on-chip microscopy | |
CN114241072A (en) | Laminated imaging reconstruction method and system | |
Chen et al. | High-resolution reconstruction of spectrum-overlapped off-axis holography by deflecting reference beam of Gaussian symmetry | |
WO2021003380A1 (en) | Calibration-free phase shifting procedure for self-interference holography | |
CN117053716A (en) | Automatic detection method for outline parameters of circular aperture interferogram | |
Zhang et al. | Deep-learning-enhanced digital holographic autofocus imaging | |
US12112452B2 (en) | Holographic ultra resolution imaging | |
KR100717414B1 (en) | Off-axis illumination direct-to-digital holography | |
Abdelsalam et al. | Digital holographic shape measurement using Fizeau microscopy | |
CN110514137B (en) | Phase unwrapping method, device, system, computer equipment and storage medium | |
Buitrago-Duque et al. | Speckle Reduction in Digital Holographic Microscopy by Physical Manipulation of the Pupil Function | |
Jin et al. | Spatial consistency calibration based on phase difference minimization for parallel slightly off-axis digital holographic microscopy | |
CN110715931B (en) | Automatic detection method and detection device for defects of transparent sample | |
Ilhan et al. | Autofocusing in digital holography | |
Sato et al. | Lens-less holographic microscope with high resolving power and no-distortion | |
Guan et al. | Fast autofocusing in off-axis digital holography based on search region segmentation and dichotomy |
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 |