Microscopy illumination error correction method with spectrum repair and adaptive optimization
Patent Information
- Application Number
- CN202611109408.1
- Authority / Receiving Office
- CN · China
- Patent Type
- Applications(China)
- Current Assignee / Owner
- Filing Date
- 2026-07-24
- Publication Date
- 2026-10-27
AI Technical Summary
[0007]针对上述问题,本发明提供一种频谱修复与自适应寻优的显微照明误差校正方法,以解决正交轴频谱伪影干扰圆形频谱特征提取、频谱修复造成非直流频谱能量变化,以及固定搜索步长难以兼顾搜索范围和定位精度的技术问题
[0060](1)本发明通过周期性平滑分解减弱图像相对边界不连续性造成的频谱伪影,并对正交轴区域的异常频谱像素进行邻域插值,能够在抑制十字型伪影的同时保留正交轴附近的频谱连续性,降低伪影对光瞳圆形特征提取的干扰。
Smart Images

Figure CN122989422A_ABST
Abstract
Description
Technical Field
[0001] This invention relates to the fields of computational imaging and microscopic imaging technology, and in particular to a method for correcting microscopic illumination errors through spectral restoration and adaptive optimization. Background Technology
[0002] Fourier Ptychographic Microscopy (FPM) is a computational imaging technique that integrates synthetic aperture and phase retrieval. This technique typically uses a programmable light source array to sequentially illuminate the sample from different spatial locations, causing high-frequency information in the sample's spectrum that is outside the objective lens passband to move into the objective lens passband under tilted illumination. After the camera acquires low-resolution images corresponding to different illumination positions, high-resolution amplitude and phase information of the sample is reconstructed through spectral stitching and iterative phase retrieval.
[0003] In Fourier layered microscopy, the spatial position of each light source unit is used to determine the corresponding illumination wave vector and the position of the object's spectral sub-region. Due to factors such as light source array installation offset, array rotation, device manufacturing errors, and assembly errors, geometric misalignment can easily occur between the actual and nominal positions of the light source units. If spectral stitching is still performed according to the nominal positions, positional deviations will occur between adjacent spectral sub-regions, leading to wrinkles, ghosting, and phase distortion in the iterative reconstruction results.
[0004] Existing illumination error correction methods include geometric calibration based on the mechanical dimensions of the imaging device, and inversion of light source position parameters based on low-resolution images or their frequency domain characteristics. Mechanical calibration relies on measuring equipment and accuracy, making it difficult to adequately compensate for overall translational, rotational, and assembly errors of the light source array. Frequency domain-based correction methods can utilize the circular pupil feature formed by tilted illumination to determine the illumination wave vector; however, discontinuities between relative boundaries in low-resolution images can create cross-shaped artifacts distributed along two orthogonal frequency axes in the discrete Fourier transform results, thus interfering with the extraction of the circular pupil feature and the calculation of light source position parameters.
[0005] While conventional periodic smoothing decomposition methods can reduce cross-shaped artifacts caused by image boundary discontinuities, they may cause non-DC spectral energy attenuation during the extraction of periodic components. Meanwhile, when using a fixed search step size to solve for the translation and rotation errors of the light source array, a larger search step size, while expanding the search range, makes it difficult to achieve high local positioning accuracy; a smaller search step size, while beneficial for fine positioning, is easily affected by initial parameters, making it difficult to balance global search range and local positioning accuracy.
[0006] Therefore, how to repair the abnormal spectrum in the orthogonal axis region while preserving the effective spectral information, compensate for the non-DC spectral energy changes caused by the spectral repair, and adjust the search step size of the translation and rotation parameters according to the optimization state are the technical problems that need to be solved in the process of microscopic illumination error correction. Summary of the Invention
[0007] To address the aforementioned problems, this invention provides a microscopic illumination error correction method based on spectrum repair and adaptive optimization, which solves the technical problems of orthogonal axis spectrum artifact interference in circular spectrum feature extraction, spectrum repair causing non-DC spectrum energy changes, and the difficulty in balancing search range and positioning accuracy with a fixed search step size.
[0008] To address the aforementioned technical problems, this invention provides a method for correcting microscopic illumination errors through spectral restoration and adaptive optimization, comprising the following steps:
[0009] S1. Spectrum Restoration and Energy Compensation: Different light source units in the light source array are controlled to sequentially illuminate the sample, and low-resolution images corresponding to each illumination position are acquired; a two-dimensional Fourier transform is performed on the low-resolution image to obtain the original spectrum; a periodic smoothing decomposition is performed on the low-resolution image to obtain periodic components, and a two-dimensional Fourier transform is performed on the periodic components to obtain the periodic spectrum; neighborhood interpolation is performed on abnormal spectrum pixels located in the orthogonal axis region of the periodic spectrum; an energy compensation coefficient is determined based on the non-DC spectrum energy of the original spectrum and the non-DC spectrum energy of the interpolated spectrum; and energy compensation is performed on the non-DC component of the interpolated spectrum while keeping the DC component unchanged to obtain a restored spectrum feature image.
[0010] S2. Circular Spectral Feature Extraction: The recovered spectral feature image is input into a pre-trained feature extraction network to obtain circular spectral features representing the pupil position.
[0011] S3. Adaptive optimization of illumination geometry parameters: Set translational search step size and rotational search step size respectively, and generate candidate geometric misalignment parameters based on the current translational search step size and the current rotational search step size; determine the corresponding simulated pupil position according to each candidate geometric misalignment parameter, and generate a simulated circular spectral feature with the same data format as the circular spectral feature obtained in step S2; calculate the feature error between the simulated circular spectral feature and the circular spectral feature obtained in step S2, update the current optimal geometric misalignment parameter according to the feature error, and decrease the translational search step size and the rotational search step size respectively when a smaller feature error is not obtained after a preset number of consecutive steps.
[0012] S4. Illumination error correction and Fourier layered microscopic reconstruction: Based on the obtained geometric misalignment parameters, the spatial position of each light source unit and the corresponding illumination wave vector are corrected, and Fourier layered microscopic reconstruction is performed based on the corrected illumination wave vector.
[0013] Further, in step S1, the periodic smooth decomposition satisfies:
[0014] ;
[0015] ;
[0016] in, For the low-resolution image, It is a periodic component. To smooth the components, The boundary term is constructed based on the boundary discontinuities of the low-resolution image. It is a two-dimensional discrete Laplace operator.
[0017] Furthermore, the two-dimensional discrete Laplace operator satisfies:
[0018] ;
[0019] in, These represent the smoothing component values on the left, right, top, and bottom sides of the current pixel, respectively.
[0020] Furthermore, the boundary item Initialize to zero, and assign values to its boundary positions according to the following formula:
[0021] ;
[0022] ;
[0023] ;
[0024] ;
[0025] in, The number of rows in the low-resolution image. The number of columns in the low-resolution image; and Used to characterize discontinuities between relative row boundaries and Used to characterize discontinuities between relative column boundaries, the symbol " "" indicates that the update is performed based on the current boundary item value.
[0026] Furthermore, performing a two-dimensional Fourier transform on both sides of the discrete Poisson equation yields:
[0027] ;
[0028] Based on the frequency domain properties of the two-dimensional discrete Laplace operator, the frequency domain representation of the smoothing component is as follows:
[0029] ;
[0030] in, , This represents a discrete index in the vertical frequency direction. Discrete indexes in the horizontal frequency direction. For smoothing components The frequency domain representation, This represents a two-dimensional Fourier transform.
[0031] The periodic spectrum satisfies:
[0032] ;
[0033] in, For periodic components The corresponding periodic spectrum.
[0034] Furthermore, neighborhood interpolation is performed on the abnormal spectrum pixels located in the orthogonal axis region of the periodic spectrum to satisfy:
[0035] ;
[0036] in, The coordinates in the interpolated spectrum are Spectral values, The coordinates of the spectrum center are used; when a spectrum pixel is located on the vertical frequency axis passing through the spectrum center but not at the spectrum center, the average interpolation is performed using two adjacent spectrum pixels along the row direction; when a spectrum pixel is located on the horizontal frequency axis passing through the spectrum center but not at the spectrum center, the average interpolation is performed using two adjacent spectrum pixels along the column direction; when a spectrum pixel is located at the spectrum center, the DC component of the spectrum center remains unchanged.
[0037] Furthermore, the non-DC spectral energy of the original spectrum and the interpolated spectrum are calculated separately:
[0038] ;
[0039] ;
[0040] in, The raw spectrum of the low-resolution image. This refers to the original non-DC spectrum energy after deducting the DC energy at the center of the spectrum. The non-DC spectrum energy is the interpolated energy after deducting the DC energy at the center of the spectrum.
[0041] The energy compensation coefficient is determined based on the original spectrum non-DC energy and the interpolated spectrum non-DC energy.
[0042] ;
[0043] in, This is the energy compensation coefficient.
[0044] Furthermore, the DC component of the interpolated spectrum is kept unchanged, and the non-DC component of the interpolated spectrum is compensated using the energy compensation coefficient:
[0045] ;
[0046] in, To recover the coordinates in the spectral feature image The spectral values.
[0047] Furthermore, the feature extraction network is a U-Net network. The U-Net network includes an encoder, a bottleneck layer, and a decoder; the number of channels in the convolutional layers of the encoder is 8, 16, 32, and 64 respectively, the bottleneck layer has 128 channels, and the number of channels in the convolutional layers of the decoder is 64, 32, 16, and 8 respectively; the U-Net network uses 3×3 convolutional kernels and downsampling is performed through 2×2 max pooling.
[0048] Furthermore, the input to the U-Net network is a single-channel recovered spectral feature image with a size of 256×256; the training samples include 300 recovered spectral feature images, which are divided into training set and validation set in a 7:3 ratio; the training process uses binary cross-entropy as the loss function, average absolute error as the evaluation index, and backpropagation is used to update the network parameters.
[0049] Furthermore, the translation search step size is initialized to 200 position offset coordinate units, the rotation search step size is initialized to 0.5 degrees, the current optimal first translation amount, the current optimal second translation amount, and the current optimal rotation amount are all initialized to zero, the local unoptimized counter and the cumulative unoptimized counter are both initialized to zero, and the historical minimum feature error is initialized to infinity.
[0050] Furthermore, based on the current translation search step size and the current rotation search step size, multiple sets of candidate first translation offset, candidate second translation offset, and candidate rotation offset are generated using uniformly distributed random numbers; the candidate first translation offset, the candidate second translation offset, and the candidate rotation offset are respectively superimposed on the current optimal first translation, the current optimal second translation, and the current optimal rotation to obtain multiple sets of candidate geometric misalignment parameters.
[0051] Further, the illumination wave vector coordinates are updated based on the candidate geometric misalignment parameters of each group, the corresponding simulated pupil position is determined according to the updated illumination wave vector coordinates, and a simulated circular spectral feature with the same image size and data format as the circular spectral feature obtained in step S2 is generated according to the simulated pupil position; the mean absolute error between the simulated circular spectral feature and the circular spectral feature obtained in step S2 is calculated.
[0052] When the mean absolute error corresponding to the candidate geometric misalignment parameter is less than the historical minimum feature error, the candidate geometric misalignment parameter is updated to the current optimal geometric misalignment parameter, the corresponding mean absolute error is updated to the historical minimum feature error, and both the local unoptimized counter and the cumulative unoptimized counter are cleared to zero.
[0053] When the mean absolute error corresponding to the candidate geometric misalignment parameter is not less than the historical minimum feature error, the current optimal geometric misalignment parameter is kept unchanged, and the count values of the local unoptimized counter and the cumulative unoptimized counter are increased respectively.
[0054] When the local unoptimized counter reaches 5, the translation search step size is reduced to 40% of its original value, the rotation search step size is reduced to 55% of its original value, and the local unoptimized counter is reset to zero. The cumulative unoptimized counter maintains its current value during the search step size reduction. Subsequently, candidate geometric misalignment parameters are generated based on the reduced translation and rotation search step sizes.
[0055] Furthermore, when the historical minimum feature error is not greater than a preset error threshold, or when the count value of the cumulative unoptimized counter reaches 50, the optimization is terminated, and the current optimal geometric misalignment parameter is output as the geometric misalignment parameter of the light source array.
[0056] Furthermore, the geometric misalignment parameters of the light source array are compensated to the nominal spatial position of each light source unit to obtain the corrected spatial position of each light source unit; the corresponding illumination wave vector is determined according to the corrected spatial position, and an adaptive phase retrieval algorithm is used to iteratively reconstruct the low-resolution image corresponding to different illumination positions to obtain the amplitude and phase information of the sample.
[0057] Further, in step S3, the geometric misalignment parameters include a first translation and a second translation generated by the light source array along a first and a second direction that are orthogonal to each other in the array plane, and a rotation generated by the light source array around the optical axis; taking the current optimal geometric misalignment parameter as the center, candidate translation parameters are generated along the positive and negative directions of the first and second directions respectively according to the current translation search step size, and candidate rotation parameters are generated along the positive and negative directions around the optical axis according to the current rotation search step size, and multiple sets of candidate geometric misalignment parameters are formed by the candidate translation parameters and the candidate rotation parameters.
[0058] Further, in step S4, the nominal spatial position of each light source unit is transformed by translation and rotation according to the first translation amount, the second translation amount, and the rotation amount to obtain the corrected spatial position of each light source unit; the illumination wave vector corresponding to each light source unit is recalculated according to the corrected spatial position, and the spectral sub-region corresponding to each low-resolution image is determined according to the recalculated illumination wave vector; phase recovery iteration is performed on each spectral sub-region in sequence to obtain the high-resolution amplitude information and phase information of the sample.
[0059] Compared with the prior art, the present invention has at least the following beneficial effects:
[0060] (1) The present invention reduces the spectral artifacts caused by the discontinuity of the relative boundary of the image by periodic smooth decomposition, and performs neighborhood interpolation on the abnormal spectral pixels in the orthogonal axis region. It can suppress cross-shaped artifacts while preserving the spectral continuity near the orthogonal axis, and reduce the interference of artifacts on the extraction of the circular feature of the pupil.
[0061] (2) The present invention determines the energy compensation coefficient based on the non-DC spectral energy of the original spectrum and the interpolated spectrum. Under the condition of keeping the DC component of the interpolated spectrum unchanged, the non-DC component of the interpolated spectrum is compensated for energy, which can reduce the effective spectral energy attenuation caused by the interpolation of the neighborhood of abnormal spectrum pixels, and make the recovered spectral feature image more accurately represent the pupil position.
[0062] (3) The present invention sets translational search step size and rotational search step size respectively, and reduces the two types of search step size according to different attenuation ratios when a smaller feature error is not obtained after a series of preset times, so that the parameter search process can gradually transition from a coarse search over a larger range to a fine search over a smaller range.
[0063] (4) This invention combines spectrum repair, pupil feature extraction, optimization of geometric misalignment parameters of light source array and Fourier stacked microscopic reconstruction, which can correct the translation error of the light source array along two mutually orthogonal directions in the array plane and the rotation error around the optical axis, thereby reducing the influence of light source position deviation on the reconstructed amplitude and phase information. Attached Figure Description
[0064] To more clearly illustrate the technical solutions in the embodiments of the present invention or the prior art, the drawings used in the description of the embodiments or the prior art will be briefly introduced below. Obviously, the drawings described below are only some embodiments of the present invention. For those skilled in the art, other drawings can be obtained based on these drawings without creative effort.
[0065] Figure 1This is a flowchart illustrating a microscopic illumination error correction method based on spectrum restoration and adaptive optimization according to the present invention.
[0066] Figure 2 This is a schematic diagram of the processing framework for spectrum restoration, circular spectrum feature extraction, adaptive optimization of illumination geometric parameters, and microscopic image reconstruction of the present invention.
[0067] Figure 3 This is a comparison chart of the spectral processing and circular spectral feature extraction results of the method described in this invention and a comparative method, wherein... Figure 3 (a) in the figure is the recovered spectrum feature map obtained after spectrum restoration and energy compensation using the method described in this invention. Figure 3 (b) in the text is based on Figure 3 (a) shows the obtained circular spectral features extracted from the data. Figure 3 (c) in the image is the spectral feature map obtained after processing using the SAA comparison method. Figure 3 (d) in the text is based on Figure 3 (c) Extracts the obtained circular spectral features.
[0068] Figure 4 This is a comparison chart showing the reconstruction results of blood smear samples using the method described in this invention and different Fourier transform microscopic reconstruction methods. Figure 4 Image (a) is a low-resolution image captured when the central illumination unit is lit. Figure 4 (a1) in the diagram represents the selected reconstruction region. Figure 4 (a2) in the image is a magnified view of the reconstructed area. Figure 4 (b1) to (g1) are amplitude reconstruction images obtained by the AA method, PC method, AA-C method, SC method, SAA method and the method described in this invention, respectively; Figure 4 In the image, (b2) to (g2) represent the phase reconstruction images obtained by the corresponding methods. Figure 4 (b3) to (g3) are magnified views of the preset local regions in the corresponding phase reconstruction images.
[0069] Figure 5 This is a schematic diagram showing low-resolution images and their spectral characteristics corresponding to different lighting locations. Figure 5 (a1) in the diagram is a schematic diagram of the vertical lighting mode of the central lighting unit. Figure 5 (a2) in the diagram is a schematic diagram of the tilted lighting mode for a non-central lighting unit. Figure 5 In the image, (b1) and (b2) are low-resolution sample images acquired under the corresponding illumination modes, respectively. Figure 5 In the diagram, (c1) and (c2) are the spectral feature maps obtained by two-dimensional Fourier transform of the local regions of the corresponding low-resolution sample images. Detailed Implementation
[0070] To make the objectives, technical solutions, and advantages of this application clearer, the following detailed description is provided in conjunction with the accompanying drawings and embodiments. It should be understood that the specific embodiments described herein are merely illustrative and not intended to limit the scope of this application. All other embodiments obtained by those skilled in the art based on the embodiments of this application without inventive effort are protected by this application.
[0071] It should be noted that the terms "comprising," "including," and "having," and any variations thereof, in the specification, claims, and accompanying drawings of this application, are intended to cover non-exclusive inclusion. For example, a process, method, terminal, product, or device that includes a series of steps or units is not limited to the listed steps or units, but may optionally include steps or units not listed, or may optionally include other steps or units inherent to these processes, methods, products, or devices. Terms such as "first" and "second" in the claims, specification, and accompanying drawings of this application, as well as relational terms, are used merely to distinguish one entity / operation / object from another entity / operation / object, and do not necessarily require or imply any such immediate relationship or order between these entities / operations / objects.
[0072] In this document, the term "embodiment" means that a particular feature, structure, or characteristic described in connection with an embodiment may be included in at least one embodiment of this application. The appearance of this phrase in various places throughout the specification does not necessarily refer to the same embodiment, nor is it a mutually exclusive, independent, or alternative embodiment. It will be explicitly and implicitly understood by those skilled in the art that the embodiments described herein can be combined with other embodiments.
[0073] In a Fourier layered microscopy system, different LED units in a programmable LED array illuminate the sample at different angles, causing different sub-regions of the sample's spectrum to move into the objective lens passband. The camera sequentially acquires low-resolution images corresponding to different illumination positions, and obtains high-resolution amplitude and phase information of the sample through spectral stitching and phase recovery. The actual spatial position of the LED unit determines the corresponding illumination wave vector and the position of the spectral sub-region. When there is an overall translational or rotational error in the LED array, reconstruction according to the nominal position will cause misalignment of the spectral sub-region, resulting in wrinkles, ghosting, and phase distortion in the reconstructed image.
[0074] Existing methods for LED position correction based on low-resolution image spectral features are susceptible to spectral cross-shaped artifacts. These artifacts penetrate or obscure the circular spectral features used to characterize the pupil position, reducing the accuracy of circular spectral feature extraction. While the Periodic Component and Smooth Component Image Decomposition (PPSID) method can mitigate these cross-shaped artifacts, it may cause non-DC spectral energy attenuation during the extraction of periodic components. Furthermore, using a fixed search step size to solve for the translation and rotation errors of the LED array makes it difficult to simultaneously achieve a large global search range and high local search accuracy.
[0075] In view of this, embodiments of this application provide a method for correcting microscopic illumination errors through spectrum restoration and adaptive optimization. On one hand, periodic smoothing decomposition reduces boundary discontinuities in the spectrum of a low-resolution image. Neighborhood interpolation is performed on abnormal spectral pixels in orthogonal axis regions, and energy compensation coefficients are determined based on the non-DC spectral energy of the original and interpolated spectra. Energy compensation is then performed on the non-DC component of the interpolated spectrum while maintaining the DC component unchanged. On the other hand, the restored spectral features are input into a pre-trained U-Net network to extract circular spectral features associated with the LED illumination location. Furthermore, translational and rotational search step sizes are set separately, and when a smaller error is not obtained after multiple iterations, the two types of search step sizes are reduced according to different attenuation ratios, gradually transitioning the parameter search process from a global coarse search to a local fine search. Finally, the obtained LED array geometric misalignment parameters are compensated in the Fourier layered microscopic reconstruction process to correct the illumination wave vector corresponding to different illumination locations.
[0076] To illustrate the technical solution of this invention, the following is combined with... Figures 1 to 5 Please provide a detailed explanation.
[0077] Please refer to Figures 1 to 4 This application provides a method for correcting microscopic illumination errors through spectrum restoration and adaptive optimization.
[0078] Specifically, the above-mentioned spectrum repair and adaptive optimization microscopic illumination error correction method includes the following steps S1 to S4.
[0079] Step S1: Spectrum repair and energy compensation.
[0080] Different light source units in the light source array are controlled to sequentially illuminate the sample, acquiring low-resolution images corresponding to each illumination position. A two-dimensional Fourier transform is performed on the low-resolution image to obtain the original spectrum. The low-resolution image is then subjected to periodic smoothing decomposition to obtain periodic and smoothed components. A two-dimensional Fourier transform is performed on the periodic components to obtain the periodic spectrum. Abnormal spectral pixels located in the orthogonal axis region of the periodic spectrum are interpolated in their neighborhood, and energy compensation is performed on the non-DC spectral components after neighborhood interpolation to obtain a recovered spectral feature image. The periodic smoothing decomposition of the low-resolution image satisfies the following conditions:
[0081] ;
[0082] in, For the low-resolution image, It is a periodic component. For smooth components. For periodic components. Satisfying periodic boundary conditions, smooth components Used to characterize discontinuities between relative boundaries in low-resolution images, and smoothing components The mean is zero. Smooth components. Satisfies the discrete Poisson equation:
[0083] ;
[0084] in, For a two-dimensional discrete Laplace operator, This is the boundary term constructed based on the boundary discontinuities of the low-resolution image. The two-dimensional discrete Laplacian operator is represented as:
[0085] ;
[0086] in, These represent the smoothing component values on the left, right, top, and bottom sides of the current pixel, respectively.
[0087] In this embodiment, the boundary terms are first... Initialize to zero, then update the boundary terms based on the pixel differences between the relative boundaries of the low-resolution image:
[0088] ;
[0089] ;
[0090] ;
[0091] ;
[0092] in, The number of rows in the low-resolution image. The number of columns in the low-resolution image; and Used to characterize discontinuities between relative row boundaries and Used to characterize discontinuities between relative column boundaries, the symbol " "" indicates that the update is performed based on the current boundary term value. A two-dimensional Fourier transform of both sides of equation (2) yields:
[0093] ;
[0094] Using the frequency domain properties of the two-dimensional discrete Laplace operator, we obtain:
[0095] ;
[0096] in: , This represents a discrete index in the vertical frequency direction. Discrete indexes in the horizontal frequency direction. For smoothing components The frequency domain representation. When When the smooth component is in frequency domain, it is represented as:
[0097] ;
[0098] in, These are the coordinates of the spectrum center. Due to the smoothing component... Since the mean of is zero, its DC component at the center of the spectrum is set as follows:
[0099] ;
[0100] The periodic spectrum corresponding to the periodic component is:
[0101] ;
[0102] in, The coordinates in the periodic spectrum are The spectral values.
[0103] Periodic smoothing decomposition can reduce spectral cross-shaped artifacts caused by boundary discontinuities. However, anomalous spectral pixels or spectral energy attenuation may still exist in the orthogonal axis region of the periodic spectrum. To maintain local spectral continuity, this embodiment uses the mean interpolation of adjacent spectral pixels with effective physical amplitudes on both sides of the orthogonal axis to interpolate the anomalous spectral pixels:
[0104] ;
[0105] in, The coordinates in the interpolated spectrum are Spectral values, The coordinates of the spectrum center are used; when a spectrum pixel is located on the vertical frequency axis passing through the spectrum center but not at the spectrum center, the average interpolation is performed using the two adjacent spectrum pixels along the row direction; when a spectrum pixel is located on the horizontal frequency axis passing through the spectrum center but not at the spectrum center, the average interpolation is performed using the two adjacent spectrum pixels along the column direction; for the coordinates of the spectrum center... For the spectral pixels, no neighborhood interpolation is performed to keep the DC component at the center of the spectrum unchanged.
[0106] Considering that neighborhood interpolation of anomalous spectral pixels may cause non-DC spectral energy loss, this embodiment calculates the non-DC spectral energy of the original spectrum and the interpolated spectrum respectively:
[0107] ;
[0108] ;
[0109] in, The raw spectrum of the low-resolution image. This refers to the original non-DC spectrum energy after deducting the DC energy at the center of the spectrum. The non-DC spectrum energy is the interpolated energy after deducting the DC energy at the center of the spectrum.
[0110] The energy compensation coefficient is determined based on the original spectrum non-DC energy and the interpolated spectrum non-DC energy.
[0111] ;
[0112] in, This is the energy compensation coefficient.
[0113] Keeping the DC component of the interpolated spectrum unchanged, and using the energy compensation coefficient Scaling the non-DC component of the interpolated spectrum yields the recovered spectral feature image, satisfying:
[0114] ;
[0115] in, To recover the coordinates in the spectral feature image Spectral values, The coordinates in the interpolated spectrum are Spectral values, This is the energy compensation coefficient. The coordinates are the center of the spectrum.
[0116] By keeping the DC component at the center of the spectrum unchanged and compensating for the energy of the non-DC components, the effective high-frequency energy attenuation caused by spectrum restoration can be reduced while preserving the macroscopic lighting background.
[0117] Please refer to Figure 3 , Figure 3 Image (a) shows the recovered spectrum feature map obtained after periodic smoothing decomposition, orthogonal axis region neighborhood interpolation, and non-DC component energy compensation using the method of this application. Figure 3 (b) in the diagram represents the pre-trained feature extraction network based on... Figure 3 (a) shows the obtained circular spectral features extracted from the data. Figure 3 (c) in the image is the spectral feature map obtained after processing using the SAA comparison method. Figure 3 (d) in the text is based on Figure 3 The circular spectral features obtained by extracting (c) in the diagram. Figure 3 As can be seen, compared with the SAA comparison method, the cross-shaped artifacts in the orthogonal axis region of the spectrum processed by the method of this application are further suppressed, and the extracted circular spectral features have a more complete outline and a more continuous region.
[0118] Step S2, circular spectral feature extraction.
[0119] It should be noted that the training process of the pre-trained feature extraction network is a pre-completed network construction stage used to obtain the pre-trained feature extraction network, and does not need to be repeated every time microscopic illumination error correction is performed. When actually performing microscopic illumination error correction, the recovered spectral feature image obtained in step S1 is input into the pre-trained feature extraction network to obtain circular spectral features representing the pupil position.
[0120] In the process of obtaining the pre-trained feature extraction network, the recovered spectral feature image obtained in step S1 can be used as training data input into the convolutional neural network. By comparing the difference between the segmentation result output by the network and the manually labeled result, the network parameters are updated using backpropagation, and the trained convolutional neural network is used as the pre-trained feature extraction network.
[0121] In this embodiment, the convolutional neural network used is the U-Net network. The dataset used to train the U-Net network includes 300 restored spectral feature images obtained after processing in step S1, which are divided into training and validation sets in a 7:3 ratio. The number of channels in the convolutional layers of the encoder part of the U-Net network is 8, 16, 32, and 64, respectively, and the bottleneck layer is set to 128 channels. The number of channels in the convolutional layers of the decoder part is 64, 32, 16, and 8, respectively. Each convolutional operation uses a 3×3 convolutional kernel, and feature map downsampling is performed through a 2×2 max pooling operation.
[0122] During training, the learning rate is set to... -4 The batch size was set to 8, the number of iterations to 500, and the Adam optimizer was used. Binary cross-entropy was used as the loss function, and the network parameters were updated via backpropagation by comparing the circular spectrum segmentation results output by the network with the manually labeled circular spectrum regions. Mean absolute error was used as the evaluation metric to characterize the deviation between the network output and the manually labeled results.
[0123] After training, network parameters with smaller validation set errors are selected as the pre-trained feature extraction network. The recovered spectral feature image corresponding to the low-resolution image to be corrected is input into the pre-trained feature extraction network, which outputs segmented circular spectral features. The circular spectral features are used to characterize the pupil position under actual lighting conditions and serve as physical prior information for subsequent optimization of LED array geometric misalignment parameters.
[0124] In this embodiment, the network training process takes approximately 2 minutes, and the inference time of the pre-trained feature extraction network for a single recovered spectral feature image is approximately 5ms to 10ms. The above training and inference times are related to the computing device used and do not constitute a limitation on the scope of protection of this application.
[0125] Step S3: Adaptive optimization of lighting geometry parameters.
[0126] Using the circular spectral features extracted in step S2 as the target circular spectral features, the translational search step size and rotational search step size of the LED array geometric misalignment parameters are set; multiple sets of candidate geometric misalignment parameters are generated based on the current translational search step size and the current rotational search step size; the corresponding simulated pupil position is determined according to each candidate geometric misalignment parameter, and a simulated circular spectral feature with the same image size, coordinate system, and data format as the target circular spectral feature is generated; the mean absolute error between the simulated circular spectral feature and the target circular spectral feature is calculated, and the current optimal geometric misalignment parameters are updated according to the mean absolute error.
[0127] In this embodiment, the translation search step size is initialized to 200 position offset coordinate units, the rotation search step size is initialized to 0.5 degrees, the current optimal first translation amount, the current optimal second translation amount, and the current optimal rotation amount are all initialized to zero, the local unoptimized counter and the cumulative unoptimized counter are both initialized to zero, and the historical minimum mean absolute error is initialized to infinity.
[0128] In each iteration, with the current optimal geometric misalignment parameter as the center, and based on the current translation search step size, candidate first translation offset and candidate second translation offset are generated along the positive and negative directions of the first and second mutually orthogonal directions in the plane where the programmable LED array is located, respectively. Based on the current rotation search step size, candidate rotation offset is generated along the positive and negative directions around the optical axis.
[0129] Each candidate first translation offset, candidate second translation offset, and candidate rotation offset is superimposed on the current optimal first translation, current optimal second translation, and current optimal rotation, respectively, to form multiple sets of candidate geometric misalignment parameters corresponding to the current iteration.
[0130] For each set of candidate geometric misalignment parameters, the nominal spatial position of each LED unit in the programmable LED array is translated and rotated according to the candidate geometric misalignment parameters to obtain the candidate spatial position corresponding to each LED unit; the candidate illumination wave vector corresponding to each LED unit is calculated according to the candidate spatial position, and the simulated pupil position in the frequency domain is determined according to the candidate illumination wave vector.
[0131] A simulated circular spectral feature is generated at the simulated pupil location using the same image size, frequency domain coordinate system, and pixel representation as the target circular spectral feature output in step S2.
[0132] The mean absolute error (MAE) between the simulated circular spectral features corresponding to each candidate geometric misalignment parameter and the target circular spectral features extracted in step S2 is calculated, and this MAE is used as an evaluation index for the corresponding candidate geometric misalignment parameter. The smaller the MAE, the closer the simulated pupil position determined by the corresponding candidate geometric misalignment parameter is to the pupil position under actual lighting conditions.
[0133] From the multiple candidate geometric misalignment parameters corresponding to the current iteration, determine the candidate geometric misalignment parameter with the smallest mean absolute error, and determine its corresponding mean absolute error as the minimum mean absolute error of this round.
[0134] When the current minimum average absolute error is less than the historical minimum average absolute error, the candidate geometric misalignment parameter corresponding to the current minimum average absolute error is accepted, updated to the current optimal geometric misalignment parameter, the current minimum average absolute error is updated to the historical minimum average absolute error, and both the local unoptimized counter and the cumulative unoptimized counter are cleared to zero.
[0135] When the minimum mean absolute error of the current round is not less than the historical minimum mean absolute error, keep the current optimal geometric misalignment parameter and the historical minimum mean absolute error unchanged, and increase the count values of the local unoptimized counter and the cumulative unoptimized counter by one respectively.
[0136] After each iteration, it is determined whether the historical minimum mean absolute error is not greater than a preset error threshold, or whether the count value of the cumulative unoptimized counter reaches 50. When either of the aforementioned conditions is met, the optimization is terminated, and the current optimal geometric misalignment parameters are output as the geometric misalignment parameters of the programmable LED array.
[0137] When the aforementioned termination condition is not met, and the count value of the local unoptimized counter reaches 5, the translation search step size is reduced to 40% of its original value, the rotation search step size is reduced to 55% of its original value, and the local unoptimized counter is reset to zero. The cumulative unoptimized counter maintains its current count value during the search step size reduction. Subsequently, the next iteration continues using the reduced translation and rotation search step sizes.
[0138] By setting the translation search step size and the rotation search step size respectively, and adjusting the translation search step size and the rotation search step size according to different attenuation ratios, the optimization process of the geometric misalignment parameters of the programmable LED array can be gradually transitioned from a coarse search over a larger range to a fine search over a smaller range.
[0139] The initial values of the translation search step size and rotation search step size, the first preset number of times corresponding to the local unoptimized counter, the second preset number of times corresponding to the cumulative unoptimized counter, the attenuation ratio of the translation search step size and rotation search step size, and the preset error threshold are all preferred parameters used in this embodiment. Those skilled in the art can adjust the corresponding parameters according to the size of the programmable LED array, the LED unit spacing, the spectral image size, and the expected correction accuracy.
[0140] Step S4: Illumination error correction and Fourier layered microscopic reconstruction.
[0141] The geometric misalignment parameters of the LED array obtained in step S3 are used as physical prior information and compensated into the program in the Fourier stacked microscopic reconstruction algorithm used to generate the spatial position and illumination wave vector of the LEDs, so as to obtain the corrected spatial position and illumination wave vector of each LED unit; based on the corrected illumination wave vector, phase recovery iterative reconstruction is performed on the low-resolution images corresponding to different illumination positions to obtain the amplitude information and phase information of the sample.
[0142] In this embodiment, the AA algorithm, or Adaptive Algorithm, is used as the Fourier stacked microscopic phase retrieval strategy. Based on the first translation, second translation, and rotation values output in step S3, the nominal spatial positions of each LED unit in the programmable LED array are transformed by translation and rotation to obtain the corrected spatial positions of each LED unit, and the corresponding illumination wave vector is recalculated based on the corrected spatial positions.
[0143] In each iteration, the spectral sub-region corresponding to the current low-resolution image is determined based on the corrected illumination wave vector. The corresponding sub-region in the current high-resolution complex amplitude spectrum is then extracted using the pupil function and transformed to the spatial domain. The spatial domain amplitude is updated based on the low-resolution image acquired by the camera, and the updated spatial domain result is then transformed to the frequency domain to update the high-resolution complex amplitude spectrum. This process is repeated for updating the spectrum of the low-resolution image corresponding to each illumination position. After a preset number of iterations, the high-resolution amplitude and phase information of the sample is obtained.
[0144] Please refer to Figure 4 , Figure 4 (a) in the image is a low-resolution image captured when the central LED unit is lit. Figure 4 (a1) in the diagram represents the selected reconstruction region. Figure 4 (a2) in the image is a magnified view of the reconstructed area. Figure 4 (b1) to (g1) in the figure are amplitude reconstruction images obtained by the AA method, PC method, AA-C method, SC method, SAA method and the method of this application, respectively; Figure 4 In the image, (b2) to (g2) represent the phase reconstruction images obtained by the corresponding methods. Figure 4 (b3) to (g3) are magnified views of the preset local regions in the corresponding phase reconstruction images.
[0145] Depend on Figure 4 The reconstruction results show that, in the presence of LED array geometric misalignment, the amplitude and phase reconstruction images obtained by the AA method exhibit significant wrinkles and background artifacts. The PC, AA-C, SC, and SAA methods can reduce these artifacts to varying degrees, but structural distortion or phase discontinuities still exist in localized areas. By using the method of this application to repair the spectrum and compensate for energy, and obtaining the LED array geometric misalignment parameters, the spectral sub-regions corresponding to each low-resolution image can be determined based on the corrected illumination wave vector, thereby reducing the impact of LED position errors on the amplitude and phase reconstruction results.
[0146] Figure 3 and Figure 4 The "SAAS" appearing in this document indicates the experimental identifier used in the experimental verification phase of the method in this application. This experimental identifier is not intended to limit the technical composition of the method in this application; the technical content of the method in this application is based on spectrum restoration, non-DC spectrum energy compensation, adaptive optimization, and lighting error correction.
[0147] The beneficial effects of this application embodiment compared with the prior art are as follows: By combining periodic smoothing decomposition, orthogonal axis anomalous spectrum pixel neighborhood interpolation, and non-DC spectrum energy compensation, the local continuity of the orthogonal axis region is restored while suppressing spectral cross-shaped artifacts, and the non-DC spectrum energy attenuation caused by spectrum restoration is reduced; inputting the restored spectrum feature image into the pre-trained U-Net network to extract circular spectrum features can reduce the interference of spectral artifacts on pupil position recognition; setting translation search step size and rotation search step size respectively, and adjusting the two types of search step sizes according to different attenuation ratios when the local unoptimized counter reaches the first preset number, so that the illumination geometric parameter search process gradually transitions from global coarse search to local fine search; using the obtained LED array geometric misalignment parameters to correct the illumination wave vector can reduce the influence of LED array installation error on the Fourier stacked microscopic amplitude reconstruction results and phase reconstruction results.
[0148] Please refer to Figure 5 , Figure 5 (a1) in the diagram represents the vertical lighting mode of the central lighting unit. Figure 5 (a2) in the diagram represents the tilted lighting mode for non-central lighting units; Figure 5 In the image, (b1) and (b2) are low-resolution sample images acquired under the corresponding illumination modes, respectively. Figure 5 In the diagram, (c1) and (c2) are the spectral feature maps obtained by two-dimensional Fourier transform of the local regions of the corresponding low-resolution sample images.
[0149] When the central illumination unit provides vertical illumination, the circular spectral feature is close to the center of the spectrum; when the non-central illumination unit provides tilted illumination, the circular spectral feature shifts relative to the center of the spectrum. Based on the shift position of the circular spectral feature, the deviation between the actual illumination wave vector and the nominal illumination wave vector can be determined, and the first and second translations generated by the light source array along the first and second mutually orthogonal directions in the array plane, as well as the rotation of the light source array around the optical axis, can be further calculated.
[0150] The above embodiments are only used to illustrate the technical solutions of this application, and are not intended to limit it. Although this application has been described in detail with reference to the above embodiments, those skilled in the art should understand that modifications can still be made to the technical solutions described in the above embodiments, or equivalent substitutions can be made to some of the technical features. The above modifications or equivalent substitutions do not cause the substance of the corresponding technical solutions to deviate from the spirit and scope of the technical solutions of the various embodiments of this application, and should all be included within the protection scope of this application.
Claims
1. A method for correcting microscopic illumination errors through spectral restoration and adaptive optimization, characterized in that, Includes the following steps: S1. Spectrum restoration and energy compensation: Control the different light source units in the light source array to illuminate the sample in sequence, and obtain the low-resolution image corresponding to each illumination position; Perform a two-dimensional Fourier transform on the low-resolution image to obtain the original spectrum; The low-resolution image is periodically smoothed to obtain periodic components, and a two-dimensional Fourier transform is performed on the periodic components to obtain a periodic spectrum. Abnormal spectrum pixels located in the orthogonal axis region of the periodic spectrum are interpolated in a neighborhood. An energy compensation coefficient is determined based on the non-DC spectrum energy of the original spectrum and the non-DC spectrum energy of the interpolated spectrum. Energy compensation is then performed on the non-DC component of the interpolated spectrum while keeping the DC component unchanged to obtain a recovered spectral feature image. S2. Circular Spectral Feature Extraction: Input the recovered spectral feature image into a pre-trained feature extraction network to obtain circular spectral features representing the pupil position; S3. Adaptive optimization of illumination geometry parameters: Set translational search step size and rotational search step size respectively, and generate candidate geometric misalignment parameters based on the current translational search step size and the current rotational search step size; determine the corresponding simulated pupil position according to each candidate geometric misalignment parameter, and generate a simulated circular spectral feature with the same data format as the circular spectral feature obtained in step S2; calculate the feature error between the simulated circular spectral feature and the circular spectral feature obtained in step S2, update the current optimal geometric misalignment parameter according to the feature error, and decrease the translational search step size and the rotational search step size respectively when a smaller feature error is not obtained after a preset number of consecutive steps. S4. Illumination error correction and Fourier layered microscopic reconstruction: Based on the geometric misalignment parameters obtained in step S3, correct the spatial position of each light source unit and the corresponding illumination wave vector, and perform Fourier layered microscopic reconstruction on the low-resolution images corresponding to different illumination positions based on the corrected illumination wave vector.
2. The method for correcting microscopic illumination errors through spectral restoration and adaptive optimization according to claim 1, characterized in that, In step S1, the periodic smooth decomposition satisfies: ; ; in, For the low-resolution image, It is a periodic component. To smooth the components, The boundary term is constructed based on the boundary discontinuities of the low-resolution image. It is a two-dimensional discrete Laplace operator.
3. The method for correcting microscopic illumination errors through spectral restoration and adaptive optimization according to claim 2, characterized in that, The boundary item Initialize to zero, and update the boundary position of the boundary item according to the following formula: ; ; ; ; in, The number of rows in the low-resolution image. The number of columns in the low-resolution image is given; a two-dimensional Fourier transform is performed on the boundary terms; when the frequency coordinates are not located at the center of the spectrum, the frequency domain representation of the smoothing component is determined according to the frequency domain expression of the two-dimensional discrete Laplacian operator; when the frequency coordinates are located at the center of the spectrum, the frequency domain value of the smoothing component is set to zero; the frequency domain representation of the smoothing component is subtracted from the two-dimensional Fourier transform result of the low-resolution image to obtain the periodic spectrum.
4. The method for correcting microscopic illumination errors through spectral restoration and adaptive optimization according to claim 3, characterized in that, At that time, adopt and Replace the corresponding abnormal spectrum pixels with the average value; when At that time, adopt and Replace the corresponding abnormal spectrum pixels with the average value; when and At that time, the DC component at the center of the spectrum remains unchanged; among which, The coordinates in the periodic spectrum are Spectral values, The coordinates are the center of the spectrum.
5. The method for correcting microscopic illumination errors through spectral restoration and adaptive optimization according to claim 4, characterized in that, After completing the neighborhood interpolation of the abnormal spectrum pixels, the energy compensation coefficient is determined based on the non-DC spectral energy of the original spectrum of the low-resolution image and the non-DC spectral energy of the interpolated spectrum. The energy compensation coefficient satisfy: ; ; ; in, The coordinates in the original spectrum obtained by the two-dimensional Fourier transform of the low-resolution image are: Spectral values, This refers to the original non-DC energy of the spectrum after deducting the DC energy at the center of the spectrum. The non-DC energy of the interpolated spectrum after deducting the DC energy at the center of the spectrum.
6. The method for correcting microscopic illumination errors by spectral restoration and adaptive optimization according to claim 5, characterized in that, While keeping the DC component of the interpolated spectrum unchanged, the energy compensation coefficient is used to compensate the non-DC component of the interpolated spectrum to obtain the recovered spectral feature image, satisfying the following: ; in, The coordinates in the recovered spectral feature image are Spectral values, The coordinates in the interpolated spectrum after neighborhood interpolation are Spectral values, The energy compensation coefficient is... The coordinates of the center of the spectrum.
7. The method for correcting microscopic illumination errors by spectral restoration and adaptive optimization according to claim 6, characterized in that, The pre-trained feature extraction network is a U-Net network. The U-Net network is trained using the recovered spectral feature image and its corresponding manually labeled circular spectral regions as training samples. The U-Net network includes an encoder, a bottleneck layer, and a decoder. The number of channels in the convolutional layers of the encoder is 8, 16, 32, and 64, respectively. The bottleneck layer has 128 channels. The number of channels in the convolutional layers of the decoder is 64, 32, 16, and 8, respectively. Each convolutional layer uses a 3×3 convolutional kernel and is downsampled using a 2×2 max pooling operation. The U-Net network outputs the circular spectral features corresponding to the recovered spectral feature image.
8. The method for correcting microscopic illumination errors by spectral restoration and adaptive optimization according to claim 1, characterized in that, In step S3, a local unoptimized counter and a cumulative unoptimized counter are set; Calculate the feature error corresponding to each candidate geometric misalignment parameter in the current iteration, and determine the minimum feature error of this round and the candidate geometric misalignment parameter corresponding to the minimum feature error of this round; when the minimum feature error of this round is less than the historical minimum feature error, update the candidate geometric misalignment parameter corresponding to the minimum feature error of this round to the current optimal geometric misalignment parameter, update the minimum feature error of this round to the historical minimum feature error, and clear the local unoptimized counter and the cumulative unoptimized counter; when the minimum feature error of this round is not less than the historical minimum feature error, keep the current optimal geometric misalignment parameter and the historical minimum feature error unchanged, and increase the count values of the local unoptimized counter and the cumulative unoptimized counter respectively; when the count value of the local unoptimized counter reaches a first preset number, reduce the translation search step size to a first proportion before reduction, reduce the rotation search step size to a second proportion before reduction, and clear the local unoptimized counter; when the historical minimum feature error is not greater than a preset error threshold, or when the count value of the cumulative unoptimized counter reaches a second preset number, terminate the optimization.
9. The method for correcting microscopic illumination errors by spectral restoration and adaptive optimization according to claim 8, characterized in that, In step S3, the geometric misalignment parameters include a first translation and a second translation of the light source array along a first and a second direction that are orthogonal to each other in the array plane, and a rotation of the light source array about the optical axis. Centered on the current optimal geometric misalignment parameter, candidate translation parameters are generated along the positive and negative directions of the first and second directions respectively according to the current translation search step size, and candidate rotation parameters are generated along the positive and negative directions around the optical axis according to the current rotation search step size. Multiple sets of candidate geometric misalignment parameters are formed by the candidate translation parameters and the candidate rotation parameters.
10. The method for correcting microscopic illumination errors by spectral restoration and adaptive optimization according to claim 9, characterized in that, In step S4, the nominal spatial position of each light source unit is transformed by translation and rotation according to the first translation amount, the second translation amount, and the rotation amount to obtain the corrected spatial position of each light source unit; the illumination wave vector corresponding to each light source unit is recalculated according to the corrected spatial position, and the spectral sub-region corresponding to each low-resolution image is determined according to the recalculated illumination wave vector; phase recovery iteration is performed on each spectral sub-region in sequence to obtain the high-resolution amplitude information and phase information of the sample.